Parameter identification for MTPA control based on a nonlinear d-q dynamic IPMSM model

Abstract

Purpose

Small highly saturated interior permanent magnet- synchronous machines (IPMSMs) show a very nonlinear behaviour. Such machines are mostly controlled with a closed-loop cascade control, which is based on a d-q two-axis dynamic model with constant concentrated parameters to calculate the control parameters. This paper aims to present the identification of a complete current- and rotor position-dependent d-q dynamic model, which is derived by using a finite element method (FEM) simulation. The machine’s constant parameters are determined for an operation on the maximum torque per ampere (MTPA) curve. The obtained MTPA control performance was evaluated on the complete FEM-based nonlinear d-q model.

Design/methodology/approach

A FEM model was used to determine the nonlinear properties of the complete d-q dynamic model of the IPMSM. Furthermore, a fitting procedure based on the nonlinear MTPA curve is proposed to determine adequate constant parameters for MTPA operation of the IPMSM.

Findings

The current-dependent d-q dynamic model of the machine models the relevant dynamic behaviour of the complete current- and rotor position-dependent FEM-based d-q dynamic model. The most adequate control response was achieved while using the constant parameters fitted to the nonlinear MTPA curve by using the proposed method.

Originality/value

The effect on the motor’s steady-state and dynamic behaviour of differently complex d-q dynamic models was evaluated. A workflow to obtain constant set of parameters for the decoupled operation in the MTPA region was developed and their effect on the control response was analysed.

Citation

Garmut, M., Steentjes, S. and Petrun, M. (2023), "Parameter identification for MTPA control based on a nonlinear d-q dynamic IPMSM model", COMPEL - The international journal for computation and mathematics in electrical and electronic engineering, Vol. 42 No. 4, pp. 846-860. https://doi.org/10.1108/COMPEL-09-2022-0331

1. Introduction

Small interior permanent magnet- synchronous machines (IPMSMs) with high power-to-weight ratios behave very nonlinear, due to the slotting effect, permanent magnets, cross coupling, cross saturation and very saturated magnetically nonlinear iron core (Hadžiselimović et al., 2007; Liu et al., 2016). As such machines are generally controlled by using linear closed-loop cascade control, which is based on the constant concentrated parameter d-q dynamic model, adequate control parameters have to be identified (Rafaq and Jung, 2020; Zhu et al., 2021; Lu et al., 2022). Different complexity levels with different dependencies of the d-q models were evaluated by Pouramin et al. (2015), Stumberger et al. (2003), Michalski et al. (2019) and Weidenholzer et al. (2013).

The discussed effects can be simulated by using a finite element method (FEM) model, which was the basis for developing a complete current- and position-dependent nonlinear d-q dynamic model (Fontana and Bianchi, 2020; Müller et al., 2022). This model was composed of various terms and dependencies, which model different nonlinear effects. By evaluating each of those effects, we analysed which of the terms and dependencies have a significant influence on the machine’s behaviour in respect to speed control. The simplest model that models the essential dynamic behaviour of the complete d-q dynamic model was chosen for the parametrization of the maximum torque per ampere (MTPA) cascade control (Li and Wang, 2019; Dianov et al., 2022; Fontana and Bianchi, 2020).

In the presented analysis, the IPMSM was operated on the MTPA curve below base speed, where the MTPA curve obtained from the full nonlinear model was chosen to be the reference for control parameter identification. The proposed control parameter set was extracted based on the inductance maps by using the constant parameter MTPA equation, where curve fitting to the reference MTPA curve was performed. This set of parameters was tested using a MTPA control implementation and compared with two different sets of constant parameters that were chosen arbitrarily. Also, an array set of inductances that were extracted from the nonlinear MTPA line and were dependent on the current was used for the control. The dynamic performance of the four parameter sets was compared and analysed.

