Bridging Scales With the Shift Frequency: Frequency-adaptive simulation of multiscale transients in power systems

©SHUTTERSTOCK.COM/LUMPPINI

One of the most remarkable features of power and energy systems pertains to their extremely far-reaching scales, which are unique in the field of engineering. In continental Europe, for example, the synchronous ac power system integrates capitals from Lisbon to Warsaw and from Athens to Copenhagen. Reaching farther out to the north, high-voltage dc (HVdc) links connect to Sweden and the United Kingdom. And as in other regions of the world, those grids are being further developed to integrate more renewables. This only adds to the diversity of technologies in the power grids and contributes to the wide range of timescales involved. Those range from electrothermal interactions in the range of minutes over electromechanical transients within seconds, down to microseconds or even faster for the electromagnetic transients of traveling waves.

Although the coverage of all those scales in an efficient and accurate manner would be very desirable for power system analysis, two distinct types of tools have emerged to address distinct timescales. There are tools based on phasors, which are complex functions to represent amplitude, frequency, and phase. Such tools treat the natural sinusoidal ac voltages and currents by tracking their envelopes. From Figure 1, it can be appreciated that tracking the red-lined envelope requires significantly fewer sampling points compared with those that would be needed to accurately track the naturally occurring blue-lined sinusoidal waveform of carrier frequency ${f}_{c}$. Reducing the number of sampling points leads to a reduced computational effort in the discrete time domain. This makes phasors popular for studying phenomena such as electromechanical transients, where observation of the envelopes is illustrative.

Figure 1. Envelope versus natural waveform tracking.

With ac waveforms missing, phasors are not applicable to the modeling of dc power systems. Thus, in tools for the accurate representation of the natural waveforms occurring in dc systems, instantaneous signals are the basis. The same applies when tracking the natural waveforms in ac systems. This approach is followed in programs of type EMT (electromagnetic transients). Those tools were, as the name readily suggests, developed with a focus on electromagnetic transients.

For users looking for a holistic viewpoint, it is desirable to have the opportunity to analyze diverse transients without the need to switch among different tools. However, is it possible to unify and integrate the established and distinctly different approaches, based on phasors versus instantaneous signals, within one modeling methodology? The answer is “yes,” and the key to the solution of a multiscale simulation is the introduction of the shift frequency. In the frequency-adaptive simulation of transients (FAST) method, the shift frequency sets the reference for other frequencies to relate to. Waveforms are then dealt with at the resulting frequency differences. Thus, for a zero shift, a 50-Hz sinusoid oscillates at 50 Hz. But the sinusoid appears at standstill for a shift frequency of 50 Hz, and this is because 50 Hz - 50 Hz = 0 Hz. As such, FAST merges the modeling approaches based on phasors and instantaneous signals.

Evolving Methods and Needs

The value of multiscale simulation is best understood in the context of emerging, existing tools and needs. In the 1960s, the development of the EMT-type program began with the purpose of finding time responses of electromagnetic transients in arbitrary single- or multiphase networks with lumped and distributed parameters (as described in Dommel 1969 and Woodford et al. 1983). Key to the success of EMT-type programs has been the arrangement of the difference equations obtained through numerical integration of the differential equations that describe lumped inductive and capacitive network elements. The difference equations are arranged in such a way that lumped inductive and capacitive network elements are modeled through so-called companion models that involve only resistive elements and current sources. Connecting the companion models then gives the overall network model of the discrete time domain. A nodal analysis is performed to obtain the nodal conductance matrix and affiliated vectors, with the purpose of representing the overall network at discrete time steps.

As opposed to the nodal analysis used in EMT-type programs, modeling through dynamic-state equations leading to the state-space form has been the most popular method when it comes to the use of phasors that represent voltages and currents. The state-space form lends itself to stability analysis as a main application of such phasor-based models (as explained in Kundur 1994 and Maksimovic et al. 2001). In this context, it is important to distinguish between quasi-static and dynamic phasor calculus. In the quasi-static formulation, the behavior of network lines, cables, and transformers is described by a set of algebraic equations. As differential equations are omitted, the electromagnetic transients of the network itself are neglected. Such tools are sometimes referred to as transient stability (TS) simulators.

Experience has shown (as discussed in Vega-Herrera et al. 2021) that neglecting electromagnetic transients in the network is not appropriate when considering the opportunities brought about by fast power-electronic controls in inverter-dominated grids. As a remedy, differential equations are considered for the modeling of network elements in dynamic phasor calculus. This leads to more accurate results and was shown to be appropriate for the study of controls of inverter-dominated grids. However, an efficient solution for the description of the propagation of electromagnetic waves along transmission lines has not been available in simulators based on dynamic phasors. Obviously, it is possible to approximate the distributed parameters of a transmission line by a sequence of connected ${\pi}$-sections of lumped elements to capture essential characteristics, such as the inductances and capacitances per-unit length. The accurate and efficient implementation of the solution to the wave equation first proposed by d’Alembert has, however, been available through companion models in EMT-type programs.

The Shift Frequency as a Simulation Parameter

The high accuracy of EMT-type programs is achieved thanks to a combination of the application of proper numerical integration in the development of companion models and the use of algorithms that allow for an exact capturing of switching events. As integration methods, the trapezoidal method or a weight-averaged method combining trapezoidal and backward-Euler methods have been popular. Despite not originally developed for covering electromechanical transients, EMT-type programs do allow for the description of such slower transients too. In general, a larger time-step size may be used for numerical integration when studying waveforms that change at a lower rate.

Although a practical upper bound of the time-step size is given by the need to maintain a satisfactory accuracy of the numerical integration, it is also insightful to consider the time-step size limits given by Shannon’s sampling theorem. According to this theorem, a waveform is sampled without distortion due to aliasing as long as the waveform bandwidth ${f}_{\max}$ is lower or equal to the Nyquist frequency ${f}_{\text{Ny}}$. The Nyquist frequency is defined as half the sampling rate, i.e., one over twice the time-step size, where the time-step size coincides with the interval between two sampling points, as indicated in Figure 1. Thus, the natural waveform of a sinusoid of 50 Hz may be sampled without aliasing at a time-step size of up to 10 ms, i.e., at a minimum sampling frequency of 100 Hz. In EMT simulation practice, the time-step size used to track the 50-Hz ac sinusoid will generally need to be much smaller than that, however, due to requirements of accuracy of the numerical integration involved. Even if a steady state is simulated, a sufficiently small time-step size must be chosen to allow for representation of the natural waveforms at the carrier frequency, for which 50 and 60 Hz are the typical values in power systems. To overcome the constraints regarding the desirable increase of the time-step size in the presence of an ac carrier waveform, the shift frequency was introduced.

2D Setting of the Shift Frequency and Time-Step Size

The shift frequency complements time-step size and constitutes a novel simulation parameter. Thanks to the addition of the shift frequency, a 2D setting of simulation parameters is made possible, as indicated in Figure 2.

Figure 2. Toward 2D parameters in time and frequency domains.

The shift frequency ${f}_{s}$ serves as the reference from which other frequencies are measured. In general, it is now the frequency difference with respect to the shift frequency, i.e., ${f}{-}{f}_{s}$, that is applicable when considering Shannon’s sampling theorem. In this sense, there is no aliasing up to ${f}_{\max}{-}{f}_{s} = {f}_{\text{Ny}}$. For ${f}_{s} = {50}{\text{ Hz}}$, a 50-Hz carrier may thus be sampled at a theoretically infinite time-step size without any aliasing appearing. This is because the 50-Hz carrier appears to be at a standstill for this setting.

Realizing Frequency Shifting

In a continuous time domain, the ac sinusoidal carrier may be represented by a signal of ${A}\,{\cos}\,{\left({2}{\pi}{f}_{c}{t} + {\Phi}\right)}$. Using the complex exponential function, an alternative representation is given by ${A} / {2}\,{\exp}{\left({j}{\left({2}{\pi}{f}_{c}{t} + {\Phi}\right)}\right)} + {A} / {2}\,{\exp}\,{\left({-}{j}{\left({2}{\pi}{f}_{c}{t} + {\Phi}\right)}\right)}$. As this alternative readily illustrates, the real oscillatory waveform is obtained by the superposition of a phasor rotating counterclockwise at a positive frequency ${f}_{c}$ and a phasor rotating clockwise at a negative frequency ${-}{f}_{c}$.

However, shifting by ${f}_{s}$ toward decreasing frequencies can be meaningful only if the signal to be shifted contains just positive frequencies. For the ac sinusoidal carrier, this would imply a complex extension of the carrier through an imaginary quadrature part to yield ${A}\,{\cos}\,{\left({2}{\pi}{f}_{c}{t} + {\Phi}\right)} + {j}{A}\,{\sin}\,{\left({2}{\pi}{f}_{c}{t} + {\Phi}\right)}$, i.e., ${A}\,{\exp}\,{\left({j}{\left({2}{\pi}{f}_{c}{t} + {\Phi}\right)}\right)}$. The resulting complex exponential function describes a phasor rotating counterclockwise at positive frequency ${f}_{c}$. Now, frequency shifting of the carrier by ${f}_{s} = {f}_{c}$ yields ${f}_{c}{-}{f}_{s} = {0}{\text{ Hz}}$ and so eliminates the oscillation. In general, to support frequency shifting in multiscale simulation of ac power systems, all waveforms are modeled through analytic signals, which are composed of the original real signals and imaginary quadrature components. The imaginary part of an analytic signal relates to the real part through the Hilbert transform the application of which is illustrated in Figure 3. Shifting of the analytic signal is then possible, as shown in Figure 4.

Figure 3. The modification of the frequency spectrum by including the imaginary quadrature part.

Figure 4. The effect on the frequency spectrum when the shift frequency equals carrier frequency.

Just as with the time-step size, the shift frequency may be constant or kept variable during a simulation. A shift frequency of 0 Hz is recommended when dc transients are observed. For slow transients visibly modulating the ac carrier as well as for electromagnetic transients above the ac carrier frequency, a shift-frequency setting of ${f}_{c}$ is recommended. Depending on the situation observed, locational and temporal modifications of the shift frequency are possible. The observation of events, such as contingencies or the monitoring of spectra, could be used as input of methods aimed at setting the shift frequency in an adaptive manner.

Recent research on power-electronic dominated ac microgrids has shown a time-step size of 50 ms at ${f}_{s} = {f}_{c}$ to be practical during transients such as those triggered by typical load changes. Instead, a time-step size of about 50 ns is commonly needed for corresponding studies when using instantaneous signals. As such, the number of time-steps is 1,000 times more than it is without the shift frequency being available. This is indicative of the efficiency gains made possible through the application of frequency shifting.

Network Modeling

As the size of the power grid to be modeled can vary greatly, the capability to efficiently and accurately represent both lumped and distributed parameters is important. As such, a variety of local phenomena, from slow-changing transients to electromagnetic waves traveling along transmission lines, are of interest in the consideration of diverse scales. Because companion models used in EMT-type programs have been shown to be effective in representing lumped elements and in offering solutions to the wave equation, this approach is also compatible with FAST-type multiscale simulation.

Lumped Parameters

In EMT-type programs, the difference equations derived from differential equations describing inductive and capacitive elements are arranged to describe resistive elements in parallel with current sources. As mentioned in the previous section, it has been a success factor of the EMT-type simulation that, at each time step, the solution process is only concerned with resistive networks and sources. With the introduction of FAST, the same general principle of using companion models is followed. Instead of just resistances, the lumped elements are now described by impedances, and the sources involve complex instead of real signals (as documented in Strunz et al. 2006 and Zhang et al. 2010).

The emergence of a companion model for an inductor is illustrated in Figure 5. The FAST-type companion model shows a complex admittance ${\underline{G}}_{L}$ and a complex history current source ${\underline{\eta}}_{L}$, which refers to the past time step. The EMT-type companion model is a special case of the more generic FAST-type formulation and is obtained for a shift-frequency setting of 0 Hz. In this special case, the imaginary parts in the FAST-type companion models disappear. Thus, the complex admittance ${\underline{G}}_{L}$ becomes a real-valued conductance ${G}_{\text{L}}$.

Figure 5. Companion models of inductive lumped element in EMT- and FAST-type simulators.

Distributed Parameters

The parameters of transmission systems are given in per-unit length because the distributed effects of the medium are noticeable. For the purpose of illustrating basic principles, the modeling of a lossless single-phase transmission line with constant parameters is considered. As indicated in Figure 6, the line is assumed to be of length l, and the inductance and capacitance per-unit length are given by L’ and C’, respectively. For such a line, wave equations can be readily formulated.

Figure 6. Companion models of lossless line in EMT- and FAST-type simulators.

The solution to the wave equations describes propagations of traveling waves along the transmission line. The EMT-type solution is achieved using the companion model shown on the left. The companion model is applicable as long as the time-step size is smaller than the time interval ${T}_{\text{wp}}$ which the waves need to travel between both line ends. At both line ends, the currents are calculated by the division through the wave impedance plus a current source component that considers the traveling wave coming from the other side. Because ${T}_{\text{wp}}$ exceeds the time-step size, both line ends of the companion model are not topologically coupled for the length of the time-step size.

The FAST-type multiscale model is valid for any time-step size. If the traveling time ${T}_{\text{wp}}$ becomes lower or equal to the time-step size, then a ${\pi}$-section model is seamlessly inserted to establish the topological coupling during the time-step interval. The EMT-type model is a special case of the FAST-type model for a shift frequency of ${f}_{s} = {0}{\text{ Hz}}$ and ${T}_{\text{wp}}$ exceeding the time-step size. Thus, as for the modeling of lumped parameters, FAST offers a generalization of the EMT-type companion model. The shift frequency, as a novel simulation parameter, enables the flexibility to cover multiscale transients, which includes the EMT-type solution as a special case.

Application

Multiscale modeling and simulation integrates the virtues of dynamic phasor calculus and EMT-type programs within one unified framework and also becomes available through tools such as CloudPSS or when combined with simulator PSCAD (with a further description given in Rupasinghe et al. 2023). As such, multiscale modeling also bridges the scopes. This becomes evident from studies elaborated upon in the following.

Wind Energy Conversion System

According to the Global Wind Energy Council, 2023 was the first year for which the addition of global wind power capacity exceeded 100 GW. Even higher annual numbers of installation are predicted over the coming years. Contributions from wind energy conversion systems (WECSs) are essential when it comes to enhancing the share of renewable sources in electric power generation. WECS technology has become ever more sophisticated.

The development of WECSs relies on the integration of diverse and complementary technologies. The corresponding subsystems are exemplarily illustrated in Figure 7 for a WECS where electric power generation is facilitated via a doubly fed induction generator (DFIG). As a part of the function of the mechanical subsystem, the translational kinetic energy associated with the wind is partially converted into rotational kinetic energy through the action of the turbine, thus resulting in a mechanical torque on the drivetrain. The latter consists of a low-speed shaft, high-speed shaft, and gearbox in between.

Figure 7. A WECS based on doubly fed induction generator (DFIG) technology. PCC: point of common coupling; PWM: pulsewidth modulation; VSC: voltage-sourced converter.

Via the high-speed shaft, a mechanical torque is applied to the rotor of the DFIG, which serves to convert rotational kinetic energy into electric energy to be processed in the electrical subsystem. In this configuration, the DFIG is a wound-rotor induction machine where the stator windings are directly connected to the grid via a transformer, while the rotor windings are connected to the grid via two back-to-back voltage-sourced converters (VSCs) with a dc link capacitor in between, a filter, and the transformer. The crowbar circuit involves resistors that may be inserted by controlling thyristors in situations where it is necessary to limit rotor currents for the purpose of protection.

Due to the ability of the back-to-back VSCs to control the rotor-side voltages and currents at frequencies other than the grid frequency, the DFIG-based WECS is suitable for variable-speed wind turbine application. Also, thanks to the back-to-back VSCs, controlled reactive power exchanges over the point of common coupling (PCC) with the grid are possible.

The objective of the rotor-side control scheme is to regulate both the active power and the reactive power on the stator side. In this context, the vector control approach is widely used for the rotor-side converter. It is performed in a stator flux dq-reference frame, leading to a decoupled consideration of active power and reactive power via the rotor current ${\boldsymbol{i}}_{\text{abcr}}$. The stator flux vector may be calculated using stator current ${\boldsymbol{i}}_{\text{abcs}}$ and stator voltage ${\boldsymbol{v}}_{\text{abcs}}$. The active power on the stator side is controlled by adjusting the rotor speed ${\omega}_{\text{r}}$, for example, by following an optimal value so that the wind turbine operates around a maximum power point. The reactive power on the stator side is controlled by means of rotor current regulation. The grid-side control scheme maintains the dc bus voltage ${v}_{dc}$ and adjusts the reactive power flowing through the filter toward the grid. The grid-side converter current ${\boldsymbol{i}}_{\text{abcg}}$ is measured for this purpose. Outputs of the rotor- and grid-side controls are the pulsewidth modulation signals for the respective VSCs.

To avoid disconnection of the DFIG during grid faults, crowbar protection is widely used. When the dc bus voltage or the rotor current exceeds threshold values, the rotor-side converter is blocked and the crowbar circuit is activated. The DFIG rotor currents flow through the crowbar instead of the rotor-side converter. When the grid fault is cleared and the rotor current decays to a safe value, the crowbar is removed. The rotor-side converter resumes its operation. The key to this protection technique is limiting the high currents and providing a bypass in the rotor circuit via the crowbar.

Along with the growth of wind power capacity available in the grid, WECSs are known to influence power system transients over a wide range of frequencies. Both low-frequency transients, such as subsynchronous oscillations, and high-frequency transients, e.g., those caused by short circuits, can appear. EMT-type programs are well suited for studying the behavior of high-frequency electromagnetic transients related to WECSs. If it is of interest to study both low- and high-frequency transients within the same study of WECSs, FAST offers an efficient and accurate solution, as discussed hereafter (with simulation data taken from Xia et al. 2020).

Initially, a WECS is assumed to operate close to steady state with a generated active power of roughly 1.5 MW. The time-step size ${\tau}$ is set to 8 ms, and the shift frequency ${f}_{s}$ is at 60 Hz in the ac parts of the circuit. The envelope of the ac current is used to represent the steady state. The phase-a current flowing through the rotor-side converter is shown in Figure 8, while the active and reactive power flows into the grid are given in Figure 9.

Figure 8. The phase-a current of the rotor-side converter. (a) EMT-type model. (b) Natural and envelope waveforms in FAST-type model. The nonbold line depicts a natural waveform, while the bold lines represent envelope waveforms.

Figure 9. The active power and reactive power at the PCC. (a) EMT-type model. The solid lines represent active power, while the dashed lines depict reactive power. (b) FAST-type model. The solid lines represent active power, while the dashed lines depict reactive power.

At ${t} = {0.25}{\text{ s}}$, a three-phase-to-ground fault occurs at the PCC. Electromagnetic transients are triggered. As a consequence of the fault, the rotor-side VSC is blocked, and the crowbar circuit is activated so that the DFIG rotor currents flow through the crowbar instead of the rotor-side VSC.

To focus on those details, natural waveforms are tracked at ${\tau} = {10}$ ns and ${f}_{s} = {0}{\text{ Hz}}$. The fault is assumed to be cleared at ${t} = {0.4}{\text{ s}}$. With the rotor currents quickly returning to tolerable values thereafter, the crowbar resistors are disconnected and the rotor-side VSC resumes control of active and reactive power flows. After approximately ${t} = {0.8}{\text{ s}}$, electromechanical transients become increasingly dominant. As such, FAST returns to envelope tracking at ${f}_{s} = {60}{\text{ Hz}}$ in the ac parts of the circuit, and the time-step size ${\tau}$ is set to 8 ms again.

Figures 8 and 9 show the curves from EMT- and FAST-type simulation results. No differences can be seen during periods of tracking natural waveforms. Only when FAST employs envelope tracking at large time-step sizes are there obvious differences in how the information is represented because the EMT-type simulation always follows natural waveforms.

Wide-Area Energy System

Beyond WECSs, it is important to understand multiscale transients of other power and energy resources as well as their systemic interactions through the power grid. The following case study illustrates the value of such modeling and simulation for a situation in China with regard to grid development for increasing the share of renewable power. In this context, the CloudPSS program makes use of the FAST method, and it is widely used for such purposes and beyond.

The layout of a power system model describing a fictitious extended West–East China regional interconnected system is depicted in Figure 10. The test system consists of 12,609 buses with 248 synchronous generators, 1,864 transmission lines, 778 transformers, and 571 loads. The connection of the Western and Eastern 50-Hz ac subsystems relies on the ±320-kV HVdc system shown in the center. The HVdc converter stations are based on VSC technology.

Figure 10. West–East China regional interconnected ac/dc system.

The three large-scale renewable power stations, denoted by G1, G2, and G3, involve multiple wind and solar parks. In steady state, those three stations, G1, G2, and G3, are able to provide a real power of the order of 1200, 800, and 600 MW, respectively. In the observed scenario, a three-phase-to-ground fault happens at bus NJ at time point ${t} = {0.5}{\text{ s}}$, following an initial stage of steady state. Such a fault triggers electromagnetic transients. The fault is cleared after five cycles, followed by a recovery process. At roughly ${t} = {1}{\text{ s}}$, fast electromagnetic transients have largely damped out, making it practical to track envelope waveforms thereafter.

An EMT-type program may be used for analyzing the fault and the recovery. Because envelope waveforms are not available, the time-step size to be used in the EMT-type program would need to be rather small. This is necessary even throughout periods of transients of low frequencies; for example, in an advanced stage of recovery. The resulting computational effort reduces the scope of practical application in a wide-area energy system analysis. The issue is addressed by the flexibility of FAST in that both the tracking of envelope and of natural waveforms is supported. During stages of low-frequency electromechanical transients or at close to the steady state, the time-step size is chosen commensurate with envelope tracking. During stages of expected fast electromagnetic transients following faults and during the initial stage of the recovery, natural waveform tracking is pursued with a time-step size as in EMT-type programs.

The simulation results are compared in Figures 11 and 12. There are no differences observed in accuracy between the multiscale simulation and the corresponding EMT-type simulation for both currents and voltages. With the availability of FAST, an efficient and accurate multiscale transients simulation of a wide-area system is possible. The system scales may no longer be considered as a limitation to practical application.

Figure 11. The phase-a current of ST-NJ transmission line. (a) EMT-type model. (b) Natural and envelope waveforms of FAST-type model. The nonbold lines depict a natural waveform, while the bold lines represent an envelope waveform.

Figure 12. The phase-a voltage of the ST side of ST-NJ transmission line. (a) EMT-type model. (b) Natural and envelope waveforms of FAST-type model. The nonbold lines depict a natural waveform, while the bold lines represent an envelope waveform.

Conclusions

The past decades have seen the transformation of power and energy systems into a heterogeneous entity involving increasingly diverse technologies. This process is largely driven by the objective of moving toward a zero-carbon society and by the need for energy security. The transformation is still very much ongoing and is set to continue for decades to come, further increasing the multiscale character of power and energy systems. With FAST, there is a method at hand for bridging those scales in modeling and simulation. As the term FAST suggests, additional flexibility is attained by a complementary focus on the frequency domain. The shift frequency, as a key modeling parameter, determines the reference to which other frequencies refer. It may be set equal to the ac carrier frequency to effectively support dynamic phasor calculus. It may also be set to zero when dc transients are of interest. In the case of a shift frequency equal to zero, EMT-type modeling, with its particular virtues, appears as a special case. EMT-type modeling, of course, supports the adaptive setting of the time-step size regarding the sampling of waveforms in the time domain. However, a true multiscale character, with all its benefits in modeling and simulation, is only achieved by the simultaneous consideration of parameter settings in the time and frequency domains. It is this opportunity of 2D parameter setting that is at the heart of FAST, and that is effective and practical in the analysis of heterogeneous resources and wide-area power and energy systems. As such, the ability of bridging the scales also supports the ability of supporting the ongoing process of analyzing and enhancing the capabilities of power and energy systems.

Acknowledgment

Supported by Deutsche Forschungsgemeinschaft (DFG) under Grant 410829522.

For Further Reading

H. W. Dommel, “Digital computer solution of electromagnetic transients in single-and multiphase networks,” IEEE Trans. Power App. Syst., vol. PAS-88, no. 4, pp. 388–399, Apr. 1969, doi: 10.1109/TPAS.1969.292459.

D. A. Woodford, A. M. Gole, and R. W. Menzies, “Digital simulation of DC links and AC machines,” IEEE Trans. Power App. Syst., vol. PAS-102, no. 6, pp. 1616–1623, Jun. 1983, doi: 10.1109/TPAS.1983.317891.

P. Kundur, Power System Stability and Control. New York, NY, USA: McGraw-Hill, 1994.

D. Maksimovic, A. M. Stankovic, V. J. Thottuvelil, and G. C. Verghese, “Modeling and simulation of power electronic converters,” Proc. IEEE, vol. 89, no. 6, pp. 898–912, Jun. 2001, doi: 10.1109/5.931486.

