An Advanced Angular Velocity Error Prediction Horizon Self-Tuning Nonlinear Model Predictive Speed Control Strategy for PMSM System

In nonlinear model predictive control (NMPC), higher accuracy can be obtained with a shorter prediction horizon in steady-state, better dynamics can be obtained with a longer prediction horizon in a transient state, and calculation burden is proportional to the prediction horizon which is usually pre-selected as a constant according to dynamics of the system with NMPC. The minimum calculation and prediction accuracy are hard to ensure for all operating states. This can be improved by an online changing prediction horizon. A nonlinear model predictive speed control (NMPSC) with advanced angular velocity error (AAVE) prediction horizon self-tuning method has been proposed in which the prediction horizon is improved as a discrete-time integer variable and can be adjusted during each sampling period. A permanent magnet synchronous motor (PMSM) rotor position control system with the proposed strategy is accomplished. Tracking performances including rotor position Integral of Time-weighted Absolute value of the Error (ITAE), the maximal delay time, and static error are improved about 15.033%, 23.077%, and 10.294% respectively comparing with the conventional NMPSC strategy with a certain prediction horizon. Better disturbance resisting performance, lower weighting factor sensitivities, and higher servo stiffness are achieved. Simulation and experimental results are given to demonstrate the effectiveness and correctness.


Introduction
Permanent magnet synchronous motor (PMSM) is frequently used in industry due to high level of power density, efficiency, and torque-ampere ratio with low weight and volume [1].The proportional-integral (PI) controller is usually used in the PMSM control system to control motor speed and stator currents to obtain desired dynamic and performance.Many improved control strategies are presented to satisfy different criteria, such as 2-degree-of-freedom (2-DoF) control, internal model control (IMC), active disturbance rejection control (ADRC), and model predictive control (MPC).Two controllers are adopted in the 2-DoF control to adjust and decouple dynamic and disturbance resisting performance [2].A model of the controlled plant is added into the structure of IMC to make the system operate in an open-loop manner to resist disturbance [3].An extended state observer (ESO) is introduced into the ADRC strategy, and all disturbances are estimated and compensated [4].The above control strategies can meet the basic requirements of the PMSM system.MPC strategy has been highlighted in recent years due to better dynamics, easier understanding, and handling, and has been widely applied in motor driving and power electronic realm [5].MPC can be divided as continuous-control-set MPC (CCS-MPC) and finite-control-set MPC (FCS-MPC) [6].Since modulation is adopted in the CCS-MPC, fewer ripple and switching losses can be obtained compared with FCS-MPC with variable switching frequency [7].But the prediction horizon in CCS-MPC is hard to be adjusted due to the limitation of the voltage function generation process.MPC belongs to a short prediction horizon method covering only one sampling period.The MPC in the PMSM system can be divided into model predictive current control (MPCC), model predictive speed control (MPSC), and model predictive torque control (MPTC) according to the primary control objective [8][9][10].A voltage smoother in [11] and an integral sliding mode disturbance observer in [12] are combined with the MPSC respectively to enhance control performances of speed or angular velocity.
Aiming at high stability and tracking performances, some improvements such as Sequential MPC (SMPC), Parallel MPC (PMPC), and NMPC are proposed to prolong the prediction horizon to several sampling periods [13,14], and the prediction horizon of Nonlinear predictive control (NPC) is extended to thousands of sampling periods which are within the Shannon upper bound [15].These are called long prediction horizon methods.
NPC is a generalized nonlinear method linearized by Lie derivation.It is usually combined with an advanced algorithm to build a nonlinear object model and to design the controller [16].
In SMPC, the cost function is split into several parts according to the multiple controlled objectives or variables [17].In PMPC, a switching vector intersection method is inserted into the MPTC structure to select the optimal vector to obtain better robustness and dynamics [18].A compensatory block of the three-phase inverter is inserted before MPC processes, the input state variables of SMPC and PMPC are at time k + 1 [19].As a result, the predictive variables in the discrete-time model are at time k + 2 [20].Due to the divided cost functions, the calculation time is increased because of the calculation of the cost function values and the selection of vectors in each step, and the number of controlled objectives is restricted by the number of vectors of the switching table.This can be solved by adding a group of extra vectors [21].
SMPC and PMPC cover two sampling periods, which improves system performances finitely.An accumulating type nonlinear cost function is presented and used in NMPC to prolong the prediction horizon [22].Some extra weighting factors are added to the cost function to enhance convergence of the system [23], and an improved extended prediction horizon method for NMPC is proposed under the condition of limited calculation burdens [24].The stability of the extended prediction horizon is verified in theory [25].A kind of adaptive dynamic programming (ADP) algorithm is combined with the NMPC to improve the robustness and dynamics [26].NMPC is used in many applications of power systems and power electronic realms.The extended prediction horizons in the traditional methods are fixed.It does not benefit all operating states, for example, the calculation burden and prediction error are increased heavily in steady-state, and hence the performances are also affected.An extended prediction self-adaptive control (EPSAC) [27] and Nonlinear EPSAC (NEPSAC) [28] to NMPC strategy are proposed in [29] to auto-tune the sampling period and the maximum prediction horizon to enhance robustness and dynamics of different systems.
MPSC strategy and NMPC structure are combined as so-called NMPSC.An NMPSC with AAVE prediction horizon self-tuning method is proposed in this paper and applied to the PMSM rotor position system.The prediction horizon is improved as a discrete-time integral variable which can be adjusted in every sampling period to improve tracking accuracy.Performance comparisons including rotor position error, servo stiffness, and weighting factor sensitivities are conducted between the proposed method and the conventional NMPSC without any prediction horizon self-tuning method.According to the simulation and experimental results, servo stiffness, weighting factor sensitivities, and ITAEs are improved effectively.