The aim of this paper was a systematic analysis of discussed influences through step-by-step reduction of the complexity of the dynamic model. This gives a deeper understanding of the behaviour of the machine and highlights the effects that are important to consider for speed control of the IPMSM. The paper presents how to determine adequate constant concentrated parameters of the d-q dynamic model for control in the MTPA operation and highlights the effect of different constant parameter sets on the dynamic behaviour of the control.

2. Theoretical background

2.1 IPMSM d-q dynamic model

The complete current-and position-dependent d-q dynamic model is presented in this section. The motion is described by:

(1) Jd2θmdt2=Te−Tl−kfdθmdt−kC(dθmdt)2,

where J is the moment of inertia, k_f is the viscous friction coefficient, k_C is the ventilation coefficient, T_e is the electromagnetic torque and T₁ is the load torque. The electrical position (angle) is θ = pθ_m, where p is the number of pole pairs and θ_m is the mechanical rotor position. The voltage-balance equation in the d-q reference frame is given by:

(2) [uduq]=[R00R][idiq]+[Ld,iLdq,iLqd,iLq,i]ddt[idiq]+dθdt([∂Ψd∂θ∂Ψq∂θ]+[0−Lq,aLd,a0][idiq]+[0Ψmd]),

Li=[∂Ψd∂id∂Ψd∂iq∂Ψq∂id∂Ψq∂iq],La=[0−ΨqiqΨdid0] and  ∂Ψ∂θ=[∂Ψd∂θ∂Ψq∂θ],

where i_d, i_a and u_d, u_q are the d-q reference frame voltages and currents, R is the resistance and Ψ_q is the flux linkage in the q-axis due to the current excitation. The total flux linkage in the d-axis is defined as Ψd*=Ψd+Ψmd, where Ψ_d is the flux linkage in the d-axis due to the current excitation and Ψ_md is the flux linkage in the d-axis due to the permanent magnet. L_d,_i, L_dq,_i, L_qd,_i and L_q,_i are the incremental inductances, as denoted by the subscript i, where L_d,_a and L_q,_a are the apparent inductances denoted by subscript a.

2.2 Reduced model’s parameter identification

Three different parameter dependencies were constructed. Firstly, the dependencies of Ψd*, Ψ_q and T_e from the position θ and the currents i_d, i_q were obtained by using an automated process. A MATLAB script was used to simulate different i_d and i_q combinations in a two-dimensional (2D) transient FEM model using Ansys Maxwell 2D. The calculated torque T_e (i_d, i_q, θ) was stored in a lookup table (LUT). Furthermore, the matrices of incremental inductances L_i (i_d, i_q, θ), apparent inductances L_a (i_d, i_q, θ) and position derivatives of the flux linkages ∂Ψ∂θ(id, iq,θ) were calculated. Figure 1 shows the workflow for the parameter identification approach.

In the second step, the average values of Ψ_d*, Ψ_q and T_e with respect to θ were considered to calculate current-dependent incremental inductances L_i (i_d, i_q) and apparent inductances L_a (i_d, i_q). The torque was defined as T_e (i_d, i_q). With this data, 2D and 3D LUTs were generated. In the third approach, L_i and L_a were considered current and position independent, whereas T_e was calculated according to:

(3) Te=32p(Ψmdiq+idiq(Ld,a−Lq,a)).

All other parameter values, such as R, J, k_f, k_C, p and Ψ_md, were assumed constant. A simulation was constructed in Simulink, by using equations (1) and (2) in combination with adequate LUTs. The model’s inputs were u_d, u_q and T₁, whereas the outputs were i_d, i_q, T_e and θ and the phase currents i_a, i_b and i_c. In the third approach, equation (3) was used to determine the torque instead of the LUTs.

2.3 Overview of reduced models

By using the three presented approaches, five models were developed, equipped with specific features, terms and dependencies. Table 1 presents all five models, where M1 was the complete current- and rotor position-dependent d-q FEM-based reduced-order dynamic model. The complexity of reduced models was decreasing, where M5 was the simplest dynamic model with the constant parameters.