J. Vega-Herrera, C. Rahmann, F. Valencia, and K. Strunz, “Analysis and application of quasi-static and dynamic phasor calculus for stability assessment of integrated power electric and electronic system,” IEEE Trans. Power Syst., vol. 36, no. 3, pp. 1750–1760, May 2021, doi: 10.1109/TPWRS.2020.3030225.

K. Strunz, R. Shintaku, and F. Gao, “Frequency-adaptive network modeling for integrative simulation of natural and envelope waveforms in power systems and circuits,” IEEE Trans. Circuits Syst. I, Reg. Papers, vol. 53, no. 12, pp. 2788–2803, Dec. 2006, doi: 10.1109/TCSI.2006.883864.

P. Zhang, J. R. Marti, and H. W. Dommel, “Shifted-frequency analysis for EMTP simulation of power-system dynamics,” IEEE Trans. Circuits Syst. I, Reg. Papers, vol. 57, no. 9, pp. 2564–2574, Sep. 2010, doi: 10.1109/TCSI.2010.2043992.

J. Rupasinghe, S. Filizadeh, D. Muthumuni, and R. Pavari, “A multi-solver framework for co-simulation of transients in modern power systems,” Electric Power Syst. Res., vol. 223, Oct. 2023, Art. no. 109659, doi: 10.1016/j.epsr.2023.109659.

Y. Xia, Y. Chen, Y. Song, and K. Strunz, “Multi-scale modeling and simulation of DFIG-based wind energy conversion system,” IEEE Trans. Energy Convers., vol. 35, no. 1, pp. 560–572, Mar. 2020, doi: 10.1109/TEC.2019.2953893.

Biographies

Kai Strunz (kai.strunz@tu-berlin.de) is with Technische Universität Berlin, 10587 Berlin, Germany.

Ying Chen (chen_ying@tsinghua.edu.cn) is with Tsinghua University, Beijing 100084, China.

Yue Xia (yue.xia@cau.edu.cn) is with China Agricultural University, Beijing 100083, China.

Digital Object Identifier 10.1109/MELE.2023.3320487

Electrical Machines in Electromagnetic Transient Simulations: Focusing on efficient and accurate models

Erfan Mostajeran, Navid Amiri, Seyyedmilad Ebrahimi, Juri Jatskevich

©SHUTTERSTOCK.COM/PHOTO SMILE

Electrical machines are extensively used in our everyday life. On the one hand, this can be seen in the rapid growth of generation from renewable energy sources such as wind, hydropower, tidal, etc. On the other hand, at the energy utilization and consumption end visible to most people, we are also witnessing revolutionary changes in many sectors, such as the electrification of all types of transportation, i.e., electric vehicles, industry-wide initiatives for more-electric aircraft and more-electric ships, military vehicles and defense systems, industrial automation, industrial and personal robots, medical devices and instruments, flying drones, electronic toys, and the multitude of household appliances and devices, all of which are designed and built with electric motors of various types and sizes.

Introduction

In many applications, electrical motors and generators are controlled using switching power electronic converters that enable efficient and flexible operation and optimize energy conversion. Engineers and researchers develop computer models of electrical machines and conduct countless offline and real-time simulations to design and operate those systems and devices. Such simulations may be run on a small system with just a few machines and their converters and controllers or on a large system together with many other components and subsystems. Just as electrical motors and generators enable many products and technologies that improve human life, models of such electrical machines for studies in electromagnetic transient (EMT) simulators also make it possible to envision, design, and analyze/operate such products and devices effectively.

The models of electrical machines may be broadly categorized depending on their underlying purpose and level of fidelity. The first group includes the finite-element analysis (FEA)- or magnetic equivalent circuit (MEC)-based models, where the dimensions of equations may reach tens to hundreds of thousands, which are determined by the discretization of the machine’s geometry and the flux paths into the finite-element mesh. For the MEC models, this would be the reluctance segments whose number may still be significant on the order of several hundreds or thousands. Although these models can provide excellent accuracy necessary for design-oriented studies and optimization on an individual machine while considering the geometry and material properties to a significant level of detail, they are computationally costly.

The second group includes the so-called lumped-parameter coupled-circuit models, where the number of electrical windings in the machine determines the size of the equations. Since the total number of windings is typically small (e.g., one per phase on the stator/rotor, and if needed, a field winding plus damper windings, etc.), the resulting models have a much smaller number of equations and may be efficiently implemented in most EMT simulators for simultaneous studies of machines with power electronic converters as well as power systems with many generators and motors. It may be estimated that for one FEA-type model, thousands of lumped-parameter coupled-circuit models are implemented and simulated every day worldwide using various offline and real-time EMT simulators. Therefore, this article focuses on the lumped-parameter coupled-circuit models of electrical machines and their interfacing methods in EMT simulators.

Electrical Machine Characteristics and Modeling

The energy conversion process in electrical machines is one of the purest, cleanest, and most efficient known to humankind and is based on Maxwell’s equations for electromagnetic fields and Lorentz force. Electrical machines perform energy conversion based on the interaction between the magnetic fields and current-carrying conductors to generate torque. Generally, an electrical machine consists of a stator and a rotor circuit, each generating its own magnetic field through external sources, such as excitation, permanent magnets (PMs), magnetization through an external magnetic field, or induction. Different types of electrical machines exist that differ depending on the method by which the stator and rotor generate their magnetic fields. These electrical machines are commonly rotating devices (i.e., their rotor is in a circular motion), although linear motors/generators also exist where the motion is in a linear direction. Nevertheless, the modeling and analysis of both rotating and linear electrical machines are greatly similar, where one is carried out using circular coordinates and angles, and the other is carried out in linear coordinates and displacements. Electrical machines that are responsible for generating most of the electricity today and converting energy to mechanical work in various applications mainly include rotating synchronous and induction machines. Therefore, this article will mainly discuss the modeling of these machines in EMT simulation programs.

Coupled-Circuit-Phase-Domain Models

The lumped-parameter coupled-circuit models of electrical machines may be obtained directly by considering the magnetically coupled windings as rotor-position-dependent inductors and Kirchhoff’s circuit laws that consider the actual connection of all coils and conductors. This approach gave rise to the so-called coupled-circuit-phase-domain (CCPD) models that consider each stator (and rotor) phase winding in their respective physical variables and phase coordinates. Figure 1 illustrates examples of three-phase induction, field-wound synchronous, and PM synchronous machines (PMSMs) and their simplified coupled-circuit representations. According to this approach, a typical three-phase induction machine may be represented by three inductors for the stator and three inductors for the rotor, all of which are magnetically coupled, and the coupling between the stator and rotor windings also changes with the rotor position. The squirrel-cage induction motor may be easily represented by shorting the rotor windings in Figure 1(a), and the doubly-fed induction machine may be implemented if the rotor circuit terminals are not shorted and made available instead. The number of equivalent rotor windings in Figure 1(a) may be increased accordingly to represent more complicated rotors, such as squirrel-cage deep-bar conductors or two layers of conductors.

Figure 1. Typical three-phase rotating ac machines and their lumped-parameter coupled-circuit representation. (a) An induction machine. (b) A field-wound synchronous machine (SM) with equivalent damper windings on the rotor. (c) A PM synchronous machine (PMSM).

The conventional wound rotor SM may be represented similarly by using three coupled inductors for the stator windings and a field winding for the rotor excitation, as depicted in Figure 1(b). Additionally, the rotor may have a specially designed damper cage, which may be similar in construction to the squirrel-cage rotor conductors of an induction motor and act to dampen the electromechanical oscillations of the rotor with respect to the stator field. This damper structure may also be represented by equivalent windings placed on the direct and quadrature rotor axes (i.e., the d- and q-axes) of the SM. Depending on the modeling accuracy and fidelity, one or several equivalent damper windings may be assumed in the CCPD model of the SM, as depicted in Figure 1(b). Generally, one or two damper windings are considered sufficient for typical transient studies, and more equivalent damper windings may be included to better represent the stator electrical transients in the higher frequency range.

PMSMs are also gaining wide acceptance in power applications such as motors and generators ranging from several hundred watts to tens of megawatts. The CCPD model for such machines may be obtained by considering the stator phases as coupled inductors with a PM rotor, as illustrated in Figure 1(c). Similar to the conventional SMs of Figure 1(b), the rotor structural/magnetic saliency (as well as the damping effect of the rotor iron structure) may be considered by assigning appropriate values to the inductances and possibly adding some damper windings to the equivalent circuit representation.

The total number of equations describing the electrical subsystem of the machine is determined by the total number of equivalent windings considered in the coupled-circuit model. The mechanical subsystem is typically characterized by Newton’s laws of motion relating to the rotor’s mass or inertia, the angular position and speed, and all torques acting on the machine. Single-mass or multi-mass systems may be assumed depending on the modeling needs. The common properties of the CCPD models of the rotating machines depicted in Figure 1 are that their inductances will be variable and dependent on the rotor position and magnetic saturation phenomena if it is modeled. Thus, if implemented directly as it is defined, these models will generally require the recalculation of all the relevant inductances at every time step of the simulation, which adds to the numerical cost when considering these models in various EMT simulators. For a long time, this computational challenge prevented the wide use of the CCPD models directly and led to the development of the two-axis models as an alternative. However, with the development of more powerful computational hardware and efficient numerical techniques, the CCPD models are now being considered again in some EMT simulators.

qd0 Models

The structure of the rotors of typical field-wound SMs determines the direct magnetic axis (the d-axis) and the orthogonal quadrature axis (the q-axis). This distinct geometrical rotor property (as well as the computational challenge with the variable inductances in the CCPD models) led the analysis of electric machinery from the early 20th century to express all machines’ electrical equations using transformed variables projected into the rotating qd system of coordinates. This is known as R. H. Park’s transformation and two-reaction theory, which transforms the original physical variables and equations of the CCPD models into the rotating qd reference frame, which is naturally fixed on the rotor for the SMs. Depending on the convention used, i.e., which axis leads the other, and the assumed direction of rotation, one may find qd- and dq-axes used equivalently in the literature. In the 1960/1970s, P.C. Krause generalized the reference frame theory for electrical machinery by relating all commonly used reference frames, showing their unique properties. This powerful tool allowed the expression of most ac machine models in the qd reference frame for analysis, control, and simulation purposes. The zero sequence is also included to uniquely transform balanced and unbalanced three-phase quantities, and thus, qd0 coordinates. The transformed qd0 equivalent circuits of the typical induction machines, field-wound SMs, and PMSMs are shown in Figure 2.

Figure 2. The qd0 equivalent circuit models of typical three-phase rotating ac machines. (a) An induction machine. (b) A field-wound SM with equivalent damper windings on the rotor. (c) A PMSM.

Figure 2(a) depicts the qd0 model of induction machines that assumes a shorted squirrel-cage type rotor. The zero-sequence circuit is also shown on the stator side to account for possible zero-sequence current if the machine stator windings are assumed to be Y-connected with an accessible neutral point. If the neutral point is not used, which is often the case, then the zero-sequence circuit can be removed from consideration. The doubly-fed induction machine model is obtained by opening up the rotor circuits and making their terminals available for external connections. The SM model is shown in Figure 2(b), where one field winding (on the d-axis) and three equivalent damper windings (one on the d-axis and two on the q-axis) are assumed on the rotor side. The stator sides of the induction and SM qd0 models look similar, except for the point that induction machines are typically round and result in equal magnetizing inductances in the q- and d-axes, but those inductances may be different for the SMs depending on their rotor construction (i.e., round or salient-pole rotor type).

The PMSM qd model may be obtained by simplifying the equivalent circuit to that in Figure 2(c), where the rotor is represented by a source to account for the PM rotor flux. Similar to the conventional SMs of Figure 2(b), the rotor structural/magnetic saliency (as well as the damping effect of the rotor iron structure) may be considered by appropriate values for the inductances (and possibly adding some damper windings to the equivalent circuit).

An essential property of the transformed equivalent circuits in Figure 2 is that since the q- and d-axes are perpendicular, the magnetic coupling is removed and replaced by the so-called coupled speed-voltage sources. This, in turn, results in constant (i.e., rotor-position-independent) inductances that describe the magnetically linear qd0 machine models. This property also greatly simplifies the computations needed for executing such models in computer simulations, which made the qd0 models by far the most widely used even today. As the three-phase ac side is typically symmetric and balanced, the overall number of equations may also be reduced from a typical three per stator side (as in the CCPD model) to two per stator side (one for q-axis and one for d-axis), assuming that the ac side does not have a zero sequence. These features make the qd0 models very compact and computationally efficient as standalone components.

However, when interfacing with external systems, the qd0 models require special consideration since their internal calculations are carried out in transformed qd variables, but all the external systems are typically represented in the actual physical variables and abc phase coordinates. The challenge introduced by this interface may also differ depending on the EMT simulation program considered.

EMT Simulators

EMT simulation programs are typically the preferred choice when studying electrical power systems in the time domain with high accuracy, e.g., detailed short circuit transient analysis, overvoltage transients, controller tuning and optimization, motor startup transients, generator synchronization, the coordination and interaction of converter-based resources, etc. The EMT simulation programs commonly used in the power industry may be broadly categorized based on their underlying solution methods into state-variable-based (SVB) and nodal-analysis- (or modified-nodal-analysis)-based EMT programs (EMTP type). Each category also includes offline and real-time industry-grade simulators. Examples of widely used SVB programs include MATLAB/Simulink, Simscape Electrical toolbox, and PLECS for offline simulations as well as Opal-RT and Typhoon HIL for real-time simulations. Examples of EMTP packages include EMTP-RV, ATP-EMTP, MicroTran, and the widely used PSCAD/EMTDC for offline simulations and RTDS for real-time simulations.

Some programs include modules based on both solution approaches and allow their interfacing, such as DIgSILENT. Both families of SVB and EMTP simulators have their advantages and limitations (due to solving the electrical circuits differently). These programs are widely used in industry and academia. Generally, the models of electrical machines are crucial components for most power systems studies.

SVB EMT Simulators

The SVB EMT simulators represent electrical system dynamics based on ordinary differential equations (ODEs). The program formulates the ODEs in state-space form or as a collection of state equations based on the state variables of the system, e.g., the currents of inductors and voltages of capacitors as the states of the energy storage circuit components. Electrical systems often contain linear and nonlinear subsystems (e.g., circuit elements, controllers, etc.). Models of such subsystems may be separated into linear and nonlinear blocks, as illustrated in Figure 3. Each subsystem may accept inputs from external sources and exchange the coupling variables with other subsystems.

Figure 3. System partitioning into linear and nonlinear models in SVB EMT simulators.

Once the network topology is determined, a state-space model may be formed for the entire linear part of the circuit (for which the standard A, B, C, and D matrices may be created). However, the nonlinear or varying parts are typically solved in their general form without forming fixed matrices to minimize additional calculations. The time-domain solution is then calculated numerically by integrating the state-space equations using either fixed- or variable-time-step ODE solvers embedded in the program. The SVB EMT simulators commonly offer a variety of ODE solvers and user-defined settings to optimize the simulation performance and accuracy. However, due to computational cost and scalability challenges, such programs are typically used for simulations of smaller-scale electrical power systems, e.g., for control designs of electric drives, interconnection studies of small-scale renewable generation with few electrical machines, etc.

Implementing models of electrical machines in SVB EMT simulators is based on numerically integrating the corresponding voltage equations (either in the qd0 form or in the CCPD form), assuming either winding currents or fluxes as the state variables. Due to their inductive nature, when formulating a proper state model that is solved by numerical integration (as opposed to differentiation), these models require the machine’s terminal voltages as the inputs and compute the windings’ currents as the outputs. Thus, the natural interface of the machine model is formed by a voltage-controlled current source. This voltage-controlled current source requirement also brings challenges when interfacing the classical qd0 machine models with the external circuits. In most EMT simulators, the components and systems outside of the electrical machines are modeled and implemented in the actual physical variables (i.e., in the abc phase coordinates for the conventional three-phase systems). Depending on the composition of the external circuit to which the machine is connected, the interface may be compatible or incompatible.

Figure 4(a) shows an example of the interface, where the qd0 machine model is interfaced with an external network that receives currents as input and computes voltages as output of its submodel. The terminal voltages ${v}_{abc}$ are transformed from the abc phase coordinates into the qd0 coordinates and fed to the qd0 model as ${v}_{qd0}$, together with mechanical torque ${T}_{m}$, speed ${\omega}_{r}$, and position ${\theta}_{r}$ from the mechanical subsystem. The machine model computes the phase currents ${i}_{qd0}$ that are then transferred to the abc physical coordinates and fed back to the external subsystem as ${i}_{abc}.$ This interface is compatible since the machine model requires phase terminal voltages as input and computes winding currents as output. Examples of networks that receive current and output voltage are circuits composed of capacitors, resistors, voltage sources, etc.

Figure 4. Interfacing of the qd0 model of ac electrical machines in SVB EMT simulators. (a) A compatible interface with the external capacitive circuit. (b) An incompatible interface with the external inductive circuit. (c) A compatible interface with the external inductive circuit using artificial snubber resistors.

However, in most practical applications, electrical machines are connected to other network elements that are also inductive, such as lines, cables, transformers, or power electronic converters with a topology that changes depending on the state of its switching devices. The typical numerical solution of such networks also requires voltage as the input, and it computes current as the output, just like the qd0 machine model. Examples of networks that require input voltages and output currents are inductive circuits, nonlinear models/elements, current sources, or qd0 models of other electrical machines. An example of such an interface is shown in Figure 4(b), where the qd0 model is assumed to be connected to an inductive external network. This interface with the qd0 machine model is incompatible because of the input–output variable mismatch. Specifically, the voltages ${v}_{abc}$, which are required as the input to the machine model and the external circuit, remain unknown.

Sometimes, this problem may be solved by invoking iterations on the subsystems involved, which increases the computational cost and may have a convergence problem. The iterative approach may be used only for offline EMT simulators because such iterations are not possible in real time. However, in most SVB EMT simulators, this problem is solved using artificial snubber circuits placed at the interfacing terminals, as illustrated in Figure 4(c). Therein, with the addition of the snubber resistors ${r}_{x}$, the output currents of the ac machine (denoted as ${i}_{abc1})$ and the external network (denoted as ${i}_{abc2})$ can flow into the snubber. This enables establishing the interfacing voltages ${v}_{abc}$ required by both subsystems to make the interface compatible. Such snubbers may be realized using simple resistors or capacitors as their equations permit calculating the terminal voltages based on the input currents.

EMTP Simulators

In nodal-analysis-based EMT programs (i.e., EMTP-type programs), the equations of individual circuit elements are discretized at the component level using the trapezoidal integration rule with a fixed simulation time step $\Delta{t}$. For example, the voltage-current relationships for basic RLC components in continuous and discrete-time representation are demonstrated in Figure 5.

Figure 5. A representation of basic RLC components for the EMTP nodal solution. (a) Voltage-current relationships in the continuous-time domain. (b) Voltage-current relationships in the discrete-time domain based on trapezoidal integration rule. (c) Equivalent circuits in the discrete-time domain.

After discretization, each component is represented as a conductance in parallel with a so-called history current source (which is equal to zero for resistors). It is noted that the history terms depend on the solution at the previous time step (i.e., ${t}{-}{\Delta}{t}$). Nonlinear and varying components can also be modeled using different techniques, e.g., piecewise-linear representations, current injections with a time-step delay, compensation-based methods, etc. All the resultant voltage-current relationships are organized in the nodal form as ${\bf{Gv}}_{n} = {\bf{i}}_{h}$. Here, ${\bf{G}}$ is the conductance matrix of the entire network, ${\bf{v}}_{n}$ is the vector of nodal voltages, and ${\bf{i}}_{h}$ is the vector containing the sum of the history and independent current sources injected at the nodes. In this formulation, each component is represented by the entries in the network conductance matrix ${\bf{G}}$ and in the voltage and current vectors ${\bf{v}}_{n}$ and ${\bf{i}}_{h}$.

Figure 6(a) illustrates the interface of a three-phase ac machine at the nodes a, b, and c to the network in EMTP. The discretized ac machine model is interfaced using the Norton equivalent circuit shown in Figure 6(b), with ${\bf{G}}_{abc}$ (the equivalent conductance matrix of the subject electrical machine) as the interfacing submatrix. Connecting the machine model to the external network at nodes a, b, and c is visualized in Figure 6(c) by writing the values of the machine submatrix ${\bf{G}}_{abc}$ to the appropriate entries in the overall ${\bf{G}}$ matrix as well as the ${\bf{v}}_{n}$ and ${\bf{i}}_{h}$ vectors.

Figure 6. The interfacing of ac machine models in EMTP programs. (a) A three-phase machine is connected to the external network at interfacing nodes a, b, and c. (b) The Norton equivalent interfacing circuit of the machine discretized model. (c) The nodal equation of the discretized network and addition of equivalent submatrix/subvectors of the ac machine model. (d) The modified augmented nodal-analysis formulation.

The nodal equation ${\bf{Gv}}_{n} = {\bf{i}}_{h}$ is solved at every time step ${\Delta}{t}$ for the unknown voltages in ${\bf{v}}_{n}$ using efficient numerical techniques (which avoid the direct matrix inversion), after which the currents of branches are calculated. This is known as H. Dommel’s algorithm, which enabled the efficient and scalable EMT simulation of power networks. It is important to note that for a network with time-invariant components, the ${\bf{G}}$ matrix will be constant, and the program needs to calculate it only once. Otherwise, the ${\bf{G}}$ matrix must be updated whenever the values of the circuit components change; thus, solving ${\bf{Gv}}_{n} = {\bf{i}}_{h}$ will require more computations.

Some offline EMT simulation programs (e.g., EMTP-RV) use the so-called modified augmented nodal analysis (MANA). Therein, the network equations are organized in the general form shown in Figure 6(d), where the vectors ${\bf{v}}_{n}$ and ${\bf{i}}_{h}$ as well as the ${\bf{G}}$ matrix correspond to their classical nodal-analysis counterparts, while there are several other matrices and vectors that are augmented to include circuit components that cannot be modeled as a conductance, e.g., ideal voltage and current sources, current-dependent circuit elements, and/or ideal switches and transformers. In MANA, the voltage/current relations of such elements are considered by the row, column, and diagonal block submatrices ${\bf{A}}_{r}$, ${\bf{A}}_{c}$, and ${\bf{A}}_{d}$ and by the vectors of unknowns ${\bf{x}}$ and knowns ${\bf{b}}$, respectively, as shown in Figure 6(d). For example, for ideal voltage (or current) sources, ${\bf{A}}_{r}$ and ${\bf{A}}_{c}$ would be incidence matrices having either 0 or ${\pm}{1}$ values possessing information about the connection of sources between nodes, ${\bf{b}}$ would include the values of the sources, and ${\bf{x}}$ is the currents (or voltages) of the sources.

For the ideal transformers, these submatrices would also include the turn ratios as the scaling factors. Ideal and nonideal switches can also be included using the new submatrices ${\bf{A}}_{r}$, ${\bf{A}}_{c}$, and ${\bf{A}}_{d}$ without changing the rank of the overall matrix. Although the MANA approach presents advantages compared to the classical nodal analysis, it may require iterations for the solution of nonlinear devices and electrical machines. Generally, it is preferred to solve the nodal equations noniteratively in real-time EMT simulators.

For implementation in EMT simulators, the models of ac machines are also discretized using the trapezoidal rule with a time step ${\Delta}{t}$. Moreover, since the classical qd0 models are derived in transformed variables and coordinates, obtaining their equivalent discretized Norton equivalent interfacing circuit with ${\bf{G}}_{abc}$ is not straightforward, and several approaches are possible, generally requiring some form of approximation. Different solutions have been proposed in the literature to overcome the incompatibility issues of the qd0 models in the EMTP solution. One approach to achieve the interfacing is using the stator voltages from the previous time step (i.e., at ${t}{-}{\Delta}{t}$) and modeling the machine as a Norton equivalent current source in parallel with artificial snubbers, as illustrated in Figure 7(a). This approach results in a fixed interfacing conductance ${\bf{G}}_{snub}$, and it is used in programs such as PSCAD. Another approach may be to average the unequal parameters of the discretized qd0 model to obtain a constant ${\bf{G}}_{abc}$ in abc coordinates and to use some form of prediction to compensate for the error. This approach is illustrated in Figure 7(b) and is used in MicroTran.

Figure 7. The interfacing of qd0 models of electrical machines in the EMTP solution. (a) Using snubbers and stator voltages from the previous time step. (b) Using the equivalent conductance matrix and some approximations.

Since the transformations between the qd0 and abc interfacing variables are nonlinear, it is also possible to iteratively solve the machine-network interfacing within the MANA formulation in offline EMT simulators. The so-called compensation-based approach may also be applied, where a Thevenin equivalent of the external network in qd0 coordinates is obtained, and the machine is interfaced and solved with the external grid in qd0 coordinates. This technique can be found in the universal-machine model in ATP. All these methods achieve the desired interfacing of the qd0 models in EMT simulators but possess different numerical properties, and their accuracy generally quickly degrades when the simulation time step is increased. Because of that, some researchers turned to discretizing the CCPD machine models and using them in EMT simulators directly. As it turned out, the discretized CCPD model naturally has a constant ${\bf{G}}_{abc}$ matrix for a magnetically linear induction machine. For SMs, however, this matrix is rotor-position-dependent.

