Adaptive Cost Function FCSMPC for 6-Phase IMs

energies Article Adaptive Cost Function FCSMPC for 6-Phase IMs Manuel R. Arahal 1,* , Manuel G. Satué 1, Federico Barrero 2and Manuel G. Ortega 1   Citation: Arahal, M.R.; Satué, M.G.; Barrero, F.; Ortega, M.G. Adaptive Cost Function FCSMPC for 6-Phase IM. Energies 2021,14, 5222. https:// doi.org/10.3390/en14175222 Academic Editor: Adolfo Dannier Received: 19 July 2021 Accepted: 19 August 2021 Published: 24 August 2021 Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations. Copyright: © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/). 1Systems Engineering and Automation Department, University of Seville, 41092 Seville, Spain; [email protected] (M.G.S.); [email protected] (M.G.O.) 2Electronic Engineering Department, University of Seville, 41092 Seville, Spain; fbarrer[email protected] *Correspondence: [email protected]; Tel.: +34-954-48-73-43 Abstract: In this paper, an adaptive cost function FCSMPC is derived from newly obtained results concerning the distribution of figures of merits used for the assessment of stator current model-based control of multi-phase induction machines. A parameter analysis of FCSMPC is carried out for the case of a six-phase motor. After extensive simulation and Pareto screening, a new structure has been discovered linking several figures of merit. This structure provides an simple explanation for previously reported results concerning the difficulty of cost function tuning for FCSMPC. In addition, the newly discovered link among figures of merit provides valuable insight that can be used for control design. As an application, a new cost function design scheme is derived and tested. This new method avoids the usual and cumbersome procedure of testing many different controller parameters. Keywords: cost function design; model-based predictive control; multi-phase induction machine; pareto dominance; titeica surface 1. Introduction Multi-phase systems are preferred in some applications due to lower torque oscillations and harmonic content compared with conventional systems. Model-based control has found a new niche in multi-phase drives, under different schemes. Among them, Finite Control Set Model Predictive Control (FCSMPC) for stator current tracking is a popular one. In FCSMPC current control stator current control is the objective, whereas additional controllers are needed for flux and speed regulation [ 1 ]. The model-based controller computes the best Voltage Source Inverter (VSI) state by minimizing a function referred to as cost function. The multi-phase VSI holds more configurations (switching states) than three-phase ones. Usually, a decomposition into ( α−β ) and ( x−y ) subspaces is considered to ease control design as these subspaces are related to different output variables such as torque production and losses [2]. In addition to α−β current tracking and x−y current rejection, the FCSMPC for current control can accommodate other objectives by an appropriate choice of the objective function. For multi-phase IM, arguably the most sought-after trait is that of low switching frequency. Such trait is, however, not easy to acquire as the tuning of the objective function is not trivial [3]. In this paper and continuing with the seminal work of [ 4 ], a selected number of variables for assessment (performance criteria) is used. The locus of these variables for different tunings of the FCSMPC reveals a not previously reported structure that links them. This finding means that a certain conservation law is in place, preventing the control system from simultaneously improving more than one performance criteria. This has implications for cost function design as will be discussed later. The methodology followed in the paper is to produce a large number of experimental values for the performance criteria obtained from different FCSMPC tunings; then, the data is fitted to a relatively simple mathematical expression. This structure provides an simple explanation for reported results. In addition, Energies 2021,14, 5222. https://doi.org/10.3390/en14175222 https://www.mdpi.com/journal/energies Energies 2021,14, 5222 2 of 14 the newly discovered links among figures of merit provides valuable insight for the task of control design. An asymmetrical six-phase drive is considered for this study although the analysis can be carried out for other topologies. The next section provides a background on previous works that are related to this one. In Section 3, the FCSMPC for current control is summarized to provide the basis for the analysis of the figures of merit. The newly discovered link between figures of merit is presented in Section 4, where Pareto analysis and cubic titeica approximation are shown. The discovery is then applied, in Section 5, to derive a new FCSMPC cost function design scheme that is experimentally tested. 2. Related Works Control systems assessment is an important part of the engineering practice; however, it seldom appears in most journals. In the case of FCSMPC for the broader field of electric machines, assessment is reduced to a few operating points and a parameter analysis is seldom performed. For instance, ref. [ 5 ] presents the surfaces of some figures of merit for a distribution of operating points covering the speed and load range. Similarly, refs. [ 1 , 6 ] presents a large set laboratory and simulated experiments enabling the assessment of the FCSMPC using a five-phase IM. In [ 7 ], four current controllers for the six-phase IM are thoroughly compared using the Root Mean Squared (RMS) value of tracking error and Total Harmonic Distortion (THD) of the stator currents as performance indices. More recently, ref. [ 8 ] presents an analysis of predictive current control for six-phase IM with alternate winding configurations. Another related topic is that of cost function design for FCSMPC systems. This topic is mainly concerned with finding adequate cost function parameters so that the closed-loop IM behavior meets requirements. This task is not easy as the cost function parameters affect every aspect of the IM behavior. In particular, the design for current controllers must face conflicting objectives. As an example consider that tracking in α−β leads to some x−y content due to the use of voltage vectors with x−y projection. Of course, a high x−y content is undesirable as copper losses increase. But, on the other hand, a poor α−β tracking can produce speed ripple. This has been studied in the seminal paper [ 4 ] where trade-offs are discussed and recently in [9] where Pareto analysis is used. Instead of the cumbersome trial and error process, usually found in the above cited works, other methods have been proposed such as neural approximations, changes in the cost function to avoid weighting factors and replacement of the cost function by decision making schemes such as fuzzy inference. For instance, in [ 2 ] an automated method is used to select the weights of the cost function for the control of a Shunt Filter of Active Power. In [ 3 ] the MPC formulation is used to include soft constraints for a nine-phase IM drive and the effect of cost function parameters is discussed. In [ 10 ] a fuzzy adaptive speed controller and adaptive weighting factors are used to reduce the speed, torque and flux ripples. Other proposals can be found in the review of [11]. Regarding weighting factor elimination, in [ 12 ] a predictive direct torque control without cost function weighting factors is presented. A multi-objective ranking is used to decide the voltage vector to be applied from a limited set. AS a result the switching frequency and computational effort can be ameliorated. In a different approach, ref. [ 13 ] presents a dynamic virtual voltage vectors strategy designed to attain zero (on average) x−y voltage production; as a result the MPC cost function can be simplified and the λxy parameter eliminated. In this paper the opposite approach is taken following the hypothesis that cost function parameters provide some flexibility that is lost otherwise. 3. Predictive Current Control of Multi-Phase IM Multi-phase IM speed regulation can be done using a classical control loop. For the stator current control, accurate tracking of α−β components and rejection of x−y currents are needed. The electro-mechanical power conversion is due solely to α−β components. The x−yplane does not contribute to torque and its components must be minimized. Energies 2021,14, 5222 3 of 14 Figure 1illustrates FCSMPC for stator current of a six-phase motor. It works by providing the control signal u that sets the state of the VSI for a whole sampling period. This is done computing, at discrete time k , the optimal switch configuration u(k+ 1 ) for the k+ 1 period. The voltage v(k+ 1 ) is impressed on the IM with the objective of producing currents is(k+2)following a reference trajectory i∗ s(k+2). VSI model A, B1, B2 IM model J Cost function min FCSMPC i*(k+2) u(k+1) M ~ 6-phase IM 6-phase VSI 6 is(k) v(k+1) Figure 1. Diagram of FCSMPC of a six-phase IM. The objective function considers predictions for (k+ 2 ) since u(k+ 1 ) will not affect the measurements made until that time. The predictions are supplied by a mathematical model of the drive which is derived from the IM equations. After some manipulation, the predictive model takes the form ˆ is(k+2|k) = A(ω)is(k) + B1u(k) + B2u(k+1) + G(k)(1) In (1), A(ω) , B1 and B2 result from discretization of the continuous time dynamics [ 1 ]. The actual angular speed ω must be measured along with the components of the state space vector is(k) = isα,isβ,isx,isy>(k) . Also, in (1) the quantity G is due to rotor currents (unmeasurable in most cases). This term is thus estimated using a backtracking method [ 1 ] as follows G(k) = (is(k)−is(k−1))T−1 s−A(ω)is(k−1)−B1u(k−1)(2) In the above equations, Ts is the sampling time and u is a vector of VSI switch states u=(K1,K2,· · · ,Kn)>. The actual control action u(k+ 1 ) is computed at discrete time k minimizing a cost function J . In this paper J is made up of terms corresponding to the penalization of tracking errors in α−β , rejection of x−y components and a term penalizing the number of switch changes in the VSI. This number is computed as SC(k) = n ∑ h=1 |uh(k+1)−uh(k)|(3) Then the objective function can be expressed as J=kˆ eαβk2+λxykˆ exyk2+λncSC (4) where kˆ eαβk2 is the quadratic deviation of predictions from reference in α−β plane and kˆ exyk2 is the quadratic deviation of predictions from reference in x−y plane. The temporal index khas been omitted for clarity. Energies 2021,14, 5222 4 of 14 3.1. Cost Function Tuning Combining previous Equations (1)–(4) it is possible to derive u(k+ 1 ) as the VSI state that minimizes J . In FCSMPC and related approaches this is usually done by exhaustive (or restrained) exploration of the possible VSI states. In this scenario, tuning of the FCSMPC consists on selecting the values for parameters λxy and λnc . This is usually done by trial and error, having as guidance some performance indices or figures of merit as the ones presented later on in the paper. The problem of cost function tuning (selection of cost functions weights) has appeared in the literature associated with stator current control of multi-phase machines, but it also appears in connection with Predictive Torque Control and in other applications where conflicting criteria must be dealt with. In fact, it has received considerable attention in recent works. The approach taken in this paper is that cost function tuning dictates the future behavior of the IM drive. Figures of merit are usually used, in the broader context of automatic control and also in the particular case of FCSMPC, to quantify such behavior in a useful way. Having this in mind, the problem of cost function tuning is viewed as the selection of weights to produce adequate figures of merit or performance indices such as the ones presented in the following. 3.2. Performance Indices Controllers for drives can be assessed by different figures of merit such as harmonic content [ 14 – 16 ], speed of response [ 17 ], steady state torque ripple [ 18 , 19 ], current imbalance [ 20 ], commutation losses [ 21 , 22 ], and total losses [ 23 ]. Of course these figures of merit are not independent. Consider as an example the usual case of sinusoidal references for α−β stator currents. The quality of tracking is then inversely proportional to Total Harmonic Distortion (THD) so both quantities are somehow linked. Another example is current content in the x−y being linked to copper losses and current THD. Oscillations in drive speed are produced by ripples in torque that are produced by imperfect tracking in α−β plane. Finally, average switching frequency is larger for tunings where current tracking requirements are more stringent, so better tracking usually requires more commutations. Continuing with the work [ 4 ], three measures ( Γ1 , Γ2 and Γ3 ) will be used as they are directly related to most of the issues above presented. Their values are computed from experiments where references for α−β are sinusoidal and zero for x−y . After the experiment is performed, data from N sampling periods is used to compute the Γ values as follows Γ1=v u u t 1 N N ∑ k=1 e2 αβ(k)(5) Γ2=v u u t 1 N N ∑ k=1 e2 xy(k)(6) Γ3=10−3 N·Ts·fe N ∑ k=1 SC(k)(7) From the definition it is clear that Γ1 is the RMS value of the stator current tracking error in α−β . From the considerations made above, this value is directly related to stator current THD and to torque ripple, so, a low Γ1 provides not only better tracking, but also less THD and less torque ripple. Similarly, Γ2 is the RMS value of the stator current tracking error in x−y plane and is related to THD and copper losses. Clearly, the lower Γ2 the better from an energy efficiency stand-point. Finally, Γ3 is the average commutation frequency for the VSI. Please notice that in FCSMPC there is not fixed commutation rate, so one must rely on averaged values such as Γ3 . This commutation frequency must be kept low due to VSI limits and due to energy efficiency considerations as commutation losses are related to this quantity. Energies 2021,14, 5222 5 of 14 In any application one can think of, the FCSMPC designer would try to minimize all Γ1 , Γ2 and Γ3 . A good predictive model and low sampling period are key factors to achieve this objective. This route has been explored in many papers dealing with identification [ 6 , 24 , 25 ], different kinds of models (including models with observers) and schemes for diminishing the application period such as in-sample or space-vector modulation [ 5 , 26 – 28 ] and virtual voltage modulation [13,29]. Tuning of the cost function J also plays a crucial role in the observed values of the Γ indices. This paper analyzes the links between the Γ values when λxy and λnc are varied over a wide range of values. The distribution of Γ values shows a remarkable pattern never reported before that links the three figures of merit. As an application, the paper also proposes an adaptive scheme for automatic cost function tuning. 3.3. Six-Phase Drive The multi-phase IM used in this paper is a motor with two three-phase sets of windings in an asymmetrical configuration. Two three-phase VSI connected as depicted in Figure 2 are used to supply voltage to the motor. Figure 2. Six-phase voltage constellation in α−βand x−yplanes. Table 1presents the main data from the name-plate and the identified parameters for the motor. The VSI are SKS21F (Semikron) inverters commanded using the TMS320F28335 Texas Instruments DSP on a MSK28335 board. Two neutral points are used and the stator currents measured with four hall-effect sensors. The speed is measured using quadrature encoders. Load can be applied independently thanks to a DC-motor sharing the shaft. Figure 3depicts the main components of the experimental setup. Table 1. Electrical parameters for the experimental motor drive obtained by identification. Symbol Parameter Value Rs Stator resistance (Ω) 1.63 Ls Stator inductance (H) 0.2792 Rr Rotor resistance (Ω) 1.08 Lr Rotor inductance (H) 0.2886 Lm Mutual inductance (H) 0.2602 Lls Stator leakage inductance (H) 0.0189 feNominal frequency (Hz) 50 Energies 2021,14, 5222 6 of 14 uv Vdc is, ωm DC SUPPLY CONTROL BOARD DUAL 3-PHASE VSI DUAL 3-PHASE IM Figure 3. Diagram of the experimental setup including photographs of various elements. 4. Cubic Titeica Approximation The aim of this section is showing that the already reported links among figures of merit (the Γ values) follow a simple formula that has been experimentally derived. This discovery is by itself interesting as it show that some sort of conservation law is in place, meaning that it is impossible to obtain better results in terms of Γ1 without altering the other two figures of merit Γ2 and Γ3 . This is important as in any application one can think of, all three indices are important, so a compromise or trade-off solution must be sought. In this section, the procedure followed for data-gathering and processing will be presented, this will lead to the discovery that the Γ values lie in a surface and thus are tied to one another. The particular surface will be identified leading the way for an application for FCSMPC design. 4.1. Data Gathering The first step for the analysis is collecting data from the operation of the 6-phase motor driven by the FCSMPC in current control mode. For the analysis, thousands of data points have to be collected corresponding to various tunings of the cost function (i.e., to different combinations of (λxy , λnc) . To do so, an adequate model must be derived to be used in simulation. Sinusoidal excitation methods have been used to identify the parameters of the IM [ 24 ]. A Runge-Kutta method for the numerical integration of the continuous-time differential equations has then be deployed. In the simulations, the controller runs at Ts= 100 ( µ s) and it is treated as a discrete-time subsystem. The reference signal for stator current traking i∗(k)uses an amplitude I∗ s(A) and a frequency fe(Hz); thus i∗ s(t) =     i∗ sα(t) i∗ sβ(t) i∗ sx(t) i∗ sy(t)     =     I∗ ssin 2πfet I∗ scos 2πfet 0 0     (8) For each (λxy , λnc) the simulation is run for a number of sampling periods N including some electrical cycles enough to compute the quantities Γ1 to Γ3 as averages according to their definition. The fit of the simulations to actual measurements has been checked using the figures of merit as shown in [ 4 ]. Figure 4shows the experimental waveforms for two choices of (λxy , λnc) . Notice that the vertical axis for the isx curve has been scaled to provide a 5×zoom for readability. It can be seen that the effect of increasing λnc is noticeable in stator current tracking, mainly (for this particular example) in x−y components. This kind of representation is useful to compare few cases and to gain insight into the Energies 2021,14, 5222 7 of 14 problem. However, to advance in the analysis one needs to resort to some quantification as provided by the Γ indices. Then new information arises, as shown in Figure 5, that will be discussed later. 0 2 4 6 8 10 12 14 16 −6 −4 −2 0 2 4 6 t (ms) i (A) isα 5 isx 0 2 4 6 8 10 12 14 16 −6 −4 −2 0 2 4 6 t (ms) i (A) isα 5 isx Figure 4. Experimental results for (λxy = 0.5, λnc = 5 × 10 −4) ( left ) and for (λxy = 0.5, λnc =15 ×10−4)(right). 0 0.05 0.1 0.15 0.2 Γ1 0 0.1 0.2 0.3 Γ2 0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0 5 10 Γ3 λnc λxy = 1 λxy = 0.5 λxy = 0.1 Figure 5. Performance indices Γ found for some combinations of controller parameters (λxy , λnc) . The legend and horizontal rule apply to all graphs. 4.2. Pareto Optimality The values (Γ1 , Γ2 , Γ3) found for a combination (λxy , λnc) might have different merits for different applications. It is clear, however, that the lower their value, the better. What one finds in many cases is that, for some tunings, one of the Γ values decreases at the expense of an increase in the other two. The concept of Pareto optimality allows us to derive a reduced set of combinations that are of interest for design in many fields. Recall that a combination is deemed Pareto-optimal if it is not dominated by other combination. A combination Γj is said to be dominated by Γk if Γk i≤Γj i for all i and for a particular m , Γk m<Γj m holds. These conditions indicate that Γk is at least as good as Γj regarding all of their indices i , and at least for some index ( m ) the combination Γkis actually better. The Pareto frontier is the set of all combinations ΓP that are Pareto-optimal (i.e., not dominated). In this paper, Pareto-optimality has been used to exclude combinations that are not optimal regardless of the application. From a set of 4923 combinations of (λxy , λnc) Energies 2021,14, 5222 8 of 14 just 1415 turned out to be on the Pareto frontier. This screening not only reduces the number of data points to be considered but also allow for a more focused search for links among figures of merit. 4.3. Data Surface The relationship between Γ= (Γ1 , Γ2 , Γ3) and Λ= (λxy , λnc) is not straightforward. It must be noted that the controller aims to minimize the cost function that includes three terms that are related to the three components of Γ ; however, the point to point minimization of J (4) does not imply a similar minimization of the three components of Γ . This can be checked in Figure 5, where the various performance indices found in the simulations are plotted against some of the (λxy , λnc) combinations used in said simulations. Some general tendencies can be seen, but not much more; in fact, the data points seem to follow intricate paths. A low order approximation (e.g., linear, quadratic) seems useless. Please note that there are gaps in the distribution. This is due to the pruning of data points performed by the Pareto screening. In particular, the line of points for λxy = 1 ends abruptly at around λnc = 0.012, this is not a mistake, is just a reflection of the fact that not all tunings yield Pareto-optimal results. On closer inspection, the 3D distribution of the Γ values seem to fit a smooth surface as shown in Figure 6. A total of 671 points have been chosen from the Pareto optimal set in order to provide a graphical representation that is not too dense or too sparse. The points are represented with marks with a different color for each value of Γ3 , in this way it is easier to relate points in the projections. It can be seen that the distribution seems to be asymptotic for large values of Γ , indicating a 3D generalization of the hyperbola. Such generalization has the amazingly simple expression x·y·z=K(9) where (x , y , z) represent the coordinates of a point in 3D space and K denotes a constant. The cubic Equation (9) yields a titeica surface studied by G. Tzitzéica [ 30 ]. Figure 7presents the titeica surface; several cuts have been performed by horizontal planes for some z=zc showing the resulting hyperbolas of the form y=1/(xzc). 0 0.1 0.2 00.05 0.1 0.15 0.2 0.25 0.3 0 2 4 6 8 10 Γ1 Γ2 Γ3 0 0.05 0.1 0.15 0.2 0 0.05 0.1 0.15 0.2 0.25 Γ1 Γ2 0 0.05 0.1 0.15 0.2 0.25 0.3 0 1 2 3 4 5 6 7 8 9 10 Γ2 Γ3 Figure 6. 3D distribution of pareto optimal Γ values ( top ) and projections ( botttom ). Color indicates height (in the Γ3direction) and has been added to enhance perception. Energies 2021,14, 5222 9 of 14 Y X Z Figure 7. Cubic surface x·y·z=K ( left ) and its projection in the XY plane ( right ) including some horizontal cuts. It is possible to fit the Γ values to the model (9) by allowing some translation of the axis. Thus, if one takes x=Γ1−∆1 , y=Γ2−∆2 , and z=Γ3−∆3 then the observed value for K is ˆ K= 0.0051 obtained for (∆1 , ∆2 , ∆3) = (− 0.0015, − 0.0157, − 0.585 ) . This is a remarkable finding as it links the performance indices to one another via a simple expression. (Γ1−∆1)·(Γ2−∆2)·(Γ3−∆3)=ˆ K(10) This clearly shows that it is impossible to improve one index without degrading at least one of the other two. Further implications will be discussed later. Figure 8shows the titeica surface fitted to data represented by Equation (10). Please notice that the upper part of the surface has been trimmed for ease of presentation. In addition, some Γ values are presented to further study the distribution. These values are commented in the following. 1. Line A is made up of Γ values from simulations where λnc = 0. The Γ values gather close to a line placed on the surface in the upper region. In fact this line is the upper limit of Γ3 . The line is not an hyperbola (although its general shape resembles one, specially viewed from above) as it contains some wiggles or bumps specially for low values of (Γ1 , Γ2) . λxy increases from right to left producing increasing values of Γ1 and decreasing values of Γ2according to (10). 2. Line B corresponds to tunings where λxy = 0.1. The value of λnc goes from λnc = 0 at the intersection with line A to a value of λnc = 0.01 near the bottom. The switching frequency goes down as both (Γ1,Γ2)increase. 3. Line C corresponds to tunings where λxy = 0.5. The line found is similar to that of the previous case. Please note that both λxy = 0.1 and λxy = 0.5 have been proposed in the literature combined with low values of λnc . This tuning produces Γ values that are close to the Z axis, meaning that (Γ1 , Γ2) have both a low value even with a decrease in Γ3 as the surface is very steep at that region. These two choices are adequate in cases where current tracking is the main concern leaving little allowance for lowering losses. 4. Line D, finally, is made up of Γ values from simulations where λxy ·λnc = 0.006. The line is almost parallel to the XZ plane. It provides tunings where Γ2 changes little at the expense of pronounced changes in the other two indices. The cases commented above are easy to implement as they just require selecting values of the design parameters λxy and λnc . These cases provide insight into the cubic surface and the underlying relationship of Equation (10). In particular, they show that the relationship between Λand Γis far from trivial.