Applied Science, Faculty of
Engineering, School of (Okanagan)
DSpace
UBCO
Lofroth, Nicole Marie
Master of Applied Science - MASc
University of British Columbia
Electromagnetic transient simulation has become a prominent method for intricate power system design and analyses. Current research emphasizes the challenges that arise in developing a computational algorithm for efficient and accurate simulations of power systems due to the conjoined behaviour of continuous states and discrete events. Common methodologies, known as discrete state event driven (DSED) approaches, aim to overcome the modeling hindrances through employing a fixed or variable simulation step-size and integration order that reduces the number of computational points required to achieve high accuracy with accelerated transient simulation. The exponential-integrator-based simulation algorithm proposed is a novel variable-step fixed-order (VSFO) method that performs sequential intermediate integration calculations between the discrete events throughout the simulation. The matrix exponential integrator portion of the algorithm offers a competitive approach for an accurate, numerically efficient, and stable technique to accelerate transient simulation. Furthermore, the algorithm is easily adapted to accommodate stiff or nonstiff circuit systems which complement the desired topology and parameters. This thesis conducts a comprehensive analysis of the foundational aspects of the theoretical model for an exponential integrator algorithm and its practical applications for simulation of power converter systems. In addition, the simulation based case studies are presented for the validation of the algorithm, and provide insightful perspective and approach to converter control modeling and operating principles. A three-phase two-level AC/DC converter, a three-phase two-level AC/DC/AC converter, and a six-pulse rectifier are simulated by the proposed matrix exponential integrator algorithm. These case studies are typical and commonly used converter typologies with active and passive switching events. They demonstrate the proposed algorithm’s ability to accurately capture the hybrid nature of events within the system before each execution of a simulation cycle. The comparison of CPU times of the proposed exponential integrator algorithm to the models built by Simulink/Simscape Electrical toolbox show the simulation efficiency improvements are 2.3749913 s to 15.8125000 s for Case Study One, 8.5514098 s to 15.6562500 s for Case Study Two, and 0.2655774 s to 2.7656250 s for Case Study Three.
Matrix Exponential IntegratorAlgorithm for Accurate and EfficientSimulation of Power ConverterSystemsbyNicole Marie LofrothBSc. Hons., University of Northern British Columbia, 2022A THESIS SUBMITTED IN PARTIAL FULFILLMENT OFTHE REQUIREMENTS FOR THE DEGREE OFMASTER OF APPLIED SCIENCEinTHE COLLEGE OF GRADUATE STUDIES(Electrical Engineering)THE UNIVERSITY OF BRITISH COLUMBIA(Okanagan)July 2024© Nicole Marie Lofroth 2024The following individuals certify that they have read, and recommendto the College of Graduate Studies for acceptance, the thesis entitled:Matrix Exponential Integrator Algorithm for Accurate and Efficient Simulationof Power Converter Systemssubmitted by Nicole Marie Lofroth in partial fulfillment of the requirementsfor the degree of Master of Applied Science.Dr. Liwei Wang, School of EngineeringSupervisorDr. Wilson Eberle, School of EngineeringSupervisory Committee MemberDr. Zheng Liu, School of EngineeringSupervisory Committee MemberDr. Sina Kheirkhah, School of EngineeringUniversity ExamineriiAbstractElectromagnetic transient simulation has become a prominent method forintricate power system design and analyses. Current research emphasizes thechallenges that arise in developing a computational algorithm for efficientand accurate simulations of power systems due to the conjoined behaviourof continuous states and discrete events. Common methodologies, known asdiscrete state event driven (DSED) approaches, aim to overcome the mod-eling hindrances through employing a fixed or variable simulation step-sizeand integration order that reduces the number of computational points re-quired to achieve high accuracy with accelerated transient simulation.The exponential-integrator-based simulation algorithm proposed is a novelvariable-step fixed-order (VSFO) method that performs sequential interme-diate integration calculations between the discrete events throughout thesimulation. The matrix exponential integrator portion of the algorithmoffers a competitive approach for an accurate, numerically efficient, andstable technique to accelerate transient simulation. Furthermore, the al-gorithm is easily adapted to accommodate stiff or nonstiff circuit systemswhich complement the desired topology and parameters. This thesis con-ducts a comprehensive analysis of the foundational aspects of the theoret-ical model for an exponential integrator algorithm and its practical appli-cations for simulation of power converter systems. In addition, the simu-lation based case studies are presented for the validation of the algorithm,and provide insightful perspective and approach to converter control mod-eling and operating principles. A three-phase two-level AC/DC converter,a three-phase two-level AC/DC/AC converter, and a six-pulse rectifier aresimulated by the proposed matrix exponential integrator algorithm. Theseiiicase studies are typical and commonly used converter typologies with activeand passive switching events. They demonstrate the proposed algorithm’sability to accurately capture the hybrid nature of events within the sys-tem before each execution of a simulation cycle. The comparison of CPUtimes of the proposed exponential integrator algorithm to the models builtby Simulink/Simscape Electrical toolbox show the simulation efficiency im-provements are 2.3749913 s to 15.8125000 s for Case Study One, 8.5514098s to 15.6562500 s for Case Study Two, and 0.2655774 s to 2.7656250 s forCase Study Three.ivLay SummaryPower converter systems containing continuous states and discrete eventsrequire a comprehensive approach to target events for efficient and accuratesimulation results. The use of digital controllers in power converter systemsprovides the ability to predict time instances at which the semiconductorswitching events occur. This thesis proposes a novel exponential integratoralgorithm for numerically efficient simulation of power electronic converters.The computational process is variable-step fixed-order but is adaptable forfixed step-sizes, and aims to pre-compute the timing of chosen discretizationsteps, which are then followed sequentially by intermediate integrations. Atheoretical analysis is introduced on the proposed simulation algorithm, andthe reasoning for the selection of power converters and the associated con-trollers is provided. Each simulated converter system is examined with thedetailed operating scheme, either with a power and voltage controller forthe detection of active switching events, or through threshold value compar-isons for detection of passive switch events. Simulation results are furtheranalyzed and discussed, which show the improved simulation efficiency of2.3749913 s to 15.8125000 s for Case Study One, 8.5514098 s to 15.6562500s for Case Study Two, and 0.2655774 s to 2.7656250 s for Case Study Three.vPrefaceThe work presented in this thesis was completed at the University of BritishColumbia - Okanagan, between September 2022 to May 2024. I assumedresponsibility for devising the proposed modeling approach, implementingthe algorithm of various test cases, and validating event detection simula-tion. The completed manuscript was authored exclusively by myself. Mysupervisor, Dr. Liwei Wang, guided me through my research and offeredsupport and thoughtful recommendations. Jared Paull, a fellow lab mem-ber, provided guidance through the simulation algorithm development andvalidation of my models. The work in Chapter 4 is based on the analysis of[1], which has been submitted to IEEE Transactions on Power Delivery.Generative AI was not used for the completion of this research.viTable of ContentsAbstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iiiLay Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vPreface . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viTable of Contents . . . . . . . . . . . . . . . . . . . . . . . . . . . . viiList of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xList of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiAbbreviations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xivCommon Symbols . . . . . . . . . . . . . . . . . . . . . . . . . . . xvAcknowledgements . . . . . . . . . . . . . . . . . . . . . . . . . . . xvi1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.1 Research Objectives . . . . . . . . . . . . . . . . . . . . . . . 21.2 Contributions . . . . . . . . . . . . . . . . . . . . . . . . . . 31.3 Thesis Organization . . . . . . . . . . . . . . . . . . . . . . . 42 Literature Review . . . . . . . . . . . . . . . . . . . . . . . . . 52.1 Overview . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52.2 Simulation Framework . . . . . . . . . . . . . . . . . . . . . 72.3 Converter Controls and Discrete Event Driven Modeling . . 82.4 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12vii3 Background on Converter Topologies and Operating Prin-ciples . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 143.1 Three-Phase Two-Level Converters . . . . . . . . . . . . . . 143.1.1 Topology . . . . . . . . . . . . . . . . . . . . . . . . . 143.1.2 Operating Principle . . . . . . . . . . . . . . . . . . . 163.1.3 Converter Switching Model . . . . . . . . . . . . . . . 183.1.4 Park Transformation . . . . . . . . . . . . . . . . . . 183.2 Three-Phase Diode Rectifier . . . . . . . . . . . . . . . . . . 203.2.1 Topology . . . . . . . . . . . . . . . . . . . . . . . . . 203.2.2 Operating Principle . . . . . . . . . . . . . . . . . . . 223.2.3 Diode Model . . . . . . . . . . . . . . . . . . . . . . . 224 Matrix Exponential Integrator Algorithm . . . . . . . . . . 244.1 State Space Equation Modeling . . . . . . . . . . . . . . . . 254.2 Algorithm Functions . . . . . . . . . . . . . . . . . . . . . . . 254.2.1 Matrix Exponential Derivation . . . . . . . . . . . . . 254.2.2 Phi-Functions and Taylor Series Expansion . . . . . . 264.2.3 Pade´ Approximation . . . . . . . . . . . . . . . . . . 274.3 Simulation Framework . . . . . . . . . . . . . . . . . . . . . 284.4 Numerical Stability . . . . . . . . . . . . . . . . . . . . . . . 335 Matrix Exponential Integrator Based Converter Modeling 355.1 Discrete Event Detection . . . . . . . . . . . . . . . . . . . . 365.1.1 Discrete Control Models . . . . . . . . . . . . . . . . 365.1.2 Control in dq-Frame . . . . . . . . . . . . . . . . . . . 365.1.3 Mathematical Overview . . . . . . . . . . . . . . . . . 375.1.4 DC Voltage Control . . . . . . . . . . . . . . . . . . . 395.1.5 Zero-Order Hold Implementation . . . . . . . . . . . 425.1.6 Phase-Locked Loop Discretization . . . . . . . . . . . 445.1.7 Controller Discretization . . . . . . . . . . . . . . . . 475.2 Passive Event Models . . . . . . . . . . . . . . . . . . . . . . 495.2.1 Iterative Zero Crossing Detection . . . . . . . . . . . 505.3 Simulation and Validation . . . . . . . . . . . . . . . . . . . 51viii5.3.1 Case Study One . . . . . . . . . . . . . . . . . . . . . 525.3.2 Case Study Two . . . . . . . . . . . . . . . . . . . . . 625.3.3 Case Study Three . . . . . . . . . . . . . . . . . . . . 716 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 806.1 Thesis Summary . . . . . . . . . . . . . . . . . . . . . . . . . 816.2 Future Work . . . . . . . . . . . . . . . . . . . . . . . . . . . 82Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84ixList of Tables5.1 AC/DC Three-Phase Two-Level Converter . . . . . . . . . . . 535.2 Computational Efficiency Comparison for Three-Phase Two-Level AC/DC VSC . . . . . . . . . . . . . . . . . . . . . . . . 615.3 AC/DC/AC Three-Phase Two-Level Converter . . . . . . . . 635.4 Computational Efficiency Comparison Three-Phase Two-LevelAC/DC/AC VSC . . . . . . . . . . . . . . . . . . . . . . . . . 715.5 Six Pulse Rectifier . . . . . . . . . . . . . . . . . . . . . . . . 725.6 Computational Efficiency Comparison Six Pulse Rectifier . . 79xList of Figures2.1 DSED Control Mechanism for VSFO Algorithms . . . . . . . 113.1 Topology of an AC/DC Single-Phase Half-Bridge Converter . 153.2 Topology of a Three-Phase Two-Level VSC . . . . . . . . . . 163.3 PWM Scheme . . . . . . . . . . . . . . . . . . . . . . . . . . 173.4 Park Transformation . . . . . . . . . . . . . . . . . . . . . . . 193.5 AC/DC/AC Three-Phase Two-Level VSC . . . . . . . . . . . 203.6 Full Wave Bridge Rectifier . . . . . . . . . . . . . . . . . . . . 213.7 Six-Pulse Rectifier . . . . . . . . . . . . . . . . . . . . . . . . 213.8 Ideal Diode Switch Behaviour (adapted from [51]) . . . . . . 234.1 Flow Chart of Matrix Exponential Integrator AlgorithmWithVariable-Step and Variable Number of φ-Function Terms (adaptedfrom [1]) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 294.2 Flow Chart of AGEI AlgorithmWith Pre-Computations (adaptedfrom [1]) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 325.1 Current-Controlled Real/Reactive Power Controller in dq-Frame . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 375.2 DC Voltage Control Scheme . . . . . . . . . . . . . . . . . . 395.3 Block Diagram id . . . . . . . . . . . . . . . . . . . . . . . . . 415.4 Block Diagram iq . . . . . . . . . . . . . . . . . . . . . . . . . 415.5 Comparison of an Iterative and Pre-Scheduling Approach forSwitch Event Detection (adapted from [3]) . . . . . . . . . . . 435.6 PLL Control Block Diagram . . . . . . . . . . . . . . . . . . . 455.7 Three-Phase AC Voltage Source Input . . . . . . . . . . . . . 54xi5.8 Real and Reactive Power Closed-Loop Control . . . . . . . . 545.9 Real and Reactive Power Closed-Loop Control Comparison . 555.10 AC/DC Three-Phase Two-Level VSC 5µs Proposed Algo-rithm and 5µs Simulink Comparison - Three Phase InductorCurrents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 565.11 AC/DC Three-Phase Two-Level VSC 5µs Proposed Algo-rithm and 5µs Simulink Comparison - Combined DC Capac-itor Voltage . . . . . . . . . . . . . . . . . . . . . . . . . . . . 575.12 AC/DC Three-Phase Two-Level VSC 1µs Proposed Algo-rithm and 1µs Simulink Comparison - Three Phase InductorCurrents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 585.13 AC/DC Three-Phase Two-Level VSC 1µs Proposed Algo-rithm and 1µs Simulink Comparison - Combined DC Capac-itor Voltage . . . . . . . . . . . . . . . . . . . . . . . . . . . . 595.14 Zoomed-In Waveform AC/DC Three-Phase Two-Level VSCProposed Algorithm and Simulink Comparison - Three PhaseInductor Currents . . . . . . . . . . . . . . . . . . . . . . . . 605.15 Three-Phase AC Voltage Source Input . . . . . . . . . . . . . 645.16 Real/Reactive Power and DC Voltage Closed-Loop ControlComparison . . . . . . . . . . . . . . . . . . . . . . . . . . . . 655.17 AC/DC/AC Three-Phase Two-Level VSC 1µs Proposed Al-gorithm and 1µs Simulink Comparison - Three Phase Induc-tor Currents AC Side One . . . . . . . . . . . . . . . . . . . . 665.18 AC/DC/AC Three-Phase Two-Level VSC 1µs Proposed Al-gorithm and 1µs Simulink Comparison - Three Phase Induc-tor Currents AC Side Two . . . . . . . . . . . . . . . . . . . . 675.19 AC/DC/AC Three-Phase Two-Level VSC 1µs Proposed Al-gorithm and 1µs Simulink Comparison - Combined DC Ca-pacitor Voltage . . . . . . . . . . . . . . . . . . . . . . . . . . 685.20 Zoomed-In Waveform AC/DC/AC Three-Phase Two-LevelVSC Proposed Algorithm and Simulink Comparison - ThreePhase Inductor Currents AC Side One . . . . . . . . . . . . . 69xii5.21 Zoomed-In Waveform AC/DC/AC Three-Phase Two-LevelVSC Proposed Algorithm and Simulink Comparison - ThreePhase Inductor Currents AC Side Two . . . . . . . . . . . . . 705.22 Three-Phase AC Voltage Source Input . . . . . . . . . . . . . 745.23 Six Pulse Rectifier Proposed Algorithm and Simulink Com-parison 1µs - Three Phase Inductor Currents . . . . . . . . . 755.24 Six Plulse Rectifier Proposed Algorithm and Simulink Com-parison 1µs - Capacitor Voltage . . . . . . . . . . . . . . . . . 765.25 Six Pulse Rectifier Proposed Algorithm - Diode Currents . . . 775.26 Zoomed-In Waveform Six Pulse Rectifier Proposed Algorithmand Simulink Comparison 1µs - Three Phase Inductor Currents 78xiiiAbbreviationsAC Alternating CurrentAGEI Adaptive-Grained Exponential IntegratorDC Direct Currentdq-Frame Direct-Quadrature-FrameDSED Discrete State Event DrivenFSFO Fixed-Step Fixed-OrderFSVO Fixed-Step Variable-OrderIGBT Insulated Gate Bipolar TransistorLTE Linear Trunacation ErrorODE Ordinary Differential EquationPCC Point of Common CouplingPI Proportional IntegratorPID Proportional Integrator DerivativePLL Phase-Locked LoopPWM Pulse Width ModulationVSC Voltage-Sourced ConverterVSFO Variable-Step Fixed-OrderVSVO Variable-Step Variable-OrderZOH Zero-Order HoldxivCommon SymbolsC CapacitanceCDC DC converter capacitanceC1 Upper DC converter capacitanceC2 Lower DC converter capacitanceFsw Switching frequencyIL Inductor currentki Integrator constantkp Proportional constantL InductanceP Real power transferQ Reactive power transferR ResistanceRoff Diode off-state resistanceRon Diode on-state resistanceron IGBT on-state resistanceTs Converter step-timeTsw Switching periodVAC AC voltage sourceVcth Condition threshold voltageVDC DC line-to-line voltageVDC1 Upper DC line-to-line voltageVDC2 Lower DC line-to-line voltageVf Forward diode voltageVpk AC voltage source peak voltagew Radial frequencyϕ Phase anglexvAcknowledgementsI wish to thank my supervisor, Dr. Liwei Wang. Your continuous supportand encouragement over the past two years has been invaluable. I havegreatly appreciated the time and effort you have spent ensuring my success,and more specifically I appreciate the patience shown me during the researchprocess. Your expertise in this field of research is to be admired.I would also like to express my gratitude to fellow lab member JaredPaull. Thank you for taking the time to work with me on research prob-lems, for all of your advice, and for your willingness to help in any situation.I have been extremely fortunate for the opportunity to work alongside youduring my Master’s degree.Lastly, I would like to thank my brother and parents for their endlesssupport in making my achievements possible.xviChapter 1IntroductionElectromagnetic transient simulations of power converter systems are chal-lenged by the frequently changing circuit topologies due to semiconductorswitching events. While certain switching devices have foreseeable switchingtimes, such as fully controlled devices, other switching devices, e.g., diodes,may have unanticipated switching instances. The complexity of these dis-crete switch events, along with the continuous states of power systems, createimpediments for accurate representations of rapidly changing circuit topolo-gies. Thus, as power converter systems continue to grow more complex,there is a demand for numerically efficient and accurate simulation algo-rithms.Fast transients introduce system state change and semiconductor switch-ing, governed by converter’s state space equations. Order selection is chosenbased on achieving high accuracy or high numerical efficiency for the desiredsimulation of a power converter system. Higher order methods are generallyimplemented to reduce error propagation which is proportional to step-size.Lower order methods may not have a large effect when reducing error, andare usually implemented to enhance numerical efficiency. Low order im-plicit ordinary differential equation solvers, such as Backward Euler methodor Trapezoidal Rule, have been used to ensure numerical accuracy and sta-bility. However, a small integration step-size must be used in conjunctionto compensate for the relatively low order. Consequently, fixed step-sizesolvers often acquire a large number of computational points that accumu-late over a simulation period, resulting in an inefficient numerical solver.Variable step-size solvers serve as alternatives to fixed-step solvers, whichaim to optimize the time at which state matrices are updated, based on the1detection of switching events. The variable-sized modality reduces the step-size when transients occur and increases the step-size once the steady stateis reached. However while the adaptation of step-size is beneficial, the sizeof the step-size is not always an obvious choice. In order to achieve optimaltransient simulation with high numerical efficiency, variable step-size solverscan be used to minimize computational time at a given numerical accuracy.Combining variable step-sizes with state space equation representations ofa power converter system introduces the opportunity for a high-order inte-gration method.Exponential integrators are a preferred numerical method to both loworder discretization methods and Taylor-series based approximations. Thehigh-order exponential method may implement a Pade´ approximation, whichincreases discretization order to increase accuracy, ensure stability, andeliminate unwanted numerical oscillations. In converter systems with highswitching frequencies, an exponential integration technique is more likelyto accurately detect these events and determine the step-size accordinglyin order to solve ordinary differential equation (ODE) systems. The algo-rithm seeks to reconstruct the exponential integrator algorithm techniqueby introducing a variable-step high-order approach with the ability to pre-dict the timing of switching events before they occur. This switching eventdetection provides a numerically efficient simulation of power converter sys-tems. The exponential integrator is a vital mathematical tool to overcomethe issue of low order integration methods that have commonly hindered theperformance of many simulation algorithms.1.1 Research ObjectivesThe objectives of the research conducted and the case studies presentedare to provide a numerically efficient and accurate simulation algorithm ofpower converter systems and the use of digital controller’s event triggeringmechanism to achieve accelerated simulation. Accurate and numerically sta-ble simulation require the state dynamics to be modeled through a matrix2exponential integrator algorithm which is introduced and broken down intoeach of its elements to offer a thorough insight into the computational pro-cess. A presentation and discussion of simulation results for diverse powerconverter topologies is included to verify the effectiveness of the algorithmfor power converter systems. The analysis provides a better understandingof the design and implementation of the algorithm. Additional power con-verter design and analysis methods required for modeling and simulationpurposes are reviewed. They provide justification of parameter selection ofcontrol methodologies, and enhance the credibility of the proposed modelingframework.1.2 ContributionsThis thesis proposes a matrix exponential integration technique which is aninnovative approach to simulate power converter systems. The prominentcontributions presented within regarding the applications of the algorithmare as follows:1. The proposed matrix exponential integrator algorithm is easily adapt-able for both fixed-step or variable-step solvers with a high-order dis-cretization strategy. This allows for high discretization accuracy ofstate space equations to accelerate electromagnetic transient simula-tion of power converter systems.2. Timing of active switching events that occur throughout a simulationperiod are pre-determined to avoid iterative calculations, speeding upthe computational process through discrete event detection control.3. The variable-step high-order algorithm does not require additional cal-culations between discrete events when achieving the desired simula-tion accuracy.4. The number of forcing function φ terms in the proposed matrix expo-nential integration algorithm can be varied in order to meet the desirederror thresholds.35. The proposed matrix exponential integration algorithm is L-stable,which makes it suitable for a large variety of power converter topolo-gies, and is flexible against diverse system parameters.1.3 Thesis OrganizationThe chapters of this thesis are predicated upon the modeling and simula-tion framework of power converter systems through a novel numerical inte-gration algorithm. Chapter 2 reviews current literature concerning circuitsimulation algorithms, converter control techniques, and potential applica-tions the research may serve. Chapter 3 highlights background on powerconverter systems and their operating principles, in addition to modelingtechniques for replicating their characteristics in steady state and dynamictransients. Chapters 4 and 5 introduce a detailed description of the matrixexponential integrator algorithm used for the simulation of power convertersystems found in Chapter 3. Three typical case studies are used to verifythe improved numerical efficiency of the converter models implemented bythe proposed exponential integration algorithm. Finally, Chapter 6 offersthe conclusion and presents potential avenues for future research.4Chapter 2Literature Review2.1 OverviewA literature review is conducted to provide thorough analyses on the mod-eling and control of power converter systems. Power electronic converters,that possess switching semiconductor devices, operate with a coexisting na-ture of continuous states and discrete events [2, 3]. Traditionally, continuousstates include capacitor voltages and inductor currents, while discrete eventsrefer to the switching devices possessing the potential to change a circuit’stopology and lead to the transition from one operating mode to another[3]. The complex nature of switching devices still remains one of the largestchallenges in electric circuit simulation due to its rapid switching behaviour.While the majority of simulation modalities have operated with an idealswitch model, two predominant problems still remain: the estimation of ex-act switching times, and the calculation of initial conditions whose valuesdepend on the updates of system states throughout the simulation [4]. Al-gorithms for accurate and numerically stable simulation of power convertersystems focus on improving simulation efficiency of power converters withhigh switching frequency and on the modeling of device intricacies [5].Classical numerical integration techniques have become less efficient incapturing switching states that possess high switching frequencies [6]. Nu-merical techniques such as Forward Euler, Backward Euler, TrapezoidalRule, and Runge-Kutta have been the most common ODE solving methodswhen attempting to conquer the challenge [4, 7]. These methods are basedon discretization of differential equations with time where values of statevariables are updated simultaneously at each step-size, providing a solution5based on a difference equation [8]. Consequently, a nonlinear solver, suchas one listed above, may contribute to serious drawbacks of the simulationwhen used for integration in nonsmooth circuits, affecting the speed and ef-ficiency of the simulation results [7, 9]. Outputs of state values throughoutthe simulation, if simulated incorrectly, will contribute to the inefficiencyand may induce an impulse response in the circuit [4]. The variety of nu-merical integration techniques proposed in many research works originatedfrom the idea that these approaches guarantee preservation of performanceobtained in continuous time domain, and under an acceptable quantity ofdiscrete time approximations [10]. Thus, the validity of these methods arecircuit dependent.Current simulation methods have sought to develop an alternative nu-merical approach when calculating switching events. A large portion ofresearch has reconstructed classical numerical integration techniques suchas Runge-Kutta, and have developed innovative methods for system staterepresentations. In [11], a diagonal Runge-Kutta method for stiff ODEs wasadapted into five different new first-stage, singly diagonally-implicit Runge-Kutta methods, for comparison to other implicit-explicit analogs. The meth-ods presented stiff-accuracy, and internal L-stability, but order reductionvaried through the five methods and remained a problem throughout [11].[12] proposed a numerical integration technique by 2-stage diagonally im-plicit Runge-Kutta, whose accuracy and stability performance are similar tothose of the trapezoidal method, but eliminates numerical oscillations whensudden changes of inductor or capacitor values occurred. [6] and [7] workedwith a class of numerical techniques, known as quantized state systems, thatreplaced time discretization used by classical numerical integration meth-ods, by the quantization of state variables. Their results demonstrated alinearly implicit quantized system state algorithm that was noticeably moreaccurate than differential algebraic system solvers. However, under realisticconditions, these models exhibited frequency discontinuities [6]. Thus, therehas still been a need for an accurate and consistent numerical integrationtechnique for simulation power systems.62.2 Simulation FrameworkHigh frequency applications with pulse-width modulation (PWM) mecha-nisms and their controls place high demands on achieving precision in simu-lations, as desirable performance is not always guaranteed without detailedefforts in electromagnetic transient model design [13]. A large number ofsimulation frameworks operate under a state matrix equation model, rep-resenting the continuous states of a dynamic system that coexist with thediscrete events and co-determine the mode of operation. State space mod-eling is a more efficient approach than conventional nodal analysis methodwhen simulating power electronic converters. It follows the solution steps,by solving the differential equations describing the system at each switchinterval, and analyzes the conditions that determine the change of the cir-cuit topology [14]. The piece-wise linear modeling technique corresponds topiece-wise linear device segments, fostering the constant values of the ma-trices that change according to the system events [2].State space simulation models can be predominantly classified into fourcategories based on flexibility: fixed-step fixed-order (FSFO) algorithms,fixed-step variable-order (FSVO) algorithms, variable-step fixed-order (VSFO)algorithms, and variable-step variable-order(VSVO) algorithms [3]. Achiev-ing high computational efficiency, accuracy, and stability necessitates a flex-ible framework that can accommodate the changing states of the system byadjusting both the step-size and the order of the solver [2, 3]. At each sim-ulation step-size, it is beneficial that the solver detects the current or nextdiscrete event and determines the initial estimate of the step-size as the timeinterval between the current time instant and the next event [15]. Variablestep-size methods are advantageous when compared to fixed step-size meth-ods due to their adaptability in use of a large step-size for high numericalefficiency [2]-[16]. In addition, variable-order solvers assist in the process ofattaining high numerical accuracy, but have shown an inability to convergewhen dealing with stiff systems [2]. Therefore, a variable-order method mayhave limited application. A fixed-order method would suffice should it guar-7antee accuracy with a small step-size when following the fastest transient.Recent literature has proposed power electronics system simulation ap-plications based on a Discrete Switching Event Driven (DSED) framework.The flexible adaptive DSED framework contributes to the improvement ofnumerical accuracy and computational efficiency by employing a VSVO al-gorithm [3, 17]. Its explicit nature is desirable for numerically efficient sim-ulation of power electronics systems, but it faces challenges for numericallystiff state equations of power converters with parasitic elements. Additionalstudies have employed this framework for high frequency power converters,backward DSED algorithms, quantized state system representation, and de-coupling strategies [2, 7, 8, 15, 17]. The results from these studies confirmedthe effectiveness of exploiting the recursive property of state space equationsto obtain time derivatives of state variables. The outcome of the VSVOalgorithm in [3, 17] was based on adaptive-order selection and step-size ad-justment embedded into the integration process [3]. Further comparison toSimulink models confirmed significantly improved numerical efficiency of themethod and the possibility for future implementations of control strategiesembedded into the DSED framework [3, 15, 18].2.3 Converter Controls and Discrete EventDriven ModelingState space representation of power electronic converters can be particularlyuseful for controller design methods [19, 20]. Control theory has been stud-ied for both linear and nonlinear models of power converters. For voltage-sourced converters (VSCs), techniques are generally implemented throughPWM for current mode control and/or voltage mode control [21]. Nonlinearmodels have been extensively used for the design of feedback controllers forthree-phase PWM VSC systems, particularly for input-output linearizabilityof the system model in the dq-reference frame [22]. The abc-dq coordinatetransformation for three-phase circuits have been the basis for control and8instantaneous power theory. Back to back AC/DC/AC power convertertopologies are among the most popular converters to enable reliable bidirec-tional power flow and power management when operating under feedbackcontrol systems [24]. Control of the AC and DC components are easily inte-grated into the circuit for regulation of sinusoidal currents and DC compo-nents to achieve ideal characteristics of the converter [25]. However, currentand voltage mode control inherit problems of instability and sub-harmonicoscillations when they operate under constant frequency PWM. A compen-sator network is usually required [21]. Proportional integrator (PI) typecompensators are the more popular candidates among the possible compen-sator techniques for regulation of error and dynamic stability performance[20]-[29].Modern day research has used control strategies based on analog con-troller implementation, which hinders the practical advantages of advancedcontrol schemes. As an alternative to analog control implementation, digitalimplementation in modern digital signal controllers offers greater flexibility,simpler implementation of complex control algorithms, and a higher tol-erance to noise [26, 34]. Digital control strategies have increasingly beenused for power electronic simulation for control of generation, transmission,distribution, and utilization levels [35]. State space digital control method-ologies are able to reconstruct the system states via sampling of the plantoutput where the states are calculated through a numerical integration al-gorithm [36, 37]. In general, it is well accepted that digital control offersmore functionality to the system compared to analog implementations [34].Power system simulation with digital control regularly operates with azero-order hold (ZOH) discretization step which ”freezes” the control sig-nal over a sampling period. A ZOH has been frequently used in practice,particularly in feedback control applications and switching device behaviourprediction [3, 36, 37, 39]. The discretization method is crucial for the im-plementation of fast digital control loops in multi-level converters and toavoiding additional computations over a simulation cycle [38, 40].9Current digital control schemes focus on predictive current controlledmodes [39, 41, 42], damping control [43, 44], or multi-phase switching con-verter control [37, 43, 45]. Predictive control techniques calculate and pre-determine the duty ratio in each switching cycle based on the previousswitching status [39]. This method however normally acquires one switch-ing cycle delay which significantly degrades the performance of digital con-trol devices in switching converters [39]. Damping proportional integratorderivative (PID) control methods require additional algorithm modificationsto attain a desired transient and steady state performance and may be sen-sitive to variations in voltage magnitudes [43, 44]. Studies conducted inrespect to control of multi-phase VSCs have the ability to preserve currentor voltage control depending on importance, and have proven their abilityto eliminate sub-harmonic oscillations [37, 39, 45].DSED literature, [3], has proposed an efficient event locating mechanism.The focus is on employing a flexible-adaptive VSVO algorithm with digitalcontrol techniques to pre-schedule events, and to avoid iterative calculations.Unlike regular-sampled PWM strategies, the DSED approach implements acontrol scheduler that directly provides the occurring time of the next con-trol event [3]. Timing of sampling events are determined by the user-definedsampling frequency and active switching instances can be pre-calculated ex-plicitly at the beginning of each cycle.Additionally, [3] proposed that VSVO algorithms under digital controlperform with optimized computational efficiency compared to FSFO andVSFO algorithms. When the time interval between consecutive events isshort, it is expected that VSVO algorithms decrease the order to reducecomputational cost, while the LTE is restricted by the error tolerance. Sim-ilarly when the time interval is large, the order can be increased to satisfythe accuracy requirements for an enlarged step-size. The trade off betweenhigh-order and less calculation points, to low-order with an increase in cal-culation points, often produces a similar efficiency outcome [3]. VSFO algo-10rithms are optimized by determining the largest step-size value that keepsthe LTE held below the specified tolerance. Since accuracy requirementsdepend upon both the algorithm and discrete events, VSFO methodologiesare viewed to be inefficient and perform redundant calculations [3]. Thismotivates the development of a variable-step high-order integration algo-rithm. A high-order method will facilitate accuracy requirements, througha ZOH specifically, when step-sizes are large and decrease computationaltime compared to small fixed-step algorithms or low fixed-order algorithms.High-order solvers will reach the desired accuracy by only updating the sys-tem state equations at switching events with no additional computations inbetween. The modulating signal will be held over the sampling period tocontribute to the decrease in computational CPU time, as shown in Fig-ure 2.1. If a VSFO method under digital control is optimized properly, themethod has the potential to outperform VSVO algorithms.Figure 2.1: DSED Control Mechanism for VSFO Algorithms11In general, digital control strategies and design procedures are amongthe more reliable methods for time descretization equations for simulationpurposes [45]. Specifically for control of three-phase systems, single-phaserectifiers, and DC-DC converters, digital techniques have been highly opti-mized for control of conventional current and voltage control [39].2.4 SummarySimulation frameworks that utilize the advantageous components of con-trollers, specifically digital controllers for discrete time simulation, have ahigher likelihood of accurate switching interpretations and greater stabilityin steady state operation. These methodologies facilitate all FSFO, FSVO,VSFO, and VSVO implementation possibilities for network equations andmatrix modifications [46]. The overall purposed simulation methods are topredict switch events and to modify simulation step-size to fall exactly onthe time instant required for turn on and turn off semiconductor switches,as time intervals between events are not identical [46, 47]. Transient anal-yses and numerical oscillation control are specific targets of interest for thenumerically accurate and stable simulation of power electronic converterswith digital control.Numerical integration algorithms, paired with simulation frameworks,have proven their improved numerical accuracy and efficiency through mul-tiple simulation studies [2, 3, 7], [15]-[17]. With such flexible algorithms,accuracy requirements can be satisfied under periods of slow change whileoperating with a minimal computational process for numerical integration[3]. Comparison to system-tool methods, such as Matlab/Simulink, has con-firmed the precision of event definition, state space modeling, event-drivenmechanisms, and numerical computation that the individual frameworkspossess. Future applications of the framework may help realize high speedsimulations and facilitate analyses of effects caused by electromagnetic infer-ence for the evaluation of overall system stability [2, 15]. Further improve-ments require algorithm development for low frequency power systems, for12switching loss interpretation, and for the simultaneous simulation of the sys-tem under practical control strategies [15]. If successfully applied to largepower electronic systems, the algorithms have the potential to prove them-selves, under various circuit topologies and input sources, providing for abroader range of applications with fewer limitations [13].13Chapter 3Background on ConverterTopologies and OperatingPrinciplesChapter 3 focuses on the topology and operating principles of AC/DC two-level converters and AC/DC full wave rectifiers. These converters are widelyused in power industry and utility and act as the fundamental buildingblocks of more complicated power electronic converters. More specifically,an AC/DC three-phase two-level VSC, an AC/DC/AC three-phase two-levelVSC, and a six-pulse diode rectifier are used in the case studies to verifythe improved numerical efficiency of the proposed exponential integratorsimulation framework. The following sections describe a detailed switchingmodel of the single-phase half-bridge AC/DC converter and a three-phasetwo-level converter. The Park transformation of a three-phase two-levelconverter into the dq-frame, needed for control purposes, is introduced andapplied to the AC/DC converter’s closed-loop controls in the later chapters.Additionally, detailed analysis is performed on single-phase and three-phasefull-wave rectifiers for the development of passive switching device modelingin the proposed exponential integration simulation framework.3.1 Three-Phase Two-Level Converters3.1.1 TopologyThe topology of a half-bridge AC/DC converter is shown in Figure 3.1. Thebasic converter presents a single-phase converter with two switch cells, also14known as IGBTs, with the upper switch as S1 and the lower switch as S2.Each individual transistor has positive current flowing from the collector tothe emitter, and each individual diode inside the IGBT module has positivecurrent flow from the anode to the cathode. The DC side capacitors, CDC ,are represented by C1 and C2. Equivalently, the AC voltage input side isshown by VAC .The half-bridge converter is connected to two DC voltage sources havingan equivalent voltage of VDC/2. In between these two sources, the voltagereference node is found, labeled as g. From the AC side, a series RL branchis an interfacing reactor that acts as a filter to ensure a low ripple for theAC side current. The interfacing reactor’s resistance and inductance arerepresented by R and L respectively.Figure 3.1: Topology of an AC/DC Single-Phase Half-Bridge ConverterBuilding off of Figure 3.1, Figure 3.2 illustrates the topology of a three-phase two-level VSC. The converter is composed of three identical half-bridge converters with each phase denoted by subscripts a, b, and c. ThisVSC is known as a two-level converter as each of the AC side terminalsassumes the voltage levels as -VDC/2 or VDC/2. All are connected to athree-phase AC voltage source. The DC side terminal follows the samestructure as the single-phase half-bridge AC/DC converter, while providing15bi-directional power flow to the three-phase AC system.Figure 3.2: Topology of a Three-Phase Two-Level VSC3.1.2 Operating PrincipleSingle-phase half-bridge converters rely on an operating principle knownas Pulse Width Modulation (PWM) strategy [48]. A PWM strategy isemployed to operate the two complementary switches of a phase in an ONor OFF state, which attain logical values 1 or 0 respectively. This mechanismis accomplished by the comparison of a triangular waveform with a periodof Tsw depicted in Figure 3.3, known as the carrier signal, with a slowlyvarying waveform called the modulating or reference signal. The intersectionbetween the two signals determines the state of the switches with values 1or 0, such that the addition of the states is equal to one. For example, S1 ina half-bridge converter obtains the value of 1, when the modulating signalhas a value larger than or equal to the carrier wave; S1 is equal to zero whenthe modulating signal has a value less than the carrier wave, respectively.The nature of the switching states for a half-bridge converter determines thewaveforms with the following assumptions [48]:• Each transistor and diode act as a short circuit with a small on-stateresistance ron in its conduction state, and act as an open circuit in itsblocking state16• Each transistor has no turn-off tailing current• Each diode has no turn-off reverse recovery current• Transitions from a conduction to a blocking state are instantaneousFigure 3.3: PWM SchemeA three-phase two-level AC/DC converter may also employ the PWMscheme with the assumption that the modulating signal, m(t), is a sinusoidalfunction with the required amplitude and frequency [48]. The AC voltage isoften then controlled by a closed-loop scheme that attempts to restrict theAC side current to directly follow the sinusoidal command. Equation 3.1demonstrates the controlled AC side voltage value, where ron is the on stateresistance of the IGBT.Vt(t) = m(t)VDC/2− roni(t) (3.1)173.1.3 Converter Switching ModelA simple switch model of an AC/DC half-bridge converter is easily trans-ferable to replicate the idea for a three-phase two-level AC/DC VSC. WhenS1 = 1 and the upper switch is closed, the AC terminal voltage is positiveand is equal to VDC/2. Alternatively, when S2 = 1, the lower switched isclosed and AC terminal voltage is equal to -VDC/2. It then follows thatthe AC/DC half-bridge converter can be mathematically characterized byEquations 3.2 and 3.3.S1(t) + S2(t) ≡ 1 (3.2)Vt(t) = (VDC/2)S1(t)− (VDC/2)S2(t) (3.3)A three-phase two-level VSC follows a similar principle, however the ACterminal side voltage for each phase is dependent on its modulation signaland phase current. In order to obtain a balanced three-phase AC side voltageand line current, ma(t), mb(t), mc(t) must establish a balanced three-phasesignal, commonly implemented by a closed-loop controller. Further analysisof the three phase modulation signals can be found in Chapter 5. Equations3.4-3.6 demonstrate the mathematical model of a three-phase two-level VSC.Vta(t) = ma(t)VDC/2− ronia(t) (3.4)Vtb(t) = mb(t)VDC/2− ronib(t) (3.5)Vtc(t) = mc(t)VDC/2− ronic(t) (3.6)3.1.4 Park TransformationThe Park Transformation is an often used mathematical technique to sim-plify control of three-phase systems. The main concept is to regulate ACwaveforms in a rotating reference frame, specifically a dq-frame, to ensurethat the steady state signals are DC. This in turn allows for compensatorswith simple structures, for lower dynamic orders, and for the possibility toachieve zero steady state tracking error with the inclusion of integral terms[50].18The conversion from a three-phase abc-frame, to a dq-frame, undergoesa transformation represented by Equation 3.7. The inverse transformationis given by Equation 3.8. This represents an a-phase to q-axis alignment.fd and fq are the projections of the rotating vector, fs, onto the orthogonaldq-axis, and it is assumed that for a balanced three-phase system, f0 mustequal zero. Figure 3.4 depicts the relation between the coordinate systems.fdfqf0 = 23sin(θ) sin(θ − 2π3 ) sin(θ + 2π3 )cos(θ) cos(θ − 2π3 ) cos(θ + 2π3 )121212fafbfc (3.7)fafbfc = sin(θ) cos(θ) 1sin(θ − 2π3 ) cos(θ − 2π3 ) 1sin(θ + 2π3 ) cos(θ + 2π3 ) 1fdfqf0 (3.8)Figure 3.4: Park Transformation19The remaining chapters of this thesis related to half-bridge convert-ers focus on the modeling and control of a three-phase two-level AC/DCand AC/DC/AC VSC. Figure 3.5 illustrates the circuit progression from anAC/DC half-bridge converter, to back-to-back three-phase two-level VSCsconnected by a DC voltage source (two series-connected DC-link capaci-tors). The fundamental operation of the AC/DC/AC VSC relies on thePWM procedure using a triangular waveform to generate twelve switchingsignals for each side of the half-bridge converter. The reference voltage isfound to be similar to that of a three-phase two-level VSC, where a nodeor ground, is located between the two DC-link capacitors and at the ACvoltage sources. All other operating mechanisms are identical to that of ahalf-bridge converter.Figure 3.5: AC/DC/AC Three-Phase Two-Level VSC3.2 Three-Phase Diode Rectifier3.2.1 TopologyThe main objective eof a full wave rectifier is to produce either a voltageor current that is purely DC, or to maintain a voltage or current valuewith a specified DC component [51]. A simplistic topology is the bridgerectifier with a RLC source load depicted in Figure 3.6. The rectifier containscomplementary diode pairs, such that D1, the upper left diode, and D4, thelower right diode, conduct simultaneously, as well as D3, the upper rightdiode, and D2, the lower left diode, which also conduct simultaneously.20Figure 3.6: Full Wave Bridge RectifierSimilar to the three-phase two-level half-bridge VSC, a full wave rectifiercan be extended to a three-phase rectifier, known as a six-pulse rectifier [51].For the purpose of this thesis, the analysis assumes an RL branch with aDC capacitor source load in the six-pulse rectifier. Voltage laws around anypath show that only one diode on the top bridge, D1 and D3 and D5, mayconduct at one time, similarly with the diodes on the bottom bridge, D2and D4 and D6. Consequently, D1 and D2, D3 and D4, and D5 and D6cannot conduct simultaneously.Figure 3.7: Six-Pulse Rectifier213.2.2 Operating PrincipleThree-phase rectifiers serve the common purpose of producing a DC voltageand current output for large source loads. The six diode currents and three-phase currents assist in the six combinations of line-to-line voltages thattransition every 60◦ with an operating fundamental frequency of the three-phase source equal to 6ω. Due to the topology of the system, ia, ib, andic represent the difference in the diode currents, iD, as demonstrated inEquations 3.9-3.11. The six-pulse rectifier conducts each diode current onethird of the time resulting in an average and rms diode current value, and acurrent source rms value, Equations 3.12-3.14. Additionally, the maximumreverse voltage across a diode is the peak line-to-line voltage.ia = iD1 − iD2 (3.9)ib = iD3 − iD3 (3.10)ic = iD5 − iD6 (3.11)ID,avg =13Io,avg (3.12)ID,rms =1√3Io,rms (3.13)Is,rms =√23Io,rms (3.14)The fundamental advantage of a six-pulse rectifier is its ability to producean average current from the AC side source equal to zero, thus eliminatingany complications associated with nonzero average source currents [51]. Asa result, a full wave rectifier produces inherently less ripple while creating aload current that is essentially DC.3.2.3 Diode ModelModeling the diode’s behaviour requires obtaining the value of the currentflowing from the anode to the cathode, and in turn, calculating the voltage22across the diode [52]. The device turn-on and turn-off conditions dependupon forward-biased and reverse-biased currents, such that the diode turnsoff when the current is less than or equal to zero [51]. The threshold voltagespecified by the user is the deciding factor for the turn on condition whenthe voltage across the diode is larger than or equal to this value. Ideally,the diode is a short circuit when it is forward-biased and is an open circuitwhen reversed-biased [51]. Figure 3.8 shows the graphical representation ofits behaviour.Figure 3.8: Ideal Diode Switch Behaviour (adapted from [51])The diode model characteristics in Figure 3.8 assume a piece-wise linearsystem. In high frequency applications, an ideal diode current will reach azero value without the dynamic characteristic of a negative reverse recov-ery current, thus avoiding an increase in the computational time requiredto reach the turn off state. Fast recovery diodes attain minuscule reverserecovery times of approximately less than 1 µs [51].23Chapter 4Matrix ExponentialIntegrator AlgorithmIn this chapter, a comprehensive analysis of a matrix exponential integratoralgorithm is proposed for the simulation of dynamic systems developed in [1].The algorithm relies on integration order, step-size, and system state equa-tions to capture dynamic changes of periodic excitation, switching events,sudden changes of loads and sources, and more [3]. Individual functionswithin the algorithm are rationalized for further insight into the character-istic behaviour that is expected within the computational process.When the exponential integrator algorithm is compared to other DSEDframeworks, specifically one found in [3], it has shown to simulate morerapidly due to pre-computed terms that reduce the number of calculationsrequired [1]. These results remain true when the exponential integrator al-gorithm does not utilize the adaptive-grained approach and strictly operateswith a matrix exponential integration scheme. Another advantageous com-ponent in the proposed method is its sequential intermediate integrationsteps between the discrete events which reduce the set of possible discretiza-tion step-sizes [1]. The exponential integrator algorithm is proficient inhandling both AC and DC sources and remains an optimal method acrossvarious circuit topologies. The efficacy of the proposed algorithm is laterverified through numerical examples of power system converters.244.1 State Space Equation ModelingPower converter systems are assembled through passive and active switchelements, controllers, load sources, and other elements that can be expressedin the form of a state equation. The piece-wise linear system expresses thetime derivatives of the independent states through matrix functions of theindependent states themselves, the input sources, and the set of outputsof the system. Four matrices, A, B, C, and D, are the state matrices thatdescribe the relationships between independent state and switching eventswithin the power converter system.x˙(t) = Ax(t) +Bu(t) (4.1)y(t) = Cx(t) +Du(t) (4.2)Equation 4.1 and 4.2 are unique, depending on the parameters of thepower converter systems, and are the initial building blocks required for thederivation of the matrix exponential algorithm.4.2 Algorithm FunctionsThe subsections below are a step-by-step approach to construct the proposedmatrix exponential integration algorithm and provide background infor-mation on the functions required for the complete exponential-integration-based circuit solver.4.2.1 Matrix Exponential DerivationSystem state matrices are used to construct the matrix exponential term forthe algorithm. The solution hinges on finding the initial state and the inputfor t ≥ 0. Assume,ddteAt = AeAt = eAtA (4.3)25Pre-multiplying e−At on both sides of Equation 4.1 yields,e−Atx˙(t) = e−AtAx(t) + e−AtBu(t) (4.4)e−Atx˙(t)− e−AtAx(t) = e−AtBu(t) (4.5)This implies,ddt(e−Atx(t)) = e−AtBu(t) (4.6)Integrating the above equation from 0 to t, the resulting exponential isgiven by,x(t) = eAtx(0) +∫ t0eA(t−τ)Bu(τ)dτ (4.7)If Bu(t) is a nonlinear function, it can be substituted for a generalizedterm, g(τ, u(τ)) [53, 54].x(t) = eAtx(0) +∫ t0eA(t−τ)g(τ, u(τ))dτ (4.8)For simulation purpose, Equation 4.8 can be expressed as an iterativefashion, where h is the step-size of the simulation.xn+1 = eAhxn +∫ t+hteA(t−τ)Bu(τ)dτ (4.9)4.2.2 Phi-Functions and Taylor Series ExpansionThe second term in Equation 4.9, the integral term, contains the Phi-function solution of the exponential integrator function [53, 54]. In [55],Lemma 5.1 states that a Taylor series expansion can be applied to a non-linear function, such as g(τ ,u(τ), to re-represent an initial value problem.Substituting a Taylor series expansion of g(τ ,u(τ) into the variation of con-stants formula and using the expression for φ-functions gives [55],xn+1 = eAhxn +∫ tn+htne(tn+h−τ)AB(u(τ))dτ (4.10)26The numerical approximation for simulation gives the solution,xn+1 = eAhxn +∞∑k=1φk(Ah)hkBd(k−1)undt(k−1)(4.11)The φ-function can be independently expressed by Equation 4.12.φk(Ah) =1(k − 1)!∫ 10e(1−θ)(Ah)θk−1dθ, k ≥ 1 (4.12)Computing the φ-functions is a challenging computation due to the ex-ponential form. [55] derives a recursive formula or Taylor series alternativeequation that can be employed when determining the outcome.φk(Ah) = Ahφk+1(Ah) +1k!, φ0(Ah) = eAh (4.13)φk(Ah) =l∑j=0(Ah)j(k + j)!(4.14)The recursive formula is a suitable option for large eigenvalues of A,however it begins to fail if the determinant of A approaches zero. TheTaylor series approach is accurate when the eigenvalues of A are small, butrequires a large amount of approximation terms for stiff systems.4.2.3 Pade´ ApproximationNow that the composition of the matrix exponential integrator algorithmis complete, the computation of the matrix exponential term, or naturalresponse term, eAh, requires the use of a numerically efficient and accuratemethod. Numerous algorithms have been proposed, however the majorityare of low numerical quality. Systems involving detailed switching transientsare best suited to use a Pade´ approximation based on a method derived in[56]. A Pade´ approximation is rational function representing the ratio oftwo power series, with the order of the method equal to the order of thenumerator and denominator in the series, (p,q) respectively. The rationalfunction is demonstrated in Equation 4.15 with the numerator Npq and the27denominator Dpq.Rpq = Dpq(Ah)−1Npq(Ah) (4.15)Npq =p∑k=0(p+ q − k)!p!(p+ q)k!(p− k)! (Ah)k (4.16)Dpq =q∑k=0(p+ q − k)!q!(p+ q)k!(q − k)! (−Ah)k (4.17)An example of a (3,3) order system to solve eAh is demonstrated inEquation 4.18.eAh =I +Ah/2 + (Ah)2/10 + (Ah)3/120I −Ah/2 + (Ah)2/10− (Ah)3/120 (4.18)A Pade´ approximation solver is considered a high-order solver specificallyfor power electronic converter transient simulation.4.3 Simulation FrameworkThe general framework of the matrix exponential integrator algorithm usesa variable step-size high-order solver based on the adaptation of computedmatrix and φ-function terms. Figure 4.1 demonstrates a flow chart of thefoundational aspects of the algorithm, using the matrix exponential integra-tion technique.28EnterEstimate step-size based on switching eventInitialize solution to homogeneous responseusing high-order Pade´ approximationCalculate φ1i=1Include higher order terms in estimatexn+1 = eAhxn +∑∞k=1 φk(Ah)hkB d(k−1)undt(k−1)i=i+1Calculate φiAccuracy threshold met?φi(Ah)hiB d(i−1)undt(i−1) < ϵtolMaximum number of terms reached?Decrease step-sizeEndYesNoYesNoFigure 4.1: Flow Chart of Matrix Exponential Integrator Algorithm WithVariable-Step and Variable Number of φ-Function Terms (adapted from [1])29The framework in Figure 4.1 demonstrates a variable step-size incorpo-rating accurate timing to update state variables using fixed discretization-order and varying φ terms up to the maximum order value. The algorithmoperates with a fixed-order Pade´ approximation solver, however the vari-able step-size for simulations is predicated for Equation 4.11, and on thoseinstances when switching events occur. The solver is a variable-step fixedhigh-order system, which is the better method compared to fixed-step low-order ODE solvers currently used in research similar to that found in [3].A variable-step fixed-order method is advantageous as switching events ofpower converter systems are only simulated at the exact points, decreasingthe total number of simulated points. Additional calculations between theseevents are considered redundant.The solver begins with an initial variable step-size, then truncates theanalytical series in Equation 4.11 up to the pth term, such that the localtruncation error (LTE) is less than the desired error tolerance ϵtol. TheLTE is the remaining term after and including (p+1 ) [1]. The LTE can beapproximated by Equation 4.19. If the maximum error has not convergedbased on the maximum number of terms, the step-size must be reducedto compute the φ terms again. The input, u, can be either an AC or DCquantity, however it is normal for the input to have a lower frequency (50Hz or 60 Hz) to facilitate a small number of φ terms, and to accuratelycapture semiconductor switching events. Since the algorithm uses a high-order Pade´ approximation solver to calculate the matrix exponential terms,a user specified tolerance can be carefully selected to reduce the error.ϵ =∞∑k=p+1φk(Ah)hkBdk−1undtk−1≈ φp+1(Ah)hp+1Bd(p)undt(p)(4.19)In addition to the variable-step requirement of the algorithm, Figure 4.2demonstrates the potential addition of an adaptive-grained pre-computationalapproach for the AGEI algorithm. If the step-size and switch event time ispre-determined, the power circuit topologies can be updated and stored be-30fore simulation execution. All the circuit topologies are found in advancewhen the the switch states are calculated and cycled through to updatethe system matrix. While there is a large amount of circuit topologies de-pendent on the number of switches, such that n switches has 2n topologies,complementary switches reduce the number of unique circuit configurations.31BeginPrecompute φ termsStart simulationEstimate step-size based on switch eventsReformat step-sizeSet intermediate step-size to largest non-zero orderApply variable-step, high-order variable φ termexponential integrator using pre-computed φ termsError tolerance met withinmaximum number of terms?Decrease intermediate step from remaining step-sizeh = h − hinterDecrease step-sizeSwitching event reached?h == 0Simulation end met?End simulationYesNoYesYesNoNoFigure 4.2: Flow Chart of AGEI Algorithm With Pre-Computations(adapted from [1])32Step-size pre-computations will contain a minimum value specified bythe user, and a maximum value determined by the switching frequency. Theintermediate points are simply different multiples of the step-size. If a cir-cuit operates with 2n topologies and m unique discretization step-sizes, thetotal pre-computed matrix exponential terms will be 2nm, and 2nmqmaxpre-computed φ terms. Some cases require large n and m values, conse-quently creating a large number of pre-computation of terms, motivatingthe addition of the adaptive-grained approach. The adaptive-grained algo-rithm immensely decreases the computational process through decreasingthe number of consecutive pre-defined step-sizes. The algorithm then alter-natively decomposes a single large step-size into multiple small sizes withvarious lengths. The computations between the minimum and maximumstep-size have the potential to be reduced by scaling factors 1 through 9, aseach magnitude step-size is now used for pre-computing the matrix expo-nential and φ terms.4.4 Numerical StabilityNumerical stability of the proposed φ-function exponential integrator isstudies in this subsection. Conventional explicit ODE solvers such as Runge-Kutta methods do not possess absolute stability. The following numericalstability analysis shows the exponential integrator solver does provide abso-lute stability (A-stability). Given the following ODE:x˙ = λx (4.20)Let R be the stability region and z be the set of all complex numberssuch that,z = hλ (4.21)|xn+1| < |xn| (4.22)If the homogeneous solution converges to zero, it is said to be stable,specifically A-stable if the stability region encompasses the entire left-hand33plane. Furthermore, the algorithm is L-stable if it is A-stable and if Re(λ)approaches negative infinity as the function approaches zero. The stabilityof the exponential integrator is evaluated by applying the homogeneous so-lution of Equation 4.11 and testing Equation 4.20.First,xn+1 = ehλxn (4.23)where |xn+1]|xn| = |ehλ| (4.24)A-stability analysis then requires that if for all Re(λ) < 0 in Equation4.24, the result must give|xn+1]|xn| < 1 (4.25)or,|xn+1| < |xn| (4.26)This confirms that the exponential integrator is indeed A-stable.Knowing that the exponential integrator is now A-stable, Equation 4.24can be used to assess the L-stability.Re(λ)→ −∞ : |ehλ| → 0 (4.27)Therefore, the exponential integrator is L-stable.It is desirable in power converter system simulation to use an A-stablemethod, as this guarantees convergence of the solution for either a fixed-step or variable-step approach. A disadvantage of A-stability however isthe presence of numerical oscillations due to switching instances. L-stablesimulation algorithms exhibit similar behaviour without the addition of nu-merical oscillations, thus making them more desirable solvers for simulationsassociated with fast transients and high switching frequencies.34Chapter 5Matrix ExponentialIntegrator Based ConverterModelingThis chapter focuses on the control modeling and implementation of AC/DCand AC/DC/AC three-phase two-level converters, as well as a six-pulse rec-tifier, being subject to the matrix exponential integrator algorithm. Thecontrol methods follow that found in [58, 59], focusing on PQ power andDC voltage control in the dq-frame. The technical aspects of the operat-ing principles are analyzed and further developed through a discretizationapproach similar to that in [60], generating the three-phase fundamental ref-erences and reference waveforms needed for the implementation of a ZOHmethodology. Passive switch event modeling is predicated on a similar modelin [60] to detect zero crossing events and to determine the numerical meth-ods appropriate for linear interpolation.Lastly, the theoretical operation is verified through simulations withinthe Matlab/Simulink environment. Three simulation case studies reviewthe numerical accuracy and efficiency of the models and control. Table 5.1,5.3, and 5.5 summarize the circuit and control parameters required for theindividual studies.355.1 Discrete Event Detection5.1.1 Discrete Control ModelsThe control system implemented for the simulation of three-phase two-levelconverters consists of a closed-loop power and DC voltage controller, whichgenerates the switching signals through a ZOH model. The PQ power andDC voltage controller follow the methodology outlined in [59], derived from adq-frame current mode controller to output the three-phase signals necessaryto construct the reference waveforms. These waveforms are compared witha carrier triangular waveform to determine the time of the switching eventsto achieve zero steady state error. The following sections provide a detailedunderstanding on the mathematical strategies and implementations of thecontrollers for the simulation and validation of the proposed algorithm.5.1.2 Control in dq-FrameThe two primary methods of real/reactive power control under the transfor-mation in dq-coordinates are voltage-mode control or current-mode control.A grid-imposed frequency VSC system, such as an AC/DC and AC/DC/ACthree-phase two-level converter, can be employed with a real/reactive powercontroller using current-mode control [59].Strict regulation of three-phase line currents is achieved by means of acontrol system utilizing the AC side terminal voltage. The real and reac-tive power is controlled by the phase angle and amplitude of the VSC linecurrent with respect to the PCC [59]. Figure 5.1 illustrates a current-modePQ power controller performed in the dq-frame using Park Transformation.The control methodology depends upon the simplicity of minimal controlloops, however there is no closed-loop control of the internal currents, cre-ating a higher possibility for control signal transients [59]. In addition, thesimplicity of the method facilitates the attainment of zero steady state errorwhen controlling variables that are DC quantities [59]. The synchronizationmechanism of the controller is achieved through a phase-locked loop (PLL).36Reference Compensators dq → abcabc→ dqVdq idqidq ref md mq mabcVsabc iLabcPQFigure 5.1: Current-Controlled Real/Reactive Power Controller in dq-Frame5.1.3 Mathematical OverviewFigure 5.1 presents the control design used for the simulation of an AC/DCthree-phase two-level VSC. The following section provides a technical overviewand detailed analysis of the components and variables in the schematic dia-gram.The mathematical derivation assumes an infinitely stiff AC system mod-eled by an ideal three-phase voltage source that is balanced, sinusoidal, andmaintains a relatively constant frequency [59]. P and Q are controlled bythe line current components, id and iq, where the feedback and feed-forwardsignals are first transformed to the dq-frame and processed by the compen-sators to produce control signals in the dq-frame. The controller is then ableto transform the signals in the abc-frame and feed them into the VSC.Equation 5.1 and 5.2 demonstrate the expressions for PQ power compo-nents in the dq-frame, measured at the point of common coupling (PCC),where Vdq and idq represent the voltage and current variables in the dq-framerespectively [59].P (t) =32[Vd(t)id(t) + Vq(t)iq(t)] (5.1)Q(t) =32[−Vd(t)iq(t) + Vq(t)id(t)] (5.2)37The dynamic model assumes that the three-phase AC system voltagecan be expressed as,Vsa(t) = Vpkcos(ω0t+ θ0) (5.3)Vsb(t) = Vpkcos(ω0t+ θ0 − ϕ) (5.4)Vsc(t) = Vpkcos(ω0t+ θ0 + ϕ) (5.5)where Vpk is the line-to-line neutral voltage (AC source peak voltage) andθ0 is zero for the purpose of this thesis. Applying a space phasor calculationto the source voltages, an equivalent equation for all three-phases is denotedas,V⃗s(t) = Vpkej(ω0t+θ0) (5.6)Using the voltage space phasor representation, the dynamics of the ACside are equivalent to,Ld⃗idt= −(R+ ron)⃗i+ V⃗t − V⃗s(t) (5.7)with i⃗t and V⃗t equal to idqejρ and Vdqejρ respectively. Combining Equa-tion 5.6 and 5.7 where fdq = fd + jfq, and decomposing the joint equationinto real and imaginary parts, it can be deduced thatLdiddt= (Ldρdt)iq − (R+ ron)id + V⃗d − Vˆscos(ω0t+ θ0 − ρ) (5.8)Ldiqdt= −(Ldρdt)id − (R+ ron)iq + V⃗q − Vˆssin(ω0t+ θ0 − ρ) (5.9)Introducing the variable ω = dρdt and allowing ρ to adopt a zero initialcondition, ω(t) is then equivalent to zero while ρ remains zero at all times.Equation 5.10 and 5.11 demonstrate the new conditions, resulting in two,decoupled, first order systems [59].Ldiddt= −(R+ ron)iq + V⃗d − Vˆscos(ω0t+ θ0) (5.10)Ldiqdt= −(R+ ron)iq + V⃗q − Vˆssin(ω0t+ θ0) (5.11)38One common case that regularly occurs is such when ω = ω0, and conse-quently ρ(t) = ω0t + θ0. The proper selection of ω and ρ(t) in the dq-framethen allow Equation 5.10 and 5.11 to take the form below.Ldiddt= Lω0iq − (R+ ron)id + V⃗d − Vˆs (5.12)Ldiqdt= Lω0id − (R+ ron)iq + V⃗q − Vˆs (5.13)5.1.4 DC Voltage ControlWith reference to Figure 5.1, Figure 5.2 shows the addition of DC volt-age regulation in current-mode control for the simulation of an AC/DC/ACthree-phase two-level VSC. If the PLL is in a steady state, the real andreactive power strictly depend upon the nature of Vd, Vq, id, and iq, and arecontrolled by these system states [59]. The following control scheme illus-trates the new, real power reference control signal with VDC input control.(·)2(·)2-1 k(s) Pmax/PminVDC - eV Pref+VDCrefPextFigure 5.2: DC Voltage Control SchemeThe introduction of idref and iqref is essential for the control system tofacilitate reference tracking at a rapid speed such that id ≈ idref and iq ≈iqref , and P ≈ Pref and Q ≈ Qref [59]. As a result, the real and reactivepower are independently controlled by their reference commands. Operatingin the steady state, the reference values are of DC quantity, translating to a39PQ/DC voltage controller possessing constant reference signals in dq-frame.idref =23VdPref (5.14)iqref = − 23VdQref (5.15)Equation 5.16 and 5.17 demonstrate the operating principle of the con-troller in a steady state condition with ω(t) = ω0.Ldiddt= Lω0iq(R+ ron)id + Vd − Vsd (5.16)Ldiqdt= −Lω0id(R+ ron)iq + Vq − Vsq (5.17)Focusing on the control of the DC voltage, Vtd and Vtq are contingenton,Vd =VDC2md(t) (5.18)Vq =VDC2mq(t) (5.19)where md(t) and mq(t) are the modulating signals produced in dq-frame.Due to the presence of Lω0, the dynamic behaviour of id and iq are coupled,requiring the result of the modulating signals to decouple the system whereud and uq will be the new inputs.md =2VDC(ud − Lω0iq + Vsd) (5.20)mq =2VDC(uq + Lω0iq + Vsq) (5.21)As the equations for diddt anddiqdt are now decoupled, the control plants inboth d and q-axis are identical, therefore the corresponding compensatorswill also be identical [59]. A compensator that is simple to optimize and oflow dynamic order can then operate with kdq(s) equal to a PI compensator.This enables tracking of a DC reference command as shown in the two blockdiagrams below.40kd(s)1Ls+(R+ron)idref ed ud id−Figure 5.3: Block Diagram idkq(s)1Ls+(R+ron)iqref eq uq iq−Figure 5.4: Block Diagram iqLetting kdq equal,kdq =kps+ kis(5.22)with proportional and integrator gains kp and ks, the loop gain is thenequivalent tol(s) =kpLss+ kikps+ R+ronL(5.23)The plant pole at s = -(R+ron)/L resides close to the origin and createsa narrow frequency window before the magnitude and phase begin to drop.The ripple effect is evident in the plant pole when it is first cancelled bythe compensator zero s = − kikp and the loop gain takes the form l(s) =kpLs .If τi is the time constant of the closed-loop system, and if kp = L/τi andki = (R + ron)/τi, then the closed-loop transfer function, l(s)/(1+l(s)),becomesid(s)idref= Gi(s) =1τis+ 1(5.24)Equation 5.22 and 5.23 indicate that the selection of ki and kp will in-fluence the response of id(t) and iq(t). Ideally τi is small for a fast current-control response but large enough such that the bandwidth of the closed-loop41control system is minuscule compared to the frequency of the VSC [59].5.1.5 Zero-Order Hold ImplementationThe pre-scheduling of discrete events requires the implementation of a ZOHwithin the controller for matrix exponential integrator simulation purposes.When a ZOH is applied to digital applications, the controller signal is fixedand constant during a switching period [48]. In each control cycle, sampling,turn on, and turn off actions are calculated based on the step-size and carrierwaveform, and the events are performed successively until the next controlcycle [3, 7]. In simulation of power converter systems, due to the presenceof discrete events, the step-size is simultaneously dependent upon the inte-gration of continuous states and the occurrences of discrete events [3]. Thestep-size previously computed by the algorithm, is now largely dependenton the time interval between two consecutive discrete events [3].Figure 5.5 demonstrates the variability of locating events using the ma-trix exponential framework with a PWM strategy for the control of the simu-lated VSCs. A detection approach samples the modulating signal and carrierwaveform at the current step-size and calculates the switching instances atthe beginning of each control cycle. The time at which the sampling eventsoccur can be detected prior, according to the user-specified sampling fre-quency. Pre-scheduling of discrete switch events can be calculated at thebeginning of each control cycle without the use of an iterative search method[3].42ZOHReference+Carrier SignalGenerated switching signal−Discrete Event Control TimingReference and Carrier Signal Scheduled timing of eventsFigure 5.5: Comparison of an Iterative and Pre-Scheduling Approach forSwitch Event Detection (adapted from [3])43For example, when a ZOH is used with a current-controlled PQ/DC volt-age controller, a three-phase two-level VSC contains six switching instancesper cycle. Based on Tsw, two switching times, t1 and t2, are computedthrough the peak value of the corresponding modulating signal within halfthe switching time. Each pre-scheduled switch event for the three-phasesystem is independently calculated without the influence of other modulat-ing phase signals. Equations 5.25-5.30 represent the general switching timesof three-phase two-level VSCs, later used for the simulation in Case StudyOne.t1a =Tsw2ma + 12(5.25)t2a = Tsw − t1a (5.26)t1b =Tsw2mb + 12(5.27)t2b = Tsw − t2b (5.28)t1c =Tsw2mc + 12(5.29)t2c = Tsw − t1c (5.30)5.1.6 Phase-Locked Loop DiscretizationA current-controlled PQ/DC voltage-controller in dq-coordinates requiresa synchronized mechanism modeling a nonlinear dynamic system, achievedthrough the PLL displayed in Figure 5.6.44H(s) ωmax/ωmin∫Vsabc/dqρFigure 5.6: PLL Control Block DiagramA PLL serves the purpose of regulating ρ at ω0t + θ0, however in view ofits nonlinear characteristics and under certain conditions, it may facilitatesmall perturbations about zero and will fall into a limit cycle [59]. Equation5.31 exhibits the nonlinear dynamic model of a PLL, where H(ρ) is a lineartransfer function.Ldρdt= H(ρ)Vˆssin(ω0t+ θ0 − ρ) (5.31)The prevention of a limit cycle can be solidified through a control lawplacing a limit on the maximum and minimum values of ω while assigningan initial value of ω(0) = ω0. If ωmax and ωmin attain values close to thatof ω0, the narrow variation of allowed ω values ensures a higher accuracy ofthe control in the dq-frame.The dynamic performance of a PLL is highly dependent upon the com-pensator, H(s), in order to achieve zero steady state error. This conditionforces the loop gain to include at least two integrators, and H(s) to consistof at least one integral term, having a pole at s=0. The remaining poles andzeros are determined solely on the basis of the PLL closed-loop bandwidth[59]. The compensator used for the simulation and validation of this thesisfollows that found it [59], given by Equation 5.32. The specifics of the designdepend on the values of Vpk, ω0, θ0, and ϕ.H(s) =685.42(s2 + 568, 516)(s2 + 166s+ 6889)s(s2 + 1508s+ 568516)(s2 + 964s+ 232, 324)[(rad/s)/V ] (5.32)45The simulation of the PLL for Case Study One and Case Study Twoin the Matlab environment requires a state matrix equation of the com-pensator. [60] derives a method for the formulation of state equations fordynamic systems and digital controls using what is called Successive Differ-entiation. First, the numerator and denominator of H(s) are expanded. Thecoefficient of the highest order in the denominator is then divided throughthe transfer function.H(s) =685.42s4 + 113780s3 + 3.94394e8s2 + 6.46856e10s+ 2.68445e12s5 + 2472s4 + 2254552s3 + 898394016s2 + 13079911184s(5.33)The new transfer function can be assigned the general form with numer-ator coefficients An and denominator coefficients Bn, where n = 5 [60].H(s) =An−1sn−1 +An−2sn−2 +An− 3sn−3 +An−4sn−4 +A0Bnsn +Bn−1sn−1 +Bn−2sn−2 +Bn−3sn−3 +Bn−4sn−4(5.34)An and Bn can be related through a matrix expression with mediatingvariables Cn given by Equation 5.35.1 0 0 0 0 0Bn−1 1 0 0 0 0Bn−2 Bn−1 1 0 0 0Bn−3 Bn−2 Bn−1 1 0 0Bn−4 Bn−3 Bn−2 Bn−1 1 0B0 Bn−4 Bn−3 Bn−2 Bn−1 1C0C1C2C3C4C5=AnAn−1An−2An−3An−4A0(5.35)The values for Cn are determined from the above matrix. By transferringin the function values, the coefficients are equivalent to,C0 = 0 (5.36)C1 = 685.42− 0C0 (5.37)C2 = 113780− 2472C1 − 2254552C0 (5.38)C3 = 3.94394e8− 2472C2 − 2254552C1 − 898394016C0 (5.39)46C3 = 3.94394e8− 2472C2 − 2254552C1 − 898394016C0 (5.40)C4 = 6.46856e10− 2472C3 − 2254552C2 − 898394016C1 − 132079911184C0(5.41)C5 = 2.68445e12−2472C4−2254552C3−898394016C2−132079911184C1−0C0(5.42)The final state matrix representation of the PLL compensator acquiresthe form of x˙ = Ax+Bu, y = Cx+Du.x˙1x˙2x˙3x˙4x˙4 =0 1 0 0 00 0 1 0 00 0 0 1 00 0 0 0 1−1 −123079911184 −898394016 −2254552 −2472x1x2x3x4x5+C1C2C3C4C5u(5.43)y =[1 0 0 0 0]x+[0]u (5.44)Completing the simulation in a Matlab environment requires a numericalmethod to calculate the output of the PLL, ρ, with either a variable or fixedstep-size. Methods in [52] have been used to calculate a solution at thepresent step-size using only those values obtained in the previous step-size[7, 12]. Forward Euler method is employed to simulate ρ which is requiredfor one of the inputs of the PQ/DC voltage controller in dq-coordinates.The general formation of Forward Euler is presented in Equation 5.28.yn+1 = yn + hf(tn, yn) (5.45)5.1.7 Controller DiscretizationThe current-controlled PQ/DC voltage controller in dq-frame is composedof four transfer functions. The simulation then requires the transformationusing the Successive Differentiation method found in [60], identical to thatused for the PLL transformation. However, both the d and q-axis utilizethe same compensators, so there will be a parallel transformation into theirstate space matrices. The first transfer function that undergoes Successive47Differentiation is kd(s) and kq(s).kd/q(s) =0.05s+ 0.815s=kps+ kis(5.46)Let A1 = 0.05, A0 = 0.815, B1 = 1, and B0 = 0 such that,kd/q(s) =A1s+A0B1s(5.47)The transfer function coefficients can then be related through a matrixform, [1 0Bn−1 1][C0C1]=[AnA0](5.48)whereC0 = A1 = 0.05 (5.49)C1 = A0 = 0.815 (5.50)The state space matrices now take the form,[x˙]=[0] [x]+[C1]u (5.51)y =[1]x+[0]u (5.52)The second compensator in the controller can be discretized by dividingout the denominator coefficients, leaving B1 = 1,l(s) =1/8e− 6s+ 1/8e− 6 (5.53)with A0 = 1/8e-6, A1 = 0, B0 = 1/8e-6, and B1 = 1.l(s) =A0B1s+B0(5.54)48The matrix representation is identical to that above,[1 0Bn−1 1][C0C1]=[AnA0](5.55)withC0 = An = A1 = 0 (5.56)C1 = A0 −Bn−1C0 = 1/8e− 6 (5.57)The final state matrix equations are,[x˙]=[−1/8e− 6] [x]+[1/8e− 6]u (5.58)y =[1]x+[0]u (5.59)The completed discretization of the controller transfer function is re-quired for simulation purposes for both Case Study One and Case StudyTwo. Forward Euler method is used to create the output signal necessaryto determine the relationships to current references generated by the PQand DC voltage control. The explicit method can be designed for stiff prob-lems with smaller stability error [3]. After the signals have been computed,the remainder of the controller produces two modulating signals in the dq-frame. Equation 3.8 reverts md(t) and mq(t) back into abc-coordinates,and are compared to the carrier signal of the VSC for simulation switchingpurposes.5.2 Passive Event ModelsPassive event modeling is an uncontrollable process which can only be de-tected by means of a system state crossing a user-defined threshold value.Specifically in diode modeling, the voltage across the device and the currentpassing through are indicative of a passive switching event. The voltagethreshold, Vcth, is the deciding factor for turn-on of the diode, such that Vf≥ Vcth. If the magnitude of the device current falls below zero, the diode49turns off [47]. When a passive event is detected at the next step-size, thesystem matrices can be updated or an iterative search can be employed tolocate the exact time of the zero crossing event [3]. The simulation of a six-pulse rectifier illustrates the diode modeling procedure through the followequations, where the diode state, ds, obtains logical values of 1 or 0 for aclosed and open state respectively. RD is the resistance across the diode,VD is the voltage across the diode, and iD is the diode current.RD = (ds)Ron + (1− ds)Roff (5.60)VD = (ds)Vf +RDiD (5.61)The simulation of the six-pulse rectifier is modeled based on the threeupper diode switches, D1, D3, and D5, with corresponding phases a, b, and crespectively. As the matrix exponential integrator algorithm has computedthe state currents and DC voltage at each simulation point in time, with theinitial diode resistance value, the three diode currents can be computed.iD1 =RD2RD1 +RD2ia − ds1 + ds2RD1 +RD2Vf − 1RD1 +RD2VDC (5.62)iD3 =RD4RD3 +RD4ib − ds3 + ds4RD3 +RD4Vf − 1RD3 +RD4VDC (5.63)iD5 =RD6RD5 +RD6ic − ds5 + ds6RD5 +RD6Vf − 1RD5 +RD6VDC (5.64)5.2.1 Iterative Zero Crossing DetectionAs introduced in the above section, the occurrences of passive events canbe described mathematically by a set of threshold equations to determinewhen a passive switching device has turned on or off. Since the diode voltageequation is an explicit function of the state currents, it is an explicit functionof the time variable. Therefore the event-locating of the passive event cansimplify to a root finding problem if desired [3]. The proposed iterativeroot finding algorithm for zero crossing detection of a six-pulse rectifier isthe Secant method. Considering the trade off between computational costs50required at each step-size and the rate of convergence, the Secant methodis more computationally reliable as it only requires the analysis of eitherthe diode voltage or the current threshold detection [3]. Equation 5.65demonstrates the premise of the method [61].xi+1 = xi − f(xi)(xi−1 − xi)f(xi−1)− f(xi) (5.65)At the initial simulation time interval the Secant method is performediteratively to detect whether the threshold conditions are met within theuser defined error tolerance. If more then one passive event is detectedduring an interval, all events will be located and stored. The system statematrices will be updated in order of the most urgent event to the least.5.3 Simulation and ValidationIn this section, case studies are presented to illustrate the effectiveness of thematrix exponential integration algorithm. The numerical accuracy and effi-ciency are compared to fixed and variable step-size Simulink/Simscape Elec-trical solvers with a step-size of 1 µs for the AC/DC and AC/DC/AC con-verter cases. The diode modeling case follows a similar procedure using thesame step-size of 1 µs. A small step-size of 1 µs is used by Simulink/SimscapeElectrical models to produce reference solutions for the case studies to val-idate the exponential model. Simulink/Simscape Electrical ode4 fixed-stepfixed-order solver provides a balance between order selection and accuracyfor a fair comparison between simulation algorithms. Each system is em-ployed with an absolute and relative error tolerance value of 10−4, Pade´ pand q orders of 3, and a φ function of order 4. The Simulink based modelsuse the Simulink/Simscape Electrical toolbox while the matrix exponentialintegrator algorithm is coded using Matlab script.The case studies demonstrate the matrix exponential integrator algo-rithm’s numerical efficiency for various circuit topologies with an AC inputsource, and provide general efficiency comparisons between the solvers. All51circuit topologies are simulated from t = 0 to 0.1 s. Confirmation of thecircuit’s parameters can be found in Table 5.1, 5.3, and 5.5. The numericalefficiency results are presented in the following case studies.5.3.1 Case Study OneCase Study One simulates an AC/DC three-phase two-level VSC, Figure3.2, with PQ control in the dq-frame, implementing controller discretiza-tion and ZOH-based digital control. The circuit and control parametersare presented in Table 5.1. The matrix exponential integration algorithm isemployed such that the determination of the timing for each discrete eventis computed only at multiples of the switching period, Tsw, operating withthree-phase abc inputs, Equations 5.66-5.68. State space formulation is de-rived in Equations 5.69-5.73 where the equations can be combined into thestate matrices, A and B. The inductor current and capacitor voltage valuesare subject to the position of the upper-position IGBTs, S1, S3, and S5. Val-idation of the algorithm’s accuracy is illustrated in Figure 5.10-5.14, wherethe model compares the matrix exponential algorithm to a fixed-step refer-ence Simulink-based model using a step-size of 5 µs and 1 µs. As shown inthe following figures, the waveforms produced by the exponential integratorsolver maintain a higher accuracy compared to the 5 µs and 1 µ reference so-lution, thus verifying the improvements of the proposed algorithm comparedto conventional solvers.52Table 5.1: AC/DC Three-Phase Two-Level ConverterCircuit Parameter Circuit Value Control Parameter Control ValueTs 1 µs f 60 HzFsw 10 kHz w0 377 rad/sTsw 1/Fsw w 2πfϕ 2π/3 wmin 376 rad/sCDC 1000 F wmax 378 rad/sC1 CDC kp 0.05C2 CDC ki 0.815IL0 0 A Pref1 2.5 MWL 0.1 mH Pref2 -2.5 MWR 7.5 mΩ Qref 1 MVARron 8.8 mΩVDC 1250 VVDC1 VDC/2VDC2 VDC/2Ua = Vpksin(ωt) (5.66)Ub = Vpksin(ωt− ϕ) (5.67)Uc = Vpksin(ωt+ ϕ) (5.68)dVDC1dt= −S1C1ia − S3C1ib − S5C1ic (5.69)dVDC2dt=(1− S1)C2ia +(1− S3)C2ib +(1− S5)C2ic (5.70)diadt=S1LVDC1 − (1− S1)LVDC2 +RLia − 1LUa (5.71)dibdt=S3LVDC1 − (1− S3)LVDC2 +RLib − 1LUb (5.72)dicdt=S5LVDC1 − (1− S5)LVDC2 +RLic − 1LUc (5.73)53Figure 5.7: Three-Phase AC Voltage Source InputFigure 5.8: Real and Reactive Power Closed-Loop Control54Figure 5.9: Real and Reactive Power Closed-Loop Control Comparison55Figure 5.10: AC/DC Three-Phase Two-Level VSC 5µs Proposed Algorithmand 5µs Simulink Comparison - Three Phase Inductor Currents56Figure 5.11: AC/DC Three-Phase Two-Level VSC 5µs Proposed Algorithmand 5µs Simulink Comparison - Combined DC Capacitor Voltage57Figure 5.12: AC/DC Three-Phase Two-Level VSC 1µs Proposed Algorithmand 1µs Simulink Comparison - Three Phase Inductor Currents58Figure 5.13: AC/DC Three-Phase Two-Level VSC 1µs Proposed Algorithmand 1µs Simulink Comparison - Combined DC Capacitor Voltage59Figure 5.14: Zoomed-In Waveform AC/DC Three-Phase Two-Level VSCProposed Algorithm and Simulink Comparison - Three Phase Inductor Cur-rents60The waveform comparisons confirm the numerical accuracy of the ma-trix exponential algorithm against conventional Simulink simulation solvers.The chosen discretization order of the Pade´ approximation and integrationterms supports the optimal step-size values over a short period of time in or-der to capture high frequency converter dynamics. Throughout the changesof switching states, the system remains nonstiff, however even with initialtransients, numerical stability is maintained.Numerical efficiency of the matrix exponential algorithm results aresummarized in Table 5.2. The algorithm is compared against industry-standard algorithms most commonly used to simulate power electronic sys-tems, thus validating potential contributions. An AC input source requiresthe solver to vary the number of integration terms for each switching cy-cle in order to meet the desired threshold values. The CPU efficiency timetherefore becomes incredibly significant. As shown below, Simulink fixed-step fixed-order solvers, ode4 (Classical Runge-Kutta), variable-step fixed-order solvers, ode45 (Dormand-Prince) and ode23t (mod.stiff/Trapezoidal),and variable-step variable-order solvers, ode15s (stiff/NDF), perform witha lower numerically efficiency to that of the variable-step matrix exponen-tial algorithm. The summarized results illustrates the algorithm’s flexibilitygiven various circuit parameters and effectiveness across complex topologies.Table 5.2: Computational Efficiency Comparison for Three-Phase Two-Level AC/DC VSCAlgorithm CPU Computation Time (s)Matrix Exponential Integrator 2.3749913 sSimulink ode4 1 µs 15.8125000 sSimulink ode45 32.2968750 sSimulink ode23t 32.4687500 sSimulink ode15s 48.7343750 s615.3.2 Case Study TwoThis study analyzes the performance of the matrix exponential algorithm foran AC/DC/AC three-phase two-level converter with PQ/DC voltage controlin the dq-frame. The simulation is subject to a voltage change at t = 0.02s, a real power change at t = 0.06 s, and a reactive power change at t =0.07 s. The state space model of the three-phase abc converter system forboth AC sides are given by Equations 5.74-5.76. The state equations arederived from Figure 3.5, and are used for calculations within the matrix ex-ponential algorithm. They can be rearranged to form state matrices A andB. The complete circuit and control parameters are summarized in Table 5.3.The benefits of the proposed algorithm, for a complex system with a largeA system state matrix, is dependent upon switch states using a step-size of 1µs. With large state space equations, it is expected that the accumulation ofcalculations each time a discrete event occurs, directly corresponds to a largeincrease in CPU time. However, due to the overall stability and efficiencyfor a nonstiff system, the matrix exponential algorithm achieves consistentlyaccurate results at a faster computational speed. Compared to Case StudyOne, an AC/DC/AC three-phase two-level converter operates with twelveswitches and 4096 (212) possible system topologies, creating a state A ma-trix of size 8 x 8 opposed to 5 x 5. The increase of eigenvalues imposes thepossibility of greater distribution between values, and whether they affectthe selection of order and step-size required for accurate results. It is as-sumed that the Pade´ approximation would be of high-order to guaranteeconvergence, yet the overall framework of the matrix exponential integratoralgorithm converges using the same order for both Case Study One and CaseStudy Two, Figures 5.17-5.19.62Table 5.3: AC/DC/AC Three-Phase Two-Level ConverterCircuit Parameter Circuit Value Control Parameter Control ValueTs 1 µs f 60 HzFsw 10 kHz w0 377 rad/sTsw 1/Fsw w 2πfϕ 2π/3 wmin 376 rad/sCDC 1000 F wmax 378 rad/sC1 CDC kp 0.05C2 CDC ki 0.815IL0 0 A Pref1 2.5 MWL1 0.1 mH Pref2 -2.5 MWL2 0.5 mH Qref 1 MVARR1 7.5 mΩ Vref1 350 VR2 7.0mΩ Vref2 VDCron 8.8 mΩVDC 1250 VVDC1 VDC/2VDC2 VDC/2Ua = Vpksin(ωt) (5.74)Ub = Vpksin(ωt− ϕ) (5.75)Uc = Vpksin(ωt+ ϕ) (5.76)dVDC1dt= −S11 + S21C1ia − S13 + S23C1ib − S15 + S25C1ic (5.77)dVDC2dt=(1− S11) + (1− S21)C2ia+(1− S13) + (1− S25)C2ib+(1− S15) + (1− S25)C2ic(5.78)dia1dt=S11L1VDC1 − (1− S11)L1VDC2 +R1L1ia1 − 1L1Ua (5.79)dib1dt=S13L1VDC1 − (1− S13)L1VDC2 +R1L1ib1 − 1L1Ub (5.80)dic1dt=S15L1VDC1 − (1− S15)L1VDC2 +R1L1ic1 − 1L1Uc (5.81)63dia2dt=S21L2VDC1 − (1− S21)L2VDC2 +R2L2ia2 − 1L2Ua (5.82)dib2dt=S23L2VDC1 − (1− S23)L2VDC2 +R2L2ib2 − 1L2Ub (5.83)dic2dt=S25L2VDC1 − (1− S25)L2VDC2 +R2L2ic2 − 1L2Uc (5.84)Figure 5.15: Three-Phase AC Voltage Source Input64Figure 5.16: Real/Reactive Power and DC Voltage Closed-Loop ControlComparison65Figure 5.17: AC/DC/AC Three-Phase Two-Level VSC 1µs Proposed Algo-rithm and 1µs Simulink Comparison - Three Phase Inductor Currents ACSide One 66Figure 5.18: AC/DC/AC Three-Phase Two-Level VSC 1µs Proposed Algo-rithm and 1µs Simulink Comparison - Three Phase Inductor Currents ACSide Two 67Figure 5.19: AC/DC/AC Three-Phase Two-Level VSC 1µs Proposed Algo-rithm and 1µs Simulink Comparison - Combined DC Capacitor Voltage68Figure 5.20: Zoomed-In Waveform AC/DC/AC Three-Phase Two-LevelVSC Proposed Algorithm and Simulink Comparison - Three Phase InductorCurrents AC Side One69Figure 5.21: Zoomed-In Waveform AC/DC/AC Three-Phase Two-LevelVSC Proposed Algorithm and Simulink Comparison - Three Phase InductorCurrents AC Side Two70Table 5.4 encapsulates the significant efficiency improvement of the ma-trix exponential algorithm. While the large system state equations do con-tribute to an increase in CPU computation time compared to an AC/DCthree-phase two-level VSC, the model still maintains a faster overall simula-tion time when evaluated against Simulink ode4, ode45, ode23t, and ode15s.Table 5.4: Computational Efficiency Comparison Three-Phase Two-LevelAC/DC/AC VSCAlgorithm CPU Computation Time (s)Matrix Exponential Integrator 8.5514098 sSimulink ode4 1 µs 15.6562500 sSimulink ode45 31.6875000 sSimulink ode23t 32.8750000 sSimulink ode15s 45.5000000 s5.3.3 Case Study ThreeCase Study Three introduces passive event detection modeling. The pur-pose of this study is to highlight the effectiveness of the matrix exponentialalgorithm when forced commutation (IGBT) switching and a ZOH digitalcontroller are not used in a power converter. The foundational model re-quired for a variable step-size iterative approach, to detect the exact timingof zero crossing events using the Secant method, is discussed in Section 5.2.1.The three-phase inputs are illustrated in Equations 5.85-5.87. State spacemodeling is derived in Equations 5.88-5.91 which are needed to form the statespace matrices. Figures 5.23 and 5.24 illustrate the accuracy of passive eventdetection against Simulink’s solver, ode4. Without the implementation ofiterative detection, the matrix exponential integrator algorithm maintainshigh accuracy and produces the expected results. The parameters foundin Table 5.5 used for this simulation produce a nonstiff system to furthervalidate the effectiveness of the proposed algorithm. Under these conditionsthe simulation framework is highly reliable. Figure 5.26 details a zero cross-ing detection without having employed an iterative solution from the matrix71exponential integrator algorithm. This method may delay the time at whichthe system states are updated when a threshold value is reached for eithercapacitor voltage or inductor currents. The current simulation techniquewill allow the inductor currents to drop below zero before turning off thediode and the capacitor voltage to cross the threshold value before turningon the diode. If an iterative solution is implemented, a smooth transitionbetween diode states would occur, similar to that shown in the Simulinkmodel. Chapter 6 proposes to improve the outcome of passive detectionby using the matrix exponential algorithm to achieve an optimal mode ofoperation.Table 5.5: Six Pulse RectifierCircuit Parameter Circuit ValueTs 1 µsf 60 Hzω 2πfCDC 0.05 FIL0 0 AL 10 mHR 5 ΩRoff 1 MΩRon 0.1 mΩVC0 20VVcth 0.3 VVf 0.3 VVpk 50 VUa = Vpksin(ωt) (5.85)Ub = Vpksin(ωt− ϕ) (5.86)Uc = Vpksin(ωt+ ϕ) (5.87)72dVDCdt=(RD2RD1+RD2ia − ds1+ds2RD1+RD2Vf − 1RD1+RD2VDC)C+(RD4RD3+RD4ib − ds3+ds4RD3+RD4Vf − 1RD3+RD4VDC)C+(RD6RD5+RD6ic − ds5+ds6RD5+RD6Vf − 1RD5+RD6VDC)C(5.88)diadt=(2Ua3L− Ub3L− Uc3L)+(−23LRD1RD2RD1 +RD2− 2R3L)ia+(13LRD3RD4RD3 +RD4)ib +(13LRD5RD6RD5 +RD6)ic+(23LRD1ds1 + ds2RD1 +RD2− 13LRD3ds3 + ds4RD3 +RD4− 13LRD5ds5 + ds6RD5 +RD6− 23Lds1 +13Lds3 +13Lds5)Vf+(23LRD1RD1 +RD2− 13LRD3RD3 +RD4− 13LRD5RD5 +RD6)VDC(5.89)dibdt=(2Ub3L− Ua3L− Uc3L)+(−23LRD3RD4RD3 +RD4− 2R3L)ib+(13LRD1RD2RD1 +RD2)ia +(13LRD5RD6RD5 +RD6)ic+(23LRD3ds3 + ds4RD3 +RD4− 13LRD1ds1 + ds2RD1 +RD2− 13LRD5ds5 + ds6RD5 +RD6− 23Lds5 +13Lds1 +13Lds3)Vf+(23LRD3RD3 +RD4− 13LRD1RD1 +RD2− 13LRD5RD5 +RD6)VDC(5.90)73dicdt=(2Uc3L− Ua3L− Ub3L)+(−23LRD5RD6RD5 +RD6− 2R3L)ic+(13LRD3RD4RD3 +RD4)ib +(13LRD1RD2RD1 +RD2)ia+(23LRD5ds5 + ds6RD5 +RD6− 13LRD1ds1 + ds2RD1 +RD2− 13LRD3ds3 + ds4RD3 +RD4− 23Lds5 +13Lds1 +13Lds3)Vf+(23LRD5RD5 +RD6− 13LRD1RD1 +RD2− 13LRD3RD3 +RD4)VDC(5.91)Figure 5.22: Three-Phase AC Voltage Source Input74Figure 5.23: Six Pulse Rectifier Proposed Algorithm and Simulink Compar-ison 1µs - Three Phase Inductor Currents75Figure 5.24: Six Plulse Rectifier Proposed Algorithm and Simulink Com-parison 1µs - Capacitor Voltage76Figure 5.25: Six Pulse Rectifier Proposed Algorithm - Diode Currents77Figure 5.26: Zoomed-In Waveform Six Pulse Rectifier Proposed Algorithmand Simulink Comparison 1µs - Three Phase Inductor Currents78Efficiency results presented in Table 5.6 compare the same Simulinksolvers used in Case Study One and Case Study Two. The proposed sim-ulation framework of the matrix exponential integrator algorithm achievesthe highest efficiency over fixed and variable step-size methodologies. Com-pared to the variable step-size solvers, the proposed algorithm operates witha CPU computational time on the approximate order of 7 times faster. Com-pared to fixed-step solvers, the algorithm operates on the order of 10 timesfaster.Table 5.6: Computational Efficiency Comparison Six Pulse RectifierAlgorithm CPU Computation Time (s)Matrix Exponential Integrator 0.2655774 sSimulink ode45 1.8437500 sSimulink ode15s 1.8593750 sSimulink ode23t 1.9531250 sSimulink ode4 1 µs 2.7656250 s79Chapter 6ConclusionDetailed semiconductor switching-based numerical simulation is crucial tothe analysis, design, and control of power electronic converter systems. Clas-sical numerical integration techniques commonly impose difficulties whensimulating complex power electronic systems. The need for locating switch-ing events accurately and varying simulation time steps accordingly in-fluences the speed and efficiency at which these methods perform due tothe consecutive step-size calculations required to update the system states.Thus, for numerically efficient simulation methods to be competitive inpower converter research, they must accurately capture the coexisting con-tinuous states and discrete events that determine the mode of operation.Fast transients that occur in the simulation is conventionally captured throughsmall fix time-step with low order integration, which causes heavy simulationburden and low simulation efficiency. Variable time-step simulation frame-works are regarded as a more efficient simulation methodology for dealingwith non-linearity, stiffness, and step-size adaptability.A novel matrix exponential integrator based algorithm presented in thisthesis provides a foundational framework for simulations of power convertersystems. The unique design allows for internal adaptations of the algorithmbest suited for computational efficiency and for system state representation,with discrete switching event detection to determine step-size length whilesimultaneously retaining fast simulation speed. Flexibility, φ-function orderselection, variable step-size, and L-stability have proven the algorithm to bewell-suited for power electronic systems, stiff or nonstiff [1]. The baseline ofthe proposed algorithm performs with a much faster simulation compared toconventional Simulink-based models. CPU times have been achieved up to8020 times faster, while still preserving accuracy over the simulation period. Itsadvantageous numerical properties make the matrix exponential integratoralgorithm an optimal simulation technique employable for various convertertopologies, and augments its range of applicability.6.1 Thesis Summary1. Chapter 1 provided an introduction and overview of the influentialfactors in numerical simulation that motivated this research. The ma-trix exponential algorithm was then introduced providing for a moreefficient and accurate methodology for electromagnetic transient simu-lation of power electronic converter systems. Research objectives wereoutlined, followed by the major contributions and impact that thisresearch would have on power converter simulation techniques. Theorganization of the thesis was presented and sequenced in a mannerto facilitate understanding.2. Chapter 2 provided a comprehensive literature review and focused onimportant topics, such as simulation frameworks, numerical solvers,digital control and discrete switch event detection. Further analysesdove into practical implementation ideas, along with future researchdirections and improvements that the matrix exponential algorithmsought to obtain. The motivation for a more stable, accurate, andefficient simulation algorithm of power electronic converter system wasadditionally enhanced.3. Chapter 3 introduced three-phase converter topology and their oper-ating principle. The structure, function, and purpose of the individualconverter components were addressed with a derivation of their math-ematical model. Park Transformation was established for the dq-framecontrol required for the simulation and validation of the topology.4. Chapter 4 presented the proposed matrix exponential integrator algo-rithm with the addition of an adaptive-grained implementation tech-81nique. A detailed mathematical description of the individual functionsthat contributed to the effectiveness of the algorithm was provided.Additionally, the numerical A-stability and L-stability of the algorithmwas illustrated.5. Chapter 5 provided mathematical controller discretization, using Suc-cessive Differentiation along with Matlab/Simulink-based case studies,for further validation and insight into the theoretical methodology.Reasoning for parameter selection and compensator design was in-cluded due to the influence on closed-loop PQ/DC voltage control forthe generation of switching signals. The performance of the algorithmwas evaluated through assessment of accuracy and stability, with acorresponding comparison to Simulink/Simscape Electrical solvers.6. Chapter 6 has concluded with a general summary. Recommenda-tions have been provided to further improve future applications ofthe methodology.6.2 Future WorkThis work has investigated the numerical capabilities possessed by the ma-trix exponential integration algorithm. Numerical efficiency analysis hasvalidated the superiority of the proposed algorithm. An adaptive selectionof time steps would facilitate the optimization of the given method throughpre-computations of matrix exponential and φ-functions demonstrated inFigure 4.2. The proposed addition maintains the integrator’s numerical sta-bility while simultaneously decreasing the required computation time whenoperating under variable step-size and φ-function terms. Three-phase VSCswould benefit from this method as accuracy is neither affected by high/lowswitching frequencies nor by high/low fundamental frequency. The expo-nential integrator then retains the ability to effectively capture the corre-sponding circuit topology matrix by cycling through the pre-calculated 2nswitching topologies when discrete switching events occur. Additionally, asthe range of possible step-sizes are known, all combinations of Ah terms82can be used to pre-compute the matrix exponential and φ terms before thesimulation loop starts. Additional work on passive event detection wouldemploy the Secant method for an iterative search or an oversampling methodto compensate for small delays currently produced by the matrix exponen-tial integrator algorithm. Since event detection cannot be obtained as anexplicit function of time, the integration process can be computationallycomplex. Oversampling methods have the potential to detect the exact tim-ing of zero crossing events through parallel computation. Optimization ofthe exponential integrator with an iterative search would require further ex-amination.Other recommendations for future work include the investigation of thelimitations of matrix exponential integration techniques for controller re-finement, and the simultaneous integration of passive and discrete switchingdevices in a circuit topology. The employment of a discrete state eventdriven approach assumes linearity between switching events, however thismay not always be the case. Quadratic or exponential curves may be presentbetween points at which the system states are updated, and further advance-ments of the matrix exponential integrator algorithm would be required toaccurately capture this behaviour. Through the exploration of controller be-havior, adaptations of the proposed algorithm can strengthen understandingand aid in the optimization of the simulation framework.83Bibliography[1] J. Paull, C. Knickle, N. Lofroth, L. Wang, and W. Li, “Adaptive-GrainedExponential Integrator Algorithm for Efficient Simulation of Power Con-verter Systems,” IEEE Transactions on Power Delivery., Submitted March16, 2024.[2] J. Ju, B. Shi, Z. Yu, Y. Zhu, and Z. Zhao, “Backward Discrete StateEvent-Driven Approach For Simulation of Stiff Power Electronic Systems,”IEEE Access., vol. 9, pp. 28573-28581, 2021.[3] Y. Zhu, Z. Zhao, B. Shi, and Z. Yu, “Discrete State Event-Driven Frame-work With a Flexible Adaptive Algorithm for Simulation of Power ElectronicSystems,” IEEE Transactions on Power Electronics., vol. 34, no. 12, pp.11692-11705, March 2019.[4] B. De Keppler, L.A. Dessaint, K. Al-Haddad, and H. Nakra, “A com-prehensive approach to fixed-step simulation of switched circuits,” IEEETransactions on Power Electronics., vol. 17, no. 2, pp. 216-224, March2002.[5] M.O. Faruque, V. Dinavahi, and W. Xu, “Algorithms for the account-ing of multiple switching events in digital simulation of power-electronicsystems,” IEEE Transactions on Power Delivery., vol. 20, no. 2, pp. 1157-1167, April 2005.[6] G. Migoni, E. Kofman, F. Bergero, and J. Fernandex, “Quantization-based simulation of switched mode power supplies,” Simulation: Transac-84tions of the Society for Modeling and Simulation International., vol. 91, no.4, pp. 320-336, 2015.[7] B. Li, Z. Zhao, Y. Yang, Y. Zhu, and Z. Yu, “A Novel Simulation Methodfor Power Electronics: Discrete State Event Driven Method,” CES Transac-tions on Electrical Machines and Systems., vol. 1, no. 3, pp. 273-282, Sept.2017.[8] E. Kofman, “Discrete event simulation of hybrid systems,” SIAM Jour-nal on Scientific Computing., vol. 25, no. 5, pp. 1771-1797, 2004.[9] V. Acary, O. Bonnefon, and B. Brogliato, “Step-sizing numerical sim-ulation of switched circuits within the nonsmooth dynamical systems ap-proach,” IEEE Transactions on Computer-Aided Design of Integrated Cir-cuits and Systems., vol. 29, no. 7, pp. 1042-1055, June 2010.[10] D. Efimov, A. Polyakov, and A. Aleksandrov, “Discretization of homo-geneous systems using Euler method with a state-dependent step,” Auto-matica., vol. 109, no. 108546, July 2019.[11] C.A. Kennedy and M.H. Carpenter, “Diagonally Implicit Runge–KuttaMethods for Stiff ODEs,” Applied Numerical Mathematics., vol. 146, pp.221-224, 2019.[12] T. Noda, K. Takenaka, and T. Inoue, “Numerical integration by the2-stage diagonally implicit Runge-Kutta method for electromagnetic tran-sient simulations,” IEEE Transactions on Power Delivery., vol. 24, no. 1,pp. 390-399, Jan. 2009.[13] G.D. Irwin, D.A. Woodford, and A. Gole, “Precision simulation of PWMcontrollers,” Proc. Of International Power Systems Transients Conference.,2001.85[14] A.M. Luciano and A.G.M. Strollo, “A fast time-domain algorithm forthe simulation of switching power converters,” IEEE Transactions on PowerElectronics., vol. 5, no. 3, pp. 363-370, July 1995.[15] B. Shi, Y. Chen, K. Chen, J. Ku, and Z. Zhao, “Event-driven approachwith time-scale hierarchical automation for switching transient simulationof SIC-based high-frequency converter,” IEEE Transactions on Circuits andSystems I: Regular Papers., vol. 68, no. 11, pp. 4746-4759, Nov. 2021.[16] J. Zheng, Z. Zhao, Y. Zeng, B. Shi, and Z. Yu, “An event-driven real-time simulation for power electronics systems based on discrete hybrid step-size algorithm,” IEEE Transactions on Industrial Electronics., vol. 70, no.5, pp. 4809-4819, May 2023.[17] B. Shi, Z. Zhao, Y. Zhu, Z. Yu, and J. Ju, “Discrete state event-drivensimulation approach with a state-variable-interfaced decoupling strategy forlarge-scale power electronics systems,” IEEE Transactions on IndustrialElectronics., vol. 68, no. 12, pp. 11673-11683, Dec. 2021.[18] E. Kofman and S. Junco, “Quantized-state systems: A DEVS approachfor continuous system simulation,” Transactions of the Society of ComputerSimulation International., vol. 18, no. 3, pp. 123-132, 2001.[19] C. Dufour, J. Mahseredjian, and J. Belanger, “A combined state spacen-odal method for the simulation of power system transients,” IEEE Powerand Energy Society General Meeting., 2011.[20] T. Geyer, G. Papafotiou, and M. Morari, “On the Optimal Switch-Mode DC-DC Converters,” Springer-Verlag., Berlin, Heidelberg, pp. 342-356, 2004.[21] C. Sreekuma and V. Argwal, “A hybrid control algorithm for voltageregulation in DC–DC boost converter,” IEEE Transactions on Industrial86Electronics., vol. 55, no. 6, pp. 2530-2538, May 2008.[22] J. Alcala, V. Cardenas, and N. Visairo, “Input-output linearization con-trol for a three-phase AC/DC/AC voltage source converter for bidirectionalpower flow,” 12th IEEE International Power Electronics Congress., 2010.[23] T. Mingxing, Y. Hong, and Y.D. Sheng, “Research on coordinate trans-formation of the three-phase circuit,” Telkomnika Indonesian Journal ofElectrical Engineering., vol. 11, no. 8, 2013.[24] S.S. Lee, Y.E. Heng, and M.A. Roslan, “Finite control set model pre-dictive control of nine-switch AC/DC/AC Converter,” IEEE InternationalConference on Power and Energy., Nov. 2016.[25] C. Pan and J. Shieh, “New space-vector control strategies for three-phase step-up/down AC/DC converter,” IEEE Transactions of IndustrialElectronics., vol. 47, no. 1, pp. 25-35, Feb. 2000.[26] C. Restrepo, T. Konjedic, F. Flores-Bahamonde, E. Vidal-Idiarte, J.Calvente, and R. Giral, “Multisampled digital average current controls ofthe versatile buck–boost converter,” IEEE Journal of Emerging and SelectedTopics in Power Electronics., vol. 7, no. 2, pp. 879-890, Dec. 2019.[27] Z. Yang, Y. Zhong-dong, H. Tao, and H. Qiu, “A novel voltage bal-ance control strategy of AC-DC-AC Converter,” International Conferenceon Energy and Environment Technology., 2009.[28] M.N. Peretz and S. Ben-Yaakov, “Time-domain design of digital com-pensators for PWMDC-DC Converters,” IEEE Transactions on Power Elec-tronics., vol. 27, no. 1, pp. 284-293, Jan. 2012.[29] A.G. Perry, F. Guang, L. Yan-Fei, and P.C. Sen, “A design method forPI-like fuzzy logic controllers for DC–DC converter,” IEEE Transactions on87Industrial Electronics., vol. 54, no. 5, pp. 2688-2696, Oct. 2007.[30] L. Qu, D. Zhang, and Z. Bao, “Output current-differential controlscheme for input-series–output-parallel-connected modular DC–DC convert-ers,” IEEE Transactions on Power Electronics., vol. 32, no. 7, pp. 5699-5711, Sept. 2017.[31] Y. Matsuda, S. Kato, Y. Wakasa, and R. Adachi, “State-feedback event-triggered control using data-driven method,” 61st Annual Conference of theSociety of Instrument and Control Engineers., Oct. 2022.[32] Y. He and F.L. Luo, “Sliding-mode control for DC–DC converters withconstant switching frequency,” IEEE Proceedings–Control Theory and Ap-plications., vol. 153, no. 1, pp. 37-45, July 2011.[33] K. Vrdoljak, N. Peric, and I. Petrovic,“Sliding mode based load-frequencycontrol in Power Systems,” Power Systems Research., vol. 80, no. 5, pp.514-527, 2010.[34] S. Chattopadhyay and S. Das, “A digital current-mode control tech-nique for DC–DC converters,” IEEE Transactions on Power Electronics.,vol. 21, no. 6, pp. 1718-1726, Nov. 2006.[35] V.R. Dinavahi, M.R. Iravani, and R. Bonert, “Real-time digital simula-tion of power electronic apparatus interfaced with digital controllers,” IEEEPower Engineering Review., vol. 21, no. 8, pp. 63-64, Oct. 2001.[36] W.J. Cao, Q. Bi, X. Liu, C. Lim, and Y.C. Soh, “Frequency-domaintransfer function of digital multirate controller with current estimator,”IEEE Transactions on Control Systems Technology., vol. 13, no. 1, pp.131-137, Jan. 2005.[37] G. Zhou and J. Xu, “Digital average current controlled switching DC–DC88converters with single-edge modulation,” IEEE Transactions on Power Elec-tronics., vol. 25, no. 3, pp. 786-793, Sept. 2010.[38] B. Wang, X. Yu, and X. Li, “ZOH discretization effect on higher-ordersliding-mode control systems,” IEEE Transactions on Industrial Electron-ics., vol. 55, no. 11, pp. 4055-4064, Nov. 2008.[39] J. Xu, G. Zhou, and M. He, “Improved Digital Peak Voltage PredictiveControl for Switching DC-DC Converters,” IEEE Transactions on Indus-trial Electronics., vol. 56, no. 8, pp. 3222-3229, Aug. 2009.[40] G. Bonanno and L. Corradini, “Digital Predictive Current-Mode Con-trol of Three-Level Flying Capacitor Buck Converter,” IEEE Transactionson Power Electronics., vol. 36, no. 4, pp. 4697-4710, April 2021.[41] S. He, Y.J. Hung, and R.M. Nelms, “A digital predictive current modecontroller using average inductor current,” IEEE Energy Conversion Congressand Exposition (ECCE)., 2014.[42] H.S. Youn, J.S. Park, K.B. Park, J. Baek, and G.W. Mon, “A DigitalPredictive Peak Current Control for Power Factor Correction With Low-Input Current Distortion,” IEEE Transactions on Power Electronics., vol.31, no. 1, pp. 900-912, Jan. 2016.[43] N. Johansson, L. Angquist, and H.P. Nee, “An Adaptive Controller forPower System Stability Improvement and Power Flow Control by Meansof a Thyristor Switched Series Capacitor (TSSC),” IEEE Transactions onPower Systems., vol. 25, no. 1, pp. 38-391, Feb. 2010.[44] Q. Fen, R.M. Relms, and J.Y. Hung, “Posicast-based digital control ofthe buck converter,” IEEE Transactions on Industrial Electronics., vol. 53,no. 3, pp. 759-767, May 2006.89[45] L. Roessetoo and S. Buso, “Digitally-Controlled Single-Phase Single-Stage AC/DC PWM Converter,” IEEE Transactions on Power Electronics.,vol. 18, no. 1, pp. 326-333, Jan. 2003.[46] N. Watson and J. Arrillagag, “State Variable Analysis,” in Power Sys-tems Electromagnetic Transients Simulation. IET, 2003. ch. 3, pp. 42-53.[47] H.W. Dommel and S. Battacharya, “Transient analysis of control sys-tems,” in Electromagnetic Transients Program (EMTP): Theory Book. Bon-neville Power Administration, 2003, ch. 8, pp. 13.1–13.3.[48] A. Yazdani and R. Iravani, “DC/AC Half-Bridge Converter,” in Voltage-Sourced Converters in Power Systems Modeling, Control, and Applications.Wiley-Blackwell, Hoboken, 2010, ch. 2, pp. 23-25.[49] A. Yazdani and R. Iravani, “Two-level, Three-Phase Voltage SourcedConverter,” inVoltage-Sourced Converters in Power Systems Modeling, Con-trol, and Applications. Wiley-Blackwell, Hoboken, 2010, ch. 5, pp. 115-125.[50] A. Yazdani and R. Iravani, “Space Phasors and Two Dimensional Frames,”in Voltage-Sourced Converters in Power Systems Modeling, Control, and Ap-plications. Wiley-Blackwell, Hoboken, 2010, ch. 4, pp. 69-107.[51] D. Hart, “Full Wave Rectifiers,” in Power Electronics. The McGraw-Hill Companies, Inc., 2010, ch. 4, pp. 111-170.[52] H.W. Dommel and S. Battacharya, “Switches,” in Electromagnetic Tran-sients Program (EMTP): Theory Book. Bonneville Power Administration,2003, ch. 10, pp. 10.1-10.4.[53] A.H. Al-Mohy and N.J. Higman, “Computing the action of the matrixexponential, with an application to exponential integrators,” SIAM Journalon Scientific Computing., vol. 33, no. 2, pp. 488-511, 2011.90[54] C. Wang, X. Fu, P. Li, and J. Wu, “Accurate dense output formula forexponential integrators using the scaling and squaring method,” AppliedMathematics Letters., vol. 43, pp. 101-107, 2015.[55] B.V. Minchev and W.M. Wright, “A review of exponential integratorsfor first order semi-linear problems,” Numerics., no. 2, 2005.[56] B. Skaflestad and W.M. Wright, “The scaling and modified squaringmethod for matrix functions related to the exponential,” Appl. Numer.Math., vol. 59, no. 3–4, pp. 783–799, Mar. 2009.[57] U.M. Ascher and C. Greif, “A First Course in Numerical Methods,”Philadelphia, PA, USA: SIAM, 2011.[58] A. Yazdani and R. Iravani, “Control of Half-Bridge Converter,” inVoltage-Sourced Converters in Power Systems Modeling, Control, and Ap-plications. Wiley-Blackwell, Hoboken, 2010, ch. 3, pp. 48-68.[59] A. Yazdani and R. Iravani, “Grid-Imposed Frequency VSC System:Control in dq-Frame,” in Voltage-Sourced Converters in Power SystemsModeling, Control, and Applications. Wiley-Blackwell, Hoboken, 2010, ch.8, pp. 204-244.[60] N. Watson and J. Arrillagag, “Analysis of continuous and discrete sys-tems, in Power Systems Electromagnetic Transients Simulation. IET, 2003.ch. 2, pp. 11-34.[61] S.C. Chapra and R.P. Canale, “Open Methods,” in Numerical Methodsfor Engineers. McGraw-Hill Companies, Inc., 2010, ch. 6, pp. 142-173.91
Thesis/Dissertation
eng
Electrical Engineering
Vancouver : University of British Columbia Library
Attribution-NonCommercial-NoDerivatives 4.0 International
Graduate
Wang, Liwei
Matrix exponential integrator algorithm for accurate and efficient simulation of power converter systems
Text