Voltage-Behind-Reactance Models

To address the challenges of interfacing the traditional qd0 models of ac machines in EMT simulation programs, the voltage-behind-reactance (VBR) models have been developed more recently. The basic idea of these types of models is to separate the stator equivalent interfacing circuit (in physical variables and abc phase coordinates) from the rotor subsystem (which uses the transformed qd variables). S. D. Pekarek first pioneered this idea in the late 1990s. This approach allows interfacing the machine model using the VBR equivalent circuit, as shown in Figure 8. Depending on the considered ac machine and the VBR model type, the interfacing inductances ${\bf{L}}^{''}_{abc}$ may be variable, i.e., dependent on the rotor position and magnetic saturation, or made constant (with some approximations), resulting in constant-parameter VBR (CPVBR) models that are particularly advantageous for EMT simulation programs.

Figure 8. The interfacing of the VBR models of ac electrical machines in SVB EMT simulators.

Example of VBR Models in SVB EMT

Due to its VBR structure, in SVB EMT simulators, this model is interfaced using current-controlled voltage sources and resistor and inductor elements, as depicted in Figure 8. The interfacing VBR circuit is solved together with the external circuit. The machine’s stator phase currents ${\bf{i}}_{abc}$ (which are also the outputs of the external circuit in this system of Figure 8) are inputs to the rotor subsystem model of the ac machine. Other inputs include mechanical torque ${T}_{m}$, rotor speed ${\omega}_{r}$, and position ${\theta}_{r}$, all of which come from the mechanical subsystem. The VBR model computes the so-called subtransient back electromotive force (emf) voltages ${\bf{e}}^{''}_{qd0}$, which are then transformed to the abc coordinates as ${\bf{e}}^{''}_{abc}$ and used as inputs for the controlled voltage sources behind the reactance in the interfacing circuit, as depicted in Figure 8. This allows connecting the VBR models to arbitrary configurations of external circuits directly without artificial snubbers, as was needed in Figure 4(c).

For numerically efficient simulations, it is desirable to have constant (i.e., time-invariant) series impedances in the VBR models. Otherwise, the program must calculate the values of the interfacing branch parameters at every time step of the simulation, which can significantly increase the computational cost. For this purpose, a family of CPVBR models has been developed for induction and synchronous machines. Figure 9 shows a unified constant-parameter decoupled VBR RL-branch equivalent circuit for interfacing ac induction and synchronous machine models with external electrical networks. This circuit includes an optional zero-sequence branch, which is needed if the neutral point of the Y-connected stator winding is actually used in the studies; otherwise, this branch is simply removed. For round and symmetrical induction machines, this circuit is derived without approximations and is composed of constant and decoupled inductors ${L}_{D}$ and resistors ${r}_{D}.$ For SMs, some approximations may be required to arrive at constant-parameter and decoupled interfacing elements with explicit noniterative solutions, but the same interfacing circuit shown in Figure 9 is achieved.

Figure 9. A unified CPVBR RL-branch equivalent circuit for interfacing induction and synchronous machine models with external electrical network circuits.

Induction Machine qd0 Versus VBR Models (Example Case 1)

A simple system composed of an induction motor fed from a source with a small inductive impedance is considered in Figure 10 to give the reader an idea about the impact of the interfacing of different models. Two models are considered to simulate a startup acceleration transient of the energized motor. The traditional qd0 model is implemented based on the equivalent circuits presented in Figure 2(a) and is interfaced with the inductive source according to Figure 4(c) using large resistive snubbers. The VBR model is implemented using the interfacing circuit of Figure 8 and is directly connected to the inductive source with no snubbers. Choosing the values for the snubbers is always a compromise decision between the modeling accuracy and the numerical performance of the simulation. Generally, one wants to select the snubbers as large as possible to make the snubber current very small, and thus, to make the model more accurate. However, this increased accuracy comes at the cost of making the model numerically stiff and possibly leading to smaller/shorter time step sizes (more time steps) and longer computation times for offline simulations. Snubbers may also limit the largest permissible time-step size for the real-time simulators.

Figure 10. A single induction motor fed from a voltage source with an inductive impedance emulating line cable. Snubber resistors are used to interface the conventional qd0 model with an inductive source.

Figure 11(a) shows the rotor phase a current obtained with the qd0 and VBR models compared to the reference solution. It is more visible in the magnified view in Figure 11(b) in which, even with the appropriately selected snubber resistors, the qd0 model has a visibly higher error and also requires more simulation time steps, making this model computationally more expensive. At the same time, the directly-interfaced VBR model uses much larger time-step sizes (i.e., fewer steps) and produces the solution points that land precisely on the reference trajectory. So, the demonstrated benefits of the VBR model include achieving better accuracy and permitting large time steps.

Figure 11. (a) The rotor current during startup acceleration transient as predicted by the qd0 and VBR models. (b) Its magnified view showing the difference between the models. Ref.: reference solution.

Modeling Magnetic Saturation

Including magnetic saturation is often considered to improve the fidelity of ac machine models in commercial EMT simulation packages. The saturation characteristic can be easily obtained from no-load tests, and it captures the nonlinear relationship between the magnetizing current ${i}_{m}$ and the main magnetizing flux magnitude ${\lambda}_{m}$ (typically the fundamental component). This nonlinear characteristic is depicted in Figure 12(a), and in the EMT simulators, it is typically represented using several piecewise linear segments. For induction machines, the saturation is typically applied to the main flux path, and since they are round and symmetrical, the saturation is incorporated into the qd0 models by considering the projections of the saturated flux vector ${\lambda}_{m}$ onto the q- and d-axes (therefore taking into account the so-called qd-axes cross-saturation), as depicted in Figure 12(b).

Figure 12. Modeling saturation in three-phase ac machines. (a) A typical nonlinear magnetizing flux–current relationship. (b) Projections of saturated magnetizing flux onto q- and d-axes in a round rotor machine. (c) Using the saliency factor to convert a salient-rotor machine to an equivalent round machine for representing magnetizing flux saturation.

For SMs, depending on the rotor construction, the saturation is sometimes applied to the d-axis only, to the q- and d-axes independently, or to the main magnetizing flux ${\lambda}_{m}$. Experiments and FEA suggest that saturation of the main flux magnitude results in better accuracy. For the qd0 model of the salient-pole SMs, the single-saturation factor or the two-saturation factor approaches can be utilized depending on the availability of the q-axis and intermediate directions of the saturation characteristic. However, in practice, only the d-axis saturation characteristic can be easily obtained experimentally from the open circuit tests. Therefore, a single saliency factor approach is often considered sufficient for the general-purpose models in EMT simulators. The saliency factor is defined as the square root of the ratio of the q- and d-axes magnetizing inductances ${L}_{mq}$ and ${L}_{md}$, as ${SF} = {\sqrt{{L}_{mq} / {L}_{md}}}$, and it is assumed to be constant for unsaturated and saturated conditions. Using this approach, the saturation is then implemented in the qd0 models by rescaling the magnetizing flux and current projections of the q-axis by the factor ${SF}$, as illustrated in Figure 12(c), which effectively converts the original salient rotor into an equivalent round rotor.

The described saturation approach assumes that the magnetizing flux wave is sinusoidally distributed in the machine’s air-gap and that this entire sinusoidal wave is being reduced due to saturation. The magnitude ${\lambda}_{m}$ of this flux wave follows the saturation characteristic presented in Figure 12(a). However, it has been shown in the literature that the flux distribution can become considerably distorted as a result of different saturation levels in different parts/zones of the machine, which creates spatial harmonics in the distribution of the magnetizing flux. A common cause is a deeper saturation in stator/rotor teeth/slots along the magnetization flux path compared to the sides. As a result, the magnetic flux density under saturation will have a flat-top waveform containing odd harmonics, as depicted in Figure 13(a). Therefore, the air-gap flux spatial harmonics can also be considered for higher modeling fidelity. Figure 13(b) demonstrates the fundamental component ${\lambda}_{m1}$ and the third harmonic ${\lambda}_{m3}$ for a typical induction machine as a function of the magnetizing current ${\boldsymbol{i}}_{m}$. Thus, it is possible to consider the air-gap flux harmonics to increase the accuracy of representing saturation in ac machine models in EMT simulators.

Figure 13. Considering air-gap flux harmonics when modeling saturation in ac machines. (a) A flat-top magnetizing flux waveform and its fundamental and third harmonic components. (b) The saturation characteristic of the magnitude of the fundamental and third harmonics flux versus the magnetizing current.

Air-Gap Flux Harmonics in Induction Machines (Example Case 2)

An induction machine example is considered here to illustrate the representation of saturation. Recent developments in advanced machine drives have given rise to robust sensorless control strategies based on monitoring certain harmonics in the machine’s voltages and currents resulting from magnetic saturation to pinpoint the magnetic flux vector position. For simulations of such drive systems in EMT simulators, it may be important to consider saturation with air-gap flux harmonics. While a wide range of harmonics can be present, the third harmonic is typically considered due to its dominant effect. Magnetic saturation may be implemented in most models, including CCPD, qd0, and variable- and constant-parameter VBR models, as illustrated in Figure 14. In all the considered models, the third harmonic ${\lambda}_{m3}$ may also be considered as an option, on top of the saturation of the fundamental magnetizing flux ${\lambda}_{m1}$. Each model considered in Figure 14 has unique interfacing properties, impacting the computational performance and numerical accuracy.

Figure 14. The implementation of the saturable (a) CCPD, (b) qd0, (c) VBR, and (d) CPVBR models of induction machines and their interface with an external network.

In the CCPD model, all the stator and rotor windings are represented as rotor-position-dependent and coupled inductive branches, as illustrated in Figure 14(a), that are solved together with the external network. The saturation is considered appropriately by defining all inductances in terms of the magnetizing current ${i}_{m}$ and the rotor position ${\theta}_{r}$. The computational cost of this model will be high due to the variable nature of all inductances.

In the qd0 model, the saturation may be represented using a flux correction method that follows the saturation characteristic of Figure 12(a) and (b) for the fundamental magnetizing flux component ${\lambda}_{m1}$. The implementation and interfacing of this model are illustrated in Figure 14(b). If the air-gap flux harmonic is considered, additional subcircuits, including the zero-sequence circuit, have to be added on top of the traditional equivalent circuit depicted in Figure 2(a) to include the effect of ${\lambda}_{m3}$. Cross-coupling between the fundamental component ${\lambda}_{m1}$ and its third harmonic ${\lambda}_{m3}$ is also considered. Snubber resistors are typically required for interfacing with this model, as depicted in Figure 14(b).

The VBR model may also be implemented with the saturation of the main magnetizing flux ${\lambda}_{m1}$ and an option of including the third harmonic ${\lambda}_{m3}$ of the air-gap flux. A direct extension from the saturable qd0 model results in a VBR model with saturation-level-dependent interfacing inductances. This property increases the computational cost, similar to the saturable CCPD model, although with only three variable coupled inductors instead of six. With some approximations, a CPVBR model may also be derived. These models and their interface are depicted in Figure 14(c) and (d).

An experimental study with a 0.25-hp squirrel-cage induction motor is presented here to demonstrate the performance of the saturable models. The experimental setup is illustrated in Figure 15(a). The measured and simulated steady-state voltages and currents for different levels of applied voltage are shown in Figure 15(b). It can be observed in Figure 15(b) that, as the voltage magnitude increases from 1 per unit (p.u.) to 1.2 p.u. and then to 1.4 p.u., the waveforms become more distorted, which is captured by the considered models. Here, all the considered models use small time steps and produce results consistent with the experiments.

Figure 15. Validation of induction machine models with saturation. (a) The experimental setup. (b) Steady-state voltage and current waveforms under different levels of saturation due to the increased voltage level.

A startup transient is simulated next, with the results presented in Figure 16. As can be seen in Figures 15 and 16, the currents are visibly distorted due to harmonics, which all models capture well along with the mechanical variables. However, Figure 16(b) shows the difference among the solutions produced by the considered models, particularly that the qd0 model requires very small time steps (i.e., a large number of steps) while being slightly off in terms of accuracy due to its snubber interface. It was also found that the CPVBR model may be about three times more efficient per step compared to the variable VBR model while achieving similar accuracy.

Figure 16. The induction machine startup transient predicted by the qd0, VBR, and CPVBR models under saturation. (a) The stator phase current and (d) its magnified view. (b) The electromagnetic torque and (e) its magnified view. (c) The mechanical speed.

Example of SM in EMTP

The role of the interfacing of ac machine models in EMTP may be illustrated by considering a single-machine infinite-bus system with a large (835-MVA) steam turbine generator subjected to a three-phase fault at the terminals. The saturation is not considered here to focus on the interface only. The considered models include variable-parameter CCPD and VBR models and two versions of the qd0 model (qd0_1 is iterative using EMTP-RV, and the second model is qd0_2 with prediction using MicroTran). It is important to mention that when using time steps on the order of 50 µs to 100 µs, all models produce reasonable results visually close to each other, with up to 2% numerical error with respect to the reference solution. Therefore, the simulation step size is increased to 1 ms to exemplify the difference among the considered models.

The result of this transient in the stator phase a current computed by the considered models is shown in Figure 17. As seen in this figure, the results of different models are visibly different. The CCPD and VBR models achieve a direct simultaneous solution, and their results are closest to the reference even at such a large time step. However, the two qd0 models produce solution points that significantly deviate. It was also confirmed that the qd0 model interfaced using snubbers and delayed stator voltages is not numerically stable with this large time step and is therefore not included here. It is concluded that the variable-parameter CCPD and VBR models possess better accuracy, although their computational cost per step may be high due to the variable equivalent conductance matrix.

Figure 17. The three-phase fault transient study of an 835-MVA steam turbine generator computed by several models. (a) The stator phase a current. (b) A magnified view showing the differences among the models.

The advantageous interfacing properties of the CCPD and VBR models in EMTP brought more attention to these modeling approaches. However, when these models are derived for synchronous machines, their interfacing conductance matrix may become rotor-position-dependent. Moreover, considering the nonlinear saturation characteristic as a piecewise-linear function additionally makes the interfacing conductance matrix saturation-segment/level-dependent.

This property inevitably adds computational cost to implementing the CCPD and VBR models and hinders their application, particularly in real-time simulators where it is highly desirable to have a constant overall network conductance matrix for efficient and noniterative solutions. Fortunately, it is also possible to derive a constant interfacing conductance matrix using careful approximations. The procedure for developing constant-parameter CCPD and VBR machine models for EMTP is summarized in steps illustrated in Figure 18. Starting from the discretized Thevenin equivalent circuit, the first step is to separate the interfacing resistance matrix of the model into two parts—one constant part and one varying (rotor-position- and saturation-dependent) part. Next, the varying part is moved into the history voltage sources, which also requires predicting/approximating the stator currents. Finally, the new Thevenin equivalent circuit is converted to the Norton equivalent with a constant conductance matrix and appropriately modified history terms.

Figure 18. Steps for obtaining the CPVBR and CCPD models for EMTP. (I) Decomposing the original resistance matrix into constant and varying parts. (II) Moving the varying part into the history terms using predicted stator currents. (III) Obtaining the constant-parameter Norton equivalent circuit for interfacing in EMTP simulators.

The merits of possessing a constant-parameter interfacing circuit are demonstrated on an IEEE 39-bus system consisting of 18 synchronous generators, as shown in Figure 19. The system has been implemented in PSCAD, where all generators have been implemented using CCPD, VBR, and CPVBR models with the saturation of the main magnetizing flux, as presented in Figure 12(c). A transient study of a three-phase-to-ground fault at bus 25 is considered. The 2-norm error of the fault current contributed by the generators G17 and G18 connected to bus 37 is considered and reported in Table 1. The simulation times for conducting a 0.5-s study in PSCAD are provided in Table 2. Table 1 shows that all the considered models produce very similar accuracy for time steps of 100 µs and 500 µs, which confirms the approximation steps in Figure 18. However, Table 2 shows that the CPVBR model significantly outperforms the variable-parameter CCPD and VBR models by a factor of 4.7 in terms of computational efficiency. This computational gain is demonstrated due to the constant overall conductance matrix of the multimachine system, testifying to the advantage of having constant-parameter machine models.

Figure 19. The single-line diagram of the modified IEEE 39-bus system comprising 18 synchronous generators.

Table 1. Average 2-norm error of fault current contributed by synchronous generators connected to bus 37 for different time-step sizes expressed in percentage (%).

Table 2. CPU time required in PSCAD for simulation of 0.5-s short circuit transient study.

Concluding Remarks and Future Work

For EMT simulators, it is crucial to have computationally efficient and accurate models of ac electrical machines. Such models should be capable of representing essential electromagnetic and mechanical transient phenomena for studies of individual machines with their controllers and power electronic converters and for system-level studies with many machines participating as generators and motors. The typical machine models include a representation of magnetic saturation, which may be important for a wide range of studies. In addition to the traditional two-axis qd0 models, the CCPD and VBR models may offer interfacing advantages, translating into numerical robustness and good accuracy even with large time-step sizes. When such models are implemented in EMT simulators, it is also important that they are realized in the most efficient way, requiring the smallest possible number of calculations. Then, such ac machine models may enable extending the application of modern offline and real-time EMT simulators to larger power grids while optimizing the use of computational resources.

Electrical machines with more than three (e.g., five, six, nine, twelve, etc.) electrical phases are becoming attractive for some applications, offering increased reliability, fault-tolerant operation, and reduced losses. Models of such machines in EMT simulations are also getting attention. Further improvement of the models’ fidelity, numerical accuracy, and computational efficiency represents ongoing research with significant potential impact and benefit to many simulation software vendors, utilities, and consulting companies.

For Further Reading

H. W. Dommel, EMTP Theory Book. Vancouver, BC, Canada: MicroTran Power System Analysis, May 1992.

P. C. Krause, O. Wasynczuk, S. D. Sudhoff, and S. Pekarek, Analysis of Electric Machinery and Drive Systems, 3rd ed. Piscataway, NJ, USA: IEEE Press, 2013.

L. Wang et al., “Methods of interfacing rotating machine models in transient simulation programs,” IEEE Trans. Power Del., vol. 25, no. 2, pp. 891–903, Apr. 2010, doi: 10.1109/TPWRD.2009.2039809.

E. Levi, “Saturation modelling in d-q axis models of salient pole synchronous machines,” IEEE Trans. Energy Convers., vol. 14, no. 1, pp. 44–50, Mar. 1999, doi: 10.1109/60.749146.

S. D. Pekarek, O. Wasynczuk, and H. J. Hegner, “An efficient and accurate model for the simulation and analysis of synchronous machine/converter systems,” IEEE Trans. Energy Convers., vol. 13, no. 1, pp. 42–48, Mar. 1998, doi: 10.1109/60.658202.

Y. Xia et al., “An efficient phase domain synchronous machine model with constant equivalent admittance matrix,” IEEE Trans. Power Del., vol. 34, no. 3, pp. 929–940, Jun. 2019, doi: 10.1109/TPWRD.2019.2897612.

N. Amiri, S. Ebrahimi, and J. Jatskevich, “Saturable voltage-behind-reactance models of induction machines including air-gap flux harmonics,” IEEE Trans. Energy Convers., vol. 38, no. 3, pp. 2022–2033, Sep. 2023, doi: 10.1109/TEC.2023.3259385.

Biographies

Erfan Mostajeran (emostajeran@epeconsulting.com) is with Electric Power Engineers LLC, Austin, TX 78738 USA.

Navid Amiri (navida@ece.ubc.ca) is with the Electrical and Computer Engineering Department, University of British Columbia, Vancouver V6T 1Z4, Canada.

Seyyedmilad Ebrahimi (ebrahimi@ece.ubc.ca) is with the Electrical and Computer Engineering Department, University of British Columbia, Vancouver V6T 1Z4, Canada.

Juri Jatskevich (jurij@ece.ubc.ca) is with the Electrical and Computer Engineering Department, University of British Columbia, Vancouver V6T 1Z4, Canada.

Digital Object Identifier 10.1109/MELE.2023.3320508

Rotating Frame, Average Value Converter Modeling: Basic principles in building analytical models

Dragan Jovcic

©SHUTTERSTOCK.COM/PPR109103

Electrical power converters have become indispensable units in modern power systems, with an increasing number of installations, application fields, and ratings. They are used for power conditioning; for controlling network, regulating generation, or load power; for integrating storage systems; and importantly for integrating most renewable energy plants. As a typical example, a modern wind generator will have two ac–dc converters back-to-back connected, facilitating optimal energy extraction from variable speed wind turbines and ensuring that the power plant meets network connection requirements.

In most cases converters are controllable, incorporating a complex control system that facilitates regulation of some variables and ensures system stability under disturbances. The inherent high controllability and increasing penetration of converters bring challenges to system designers and operators because of the impact that converters have on system stability and more generally on security of power exchange.

By its nature, converters will facilitate electrical power exchange between two different systems, and a typical example will be power conditioning between dc and three-phase ac circuits. These two electrical systems have fundamentally different electrical variables (oscillating versus nonoscillating), causing complications in joint system design and operation.

Converters are built using semiconductor power switches, typically from transistor or thyristor families. In power electronic applications, these semiconductors are exclusively used as switches, i.e., they are in either in an ON (conducting) or OFF (nonconducting) state. The switching between states occurs at high frequencies (many times per second) and can be controlled. However, the current flowthrough of a semiconductor while in the ON state cannot be controlled. The current control through a converter (at the system level) is only indirectly controlled by modulating turn ON and OFF instants. The changes in switch states are abrupt and lead to highly nonlinear waveshapes of converter variables.

The above aspects bring a range of challenges in converter modeling. They need substantially different models compared with traditional electromechanical power system elements, like transformers or rotating generators.

Power systems engineers require converter models for numerous system studies, including component dimensioning, power flow, stability, protection, or analysis of contingencies during all stages of system design and operation.

Detailed Versus Analytical Modeling of Converters

Thanks to the rapid development of computing power and software tools in the past 30 to 40 years, engineers are nowadays able to perform quite detailed and accurate converter modeling and simulation. Converter simulation in the time domain is often employed using detailed modeling on electromagnetic transient (EMT) platforms with several commercial suppliers established worldwide. These platforms perform calculations of all system variables in each step and simulation progresses step-by-step in the time domain using small steps of a few microseconds. Mathematically, the simulation methodology is simple, and it achieves good convergence and robustness. Their major advantage is in flexibility, since users can enter any circuit, but changes in topology are also allowed while the simulation is running (like switching of semiconductors). The models are accurate, flexible, with good graphical interface, and they are widely adopted in the professional community. The EMT modeling and most converter integration studies can be done by practicing engineers without advanced skills or research.

However, EMT simulation has some limitations and challenges, including:

Simulation time may become prohibitively long for complex systems. As an example, a high-voltage converter used with high-voltage dc (HVdc) transmission, may have over 300 submodules in each of six semiconductor valves, and each submodule has two basic transistors, which are individually controlled. It has been reported that simulation may take over an hour for 1 s of real time, when a system with multiple HVdc converters is studied.
Only time-domain analysis is possible, which means that predominantly trial-and-error design and analysis approaches are used. Repeated simulation runs are always required. The performance indicators are largely based on minute changes in time-domain waveshapes of variables. With nonlinear detailed digital simulation, it is very difficult to perform qualitative system analysis or parametric studies.

Scientists and engineers have always strived to develop models that capture the essence, key properties, and dependencies of the system under examination, while being simple and convenient for the studies. Some ingenious scientific principles and engineering modeling methods have been recorded in history, and a few that have made notable impact in electrical power engineering will be mentioned.

Park’s transformation, introduced in late 1920s by Robert H. Park, simplified analyses of three-phase ac circuits, and in particular helped understanding ac rotating machines. This method introduces a rotating coordinate frame that eliminates time-varying parameters and translates oscillating variables into dc variables, which are more convenient for analysis. In turn, Park’s transformation builds on the more general earlier coordinate frame transformation theories, like the theory of moving observer used by Einstein in his theory of relativity.

Another very important principle is Fourier series expansion, introduced by Joseph Fourier in the 18th century, that enables representation of any signal as a sum of simple oscillating signals at a range of frequencies. This method is widely used in engineering for studying periodic signals, and is especially important with ac power circuits, where all variables are oscillating at the operating frequency (typically 50 or 60 Hz). Thanks to the Fourier series, each ac variable can be expressed as a simple sine signal at fundamental frequency and numerous other sine signals at harmonic frequencies.

The third essential method is Taylor series expansion, introduced by Brook Taylor in the 18th century. It facilitates representing a nonlinear function/signal as a sum of simple linear terms dependent only on the input variable and its derivatives. Considering the highly nonlinear nature of converter behavior, linearization plays a key role in converter studies.

