Comparative Study of Sphere Decoding Algorithm and FCS-MPC for PMSMs in Aircraft Application

1. Introduction

Advances in power electronics and modern control techniques are facilitating the integration of more electrical subsystems into aircrafts, reinforcing flight control and improving passenger and cargo comfort [1]. The transition from hydraulically, pneumatically, and mechanically powered subsystems to electrically driven systems—known as the “More Electric Aircraft” (MEA) concept—has gained traction due to its potential to reduce weight, maintenance costs, and carbon emissions while improving fuel efficiency. As the MEA concept continues to evolve, onboard power demand and generator capacity have risen significantly, increasing from a few hundred kilovolt-amperes before World War II to megavolt-amperes in recent decades [2].

Increasing the power capacity of aircraft electrical systems inevitably adds weight and volume due to larger machines, drives, and thicker transmission lines. Even when inactive, this added weight increases fuel consumption. Therefore, minimizing size and weight while maintaining efficient power generation is a key priority. Permanent magnet (PM) machines are widely preferred in aerospace due to their high efficiency, reduced losses, and compact design. Unlike electrically excited machines, PMs eliminate excitation losses , simplifying thermal management and improving overall performance [3]. They also eliminate mechanical components like brushes and commutators, reducing maintenance and sparking risks [4,5]. PMSMs are categorized by back-EMF waveform into sinusoidal and trapezoidal types, where sinusoidal PMSMs offer smooth torque when harmonics are well-filtered, but they require of high switching frequencies which in turn result in increased losses. Trapezoidal PMSMs are cost-effective but suffer from core losses at high speeds. Based on magnet placement, PMSMs are classified as surface-mounted (SPMSM), insert-type, or interior-mounted. SPMSMs are particularly suited for aerospace applications due to their superior power density, efficiency, and reliability [2]. Their high torque-to-weight ratio is crucial for minimizing aircraft weight while maintaining performance. The absence of excitation windings enhances efficiency and thermal management, which is critical in aerospace environments with limited cooling options. Additionally, their fast dynamic response, smooth torque, and brushless design ensure mechanical stability, reduced maintenance, and long-term reliability under extreme conditions such as high vibrations and temperature variations.

Various drive control strategies as shown in Figure 1 have been proposed for SPMSM drives, including Field-Oriented Control (FOC), Direct Torque Control (DTC), constant volt-per-frequency control (V/F) and constant current-per-frequency control (I/F). Among these, FOC is the most commonly employed due to its high state performance. FOC decouples torque and flux, enabling independent control of these parameters. This approach enables smooth torque and speed regulations by using Proportional Integral(PI) controllers to manage the currents in the d- and q- axis frame [6], where Pulse Width Modulation is typically employed to implement the control. Direct Torque Control (DTC) is a control strategy that directly regulates the motor’s torque and flux by selecting appropriate voltage vectors. It employs hysteresis controllers to ensure that torque and flux remain within a predefined bound. By continuously measuring the stator flux and torque, the system selects the optimal voltage vector from a predefined lookup table to restore the torque and flux to their desired range. Unlike the FOC, DTC does not rely on current control loops, simplifying its computational requirements. However, a drawback of DTC is its tendency to produce a higher torque ripple[7].

MPC has recently emerged as a powerful strategy in motor drives systems, enabled by the advancement in faster microprocessors with enhanced computational capabilities. Its key features such as simplicity in design, ease of implementation, explicit handling of constraints and capacity to manage complex, non-linear dynamics have facilitated its widespread application particularly in power electronics.

MPC can be classified into Continuous control set MPC (CCS-MPC) and finite Control set MPC (FCS-MPC) as shown in Figure 2. The CCS-MPC computes a continuous control signal, which is then processed through modulation to generate the desired output voltage for the power converter. The modulation strategy choice depends on the converter topology, a notable advantage of CCS-MPC is its ability to maintain a fixed switching frequency, a critical requirement in power electronic applications. Two widely used CCS-MPC strategies are the generalised predictive control (GPC) and Explicit Model Predictive Control (EMPC). The former is well suited for linear and unconstrained problems while the latter is designed for systems with nonlinear dynamics and constraints. However, both face challenges in their formulation particularly when applied to power converters due to the inherent complexity of these systems [8]. FCS-MPC, also known as direct MPC, simplifies the control process by leveraging the discrete nature of power converters and eliminating the need for an external modulator [9]. This approach directly incorporates the system’s hybrid nature making it conceptually straightforward. Like CCS-MPC, FCS-MPC is also categorized into two sets which are the Optimal Switching Vector (OSV-MPC) and Optimal Switching Sequence MPC (OSS-MPC). The former is the most widely used method in power electronics, OSV-MPC employs the discrete output voltage vectors of the converter as its control set. The optimization process is reduced to a simple enumerated search, where one voltage vector is applied per switching period, but this leads to variable switching frequency unless additional constraints are imposed. To address the variable switching frequency limitation of OSV-MPC, OSS-MPC evaluates a limited switching sequence within each control period. By treating switching times as additional decision variables, OSS-MPC enhances control performance and mitigates the drawback associated with OSV-MPC [10].

MPC optimizes motor performance by predicting future system behavior. It defines control actions by minimizing a cost function that represents the desired system performance, offering precise control and improved flexibilty for modern drive applications. MPC has been applied with FOC to control the internal stator current (MPCC) and also in DTC to directly control the machine’s torque and flux (MPTC). MPC can also be classified into single-step MPC and multi-step MPC / Long Prediction Horizon (LPH) as shown in Figure 3. FCS-MPC has been widely applied to drives, demonstrating substantial benefits from incorporating a long prediction horizon (LPH) in terms of steady state operation and reduced current distortions. In this sense, Li et.al., [11] confirmed that LPH enhances steady state performance of three phase drive systems, while Dorfling et.al.,[12] validated the advantages of LPH-FCS-MPC for three level neutral point clamped converters with RL loads. Moreover, LPH-FCS-MPC significantly improves performance in higher-order systems, particularly by reducing harmonic distortions [13]. While the multi-step offers superior performance compared to the single-step approach, it comes with a significantly higher computational burden. To address this a modified sphere decoding algorithm was proposed in [14] significantly reducing computational load but introducing unpredictability in execution time for each sampling interval. To overcome this unpredictability the K-best SDA was introduced in [15] providing fixed computational burden per sampling instant. A comparative analysis in [16] demonstrates that the K-best approach provides fixed computational burden per sampling instant with tuneable Kb value but the modified sphere decoding algorithm provides better performance providing lowest THD.

