Voltage-Behind-Reactance Non-Excited Synchronous Machine

[Beat_Arnet]

This document describes the structure of the Non-Excited Synchronous Machine model implemented in PLECS. The machine model is applicable to three-phase permanent magnet (PM) and synchronous-reluctance machines. The machine is parameterized using non-linear flux-maps that are typically exported from a finite element analysis (FEA) tool or from lab measurements. The Voltage-Behind-Reactance (VBR) machine formulation allows the machine to interface with arbitrary electrical networks, making it suitable for simulating in converters with switching dead-time or in various failure modes.

Notation

Space vectors are written as real 2x1 column vectors, and rotation operators as 2x2 matrices. The key operators are the 90° rotation matrix \mathbf{J} (the matrix counterpart of multiplication by j in complex space-vector notation) and the frame rotation matrix \mathbf{T}(\theta). Both are defined in Appendix A . Bold symbols denote 2x2 matrices; non-bold symbols denote scalars or 2-vectors.

A space vector \vec{x} is the Clarke projection of a three-phase quantity onto two orthogonal axes, written as \vec{x} = [x_\alpha,\ x_\beta]^T in the fixed frame or \vec{x} = [x_d,\ x_q]^T in the synchronous (dq) frame. The variables used in this document are:

• \vec{v}_s, \vec{i}_s — stator voltage and current vectors in the synchronous (dq) frame.
• \vec{\psi}_{m,s} — stator magnetizing flux vector in the synchronous frame; output of the flux-map LUTs (see Flux-map lookup). The magnetizing qualifier is essential: the flux-map LUTs tabulate the airgap-coupled flux only, not the leakage component. See the convention paragraph below.
• \vec{e}_s — stator back-EMF vector in the synchronous frame; computed from \vec{\psi}_{m,s} and \vec{i}_s (see Back-EMF construction).
• R_s — stator phase resistance (scalar).
• L_{\sigma s} — stator leakage inductance scalar, assumed equal across all three phases. Treated as a separate model parameter, not included in the flux maps. See Leakage handling for the implementation, including the singularity argument that requires L_{\sigma s} > 0.
• \mathbf{L}_{mi} — stator incremental inductance, a 2x2 matrix that depends on the operating point; defined and computed in Flux-map lookup.
• \omega = d\theta/dt — electrical rotor speed; \theta is the rotor electrical angle.
• p — pole pair count (appears in the torque equation).

Convention: flux maps are magnetizing-only; leakage is separate

The flux-map LUTs \vec{\psi}_{m,s}(i_d, i_q) tabulate the magnetizing flux linkage — flux that crosses the airgap and couples to the rotor. They do not include leakage flux that links a phase winding without crossing the airgap (slot leakage, end-winding leakage). The leakage inductance L_{\sigma s} is treated as a separate scalar parameter, assumed equal on all three phases, and added to the model along the two paths described in Leakage handling.

Architecture

The machine is modeled as a three-phase stator winding with a magnetically anisotropic, possibly magnetized rotor. Any PM contribution is absorbed into the flux-map lookups. A pure synchronous-reluctance machine is the special case with no PM contribution. The circuit is configured with a floating neutral: the star point of the three-phase winding is not hardwired to ground. No assumption is made that the three phase currents sum to zero. The model includes a zero-sequence inductance component (see Leakage handling in Transformation to the abc frame), so the neutral point can be connected to an external electrical circuit and zero-sequence currents will flow correctly.

Appendix A derives the synchronous-frame stator voltage equation:

\vec{v}_s = R_s\,\vec{i}_s + (L_{\sigma s}\mathbf{I} + \mathbf{L}_{mi})\,\frac{d\vec{i}_s}{dt} + \omega\mathbf{J}\,(L_{\sigma s}\vec{i}_s + \vec{\psi}_{m,s})

The first two terms are the resistive drop and the inductive drop on the time derivative of the current. The third is the rotational EMF from the total stator flux rotating at speed \omega, with one contribution from the leakage flux L_{\sigma s}\vec{i}s and one from the magnetizing flux \vec{\psi}_{m,s}.

This form is correct, but the inductive content is split across two terms. The \frac{d\vec{i}_s}{dt} piece and the \omega\mathbf{J}\,L_{\sigma s}\vec{i}_s piece sit separately, and the same is true for the \mathbf{L}_{mi} contributions. The VBR rearrangement regroups them so that all the inductive content is acted on by a single \left(\frac{d}{dt} + \omega\mathbf{J}\right) operator.