This article describes analytical modeling principles for power converters. Development of a converter analytical model is generally much more challenging than EMT modeling, since it requires in-depth converter knowledge and the use of advanced mathematical concepts, like coordinate frame transformation and Fourier or Taylor series expansion. It is not routinely used in industry, although most operators and manufacturers will use some analytical models at least at some stages of substantial studies. The effort in developing an analytical model may be justified when new topology is evaluated, repeated (control) issues are observed in practice, or a highly complex system is brought for feasibility evaluation.

Analytical modeling inevitably involves simplifications, which require understanding the tradeoff between accuracy and convenience, and therefore necessitate careful verification of the model accuracy. It may not be easy to understand modeling assumptions, and some practicing engineers will have misgivings about the validity of analytical models and may be concerned about model robustness and flexibility. However, many good analytical models have been reported in the professional community for a range of converter topologies. Numerous research teams have provided invaluable insights into converter design and integration based on analytical models. In particular, analytical models facilitate parametric studies, which lead to qualitative conclusions without repetitive time-domain simulations. One practical example is a system stability study, which can be performed either in time domain by observing the evolution of system variables following a perturbation, or using parametric domain methods from control system theory, like eigenvalue study. Arguably, analytical models bring the most value with multidimensional design tasks, like tuning multivariable controllers or broad robustness studies.

This article will describe essential steps in building analytical models for few typical converter topologies, including averaging, coordinate frame transformation, linearization, and validation.

An illustration of analytical model use for stability analysis on a practical converter test system will be provided in the last section. This case study will include comparison with stability analysis based on detailed EMT modeling for the same test system to highlight differences and limitations of each method.

Average Value Converter Modeling

Basic Principles

The first principle with converter analytical modeling is the averaging step. Because of the use of semiconductors, converter variables (voltages and currents) will have “chopped” waveforms, which are described as nonlinear and discrete in mathematical terms. The modeling problem is made worse by the fact that typically both dc and ac chopped signals are involved. Another challenge is that the converter control signal is active (affects the system) only at the switching instants and basically has no impact between turn ON and OFF. However, switchings occur at high frequencies (hundreds of hertz or a few kilohertz) compared to dominant system dynamics, which are typically restricted to a few tens of hertz. This separation of frequency ranges justifies an averaging approach.

The averaging process derives converter variables as continuous and simple signals, which nevertheless should capture well the key converter behaviors. It follows different methodology for ac and dc variables.

Also, averaging differs between converter topologies, and here it will be illustrated through example converters from two common topologies: voltage source converters (VSC) and current source converters (CSC).

VSC Average Modeling

Figure 1 shows a typical ac–dc VSC and its key variables. This converter enables power exchange between the dc system and three-phase ac system. A typical application would be converter-enabling interconnection of a battery with an ac grid, or one of two converters in a modern wind generator. It consists of six semiconductor valves (each may contain many semiconductors) S₁–S₆, an inductor on the ac side L, capacitor on the dc side C, and a controller that sends control signals to switches m₁–m₆. There are multiple (21 in this example) control pulses per each 50-Hz cycle, and timing of the pulse is achieved according to the desired ac voltage waveform. With VSC topology, dc voltage is constant between switchings, which is facilitated by large capacitor C on the dc side.

Figure 1. VSC with detailed and average variables. (a) VSC topology. (b) Detailed and average variables.

The actual converter ac variables V_acA, V_acB, V_acC, I_acA and dc variables V_dc and I_dc are shown in blue in Figure 1. The corresponding average variables are shown in red in Figure 1, which are obtained using Fourier series expansion on the ac side and simple averaging on the dc side. The obtained average variables are much more convenient for further analysis and design, and yet they capture essential converter dynamics of interest. A designer is primarily interested in converter variables at fundamental frequency (50 Hz) on the ac side and average variables on the dc side.

The average waveforms on the ac side are obtained using Fourier series expansion and neglecting harmonics. The waveforms on the dc side are commonly obtained by equating powers on the ac and dc sides of the converter, assuming that there is no loss or storage of energy in the converter. The dc voltage can be assumed as a firm, slow-varying variable with all converters from the VSC family.

Table 1 shows the key parameters for the test system converters considered in this article.

Table 1. Key parameters of the test systems.

CSC Modeling

Figure 2 shows the test CSC and its key variables. The essential differences compared with VSC are:

The dc current is constant between switching, achieved using inductor L on the dc side.
Thyristors are used as switches because of reverse blocking capability.
Switching occurs only once per half-cycle, for each switch.

Figure 2. Line commutated converter with detailed and average variables. (a) Line commutated converter topology. (b) Detailed and average variables.

The actual converter ac variables V_acA, V_acB, V_acC, I_acA and dc variables V_dc and I_dc are shown in blue in Figure 2. The corresponding average variables are shown in red in Figure 2, which are obtained using Fourier series expansion on the ac side and equating powers for dc side variables.

The dc current can be assumed as a firm, slow-varying variable with all converters from the CSC family.

Converter Modeling in Rotating Coordinate Frame

Oscillating variables are difficult for analysis, as it is awkward to observe and measure changes in such variables. Also, controllers need feedback measurements of system variables, but they have difficulties in processing oscillating signals. The theory of Park’s transformation gives elegant proof that we can use few dc variables to analyze an oscillating signal if our new coordinate frame rotates at the same frequency of the oscillating signal. In practical terms, this means that the observer is positioned to oscillate in synchronism with the system variables at the operating frequency. The key assumptions for the use of Park’s transformation are:

the operating frequency is constant
the three-phase system is symmetrical and balanced (circuits on each of three phases are identical).

These conditions are usually satisfied with ac transmission/distribution systems, and therefore this transformation is widely applied.

An oscillating signal is characterized by three variables—frequency, magnitude, and phase angle—and if frequency is eliminated (it is not changing), then only two variables are adequate. They can be represented in one of the following two ways: using a polar coordinate frame (magnitude and angle) or rectangular coordinate frame (d and q components).

Figure 3 shows the VSC converter phase A ac voltage in both a static ABC frame and in rotating dq coordinate frames. To illustrate correlation between these two coordinate frames, two step changes are applied:

at .6 s there is a change in d component of voltage V_acAd, from V_acAd = 287 kV to V_acAd = 224 kV
at .65 s there is a change in q component of voltage V_acAq, from V_acAq = −128 kV to V_acAd = 128 kV.

Figure 3. VSC converter phase A voltage in static (a) and rotating (b) coordinate frames.

In the rotating dq frame, the same signal is illustrated using polar (V_acAm, FiV_acA) and orthogonal components (V_acAd, V_acAq) for completeness. It is seen that in the dq frame the signal change is readily observed for both excitations. However, the permanently oscillating signal nature clutters the picture and obscures important information in the static ABC frame. The converter controller establishes the position of the rotating coordinate frame by measuring some reference voltage (usually grid connection point), as shown by the signal V_acAref in the ABC frame.

The detailed EMT model and the physical system will have all ac variables as oscillating signals. However, Park’s transformation is commonly implemented in practice for online measurement of various ac variables. These dq frame signals are used for monitoring and also as feedback signals in controllers.

The analytical model has no oscillating signals. It represents all internal dynamics and dependencies using dq ac signals and dc signals only. If required, the analytical model can also produce ABC frame oscillating signals by applying an inverse Park’s transformation on the obtained dq frame variables.

Model Linearization and Format Conversion

Averaged models for most converter topologies will have some nonlinear elements: as an example, multiplication of two signals. Such a nonlinear model can be used for time-domain simulation, for power flow study, component dimensioning, and with some nonlinear controllers.

However, most stability analysis and controller design theories require linear models, and typically they should be represented in state-space format. A linearized model is obtained by applying Taylor series expansion, or similar methods like describing function, to each of the identified nonlinear elements. A very important assumption with linearization is that variables do not deviate much around the steady-state operating point. This is an important factor that has an impact on accuracy and requires careful consideration; the linearized model is also commonly called a small-signal model. The same model may not be valid for large disturbances, like faults.

The model is commonly converted to a specific format, like state-space, depending on the intended use. State-space is the widely adopted model representation in the control engineering field.

Verification of Accuracy of the Average Model

The average value model always has a cloud of uncertainty related to accuracy hanging over it. The problems with accuracy may arise because of the overly liberal use of simplifications, contravention of some assumptions with transformations and methods, linearization assumptions, or simply because of errors caused by manual modeling associated with system complexity.

Model accuracy should be examined by comparing responses against a detailed model of known accuracy and for the same excitation. Ideally, a broad range of signals should be observed (dc side, ac side, control signals) to confirm accuracy of all subunits, and testing should be repeated for different excitations.

Figure 4 shows comparisons of dc variables (dc voltage V_dc) and an ac variable (d component of ac current I_acd) between a detailed and average model, for a 5% step increase in power transfer order. We can conclude that the average model has good accuracy based on two observations:

Steady-state values show good matching, as observed before and after excitation.
Essential dynamics are in good agreement, as demonstrated in the first 60 ms after perturbation. System dynamics could be further characterized in more detail using control theory indicators for step response, like overshoot and settling time.

Figure 4. Accuracy verification of average model.

Figure 4 strengthens arguments on the benefits of an average model considering clarity of communicating dynamic phenomena. It is much easier to evaluate dynamic response (like, for example, percentage overshoot) using an average model. The detailed model variables contain lots of noise caused by converter switchings and feedback control amplification.

Stability Analysis

Relevance and Challenges

In electrical power systems there are multiple forms of stability, but preeminent is small signal stability. It describes system response to small perturbations in steady-state operation. In practice, perturbations occur continuously: as an example, small changes in customer circuits, like loading demand or circuit topology (connection/disconnection). A system is stable if all variables remain bounded following a disturbance, while it is more desired that a system is asymptotically stable where all variables are bounded but also converge to a new steady state.

Robust stability is of particular importance but may be challenging in a general case. In practice, a system should be stable when operating under different parameters: as an example, different network characteristics or controller parameters. An optimally tuned controller may turn out to bring disappointment to a system designer when instability is observed under a different but plausible operating scenario.

A system designer would normally first provide a list of all parameters that are expected to vary, and the range of variation for each. Guaranteeing robust stability is difficult, as the number of possible scenarios could be large and dependencies between parameters could be highly nonlinear and unpredictable.

One of the main reasons for using analytical models is multivariable analysis of robustness of small-signal stability. A study of robustness can be tedious with detailed time-domain simulation. Analytical models, however, facilitate parametric stability study, which can be very effective for robustness analysis, as will be elucidated below.

Illustration of Stability Analysis Using a Detailed Model

Figure 5 shows a stability analysis case with the detailed EMT model of an HVdc system consisting of two VSC converters, as given in Table 1. The system operates normally at rated power in a stable manner until 1 s. At 1 s controller gain is increased 10 times and it can be observed that the system becomes unstable. Both the dc variable (dc voltage V_dc) and ac variable (d component of ac current I_acd) show oscillatory instability. In practice, such a large variation in variables would be detected by converter self-protection and the converter would be tripped, resulting in loss in capacity. This study enables us to estimate frequency of instability at around 78 Hz, since 7.8 full cycles can be counted in any .1-s interval.

Figure 5. Simulation of instability on detailed model. Controller gain increased 10 times at 1 s.

If a designer is interested in stability at another gain value (say K = 4) they would need to adjust the gain and repeat EMT simulation through the time interval from the beginning until the system settles in the steady state. Unfortunately, each new run demands model preparation by the software platform, which involves building and conditioning model matrices. In the case that new gain value gives stable response, the difference in the time-domain response may not always become readily perceptible.

System Stability Analysis in Paramedic Domain Using an Average Model

Figure 6 shows stability analysis for the same system as in the section above, but using an analytical model. The location of system eigenvalues is shown in complex plane with real and imaginary axes. Eigenvalues can be calculated routinely using an analytical state-space model and they provide indication of system stability. To achieve stable operation, all system eigenvalues should be located in the left-half-plane (negative real component). Eigenvalues closest to the imaginary axis (closest to instability boundary) have a dominant impact on the system response. The imaginary component of eigenvalues indicates frequency of the associated oscillating mode in radians per second.

Figure 6. Instability analysis on average linear state-space model. (a) Original system eigenvalues. (b) Eigenvalue location as controller gain increases. Im: imaginary; Re: real.

As seen in Figure 6(a), the original system has all eigenvalues in the stable region, which is consistent with the conclusions of stable operation obtained with the detailed EMT model. There are in total 56 eigenvalues for this system, but only a subset of eigenvalues closest to the imaginary axis is shown. In addition to stability conclusion, a designer engineer can exploit rich control theory to analyze system performance and robustness using eigenvalues. To this end, eigenvalue sensitivity and participation factor analysis would normally be performed.

In Figure 6(b), an illustration of parametric stability analysis is given. The controller multiplying gain K is varied in the range 1 < K < 10 and the position of eigenvalues is recorded. The original position of eigenvalues (K = 1) is indicated with blue circles, the final position (K = 10) with red triangles, and the intermediate positions with black crosses [Figure 6(b)]. It is seen that by increasing gain, the position of most eigenvalues deteriorates (they are moving toward the imaginary axis). Surprisingly however, it is not the dominant eigenvalues that are mostly sensitive to the gain variation, as another pair of eigenvalues first crosses the stability boundary at approximately K = 9. We know that this conclusion is correct, as both the gain value and frequency of instability are consistent with the simulations on the detailed EMT model shown in Figure 5. On the imaginary axis we can read the frequency of the unstable eigenvalues as around 490 rad/s [or 490/(2${\pi}$) = 78 Hz], which corresponds well with the frequency observed using the detailed EMT model.

The parametric study in Figure 6(b) can usually be performed fast, even with complex systems containing numerous eigenvalues and many gain values. This study also provides a full stability picture, showing all eigenvalues, with stability margin and frequency for each eigenvalue. In more advanced design tasks, multiple controller gains can be varied simultaneously to support multidimensional decision making.

Conclusion

Analytical converter modeling complements detailed digital simulation and may expedite system studies. It usually involves multiple development steps, including averaging of all variables, conversion to rotating coordinate frame, linearization, model format conditioning, and verification.

Analytical models facilitate parametric studies, which are not possible with detailed digital models based on time-domain simulation. A good example of model use is the eigenvalue analysis, which provides expedient and comprehensive stability examination.

Work is ongoing on automating the process of developing analytical converter models, and on connecting digital simulation with analytical models.

For Further Reading

D. Jovcic, High Voltage Direct Current Transmission: Converters, Systems and DC Grids, 2nd ed. Hoboken, NJ, USA: Wiley, 2019.

I. Papicˇ, “Mathematical analysis of FACTS devices based on a voltage source converter: Part 1: Mathematical models,” Electric Power Syst. Res., vol. 56, no. 2, pp. 139–148, Nov. 2000, doi: 10.1016/S0378-7796(00)00107-3.

D. Jovcic, N. Pahalawaththa, and M. Zavahir, “Novel current controller design for elimination of dominant oscillatory mode on an HVDC line,” IEEE Trans. Power Del., vol. 14, no. 2, pp. 543–548, Apr. 1999, doi: 10.1109/61.754101.

C. M. Osauskas, D. J. Hume, and A. R. Wood, “Small signal frequency domain model of an HVDC converter,” Inst. Elect. Eng. Proc. Gener. Transmiss. Distrib., vol. 148, no. 6, pp. 573–578, Nov. 2001, doi: 10.1049/ip-gtd:20010533.

Biography

Dragan Jovcic (d.jovcic@abdn.ac.uk) is professor and director of the Aberdeen HVDC (High Voltage Direct Current) Research Centre, University of Aberdeen, AB24 3FX Aberdeen, U.K.

Digital Object Identifier 10.1109/MELE.2023.3320486

Zero-Emission Marine Vessels: Multidomain Modeling and Real-Time Hardware-in-the-Loop Emulation on Adaptive Compute Acceleration Platform: Zero-emission marine vessels: modeling and real-time emulation

Chengzhang Lyu, Venkata Dinavahi

©SHUTTERSTOCK.COM/AUDIO UND WERBUNG

Proton exchange membrane fuel cells (PEMFCs) are becoming increasingly common in modern marine electric vessels, helping achieve the sustainable development objective of lowered emissions and zero-emission marine transportation. This article introduces a hierarchical hardware-in-the-loop (HIL) emulation scheme for the zero-emission marine vessel at the system level and device level. A comprehensive computational PEMFC model is presented in electrical, thermal, and fluid domains. Electromagnetic transient program models are utilized for batteries, power converters, and electric thrusters to describe their behavior and dynamic response and analyze the influence of system components on their performance. The real-time hardware emulation on the Xilinx Versal adaptive compute acceleration platform (ACAP) provides the ability to simulate the complete system at microsecond-level time intervals, which is essential for validating the marine vessel’s dynamic behavior under various operating conditions. The results demonstrate that the multidomain PEMFC model effectively captures the complex electrical, fluid, and thermal behavior and the interaction with marine vessels. Additionally, the proposed hierarchical HIL emulation scheme is proven to be a valuable tool for the design and testing of zero-emission marine vessels, which enables comprehensive assessment and verification of vessel performance.

Introduction

In recent years, there has been an increase in demand for hydrogen as a fuel for the future due to its immense environmental benefits, particularly the reduction of carbon emissions. Hydrogen fuel is considered a zero-emission energy source, which means there is no net release of carbon dioxide or other greenhouse gas into the atmosphere during operations. The emission profile of hydrogen makes it environmentally friendly and contributes to reducing greenhouse gas emissions and addressing climate change concerns. It seems that it is only a matter of time before hydrogen permeates additional facets of society, such as marine transportation (Skjong et al. 2017).

In 2018, the initial greenhouse gas strategy was adopted by the Marine Environment Protection Committee of the International Marine Organization (IMO). Figure 1 demonstrates the planned actions from IMO to cut the greenhouse gas emissions from marine transport work. This strategy’s urgent goal is to phase hydrogen-powered vessels in as quickly as feasible among the other short- and long-term candidate initiatives for this century. In particular, the energy efficiency design index for new ships will be implemented in additional phases to reduce the carbon intensity from marine vessels. As shown in Figure 1, the goal is to decrease the carbon intensity of international shipping by reducing greenhouse gas emissions per unit of transport work. The target is to achieve a reduction of at least 40% by 2030 compared with the 2008 level. The ultimate objective is to reduce total annual emissions by at least 50% by 2050 compared with the 2008 levels (according to IMO). Efforts are being made to pursue a pathway toward zero-emission marine transport. This involves not only reducing emissions but also working toward phasing hydrogen-powered vessels out entirely. The vision is to move toward a future where marine transport operates with zero emissions.

Figure 1. A decade of action for IMO to cut emissions from marine transport work. GHG: greenhouse gas.

Furthermore, according to the IMO policy, “the emissions estimates demonstrate that gains in efficiency are crucial in containing emissions growth, but even the most considerable improvements modeled do not result in a downward trend. If fossil fuels continue to predominate, changes in the fuel mix have a minimal effect on greenhouse gas emissions compared to regulatory or market-driven efficiency improvements.”

To achieve this goal, liquid hydrogen ${(}{\text{LH}}_{2}{)}$ fuel cell technology is widely implemented on short-ocean marine vessels and passenger ferries. Meanwhile, the manufacturers are accelerating the progress toward developing zero-emission marine vessels. Figure 2(a) shows the first ${\text{LH}}_{2}$-powered ferry in the world, which is operated in Norway, making history when it began using hydrogen fuel that emits no emissions. With its 80 cubic meter hydrogen storage tank, the ferry should be able to reduce its carbon emissions annually by up to 95%. Figure 2(b) demonstrates the Switch passenger ferry in America, which has a capacity of seventy-five passengers and a maximum cruise speed of twenty-one knots thanks to the 410 kW ${\text{LH}}_{2}$ PEMFCs and lithium-ion batteries on board. Table 1 summarizes the specifications of marine vessels and passenger ferries. The Zero-V research vessel is shown in Figure 2(c), and Figure 2(d) shows the Polish ferry Unity Line, which will be operated in the Baltic Sea in the future.

Figure 2. Hydrogen fuel cell-powered marine vessels: (a) the MF Hydra in Norway; (b) the Switch passenger ferry in the United States; (c) the Zero-V research vessel in the United States; and (d) the Polish ferry Unity Line, which will be operated in the Baltic Sea.

Table 1. The marine vessel and passenger ferry specifications.

Among the rapidly developing large-scale modular technologies, PEMFCs have emerged as one of the most appealing power sources for marine vessels. This is primarily due to the following key advantages they offer:

Environmentally friendly: PEMFCs are environmentally friendly power sources as they only emit water and oxygen during their operation, which supports the achievement of sustainable and net-zero-emission marine transport. By utilizing PEMFCs, vessels can significantly reduce their carbon footprint and contribute to environmental conservation.
Highly efficient: PEMFCs exhibit high efficiency, offering an energy conversion efficiency of approximately 60%. When considering the overall fuel cell systems, their efficiency ranges from 45% to 55%.
Long distance to empty (DOE): PEMFCs are particularly suitable for medium- and long-range marine transport due to their extended DOE capability. Unlike lithium-ion batteries, PEMFCs do not require additional charging processes after operations. This advantage allows marine vessels equipped with PEMFCs to cover longer distances without interruptions for recharging, thereby improving operational efficiency.
Stability: PEMFCs have highly stable performance and a mature integrated liquid cooling system, which contribute to impressive thermal performance and reliable power generation and enhance the overall durability and longevity of the energy storage system.

There are three common models used for PEMFCs: equivalent-circuit-based models, physics-based models, and machine learning models. At the PEMFC layer level, various dynamic models have been developed to investigate the electrical behavior and thermal performance (Gao et al. 2012). The equivalent-circuit-based model has undergone improvements, particularly regarding the heat component for individual cells. To incorporate fuel cell degradation into the equivalent-circuit model, resistance and capacitance elements are combined. These models provide valuable tools for understanding and studying the performance and behavior of PEMFCs.

In medium- and high-voltage applications such as marine vessels, the device-level models play a critical role in the interconnection of power converters and energy systems. Hierarchical HIL emulation is more critical in operating a modular multilevel converter (MMC) (Dinavahi and Lin 2021), which facilitates the interaction between multiterminal energy and power and a high-power electric propulsion system. The increasing number of energy storage devices naturally brings new challenges for the HIL emulation scheme to maintain in real time.

The hierarchical real-time hardware emulation approach is now regarded as a useful tool for the power system of marine vessels during the initial design and testing phases. Comprehensive dynamic device-level models, composed of PEMFCs, batteries, and power switches, are built to describe behavioral transients during operation. The development of system-on-a-chip and multicore technologies has led to significant advancements in HIL and software-in-the-loop emulation over the past ten years. Enhancing multidomain modeling for devices in subsystems and putting in place a cutting-edge computing platform for real-time hardware emulation of zero-emission maritime vessels are crucial. The following is a list of the key expectations:

The multidomain models provide explicit representations of dynamic behavioral transients in PEMFCs through ordinary differential equations from the perspectives of electrical, fluid, and thermal domains. These ordinary differential equations mathematically couple the interactions among these domains during the solution process. In addition, device-level modeling of power switches is implemented for real-time emulation of power converters in control systems.
The hardware emulation scheme is improved for representing both device-level and system-level features of the energy storage and electric propulsion systems in marine vessels. The proposed hierarchical HIL scheme enables the development of detailed models for PEMFCs to explicitly present device-level steady-state and dynamic characteristics in real time; meanwhile, the system-level emulation is carried out at the same time.
Xilinx Versal ACAP offers exceptional compute performance with embedded artificial intelligence engines (AIEs) and programmable logic (PL) arrays. These platforms facilitate the mathematical separation and parallel solving of device-level and subsystem models.

In this article, the Zero-V marine research vessel is used as an example for HIL emulation. The comprehensive PEMFC model is validated with the empirical model from Simulink. The hardware emulation is performed on the Xilinx VCK190 board with the ACAP device XCVC1902. The results indicate that the proposed hierarchical hardware emulation effectively represents detailed information at the device and system levels during steady-state and transient conditions.

Marine Vessel Power Systems

Figure 3 shows an example of onboard power systems for marine transport, including power converters, electrical propulsion systems, and energy storage systems. The electrical propulsion system provides the main thrust and maintains the position of the ship, while the PEMFC and battery system is the only power supply and energy storage. In this section, these main subsystems are introduced and analyzed.

Figure 3. The diagram of hydrogen and battery hybrid marine vessel: (a) diagram for electrical systems of hybrid marine vessel, (b) power converter and electric propulsion system, and (c) onboard energy storage systems.

Electric Propulsion System

The marine vessel is equipped with stern thrusters in each outer hull and a retractable azimuth bow thruster to facilitate precise position-keeping during on-station science activities. To optimize the steering forces generated by the main propellers during station preservation, high-lift flap rudders are also incorporated. This propulsion system provides sufficient maneuverability and positioning capability for the ship under the required operating conditions, and the specific powertrain parameters can be found in Table 2. The 500-kW motors possess sufficient reserve power for safe operation in rough sea conditions, dynamic positioning, and meeting the diverse mission requirements. The permanent magnet motor (PMSM) has been chosen for the electric propulsion system as it offers advantageous characteristics such as effective and quiet operation when directly connected to the propeller shaft.

