Model Implementations

• 1: Reduced Order Generator Implementation
• 2: Switch and Load Implementation
• 3: Source Implementation
• 4: Injection and Compensation Implementation
• 5: SSN Component Implementation
• 6: Signal Component Implementation
• 7: DP Ph1 Averaged VSI Implementation
• 8: DP Ph3 Averaged VSI Implementation
• 9: EMT Ph3 Averaged VSI Implementation
• 10: EMT Ph3 Grid-Forming VSI Implementation

1 - Reduced Order Generator Implementation

The equations are derived under reduced order machine models. This page covers only their arrangement in code.

Class hierarchy

Base::ReducedOrderSynchronGenerator<VarType> holds everything independent of domain and of order: the per unit base values, the operational parameters, the mechanical states, the controller attachments and the discretisation coefficients. It is templated on Real for EMT and Complex for DP and SP, which is why the axis frame quantities appear twice, as mVdq0/mIdq0 in the real specialisation and mVdq/mIdq in the complex one.

Each domain then provides a ReducedOrderSynchronGeneratorVBR layer holding the frame transform, and each order a concrete class. The order is recorded in mSGOrder, which selects which coefficients are computed.

Network interface

setModelAsNortonSource chooses between the two interface forms. The default is the Norton equivalent, in which the machine contributes only to the right hand side vector and requests no virtual nodes. The Thevenin form requests two virtual nodes instead. Both represent the same model; the Norton form is cheaper because it leaves the system matrix untouched between steps [Wang2010].

Coefficients

calculateAuxiliarConstants computes the discretisation coefficients once, since they depend only on the parameters and the step size. The member names map to the symbols on the theory page as follows.

The naming looks wrong at first and is not. Zd_t is built from $L_q - L_q’$ and Zq_t from $L_d - L_d’$, because each is named for the axis whose coefficient it feeds rather than for the parameters it is assembled from. That follows the physics: the d-axis internal voltage arises from q-axis rotor flux and decays with $T_{q0}’$, so mAd_t correctly combines $L_q - L_q’$ with $T_{q0}’$ and multiplies the q-axis current.

Read a coefficient’s use rather than its assignment line before concluding an axis is swapped.

Step sequence

mnaCompPreStep runs before the network solve and does three things in order. It advances the controllers, saving mEf_prev and mMechTorque_prev first because the trapezoidal history terms need the previous values. It calls stepInPerUnit, which updates the frame transforms from mThetaMech, recomputes the axis frame state from the terminal quantities, and evaluates the history voltage into mEh_vbr. It then stamps the result into the right hand side vector.

Each concrete order implements only specificInitialization and stepInPerUnit. Everything else is inherited.

Initialization

Initialization runs from the powerflow solution, not from user supplied states. The base class computes the load angle as the phase of $V + j L_q I$, projects the terminal voltage and current onto the axis frame, and derives the field voltage from the no-load relation. Only then does specificInitialization set the order specific states, which is why a concrete class can assume mVdq and mIdq are already populated.

Attached controllers are initialized afterwards from the machine’s own initial values, so an exciter or governor never needs its own operating point.

Controllers

Excitation, governor, turbine and power system stabilizer attach through the base class and are optional, guarded by mHasExciter, mHasGovernorAndTurbine, mHasTurbineGovernor and mHasPSS. The stabilizer output feeds the exciter within the same step, and the governor output feeds the turbine, so the order of the calls in mnaCompPreStep is load bearing.

Source code

• Base class header, implementation
• Concrete orders under dpsim-models/src/{SP,DP,EMT}/ named SynchronGenerator<N>OrderVBR
• Availability per domain is listed under model availability

References

• [Wang2010] IEEE Xplore document 5411963. Cited in the machine model pages as the basis for interfacing a machine to a nodal solver through a current source that leaves the system matrix unchanged.

2 - Switch and Load Implementation

The models are derived under switches and loads. This page covers only their arrangement in code. Availability per domain is in model availability.

Switches

Switch implements Base::Ph1::Switch and stamps one admittance chosen by mIsClosed, using MNAStampUtils::stampAdmittance so the grounded-terminal cases are handled centrally. SeriesSwitch folds a series resistance into the same branch.

varResSwitch additionally implements MNAVariableCompInterface, which is what allows it to change the system matrix during a run. Its hasParameterChanged is called each step and drives the transition:

• Opening multiplies the resistance by mDeltaResOpen each step until it passes the target open value, then clamps to it and reports the transition finished.
• Closing uses mDeltaResClosed, which is 0, so the first step takes the resistance to zero, the clamp catches it and sets the closed value. Closing is therefore immediate by construction, not by a separate code path.

