# Voltage-Behind-Reactance Non-Excited Synchronous Machine

**URL:** <https://forum.plexim.com/t/voltage-behind-reactance-non-excited-synchronous-machine/3073>\
**Category:** Support Solutions\
**Tags:** plecs, electric-machine, flux-maps, vbr\
**Created:** [May 15, 2026, 3:42pm UTC](https://forum.plexim.com/t/voltage-behind-reactance-non-excited-synchronous-machine/3073 "2026-05-15T15:42:29Z")\
**Posts on this page:** 1\
**Page:** 1

<div class="post-metadata">

**Author:** ![Beat\_Arnet](https://dub1.discourse-cdn.com/flex006/user_avatar/forum.plexim.com/beat_arnet/32/9_2.png) [@Beat\_Arnet](https://forum.plexim.com/u/Beat_Arnet)\
**Post date:** [May 15, 2026, 3:42pm UTC](https://forum.plexim.com/t/voltage-behind-reactance-non-excited-synchronous-machine/3073/1 "2026-05-15T15:42:29Z")

</div>

![Mapping from the s-plane to the z-plane under z = e^{sT}](https://europe1.discourse-cdn.com/flex006/uploads/plecs/original/2X/3/31737577be318cc655ba526f4b8e9730ecbd1100.png)

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](#appendix-a-derivation-of-the-synchronous-frame-stator-voltage-equation) . 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](#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](#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](#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](#leakage-handling).

## Architecture

 ![Mapping from the s-plane to the z-plane under z = e^{sT}](https://europe1.discourse-cdn.com/flex006/uploads/plecs/original/2X/3/3581e2dd30b61852d109a9cceed4a23fe909decc.png)

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](#leakage-handling) in [Transformation to the abc frame](#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](#appendix-a-derivation-of-the-synchronous-frame-stator-voltage-equation) 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

 ![Mapping from the s-plane to the z-plane under z = e^{sT}](https://europe1.discourse-cdn.com/flex006/uploads/plecs/original/2X/f/f850934ec3557425b4d6b47780f5e80f26944ce2.png)

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](#appendix-b-initialization-script) 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:

 ![Mapping from the s-plane to the z-plane under z = e^{sT}](https://europe1.discourse-cdn.com/flex006/uploads/plecs/original/2X/0/0e470c5bd360b078ecc971ffd0ab168a9831f2f4.png)

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](#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](#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.

 ![](https://europe1.discourse-cdn.com/flex006/uploads/plecs/original/2X/b/b2aec9810bbd30278bbfad17d8318bd66c46b027.png)

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](#appendix-a-derivation-of-the-synchronous-frame-stator-voltage-equation)) 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.

 ![](https://europe1.discourse-cdn.com/flex006/uploads/plecs/original/2X/5/55189387995db5f643054299b65feb30c1d0a8ea.png)

### `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.

 ![](https://europe1.discourse-cdn.com/flex006/uploads/plecs/original/2X/8/85b01d3d4b496f68a84271cc8425aa2104f9537e.png)

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](#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.