Table 2. The Zero-V marine research vessel powertrain.

In the case study, conventional tunnel thrusters, which are based on PMSMs, were utilized for the stern thrusters (Klebanoff et al. 2018).

The Zero-V marine vessel operates at a dc bus voltage of 900 V. Power converters are required for the energy storage devices to regulate the output voltage. Furthermore, considering that the efficiency of PEMFCs is dependent on the operating voltage and current, it will be more efficient for the marine vessel to operate the PEMFC stack at its peak efficiency. Therefore, the controller is modified to consider the efficiency constraint, ensuring that the PEMFC stack operates at its optimal efficiency level.

The onboard power systems, which are made up of thousands of energy storage devices coupled with power converters, are the only energy sources for zero-emission vessels, while the electric propulsion systems are usually operated by electric machines such as induction and synchronous motors. Many models for energy storage devices only consider electrical properties. To offer electrical and thermal data for energy management techniques, it is crucial to build multidomain models for batteries and fuel cells.

Energy Storage System

Table 3 compares the PEMFC with lithium-ion batteries in terms of power density, rated current, operating temperature, and peak efficiency. The batteries have better power density, thermal performance, and peak efficiency, while a single PEMFC cell can deliver higher operating currents. In the marine vessel, the PEMFC plays the role of main energy storage, while batteries are energy buffering to improve the overall benefits of energy storage systems.

Table 3. Comparison of energy storage devices.

As shown in Figure 3(c), the energy storage system is based on the PEMFCs and batteries. PEMFCs are versatile and can be applied to both medium- and low-voltage ac and dc power systems. They can be explored in combination with batteries or other energy storage devices, allowing for various system configurations. The integration of PEMFCs into vessels can be done in different ways, including fully hydrogen-electric systems or as part of hybrid power systems. By incorporating the PEMFC power plant, marine vessels can significantly extend their running hours in zero-emission mode. Moreover, shore charging infrastructure provides even more flexibility and energy storage capabilities.

In this case study, batteries are introduced to the onboard energy storage system since batteries have better power density than PEMFC, as shown in Table 3. Batteries play an important role in energy buffering, improving the maximum output current capacity of the entire system. In addition, fuel cells have high recyclability, contributing to the overall sustainability of the onboard power system.

Thanks to the onboard fuel tank and ${\text{LH}}_{2}$ circulation systems, direct trailer refueling has been identified as the preferred strategy, where ${\text{LH}}_{2}$ trailers dock alongside the vessel and connect to a stationary fueling stanchion. This method allows for safe and efficient refueling without the need for the vessel to move during the process. While initially expensive, having a large and reliable customer for renewable hydrogen fuel would encourage suppliers to make renewable ${\text{LH}}_{2}$ more widely available nationwide. This demonstrates the potential for sustainable and renewable hydrogen as a fuel source for marine transportation applications.

Multidomain Modeling and Real-Time Emulation

In this section, the multidomain modeling for PEMFC and real-time HIL emulation are developed and implemented on a Versal ACAP. The main PEMFC modeling is introduced, and the hardware platform is also introduced.

Proton Exchange Membrane Modeling

As the main energy source in the marine vessel, the PEMFC stacks need to be precisely modeled in the electrical, fluid, and thermal domains to reveal detailed information for energy management systems. Figure 4 demonstrates the insight into a single PEMFC on the Zero-V marine vessel, where the hydrogen fuel goes through the anode channel to the cathode catalyst. Then ${\text{H}}^{+}$ ions can cross the proton exchange membrane to the cathode catalyst, where the chemical reaction occurred. The extra oxygen and water come out from the cathode channel.

Figure 4. A diagram of a PEMFC.

Electrical Domain

In general, the output voltage of PEMFC ${(}{U}_{\text{FC}}{)}$ can be calculated as follows: \[{U}_{\text{FC}} = {E}_{\text{FC}}{-}{U}_{\text{ohm}}{-}{U}_{\text{act}}{-}{U}_{\text{con}}\] where ${E}_{\text{FC}}$ is the is the ideal potential, ${U}_{\text{ohm}}$ is the ohmic loss, ${U}_{\text{act}}$ is the activation loss, and ${U}_{\text{con}}$ is the loss caused by convection and diffusion naturally caused by the electrical reaction.

Actually, the PEMFC model is a time-varying and nonlinear model that requires a nonlinear solver during the device-level emulation. To accelerate the HIL emulation and achieve real-time emulation, the advanced software architecture and hardware platform are implemented in this article.

Fluid Domain

The Maxwell–Stefan equations are widely used for description of the diffusion phenomenon, where diffusive fluxes, ${P}_{i}$, of species through a plane depend on all independent driving forces. The gas diffusion of the ith species in the layers is described by \[{\Delta}{P}_{i} = \frac{{\delta}\,{R}\,{T}}{{P}_{t}\,{S}} \mathop{\sum}\limits_{{j}\,{\ne}\,{i}} \frac{{P}_{i} \frac{{q}_{j}}{{M}_{j}}{-}{P}_{j} \frac{{q}_{i}}{{M}_{i}}}{{D}_{ij}}\] where ${\delta}$ is the gas diffusion layer thickness, ${P}_{t}$ is the total gas pressure, ${S}$ is the gas diffusion layer, ${M}_{j}$ is the gas molar mass, j stands for species other than the ith species, and ${D}_{ij}$ is the binary diffusion coefficient between the ith and jth species.

Thermal Domain

Figure 5 illustrates a thermal network for PEMFC, which incorporates four heat sources and two resistance–capacitance couples to account for convection and radiation effects. In addition, the main thermal energy of a PEMFC cell is the difference between the chemical energy of hydrogen and oxygen and the electrical energy of the stack. This thermal energy is distributed as the sensible heat of the coolant and reactants, the latent heat during the phase change of the water, and heat loss to the surroundings by convection.

Figure 5. The thermal network model for a PEMFC.

The main internal resistance of the PEMFC arises from the variable resistance of the membrane, which is noted as ${R}_{\text{mem}}$, which leads to an irreversible voltage drop. As power is generated through the electrochemical reactions in the PEMFC, heat is also produced. These phenomena can be quantitatively calculated using specific equations or mathematical models. The abovementioned phenomenon can be calculated by \[{Q}_{\text{FC}} = {I}_{\text{FC}}^{2}{R}_{\text{mem}} + {I}_{\text{FC}}\,{\cdot}\,{\left({\Delta}{U}_{\text{act}}{-}\frac{{T}_{C}\,{\cdot}\,{\Delta}{S}}{2F}\right)} + \frac{{\lambda}_{\text{mem}}\,{\cdot}\,{A}_{\text{mem}}\,{\cdot}\,{\left({T}_{\text{mem}}{-}{T}_{c}\right)}}{{d}_{\text{mem}}}\] where ${R}_{\text{mem}}$ is the equivalent resistance of the proton exchange membrane, ${\Delta}{S}$ is the entropy changes, ${\lambda}_{\text{mem}}$ is the membrane’s material thermal conductivity, ${A}_{\text{mem}}$ is the contact area between the membrane and catalyst of the cathode side, and ${d}_{\text{mem}}$ is the layer thickness. ${T}_{\text{mem}}$ and ${T}_{c}$ are the membrane and catalyst temperatures, respectively.

The electrical, fluid, and thermal domains are mathematically coupled with each other since the electrical behavior calculations require fluid and thermal information. In this case, it is important to advance the HIL scheme to achieve real-time device-level emulation.

Real-Time HIL Emulation

Figure 6 depicts the architecture of the Xilinx ACAP board, which comprises various components such as AIEs, PL, and an Advanced RISC Machine (ARM) core integrated into the chip. The AIE array serves as a higher-level hierarchy within the AIE architecture and consists of a two-dimensional array of AIE tiles. The AIE can establish connectivity with other parts of the device either through the network-on-chip or by directly interfacing with the PL via the AIE array interface.

Figure 6. The architecture of the Xilinx Versal ACAP. PL: programmable logic, RISC: reduced instruction set computer.

Each AIE tile, as shown in Figure 6, includes an engine capable of accessing at most four memory modules, a dedicated memory module itself, and a tile interconnect module. The ACAP kit used in this article has a maximum capacity of four hundred AIE tiles. The AIE processor is highly optimized to support both fixed- and floating-point precision.

Programming the AIE involves a two-stage approach within the Vitis environment, which includes kernel programming and graph programming. The communication between the AIEs and other components of the device is represented as a graph, where kernels are instantiated and connected through buffers and streams. This approach facilitates efficient and scalable programming of the AIE for implementing complex algorithms and computations.

General-purpose computing is provided by the ARM processing unit, which is based on the ARM Cortex-A72 CPU core. It is suitable for solving the multidomain PEMFC model due to its capabilities and high clock frequency. OpenCL is adopted for software programming, enabling parallel execution of multiple kernels with an initialized command queue, resulting in improved speed. The PL is a flexible structure that allows for the creation of various hardware functionalities. It consists of adjustable logic blocks, configurable RAM, and digital signal processor engines. These components can be interconnected and combined to generate a wide range of functions, including processors, functional pipeline units, peripherals, and accelerators. The PL provides the flexibility to tailor the hardware configuration to specific requirements.

In the design of the energy storage system, the global memory input/output port plays a crucial role. This port serves as a connection point for external memory, which can be mapped to or from the main memory of the system. It enables efficient and seamless data transfer between the system and external memory sources. This connection is essential for effective memory management and enables the system to handle large datasets efficiently.

To set up and establish connections among the components within the PL and ARM cores, specialized software tools such as the Vivado design suite and the Vitis integrated software platform toolset are utilized. These tools provide a comprehensive environment for specifying the configuration, connectivity, and interoperation of the PL components. They offer features and functionalities that facilitate the design, development, and optimization of the PL and ARM cores.

By utilizing these software tools, designers can define and configure the PL components, specify their connectivity, and optimize their performance. The tools enable the creation of a programmable device image that can be loaded onto the hardware to implement the desired functionality. This image represents the configured PL components and their interconnections, allowing the hardware to execute the required operations as defined in the design.

Real-Time HIL Emulation Results and Discussion

In the case study, the configuration and specifications of the Zero-V marine vessel are adopted in this article. The propulsion switchboard receives electricity from each PEMFC and battery pack via a dc–dc converter, which changes the voltage to a constant nominal voltage of 900 V dc. The propulsion switchboard supplies power through dc–ac drives to the different heavy loads, including the propulsion, thruster, and winch motors. Moreover, redundant dc–ac power converters provide ship service electrical power to the 480 V ac ship service switchboard. The ship service switchboard supplies smaller loads like lights, fans, and pumps. The time step for the HIL emulation of PEMFCs and power converters is ${0.5}{\mu}{\text{ s}}$.

Figure 7 illustrates the operations of power converters and multidomain modeling of the onboard PEMFC in the electrical and fluid domains. In Figure 7(a), the three-phase output currents of the MMC are shown using an average model under preset operating conditions. The PMSM torque and speed are generated according to the preset operating conditions of the marine vessel. Then, these currents serve as inputs for device-level emulation.

Figure 7. The PEMFC modeling and HIL emulation results: (a) the inverter current of the inverter between the PEMFC and the PMSM, (b) the PEMFC current density versus the output voltage, (c) the PEMFC current density versus the PEMFC output power, (d) the PEMFC efficiency, (e) the power consumption of the PEMFC during varying LH₂, and (f) the PEMFC power loss during varying LH₂.

Figure 7(b) presents the basic relationship between the current and the output, while Figure 7(c) depicts the current and output power curve in steady state during PEMFC operations. In Figure 7(b), the voltage output experiences a significant jump/drop when the current density falls below 0.18 or exceeds 1.03. Additionally, the output power increases with increasing current density but starts to decline when the current exceeds 1.03 times the rated current. In addition, the output power sharply decreases when the PEMFC is overloaded. This behavior is attributed to the negative effects of overcurrent on the electrical performance of PEMFCs.

Figure 7(d), (e), and (f) demonstrates the efficiency, power consumption, and power loss of the PEMFCs during varying ${\text{LH}}_{2}$ supply. In the presented scenario, the hydrogen fuel supply increases between 1 s and 4 s, leading to a decrease in efficiency and an increase in power consumption and power losses. This phenomenon occurs because the excessive supply of ${\text{LH}}_{2}$ cannot be sufficiently reacted with oxygen, resulting in increased power losses and overall power consumption. Despite the output power remaining relatively constant, the overall efficiency of the PEMFCs decreases throughout the process. Overall, the proposed HIL emulation scheme can reveal the insight into power converters and PEMFC from the electrical, fluid, and thermal domains in real time.

Conclusions

Significant innovations in zero-emission energy storage, such as lithium-ion batteries and hydrogen fuel cells, have been driving urgent research in zero-emission marine transportation. These innovations play a critical role in reducing emissions across all industry sectors. In this context, it is essential to advance device-level modeling to understand the intricacies of power systems and improve control and management objectives through real-time emulation. HIL emulation plays a crucial role in marine vessel design, and device-level modeling provides valuable information across the electrical, chemical, hydraulic, and thermal domains. The advanced computing platform, ACAP, represents the state of the art in parallel computation, offering significant potential for large-scale device-level emulation and improved computation speed for real-time emulation. Furthermore, the hierarchical nature of HIL emulation allows for mathematical decoupling and separate solving in ACAP, thereby reducing the emulation time. Future research will focus on developing multidomain models for various fuel cell types, including solid oxide fuel cells. Additionally, real-time HIL emulation will be developed for zero-emission transport sections with megawatt-level power consumption and hydrogen-based oceangoing marine vessels Furthermore, hybrid electromagnetic transient and machine learning methods will be investigated for device-level modeling, and parallel-in-time-and-space methods will be implemented for HIL emulation to improve efficiency and maintain accuracy.

For Further Reading

E. Skjong, T. A. Johansen, M. Molinas, and A. J. Sørensen, “Approaches to economic energy management in diesel-electric marine vessels,” IEEE Trans. Transport. Electrific., vol. 3, no. 1, pp. 22–35, Mar. 2017, doi: 10.1109/TTE.2017.2648178.

“IMO’s work to cut GHG emissions from ships,” International Marine Organization, London, U.K., 2023. [Online] Available: https://www.imo.org/en/MediaCentre/HotTopics/Pages/Cutting-GHG-emissions.aspx

H. Ahmadi, M. Rafiei, M. A. Igder, M. Gheisarnejad, and M.-H. Khooban, “An energy efficient solution for fuel cell heat recovery in zero-emission ferry boats: Deep deterministic policy gradient,” IEEE Trans. Veh. Technol., vol. 70, no. 8, pp. 7571–7581, Aug. 2021, doi: 10.1109/TVT.2021.3094899.

F. Gao, B. Blunier, D. Chrenko, D. Bouquain, and A. Miraoui, “Multirate fuel cell emulation with spatial reduced real-time fuel cell modeling,” IEEE Trans. Ind. Appl., vol. 48, no. 4, pp. 1127–1135, Jul./Aug. 2012, doi: 10.1109/TIA.2012.2198909.

V. Dinavahi and N. Lin, Real-Time Electromagnetic Transient Simulation of AC-DC Networks, 1st ed. Hoboken, NJ, USA: Wiley/IEEE Press, 2021.

L. Klebanoff et al., “Feasibility of the zero-V: A zero-emission, hydrogen fuel-cell, Coastal research vessel,” Sandia Nat. Lab., Livermore, CA, USA, Rep. SAND2018-4664, 2018.

E. Chemali, M. Preindl, P. Malysz, and A. Emadi, “Electrochemical and electrostatic energy storage and management systems for electric drive vehicles: State-of-the-art review and future trends,” IEEE J. Emerg. Sel. Topics Power Electron., vol. 4, no. 3, pp. 1117–1134, Sep. 2016, doi: 10.1109/JESTPE.2016.2566583.

Biographies

Chengzhang Lyu (clyu1@ualberta.ca) is with the University of Alberta, Edmonton, AB T6G 2H5, Canada.

Venkata Dinavahi (dinavahi@ualberta.ca) is with the University of Alberta, Edmonton, AB T6G 2H5, Canada.

Digital Object Identifier 10.1109/MELE.2023.3320509

The Unsung Hero of the Electric Vehicle Revolution: The role of computational electromagnetics in electric machine design and analysis

Dionysios Aliprantis, Steven Pekarek

©SHUTTERSTOCK.COM/PAUL CRAFT

The electrification of transportation hinges on a number of core technologies, such as power electronics and energy storage. In this article, the spotlight is on electric machines (EMs), a term that encompasses all electromechanical energy conversion devices, i.e., motors and generators. Computational electromagnetics, like the finite element method (FEM), play a vital role in a modern EM design workflow. Such methods have been under development for decades, however, with recent advances in computing, they are now mainstream. One of their key features is that they allow us to search a design space with unprecedented accuracy. But perhaps more importantly in today’s rapidly changing landscape, they are accelerating the pace of innovation, reducing prototyping costs and the time it takes to go from concept to product.

Let us consider, for example, a modern electric vehicle (EV). First, we note that the term EV does not apply only to a passenger car. Today, when we speak of an EV, we may be referring to an electric scooter or bicycle for micromobility, a zero-emissions city bus, a medium-duty box truck for local deliveries, a heavy-duty semitrailer truck for long-haul freight transport, a forklift employed on a factory floor or in a port, a loader used at a construction site, a hauler in a mine, an aircraft or watercraft, and the list goes on, with many products already available in all these markets. The propulsion motors found on board these vehicles are obviously not commercial off-the-shelf devices. For maximum performance and competitiveness in the marketplace, each EM has to be crafted very carefully for the specific requirements of a given vehicle, while accounting for efficiency, weight, size, cost, reliability, and the vehicle’s drive cycle and operating conditions. In addition to the main propulsion motors, we should not forget about all the other actuators that are also being electrified by virtue of having access to an underlying EV power system architecture.

The EM design problem is far from trivial. Even though EMs have been around for more than a century, modeling their behavior remains a challenge, especially when they are connected to power electronic (switched) circuits. This is due to a variety of reasons, but mainly due to their complicated structures as well as the nonlinearities of ferromagnetic materials and our imperfect modeling of their physics (e.g., hysteretic effects). It could be argued that the basic dimensions of an EM can be obtained using relatively simple analytical formulas found in textbooks dating back to the early 20th century. Nevertheless, such designs serve only as an initial point for further fine-tuning. We now rely on sophisticated computer-aided design software that can predict the electromagnetic, mechanical, and thermal behavior of an EM very accurately over its entire operating range. We can thus incorporate every geometric detail, no matter how minute; however, it should be understood that our models are also limited by uncertainties in material parameters and other assumptions.

The Finite Element Method

The FEM is an established and widely used technique in a variety of disciplines. We may think of the FEM as a numerical microscope, allowing us to examine what is transpiring inside EMs over time. For instance, the FEM can be employed to establish the mechanical properties of an EM to ensure the structural integrity of a high-speed rotor, or to predict noise, vibration, and harshness (NVH) characteristics. In what follows, we describe the electromagnetic FEM in more detail, providing a high-level explanation of what it does and how it works in the context of EMs.

Simply speaking, the electromagnetic FEM solves Maxwell’s equations, which are partial differential equations (PDEs) quantifying the variation of the electromagnetic field through space and time. EMs are relatively low-frequency devices, so we can safely ignore electromagnetic radiation effects. We thus have a simpler quasi-magnetostatic problem, where the unknown is the magnetic field throughout a finite region in space surrounding the device. It is possible to formulate the FEM in either two or three dimensions. In two dimensions, it is assumed that nothing changes along the third axis. Luckily, because of how EMs are built, a 2D model often yields a decent approximation of the fields inside the active region of the EM, which is where the ferromagnetic material resides (in contrast to the end winding region, which has to be analyzed separately). A 3D model can be more accurate, but this comes at a significantly higher computational cost.

Physics tells us that there are two sources for the magnetic field. First, we have free current sources due to the flow of electrons inside conductors, such as current in windings or induced eddy currents in conductive regions. Second, we must account for bound current sources; these are more difficult to conceptualize, but they essentially represent the magnetization of matter—and in EMs, they are the dominant source because their effect is what creates the strong magnetic field required to make them function properly. In the FEM, bound currents are not modeled explicitly; rather, their presence is captured indirectly through the magnetization (B–H) characteristics of magnetic materials.

To initiate a finite element analysis (FEA) study, a user has to input the geometry. In a 2D analysis, this corresponds to a cross section of the EM on the plane of interest (i.e., perpendicular to the shaft). In the most elementary fashion, a geometry can be described by two lists: 1) a list of nodes and their coordinates and 2) a list of linear segments interconnecting pairs of nodes. This is similar to how a connect-the-dots game looks like. For convenience, the geometry can also be described through a graphical interface (but the software then creates internally lists of nodes and segments) or even through an application programming interface (which facilitates parametric studies). Curved surfaces can be approximated by polygonal shapes. Subsequently, a user has to define the material type occupying the various regions in space (e.g., air, copper, steel, or permanent magnets) and the electromagnetic parameters of each material (e.g., using a constant magnetic permeability or a B–H curve). After the geometry definition is complete, a user has to specify the variation of the winding currents and the rotor position over time. For example, in a three-phase permanent magnet synchronous machine (PMSM) operating in steady state, the currents can be defined as a balanced set, and the rotor position could be a linear function of time. A separate FEA would then be conducted for each time instant and mechanical position. In more advanced FEM implementations, the variation of the currents and/or rotor position can be determined dynamically by interfacing the FEA with an external circuit and/or mechanical subsystem transient solver. The FEM can also be used for an analysis of linear EMs, where the mechanical degree of freedom is a linear position. Finally, boundary conditions are set, which capture the impact of external magnetic fields (if any). This information is really all that is needed to conduct an FEA study.

What happens under the hood of an FEA program? The first thing that takes place is a partitioning of the geometric domain, which is subdivided into a high number of relatively small finite elements. This action is called meshing. Typically, triangular elements are used for 2D EM FEA because they work fairly well for meshing the common internal features of EMs (e.g., teeth, slots, pockets, or curved surfaces). For a 3D FEA, the finite elements could be, for instance, tetrahedral. There could be tens of thousands of such elements. Commercial FEA software employs sophisticated meshing algorithms, which are intelligent enough to determine which regions need to be meshed with a higher density than others. In general, a decent mesh has more finite elements concentrated in regions where the magnetic field varies rapidly, while a sparser mesh can be used where the field is more constant. (Typically, we avoid meshing the domain with a uniformly small finite-element size because doing so increases the computational burden.) Various mathematical results exist that correlate the size of the finite elements with the error between the FEM solution and the actual solution. An example mesh is shown in Figure 1, depicting an elementary surface-mounted PMSM. To take advantage of the periodicity of the magnetic field, we meshed a single pole (30 mechanical degrees).

Figure 1. FEA mesh and magnetic flux lines for a surface-mounted PMSM. PM: permanent magnet.

Once a mesh has been established, the FEM proceeds with numerical calculations to obtain the magnetic field. We do not describe these details here. The key idea is that the FEM searches among a set of functions (e.g., piecewise linear functions) that are defined on the given mesh. Interestingly, we do not work directly with Maxwell’s PDEs; rather, by using vector calculus identities, we obtain an equivalent variational calculus formulation. The new objective is to identify the function that minimizes an energy-related integral over the domain. So, essentially, we are solving an optimization problem rather than the original PDE.

At the end of the day, the FEM assembles a system of equations where, typically, the unknown is the magnetic vector potential (MVP), denoted as A. These equations become nonlinear if magnetic saturation characteristics of materials are modeled. Specifically, the unknown of the underlying problem is an array of MVP values. If 2D linear triangular elements are used, where a linear variation of the MVP is assumed over each element, the unknown MVP values correspond to potentials at the vertices of the finite elements. The MVP is a physical quantity that is analogous to the electric potential, V or $\varphi$, which is perhaps a more familiar scalar field related to the distribution of electric charge around space. The sources of the MVP are free and bound currents.

Finally, as a postprocessing step, the FEA program will compute various quantities of interest and will generate plots as needed. Typical outputs are the flux linkages of the coils (see Figure 2), electromagnetic torque (see Figure 3), various losses (ohmic and magnetic), and the magnetic forces acting on the inner parts of the device (that inform a subsequent NVH analysis). The visualization aspect is important in its own right because it can help identify issues and possibilities for improving a design. Magnetic flux lines are directly obtained as contour lines of the MVP (see Figure 1). The magnetic flux density (i.e., the B-field) can be readily found as the curl (an operator from vector calculus) of the MVP; when plotted over the domain, it can highlight potential saturation or high-loss regions.

Figure 2. Variation of the a-phase flux linkage (per unit of length in the axial direction) with rotor position. MoM: method of moments.

Figure 3. Variation of electromagnetic torque (per unit of axial length) with rotor position. MoM: method of moments.

The Method of Moments