The derivation is as follows:

\begin{align} \vec{v}_s &= R_s\,\vec{i}_s + (L_{\sigma s}\mathbf{I} + \mathbf{L}_{mi})\,\frac{d\vec{i}_s}{dt} + (L_{\sigma s}\mathbf{I} + \mathbf{L}_{mi})\,\omega\mathbf{J}\,\vec{i}_s + \omega\mathbf{J}\,(L_{\sigma s}\vec{i}_s + \vec{\psi}_{m,s}) - (L_{\sigma s}\mathbf{I} + \mathbf{L}_{mi})\,\omega\mathbf{J}\,\vec{i}_s\,\\ &= R_s\,\vec{i}_s + (L_{\sigma s}\mathbf{I} + \mathbf{L}_{mi})\!\left(\frac{d\vec{i}_s}{dt} + \omega\mathbf{J}\,\vec{i}_s\right) + \omega\mathbf{J}\,(L_{\sigma s}\vec{i}_s + \vec{\psi}_{m,s}) - (L_{\sigma s}\mathbf{I} + \mathbf{L}_{mi})\,\omega\mathbf{J}\,\vec{i}_s \\ &= R_s\,\vec{i}_s + (L_{\sigma s}\mathbf{I} + \mathbf{L}_{mi})\!\left(\frac{d\vec{i}_s}{dt} + \omega\mathbf{J}\,\vec{i}_s\right) + \omega\mathbf{J}\,\vec{\psi}_{m,s} - \mathbf{L}_{mi}\,\omega\mathbf{J}\,\vec{i}_s \\ \end{align}

In the first step add and subtract (L_{\sigma s}\mathbf{I} + \mathbf{L}_{mi})\,\omega\mathbf{J}\,\vec{i}_s. The second and third terms then share (L_{\sigma s}\mathbf{I} + \mathbf{L}_{mi}) as a common left factor, which is grouped in the second step.

The remaining two \omega\mathbf{J} terms expand with the L_{\sigma s}\omega\mathbf{J}\vec{i}_s contributions cancelling (L_{\sigma s} is a scalar), leaving

\omega\mathbf{J}\,\vec{\psi}_{m,s} - \mathbf{L}_{mi}\,\omega\mathbf{J}\,\vec{i}_s =: \vec{e}_s

Putting it all together,

\boxed{\;\vec{v}_s = R_s\,\vec{i}_s + (L_{\sigma s}\mathbf{I} + \mathbf{L}_{mi})\!\left(\frac{d\vec{i}_s}{dt} + \omega\mathbf{J}\,\vec{i}_s\right) + \vec{e}_s \qquad \vec{e}_s = \omega\mathbf{J}\,\vec{\psi}_{m,s} - \mathbf{L}_{mi}\,\omega\mathbf{J}\,\vec{i}_s\;}

The stator voltage is a sum of three terms: a resistive drop R_s\,\vec{i}_s, an inductive drop (L_{\sigma s}\mathbf{I} + \mathbf{L}_{mi})\!\left(\frac{d\vec{i}_s}{dt} + \omega\mathbf{J}\,\vec{i}_s\right), and a back-EMF \vec{e}_s. The back-EMF carries the rotational EMF from the magnetizing flux (\omega\mathbf{J}\,\vec{\psi}_{m,s}) plus the inductive cross-coupling correction -\mathbf{L}_{mi}\,\omega\mathbf{J}\,\vec{i}_s that was pushed out of the inductance bracket.

This equation is in the synchronous (dq) frame. To be implemented as an electrical circuit at the abc terminals, each of the three voltage terms must be computed and then transformed back into the fixed (abc) frame. The remainder of this document describes how each term makes that journey.

Computations in the synchronous (dq) frame

The dq frame rotates with the rotor and sees a static magnetic structure. For this model, flux is a function of current alone, with no dependence on rotor angle. This is where the nonlinear magnetic behavior of the machine is modeled. Based on the measured phase currents and the rotor position and speed inputs, three quantities are computed in the dq frame: the incremental inductance matrix \mathbf{L}_{mi}, the back-EMF vector \vec{e}_s, and the electromagnetic torque T_e.