The application of Model Predictive Control (MPC) in PMSM faces several challenges with its performance heavily reliant on the accuracy of the motor’s mathematical model. Key issues include variations in motor parameters such as its internal resistance, which is influenced by temperature, its inductance, which can vary as a function of current and permanent flux linkage which in turn is affected by environmental conditions, and additionally, manufacturing errors that may result in parameter inaccuracies either at runtime or during initial measurements. These factors collectively degrade the performance of MPC systems. To address parameter uncertainties, model-free approaches have been introduced in [17] to mitigate the impact of inaccurate motor parameters. However, this method has two notable drawbacks; firstly, two current measurements are required in each sampling period for the detection of current spikes and secondly, it suffers from the issue of stagnant current variation updates, where the switching state can be unchanged for several sampling periods, adversely affecting current prediction performance. Nonetheless, these drawbacks can be eliminated while retaining the merits of a model-free predictive control approach [18]. Alternative methods have also been explored, for instance, an improved current MPC for parameter mismatch in [19] introduced a novel current update mechanism. Unlike earlier methods, where the current variation did not update within each control period (causing spikes at high speeds and impairing id and iq current tracking), this mechanism ensures that the predicted and measured stator currents are updated during each control period significantly improving performance.

In [11] the effect on parameter mismatch on PMSM was comprehensively analyzed, this study implemented a moving horizon estimator as a disturbance observer, enabling estimation of external disturbances and their integration into the direct control problem. This approach enhanced the robustness of FCS-MPC across a wide range of operating conditions. In [20], the effect of positive and negative parameter mismatch on the THD and Switching Frequency was studied on a conventional MPC and the proposed velocity model formulated by adding an integrator to the control problem without augmenting the states. In [21], the trade-off between switching frequency and THD was studied for induction motors using monte carlo simulations (10,000 random samples), but the comparison was only based on the fitted results, with the effects of sampling time and weighing factor not being fully disclosed. Likewise, in [22] another comparison for switching frequency and THD was presented, which did not consider the effect of sampling time and the weighing factor. In [13] the comparison between 1-norm and 2-norm approach was carried out, considering various prediction horizons from $N_p=1$ to $N_p=20$ for different switching frequencies for an induction machine, but a clear comparison between the sampling time and the weighing factor effect on the THD and switching frequency was missing. In this sense, to the best of our knowledge, none of the above references have looked at the direct and complete mapping from these two variables (sampling time and weighing factor) to the control performance (THD and Switching Frequency). This is the main gap that is being fulfilled in this paper, where a comprehensive comparison between the SDA and FCS-MPC for SPMSM is presented, looking at the effects of parameter mismatch on the compared models.

2. Problem Formulation for LPH PMSM Control

2.1. Performance Metrics

In this work, two primary performance metrics are employed: total harmonic distortion (THD) and switching frequency. These metrics provide a quantitative basis for evaluating both the quality of the output waveform and the dynamic behavior of the control strategy.

2.1.1. Total Harmonic Distortion (THD)

THD quantifies the extent to which a waveform deviates from its fundamental sinusoidal component. Consider a signal $x\left(t\right)$ represented as:

$$x\left(t\right)=\mathrm{A}·cos(\omega t-\mathrm{ph})+(\mathrm{higher}\mathrm{harmonics}),$$

(1)

where the amplitude of the fundamental component is determined by first computing its cosine and sine coefficients:

$$a_{cs}={\displaystyle \frac{2}{T}}{\int }_0^{T}x\left(t\right)cos\left(\omega t\right)dt,b_{sn}={\displaystyle \frac{2}{T}}{\int }_0^{T}x\left(t\right)sin\left(\omega t\right)dt.$$

(2)

The amplitude of the fundamental is then given by:

$$\mathrm{A}=\sqrt{a_{cs}^2+b_{sn}^2}.$$

(3)

The root-mean-square (RMS) value of the entire waveform is computed as:

$${\mathrm{RMS}}^2={\displaystyle \frac{2}{T}}{\int }_0^{T}x{\left(t\right)}^{2}dt,$$

(4)

and therefore, the Total Harmonic Distortion (THD) is subsequently defined by:

$$I_{THD}=\sqrt{{\displaystyle \frac{{\mathrm{RMS}}^2}{{\mathrm{amp}}^2}}-1}.$$

(5)

This metric effectively measures the relative contribution of higher harmonics in the signal.

2.1.2. Switching Frequency

The switching frequency $f_{sw}$ reflects the average number of switching events per unit time. It is defined by the formula:

$$f_{sw}=\underset{M\to \infty }{lim}{\displaystyle \frac{1}{m{c}_{k}M{T}_s}}\sum _{k=0}^{M-1}{\parallel \Delta u\left(k\right)\parallel }_1$$

(6)

where $\Delta u\left(k\right)$ represents the change in the control signal at time step k, $T_s$ is the sampling period, and $m_{ck}$ is the number of power semiconductor switches in the power converter of interest which is 3. $c_k$ is the converter dependent correction factor which is $c_k=2$ for two level converter [13]. This metric allows us to assess the dynamic behavior and efficiency of the control strategy, particularly with respect to the operation of the power converter.

2.2. Model Predictive Control Formulation for PMSM Drive

2.2.1. SPMSM Mathematical Model

In [23], several assumptions were made regarding the SPMSM were adopted in this article: magnetic saturation, cogging torque, and cross-saturation effects in [24] were neglected. This implies that the d-axis inductance is approximately equal to the q-axis inductance, as well as to the inductance in the stationary alpha and beta axes.

The mathematical model of the PMSM drive in the dq axis can be presented as:

$$V_d=R_s{I}_d+L{\displaystyle \frac{d{i}_d}{d_t}}-{\omega }_{e}L{I}_q$$

(7)

$$V_q=R_s{I}_q+L{\displaystyle \frac{d{i}_q}{d_t}}+{\omega }_{e}L{I}_d+{\omega }_e{\psi }_m$$

(8)

$$T_e=1.5p{\psi }_m{I}_q$$

(9)

$$T_e-T_l=n{\displaystyle \frac{d{w}_m}{d_t}}+B{w}_m$$

(10)

where:

• $I_d$ and $I_q$ are the currents on the d- and q-axis, respectively.
• $V_d$ and $V_q$ are the d-axis and q-axis voltage inputs.

2.2.2. Model Predictive Current Control

MPC consist of three components: a prediction model which is a mathematical or state space representation of the system, a cost function and an optimisation algorithm. The current dynamics of a Permanent Magnet Synchronous Motor (PMSM) in the d-q reference frame are given by:

$${\displaystyle \frac{d{I}_d}{dt}}=-{\displaystyle \frac{R_s}{L_d}}{I}_d+{\omega }_e{I}_q+{\displaystyle \frac{V_d}{L_d}}$$

(11)

$${\displaystyle \frac{d{I}_q}{dt}}=-{\displaystyle \frac{R_s}{L_q}}{I}_q-{\omega }_e{I}_d+{\displaystyle \frac{V_q-{\omega }_e\left(k\right){\psi }_m}{L_q}}$$

(12)

The FCS-MPC in Figure 4 follows a conventional single step prediction, which uses Equation (11) and Equation (12) as its prediction model using the forward Euler discretization method, thus the design has to be done in the discrete time domain. The forward Euler discretization is given by:

$$i_d(k+1)=i_d\left(k\right)+T_s\left(-{\displaystyle \frac{R_s}{L_d}}{i}_d\left(k\right)+{\omega }_e\left(k\right){\displaystyle \frac{L_q}{L_d}}{i}_q\left(k\right)+{\displaystyle \frac{v_d\left(k\right)}{L_d}}\right),$$

(13)

$$i_q(k+1)=i_q\left(k\right)+T_s\left(-{\displaystyle \frac{R_s}{L_q}}{i}_q\left(k\right)-{\omega }_e\left(k\right){\displaystyle \frac{L_d}{L_q}}{i}_d\left(k\right)-{\displaystyle \frac{{\psi }_m{\omega }_e\left(k\right)}{L_q}}+{\displaystyle \frac{v_q\left(k\right)}{L_q}}\right).$$

(14)

Because the mechanical time constant of the SPMSM is much larger than the electromagnetic time constant,the motor rotor speed can be considered as constant during two control periods namely ${\omega }_e\approx {\omega }_e(k+1)$. In addition, according to the first order forward euler discretization, the next instant $\theta$ can simply be obtained as $\theta (k+1)=\theta \left(k\right)+{\omega }_e\left(k\right)T_s$. At this stage, the desired control targets are incorporated into a cost function $J_N$, also known as the objective function. This function forecasts the system’s behavior over a finite prediction horizon N, expressing it in terms of tracking errors. For single-step Model Predictive Control (MPC), the prediction horizon is set to $N=1$, and the cost function is referred to as J.

$$J=\parallel I_d{(k+1)}^{*}-I_d(k+1)\parallel +\parallel I_q{(k+1)}^{*}-I_q(k+1)\parallel .$$

(15)

Many works in power electronics use the $L_1$-norm of the error vector for its computational simplicity [25,26,27,28], although recent literature suggests the $L_2$-norm for improved closed-loop performance, resulting in closer reference tracking and reduced ripple [13,29]. Another common criterion in cost function design is minimizing the control effort. This criterion recognizes that control actions not only impact system output but also come with costs, such as switching effort in power converters. A weighting factor ${\lambda }_u$ is introduced to balance this trade-off. A smaller value for ${\lambda }_u$ (close to but not zero) results in faster system dynamics. When the error term in the cost function approaches zero (i.e., when the measured output closely matches the reference), the second term, related to control effort, becomes more dominant, guiding the system towards the desired steady-state behavior.

$$J=\parallel I_d^{*}(k+1)-I_d{(k+1)\parallel }_2^2+\parallel I_q^{*}(k+1)-I_q{(k+1)\parallel }_2^2+{\lambda }_u{\parallel \Delta U\parallel }_2^2$$

(16)

This type of cost function in Equation (16), when used in FCS-MPC, can ensure system stability and improve closed-loop performance. The most common implementation of Finite Control Set MPC (FCS-MPC) considers a one-step prediction horizon, where the stage cost does not include the control effort term, as seen in several references [25,26,27,28,30,31].

The computational stage in Equation (15) is more efficient when using the $L_1$-norm instead of the $L_2$-norm. This is advantageous from a computational perspective, particularly because Model Predictive Control (MPC) must operate in real time, often within a few tens of microseconds. Consequently, much of the literature focuses on the $L_1$-norm. However, as noted in [13], using an $L_2$-norm ensures closed-loop stability (or practical stability). In contrast, the $L_1$-norm can lead to closed-loop instability unless ${\lambda }_u\ne 0$. Setting ${\lambda }_u=0$ is not advisable as it results in the converter operating at its maximum switching frequency, and then, the stated switching frequency is limited to the chosen sampling interval. Moreover, operating at maximum switching frequency does not guarantee better performance, as demonstrated in Table 2 of [16], the comparison involved two MPC formulations: one using the $L_2$-norm with ${\lambda }_u=0$, and the other an unconstrained solution where ${\lambda }_u\ne 0$.

Results using Figure 4 with Equation (15) and Equation (16) are presented in Table 1 and Table 2, respectively, showing sampling frequency and THD values for various sampling times. Results in Table 1 were operating at maximum switching frequency as ${\lambda }_u=0$, and in Table 2 control effort was set to ${\lambda }_u=5e-06$.

Simulations were conducted for Figure 4 with Equation (15) and Equation (16), respectively, with a prediction horizon corresponding to $N_p=1$ and control effort tuning ranging from 0 to 4, incremented in steps of 0.01, in a bid to control the average switching frequency. Stability was lost for the 1-norm formulation(Equation (15)) when ${\lambda }_u$ exceeded 2.53.

The cost function in Equations (15) and (16) is evaluated for all possible control input values in Figure 5. However, as the prediction horizon increases, the number of possible values grows exponentially, making this approach computationally demanding. As previously discussed, increasing the prediction horizon offers benefits such as improved closed-loop stability and reduced THD. To balance these advantages with computational efficiency, the optimization problem is reformulated as follows.

2.3. Long Prediction Horizon MPC Formulation for PMSM Drive

The Figure 6 represents the reformulated optimization problem for long prediction with the cost function block illustrated in Figure 4 incorporated in the sphere decoding algorithm block. Here, the system model is formulated as a continuous-time prediction model in a state space form by solving the set of equations:

$${\displaystyle \frac{dx\left(t\right)}{dt}}=Ax\left(t\right)+Bu\left(t\right)+E,$$

(17)

$$y\left(t\right)=Cx\left(t\right)+Du\left(t\right),,$$

(18)

where the matrices $A\in {\mathbb{R}}^{n_x\times n_x}$ and $B\in {\mathbb{R}}^{n_x\times n_u}$ are the system and input matrices in continuous-time, with dimensions of $(n_x\times n_x)$ and $(n_x\times n_u)$, respectively. The matrices A and B defined are dependent on the system and can be considered to be both time-variant or time-invariant. Moreover, the output matrix $C\in {\mathbb{R}}^{n_y\times n_x}$ is time-invariant with a dimension of $(n_y\times n_x)$, which is used to define y from the available x. Here, x stands for the plant state vector consisting of $n_x$ number of system states ${\left[\begin{array}{cc}{I}_d& I_q\end{array}\right]}^{\top }$ that fully defines the dynamics of the system.

The controlled plant-input vector, u, comprises $n_u$ number of control plant-input which are the manipulated variables that belong to an integer-valued finite control set (FCS), $U=V^{n_u}$, which is expressed as $S(a,b,c)$=$\{0,1\}$. The output vector, y, is composed of $n_y$ output variables of the system which are a function of the state variables ${\left[\begin{array}{cc}{I}_d& I_q\end{array}\right]}^{\top }$. Then, the called feed-forward matrix D is neglected or assumed to be zero, as in most power electronics control applications this direct relationship is generally absent without any loss of generality.

Thus, the PMSM system in state-space form can be represented by the following matrices:

$${\mathbf{Ac}}_{\mathbf{xy}}=\left[\begin{array}{cc}-{\displaystyle \frac{R_s}{L_d}}& {\omega }_e\left(k\right)\\ -{\omega }_e\left(k\right)& -{\displaystyle \frac{R_s}{L_q}},\end{array}\right]$$

$${\mathbf{Bc}}_{\mathbf{xy}}=\left[\begin{array}{cc}{\displaystyle \frac{1}{L_d}}& 0\\ 0& {\displaystyle \frac{1}{L_q}}\end{array}\right]*T_{dq}*T_r*V_{dc},$$

$$\mathbf{C}=\left[\begin{array}{cc}1& 0\\ 0& 1\end{array}\right],$$

$$\mathbf{D}=\left[\begin{array}{cc}0& 0\\ 0& 0\end{array}\right],$$

$${\mathbf{Ec}}_{\mathbf{xy}}=\left[\begin{array}{c}0\\ {\psi }_m{\omega }_e\left(k\right)/L_q\end{array}\right].$$

As stated earlier the MPC controller is in the discrete time-domain, therefore the above matrices are discretised using the forward Euler approximation such that:

$$A=I_2+A_{c\_xy}·T_s,$$

$$B=B_{c\_xy}·T_s,$$

$$\psi =E_{c\_xy}·T_s.$$

Matrix B is multiplied by the Clarke and Park transforms in the dq frame, where the Park transform ($T_{dq}$) is given by:

$$T_{dq}=\left[\begin{array}{cc}cos\left(\theta \right)& sin\left(\theta \right)\\ -sin\left(\theta \right)& cos\left(\theta \right),\end{array}\right]$$

and the Clarke transform ($T_r$) by:

$$T_r=\sqrt{{\displaystyle \frac{2}{3}}}\left[\begin{array}{ccc}1& -{\displaystyle \frac{1}{2}}& -{\displaystyle \frac{1}{2}}\\ 0& {\displaystyle \frac{\sqrt{3}}{2}}& -{\displaystyle \frac{\sqrt{3}}{2}}\end{array}\right].$$

By inserting Equation (17) into Equation (18), we obtain:

$$Y\left(k\right)=\Gamma x\left(k\right)+{\rm Y}U\left(k\right)+\psi ,$$

(19)

where the matrices $\Gamma$ and ${\rm Y}$ are defined as :

$$\Gamma =\left[\begin{array}{c}CA\\ C{A}^2\\ ⋮\\ C{A}^{N_p}\end{array}\right]\mathrm{and}{\rm Y}=\left[\begin{array}{cccc}CB& 0_{ny\times 3}& \cdots & 0_{ny\times 3}\\ CAB& CB& \cdots & 0_{ny\times 3}\\ ⋮& ⋮& \ddots & ⋮\\ C{A}^{N_p-1}B& C{A}^{N_p-2}B& \cdots & CB\end{array}\right].$$

Recalling that the cost function in Equation (16) has two parts, the tracking and the control effort for the tracking effort, the tracking effort is represented as $J_a$:

$$J_a={\parallel Y^{*}\left(k\right)-(\Gamma x\left(k\right)+{\rm Y}U\left(k\right)+\psi )\parallel }_2^2.$$

(20)

Similarly, the second term in the cost function (Equation (16)), the control effort can be rewritten as:

$$\begin{array}{ccc}\hfill J_b& =& \sum _{i=k}^{k+N_p-1}{\lambda }_u{\parallel \Delta u\left(i\right)\parallel }_2^2={\lambda }_u\sum _{i=k}^{k+N_p-1}{(u\left(i\right)-u(i-1))}^{\top }(u\left(i\right)-u(i-1)),\hfill \\ & =& {\lambda }_u{(SU\left(k\right)-Eu(k-1))}^{\top }(SU\left(k\right)-Eu(k-1))={\lambda }_u{\parallel SU\left(k\right)-Eu(k-1)\parallel }_2^2.\hfill \end{array}$$

(21)

Where matrices S and E are defined as:

$$S=\left[\begin{array}{cccc}{I}_3& 0_{3\times 3}& \cdots & 0_{3\times 3}\\ -I_3& I_3& \cdots & 0_{3\times 3}\\ 0_{3\times 3}& -I_3& \cdots & 0_{3\times 3}\\ ⋮& ⋮& \ddots & ⋮\\ 0_{3\times 3}& 0_{3\times 3}& \cdots & I_3\end{array}\right]\mathrm{and}E=\left[\begin{array}{c}{I}_3\\ 0_{3\times 3}\\ 0_{3\times 3}\\ ⋮\\ 0_{3\times 3}\end{array}\right].$$