Its initializeFromNodesAndTerminals carries a comment saying it is not used.

Loads

RXLoad is a CompositePowerComp. In initializeFromNodesAndTerminals it converts the powers to element values and builds sub-components:

• a resistor, only if the active power is non-zero
• an inductor if the reactance is positive, a capacitor if negative, and nothing if the reactive power is zero

PQLoadCS wraps a current source and sets its reference in updateSetPoint from conj(S / mNomVoltage). The nominal voltage, not the terminal voltage, is deliberate; the line using the terminal voltage is present but commented out. Changing it would make the component nonlinear and require an iterative solve.

Shunt takes a conductance and a susceptance directly and additionally carries per-unit attributes, since it is the form the powerflow solver consumes.

Source

• Switches: {SP,DP,EMT}_Ph{1,3}_Switch, DP_Ph3_SeriesSwitch, EMT_Ph3_SeriesSwitch, {DP,SP}_Ph1_varResSwitch under dpsim-models/src/
• Loads: DP_Ph1_RXLoad, EMT_Ph3_RXLoad, DP_Ph1_RXLoadSwitch, DP_Ph1_PQLoadCS, SP_Ph1_Load, {SP,DP}_Ph1_Shunt, EMT_Ph3_Shunt

3 - Source Implementation

The models are derived under sources. This page covers only the code.

Ideal sources

VoltageSource requests one virtual node, which carries the source current as the extra unknown, and stamps the constraint rows that fix the terminal voltage difference. CurrentSource requests none and contributes only to the right hand side.

ControlledVoltageSource and ControlledCurrentSource are the same components with their reference supplied as an attribute rather than a parameter, so another component or an interface can drive them. The reference is read during the pre-step, which is why it is the previous step’s value.

VoltageSourceNorton

Stamps directly rather than through a virtual node. mnaCompApplySystemMatrixStamp adds mConductance to both diagonal entries and subtracts it from the two off-diagonal entries, guarded by terminalNotGrounded, and mnaCompApplyRightSideVectorStamp sets the equivalent current mIntfVoltage / mResistance with opposite signs at the two terminals.

mConductance is computed in setParameters as 1 / resistance, so calling setParameters is mandatory before the run and a zero resistance is a division by zero rather than an ideal source.

VoltageSourceRamp

A composite wrapping a VoltageSource whose reference it rewrites each step in updateState(time). Three regimes: before mSwitchTime the reference is unchanged; during mRampTime the added voltage is interpolated linearly while the added frequency is blended by a raised sine 0.5 + 0.5 * sin(pi * t / T - pi/2); afterwards both are fully applied.

The two are blended differently on purpose. A linear frequency interpolation applied as a phase offset would step the phase at both ends of the ramp; the raised sine has zero derivative at both ends, so the frequency contribution enters and leaves smoothly. The consequence is that the instantaneous frequency during the ramp is not the linear interpolation between the two values, and reading mAddSrcFreq as “the frequency at the midpoint” is wrong.

Note also that the added frequency term is applied as mAddSrcFreq * time, using absolute simulation time rather than time since the switch, so the phase contribution depends on when in the run the ramp occurs.

ProfileVoltageSource

Holds a std::filesystem::path, a sample vector and an index, and reads the file in readFromFile at construction. It implements DAEInterface in addition to the MNA hooks.

The samples are stepped by index rather than interpolated against simulation time, so the profile’s sample rate and the simulation step must match for the waveform to have the intended duration. It is bound in Python and constructing it with a file that is not a readable sample list raises rather than crashing, which is covered by a test.

Source

Under dpsim-models/src/{SP,DP,EMT}/. Availability per domain is in model availability.

4 - Injection and Compensation Implementation

The models are derived under network injection and compensation. This page covers only the code.

NetworkInjection

A CompositePowerComp wrapping a single VoltageSource sub-component. It owns no equations of its own; it exists so that the external network is a named component rather than a bare source, and so that the driving waveform can be swapped without changing the network description.

setParameters is overloaded by the kind of generator wanted behind it: a constant phasor for a fixed source, a start frequency with a rate of change for a ramp, and an initial phasor with a modulation frequency for a modulated one. Which overload is called determines which SignalGenerator the sub-source is given; see signal component implementation.

Because the source is ideal, adding an impedance to represent a finite short circuit level is the caller’s job. Nothing in the component does it.

SVC