Although the FEM is perhaps the most well-known approach, there exist alternative numerical approaches to solve for magnetic fields within low-frequency electromagnetic components. These are often placed under the umbrella terms method of moments (MoM) or integral-based methods, of which there are some different flavors. Here we highlight one in particular that we refer to as an MoM, which connects nicely with the earlier description of magnetic materials and with the FEM. Specifically, one can consider that a magnetized material has bound currents that, in addition to the free currents in conductors, act to create magnetic flux. It is these bound currents that provide, for example, the remanent flux that resides in and around a material when there is no externally sourced magnetic field applied. These bound currents also change (respond) when a material is subjected to an externally generated magnetic field. One can show that if a material has uniform permeability, the bound currents reside only on the surface. If the permeability is nonuniform, bound currents also reside within the material volume.

Derivation of the considered MoM starts by asking a basic physics question: What is the value of the magnetic flux density everywhere in a system? It is then recognized, from Ampère’s law, that the flux density at any location can be determined by summing the contributions from all free and bound currents. The flux density within a magnetic material can also be obtained with knowledge of the magnetization and the relative permeability. Finally, there is a straightforward relationship between magnetization and the bound current at the surface of the region. Taken together, these relationships are used to generate an MoM model in which the free currents serve as an input, and the bound currents (or the tangent component of magnetization) on or within magnetic materials are the unknowns.

One may ask, “Where does the discretization take place?” Any region in the system that has a relative permeability greater than one has bound currents and must be discretized. When considering an EM, the stator and rotor iron (and permanent magnets) must be meshed. In the easiest of cases, if the machine is modeled in 2D and any of these materials is assumed to have uniform permeability throughout their respective regions, the mesh consists of line segments that outline the edge of the respective materials. In Figure 4, a line-segment-based MoM mesh is shown for a surface-mounted PMSM. If the materials have nonuniform permeability, which often needs to be represented if a region is saturated, the interior of the material must be meshed. This can be done using the same discretization geometries (e.g., triangles) that are used in an FEA, although other polygons can also be applied.

Figure 4. MoM mesh and magnetic flux lines for a surface-mounted PMSM. PM: permanent magnet.

Looking under the hood of an MoM-based model, one will find the critical postprocessing functions needed to compute winding flux linkage (see Figure 2) and electromagnetic torque (see Figure 3). Similar to the techniques applied in an FEA, flux linkage is determined from the value of the MVP at the winding conductors. In a 2D MoM, one can leverage a common expression for the MVP from filamentary currents in free space to predict the MVP at the stator conductors resulting from all free and bound currents. Because the MoM formulation is based on determining flux everywhere in a system and is structured to solve for bound currents, the forces acting throughout the stator and rotor and the electromagnetic torque are readily predicted by calculating the Lorentz force acting on each element.

One can observe in Figures 2 and 3 that the FEA and MoM responses are close, but they are not exact matches. The differences are attributed to several factors, including that in the MoM model, only a surface mesh was applied and the relative permeability is assumed constant, whereas the FEA model employs a nonlinear B–H curve. In addition, in the MoM-based model, stator slot current is modeled using a single filamentary conductor at the center of the slot. In contrast, within the FEM-based model, the conductors are partitioned by phase into the top and bottom regions of the slot, and the respective current density is modeled across the respective region.

Conclusions

As electrification continues to spread across land-, sea-, and air-based transportation, computational electromagnetics will undoubtedly continue to play a central role. A search of recent literature shows that a key goal is to leverage integrated multiphysics (electromagnetic, thermal, and mechanical) models to enhance machine and system performance. Nevertheless, predicting behavior in space and time (for a true 4D analysis) over different physical domains is still very demanding computationally as well as from a modeling perspective. Artificial intelligence-based techniques, such as deep neural networks that can be trained to learn the EM response surface, and next-generation (e.g., quantum) computing are likely to serve as key tools to enable widespread use of 4D models in multiobjective optimization and real-time control.

For Further Reading

P. P. Silvester and R. L. Ferrari, Finite Elements for Electrical Engineers, 3rd ed. Cambridge, U.K.: Cambridge Univ. Press, 1996.

J. P. A. Bastos and N. Sadowski, Electromagnetic Modeling by Finite Element Methods. New York, NY, USA: Marcel Dekker, 2003.

N. Bianchi, Electrical Machine Analysis Using Finite Elements. Boca Raton, FL, USA: CRC Press, 2005.

D. Aliprantis and O. Wasynczuk, Electric Machines: Theory and Analysis Using the Finite Element Method. Cambridge, U.K.: Cambridge Univ. Press, 2023.

T. C. O’Connell and P. T. Krein, “A time-harmonic three-dimensional vector boundary element model for electromechanical devices,” IEEE Trans. Energy Convers., vol. 25, no. 3, pp. 606–618, Sep. 2010, doi: 10.1109/TEC.2010.2042811.

R. Howard and S. Pekarek, “Two-dimensional Galerkin magnetostatic method of moments,” IEEE Trans. Magn., vol. 53, no. 12, pp. 1–6, Dec. 2017, doi: 10.1109/TMAG.2017.2741449.

Biographies

Dionysios Aliprantis (dionysis@purdue.edu) is with Purdue University, West Lafayette, IN 47907 USA.

Steven Pekarek (spekarek@purdue.edu) is with Purdue University, West Lafayette, IN 47907 USA.

Digital Object Identifier 10.1109/MELE.2023.3320510

Electromagnetic Transients Simulation Program: A unified simulation environment for power system engineers

Jean Mahseredjian, Mohammed Naidjate, Mehdi Ouafi, Juan Antonio Ocampo Wilches

©SHUTTERSTOCK.COM/ROMAN ZAIETS

This article introduces a comprehensive and unified environment that enables the exploration of power systems across various modeling levels, encompassing load-flow analysis and extending to intricate time-domain simulations featuring inverter-based resources (IBRs).

Initializations of IBR circuits and parallel computations are presented for practical system problems. In the second part of this article, it is demonstrated that the time-domain accuracy level can be further increased through detailed and integrated modeling of magnetic circuits in the proposed unified environment.

Unified Environment, From Load-Flow to Time Domain

The Electromagnetic Transients Program (EMTP) is designed for the simulation of electromagnetic transients. Since its computational engine is circuit based, it allows simulations for phenomenon frequencies ranging from dc to 100 MHz and more. It is possible to expand into a large spectrum of power system analysis applications.

The main concept is the usage of an expandable system of equations (ESE) that accommodates various circuit models through generic relations. The ESE is based on modified augmented nodal analysis (MANA) created for the first time for the development of EMTP code. It allows formulating a generic system of equations stated as \[{\bf{Ax}} = {\bf{b}} \tag{1} \] where the matrix ${\bf{A}}$ is the matrix of coefficients, ${\bf{x}}$ is the vector of unknowns, and ${\bf{b}}$ is the vector of known values. The matrix ${\bf{A}}$ consists of the classic nodal admittance matrix but augments it with various other model equations, such as ideal transformers, ideal or nonideal switches, independent source functions, dependent source functions, and any other type of equation derived for a model.

The time-domain solution of (1) is based on converting differential equations of models into the well-known companion circuit representation through numerical integration. For power electronics converters, the ESE allows the modeling of switching devices, such as insulated-gate bipolar transistors (IGBTs) or diodes, through ideal switch models or detailed nonlinear functions. Conversely, controlled sources can be employed to furnish the converters with average value models, depending on the modulation index. In the time-domain, the augmented part of ${\bf{A}}$ also contains the relations for nonlinear functions, representing, for example, a nonlinear inductance function for the transformer magnetization branch or IGBT nonlinear characteristic. The matrix ${\bf{A}}$ is time variant and must be iteratively updated to account for switch position changes or changes in nonlinear function slopes. In fact, ${\bf{A}}$ is the Jacobian matrix used in Newton’s method.

The ESE formulation method is used, furthermore, for performing multiphase unbalanced load-flow simulations with phasors. Such simulations allow the study of operating conditions in transmission, distribution, and microgrid networks. The ESE for load-flow can easily accommodate traditional load-flow constraints and integrate more recent constraint (such as droop functions, for example) requirements for power electronics converters. In the load-flow solution, (1) is a complex equation, with the augmented part of ${\bf{A}}$ now containing load-flow constraint equations. Once again, ${\bf{A}}$ becomes the Jacobian matrix and must be updated iteratively until convergence for a given tolerance. It has been shown that this formulation is more efficient and converges much faster than traditional mismatch-based load-flow solvers and various other methods available in the literature. Its advantages are remarkable for ac/dc systems, including microgrids in connected or islanded modes of operation. It can accommodate multifrequency conditions in islanded mode. This formulation accommodates arbitrary balanced or unbalanced multiwire networks. The time-domain models are seamlessly compatible with load-flow equations, facilitating a direct transition into the steady-state solution that serves as the initialization for the time-domain models. Contrary to the load-flow solution, the steady-state solution of (1) uses lumped circuits for all components based on initialization from load-flow solution phasors.

The multiphase load-flow solver is used to initialize IBRs for achieving a quick steady state for time-domain waveforms from which various analyses can be performed. This is an essential aspect not only for dramatically reducing computing times but also for finding the correct and stable steady-state solution. A unified environment is realized for modeling IBRs in the steady state and time-domain, with all related control systems. The computational burden of time-domain simulations with IBRs is alleviated through parallel processing methods.

The simulation of very large grids with IBRs is a challenging problem. At first, it is necessary to use advanced graphical interfaces to present system schematics and enter data for components. A schematic editor must be used to simplify the tasks of the study engineer. Hierarchical arrangements are necessary to achieve high-level visualization. Detailed component model information is encapsulated within multilevel subcircuits, providing a structured approach for data entry through masks.

A sample studied system is presented in Figure 1. The grid-forming converters of this system are initialized from the shown load-flow solution. IBRs can effectively benefit from load-flow PV or PQ constraints. It is possible to also include droop control. Such a system includes all components, with detailed control systems for conducting performance studies with traditional generation sources and modern power electronics components, such as grid-forming converters. Scenario variation tools allow comparing the performances of synchronous generators replaced by grid-forming converters. Sample simulation results are shown in Figure 2 for a three-phase-to-ground fault (see the fault location in Figure 1) appearing at 1 s and disappearing at 1.0833 s. It is apparent that the grid-forming inverter performs better than the synchronous generator for the same fault condition.

Figure 1. A sample test system integrating grid-forming converters with conventional generation.

Figure 2. A performance analysis with synchronous generators replaced by grid-forming inverters. (a) Grid-forming inverter 3, with per-unit active power, replacing synchronous generator 3. (b) Synchronous generator 3, with per-unit active power.

Enhanced unification is achieved by including all the protection relays in Figure 1. Such complete system representation is in line with the digital twin concept since all necessary components are included to create an electronic copy of the studied grid.

A much larger grid model is presented in Figure 3. It contains 10 wind parks (WPs) with aggregated models. In such a system, the WPs can be simulated using different numerical integration time-steps by decoupling them through propagation delay-based transmission line models. This approach is sufficiently accurate due to the interpolation process in the history buffer of transmission lines or cables. This concept is illustrated in Figure 4, where the network is simulated using a 50-µs time-step, with average value models for all IBRs, and only one of the WPs is simulated using a time-step of 10 µs, with detailed modeling of its converter (IGBT nonlinear models). A phase-A-to-phase-B-to-ground fault (see the fault location in Figure 3) is occurring at 2 s, with unsuccessful fault clearing followed by line tripping at 2.25 s. The results are compared to 10-µs usage for the entire system to validate accuracy. Such a setup allows reducing by five times the computing time on two CPUs.

Figure 3. A transmission grid with 10 wind parks.

Figure 4. The active power of a WP, using two-time-step simulation.

Initialization can be further improved through more detailed initialization of IBR control blocks. The basic top-level procedure is given in Figure 5 for a WP. It is applicable also to other types of IBRs. In this version, the ideal voltage source phasor is found from the load-flow solution. It is used in the steady-state solution and then transposed into the time-domain. It stays connected for a duration of ${T}_{i}$ (typically less than 200 ms). This duration allows the IBR to initialize rapidly and reach quasi-steady-state conditions when the source switch is disconnected to interface the IBR directly with its network. The IBR contents must be also initialized automatically, including the initialization of control system blocks and the collector grid.

Figure 5. The WP initialization procedure.

As demonstrated in Figure 6, all the IBRs of Figure 3 achieve steady-state conditions before 200 ms. A phase-A-to-ground fault (see the location in Figure 3) occurs at 1 s, followed by line opening for clearing the fault and then reclosing to recover the system and return to a stable steady state. The entire system is simulated with average value models for WP converters and a time-step of 50 µs.

Figure 6. The active powers of WPs in Figure 3.

The needs for modeling and simulating large-scale power systems in EMT mode are constantly increasing. Recently, the Chilean system operator (Coordinador Eléctrico Nacional) performed EMT simulations for the Chilean power grid. This work was part of a road map to enable the decommissioning of fossil fuel-based power plants and reaching a 100% share of renewable energies as of 2030. This accelerated energy transition scenario must meet the enabling conditions to prepare the grid to integrate new technologies and execute the necessary investments in renewable generation and storage. Network reliability being a key issue, the integration of IBRs requires detailed modeling and analysis for the large-scale Chilean grid. The current EMT (EMTP) model contents are detailed as follows:

WPs: 27, average value generic models
PV parks: 32, average value generic models
Transmission lines/cables: 297 distributed line models (constant parameters) and 307 PI circuits
Synchronous generators: 57, with governor and exciter controls, with saturation data when available
Transformers: 48, with nonlinear magnetization branches
Total number of nodes: 6,785; the size of the MANA (ESE) matrix is 10,664 × 10,664
Control diagram blocks: 57,708.

The above IBR models are generic in preparation for studies of the eventual integration of manufacturer models. The generic models are adequate for testing various performance scenarios and operating conditions.

The computing times for such a system can reach 120 s for a 1-s simulation interval, for a numerical integration time-step of 50 µs when using a single CPU on a laptop computer. This speed can be improved significantly by using several CPUs in parallel. With 60 CPUs, the computing time for 1 s drops to 13 s. A sample simulation waveform for a fault applied at 1 s is shown in Figure 7 to validate the parallel solution.

Figure 7. The active power at one WP at the connection point.

More recently, the ESE approach has been used to integrate the capability to model in detail magnetic circuits. It is presented in the following section.

Magnetic Circuits

Computational electromagnetics aims to develop high-performance numerical models and provides engineers with tools to conceive devices that meet technical requirements and to analyze them in various scenarios. Simulation allows testing different situations that may be challenging to achieve experimentally due to safety concerns, high costs, lack of equipment, or time constraints. Simulations can replicate extreme cases, transient states, and faulty conditions, which deepens the understanding of complex situations and helps prevent or mitigate them.

Historically, the evolution of simulation tools has been closely intertwined with the progress of computer computational capabilities. It has transitioned from basic analytical models and expert-derived formulas to very complex numerical models and, further, to multiphysics simulations. Presently, breakthroughs in cloud computing, edge computing, and artificial intelligence are igniting ambitions toward real-world simulation. This involves modeling the subject device while considering its ultimate environment. Evidence of this current trend is the rising popularity of the concept of digital twins in both academic and industrial spheres. This concept serves to shape the next generation of power grids, characterized by their minimized environmental footprint but also well suited to withstand the evolving impacts of climate changes.

In numerous instances, the device of interest operates within a large system with interactions among multiple components. Thus, it is crucial to adopt a system approach that considers the device’s environment. For instance, when studying faults affecting transformers, it is insufficient to analyze only their impacts on the transformers themselves. Effects on other power grid equipment must also be evaluated. Similarly, transients occurring in the power grid, such as lightning, ferroresonance, and geomagnetically induced currents (GICs), can significantly affect transformers and should be considered during design and operation phases.

While detailed simulators that use classical numerical methods, such as finite elements, are powerful modeling tools, they often consume considerable computing time and focus on devices individually, rendering them incompatible with complete simulations. Several issues are encountered, including very long delays to reach the steady state.

A possible solution is to use a circuit-based simulation tool, such as EMTP. EMTP is a comprehensive simulation software package. Its extensive library empowers users to generate digital replicas of the studied systems. The ESE-based kernel of EMTP incorporates the ability to model the magnetic circuits of devices, such as transformers and rotating electrical machines, using lumped circuit-based methods.

Circuit-based modeling methods offer several distinct advantages, the most important one being seamless integration within complex power systems networks. This integration allows harnessing the benefits of a system approach, enabling a more comprehensive analysis of devices within large systems. One noteworthy advantage is the combination of a magnetic circuit meshing approach with circuit methods, based on a device’s geometry. This fusion is highly valuable for dimensioning and optimization studies, facilitating precise analysis and fine-tuning of designs. Moreover, the mesh approach ensures a generic aspect of the modeling process, leading to the automatic generation of models. This automation streamlines the modeling procedure and enhances efficiency.

Furthermore, it is essential to highlight that the precision achieved by circuit-based methods is comparable to that of classical numerical methods, particularly when the meshing resolution is equivalent. This means that by appropriately adjusting the meshing parameters, circuit-based models can yield results on par with traditional methods, ensuring accuracy and reliability in the simulation process. This attribute further solidifies the value of circuit-based modeling as a viable alternative capable of providing fast, reliable, and precise results while benefiting from seamless integration with circuit-based software.

Diverse circuit-based methods for magnetic circuit modeling have been proposed in the literature, among which the partial element equivalent circuit method and the transmission line matrix method are prominent examples. This presentation adopts the equivalent circuit-based modeling approach for magnetic circuits, including all windings and interfacing with surrounding electric circuits. This choice is reinforced by compelling reasons, including robustness, simplicity of implementation, and excellent suitability for low- and midfrequency applications, which aligns perfectly with the characteristics of most electrical network components.

The following presents an application focused on transformers. Transformers exhibit nonlinear characteristics and operate at diverse power grid levels, subject to various transients. Three circuit-based methods can be selected:

The first method, known as the Hopkinson analogy, embodies the magnetic circuit through reluctances. Mutators are utilized to couple the magnetic circuit with the electrical circuit, ensuring seamless integration between the two domains.
The second approach, termed the Buntenbach analogy, employs capacitors to model the magnetic behavior. This methodology relies on the use of gyrators to accomplish the coupling between the magnetic circuit and the electrical circuit.
Finally, the duality principle employs inductors to create the magnetic circuit and ideal transformers to establish coupling with the electrical part.

The ESE of EMTP accommodates all dependent sources needed for the integration of reluctance network models. Elements such as current-controlled voltage sources (CCVSs), current derivative-controlled voltage sources (CDCVS), current controlled current sources (CCCS), and voltage-controlled voltage sources (VCVSs) can be seamlessly integrated within this framework. These elements, when combined, effectively operate as mutators, gyrators, or ideal transformers. This inherent versatility ensures a simultaneous and strong coupling between magnetic and electrical models. The ESE framework also offers flexibility by including nonlinear components solved iteratively for the highest accuracy and stability. This adaptability aids in optimizing convergence, particularly when dealing with huge numbers of nonlinear elements.

For all three analogies above, the magnetic circuit model is generated through an automated procedure. Initially, the transformer topology undergoes discretization, effectively partitioning the entire domain into a mesh structure, as exemplified in Figure 8. A three-phase and three-column transformer is modeled in this figure. Each area (domain) within the transformer can have a distinct mesh size, contributing to the establishment of an adaptive and tailored mesh, as in traditional finite-element modeling (FEM). Areas are determined based on two physical considerations. First, there is the inherent nature of the materials involved (ferromagnetic, dielectric, conductive, and so on). Second, the presence or absence of magnetomotive forces further influences the demarcation of these regions.

Figure 8. The transformer geometry mesh and examples of cells used to represent an element of the mesh for each of the three circuit-based methods.

Subsequently, to formulate the circuit-based equivalent model corresponding to the meshed system, elementary circuit “cells” are systematically associated with each individual mesh element; the code dynamically selects the appropriate type of cell, depending on the diverse characteristics exhibited by each element’s area. A novel library has been developed, consisting of an exhaustive array of cells to encompass all conceivable configurations that may arise in electromagnetic device modeling. The parameters of these cells are derived from the geometric and physical attributes of the individual mesh element. By effectively interconnecting these cells, a cohesive and comprehensive equivalent circuit model is established in the EMTP graphical user interface (see Figure 9).

Figure 9. A detailed magnetic circuit model of a three-phase transformer generated automatically in EMTP.

In essence, this automated approach allows for the systematic generation of models encapsulated as subcircuit devices. Upon double-clicking the device, users are granted access to a detailed representation of the model, faithfully aligned with the employed mesh configuration, as demonstrated in Figure 9.

This hierarchical structure closely mirrors the topological design of the studied device. The users can access internal magnetic circuit details in addition to surrounding complex power grid models.

The results stemming from the three circuit-based models, generated within EMTP, are juxtaposed against FEM outcomes acquired via ANSYS Electromagnetics. These evaluations pertain to a three-phase transformer (15 MVA, 25.63 to 4.53 kV). In Figure 10, the voltage current characteristics extracted from both open-circuit and short-circuit tests are presented. In the open-circuit test, the secondary is unloaded, and the primary winding is excited across varying voltage levels. The corresponding primary winding current is determined for each voltage level. Subsequently, in the short-circuit test, the secondary winding is short-circuited, and the primary winding is excited using different voltage levels. The secondary winding current is then deduced for each voltage condition.

Figure 10. The transformer V-I characteristic for (a) open-circuit and (b) short-circuit tests.

The observed maximum discrepancies with the FEM approach amount to 1.5% for the open-circuit test and 3.6% for the short-circuit. In this simulation, the circuit-based model utilizes a mesh including 2,848 elements, with the inclusion of over 2,200 nonlinear components. In contrast, the FEM mesh involves 3,292 elements. The simulations were conducted over a 32-ms interval, employing a time-step of 1 ${\mu}$s. The circuit-based model demonstrates a speed advantage, operating approximately three times faster than FEM, particularly under conditions characterized by significant nonlinearity and saturation effects. The good alignment between the outcomes of the circuit-based model and FEM results underscores the accuracy and efficiency of the proposed approach in reproducing electromagnetic behavior. One more important advantage is the capability to initialize systems dramatically faster than in the FEM approach.

Another validation scenario encompassing an electromagnetic transient is the study of inrush current arising from transformer energization. This investigation is presented in Figure 11. Once more, the maximum relative error observed between the circuit-based method and FEM remains below 3%.

Figure 11. The transformer energization transient.

This modeling approach not only facilitates the determination of global quantities but also possesses the capacity to capture localized behaviors. For instance, it can effectively provide the induction level at specific points within the transformer over time. Furthermore, the mesh data can be exploited with circuit solutions to generate field maps and animations.

To illustrate this capability, an internal fault scenario occurring in phase A of the primary winding is simulated. It represents a short circuit between internal winding turns appearing at 50 ms. The global quantities in Figure 12 enable studying the transient stemming from this fault. Complementing this, Figure 13 provides field maps that depict the distribution of induction levels at the instances of 52.5 and 60 ms. These maps offer a visual representation that readily reveals potential hot spots arising due to system malfunction.

Figure 12. The transient due to a short circuit between primary winding turns that occurs at 50 ms: the (a) primary current for the three phases and (b) active and reactive powers.

Figure 13. The distribution of induction levels, in teslas, at (a) 52.5 and (b) 60 ms for a fault case.

An additional advantage of these models lies in their ability to harness the EMTP environment for the purpose of conducting study scenarios necessitating a systemwide approach. Notably, GICs provoked by solar storms exemplify such scenarios. This phenomenon spreads over vast distances. Consequently, simulating such scenarios using conventional numerical methods becomes either prohibitively costly in computing times or practically unattainable.

The proposed circuit-based modeling approach offers a tangible solution. It not only facilitates the simulation of such complex scenarios within reasonable computing times but also provides access to internal constraints. This holds considerable value, not only for designers seeking to fortify transformer resilience but also for system operators tasked with implementing safeguards and protections to preemptively mitigate the impact of such disturbances.

Figure 14 illustrates a simplified configuration adopted to emulate a GIC scenario. The unloaded transformer is connected to the network. The GIC disturbance is injected in the system through the transformer neutral. The transformer is tested with two connectivity modes: Y-ground primary and Y-grounded secondary and Y-grounded primary and delta secondary. Figure 15 clearly demonstrates the significant differences in GIC impacts on the Y and delta configurations. This is explained by the fact that the delta connection does not allow the third harmonic. The shown reactive power oscillation exhibits slower stabilization compared to actual behavior, mainly due to the absence of iron losses in the current circuit-based model. These losses could be added as a parallel resistor to each winding for better accuracy.

Figure 14. The configuration for studying the GIC effect on a transformer.

Figure 15. The disturbance due to a GIC in the case of Yy0 and Yd11 connections: the (a) primary current and (b) reactive power.

The circuit-based model enables estimating iron losses through postprocessing. These estimations encompass both a comprehensive global value and a meticulous evaluation at the granularity of individual elements within the mesh. This process discerns the contributions of hysteresis losses and eddy current losses. This capability is demonstrated in Figure 16, where the aggregate losses over a single period, at a voltage level of 50% of the nominal supply, are estimated at 11.1 kW. Within this estimation, 23% of the losses are attributed to eddy current losses, while the remaining 77% arise from hysteresis losses. Moreover, the comprehensive distribution of these losses is pictured through the colormap.