Summing up Equation (21) and Equation (20) yields the cost function in vector notation:

$$J=\parallel Y^{*}{\left(k\right)-\Gamma x\left(k\right)-{\rm Y}U\left(k\right)-\psi \parallel }_{\tilde{Q}}^2+{\lambda }_u{\parallel SU\left(k\right)-Eu(k-1)\parallel }_2^2.$$

(22)

Then, the cost function can be reformulated according to [14,32] as:

$$J=(U{\left(k\right)}^{T}HU\left(k\right)+2{\left(\Theta \left(k\right)\right)}^{T}U\left(k\right)+\theta \left(k\right)),$$

(23)

where:

$$H={{\rm Y}}^T\tilde{Q}{\rm Y}+{\lambda }_u{S}^{T}S,$$

(24)

$${\left(\Theta \left(k\right)\right)}^T=-{(Y^{*}\left(k\right)-\Gamma x\left(k\right)-\psi )}^T\tilde{Q}{\rm Y}-{\lambda }_u{\left(E_u(k-1)\right)}^{T}S,$$

(25)

and,

$$\theta \left(k\right)=\parallel Y^{*}{\left(k\right)-\Gamma x\left(k\right)-\psi |}_{\tilde{Q}}^2+{\lambda }_u{\parallel E_u(k-1)\parallel }_2^2.$$

(26)

The term $\theta \left(k\right)$ is independent of $U\left(k\right)$ and has no effect on the optimal solution, thus can be omitted. Therefore, when substituting back into the , it gives:

$$U_{\mathrm{Opt}}\left(k\right)=argmin{\left(U\left(k\right)+H^{-1}\Theta \left(k\right)\right)}^T\times H\left(U\left(k\right)+H^{-1}\Theta \left(k\right)\right)$$

(27)

also subject to the constraint conditions $U\left(k\right)\in \mathcal{U}$, and ${\parallel \Delta u\left(l\right)\parallel }_{\infty }\le 1$.

The solution obtained from Equation (27) without considering these constraints is known as the unconstrained solution. While this solution may not be directly applicable to the converter due to disregarding the system constraints, it serves as a starting point for finding the constrained solution given by:

$$U_{\mathrm{unc}}\left(k\right)=-H^{-1}\Theta \left(k\right)$$

(28)

The matrix H being symmetric and positive definite ($\lambda >0$), ensures its Cholesky factorization exists and yields to a unique lower triangular matrix V. This property is crucial for the optimization process, as it allows for efficient computation of the constrained solution. Cholesky factorization decomposes the positive definite matrix H into the product of a lower triangular matrix V and its transpose, $V^{T}V=H$[33], i.e., with ${\overline{U}}_{\mathrm{unc}}\left(k\right)=V{U}_{\mathrm{unc}}\left(k\right)$.

The cost function can then be expressed in terms of V and ${\overline{U}}_{\mathrm{unc}}\left(k\right)$ as:

$$U_{\mathrm{Opt}}\left(k\right)=argmin{\parallel VU\left(k\right)-{\overline{U}}_{\mathrm{unc}}\left(k\right)\parallel }_2^2$$

(29)

Geometrically, the matrix V serves as a lattice generator matrix, defining a discrete space where the optimal solution resides. This transformation allows us to represent the optimization problem in a new space where the constraints are simplified and the problem becomes computationally tractable for long prediction horizons. In this transformed space, the optimization problem involves finding the control sequence $U_{\mathrm{Opt}}\left(k\right)$ that minimizes the distance J to the unconstrained solution $U_{\overline{\mathrm{unc}}}\left(k\right)$, while satisfying all system constraints. Operating in the transformed space reduces the computational complexity of the optimization problem, making it feasible for practical applications in power converters.

2.4. Sphere Decoding Algorithm Implementation for the PMSM Drive

The challenge in Equation (29) which relates to finding the shortest distance between two vectors, is known as the closest point problem, where several conventional algorithms have been already addressed in [34].

Two algorithms are highlighted, known as the Pohst and Kannan algorithms, with the former using a hypersphere for the devising of methods regarded as sphere decoders such as the SDA although it was been primarily applied in communication theory. Nonetheless, in the context of power electronics, specifically in addressing the Long Prediction Horizon (LPH)-FCS-MPC problem for power converters, a modified sphere decoding algorithm has been developed and adapted [14],this approach tailors the SDA to suit the requirements of power converter control. This problem can be interpreted as minimizing the hypersphere of radius $\rho$, which corresponds to the value of the cost function, centered at $U_{\mathrm{unc}}\left(k\right)$. The Sphere Decoding Algorithm (SDA) evaluates only those candidate solutions that lie within the computed hypersphere, selecting the one with the smallest cost function as the optimal solution. Specifically, the optimal solution $U^{*}\left(k\right)$ is the n-dimensional lattice point that has the minimum Euclidean distance from the unconstrained solution ${\overline{U}}_{\mathrm{unc}}$ as shown in Figure 5. Equivalently, this solution is the unique lattice point enclosed within the smallest hypersphere centered at ${\overline{U}}_{\mathrm{unc}}$ with radius $\rho$.

A key advantage of SDA is that it ensures optimality with significantly higher computational efficiency compared to the exhaustive enumeration used in Figure 7, coupled with the fact that it does not need to evaluate every possible control input, reducing the exponential computational burden when the prediction horizon is increased. The efficiency of SDA is largely dependent on the initial choice of $\rho$, which must be small enough to minimize the number of evaluated nodes, but not too small as to exclude all possible solutions. The initial radius is computed as proposed in [35]:

$$\rho =min\{{\rho }_a,{\rho }_b\},$$

(30)

where

$${\rho }_a={\parallel {\overline{U}}_{\mathrm{unc}}-V{U}_{\mathrm{bab}}\left(k\right)\parallel }^2/,,$$

(31)

and,

$${\rho }_b={\parallel {\overline{U}}_{\mathrm{unc}}-V{U}_{\mathrm{ed}}\left(k\right)\parallel }^2.$$

(32)

These radii are based on two different approximations: the “Babai estimate” $U_{\mathrm{bab}}$ and the “educational guess” $U_{\mathrm{ed}}$. The Babai estimate represents the rounded unconstrained solution to the nearest integer vector:

$$U_{\mathrm{bab}}\left(k\right)=\lfloor U_{\mathrm{unc}}\left(k\right)\rceil /,,$$