Not composite. It computes a susceptance each step and realises it by reconfiguring an internal reactive element, so it implements the variable-component interface and forces a refactorisation whenever the value changes.

updateSusceptance performs both lags with the trapezoidal rule, using precomputed constants Fac1 = dt / (2 Tr), Fac2 = dt Kr / (2 Tr) and Fac3 = dt / (2 Tm). The measurement lag is applied first, then the error is formed in per unit against mNomVolt, then the susceptance follows from the previous value and the present and previous error.

The result is clamped to mBMax and mBMin before use, and the internal element is only rebuilt when the value actually changed. The sign of the clamped susceptance selects which element is formed: positive gives an inductance 1 / (omega * B * mBN), negative a capacitance B * mBN / (-omega). mBN is the base susceptance, so B is per unit.

SolidStateTransformer

A CompositePowerComp that represents each side as a current source rather than as a coupled winding pair. setParameters(nomV1, nomV2, Pref, Q1ref, Q2ref) takes the two nominal voltages and three power set points; the active power is common to both sides, while the reactive powers are set per side.

Values are held in per unit internally, so the nominal voltages are the base rather than a turns ratio. There is no magnetising branch, no leakage impedance and no angle dependence, which is the representation the concept page describes and not an omission.

Source

• DP_Ph1_NetworkInjection, and the SP and EMT variants alongside it
• DP_Ph1_SVC
• SP_Ph1_SolidStateTransformer

Availability per domain is in model availability.

5 - SSN Component Implementation

The method is derived under state-space nodal components and state-space nodal. This page covers only the code.

Base hierarchy

SSNComp holds the continuous matrices, the discrete pair, the equivalent admittance mW, the history vector mYHist and the state attribute x. Two branches specialise it by which quantity is the input:

• VTypeSSNComp takes voltage in and gives current out, so it stamps as an admittance
• ITypeSSNComp is the dual

Terminal-count layers sit on top: TwoTerminalVTypeSSNComp, TwoTerminalITypeSSNComp and FourTerminalVTypeSSNComp handle the mapping from terminal quantities to the model input, and the Variable layers add re-forming of the model between steps. Every EMT::Ph3 base sets PhaseType::ABC in its constructor, which the concrete components rely on.

What a component provides

A fixed-model component only calls SSNComp::setParameters(A, B, C, D) with its chosen state, input and output. EMT::Ph3::SSN::Inductor is the whole pattern:

Matrix aMatrix = Matrix::Zero(3, 3);   // x = i_abc
Matrix bMatrix = inductance.inverse(); // u = v_abc
Matrix cMatrix = Matrix::Identity(3, 3);
Matrix dMatrix = Matrix::Zero(3, 3);
SSNComp::setParameters(aMatrix, bMatrix, cMatrix, dMatrix);

The base does the rest: recomputeDiscreteModel calls Math::calculateStateSpaceTrapezoidalMatrices and sets mW = mC * mdB + mD, calculateHistoryVector returns mC * (mdA * x + mdB * u), and the post step updates the state from the old and new input.

A varying component additionally overrides updateStateSpaceModel (a no-op for linear components) and, for the variable bases, updateComponentParameters to report whether the model changed. Only when it reports a change is the system matrix refactorised.

Domain differences

The formulation differs by domain, and so does the code path. The theory is under SSN across domains.

A component supplies the same real (A, B, C, D) in either domain. EMT::SSNComp discretises them directly. DP::SSNComp does not: buildAugmentedA(omega) assembles the real-augmented 2n x 2n matrix with A on both diagonal blocks and +wI / -wI off-diagonal, buildAugmentedB places B on both diagonal blocks, and the result goes through the same Math::calculateStateSpaceTrapezoidalMatrices helper as EMT. The discrete blocks are then folded back into complex form as topLeft + j * bottomLeft, which is the inverse of the [[P, -Q], [Q, P]] representation. mW and the history vector are complex as a result.

recomputeDiscreteModel therefore takes omega in DP and takes no argument in EMT. A component that hardcodes a frequency here rather than using the value handed to mnaCompInitialize is wrong at any other system frequency.

Frame metadata

getLocalAbcStateBlocks returns nothing by default and should be overridden only for states that genuinely form physical abc triples. It is consumed by tooling that reasons about the state vector in the phase frame, and declaring a block that is not one produces wrong groupings rather than an error.

Initialization