To evaluate the steady-state performance of the models, the motion equation (1) was neglected and the electrical angular velocity ω = dθdt was chosen as an input for all five models. A representative operation point (OP) was defined by u_d = −0.6 p. u., u_q = −0.9 and ω = 1.3535 p. u., referred to as OP1. All models were analysed in this OP. The constant parameters for M5 were obtained by using the steady-state i_d and i_q of M4 in OP1.

The higher harmonic content of M1, M2 and M3 was similar to the FEM model, which is presented in Figure 2 and Figure 3 when comparing the THD (i_a) and the torque ripple waveform of the three models versus the reference FEM results. By those metrics, M1 models the nonlinearities of the FEM model the most accurate. Furthermore, M4 and M5 had almost identical result of the current i_a and torque T_e in discussed comparison, as the steady-state i_d and i_q of M4 were used for obtaining the constant inductances of M5. The higher harmonics in M4 and M5 were very low, and there was no torque ripple, as the output was the mean torque. The difference between the rms(i_a) of FEM and M4 or M5 was 0.04%, where the difference between the avg(T_e) between FEM and M4 or M5 was 0.3%, which was in both cases negligibly small.

2.4 Operation of reduced models

The dynamic behaviour and steady-state values at different OPs was compared between M1, M4 and M5. In this analysis, the motion equation (1) was included, where OP1 was simulated by applying a load torque T₁ = 1 p. u. and voltage u_d = −0.6 p. u., u_q = −0.9 p. u. The operation was then changed to OP2 at 0.5 s with a load torque T₁ = 0.75 p. u. and to OP3 at 1.5 s with T₁ = 1.25 p. u. In Figure 4, the angular velocity ω and current i_d, i_q response to the change of the OPs are presented.

The stationary and dynamic behaviours of M1 and M4 in all three OPs were very similar for the current i_d, i_q and angular velocity ω response as presented in Figure 4. Steady-state values between the two models differ by less than 1%, as presented in Table 2. Only M5 displays very different current i_d, i_q and angular velocity ω values and response for OP2 and OP3 compared with M1 (Table 2). The current i_d was 145% higher in the case of M5 compared with M1 in OP3, which accounted for the biggest difference.

M4 includes the current dependency of the apparent L_a (i_d, i_q) and incremental inductances L_i (i_d, i_q); consequently, those dependencies play a major role when modelling essential dynamic behaviour of the machine, for synthesis of the discussed cascade control. In Figure 5, the incremental inductance maps dependent on the current i_d and i_q of M4 are presented. On the maps, the discussed OPs are marked, and the values of the inductances are presented in Table 3.

The nonlinear and significant change of the inductances in different OPs has a major impact on the machine’s behaviour and cannot be neglected and assumed constant. For example, the incremental inductances change in the q-axis from OP2 to OP3 for 26.6%. This also explains why M5 with constant parameter performed so different in OP2 and OP3 compared with the M4 and M1, as the parameters were identified based on the OP1. Therefore, the current dependency of the inductance was critical for the further development of control algorithms and must be considered. We can conclude that the position dependency, the terms ∂Ψ∂θ and cross-saturation did not influence the discussed specific dynamic behaviour of the machine with the change of the OPs and have just a direct influence on the higher harmonic components of the currents and, consequently, a higher torque ripple.

2.5 Control parameter identification for MTPA control

A closed-loop cascade speed control was implemented for controlling the machine, which is based on constant parameters in the control algorithm (Dianov et al., 2022). Proportional integral (PI) controllers were used for the speed and current control, and decoupling was added between the d- and q-axis voltage-balance equations. Feedforward MTPA control was implemented (Morimoto et al., 1994; Rahman and Dwivedi, 2019). By using the presented implementation, the machine was operated only on the MTPA curve and below base speed.