(33)

whereas the educational guess is given by the previously applied optimal solution $U^{*}(k-1)$ shifted by one time step. Once $U^{*}\left(k\right)$ is determined, following the receding horizon policy, only its first element is applied by the converter, while the remaining elements are discarded. The process then repeats at the next time step over a shifted prediction horizon, incorporating updated measurements and estimates.

After determining the initial radius, the algorithm explores the search tree, evaluating one switching element at a time, effectively reducing the problem to a one-dimensional search at each step. The distance of the first selected element from the unconstrained solution is calculated, and if it falls within the sphere’s boundary, the algorithm proceeds with that branch. Otherwise, the entire branch is discarded, as no further elements from that branch can contribute to a feasible solution within the sphere. This pruning mechanism is fundamental to the SDA and makes it well-suited for long prediction horizons.

Once a feasible solution is found, the radius is recalculated and updated, further refining the search space. This iterative process continues, progressively tightening the sphere’s boundary, thereby eliminating more suboptimal solutions. The final solution obtained at the end of this process is declared optimal.

Figure 8 illustrates the operation of the sphere decoding algorithm for $N_p=1$. At each node in the search tree, the algorithm calculates the cumulative squared deviation up to the current step, and then compares it against the squared radius. If the current node’s cost exceeds the radius (yellow triangle), the branch is pruned (grey) and the algorithm backtracks (orange arrow) to explore other nodes. This simplifies the verification process and enhances computational efficiency, ensuring that the sphere decoding algorithm can quickly and effectively prune branches of the search tree that cannot lead to a valid solution within the prescribed radius. When the algorithm reaches a bottom node and finds a candidate solution within the current sphere, it tightens the sphere around this solution. The process is repeated over and over, further tightening the sphere to fewer and fewer candidate solutions, thus refining the search. The SDA terminates when the entire search tree has been explored; only at this point can the algorithm guarantee that the found solution is optimal (green circle) having examined all possible candidate solutions within the defined sphere.

3. Simulation Results and Comparison

In this section we analyze the impact of the sampling time and a long prediction horizon on the performance of the PMSM, such as the THD, Switching Frequency and the effect of the weighing factor ${\lambda }_u$ on the switching frequency. This was also implemented for long predictions of $N_p=3$ and 5.

3.1. Effect of Sampling Time on Sphere Decoding Algorithm for PMSM

The SDA algorithm was implemented for a prediction horizon of 2, 3, and 5 over sampling times ranging from 1 to 80 microseconds. Figure 9 shows the trends for FCS-MPC and SDA for various sampling times. This comparison revealed that FCS-MPC outperformed SDA in terms of THD for sampling times less than 10 $\mu$$\mathrm{s}$ However, for sampling times greater than 10 $\mu$$\mathrm{s},$ SDA consistently exhibited lower THD than FCS-MPC, albeit with a significant issue of undershoot as illustrated in Figure 9.

The SDA introduced in Figure 10 was implemented for a two-step prediction horizon at a sampling time of 40 $\mu$$\mathrm{s}\mathrm{a}$nd 80 $\mu$$\mathrm{s}$ It is observed that the undershoot in the current dq axis becomes more pronounced with higher sampling times as shown in Figure 10. At 80 $\mu$$\mathrm{s}\mathrm{t}$here is a significant undershoot, with the THD in the current abc axis being 23.13% and an average switching frequency of 1.9 $\mathrm{KHz}$ However, at a sampling time of 40 $\mu$$\mathrm{s}\mathrm{t}$he undershoot in the dq axis diminishes but is still present having a THD of 12.32% in the current abc axis and average switching frequency of 3.94 $\mathrm{KHz}$

To clearly illustrate the impact of load transients on Total Harmonic Distortion (THD) and switching frequency, in Figure 10(c) black arrows have been used to denote abrupt torque changes. Specifically, a positive torque step of +10 was introduced at 0.04 seconds, and a subsequent negative torque step of -8 was introduced at 0.1 seconds. These deliberate disturbances simulate realistic load variations that may occur in aerospace motor operations, yet the methodology remains applicable to a broad range of high-performance drive systems. By emphasizing these transient events, we can effectively evaluate how the control algorithm responds under load changes.

At lower sampling times 40 $\mu$$\mathrm{s}\mathrm{t}$he undershoot diminishes until it was eliminated at around 10 $\mu$$\mathrm{s}\mathrm{a}$s shown in Figure 11 with a THD of 6.48%. however with an average switching frequency of 16.1 $\mathrm{KHz}$ A three step and 5 step prediction horizon was also simulated for the same parameters, the 3 step prediction had a THD of 6.28% at average switching frequency of 17.9 $\mathrm{Khz},$ while the 5-step SDA had a THD of 6.3% at a average switching frequency of 17.73 $\mathrm{Khz},$ showing that at maximum switching frequency the effect of long prediction horizon is not beneficial in terms of THD reduction. However, the average switching frequency can be controlled when using Equation (10), by tuning the weighing factor ${\lambda }_u$.

3.2. Effect of Weighing Factor ${\lambda }_u$ on THD and Average Switching Frequency $F_{sw}$

In this approach, the value of the weighting factor $\lambda$ was varied from 0.1 to 20 in increments of 0.1, as illustrated in Figure 12. The figure highlights the relationship between the SDA and the FCS-MPC (Equation (16)) across different $\lambda$ values. A more detailed view of this trend—specifically involving the Total Harmonic Distortion (THD), average switching frequency, and weighting factor—is also presented in Figure 12. The results indicate that increasing the weighting factor reduces the switching frequency but leads to higher THD.

In the comparison shown in Figure 12, each marker point represents a separate simulation run using identical parameters. The simulation corresponding to the highest $\lambda$ value from Figure 12 is further illustrated in Figure 13. In this case, the FCS-MPC shows noticeable instability in the $dq$ axis and a high THD in the $abc$ axis. In contrast, the SDA demonstrates a slight undershoot, with an $abc$ axis THD of 19.8% and an average switching frequency of 3.82 $\mathrm{kHz}$

3.3. Sphere Decoding Algorithm for $N_p>2$