calculateSteadyStateStateFromInput evaluates (jωI − A)⁻¹ B u, which requires the continuous model to be set first. Components with real control states cannot use the default initializeFromNodesAndTerminals on the mixed base; see DP Ph1 averaged VSI implementation for that case and for the requirement that the state matrix be handed over already carrier shifted.

The components

Fixed models: SSN_Full_Serial_RLC, SSN_Capacitor, SSN_Inductor, SSNTypeV2T, SSNTypeI2T. Varying models: SSN_Variable_Serial_RLC, PiecewiseLinearInductor, and the inverter models under power electronics. The Generic two- and four-terminal classes take the matrices from the caller instead of forming them, so they are the route to an SSN component without writing C++. Availability per domain is in model availability.

6 - Signal Component Implementation

The models are derived under converter control and signal processing blocks. This page covers only the code.

Base and scheduling

Signal blocks derive from SimSignalComp and take no part in the nodal solve. They contribute tasks through getTasks() rather than through the MNA hooks, and the scheduler orders them from the attribute dependencies those tasks declare. A block that reads an attribute without declaring it may still produce the right answer, by luck of ordering, and then change behaviour when an unrelated component is added; see adding tasks to a component.

Most blocks follow a two-task shape: a PreStep that copies the current values into the previous ones, and a Step that computes the new state and output. The split exists so that a value consumed by another block within the same step is unambiguous about which timestep it belongs to.

The mInputPrev / mInputCurr pattern

Blocks that integrate with the trapezoidal rule need both the present and the previous input, so they carry mInputPrev, mInputCurr, mStatePrev, mStateCurr and the matching output pair. The PreStep task performs the shift. Integrator::signalStep is the whole pattern:

**mStateCurr = **mStatePrev + mTimeStep / 2.0 * **mInputCurr
                            + mTimeStep / 2.0 * **mInputPrev;
**mOutputCurr = **mStateCurr;

VCO::signalStep deliberately does not, using state + dt * input, because it accumulates an angle rather than integrating a control signal.

Every one of these blocks needs setSimulationParameters(timestep) before the run, since the step size appears directly in the update. Blocks that expose setInitialValues must also have it called, or they start from zero rather than from the operating point.

State-space blocks

PLL is written as an explicit state-space block rather than as arithmetic, setting

mA << 0, mKi, 0, 0;
mB << 1, mKp, 0, 1;
mC << 1, 0, 0, 1;
mD << 0, 0, 0, 0;

The first input is the nominal frequency and is held constant, which is how the feed-forward term enters. Writing it this way means the block can be discretised by the same helpers as anything else rather than by hand.

FIRFilter

FIRFilter keeps a circular buffer and a write index, and step sums mFilter[i] * mSignal[...] over the filter length before advancing the index. It contributes a single Step task. Filter coefficients are supplied by the caller; nothing validates their length against the buffer or checks that they sum to a sensible gain.

Generators

SignalGenerator is the abstract base; SineWaveGenerator, DCGenerator, CosineFMGenerator and FrequencyRampGenerator are the concrete ones, and all expose their value through a sigOut attribute that a source component references.

Source

Under dpsim-models/src/Signal/. Availability is in model availability; these blocks are domain independent and appear there as a list rather than a matrix.

7 - DP Ph1 Averaged VSI Implementation

The equations are derived under DP Ph1 averaged voltage source inverter. This page covers only their arrangement in code.

Class and base

DP::Ph1::AvVoltSourceInverterStateSpace is final and derives from DP::Ph1::MixedVTypeVariableSSNComp. The mixed base is what makes the model possible in this domain: eight of the twelve states are real baseband control states and only the last four are the real and imaginary parts of the two carrier-band envelopes, so the component cannot use the plain complex SSN base.

State layout

The state order is fixed by a private StateIndex enum, which the linearization indexes directly.

The base does not impose this ordering. It is told only how many real and how many complex states there are, and sizes the packed real vector as realStateCount + 2 * complexStateCount. The three-phase model orders its states the other way round, envelopes first and controls after, and is equally valid. What the base does require is that the derived class hand it a state matrix that is already carrier shifted: the steady-state solve assumes it, and a model that supplies an unshifted matrix initializes to the wrong operating point rather than failing.

The default initializeFromNodesAndTerminals throws unless realStateCount is zero, so any model with real control states, which includes this one, must override it.

Tracking Psi rather than the raw PLL angle keeps the tracked quantity bounded. The raw angle grows without limit, which costs relinearization accuracy as a run gets longer.

Parameters

initializeFromNodesAndTerminals derives the initial state from the connected node voltage, so the operating point comes from the powerflow rather than from user supplied states.