For determining the constant inductance parameters, the M1 and M4 MTPA curves were extracted; however, they resulted in identical curves, as presented in Figure 6. To obtain the parameters for M5, curve fitting was performed on the reference M4 MTPA curve data, by using the MTPA equation:

(4) id=Ψmd2a−Ψmd24a2+iq2

where the inductance difference a = L_q,a − L_d,a was determined using the least squares method. The calculated value of a was 4.025 µH.

The next step was to find the adequate inductance difference a on reference MTPA M4 curve. The inductances L_q,a and L_d,a were extracted from the inductance maps across the reference M4 MTPA curve (presented in Figure 8) and are presented in Figure 7 with respect to the total current i=id2+iq2. The difference that matched the determined coefficient a between them was determined at a specific current, as presented in Figure 7a, and the inductances were determined at that point, as presented in Figure 7b. The values of the of inductances were L_q,a = 78.95 µH and L_d,a = 38.7 µH.

Besides the parameter set (L_q,a and L_d,a), referred to as set 3 (S3) and identified with the discussed method, two additional sets were determined for comparison reasons. The first set S1 was determined from the OP2, and the second set S2 was randomly chosen on the reference M4 MTPA curve (not with respect to the adequate a). All the chosen constant parameter sets are as follows:

• S1: apparent inductances from the OP2: L_q,a = 77.47 µH and L_d,a = 46.70 µH;

• S2: random apparent inductances from the M4 MTPA curve: L_q,a = 82.90 µH and L_d,a = 38.25 µH; and

• S3: apparent inductance fitted to the M4 MTPA curve: L_q,a = 78.95 µH and L_d,a = 38.70 µH.

Finally, a vector parameter set (S4) of the apparent inductances was evaluated. The apparent inductances according to the reference M4 MTPA curve as a function of the total current i were chosen as control parameters and are denoted as follows:

• S4: L_q,a (i) and L_d,a (i).

All four parameter sets are presented on the apparent inductance map in Figure 8, where the change of the L_q,a and L_d,a can be observed according to the current change.

3. Results

The evaluation of the control parameters was performed based on M1 and an implementation of the cascade feedforward MTPA speed control algorithm. The apparent inductances were included in the following control parameters:

• time constant of current d- and q-axis PI regulator;

• time constant of speed PI regulator;

• decoupling; and

• MTPA control algorithm under base speed.

All other control parameters were assumed constant. The analysis of implemented control was performed using the following profile: first, the machine was started up to the reference angular velocity ω_ref and load torque T₁. When a steady-state operation was achieved, the load torque T₁ was increased. The control operation was described by:

ωref=0.9425 p.u. and Tl(t)={0.75 p.u., t<0.2 s1.25 p.u., t=>0.2 s

The responses of currents i_d, i_q and the angular velocity ω were analysed. Figure 9 shows the current i_d, i_q responses, and Figure 10 shows the angular velocity ω responses of the controlled machine for all four parameter sets.

The responses were stable in all four cases, and the reference angular velocity ω_ref was achieved. When evaluating the current responses in Figure 9, different dynamic responses between the models were observed. For the further analysis, the following performance indicators were calculated:

• M_s – steady-state value;

• t_s – settle time;

• e% – relative overshoot;

• integral square error or ISE=∫0Te(t)2dt, where e (t) is the difference to the steady-state value; and

• integral absolute error or IAE=∫0T|e(t)|dt.

For all four parameter sets and for two load torques T₁ applied, the performance indicators were calculated, as presented in Table 4. The least adequate control performance was achieved with S1, where the performance indicators were the highest almost in all cases. When comparing just the constant parameter sets, S3 that was calculated by fitting the reference M4 MTPA curve outperformed S1 and S2. The settling time t_s of S2 and S3 had similar values, whereas S1 in the case of i_d and T₁ = 1.25 p. u. did not reach the steady-state value. When comparing the S4 response with the S3 response, it performed better in the case of the current response in the d-axis i_d, where the angular velocity ω and the current in the q-axis i_q responded better to S3, as presented in Table 4.

When analysing the steady-state currents i_d and i_q shown in Table 4 of the S4 current-dependent parameters, they produce at all load torques T₁, the optimal nonlinear MTPA point, as expected, as the M4 and M1 MTPA curves align perfectly, as shown in Figure 6. The biggest difference of i_d and i_q to the M4 MTPA values was observed for the parameters of S1, where in the case of i_d and T₁ = 0.75 p. u., the difference was 38%. The lowest difference at T₁ = 0.75 p. u. of i_d and i_q to the M4 MTPA values was observed for the parameters of S2. The lowest difference at T₁ = 1.25 p. u. of i_d and i_q to the M4 MTPA values was observed for the parameters of S3. This was because both parameter sets of S2 and S3 give different MTPA approximations and were closer to the nonlinear MTPA curve of model M1 in different OPs, defined by the load torque T₁. Constant parameters for control purposes must be chosen according to a criterion, as the results showed that the two arbitrarily chosen parameter sets S1 and S2 performed worse compared with S3. The criterion chosen in these analyses was to fit the nonlinear M4 MTPA curve with constant inductances. Adequate results can be also achieved with a vector set of the inductances changing with respect to the current along the MTPA curve.

4. Conclusion

The nonlinear FEM model was used to extract the inductance matrices L_i, L_a and position derivatives of the flux linkages ∂Ψ∂θ with respect to the current i_d, i_q and position θ change. A reduced-order dynamic model was developed, which very accurately modelled the FEM model. The model was step-by-step reduced to the current-dependent M4 and the constant parameter M5. M4 modelled the essential dynamic behaviour of the nonlinear FEM-based M1, where the M5 showed a big deviation compared with the M1 when not operating in the specific operation point for which its parameters were determined. It was concluded that the current-dependent M4 can be further used for control development and the θ dependency and the neglected additional terms did not affect the dynamic behaviour. Furthermore, a reference M4 MTPA nonlinear curve was extracted, where a curve fit was performed with the constant parameter MTPA equation to obtain the constant apparent inductances from the inductance maps for S3. S3 was evaluated on a cascade MTPA speed control versus two arbitrarily chosen parameter sets (S1 = OP2 and S2 = random from the M4 MTPA curve). Additionally, the inductances across the reference M4 MTPA array in dependency from the total current i were implemented. The dynamic response performance indicators showed that the least adequate parameters were the ones from S1. The most adequate parameters were from S3, as they performed better than the constant S1 and S2 parameters and gave better responses of the current in the d-axis i_d and the angular velocity ω compared with S4. The S4 parameter achieved at each load torques T₁ the optimal nonlinear MTPA point. The results show that the proposed method for extracting the constant parameters gives adequate results and can be used in applications where the parameters for MTPA operation are limited to constant values, despite the IPMSM being highly nonlinear.

Figures

Figure 1.

Workflow for the position and the current dependencies parameter identification

Figure 2.

The phase current in phase a of all five models compared with the results of the reference FEM model

Figure 3.

The electromagnetic torque of all five reduced-order models compared with the FEM model

Figure 4.

Angular velocity and current response of M1, M4 and M5 to the change of OPs

Figure 5.

Incremental inductance maps for M4 with highlighted OP1, OP2 and OP3

Figure 6.

Best fit of the constant parameter MTPA curve on the reference M4 MTPA curve data

Figure 7.

Apparent inductances and their difference according to the reference M4 MTPA curve depending on the total current

Figure 8.

Inductance maps for M4 with highlighted S1, S2, S3 and M4 MTPA curve

Figure 9.

i_d and i_q dynamic response to a load torque T₁ step for S1, S2, S3 and S4

Figure 10.

ω dynamic response to a load torque step T₁ for S1, S2, S3 and S4

Table 1.

Model complexity reduction from Model 1 to Model 5

Model	Model features
M1	Defined by equations (1) and (2), parameters were dependent on θ, id and iq. Torque LUT - Te (id, iq, θ)
M2	M1 neglecting cross-saturation terms, i.e. Ldq,i (id, iq, θ) and Lqd,i (id, iq, θ)
M3	M2 neglecting term ∂Ψ∂θ(id, iq, θ)
M4	M3 neglecting θ dependence. Torque LUT - Te (id, iq)
M5	M4 with constant parameters that were operation point specific. Torque was defined by equation (3)

Table 2.

Steady-state values and the relative difference to M1 of M4 and M5 for all three OPs

Quantity	OP1			OP2			OP3
	M1	M4	M5	M1	M4	M5	M1	M4	M5
ω (p. u.)	1.3567	1.3510	1.3535	1.8475	1.8487	1.8368	1.058	1.0605	0.9235
Rel. diff. to M1	/	0.42%	0.24%	/	−0.07%	0.58%	/	−0.24%	12.7%
id (p. u.)	−0.692	−0.689	−0.687	−1.099	−1.095	−1.164	−0.240	−0.247	0.110
Rel. diff. to M1	/	0.47%	0.83%	/	0.39%	−5.94%	/	−2.81%	145.6%
id (p. u.)	0.639	0.643	0.642	0.484	0.484	0.452	0.890	0.890	0.997
Rel. diff. to M1	/	−0.59%	−0.39%	/	0.10%	6.66%	/	−0.03%	−12.1%

Table 3.

Inductance values at OP1, OP2 and OP3 for M4

Operation points	M4
	id (p. u.)	iq (p. u.)	Ld,a (µH)	Lq,a (µH)	Ld,i (µH)	Lq,i (µH)
OP1	−0.689	0.643	44.16	81.91	47.24	68.72
OP2	−1.095	0.484	46.70	77.47	46.64	71.85
OP3	−0.247	0.890	38.57	78.04	39.50	52.71

Table 4.

Dynamic response performance indicators of S1, S2, S3 and S4

Quantity	S1		S2		S3		S4
	T1 = 0.75 p. u.	T1 = 1.25 p. u.	T1 = 0.75 p. u.	T1 = 1.25 p. u.	T1 = 0.75 p. u.	T1 = 1.25 p. u.	T1 = 0.75 p. u.	T1 = 1.25 p. u.
id [p. u.] − Ms	−0.86	/(−0.209)	−0.118	−0.276	−0.109	−0.257	−0.139	−0.254
iq [p. u.] − Ms	0.563	0.892	0.555	0.873	0.557	0.879	0.551	0.880
ω [p. u.] − Ms	0.943	0.943	0.943	0.943	0.943	0.943	0.943	0.943
id [p. u.] − ts	0.22	/	0.13	0.067	0.12	0.062	0.11	0.054
iq [p. u.] − ts	0.18	0.20	0.11	0.48	0.11	0.045	0.11	0.052
ω [p. u.] − ts	1.56	0.16	0.1	0.037	0.1	0.034	0.10	0.045
id − e%	452.07	159.65	407.69	121.67	404.12	108.26	212.72	43.42
iq − e%	112.49	18.87	118.17	15.43	117.63	12.54	122.39	18.02
ω − e%	1.95	1.72	0.74	0.14	0.76	0.10	1.20	0.47
id – ISE	10.36	2.77	13.09	1.01	10.99	0.71	4.92	0.19
iq −ISE	32.68	1.85	29.56	0.84	30.63	0.83	33.54	1.13
ω −ISE	3.56 105	642.64	3.49 105	397.8	3.47 105	359.16	3.53 105	497.90
id – IAE	1.03	0.57	0.93	0.16	0.86	0.14	0.58	0.08
iq −IAE	1.68	0.42	1.42	0.14	1.44	0.13	1.50	0.18
ω −IAE	132.26	6.67	126.53	3.39	126.29	3.09	128.43	3.52