To show the benefit of the long prediction horizon a three step and a 5 step prediction horizon SDA was implemented for the same parameter as shown in Figure 14, there was a decrease in THD from the 2 to 5 step, at the third step the THD was 16.27 % with an average switching frequency of 5.18 $\mathrm{KHz},$ at a prediction horizon of 5, the THD in the abc axis reduced to 14.49 % with an average switching frequency of 5.79 $\mathrm{KHz}$ In this case increase in the prediction horizon leads to higher switching frequency and lesser THD. To fully assess the impact of the long prediction horizon on system performance, numerous simulations were conducted using the sphere decoding algorithm (SDA) with prediction horizons of 2, 3, and 5 steps. These simulations were carried out for different sampling times from 10 - 80 $\mu$$\mathrm{s}($8 sampled values).

By extending the prediction horizon, the controller is provided with enhanced foresight over future system states, which can facilitate more effective optimization of switching decisions. This improvement in anticipation is expected to yield a reduction in THD by enabling the controller to better manage and mitigate harmonic distortions in the output waveform. Figure 15 illustrates the comparative results for the 2-, 3-, and 5-step SDA configurations. The red arrow shows the increment of the weighing factor while the blue arrow shows the increment in sampling time . The simulation outcomes demonstrate that as the prediction horizon is extended, there is a notable decrease in THD, highlighting the benefits of a longer horizon. At the same time, variations in the switching frequency are observed, reflecting the inherent trade-off between achieving lower Total harmonic distortion and managing the average switching frequency.

Figure 16 presents results from a simulation using a randomly selected weighting factor of 10. For the two-step prediction, the THD was 12.98% with an average switching frequency of 6.37 $\mathrm{kHz}$ In the three-step prediction, the THD increased to 13.36%, while the average switching frequency decreased to 5.96 $\mathrm{kHz}$ For the five-step prediction, the THD dropped to 10.64% with a higher average switching frequency of 8.69 $\mathrm{kHz}$ Notably, increasing the prediction horizon from two to three steps resulted in higher THD but a reduction in switching frequency.This highlights the complex and nonlinear relationship between the prediction horizon, switching frequency, and total harmonic distortion (THD). The arrow in the upper three sub figures highlights the decrease trend of the $I_q$ current ripple.

3.4. Computational Time $N_p>2$

In this section the computational time of the FCS and SDA were discussed. For the results the simulation was run 10 times for each algorithm and the minimum and maximum value is taken as the range of the algorithm.

Table 3. Computational Time ($\mu \mathrm{s}$)

Optimization Algorithm	$N_p=1$	$N_p=2$	$N_p=3$	$N_p=5$
FCS-MPC	5.63 - 16	-	-	-
SDA	4.5 - 11	46 - 69	50 - 57	55-57

Table 3. Computational Time ($\mu \mathrm{s}$)

Optimization Algorithm	$N_p=1$	$N_p=2$	$N_p=3$	$N_p=5$
FCS-MPC	5.63 - 16	-	-	-
SDA	4.5 - 11	46 - 69	50 - 57	55-57

The computational time analysis shows that for a short prediction horizon ($N_p=1$), the Sphere Decoding Algorithm (SDA) (4.5–11 $\mu$$\mathrm{s})$ is more efficient than Finite Control Set Model Predictive Control (FCS-MPC) (5.63–16 $\mu$$\mathrm{s})$. As shown in Figure 5, this is because SDA does not need to evaluate every possible control sequence, thereby reducing computational burden. In contrast, SDA exhibits a significant increase in computational time from $N_p=1$ to $N_p=2$ (46–69 $\mu$$\mathrm{s})$. However, beyond this point, the increase stabilizes, with computational time ranging between 50–57 $\mu$$\mathrm{s}\mathrm{f}$or $N_p=3$ and 55–57 $\mu$$\mathrm{s}\mathrm{f}$or $N_p=5$. This trend suggests that while SDA initially experiences a sharp rise in computational cost, it becomes more predictable for higher horizons. Since FCS-MPC was not run beyond $N_p=1$, its computational scalability remains uncertain. However, SDA remains feasible even for larger horizons, making it a promising option for computationally constrained applications where extended predictive control is required.

3.5. Effect of Parameter Mismatch on Long Prediction Horizon MPC

In this section, we analyze the impact of parameter mismatches on the performance of the Sphere Decoding Algorithm (SDA) applied to a Permanent Magnet Synchronous Motor (PMSM) drive, paying particular attention to the effect of variations in resistance, inductance, and flux linkage on the system’s Total Harmonic Distortion (THD), switching frequency, and reference current tracking.

3.5.1. Resistance Mismatch

The simulation results demonstrate that variations in the resistance parameter within a reasonable range (e.g., 0.2, 0.5, 5, and 10 times the nominal resistance $R_s$) have a negligible impact on the Total Harmonic Distortion (THD) and the switching frequency of the SDA-based FCS-MPC for PMSM drives. As shown in Table 4, regardless of whether the resistance is increased (5×$R_s$ or 10×$R_s$) or decreased (0.2×$R_s$), the THD remains relatively constant, and the corresponding switching frequencies vary slightly.

This minimal influence of resistance mismatch can be primarily attributed to the inherent nature of resistance in PMSM drive systems. Resistance predominantly affects the steady-state voltage drops and the ohmic losses, rather than the dynamic response of the system. Since the SDA-based FCS-MPC controller is designed to dynamically adjust control inputs based on real-time current tracking, it effectively compensates for any deviations in resistance during each control cycle. Consequently, the current waveforms remain consistent, and the phase currents in the abc frame do not exhibit significant distortions.

Furthermore, the structure of the cost function in the FCS-MPC controller primarily optimizes the voltage application to minimize the current error. Since resistance impacts mainly the steady-state power loss and does not alter the system’s dynamic characteristics (such as pole locations), its variations do not lead to substantial differences in the controller’s performance.

3.5.2. Flux Mismatch

Flux linkage mismatch has a relatively minor impact on Total Harmonic Distortion (THD) and switching frequency when compared to inductance mismatch. The THD values remain fairly stable across different flux mismatch conditions, and the switching frequency exhibits little variation as shown in Table 5. However, flux mismatch introduces noticeable variations in the reference values of the dq-axis current components as shown in Figure 17. This behavior can be attributed to the direct influence of flux linkage on the back-electromotive force in the PMSM, which in turn affects the generated torque. When the flux linkage is incorrectly estimated, the predictive model used within the SDA-based FCS-MPC either underestimates or overestimates the required voltage commands. Consequently, the computed current references deviate from their ideal values, which may manifest as apparent overshoot or undershoot in the dq current trajectories. Importantly, this perceived overshoot is not indicative of instability; rather, it reflects a shift in the reference due to the flux estimation error. The controller is still capable of tracking these altered references accurately, as the overall dynamic performance of the system remains stable.