Source and examples

• Source: header, implementation
• C++ example
• Python notebook
• The state-space nodal interfacing this model relies on is described under state-space nodal

8 - DP Ph3 Averaged VSI Implementation

The equations are derived under DP Ph3 averaged voltage source inverter. This page covers only their arrangement in code.

Class and base

DP::Ph3::AvVoltSourceInverterStateSpace is final and derives from DP::Ph1::MixedVTypeVariableSSNComp, the same mixed base as the single-phase model. Per-phase complex quantities are carried as std::array<Complex, 3>.

State layout

Twenty states by default, ordered envelopes first and controls afterwards, or twenty-two with the optional negative-sequence loop enabled.

This is the reverse of the single-phase ordering, which places controls first. The base does not care: it is given only the counts of real and complex states and sizes the packed real vector as realStateCount + 2 * complexStateCount. What it does require is a state matrix that is already carrier shifted, since the steady-state solve assumes it.

The last two states are the difference from the single-phase model beyond the per-phase filter. Three independent phase envelopes admit a negative-sequence component that a single positive-sequence envelope cannot represent, so the controller carries its own negative-sequence integrator pair.

Enabling the negative-sequence loop

The constructor takes an enableNegSeqControl flag, false by default. The two references $i_{nd,\mathrm{ref}}$ and $i_{nq,\mathrm{ref}}$ are the last two arguments of setParameters and default to zero, which makes the loop a suppressor rather than an injector. The measured $i_{rc,nd}$ and $i_{rc,nq}$ are exposed as the irc_n_d and irc_n_q attributes, and stay at zero while the loop is disabled.

The two integrators are appended after the control block rather than inserted next to the other control states, so enabling the flag leaves every envelope and positive-sequence control index unchanged. Code indexing into the state vector therefore does not need to know about the flag.

The theory behind the loop is derived under DP Ph3 averaged VSI.

Source and examples

• Source: header, implementation
• C++ example
• Python notebook
• Single-phase counterpart: DP Ph1 averaged VSI implementation

9 - EMT Ph3 Averaged VSI Implementation

The equations are derived under EMT Ph3 averaged voltage source inverter. This page covers only their arrangement in code.

Class and base

EMT::Ph3::AvVoltSourceInverterStateSpace is final and derives from EMT::Ph3::TwoTerminalVTypeVariableSSNComp. Unlike the dynamic phasor ports of this model, every state here is real, so it uses the plain variable state-space nodal base rather than the mixed one. The base sets PhaseType::ABC in its constructor, which this component relies on.

State layout

Fourteen real states, controls first and filter states afterwards.

The first state is the raw PLL angle. The dynamic phasor ports track the deviation from the nominal carrier phase instead, because there the angle is compared against a carrier and an unboundedly growing value costs relinearization accuracy. In EMT there is no carrier to drift against, so the raw angle is used directly.

Six real filter states here correspond to two complex envelopes in the single-phase dynamic phasor model and six in the three-phase one. That correspondence is the practical statement of what the envelope transform buys.

Source and examples

• Source: header, implementation
• C++ example
• Python notebook
• Dynamic phasor counterparts: DP Ph1, DP Ph3

10 - EMT Ph3 Grid-Forming VSI Implementation

The equations and the linearization mathematics are derived under EMT Ph3 grid-forming voltage source inverter. This page covers only their arrangement in code.

Class and base

EMT::Ph3::SSN_GFM is final and derives from EMT::Ph3::TwoTerminalVTypeVariableSSNComp. All seventeen states are real.

State layout

Omega and Theta being states rather than inputs is what makes this grid forming: the converter carries its own frequency and angle instead of tracking a measured one through a PLL.

Numerical linearization

The Jacobians are not written out by hand. calculateNumericalJacobians forms all four by central differences of the nonlinear state and output functions, so a change to the control equations needs no matching change to any matrix code.

The perturbation for column $j$ is absoluteStep + relativeStep * max(1, |x_j|), defaulting to 1e-8 and 1e-6 and adjustable at runtime. The max(1, ...) floor means the step is effectively absolute for small states and relative for large ones, which keeps the difference well conditioned across states whose magnitudes differ by orders.

Because the model is time varying, the state-space form and its stamp are recomputed every step rather than cached. That is the cost of this approach and the reason it is used only where the control is genuinely nonlinear.

Source and examples

• Source: header, implementation
• C++ example, IEEE 9-bus mixed machine and inverter
• Python notebook
• Grid-following counterpart: EMT Ph3 averaged VSI implementation