Problem Description
Stability and rapidity are the main indicators of the PMSM control system.High accuracy and low operation burden of the system are also required during the steadystate [13].
Time delay has to be considered for the selected switching signal, which results that the output cannot be output in time in the digital control system [10].NMPC predicts the variables at time k + j based on the sampled signals and selected prediction horizon by repetitive predictive processes.Figure 1 shows the time flows of two NMPC strategies.The prediction calculation time of the conventional NMPC is constant in solid line as shown in Figure 1a and it is adjustable in variable prediction horizon NMPC in dash-line as shown in Figure 1b.The delay time is not fixed in the variable prediction horizon NMPC.Prediction processes and cost function calculation processes are two main parts of the prediction time slot.Since the number of voltage vectors is pre-selected, the prediction time is changed with the number of prediction processes which equals the prediction horizon.

Problem Description
Stability and rapidity are the main indicators of the PMSM control system accuracy and low operation burden of the system are also required durin steady-state [13].
Time delay has to be considered for the selected switching signal, which resul the output cannot be output in time in the digital control system [10].NMPC pr the variables at time k + j based on the sampled signals and selected prediction h by repetitive predictive processes.Figure 1 shows the time flows of two NMPC gies.The prediction calculation time of the conventional NMPC is constant in sol as shown in Figure 1a and it is adjustable in variable prediction horizon NM dash-line as shown in Figure 1b.The delay time is not fixed in the variable pred horizon NMPC.Prediction processes and cost function calculation processes ar main parts of the prediction time slot.Since the number of voltage vect pre-selected, the prediction time is changed with the number of prediction pro which equals the prediction horizon.The control performance of NMPC is related to the prediction horizon.For th sient state, a longer prediction horizon is usually required to achieve better dy performances to resist disturbance.This induces a heavier calculation burden and ers accuracy due to longer calculation time and more neglected sampling points [2 for the steady-state, a shorter prediction horizon with better accuracy and a light culation burden is required to enhance steady performances and to decrease static Therefore, the prediction horizon should be adjusted for different operating states alize optimal operating performance.

Discrete-Time Model
The state functions of the PMSM model can be expressed based on the two synchronized rotation reference frame.The functions including the stator voltage v linkage ψs and electromagnetic torque Te for the Surface PMSM (SPMSM) are exp as: The control performance of NMPC is related to the prediction horizon.For the transient state, a longer prediction horizon is usually required to achieve better dynamic performances to resist disturbance.This induces a heavier calculation burden and lowers accuracy due to longer calculation time and more neglected sampling points [24].As for the steady-state, a shorter prediction horizon with better accuracy and a lighter calculation burden is required to enhance steady performances and to decrease static error.Therefore, the prediction horizon should be adjusted for different operating states to realize optimal operating performance.

Discrete-Time Model
The state functions of the PMSM model can be expressed based on the two-phase synchronized rotation reference frame.The functions including the stator voltage v s , flux linkage ψ s and electromagnetic torque T e for the Surface PMSM (SPMSM) are expressed as: and the electrical angular velocity ω r can be described by the equation as: where subscripts d and q mean that the components of variables locate in the d-axis and q-axis respectively, such as i sd and i sq are stator current components on the d-axis and q-axis respectively.Moreover, L s is stator self-inductance which is satisfying L s = L sq = L sd for SPMSM, T L is load torque, R s is stator resistance, p is the number of pole pairs, J is rotor inertia, B is friction coefficient and ψ m is rotor flux magnitude.
Combining Equations ( 1) to (4), the space state equation of PMSM can be expressed as: where The two-step Euler interpolation method including predicting step x p (k + 1) and correcting step x(k + 1) is used to predict and discretize the model based on a sampling period T s , and the expressions are [5]:

NMPSC Strategy
The typical cost function for NMPC is [29]: where N 1 and N 2 are the minimum and the maximum prediction horizons which are constant values.The basic idea of the method is to calculate a sequence of future control signals to minimize a cost function defined over a prediction horizon.For the PMSM control system, N 1 can be selected as 1 which meets the prediction horizon in MPSC, and N 2 can be selected according to the Shannon sampling principle and calculation burden of the processor.Speed reference tracking and torque ampere optimization are two main objectives in the PMSM system, therefore, the cost function for NMPSC can be determined as [5]: where λ ω and λ i are two weighting factors reflecting the importance of the objectives.The reference tracking is realized by the first term of the cost function.i sq is related to ω r according to the discrete-time model in Equations ( 4) and ( 5), as a result, it can be controlled by the first term at the same time.The torque ampere optimization objective is realized by the second term of the cost function.
To enhance convergence of the system, a group of weighting factors Q j converging to zero are added into the cost function, Equation (11), as: and the weighting factor Q j has a common tuning function:

Proposed AAVE Method
Redefining the prediction horizon N 2 in Equation ( 12) as a positive discrete-time integral variable P k+1 , the cost function is improved as: where ω r * is angular velocity reference.
There are two fundamental objectives: (1) to ensure that the angular velocity reference ω r * is tracked, and (2) to guarantee that a virtual reference ω rv * is generated and to modify the prediction horizon to the optimal value automatically.The flowchart of the proposed method is shown in Figure 2, where i and j are two cyclic variables, i is used to make sure that all voltage vectors are fully considered, and j is used to judge whether the prediction horizon accumulating processes are accomplished.
trolled by the first term at the same time.The torque ampere optimizati realized by the second term of the cost function.
To enhance convergence of the system, a group of weighting factors to zero are added into the cost function, Equation (11), as: and the weighting factor Qj has a common tuning function:

Proposed AAVE Method
Redefining the prediction horizon N2 in Equation ( 12) as a positive d tegral variable Pk+1, the cost function is improved as: where ωr* is angular velocity reference.
There are two fundamental objectives: (1) to ensure that the angular ence ωr* is tracked, and (2) to guarantee that a virtual reference ωrv* is ge modify the prediction horizon to the optimal value automatically.The fl proposed method is shown in Figure 2, where i and j are two cyclic variab make sure that all voltage vectors are fully considered, and j is used to the prediction horizon accumulating processes are accomplished.

AAVE prediction horizon
self-tuning method

Measurement and detection;
The angular velocity ω r is measured by an encoder, and the prediction horizon P k at the last sampling period is stored in the register for the calculating process at the k-th sampling period.

2.
Virtual reference and step generation; Once the current angular velocity ω r is different from that of the previous sampling period, the virtual reference ω rv needs to be redefined automatically during each sampling period so that it can converge to a new optimal prediction horizon when there is an abrupt change in the system.
The virtual reference ω rv and the step increment ∆P k+1 of the prediction horizon at k + 1-th sampling period can be calculated as: and and K 1 , K 2, and M are coefficients to ensure the same order of magnitudes of ω rv as that of ω r and the limited change of the prediction horizon, and a constraint as 2K 2 < M < ∞ should be satisfied to obtain a determined sign of the virtual reference.K 1 , K 2 , and M are selected as 0.05, 1, and 2 respectively in the strategy.The functions v and u are shown in Figure 3.The virtual reference ω rv and prediction horizon P k+1 have positive slopes when the angular velocity error e k at the k-th sampling period increases and vice versa.

Measurement and detection;
The angular velocity ωr is measured by an encoder, and the prediction horizon Pk at the last sampling period is stored in the register for the calculating process at the k-th sampling period.

Virtual reference and step generation;
Once the current angular velocity ωr is different from that of the previous sampling period, the virtual reference ωrv needs to be redefined automatically during each sampling period so that it can converge to a new optimal prediction horizon when there is an abrupt change in the system.
The virtual reference ωrv and the step increment ΔPk+1 of the prediction horizon at k + 1-th sampling period can be calculated as: and   where and K1, K2, and M are coefficients to ensure the same order of magnitudes of ωrv as that of ωr and the limited change of the prediction horizon, and a constraint as 2K2 < M < ∞ should be satisfied to obtain a determined sign of the virtual reference.K1, K2, and M are selected as 0.05, 1, and 2 respectively in the strategy.The functions v and u are shown in Figure 3.The virtual reference ωrv and prediction horizon Pk+1 have positive slopes when the angular velocity error ek at the k-th sampling period increases and vice versa.According to Equation (11), Pk+1 should be integral betwe in (21) should be handled as follow: where the function ceil (.) is a value round-up function that m tion horizon is an integral value as shown in Figure 4.

1.
Prediction horizon calculation and range judgment.
The prediction horizon at k + 1-th sampling period P k+1 can be calculated as follow: According to Equation (11), P k+1 should be integral between 1 and N 2 , and the result in (21) should be handled as follow: where the function ceil (.) is a value round-up function that makes sure that the prediction horizon is an integral value as shown in Figure 4.
The prediction horizon at k + 1-th sampling period Pk+1 can be calculated as follow: According to Equation (11), Pk+1 should be integral between 1 and N2, and the result in (21) should be handled as follow: where the function ceil (.) is a value round-up function that makes sure that the prediction horizon is an integral value as shown in Figure 4.The resulted Pk+1 from Equation ( 22) is brought into Equation ( 14) to form the cost function.After all of the voltage vectors are evaluated in the cost function, the voltage vector with the minimum cost function value can be selected as the optimal vector.
According to the conventional NMPSC cost function in Equation ( 12), usually N1 = 1, N2 is a fixed value that satisfies the stability condition in [30] and maintains constant during the operation of the system.For the proposed method, N2 is also selected from the range to ensure stability, but the prediction horizon Pk+1 is adjusted within the range to achieve optimum performances in every sampling period.

Simulation Results and Analysis
The structure diagram of NMPSC with AAVE prediction horizon self-tuning method on a PMSM rotor position system is shown in Figure 5.A simulation environment has been built in the MATLAB/Simulink software according to the structure.The weighting factors λω and λi are selected by the branch and bound algorithm as 1 and 0.5 The resulted P k+1 from Equation ( 22) is brought into Equation ( 14) to form the cost function.After all of the voltage vectors are evaluated in the cost function, the voltage vector with the minimum cost function value can be selected as the optimal vector.
According to the conventional NMPSC cost function in Equation ( 12), usually N 1 = 1, N 2 is a fixed value that satisfies the stability condition in [30] and maintains constant during the operation of the system.For the proposed method, N 2 is also selected from the range to ensure stability, but the prediction horizon P k+1 is adjusted within the range to achieve optimum performances in every sampling period.

Simulation Results and Analysis
The structure diagram of NMPSC with AAVE prediction horizon self-tuning method on a PMSM rotor position system is shown in Figure 5.A simulation environment has been built in the MATLAB/Simulink software according to the structure.The weighting factors λ ω and λ i are selected by the branch and bound algorithm as 1 and 0.5 respectively [5].The sampling frequency f s and the maximum prediction horizon N 2 are selected as 10 kHz and 100.Some important parameters of PMSM are listed in Table 1.
Electronics 2021, 10, x FOR PEER REVIEW 8 of 20 respectively [5].The sampling frequency fs and the maximum prediction horizon N2 are selected as 10 kHz and 100.Some important parameters of PMSM are listed in Table 1.
Position Controller AAVE Prediction Horizon Self-tuning Method

Tracking Performance Analysis
The rotor position reference θ* is given by a sinusoidal signal with amplitude 10, 5 Hz, and zero initial phases.The reference, rotor position, and error waveforms of the NMPSC with AAVE prediction horizon self-tuning method and the conventional NMPSC with constant prediction horizon value are shown in Figure 6.The left y-axis in Figure 6a reflects rotor position and rotor position waves are all sinusoidal, and the right y-axis in Figure 6a reflects rotor position error and the rotor position error waves have doubled frequency comparing with rotor position waves.The part enlarged waves including both rotor positions and errors are shown in Figure 6b,c.Figure 6b,c shows that the rotor position in the NMPSC with AAVE prediction horizon self-tuning method can track the reference with a delay time of about 4ms, and its maximum static error is about 1rad.As for the conventional NMPSC strategy, they are 4.2 ms and 1.2 rad.During normal operation, the maximal error occurs near the crossing zero points.
The prediction horizon waves at starting operation stage are shown in Figure 6d.Comparing with the fixed prediction horizon in the conventional strategy, the P k+1 of the proposed method is changed with the errors.Around 22 ms, the prediction horizons are converged to 1 when the actual rotor position signal tracks the reference successfully.Because the delay errors, which exist in measurement, calculation, transmission, and action processes of the digital system, are neglected in the simulation environment, the ideal errors are not large enough to trigger the AAVE prediction horizon self-tuning method during normal tracking processes.
The converter efficiency is also tested in simulations.It is about 0.9918 for the proposed method and 0.98 for conventional strategies.The increment is about 1.204%.
ITAE is a comprehensive index defining the absolute value of accumulating error during the whole operating process which can be expressed as: The rotor position ITAE, speed ITAE of the NMPSC with the proposed method, and the conventional NMPSC within 1.5s are listed in Table 2.The lowest ITAE values are obtained for the proposed method due to decreased errors at the start operation state mainly.The prediction horizon waves at starting operation stage are shown in Figure 6d.Comparing with the fixed prediction horizon in the conventional strategy, the Pk+1 of the proposed method is changed with the errors.Around 22 ms, the prediction horizons are converged to 1 when the actual rotor position signal tracks the reference successfully.Because the delay errors, which exist in measurement, calculation, transmission, and action processes of the digital system, are neglected in the simulation environment, the A 3 N•m load torque is added at 0.3 s, the torque currents (q-ax resist load torque, and the current increasing process of the proposed shorter than the conventional NMPSC.
Rotor position and speed ITAEs within 0.8 s are listed in Table be decreased by about 3.53% for the proposed strategy and increa when the speed direction is negative.

Weighting Factor Sensitivities
Weighting factors are inserted into the cost function when the cludes several control targets, variables or constraints.The branch an is a frequently used weighting factor tuning method based on lots of perimental results.To assure control performance with inaccurate w system should not be sensitive to the weighting factors and that is an make the system operate in a suitable state.
There are two weighting factors in ( 14) and the relationships A 3 N•m load torque is added at 0.3 s, the torque currents (q-axis) are increased to resist load torque, and the current increasing process of the proposed strategy is slightly shorter than the conventional NMPSC.
Rotor position and speed ITAEs within 0.8 s are listed in Table 3, ITAE values can be decreased by about 3.53% for the proposed strategy and increased by about 1.88% when the speed direction is negative.

Weighting Factor Sensitivities
Weighting factors are inserted into the cost function when the control strategy includes several control targets, variables or constraints.The branch and bound algorithm is a frequently used weighting factor tuning method based on lots of simulation and experimental results.To assure control performance with inaccurate weighting factors, the system should not be sensitive to the weighting factors and that is an important index to make the system operate in a suitable state.
There are two weighting factors in (14) and the relationships between ITAEs and different weighting factors are shown in Figure 8. ITAEs are decreasing with the increasing of λ ω and decreasing of λ i .The rotor position ITAE variations are 35.876%and 63.456% for the NMPSC with proposed method and the conventional NMPSC method respectively, and the average rotor position ITAE values are 0.0902 and 0.1112 respectively, the increment is about 18.885%.The smoother waves of the proposed strategy imply lower sensitivities of the weighting factors.

Propotional Parameter Selection
A proportional controller is adopted as the rotor position controller to resist overshoot and obtain better dynamics [31].The waveforms for rotor position ITAEs and maximal delay time with different proportional parameters are shown in Figure 9, the rotor position ITAEs and maximal delay times have a similar tendency and the waveforms decrease with increasing of proportional coefficient values.As shown in the figures, when the proportional coefficient takes a value from 150 to 400 according to the branch and bound algorithm, the rotor ITAE and delay time are within appropriate ranges.
Electronics 2021, 10, x FOR PEER REVIEW 13 of 20 A proportional controller is adopted as the rotor position controller to resist overshoot and obtain better dynamics [31].The waveforms for rotor position ITAEs and maximal delay time with different proportional parameters are shown in Figure 9, the rotor position ITAEs and maximal delay times have a similar tendency and the waveforms decrease with increasing of proportional coefficient values.As shown in the figures, when the proportional coefficient takes a value from 150 to 400 according to the branch and bound algorithm, the rotor ITAE and delay time are within appropriate ranges.

Proportional Parameter Selection
Some parameters including stator resistance R s , stator induction L s and magnet flux linkage ψ m are always changing during operating states [4].Defining R se , L se, and ψ me as model parameters, different mismatches including R se /R s , L se /L s, and ψ me /ψ m are used to test stable operation ranges.The rotor position, error, and mismatch waves are shown in Figure 10.As shown in the figures, the stable ranges of parameter mismatches are 0.1~50, 0.1~1000 and 0.1~25 for R se /R s , L se /L s and ψ me /ψ m respectively.

Servo Stiffness Analysis
Servo stiffness is an important index for the rotor position control system.This index reflects the rotor rotating angle generated by the load torque T L which is uploaded on the shaft with a constant rotor reference.The definition is: The servo stiffness waveforms with different weighting factors are shown in Figure 11.The servo stiffness K increases as the load torque T L increases.The one with λ ω = 1 and λ i = 0.2 has the highest servo stiffness, on the contrast, the one with λ ω = 0.2 and λ i = 1 has the lowest servo stiffness, which means that the former generates a less rotor angle change than the latter when a load torque is uploaded on the shaft with the same constant reference.

Servo Stiffness Analysis
Servo stiffness is an important index for the rotor p dex reflects the rotor rotating angle generated by the lo on the shaft with a constant rotor reference.The definitio The servo stiffness waveforms with different weigh 11.The servo stiffness K increases as the load torque TL and λi = 0.2 has the highest servo stiffness, on the contra 1 has the lowest servo stiffness, which means that the for change than the latter when a load torque is uploaded stant reference.The average values and mean square error values of servo stiffness with different weighting factors for the proposed method and conventional method are listed in Table 4.Under the same weighting factors, the average value and the mean square error value are larger in the proposed method than that in the conventional method.The average values and mean square error values of servo stiffness with different weighting factors for the proposed method and conventional method are listed in Table 4.Under the same weighting factors, the average value and the mean square error value are larger in the proposed method than that in the conventional method.Moreover, the ratio of the speed weighting factor λ ω and the current weighting factor λ i is also an important index for servo stiffness according to the waveforms and data.The same λ ω /λ i ratio will generate the same servo stiffness, for example, when λ ω /λ i = 1 (λ ω = λ i = 1, 0.6, 0.2), these three waveforms are the same as shown in Figure 11 and Table 4.The servo stiffness K increases as the ratio increases until the servo stiffness reaches the maximal value when the ratio equals 5, and vice versa.

Tracking Comparison
The experimental setup is shown in Figure 12.A three-phase inverter with IGBTs (FGL35N120FTD) and a 1kW PMSM (INOVANCE ISMH2-10C30CD) with an incremental encoder (INOVANCE EI34H) are combined into the main circuit where the rated input voltage of the inverter is 800 V.The control strategy is accomplished by a DSP (TMS320F2812PGFA) and a CPLD (EPM240T100I-5N).The parameters of the PMSM in the experiment are the same as that of the simulation in Table 1, and the sampling period T s is 0.1 ms.Importantly, the maximum prediction horizon N 2 is selected as 4 because of the calculation speed and the storage volume limitation of the hardware.λ ω , λ i, and the proportional parameters are 1, 0.5, and 300 respectively.The prediction horizon waveform of the proposed method is shown in Figure 13, it is adjusted continuously rather than remaining unchanged value in the conventional strategy.The prediction horizon waveform of the proposed method is shown in Figure 13, it is adjusted continuously rather than remaining unchanged value in the conventional strategy.The prediction horizon waveform of the proposed method is shown in Figure 13, it is adjusted continuously rather than remaining unchanged value in the conventional strategy.Figure 15 provides the experimental results of conventional NMPSC strategy under the same conditions.It shows that longer delay time and larger static error are obtained than that in Figure 14, and a slight additional delay time, which is circled in the red dashed line shown in Figure 15b, is generated when the load torque uploads on the shaft.
Speed and filtered torque current (q-axis) waveforms for both strategies are shown in Figure 16a,b respectively.Both speed waves are operated in sinusoidal shape, and the amplitudes of speed and current waves are different.Since the rotor position can track the reference more closely by the proposed method, the speed continuously operates in a low-speed region comparing with the conventional NMPSC strategy, and the speed ripple is increased because detecting error cannot be eliminated on time due to the limited encoder accuracy, especially in the start and stop operating processes.Figure 15 provides the experimental results of conventional NMPSC strategy under the same conditions.It shows that longer delay time and larger static error are obtained than that in Figure 14, and a slight additional delay time, which is circled in the red dashed shown in Figure 15b, is generated when the load torque uploads on the shaft.Speed and filtered torque current (q-axis) waveforms for both strategies are shown in Figure 16a,b respectively.Both speed waves are operated in sinusoidal shape, and the amplitudes of speed and current waves are different.Since the rotor position can track the reference more closely by the proposed method, the speed continuously operates in Some related experimental results of the NMPSC with the proposed method and the conventional NMPSC are listed in Table 5.The rotor position ITAE, maximal static error, and delay time within three cycles of the proposed method have been improved about 15.033%, 10.294%, and 23.077% respectively.Under the condition of the same proportional coefficient, obtained speed reference is lower than the conventional method, as a result, the maximum speed is less than the one achieved by the conventional method.Some related experimental results of the NMPSC with the proposed method and the conventional NMPSC are listed in Table 5.The rotor position ITAE, maximal static error, and delay time within three cycles of the proposed method have been improved about 15.033%, 10.294%, and 23.077% respectively.Under the condition of the same proportional coefficient, obtained speed reference is lower than the conventional method, as a result, the maximum speed is less than the one achieved by the conventional method.

Computation Burden
In each sampling period, 8 vector candidates need to be evaluated in the cost function (14) of the proposed method.The circulate time of the conventional MPSC algorithm in [5] on the available experimental platform is about 25 µs, and it is about 65 µs for the conventional NMPSC algorithm above-mentioned in which the prediction horizon equals 3.For the NMPSC with a proposed method with 4 maximum prediction horizons, the maximal and average circulate time is about 85 µs and 64.415 µs, and the sampling interval is sufficient to execute without overruns.Due to the fact that the weighting factor optimization is performed offline, there is no additional computation burden during the execution of the algorithm [9], and the average circulates time of the proposed method is decreased about 0.9% compared with the conventional NMPSC strategy.
In a steady state, the calculation burden will be decreased because of the decreasing of the prediction horizon, and in a transient state, it is not increased even though the prediction horizon reaches 4 sometimes.

Conclusions
An NMPSC with AAVE prediction horizon self-tuning method with a discrete-time integral prediction horizon is proposed in this paper.The prediction horizon is adjusted online to the optimal value for different operating states.The PMSM angular velocity error is used to generate a virtual reference to adjust the prediction horizon for each sampling period.A comparative study is carried out between the proposed method and the conventional NMPSC with a constant prediction horizon.It shows that higher efficiency and servo stiffness, and lower weighting factor sensitivities, ITAEs, static errors, and delay times are obtained.Moreover, the experimental results show that high tracking performances are achieved for all operating processes, which are not demonstrated in an ideal simulation environment without delays of digital system.
Comparing with the conventional NMPSC, the rotor position ITAE of the NMPSC with the proposed method in three sinusoidal cycles is decreased about 15.033%, and the maximal static error and delay time is decreased about 10.294% and 23.077% respectively in experiments.Besides, the average value of calculation burden is decreased about 0.9% by adjustable prediction horizon, and the sensitivities of weighting factors and servo stiffness are improved about 18.885% and 42.102% respectively.The AAVE prediction horizon self-tuning method shows a better performance comparing with the conventional method.
Based on the proposed method, the future research work is to combine deterministic artificial intelligence (D.A.I.) with the prediction horizon tuning method to further improve the online tuning process and control performance, and the experimental realization of D.A.I. and MPC needs to be considered.

Figure 2 .
Figure 2. Flowchart of NMPSC with the proposed VVAE method.

Figure 2 .
Figure 2. Flowchart of NMPSC with the proposed VVAE method.

Figure 3 .
Figure 3. Ideal functions u and v of the proposed AAVE method.(a) f 1. Prediction horizon calculation and range judgment.The prediction horizon at k + 1-th sampling period Pk+1 can

Figure 3 .
Figure 3. Ideal functions u and v of the proposed AAVE method.(a) function v; (b) function u.

Figure 3 .
Figure 3. Ideal functions u and v of the proposed AAVE method.(a) function v; (b) function u.

Figure 4 .
Figure 4. Ideal functions ceil of the proposed AAVE method.

Figure 4 .
Figure 4. Ideal functions ceil of the proposed AAVE method.

Figure 5 .
Figure 5. Structure of PMSM rotor position system with the proposed method.

Figure 5 .
Figure 5. Structure of PMSM rotor position system with the proposed method.

Figure 6 .
Figure 6.Rotor position, error, and prediction horizon waveforms of the proposed NMPSC with AAVE method and the conventional NMPSC.(a) Rotor positions and errors; (b) Part enlarged rotor positions; (c) Part enlarged errors; (d) Prediction horizon.

Figure 6 .Figure 7 .
Figure 6.Rotor position, error, and prediction horizon waveforms of the proposed NMPSC with AAVE method and the conventional NMPSC.(a) Rotor positions and errors; (b) Part enlarged rotor positions; (c) Part enlarged errors; (d) Prediction horizon.Two ramp signals with slope ±200 rad/s are adopted as rotor position references to compare performance in the same speed region.Simulation results are shown in Figure7.It is shown that rotor position waves of two strategies can track the references with similar delay times, and speeds are around ±95 rpm respectively.

Figure 7 .
Figure 7. Simulation results of ramp references of the proposed NMPSC with AAVE method and conventional NMPSC.(a) Rotor positions with a positive slope; (b) Speed and torque current (q-axis) with a positive slope; (c) Rotor position with a negative slope; (d) Speed and torque current (q-axis) with a negative slope.

Figure 8 .
Figure 8. Contour maps of rotor position ITAEs with different weighting factors for the conventional NMPSC and the NMPSC with the proposed method.(a) λω to ITAE for the NMPSC with AAVE method; (b) λω to ITAE for the conventional NMPSC; (c) λi to ITAE for the NMPSC with AAVE method; (d) λi to ITAE for the conventional NMPSC.

4. 3 .Figure 8 .
Figure 8. Contour maps of rotor position ITAEs with different weighting factors for the conventional NMPSC and the NMPSC with the proposed method.(a) λ ω to ITAE for the NMPSC with AAVE method; (b) λ ω to ITAE for the conventional NMPSC; (c) λ i to ITAE for the NMPSC with AAVE method; (d) λ i to ITAE for the conventional NMPSC.

4. 4 .
Proportional Parameter Selection operating states [4].Defining Rse, Lse, and ψme as model parameters, different mismatches including Rse/Rs, Lse/Ls, and ψme/ψm are used to test stable operation ranges.The rotor position, error, and mismatch waves are shown in Figure 10.As shown in the figures, the stable ranges of parameter mismatches are 0.1~5 0, 0.1~1000 and 0.1~25 for Rse/Rs, Lse/Ls and ψme/ψm respectively.

Figure 10 . 20 Figure 11 .
Figure 10.Parameter mismatches for the proposed method.(a) Rotor position and error with R se ; (b) different mismatch degrees for R se ; (c) Rotor position and error with L se ; (d) Different mismatch degrees for L se ; (e) Rotor position and error with ψ me ; (f) Different mismatch degrees for ψ me .Electronics 2021, 10, x FOR PEER REVIEW 15 of 20

Figure 11 .
Figure 11.Servo stiffness with different weighting factors for the proposed method.

Figure 13 .Figure 14 .
Figure 13.Rotor position and prediction horizon experimental results of the NMPSC with the proposed method.The rotor position and reference waveforms at steady-state and transient state are shown in Figure 14a,b respectively.It is shown that the actual rotor position can track the reference with a short delay time, and the disturbance can be resisted with little influence when a 2 N•m load torque is uploaded onto the shaft.

Figure 15
Figure15provides the experimental results of convent the same conditions.It shows that longer delay time and la than that in Figure14, and a slight additional delay time dashed line shown in Figure15b, is generated when the shaft.

Figure 14 .
Figure 14.Rotor position and reference experimental results of the NMPSC with the proposed method.(a) Steady-state; (b) Transient state.

Figure 14 .
Figure 14.Rotor position and reference experimental results of the NMPSC with the proposed method.(a) Steady-state; (b) Transient state.

Figure 15 .
Figure 15.Rotor position and reference experimental results of the conventional NMPSC.(a) Steady-state; (b) Transient state.

Table 2 .
Performance comparison within 1.5 s.

Table 3 .
Ramp reference performance comparisons within 0.8 s.

Table 3 .
Ramp reference performance comparisons within 0.8 s.

Table 4 .
Comparison of servo stiffness with different weighting factors.

Table 5 .
Performance comparison in three cycles.