Flux-map lookup

The stator phase currents are transformed into the synchronous frame by the inbound Park transform , yielding the dq current vector \vec{i}_s = [i_{sd},\ i_{sq}]^T. This is the input to the DQ Inverse Saturation Lookup block, shown in the image above.

The block first checks whether the operating point lies within the calibrated region of the tables. Four threshold comparators test i_{sd} and i_{sq} against the minimum and maximum breakpoints of the LUT axes and raise a warning flag if either limit is exceeded.

Within the valid region, two 2-D lookup tables return the magnetizing flux vector

\vec{\psi}_{m,s}(i_{sd}, i_{sq}) = \begin{bmatrix} \psi_{md}(i_{sd}, i_{sq}) \\ \psi_{mq}(i_{sd}, i_{sq}) \end{bmatrix}

output as the Psim signal. In parallel, four 2-D lookup tables return the entries of the incremental inductance matrix

\mathbf{L}_{mi}(i_{sd}, i_{sq}) = \begin{bmatrix} \dfrac{\partial \psi_{md}}{\partial i_{sd}} & \dfrac{\partial \psi_{md}}{\partial i_{sq}} \\[6pt] \dfrac{\partial \psi_{mq}}{\partial i_{sd}} & \dfrac{\partial \psi_{mq}}{\partial i_{sq}} \end{bmatrix} = \begin{bmatrix}L_{mi,dd}&L_{mi,dq}\\L_{mi,qd}&L_{mi,qq}\end{bmatrix}

output as the Lmi signal.

For a saturating machine, \vec{\psi}_{m,s}(\vec{i}_s) is a nonlinear function of current. The flux value at an operating point and the local slope \mathbf{L}_{mi} = \partial\vec{\psi}_{m,s}/\partial\vec{i}_s are independent: neither can be recovered from the other. The diagonal entries capture self-saturation and saliency (L_{mi,dd} \neq L_{mi,qq}). The off-diagonal entries capture cross-saturation, where d-axis current saturates the q-axis iron and vice versa. Magnetic reciprocity gives the matrix its symmetry (\partial\psi_{md}/\partial i_{sq} = \partial\psi_{mq}/\partial i_{sd}), so only three of the four entries are independent, with all four shown in the equation for clarity. At a saturating operating point, both diagonal entries vary with current, the off-diagonal entries become nonzero, and neither can be read from the flux value alone. This is why the model requires five separate lookup tables.

The initialization script lets the user supply just the flux maps and compute the incremental inductances automatically by numerical differentiation using MATLAB’s gradient function. The inductance tables are an optional override for users who prefer to supply maps from FEA or analytical expressions. See Appendix B for additional details on the initialization script.

For a linear, non-saturating machine, the flux map is linear in current and the incremental inductance reduces to a constant matrix. For a salient but non-saturating machine without cross-coupling (L_{mi,dq}=0, L_{mi,dd} \neq L_{mi,qq}), the off-diagonal entries vanish and the matrix is diagonal:

\mathbf{L}_{mi} = \begin{bmatrix} L_{mi,dd} & 0 \\ 0 & L_{mi,qq} \end{bmatrix}

For a non-salient, non-saturating machine (L_{mi,dd} = L_{mi,qq} = L_m), the matrix reduces to a scalar multiple of the identity:

\mathbf{L}_{mi} = L_m\,\mathbf{I}

These are the cases for which the LUT-based formulation reduces to the textbook constant-parameter model.

Back-EMF construction

The Es block computes the synchronous-frame back-EMF vector

\vec{e}_s = \omega\mathbf{J}\,\vec{\psi}_{m,s} - \mathbf{L}_{mi}\,\omega\mathbf{J}\,\vec{i}_s

as defined in Architecture. The block is shown below:

The top branch multiplies Psim by \omega\mathbf{J}, giving the rotational EMF from the magnetizing flux. The bottom branch multiplies is by \omega\mathbf{J}, then applies \mathbf{L}_{mi} via the M·x block. The order matters: \mathbf{L}_{mi} and \omega\mathbf{J} do not commute as 2x2 matrices. The two results are subtracted to yield \vec{e}_s, available as the Es_dq signal. The block also produces Es_uvw by applying the inverse Park transform, making the result available directly in the stationary abc frame.

Torque

The Fcn block evaluates the electromagnetic torque from the magnetizing flux and the dq current:

T_e = \frac{3}{2}\,p\,(\psi_{md}\,i_{sq} - \psi_{mq}\,i_{sd})

using Psim for the flux and is for the current.

Transformation to the abc frame

The three voltage terms from the Architecture section each require a different outbound transformation to reach the abc terminals.

The resistive drop R_s\vec{i}_s needs no transformation. Resistance is a scalar and therefore frame-independent. The abc-frame resistive drop is simply R_s\vec{i}_{s,abc}, computed directly from the phase current measurements.

The back-EMF \vec{e}_s is transformed by a standard inverse Park rotation. The Es block (Back-EMF construction) already outputs Es_uvw by applying the inverse Park transform to Es_dq, so no separate block is needed.

The inductive drop is the most involved, requiring two successive frame changes to reach the abc terminals: first from dq to αβ, then from αβ to abc. The starting point is the VBR dq-frame inductive term (L_{\sigma s}\mathbf{I}+\mathbf{L}_{mi})\!\left(\frac{d\vec{i}_s}{dt}+\omega\mathbf{J}\vec{i}_s\right) from the Architecture section; the goal is the 3×3 abc inductance matrix \mathbf{L}_{s,abc} that drives the PLECS Variable Inductor. This is the job of the Lss block, highlighted below.

The first step converts from dq to αβ. Differentiating \vec{i}_s = \mathbf{T}^{-1}(\theta)\,\vec{i}_{s,\alpha\beta} by the product rule and applying the rotation-matrix derivative \frac{d\mathbf{T}^{-1}}{dt} = -\omega\mathbf{J}\mathbf{T}^{-1}(\theta) (derived in Appendix A) gives the identity

\frac{d\vec{i}_s}{dt}+\omega\mathbf{J}\vec{i}_s = \mathbf{T}^{-1}(\theta)\,\frac{d\vec{i}_{s,\alpha\beta}}{dt}

Substituting into the VBR inductive term gives (L_{\sigma s}\mathbf{I}+\mathbf{L}_{mi})\,\mathbf{T}^{-1}(\theta)\,\frac{d\vec{i}_{s,\alpha\beta}}{dt}, which is still a dq-frame voltage. Left-multiplying by \mathbf{T}(\theta) converts it to the αβ frame:

\mathbf{L}_{s,\alpha\beta}\,\frac{d\vec{i}_{s,\alpha\beta}}{dt}, \qquad \mathbf{L}_{s,\alpha\beta} := \mathbf{T}(\theta)(L_{\sigma s}\mathbf{I} + \mathbf{L}_{mi})\mathbf{T}^{-1}(\theta)

The second step expands the 2×2 matrix \mathbf{L}_{s,\alpha\beta} to the full 3×3 abc inductance matrix \mathbf{L}_{s,abc}.

PLECS does not support matrix multiplication as a native block operation. The similarity transform and the αβ-to-abc expansion are implemented as networks of scalar f(u) function blocks, each computing one algebraic combination of the matrix entries. The two sub-blocks Mdq→Mab and Mab→Muvw carry out this calculation in two stages.

Mdq → Mab: dq to αβ

The Mdq → Mab sub-block computes the similarity transform \mathbf{T}(\theta)(L_{\sigma s}\mathbf{I} + \mathbf{L}_{mi})\mathbf{T}^{-1}(\theta) entry by entry. The block takes a 3-element input vector. The incremental inductance matrix \mathbf{L}_{mi} is symmetric by magnetic reciprocity, so its four entries reduce to three independent values: L_{mi,dd}, L_{mi,qq}, and L_{mi,dq}. Before entering the block, the stator leakage is added via the constant vector [Lls 0 Lls], which adds L_{\sigma s} to L_{mi,dd} and L_{mi,qq} while leaving L_{mi,dq} unchanged. This defines (L_{\sigma s}\mathbf{I} + \mathbf{L}_{mi}) = \lbrack L_{dd}, L_{dq} , L_{qd}, L_{qq} \rbrack with L_{dq} = L_{qd}.

The three independent entries of the symmetric \mathbf{L}_{s,\alpha\beta} matrix are:

\begin{align*} L_{\alpha\alpha} &= L_{dd}\cos^2\!\theta + L_{qq}\sin^2\!\theta - 2L_{dq}\sin\theta\cos\theta \\ L_{\beta\beta} &= L_{qq}\cos^2\!\theta + L_{dd}\sin^2\!\theta + 2L_{dq}\sin\theta\cos\theta \\ L_{\alpha\beta} = L_{\beta\alpha} &= L_{dq}(\cos^2\!\theta - \sin^2\!\theta) + (L_{dd}-L_{qq})\sin\theta\cos\theta \end{align*}

The block uses the double-angle trigonometric entities to simplify the expressions into functions of \cos(2\theta) and \sin(2\theta), avoiding runtime squaring and products. Applying the double-angle identities collapses the squared terms:

\begin{align*} L_{\alpha\alpha} &= \bar{L} + \Delta L\cos 2\theta - L_{dq}\sin 2\theta \\ L_{\beta\beta} &= \bar{L} - \Delta L\cos 2\theta + L_{dq}\sin 2\theta \\ L_{\alpha\beta} = L_{\beta\alpha} &= \Delta L\sin 2\theta + L_{dq}\cos 2\theta \end{align*}

where \bar{L} = \tfrac{L_{dd}+L_{qq}}{2} is the isotropic average and \Delta L = \tfrac{L_{dd}-L_{qq}}{2} is the half-saliency. The inductance matrix decomposes cleanly into an angle-independent part (\bar{L}) and an angle-dependent part driven by saliency (\Delta L) and cross-saturation (L_{dq}), both modulated at twice the electrical angle.

The output of the block a 3-element vector Mab = (L_{\alpha\alpha},\, L_{\alpha\beta},\, L_{\beta\beta}).

Mab → Muvw: αβ to abc

The Mab → Muvw sub-block expands the αβ matrix into the 9-element 3x3 abc inductance matrix. Since matrix operations are not available in PLECS, six scalar f(u) blocks evaluate the required algebraic combinations of the αβ entries. A separate scalar Lls input adds the zero-sequence leakage component, and a x1/3 gain is applied at the output.

The six expressions and their roles in the assembled matrix are:

f(u)	Expression	Role
f_1	2*u[1]	L_{aa} diagonal
f_2	-u[1] + sqrt(3)*u[2]	L_{ab} off-diagonal
f_3	-u[1] - sqrt(3)*u[2]	L_{ac} off-diagonal
f_4	u[1]/2 + 3*u[3]/2 - sqrt(3)*u[2]	L_{bb} diagonal
f_5	u[1]/2 - 3*u[3]/2	L_{bc} off-diagonal
f_6	u[1]/2 + 3*u[3]/2 + sqrt(3)*u[2]	L_{cc} diagonal

where u[1] = L_{\alpha\alpha}, u[2] = L_{\alpha\beta}, u[3] = L_{\beta\beta}. The block takes three independent entries, exploiting the symmetry L_{\alpha\beta} = L_{\beta\alpha} that holds because \mathbf{L}_{s,\alpha\beta} is the similarity transform of a symmetric dq matrix. The output matrix is

\mathbf{L}_{s,abc} = \frac{1}{3}\begin{bmatrix} f_1 + L_{ls} & f_2 + L_{ls} & f_3 + L_{ls} \\ f_2 + L_{ls} & f_4 + L_{ls} & f_5 + L_{ls} \\ f_3 + L_{ls} & f_5 + L_{ls} & f_6 + L_{ls} \end{bmatrix}

Leakage handling

Stator leakage flux links all three phase windings equally and independently of the rotor. In the abc frame this appears as L_{\sigma s} on every diagonal of \mathbf{L}_{s,abc} and zero on every off-diagonal, giving L_{\sigma s}\mathbf{I}_3. The dq frame is a 2x2 representation and cannot capture the full structure of a 3x3 matrix. The Clarke transform projects out the zero-sequence subspace entirely. Any quantity that is equal across all three phases has a zero-sequence component that is invisible in dq. The leakage inductance must be introduced via two separate paths, one for each subspace.

Path A (αβ component): the vector [Lls 0 Lls] adds L_{\sigma s} to L_{mi,dd} and L_{mi,qq} while leaving L_{mi,dq} unchanged, which introduces L_{\sigma s}\mathbf{I} into the 2x2 dq matrix. A scalar multiple of the identity commutes with any rotation and passes through the Mdq→Mab similarity unchanged. After the Mab→Muvw αβ-to-abc expansion and the x1/3 gain, this contributes +\tfrac{2}{3}L_{\sigma s} on each diagonal and -\tfrac{1}{3}L_{\sigma s} on each off-diagonal.