Figure 16. The estimation of the distribution of iron losses, in watts, for a healthy case.

Conclusions

The presented unified environment is circuit based and allows simulating power systems, from load-flow solutions to fast electromagnetic transients. The load-flow solution plays a pivotal role in initializing simulated networks, enabling rapid attainment of steady-state conditions, even in the presence of IBRs, such as WPs.

Furthermore, the presented environment integrates the capability to model in detail magnetic circuits, with the noteworthy advantage of combining magnetic circuit meshing based on a device’s geometry.

For Further Reading

A. Haddadi and J. Mahseredjian, “Power system test cases for EMT-type simulation studies,” CIGRE Brochure, Paris, France, Rep. WG C4.503, Aug. 2018.

N. Rashidirad, J. Mahseredjian, I. Kocar, U. Karaagac, and O. Saad, “MANA-based load-flow solution for islanded AC microgrids,” IEEE Trans. Smart Grid, vol. 14, no. 2, pp. 889–898, Mar. 2023, doi: 10.1109/TSG.2022.3199762.

J. Mahseredjian, U. Karaagac, S. Dennetiere, and H. Saad, “Simulation of electromagnetic transients with EMTP-RV,” Numer. Anal. Power Syst. Transients Dyn. Inst. Eng. Technol., vol. 78, pp. 103–134, Jan. 2015.

U. Karaagac, J. Mahseredjian, S. Jensen, R. Gagnon, M. Fecteau, and I. Kocar, “Safe operation of DFIG-based wind parks in series-compensated systems,” IEEE Trans. Power Del., vol. 33, no. 2, pp. 709–718, Apr. 2018, doi: 10.1109/TPWRD.2017.2689792.

M. Ouafi, J. Mahseredjian, J. Peralta, H. Gras, S. Dennetière, and B. Bruned, “Parallelization of EMT simulations for integration of inverter-based resources,” Electric Power Syst. Res., vol. 223, Oct. 2023, Art. no. 109641, doi: 10.1016/j.epsr.2023.109641.

W. Nzale, J. Mahseredjian, X. Fu, I. Kocar, and C. Dufour, “Improving numerical accuracy in time-domain simulation for power electronics circuits,” IEEE Open Access J. Power Energy, vol. 8, pp. 157–165, Apr. 2021, doi: 10.1109/OAJPE.2021.3072369.

H. Xue, A. Ametani, and J. Mahseredjian, “Very fast transients in a 500 kV Gas-insulated substation,” IEEE Trans. Power Del., vol. 34, no. 2, pp. 627–637, Apr. 2019, doi: 10.1109/TPWRD.2018.2874331.

M. Naïdjate, N. Bracikowski, M. Hecquet, M. Fratila, M. M. Duro, and J. P. Ducreux, “An intelligent reluctance network model for the study of large power and distribution transformers,” in Proc. 6th Int. Adv. Res. Workshop Transformers (ARWtr), Cordoba, Spain, 2019, pp. 89–92, doi: 10.23919/ARWtr.2019.8930174.

S. R. Pordanjani, M. Naïdjate, N. Bracikowski, M. Fratila, J. Mahseredjian, and A. Rezaei-Zare, “Electromagnetic modeling of transformers in EMT-type software by a circuit-based method,” IEEE Trans. Power Del., vol. 37, no. 6, pp. 5402–5413, Dec. 2022, doi: 10.1109/TPWRD.2022.3177137.

S. R. Pordanjani, J. Mahseredjian, M. Naïdjate, N. Bracikowski, M. Fratila, and A. Rezaei-Zare, “Electromagnetic modeling of inductors in EMT-type software by three circuit-based methods,” Electric Power Syst. Res., vol. 211, Oct. 2022, Art. no. 108304, doi: 10.1016/j.epsr.2022.108304.

Biographies

Jean Mahseredjian (jeanm@polymtl.ca) is with Polytechnique Montréal, Montréal, QC H3T 1J4, Canada.

Mohammed Naidjate (mohammed.naidjate@polymtl.ca) is with Polytechnique Montréal, Montréal, QC H3T 1J4, Canada.

Mehdi Ouafi (mehdi.ouafi@rte-france.com) is with RTE, 69330 Lyon, France.

Juan Antonio Ocampo Wilches (juan-antonio.ocampo-wilches@polymtl.ca) is with Polytechnique Montréal, Montréal, QC H3T 1J4, Canada.

Digital Object Identifier 10.1109/MELE.2023.3320511

Basics of Electromagnetic Transients: Underlying mathematics

Xin Ma, Xiao-Ping Zhang

Electromagnetic transient (EMT) simulation plays an important role is power system operation, control, and planning. Considering that the digital computer is unable to give a continuous history of transient phenomena, it is essential to solve the time response at discrete intervals reasonably. In the 1960s, Dommel first started an EMT solution by a nodal matrix method. To reduce solution time for nonlinear and time-varying components, Dommel proposed a compensation method and developed equivalent circuits for those components in 1971. To improve the calculation speed, Dommel introduced a new technique for solving EMTs with the implicit trapezoidal rule of integration in 1972. In 1986, Dommel published his EMT Program (EMTP) theory book, where all of the EMT component modeling and overall solution methods are detailed.

EMT simulation tools can be implemented either offline or in real time. The well-known Electromagnetic Transients Program (EMTP®) development started in 1966 by Dommel. In addition to EMTP, power system computer aided design (PSCAD) was first conceptualized by Manitoba Hydro International Ltd. in 1988 and has become a competitive alternative to network transient analyzers. Using a parallel processing structure, the real-time digital simulator (RTDS) was first introduced by RTDS Technologies Inc. in 1992 to simulate EMT in real time with a time step of 50–100 ${\mu}$s.

In these commercial programs, a library of components is available and a circuit model can be constructed by dragging, dropping, and connecting the appropriate model block on the drawing canvas. Although drawing a schematic diagram is simple, a typical student registered in undergraduate or postgraduate courses is subjected to learn exclusively from fundamental operations. Because of huge pressure and limited time, it is difficult to be familiar with a large quantity of materials in the curriculum. Nearly all EMT simulation tools start from modeling basic components. A poorly designed EMT solution can result in matrix singularity and nonconvergence, which will make tasks become unsolvable for the next time step. With the very strict requirement of accuracy, it is essential for engineers to be aware of EMT issues and how to avoid them.

The purpose of this article is to provide a practical and concise EMT modeling guideline to avoid common mistakes. The contributions can be summarized as follows:

For numerical stability, three integration methods, including the trapezoidal rule, backward Euler, and forward Euler, are explained and compared graphically.
For accurate transient behavior, linear and nonlinear components are developed in a step-by-step procedure showing block diagram and mathematical models, respectively.
To reduce solution complexity, the node voltage method, compensation method, frequency shift method, and an offline/real-time solution are detailed in theoretical and graphical analysis.

The conclusions are supported for both offline simulators and real-time simulators. It is designed for rapid test of research and development ideas, especially for use in an advanced undergraduate and graduate level.

The article is organized as follows: For numerical stability, the “Integration Methods” section compares numerical integration methods. For discrete-time modeling, the “EMT Component Models” section provides a step-by-step guide to develop ordinary differential equations for different EMT components. To speed up EMT computation, the “Network Solution” section provides efficient methods. For verifying these models and methods, the “Case Studies” section gives a four-machine 11-bus case study. The last section summarizes conclusions.

Integration Methods

Basically for EMT studies, all system components are normally represented by difference equations. EMT simulators are expected to provide fast and accurate solutions for power grids with different components. Digital techniques are unable to give a continuous calculation of the transient phenomena, but rather a sequence of discrete intervals. A continuous period is required to be discretized from the current time step to the next time step until the maximum time is reached. To meet the simulation requirement, a simulation time step needs to be extremely small for discretization. For solving ordinary differential equations, common numerical methods include back Euler, forward Euler, and the trapezoidal rule. Figure 1 details a comparison among back Euler, forward Euler, and the trapezoidal rule:

Forward Euler [Figure 1(a)] is one of the most basic explicit integration methods. Only history value at previous time step ${t}{-}\Delta{t}$ is used to approximate future value for the whole time range. Since the future is obtained directly using history values ${f}{(}{t}{-}\Delta{t}{)}$ and ${u}{(}{t}{-}\Delta{t}{)}$, forward Euler is an explicit method. Without future value consideration, it is easy to grow an error of the first order.
Backward Euler [Figure 1(c)] is an implicit integration method. By comparison, this method only uses a future value at time step t to approximate neighborhood values. To obtain f(t), u(t) is required to be solved at time ${t}{-}\Delta{t}$. Therefore, backward Euler is an implicit method. According to the Taylor series, a local truncation error is introduced.
Trapezoidal rule [Figure 1(b)] is a second-order integration method, which uses both history values and future values to approximate the whole range. Since f(t) cannot be solved directly without u(t) time ${t}{-}\Delta{t}$, the trapezoidal rule is an implicit method. By providing damping to errors until coming back to the stable condition, this method is an absolutely stable integration method discussed by Dommel in 1969. In other words, history errors will not accumulate along with calculation time.

Figure 1. (a)–(c) Integration methods comparison.

By comparison, the trapezoidal rule is simple, numerically stable, and accurate enough to integrate the ordinary differential equations in EMT. Dommel discussed that discretization using forward Euler or backward Euler causes truncation errors and can lead to numerical instability, in the EMTP Theory Book in 1986. Therefore, the trapezoidal method is the most frequently used in EMT simulators, for instance EMTP. In the “Network Solution” section, how this is implemented in EMT simulation will be elaborated.

EMT Component Models

Based on arithmetic operations, EMT component models can be divided into linear and nonlinear components:

Linear components can be modeled with linear equations only. Most common linear components include passive circuits.
Nonlinear components contain nonlinear equations, such as sine and cosine functions.

By applying the trapezoidal rule, a linear component can be converted into a Norton equivalent circuit, which has an equivalent conductance in parallel with a known history current source. This is known as the companion circuit model.

Linear Components

Linear lumped components are resistances, inductances, and capacitances (RLC) branches.

To derive an accurate model for an arbitrary RLC branch between node k and node m, we can start from a single inductor of the circuit. Any current change through an inductor creates a changing flux and induces voltage. Figure 2 details the derived procedure and equivalent circuit for the inductor.

Figure 2. (a) and (b) EMT solution for RLC branch.

According to Faraday’s law of induction, inductance between node k and m starts from an ac steady-state solution. After the trapezoidal rule is applied to the integral of the inductor, the exact differential equation is replaced by the approximate discrete equation, and can be rewritten as an equivalent conductance in parallel with a history term, already known from the solution at the previous solution point at time ${t}{-}\Delta{t}$. Therefore, once currents and voltages are obtained at time ${t}{-}\Delta{t}$, the history terms can be updated for each individual inductor at time t.

To simplify modeling, it is essential to get general representation of an RLC branch. Similar to a single inductor, the EMT solution of an arbitrary RLC branch can be obtained if the trapezoidal rule is applied to the integral. Figure 2 details the mathematical models and equivalent circuits of an RLC branch. In spite of series or parallel connections, any RLC branch can be represented by a conductance G_RLC in parallel with a history term I_RLC. For the same RLC branch between node k and node m, parameters G_RLC, H_RLC, and J_RLC remain constant for the whole calculation.

Table 1 summarizes the detailed parameters and physical representation for common RLC branches. For example, The RC branch can be treated as a two-terminal network by combining unessential buses; therefore, the proposed mathematical models and parameters can effectively reduce computational burden.

Table 1. Comparison between different RLC branches.

The transmission lines modeling can be classified into a lumped parameter model and distributed parameter model. Based on the traveling wave model, the distributed parameter model implements a three-phase distributed line model where the RLC are uniformly distributed along the line. Contrary to the distributed parameters model, a three-phase ${\pi}$-section line lumped parameter model implements a balanced three-phase transmission line model with lumped parameters in a ${\pi}$-section circuit. For transient simulation, the traveling wave model is faster and more accurate than the lumped parameter model, as shown in Dommel’s EMTP Theory Book. In 1995, Zhang proposed the theoretic analysis and guidelines on the selection of a different number of equivalent ${\pi}$ sections of the transmission line model according to the frequency range of interest. A lumped model can be well represented by combining RLC branches. Table 2 details the comparison between the lumped parameter model and distributed parameter model.

Table 2. Comparison between lumped parameter line model and distributed parameter traveling wave model.

To get an accurate model of three-phase distributed transmission lines with loss, we start from single-phase lossless transmission lines. Figure 3 details that an EMT solution for transmission lines needs to go through four steps. At the first step, the surge impedance and travel time can be obtained based on traveling wave theory. At the second step, the distributed parameter model of lossless transmission line can be represented as an equivalent resistance and history term, similar to the RLC branch. The difference is that the history term is not instantly integrable because of travel time ${\tau}$. At the third step, a transmission line with loss is modeled with lumped R/4 in the terminals of transmission lines and R/2 in the middle of the transmission line, in which great accuracy were proved in experiments in Dommel’s article, published in 1969. At the fourth step, mutual three phases are converted to uncoupled modal domains to avoid an off-diagonal conductance matrix. By removing mutual couplings, each modal domain equation can be solved independently using surge impedance and travel time.

Figure 3. (a) and (b) EMT solution for a transmission line.

Nonlinear Components

Contrary to linear components, nonlinear components cannot be well represented by a simple Norton circuit. This is because additional nonlinear calculations, such as sine wave function and limiter, are involved in modeling nonlinear components. Typical nonlinear components include the synchronous machine and control system.

Synchronous machines play a significant role in controlling power system stability. A universal model is built up to represent synchronous machines, induction machines, and dc machines. Based on electric machine theory, synchronous machines consist of magnetic rotors and stators. It is always challenging to model electrical and mechanical parts simultaneously in a large system.

The key solution is to integrate differential equations for a whole synchronous machine. Because of the coupling relationship, the electrical part and mechanical part are unable to be solved simultaneously. To overcome this, we start from predicting the mechanical part because electrical changes are much faster. Figure 4 presents the detailed calculation of synchronous machines. The first step is to predict mechanical angular speed and angle by linearization so as to eliminate numerical errors. The second step is to convert physical quantities from abc reference of frame to dq0 reference of frame using Park’s transformation. Based on the trapezoidal rule, the third step is to update currents in armature and field windings according to a voltage–current–flux relationship. The fifth step is to update electrical torque. By back substitution, the sixth step is to solve the question of the exact mechanical part model, and update history values and then go back to the first step until the maximum solution time duration is reached.

Figure 4. EMT solution for synchronous machine.

The modeling of a control system, such as an excitation system, is the combination of basic arithmetic operations, limiters, and a proportional-integral-derivative (PID) controller. To represent a wide variety of excitation systems, IEEE has standardized 12 models for the excitation system currently in use, such as ac1A, dc1A, and ac4A, etc. To obtain an accurate model of an excitation system, we start with an ac4A excitation system to provide field voltage to control cd and terminal voltage. Figure 5 presents the mathematical model and block diagram of the ac4A system. At the first step, the output can be obtained by adding the referenced stator terminal voltage, stabilization voltage, and compensated voltage. At the second step, the limiter block ensures that the output will not exceed the upper bound and lower bound. At the third step, transient gain reduction is implemented by a PID control block with transient gain reduction lead and lag time constants, which can be discretized as a history term and current term. At the fourth step, the main regulator is implemented by a nonwindup limiter. With the nonwindup limits, the output comes off the limits as soon as input re-enters the range within limits. To overcome the coupling relationship between input and output, predicted value is linearized by history values. Applying the trapezoidal rule, nonwindup limiter can be partitioned as linearization, PID control block, and limiter in the discrete time domain.

Figure 5. (a) and (b) EMT solution for ac4A excitation system.

In addition to excitation systems, prime movers and governor systems are to provide mechanical torques for synchronous machines. Steam turbines, as one of the most common prime movers, can convert high temperature and high pressure into rotating energy.

Figure 6 shows the block diagram of a steam prime mover and a governor system. The first step is adder. The second step is multiplier. The third step is subtracter. At the third step, the load setup point is obtained by a subtracter. Based on the trapezoidal rule, a low-pass filter is implemented at the fourth and fifth steps. Therefore, steam governor and prime mover can be discretized step by step.

Figure 6. (a) and (b) EMT solution for steam governor and prime mover system.

Network Solution

Nodal Voltage Solution

Node voltage solution is able to straightforwardly provide EMT solutions for a network only consisting of linear components. This is because all linear components, such as the RLC branch, can be represented with an equivalent resistance in parallel with a history term. Equivalent resistance is involved in the admittance matrix, while the history term is considered as an injected current source. For a known node m, such as ground and voltage source, row m in the Y matrix can be set zero except set 1 for diagonal components. By Gaussian elimination, the remaining unknown node voltage can be obtained by ${U} = {\text{inv}}{(}{Y}{)}\,{\ast}\,{I}$. Therefore, node voltage and branch current at time t can be updated consistently using a node voltage solution.

Compensation Method

When nonlinear components are involved, the compensation method is required to deal with nonlinearity. The key of the compensation method lies in isolating the nonlinear part from the linear parts. In the compensation method, nonlinear components are treated as current injections, which are superimposed on the rest of the linear network. This treatment was initiated by Dommel and has been detailed in his 1986 EMTP Theory Book.

Figure 7 details applying a compensation method to the simultaneous solution of nonlinear and linear equations. At the first step, nonlinear components are excluded from the original network to get a Norton circuit of the rest of the open-circuit linear network. At the second step, open-circuit voltage can be obtained by setting the current injection of the nonlinear component as zero. At the third step, equivalent conductance can be calculated by replacing the current vector with a new vector whose components are all zero, except for +1 in row k and −1 in row m. At the third step, closed-circuit terminal voltage can be obtained after nonlinear components calculation. At the fourth step, the closed-circuit voltage superimposes the response into the current vector to complete the node voltage solution. Therefore, the compensation method can effectively get the simultaneous solution of both linear equations and nonlinear equations, simplifying computation burden.

Figure 7. (a) and (b) Application of compensation method.

Frequency-Shift Method

Power system transients are characterized by containing very high-frequency components. In 1999, Henschel proposed the frequency-shift method in his Ph.D. thesis. Basically, If the system returns to steady state, the signal spectrum becomes around the system frequency, which is called a bandpass signal. Frequency shift is to analyze a time function in the shifted time domain.

The basic procedures for frequency shift as shown by Henschel in his Ph.D. work in 1999 include: At the first step, a power system signal and dynamic phasor is represented with in-phase and quadrature. At the second step, an analytic signal is defined based on a Hilbert transform. At the third step, an analytical signal is derived using Fourier transform and system frequency is shifted by angular frequency ${-}{\omega}_{0}$.

This method allows a much larger simulation step to be applied to improve simulation efficiency without loss of accuracy.

Real-Time and Offline Simulation Method

To explore the detailed behaviors of power systems, it is necessary to construct circuits and observe simulation results. Depending on execution time and the simulation step, EMT solutions consist of offline simulation methods and real-time simulation methods. Table 3 illustrates the difference between real-time simulation and offline simulation.

In real-time simulation, execution time T_e is shorter than the selected time step. A real-time simulator needs to solve model equations for one time step within the same time in a real-world clock. To make a real-time simulation relevant, all necessary operations, including driving inputs and outputs to and from externally connected devices, require faster-than-real-time calculation speed. Typical real-time simulators include RTDS, as detailed in the online reference from RTDS Technologies in 2014 and field programming gate array (FPGA). In addition to accuracy, real-time simulators also need to be verified in calculation time and resource utilization.
In offline simulation, execution time T_e is normally greater than the selected time step, in which overruns occur. Offline software aims to solve mathematical equations as fast as possible. The actual solution at each step will normally take much longer than the simulation time step. Typical offline software includes PSCAD, EMTP, and MATLAB, etc.

Table 3. Comparison between real-time and offline simulation.

Case Studies

All of the case studies are implemented using a single Virtex6 ML605 FPGA board. This board has been equipped with 768 digital signal processors (DSPs), 416 RAMs, 301,440 registers, and 150,720 look-up tables (LUTs). A 100-MHz system clock is selected for global synchronization. All of the floating-point calculations are using intellectual property cores (IP Cores). The simulation step for EMT programs is 50 ${\mu}$s.

Verification

Real-time performance is validated in three aspects, including resource utilization, accuracy, and timing.

Resource utilization: Resource utilization is to estimate occupied resources, such as LUT, registers, and DSP, on a single ML605 FPGA board. Lower resource utilization means that less hardware is required for larger power system simulations. In order to generate a bitstream file, all of the resource utilization percentage should be less than 100% to avoid overused resources.
Accuracy: Accuracy is compared with a referenced offline simulation benchmark in MATLAB 2018b in Windows 64-bits system. The simulation step is 50 ${\mu}$s in both MATLAB and FPGA. The simulation results in FPGA are obtained in a binary format from the Chipscope tool, which allows inserting a logic analyzer, system analyzer, and virtual I/O low-profile software cores to view any internal signal.
Timing schedule to ensure real-time implementation: Timing schedule is to analyze the total time schedule so as to determine whether the FPGA-based EMT model is implementable for the chosen simulation time step $\Delta{t}$ (typically 50 ${\mu}$s or smaller). Implementable means execution time should be smaller than the chosen simulation time step, avoiding timing failures.

In summary, to ensure real-time implementation, resource utilization and timing constraints should be operated. In other words, the bitstream file, executable program is unable to be generated when overused resources or failed timing constraints occur.

Four-Machine 11-Bus Network

To verify the effectiveness of EMT mathematical models, a four-machine 11-bus system given in Figure 8(a) is implemented on a single ML605 FPGA board. A single-phase grounding fault and three-phase fault are applied at bus 7 from 1 s to 1.1 s, respectively. Figure 8 presents simulation results of electric torque, mechanical torque, and current and time schedule. Figure 8(c)–(e) presents simulation results with a single-phase grounding fault, and Figure 8(f)–(h) presents simulation results with a three-phase grounding fault.

Figure 8. (a)–(g) Four-machine 11-bus topology and simulation results.

It can be seen from Figure 8 that the relative error of electric torque, currents, and mechanical torque can be maintained less than 5%. With the trapezoidal rule, this relative error can be maintained in a limited range for both single-phase fault and three-phase fault. In particular, the d axis current and electric torque simulation results illustrate the proposed models can well represent dynamics of nonlinear calculations, for instance Park’s transform. It can be demonstrated that the proposed mathematical models, such as synchronous machines and RLC branches, can well represent the dynamic behaviors under both symmetric and asymmetric grounding faults.

Figure 8(b) shows the total calculation only costs 25.4 ${\mu}$s at a simulation time step 50 ${\mu}$s. Therefore, further calculation blocks are still allowed at available time slots. Seen from Table 4, it is observed that only 35% of registers and 87% of LUTs are utilized for the four-machine 11-bus system. Therefore, a larger system is still allowed on the same FPGA board, which shows high portability and flexibility.

Table 4. Resource utilization of 11-bus four-machine system.

To trigger higher frequency transients, the standard deviation is used to measure the variation of electric torque. In the referenced model, standard deviations of electric torque are 0.3318 and 0.0892, respectively, where the impact of the single-grounding fault is much more significant than that of the three-phase grounding fault. This is solid evidence as to why a single-phase grounding fault is more preferable for EMT simulation verification than a three-phase grounding fault.

Conclusion

This article has presented EMT models and solutions of power systems with components, such as RLC branches and synchronous machines. Theoretical analysis and graphical comparison can help beginners develop their fast familiarity with EMT models and solutions. Based on theoretical analysis and simulation results, the following conclusions can be given:

The trapezoidal rule is an absolute stable integration method to limit numerical errors for EMT simulation.
Linear component models, such as RLC branches and transmission lines, can be well represented by equivalent resistance and history term in a discrete time domain. Nonlinear components, for instance synchronous machines, can be discretized by following sequential steps.
Node voltage solution and a compensation method can effectively discretize linear and nonlinear components in the same system, reducing computational burden.
The dynamic transient behaviors can be well represented by developed mathematical models, such as synchronous machines and RLC branches, with relative errors less than 5%.
For fault analysis, a single-phase grounding fault can trigger higher frequency transients compared with three-phase grounding fault.

For Further Reading

H. W. Dommel, EMTP Theory Book. Portland, OR, USA: Bonneville Power Admin., Aug. 1986.

H. W. Dommel, “Nonlinear and time-varying elements in digital simulation of electromagnetic transients,” IEEE Trans. Power App. Syst., vol. PAS-90, no. 6, pp. 2561–2567, Nov. 1971, doi: 10.1109/TPAS.1971.292905.

H. W. Dommel and N. Sato, “Fast transient stability solutions,” IEEE Trans. Power App. Syst., vol. PAS-91, no. 4, pp. 1643–1650, Jul. 1972, doi: 10.1109/TPAS.1972.293341.