Moreover, since flux linkage does not directly affect the motor’s inductive reactance, it does not significantly alter the rate of change of the currents or the system’s time constants. This explains why the switching frequency remains largely unaffected by flux mismatch. The cost function in the FCS-MPC framework focuses on minimizing the error between predicted and actual currents, and any shift in the dq references caused by flux mismatch is dynamically corrected by the controller.

In summary, while flux linkage mismatch may lead to observable shifts in the dq-axis reference values, these shifts do not compromise the stability or overall performance of the control system. The robustness of the SDA-based FCS-MPC against flux mismatch underscores that, within practical limits, the control strategy is more sensitive to dynamic parameters (such as inductance) than to static parameters like flux linkage.

3.5.3. Inductance Mismatch

Inductance mismatch is the most influential parameter affecting the performance of the SDA-based FCS-MPC in PMSM drives, as demonstrated in Table 6. When the inductance is drastically underestimated (e.g., $0.1\times L_{dq}$), the controller exhibits very THD values, around 29–30%, across all prediction horizons with a correspondingly low switching frequency of approximately 1.6 $\mathrm{KHz}$ This behavior indicates that an underestimation of inductance drastically reduces the motor’s electrical time constant, $\tau ={\displaystyle \frac{L}{R}}$, causing the current to change too abruptly in response to voltage variations. Consequently, the voltage vectors selected by the controller fail to correct these rapid changes, resulting in a distorted current waveform with high THD and minimal corrective switching.

Conversely, when the inductance is moderately underestimated (e.g., $0.5\times L_{dq}$), the performance improves significantly. For a 2-step prediction, the THD drops to 12.85% with an increased switching frequency of about 8.4 $\mathrm{KHz},$ and further improvements are observed in the 3-step and 5-step predictions, with THD values decreasing to 11.20% and 9.57%, respectively. However, this improvement comes at the cost of an increased switching frequency, reflecting a more aggressive corrective action as the controller compensates for the residual mismatch. This compensation is indicative of the controller’s attempt to maintain current tracking despite an imperfect model.

When the inductance is overestimated (e.g., $2\times L_{dq}$), the effects vary with the prediction horizon. For a 2-step prediction, the THD increases to 18.37% with a lower switching frequency of 4.41 $\mathrm{KHz},$ implying that the controller is under-compensating for the slower response induced by the higher inductance. In contrast, for the 3-step and 5-step predictions, the THD improves to 12.48% and 13.20%, with switching frequencies around 7.10 $\mathrm{KHz},$ indicating that longer prediction horizons can partially mitigate the adverse effects of an overestimated inductance by providing better foresight into system dynamics.

In the case of extreme overestimation (e.g., $5\times L_{dq}$), performance deteriorates markedly for shorter prediction horizons; a 2-step prediction yields a very high THD of 39.60% with a low switching frequency of 3.11 $\mathrm{KHz},$ while a 3-step prediction shows a THD of 22.22% at 4.33 $\mathrm{KHz}$ Interestingly, the 5-step prediction in this scenario achieves a much improved THD of 12.09% with a switching frequency of 7.38 $\mathrm{KHz},$ suggesting that extended prediction horizons can partially compensate for the mismatch by aligning control actions with the slower system response.

From a theoretical perspective, these observations are well explained by control theory. The electrical time constant, $\tau ={\displaystyle \frac{L}{R}}$, is critical in determining the system’s dynamic response. A lower L results in a shorter $\tau$, causing faster but more volatile current changes, while a higher L increases $\tau$, leading to a sluggish response that the controller attempts to counterbalance with increased switching. In the closed-loop transfer function of the current control loop, the pole positions are directly influenced by L; significant deviations in L alter these poles, resulting in either overly aggressive or insufficient corrective actions. Therefore, the controller’s performance is highly sensitive to the accuracy of the inductance parameter.

4. Conclusions

In this paper, we investigated the performance of a Sphere Decoding Algorithm (SDA)-based Finite Control Set Model Predictive Control (FCS-MPC) for PMSM drives, with a particular focus on the effects of sampling time, control effort weighting, and parameter mismatches. Our simulation results demonstrate that at longer sampling times (40 and 80 $\mu$$\mathrm{s})$, significant undershoot occurs in the current response. As the sampling time is reduced below 40 $\mu$$\mathrm{s},$ the undershoot gradually diminishes and disappears entirely at a 10 $\mu$$\mathrm{s}\mathrm{s}$ampling time. At this shorter sampling time, while the inclusion of a control effort term successfully reduces the switching frequency, a small undershoot remains, indicating that the predicted output closely follows the reference but does not perfectly align with the midpoint of the desired waveform. Furthermore, our investigation into parameter mismatches revealed that variations in resistance and flux linkage have a negligible effect on the system’s THD and switching frequency. In contrast, inductance mismatches significantly impact the controller’s performance. Depending on the magnitude of the inductance error, the system sometimes compensates by increasing the switching frequency; however, for severe mismatches, the system can become unstable. This outcome underscores the critical importance of accurate inductance estimation in ensuring effective current tracking and maintaining low THD.

From an aerospace perspective, achieving low torque ripple and high efficiency is paramount for PMSM applications in electric propulsion and actuation systems. The findings suggest that optimizing the prediction horizon and control effort weighting can significantly improve dynamic performance, reducing electromagnetic interference and enhancing system reliability in aircraft environments where precise control is essential. Looking ahead, future research should focus on enhancing the robustness of the SDA-based FCS-MPC approach in the presence of parameter mismatches. One promising direction is the integration of Recursive Least Squares (RLS) for real-time parameter estimation. RLS can continuously update key parameters such as resistance, inductance, and flux linkage, mitigating the adverse effects of mismatches and enabling the controller to maintain optimal performance despite variations due to temperature, aging, or manufacturing tolerances. Additionally, exploring hybrid adaptive estimation methods, such as combining RLS with an Extended Kalman Filter, may further improve noise rejection and estimation accuracy. Experimental validation on a physical PMSM setup under various operating conditions, including thermal and load variations, will be sought in future studies to confirm the practical viability of the proposed solutions in aerospace applications.