Path B (zero-sequence component): the zero-sequence part of L_{\sigma s}\mathbf{I}_3 is the uniform matrix \tfrac{L_{\sigma s}}{3}\mathbf{1}_{3\times3}. It has no representation in the 2x2 dq frame and cannot pass through the similarity transform. It is injected directly inside Mab→Muvw via the scalar Lls port, which adds L_{\sigma s} uniformly to all nine entries before the x1/3 gain, contributing +\tfrac{1}{3}L_{\sigma s} to every entry.

The two paths together give L_{\sigma s}\mathbf{I}_3 at the abc output:

\underbrace{\frac{L_{\sigma s}}{3}\begin{bmatrix} 2 & -1 & -1 \\ -1 & 2 & -1 \\ -1 & -1 & 2 \end{bmatrix}}_{\text{Path A: αβ component}} + \underbrace{\frac{L_{\sigma s}}{3}\begin{bmatrix} 1 & 1 & 1 \\ 1 & 1 & 1 \\ 1 & 1 & 1 \end{bmatrix}}_{\text{Path B: zero-sequence component}} = \underbrace{L_{\sigma s}\begin{bmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{bmatrix}}_{L_{\sigma s}\mathbf{I}_3}

Path B is required in two distinct cases: when the neutral point has any active electrical connection, and unconditionally for numerical reasons.

Physical case. If the neutral point is connected to an external circuit, the constraint i_a + i_b + i_c = 0 no longer holds. Zero-sequence current can flow, and the model must include the inductance that opposes it. Path B provides this inductance: it contributes \tfrac{1}{3}L_{\sigma s}\mathbf{1}_{3\times3} to \mathbf{L}_{s,abc}, which acts on the zero-sequence component of the current. Without Path B, the model would predict zero inductive opposition to zero-sequence current, which is unphysical.

Numerical case. Even with a floating neutral, where zero-sequence current is identically zero, Path B is still required. The Path A contribution alone — the αβ-only matrix shown in the underbrace above — is rank-2 in the 3x3 abc representation. The eigenvalues of \mathbf{L}_{s,abc} assembled from Path A only are \{0,\, L_d^{\text{eff}},\, L_q^{\text{eff}}\}, with the zero eigenvalue corresponding to the zero-sequence direction [1,1,1]/\sqrt{3}. This is a structural property of the αβ-to-abc expansion: the zero-sequence direction is in the null space of any matrix produced by the similarity transform of a 2x2 dq matrix.

The PLECS Variable Inductor requires a non-singular \mathbf{L}_{s,abc} at every time step to invert it for the state-space formulation d\vec{i}_{abc}/dt = \mathbf{L}_{s,abc}^{-1}\,\vec{v}_L. A singular matrix breaks the solver. Path B’s contribution fills the zero-sequence gap, producing eigenvalues \{L_{\sigma s},\, L_d^{\text{eff}},\, L_q^{\text{eff}}\} — all positive provided L_{\sigma s} > 0.

This is why the init script enforces Lls > 0 and forbids Lls == 0. The model requires a strictly positive value to keep \mathbf{L}_{s,abc} invertible, regardless of how small the physical leakage is. For floating-neutral operation with negligible physical leakage, a small numerical value (e.g., 10^{-6} H) is sufficient. For neutral-connected operation, Lls must reflect the actual physical leakage so that the predicted zero-sequence current is correct.

The trailing zero vector

The [0 0 0 0 0 0 0 0 0] constant appended to Lss_uvw completes the 18-element control signal required by the PLECS Variable Inductor. The first 9 entries are the inductance matrix entries. The second 9 are their time derivatives dL_{ij}/dt. The derivatives are hardwired to zero because the rotation-derivative contribution has already been absorbed into \vec{e}_s in the Back-EMF construction section. Supplying the actual dL/dt here would double-count it.

Appendix A: Derivation of the synchronous-frame stator voltage equation

This appendix derives the synchronous-frame stator voltage equation, starting from the stator winding voltage equation in its own (stationary) frame. Space vectors are written as real 2x1 column vectors and the rotation operators as 2x2 matrices throughout.

In the stator (stationary) frame \vec{x}_{s_s} = [x_\alpha,\ x_\beta]^T; in the rotor (synchronous) frame \vec{x}_s = [x_d,\ x_q]^T. The inner subscript s identifies a stator winding quantity and r a rotor winding quantity. An outer subscript names the frame of expression when it is not the synchronous frame. So \vec{x}_{s_s} is a stator quantity in the stator frame and \vec{x}_{r_r} is a rotor quantity in the rotor frame. Note that \vec{x}_{s_s} has both subscripts equal to s: the inner identifies the stator winding, the outer identifies the stator (stationary) frame.

The 90° rotation matrix and the frame rotation matrix are

\mathbf{J} = \begin{bmatrix} 0 & -1 \\ 1 & 0 \end{bmatrix} \qquad \mathbf{T}(\theta) = \begin{bmatrix} \cos\theta & -\sin\theta \\ \sin\theta & \cos\theta \end{bmatrix}

\mathbf{J} takes the role of j in complex notation: \mathbf{J}^2 = -\mathbf{I} and \mathbf{J} rotates a 2-vector by 90°. \mathbf{T}(\theta) takes the role of e^{j\theta}: it rotates a 2-vector by angle \theta, satisfies \mathbf{T}(-\theta) = \mathbf{T}(\theta)^{-1} = \mathbf{T}(\theta)^T, and \mathbf{T}(\theta_1)\mathbf{T}(\theta_2) = \mathbf{T}(\theta_1+\theta_2). The two are related by

\frac{d}{dt}\mathbf{T}(\theta(t)) = \dot\theta\,\mathbf{J}\,\mathbf{T}(\theta(t)) = \dot\theta\,\mathbf{T}(\theta(t))\,\mathbf{J}

where the two equivalent forms follow from the fact that \mathbf{J} and \mathbf{T}(\theta) commute with each other (both being rotations) and with any scalar. \mathbf{J} does not commute with a general 2x2 matrix such as \mathbf{L}_{mi}. The rotor electrical angle is \theta and \omega = \dot\theta.

The stator winding voltage equation in its own (stationary) frame is

\vec{v}_{s_s} = R_s\,\vec{i}_{s_s} + \frac{d\vec{\psi}_{s_s}}{dt}

Two-winding flux coupling between stator and rotor gives

\vec{\psi}_{s_s} = L_s\,\vec{i}_{s_s} + L_{sr}\,\mathbf{T}(\theta)\,\vec{i}_{r_r}

Here L_s is the stator self-inductance, L_{sr} is the stator-rotor mutual inductance, and \mathbf{T}(\theta) rotates the rotor-frame current \vec{i}_{r_r} into the stator frame. For a PM machine, \vec{i}_{r_r} is a constant fictitious current representing the magnet flux. For a pure synchronous-reluctance machine it is zero.

To express each stator quantity in the rotor frame, apply \mathbf{T}(\theta)^{-1} = \mathbf{T}(-\theta) to it. Define

\vec{v}_s := \mathbf{T}(-\theta)\,\vec{v}_{s_s} \qquad \vec{i}_s := \mathbf{T}(-\theta)\,\vec{i}_{s_s} \qquad \vec{\psi}_s := \mathbf{T}(-\theta)\,\vec{\psi}_{s_s}

These are the synchronous-frame stator quantities used in the body of the document. The flux equation becomes

\vec{\psi}_s = L_s\,\vec{i}_s + L_{sr}\,\vec{i}_{r_r}

since the inverse rotation cancels the \mathbf{T}(\theta) on the rotor term. All three quantities now live in the same (rotor) frame.

The voltage equation requires the chain rule because of the time-varying rotation. From \vec{\psi}_{s_s} = \mathbf{T}(\theta)\,\vec{\psi}_s,

\frac{d\vec{\psi}_{s_s}}{dt} = \mathbf{T}(\theta)\,\frac{d\vec{\psi}_s}{dt} + \omega\,\mathbf{J}\,\mathbf{T}(\theta)\,\vec{\psi}_s

Substituting into the stator-frame voltage equation and left-multiplying by \mathbf{T}(-\theta):

\vec{v}_s = R_s\,\vec{i}_s + \frac{d\vec{\psi}_s}{dt} + \omega\mathbf{J}\,\vec{\psi}_s

This is the standard synchronous-frame stator voltage equation. The \omega\mathbf{J}\vec{\psi}_s term is the rotational EMF that appears when transforming a stator equation into a rotating frame.

Next, split the stator self-inductance as L_s = L_{\sigma s} + L_{sr}, where L_{\sigma s} := L_s - L_{sr} is the stator leakage inductance. Substituting into the synchronous-frame flux expression:

\vec{\psi}_s = L_{\sigma s}\vec{i}_s + L_{sr}(\vec{i}_s + \vec{i}_{r_r})

Define the magnetizing current \vec{i}_m := \vec{i}_s + \vec{i}_{r_r} and the magnetizing inductance L_m := L_{sr}, so that \vec{\psi}_s = L_{\sigma s}\vec{i}_s + L_m\vec{i}_m. The magnetizing flux is \vec{\psi}_{m,s} := L_m\,\vec{i}_m. The stator flux then decomposes as

\vec{\psi}_s = L_{\sigma s}\,\vec{i}_s + \vec{\psi}_{m,s}

In the linear two-winding model, \vec{\psi}_{m,s} = L_m(\vec{i}_s + \vec{i}_{r_r}) is linear in \vec{i}_s for fixed \vec{i}_{r_r}. In the saturating machine this model implements, \vec{\psi}_{m,s} is replaced by a nonlinear flux map of \vec{i}_s alone. The rotor current is fictitious and constant, or absent for a synchronous-reluctance machine. The saturating extension replaces L_m\,\vec{i}_m with the tabulated function \vec{\psi}_{m,s}(\vec{i}_s); see Flux-map lookup.

Substituting the flux decomposition into the synchronous-frame voltage equation:

\vec{v}_s = R_s\vec{i}_s + \frac{d}{dt}(L_{\sigma s}\vec{i}_s + \vec{\psi}_{m,s}) + \omega\mathbf{J}(L_{\sigma s}\vec{i}_s + \vec{\psi}_{m,s})

L_{\sigma s} is a constant parameter, so it commutes with d/dt. The magnetizing flux depends on \vec{i}_s alone, so its time derivative is

\frac{d\vec{\psi}_{m,s}}{dt} = \frac{\partial\vec{\psi}_{m,s}}{\partial\vec{i}_s}\,\frac{d\vec{i}_s}{dt} = \mathbf{L}_{mi}\,\frac{d\vec{i}_s}{dt} \qquad \mathbf{L}_{mi} := \frac{\partial\vec{\psi}_{m,s}}{\partial\vec{i}_s}

This defines \mathbf{L}_{mi} as the local 2x2 Jacobian of the flux map. Substituting:

\boxed{\;\vec{v}_s = R_s\,\vec{i}_s + (L_{\sigma s}\mathbf{I} + \mathbf{L}_{mi})\,\frac{d\vec{i}_s}{dt} + \omega\mathbf{J}\,(L_{\sigma s}\vec{i}_s + \vec{\psi}_{m,s})\;}

This is the synchronous-frame stator voltage equation in matrix form.

Appendix B: Initialization Script

The model’s mask initialization script accepts the following user inputs.

Required inputs:

• Id_r, Iq_r — row vectors of the dq current axis breakpoints (at least 2 elements each)
• Psid_r, Psiq_r — 2-D matrices of flux linkage values, dimensioned [length(Id_r), length(Iq_r)]
• Lls — stator leakage inductance (scalar, must be nonzero)

Optional inputs:

• Lmidd_r, Lmidq_r, Lmiqq_r — 2-D incremental inductance maps, same dimensions as the flux maps
• Nxdq — desired interpolated table size

Behavior:

• If the user supplies the incremental inductance maps, they are validated for dimension consistency and used directly.
• If the incremental inductance maps are not supplied, the script computes them by numerical differentiation of the flux maps using MATLAB’s gradient function.
• If Nxdq is supplied, the flux maps are re-interpolated onto an equally-spaced grid of the requested size. This is also done automatically when the incremental inductance maps need to be generated, since gradient requires uniform spacing.

The user-facing API requires only the flux maps. The dual-table option is for users who want to override the auto-generated incremental inductance.

The Lls == 0 guard:

The script does not allow zero stator leakage. The variable-inductance circuit element requires a nonzero diagonal to remain non-singular under all operating conditions. Users with negligible physical leakage should specify a small positive value rather than zero.