M. Gole, O. B. Nayak, T. S. Sidhu, and M. S. Sachdev, “A graphical electromagnetic simulation laboratory for power systems engineering programs,” IEEE Trans. Power App. Syst., vol. 11, no. 2, pp. 599–606, May 1996, doi: 10.1109/59.496082.

RTDS Technology Inc., Winnipeg, MB, Canada. Real-Time Simulation. (2014). [Online] . Available: http://www.rtds.com

X.-P. Zhang and H. Chen, “Analysis and selection of transmission line models used in power system transient simulations,” Int. J. Electr. Power Energy Syst., vol. 17, no. 4, pp. 239–246, Aug. 1995, doi: 10.1016/0142-0615(95)00039-S.

S. Henschel, “Analysis of electromagnetic and electromechanical power system transients with dynamic phasors,” Ph.D. dissertation, Dept. Elect. Eng., Univ. British Columbia, Vancouver, BC, Canada, 1999.

Biographies

Xin Ma (x.ma.1@bham.ac.uk) is with the Department of Electronic, Electrical, and Systems Engineering, University of Birmingham, B15 2TT Birmingham, U.K.

Xiao-Ping Zhang (x.p.zhang@bham.ac.uk) is with the University of Birmingham, B15 2TT Birmingham, U.K.; the Birmingham Energy Institute, Birmingham B25 8FJ, U.K.; and the Birmingham Center for Energy Storage, University of Birmingham, B15 2TT Birmingham, U.K.

Digital Object Identifier 10.1109/MELE.2023.3320485

Electromagnetic Transient Simulation: Moving to the Mainstream

Brian K. Johnson

Electromagnetic transient (EMT) simulation has moved from a tool used for a few specialist applications, such as insulation coordination, to becoming a common tool for interconnection studies for inverter-based resources (IBRs). Historically, many utilities had few, if any, engineers with a background to perform EMT studies. Today, the proliferation of inverter-based generation and storage resources makes EMT simulation a critical tool for protection and planning studies.

In early power systems, engineers performed transient analysis using hand calculations to solve differential equations. This approach was limited to small-scale problems. Applications included evaluating transient overvoltages due to switching transients such as capacitor switching or opening breakers. Some utilities applied transient network analyzers, which were scale model power systems that represented transmission or distribution lines as a series of modules with shunt capacitors and series resistor/inductor pairs. Adding more model segments improved the transient transmission line model to come closer to the real system. The transient network analyzers took a long time to set up but always ran in real time. Some utilities and equipment vendors used transient network analyzers for hardware-in-the-loop (HIL) testing of protective relays and of controls for high-voltage direct current transmission and static volt-amperes reactive compensators. Transient network analyzers continued to see use into the 1990s but were complicated to set up and had limits to scalability.

Hermann W. Dommel and W. Scott Meyer developed an EMT program for the computer simulation of EMTs at the Bonneville Power Administration in the late 1960s. The underlying algorithm is used in many of the EMT simulation tools in common use today. Initially, engineers simulated switching transients and performed insulation studies. Due to limitations in computer capabilities and difficulties in building and verifying models, EMT simulations were limited to cases with a limited number of nodes. A small number of engineers used transient simulation, and many utilities often had few, if any, engineers familiar with EMT simulation.

Protection engineers were among the earlier adopters of EMT simulation tools, especially for classes of protection studies not amenable to using phasor-based programs. One significant example is protection studies for series compensated lines, where transient low-frequency oscillations in the voltage and current pass through the digital filtering and cause relay overreach. Relay engineers started out using offline EMT studies with simulation results played to relays using amplifiers or directly to relays as low-level analog outputs to test protection settings.

The growing importance of power electronic controls in power systems led to another area where EMT simulation had a growing importance over the years. Simple power electronic switch models and control system models provided the ability to model the Pacific High Voltage dc (HVdc) intertie in the Bonneville Power Administration EMT program. The Manitoba HVdc Research Centre developed an EMT program based on Dommel’s paper to simulate the Nelson River HVdc system. An early EMT program application was to identify the impact of a sensor’s time constant on the control response.

EMT simulation soon moved on to the modeling of other power converter applications, including power quality studies. Researchers started modeling voltage source converters for photovoltaic inverters, flexible ac transmission systems (FACTS) devices, and many other applications. One of the challenges with modeling power converter systems is access to accurate models of the converter controls. Unlike generator exciters and governors, standard models are not available. Accurate modeling of protective relays faces a similar challenge. It is possible to build generic models, but obtaining sufficiently accurate models of a specific vendor’s control is an ongoing challenge.

Real-time EMT simulation grew out of the need to replicate the transient network analyzer capabilities to test the closed response of commercial protection and control devices for power system conditions. These simulators run Dommel’s EMT simulation algorithm in real time with high-speed input and output to external devices, enabling HIL simulation. A key challenge was dividing the system into small parts solvable in parallel without compromising the simulation accuracy. The solution takes advantage of the propagation delay for transients moving over transmission lines to divide the problem into pieces that can be solved in parallel. The approach requires that the propagation delay is longer than the signal delay for exchanging information between processors. There are other options to divide the problem as well. Protection engineers were early adopters of real-time EMT simulations, with commercial relays connected to simulators to test relay responses in an HIL environment. HIL testing is commonly used in the design and testing of controls for HVdc and FACTS systems. Asset owners are increasingly purchasing replica copies of their converter controls to use in simulator studies, both to test changes to controls and to test firmware updates from vendors.

The capabilities of real-time simulators have evolved as processing power has improved. User requirements for real-time simulation have increased as well. The increased use of automation schemes and the advent of digital substations introduce an arena where HIL simulation supports design and testing of automation schemes, protection schemes, and the underlying communication systems. Consulting firms, utilities, and equipment vendors are designing and testing their automation schemes, including remedial action schemes, through HIL simulation. These applications are moving beyond power engineering to support industrial control systems in other industries. HIL simulation is also a growing platform for some areas of cybersecurity research.

The capabilities of conventional offline EMT simulation tools and available models have evolved to meet the demands of users, taking advantage of improvements in computing platforms. Improved graphical interfaces for EMT programs eliminate a significant source of errors in earlier programs. EMT programs can now simulate much larger systems. However, building and validating a model of a truly large system is challenging. Several simulation programs support black box models from equipment vendors, such as controls for HVdc converters and IBRs, allowing engineers to perform studies with more accurate models of the device response.

Modern grids with large-scale installations of IBRs, such as wind turbines, photovoltaic generation, and energy storage systems, are increasing the need for applying EMT simulation in planning studies and in protection studies. The impacts of fast controls on the fault current responses of voltage source converter-based IBRs, HVdc converters, and FACTS devices is difficult to capture in phasor-based fault programs. In addition, aspects of the dynamic responses of these resources are difficult to capture in other standard tools. EMT simulation is a valuable tool for a growing number of types of system studies.

There are challenges with the increased use of EMT simulation. Acquiring adequate model data remains a challenge, and verifying models is an increasingly significant challenge, especially for systems with high levels of IBR generation from different vendors. Utilities build these models through a combination of in-house expertise, external consultants, or equipment suppliers. There is a need to help engineers learn to perform all the steps of the EMT simulation process, requiring knowledge of the applications as well as the specific simulation tools they will use. While user interfaces have made EMT simulation much easier to set up, there is still a need for specialized training.

EMT simulation is now an indispensable tool for power engineers to ensure resilient operation of modern power systems.

Biography

Brian K. Johnson (bjohnson@uidaho.edu) is with the University of Idaho, Moscow, ID 83844 USA.

Digital Object Identifier 10.1109/MELE.2023.3320482

Computer Simulation Tools for Electrification

Caisheng Wang, Lingling Fan

Analysis, design, implementation, and control of an engineering system always start with modeling. On the one hand, modeling is always a process of simplification. For example, we can simply model a battery by an ideal constant voltage source and a light bulb by a resistor. Then, we can analyze the circuit and calculate circuit quantities (current, voltage, power, etc.). On the other hand, depending on the scope of study, we need to have enough details or a certain level of fidelity in the models so that meaningful studies can be carried out. For example, copper polygons in low-frequency (e.g., ≤10 kHz) printed circuit boards (PCBs): they can be simply treated as copper conductors using lump parameters, while we should use transmission line theory when handling polygons in high-frequency (e.g., over gigahertz) PCBs.

From a microgrid with a handful of distributed energy resources, such as solar panels and energy storage, and a group of local loads to a large, interconnected power system with thousands of components and even millions of subcomponents and systems, we need a variety of tools to help us to analyze, design, implement, and operate such systems. Electromagnetic transient (EMT) analysis is one of the important tools. In a world increasingly driven by sustainable energy development and green transportation, the role of EMT analysis cannot be overstated. It is the hidden force behind the scenes, the technology that powers the clean energy revolution, enables efficient transportation electrification, paves the way for carbon neutrality, and propels the electric vehicle industry into the future.

This issue is on EMT analysis fundamentals, methods, tools, data preparation, and applications. Eight remarkable feature articles, along with two columns authored by leading experts from around the globe, shed light on this important field, encompassing a wide spectrum of topics. From system and component modeling, multidomain and multiphysics modeling, to advanced simulation and emulation methods, such as hardware-in-the-loop (HIL), and applications in the design and implementation of zero-emission energy systems, this issue takes us on a journey through the heart of EMT analysis.

Brian K. Johnson from the University of Idaho contributes to the “Technology Leaders” column with “Electromagnetic Transient Simulation: Moving to the Mainstream.” He provides his insightful perspective on the evolution of EMT analysis from a highly specialized tool for a small group of people at the beginning to a common tool today for interconnection studies for inverter-based resources (IBRs) and design and implementation of transportation electrification systems. The power systems around the world have gone through significant changes during the same period of time. When the first author first saw a 500-kV transmission line and tower in the 1980s, he was really fascinated as an elementary school student. Now, there are transmission lines of over 1,000 kV, and the global electricity generation from wind and solar has gone up from almost zero in 1982 to over 3 trillion kWh in 2022. EMT tools have a unique contribution to this fantastic development.

With the promise that computer simulation tools can offer, Peter Wung penned the “Viewpoint” column titled “To Test or Not To Test?” to remind the reader that computer simulation tools are based on physics we know, and in the territories that physics are still to be learned, the trustworthiness of computer simulation results remains a question. Fortunately, for electromagnetism, power electronic circuits, and controls, we have the full confidence of physics and therefore computer simulation tools have been vastly deployed.

The first article is contributed by Xin Ma and Xiao-Ping Zhang from the University of Birmingham, providing a brief history of EMT development, and reviewing the mathematical fundamentals (such as trapezoidal rule in numerical analysis) of EMT models and solutions of power systems with components, such as RLC branches and synchronous machines. The second article is written by Dragan Jovcic, an expert on high-voltage dc from the University of Aberdeen. This article provides a great introduction to modeling power electronic converters, particularly via the average value approach based on a rotating frame. Power converters have become indispensable units in modern power systems. Power converters are the core equipment for most renewable and clean alternative energy systems, such as wind, solar, and fuel cell generation systems, as well as for transportation electrification. It has been estimated that 80% of utilized electricity will flow through power electronics devices by 2030. Analytical converter modeling complements detailed digital simulation and can expedite system studies.

The third feature article is the outcome of a collaboration between Kai Strunz from Technische Universität Berlin, Ying Chen of Tsinghua University, and Yue Xia from China Agricultural University. The article introduces the shift frequency method for frequency-adaptive simulation of wide-scale transients in power systems. The so-called frequency-adaptive simulation of transients method can handle power system transients in a vast range, from minutes down to microseconds or even faster for the EMTs of traveling waves in an efficient and accurate manner.

Electrical machines are the core equipment in power systems and are used in almost every aspect of our daily life. Modeling and simulation of electric machines play a critical role in EMT studies. The fourth article, coauthored by industry expert Erfan Mostajeran (from Electric Power Engineers LLC.) and academic researchers Navid Amiri, Seyyedmilad Ebrahimi, and Juri Jatskevich (from The University of British Columbia), provides a comprehensive review of this important topic of electrical machine modeling. Typical electrical machine models introduced with appropriate details in the article represent essential electromagnetic and mechanical transient phenomena. Electrical machines with more than three (e.g., five, six, nine, 12, etc.) electrical phases are also introduced in the article for applications with some specific requirements, such as increased reliability and fault-tolerant operation. Discussions on further improvement of the electrical machine model’s fidelity, numerical accuracy, and computational efficiency are also given in the article. Electrical machines are the core equipment in power systems and are used in almost every aspect of our daily life. Modeling and simulation of electric machines play a critical role in EMT studies.

The fifth article from the University of Alberta by Chengzhang Lyu and Venkata Dinavahi covers another important aspect of EMT studies: multidomain modeling of EMTs. Zero-emission marine vessels powered by fuel cells with battery storage are presented in the article as an interesting example of EMT modeling across electrical, chemical, hydraulic, and thermal domains. HIL emulation has often been considered for expediting the simulation process and for obtaining more accurate results. An FPGA-based, real-time HIL emulation system, called adaptive compute acceleration platform, is introduced for emulation studies of a zero-emission marine vessel.

The sixth article, contributed by Dionysios Aliprantis and Steven Pekarek from Purdue University, continues the discussions on computational electromagnetics for electric machine analysis, modeling, and design. Computational electromagnetics has played and will continue to play an even more important role in transportation electrification thanks to the recent advances in computing. After introducing the finite element method and the method of moments in the article, the authors provide their perspective on multiphysics modeling and 4D (space and time) simulation studies of electric machines.

The seventh article is another great piece of international collaboration work between Jean Mahseredjian, Mohammed Naidjate, and Mehdi Ouafi from Polytechnique Montreal and Juan Antonio Ocampo Wilches from Réseau de Transport d’Électricité in France. This article introduces a comprehensive and unified environment that enables the exploration of power systems across various modeling levels, encompassing load-flow analysis and extending to intricate time-domain simulations featuring IBRs. The authors first start with the initialization of IBR circuits and parallel computations for practical system problems. Then, they move on to demonstrate that the time-domain accuracy level can be further increased through detailed and integrated modeling of magnetic circuits in the proposed unified environment.

Last but not least, the eighth feature article in this issue is on data preparation for EMT simulation studies. The authors, Taku Noda and Tomo Tadokoro from Grid Innovation Research Laboratory, Central Research Institute of Electric Power Industry (CRIEPI), Japan and Takashi Dozaki from National Renewable Research Laboratory, discuss this important topic, and present application examples of the developed automatic data generation system in Japan. Simulation data preparation, although not directly focused on EMT modeling, is an indispensable segment and has become the most time-consuming part of EMT studies today.

There are also two timely newsfeeds in this issue. One newsfeed is on the past IEEE Transportation Electrification Conference and Expo (ITEC) conference held in Detroit, MI, USA, contributed by Matt Woongkul Lee from Michigan State University. A panel consisting of academia and industry experts was successfully organized by Dr. Lee in the conference to discuss high-performance rare-earth magnet-free electric machine design. The panel presentations and discussions underscore the importance of collaboration, innovation, and optimization to unlock the full potential of rare-earth magnet-free designs, which are critical to the sustainable development of transportation electrification. The other newsfeed is provided by Taku Noda from CRIEPI, Japan about the 2023 International Conference on Power Systems Transients (IPST). This year’s IPST hosted about 200 EMT experts around the world, seeking solutions for up-to-date and long-standing issues in the power system EMTs. We believe the readers will find the detailed descriptions of the conference interesting. More interestingly, beautiful pictures are also shared in this newsfeed.

Appendix: Related Articles

B. K. Johnson, “Electromagnetic transient simulation: Moving to the mainstream,” IEEE Electrific. Mag., vol. 11, no. 4, pp. 6–7, Dec. 2023, doi: 10.1109/MELE.2023.3320482.

P. Wung, “To test or not to test?” IEEE Electrific. Mag., vol. 11, no. 4, p. 96, Dec. 2023, doi: 10.1109/MELE.2023.3320529.

X. Ma and X.-P. Zhang, “Basics of electromagnetic transients: Underlying mathematics,” IEEE Electrific. Mag., vol. 11, no. 4, pp. 8–19, Dec. 2023, doi: 10.1109/MELE.2023.3320485.

D. Jovcic, “Rotating frame, average value converter modeling: Basic principles in building analytical models,” IEEE Electrific. Mag., vol. 11, no. 4, pp. 20–28, Dec. 2023, doi: 10.1109/MELE.2023.3320486.

K. Strunz, Y. Chen, and Y. Xia, “Bridging scales with the shift frequency: Frequency-adaptive simulation of multiscale transients in power systems,” IEEE Electrific. Mag., vol. 11, no. 4, pp. 29–37, Dec. 2023, doi: 10.1109/MELE.2023.3320487.

E. Mostajeran, N. Amiri, S. Ebrahimi, and J. Jatskevich, “Electrical machines in electromagnetic transient simulations: Focusing on efficient and accurate models,” IEEE Electrific. Mag., vol. 11, no. 4, pp. 38–53, Dec. 2023, doi: 10.1109/MELE.2023.3320508.

C. Lyu and V. Dinavahi, “Zero-emission marine vessels: Multi domain modeling and real-time hardware-in-the-loop emulation on adaptive compute acceleration platform,” IEEE Electrific. Mag., vol. 11, no. 4, pp. 54–63, Dec. 2023, doi: 10.1109/MELE.2023.3320509.

D. Aliprantis and S. Pekarek, “The unsung hero of the electric vehicle revolution: The role of computational electromagnetics in electric machine design and analysis,” IEEE Electrific. Mag., vol. 11, no. 4, pp. 64–68, Dec. 2023, doi: 10.1109/MELE.2023.3320510.

J. Mahseredjian, M. Naidjate, M. Ouafi, and J. A. O. Wilches, “Electromagnetic transients simulation program: A unified simulation environment for power system engineers,” IEEE Electrific. Mag., vol. 11, no. 4, pp. 69–78, Dec. 2023, doi: 10.1109/MELE.2023.3320511.

T. Noda, T. Tadokoro, and T. Dozaki, “Automatic generation of power system simulation data cases from utility databases: Introducing a new technology,” IEEE Electrific. Mag., vol. 11, no. 4, pp. 79–85, Dec. 2023, doi: 10.1109/MELE.2023.3320521.

W. Lee, “Rare earth magnet-free electric machine design: Unlocking sustainable electrification,” IEEE Electrific. Mag., vol. 11, no. 4, pp. 88–89, Dec. 2023, doi: 10.1109/MELE.2023.3320522.

T. Noda, “International conference on power systems transients 2023,” IEEE Electrific. Mag., vol. 11, no. 4, pp. 90–92, Dec. 2023, doi: 10.1109/MELE.2023.3320556.

Digital Object Identifier 10.1109/MELE.2023.3320479

Editorial board

Editor-in-Chief

Lingling Fan

University of South Florida, USA

fll@ieee.org

Associate Editors

Fabio D’Agostino

University of Genova, Italy

David Wenzhong Gao

University of Denver, USA

Wei Hua

Southeast University, China

Jahangir Khan

BC Hydro, Canada

Hiva Nasiri

K&A Engineering Consulting, USA

Ehsan Nasr

Microsoft Corporation, USA

Ioannis (John) Prousalidis

National Technical University of Athens, Greece

Farrokh Rahimi

Open Access Technology International, Inc., USA

Deepak Ramasubramanian

Electric Power Research Institute, USA

Doris Sáez Hueichapan

University of Chile, Chile

Marcelo Godoy Simões

University of Vaasa, Finland

Caisheng Wang

Wayne State University, USA

Ying Xue

South China University of Technology, China

Bo Yang

Hitachi America, USA

Xiao-Ping Zhang

University of Birmingham, UK

IEEE PERIODICALS MAGAZINES DEPARTMENT

445 Hoes Lane,

Piscataway, NJ 08854 USA

Christie Inman

Journals Production Manager

Patrick Kempf

Senior Manager, Journals Production

Gail A. Schnitzer
Associate Art Director

Theresa L. Smith

Production Coordinator

Felicia Spagnoli

Advertising Production Manager

Peter M. Tuohy

Production Director

Kevin Lisankie

Editorial Services Director

Dawn M. Melley

Senior Director, Publishing Operations

advertising sales

Aviva Rothman

Naylor Association Solutions

Phone: +1 352 333 3435

arothman@naylor.com

Mission Statement: IEEE Electrification Magazine is dedicated to disseminating information on all matters related to microgrids onboard electric vehicles, ships, trains, planes, and off-grid applications, including electrification of remote communities and decarbonization strategies. Microgrids refer to electric networks in a car, a ship, a plane, or an electric train, which have a limited number of sources and multiple types of loads as well as off-grid applications that include small scale multi-carrier energy systems supplying electricity and heat in areas away from high voltage power networks. Feature articles focus on advanced concepts, technologies, and practices associated with all aspects of electrification in the transportation and off-grid sectors from a technical perspective in synergy with nontechnical areas such as business, environmental, social and policy topics.

IEEE Electrification Magazine (ISSN 2325-5897) (IEMECM) is published quarterly by the Institute of Electrical and Electronics Engineers, Inc. Headquarters: 3 Park Avenue, 17th Floor, New York, NY 10016-5997 USA. Responsibility for the contents rests upon the authors and not upon the IEEE, the Society, or its members. IEEE Operations Center (for orders, subscriptions, address changes): 445 Hoes Lane, Piscataway, NJ 08854 USA. Telephone: +1 732 981 0060, +1 800 678 4333. Individual copies: IEEE members US$20.00 (first copy only), nonmembers US$140.00 per copy. Subscription Rates: Society members included with membership dues. Subscription rates available upon request. Copyright and reprint permissions: Abstracting is permitted with credit to the source. Libraries are permitted to photocopy beyond the limits of U.S. Copyright law for the private use of patrons 1) those post-1977 articles that carry a code at the bottom of the first page, provided the per-copy fee indicated in the code is paid through the Copyright Clearance Center, 222 Rosewood Drive, Danvers, MA 01923 USA; 2) pre-1978 articles without fee. For other copying, reprint, or republication permission, write Copyrights and Permissions Department, IEEE Operations Center, 445 Hoes Lane, Piscataway, NJ 08854 USA. Copyright © 2023 by the Institute of Electrical and Electronics Engineers, Inc. All rights reserved. Postmaster: Send address changes to IEEE Electrification Magazine, IEEE Operations Center, 445 Hoes Lane, Piscataway, NJ 08854 USA. Canadian GST #125634188 Printed in U.S.A.

Digital Object Identifier 10.1109/MELE.2023.3319893

Title

Bridging Scales With the Shift Frequency: Frequency-adaptive simulation of multiscale transients in power systems

Evolving Methods and Needs

The Shift Frequency as a Simulation Parameter

2D Setting of the Shift Frequency and Time-Step Size

Realizing Frequency Shifting

Network Modeling

Lumped Parameters

Distributed Parameters

Application

Wind Energy Conversion System

Wide-Area Energy System

Conclusions

Acknowledgment

For Further Reading

Biographies

IEEE Electrification Magazine / DECEMBER 2023

Title

Electrical Machines in Electromagnetic Transient Simulations: Focusing on efficient and accurate models

Introduction

Electrical Machine Characteristics and Modeling

Coupled-Circuit-Phase-Domain Models

qd0 Models

EMT Simulators

SVB EMT Simulators

EMTP Simulators

Voltage-Behind-Reactance Models

Example of VBR Models in SVB EMT

Induction Machine qd0 Versus VBR Models (Example Case 1)

Modeling Magnetic Saturation

Air-Gap Flux Harmonics in Induction Machines (Example Case 2)

Example of SM in EMTP

Concluding Remarks and Future Work

For Further Reading

Biographies

IEEE Electrification Magazine / DECEMBER 2023

Title

Rotating Frame, Average Value Converter Modeling: Basic principles in building analytical models

Detailed Versus Analytical Modeling of Converters

Average Value Converter Modeling

Basic Principles

VSC Average Modeling

CSC Modeling

Converter Modeling in Rotating Coordinate Frame

Model Linearization and Format Conversion

Verification of Accuracy of the Average Model

Stability Analysis

Relevance and Challenges

Illustration of Stability Analysis Using a Detailed Model

System Stability Analysis in Paramedic Domain Using an Average Model

Conclusion

For Further Reading

Biography

IEEE Electrification Magazine / DECEMBER 2023

Title

Zero-Emission Marine Vessels: Multidomain Modeling and Real-Time Hardware-in-the-Loop Emulation on Adaptive Compute Acceleration Platform: Zero-emission marine vessels: modeling and real-time emulation

Introduction

Marine Vessel Power Systems

Electric Propulsion System

Energy Storage System

Multidomain Modeling and Real-Time Emulation

Proton Exchange Membrane Modeling

Electrical Domain

Fluid Domain

Thermal Domain

Real-Time HIL Emulation

Real-Time HIL Emulation Results and Discussion

Conclusions

For Further Reading

Biographies

IEEE Electrification Magazine / DECEMBER 2023

Title

The Unsung Hero of the Electric Vehicle Revolution: The role of computational electromagnetics in electric machine design and analysis

The Finite Element Method

The Method of Moments

Conclusions

For Further Reading

Biographies

IEEE Electrification Magazine / DECEMBER 2023

Title