R Delpoux
email: romain.delpoux@ec-lille.fr

M Bodson
email: bodson@eng.utah.edu

T Floquet

Parameter estimation of permanent magnet stepper motors without mechanical sensors

Keywords: Permanent magnet stepper motors, synchronous motors, sensorless systems, parameter identification, elimination theory

The paper presents a new sensorless parameter identification method for permanent magnet stepper motors. Current sensors are assumed available, but mechanical sensors are not. Data is obtained with open-loop commands at multiple speeds. A new frame is proposed that presents advantages similar to the d -q frame, but without the need for a position sensor. The method exploits derived linear parameterizations and least-squares algorithms. In some cases, overparameterization is resolved using elimination theory. The parameters identified using the new procedure are found to be very close to those obtained with sensors. The approach is potentially applicable to other types of synchronous motors.

Introduction

Permanent Magnet Stepper Motors (PMSM's) are widely used in industry for position control, especially in manufacturing applications. PMSM's are more robust than brush DC motors and produce high torque per volume.

They are often controlled in open-loop, although the potential loss of synchronism limits operation away from resonances and from high acceleration trajectories. These problems can be resolved by using closed-loop control methods with position sensors of sufficient precision. Recent research has focussed on whether the performance of closed-loop control methods could be achieved using sensorless systems. In this case, sensorless refers to systems that do not have position sensors, although current sensors are still assumed to be available. Sensorless control is useful to reduce the cost of the application, or when there is no space for a mechanical sensor. Current sensors can reconstruct the position of the rotor through the induced back-emf voltages at nonzero speeds [START_REF] Johnson | Review of sensorless methods for brushless DC[END_REF], [START_REF] Schroedl | Sensorless control of permanent-magnet synchronous machines: An overview[END_REF], [START_REF] Tomei | Observer-Based Speed Tracking Control for Sensorless Permanent Magnet Synchronous Motors With Unknown Load Torque[END_REF] and [START_REF] Shah | Sensorless Speed Control of Nonsalient Permanent Magnet Synchronous Motors[END_REF]. For such methods to succeed, the model of the motor and its parameters have to be well known, which brings to the forefront the question of parameter identification without position or velocity sensors, and notably an initial scenario for off-line parameter identification. Indeed, the electric motor manufacturers provide the parameters of the motor itself, without load. Moreover, these parameters are nominal, and therefore uncertain. Ultimately, sensorless identification could be used to provide autotuning of a sensorless control law, real-time adaptation, and fault detection.

The estimation of PMSM parameters was studied in [START_REF] Blauch | High-speed parameter estimation of stepper motors[END_REF], [START_REF] Kim | Improved current regulators for IPM machine drives using on-line parameter estimation[END_REF], and [START_REF] Mobarakeh | On-line Identification of PMSM Electrical Parameters Based on Decoupling Control[END_REF], but with rotor position information. Position sensorless identification was applied using special signals at standstill or under load condition in [START_REF] Nee | Determination of d and q reactances of permanent-magnet synchronous motors without measure-ments of the rotor position[END_REF], but for the identification of the d and q reactances only. Other methods to identify motor parameters online include [START_REF] Bolognani | On-line parameter commissioning in sensorless PMSM drives[END_REF], [START_REF] Lee | An online identification method for both stator resistance and back-EMF coefficient of PMSMs without rotational transducers[END_REF], but [START_REF] Bolognani | On-line parameter commissioning in sensorless PMSM drives[END_REF] only provides simulation results and in [START_REF] Lee | An online identification method for both stator resistance and back-EMF coefficient of PMSMs without rotational transducers[END_REF], only the stator resistance and the back EMF constant are identified. In [START_REF] Ichikawa | Sensorless control for all types of synchronous motors using an on-line parameter identification[END_REF], [START_REF] Ichikawa | Sensorless control of permanent-magnet synchronous motors using online parameter identification based on system identification theory[END_REF], [START_REF] Yoshimi | Parameter identification for IPMSM position sensorless control based on unknown input observer[END_REF], parameter identification is realized in the d -q frame, where the position needed for the d -q transformation is estimated using identified parameters. This type of structure may be successful in practice, but guarantees of stability and convergence are absent, because parameter estimation depends on position estimation and vice-versa. This paper presents a new experimental off-line method for the identification of the parameters of a PMSM without position or velocity sensors, using open-loop command of the PMSM and assuming that the velocity is equal to the reference velocity on the average in steady-state, i.e., that the motor synchronism is kept. Based on well-known parameter identification approaches such as the least squares algorithm and elimination theory, the contribution of this paper originates from a new change of variables leading to a frame of reference which is advantageous for sensorless applications. Compared to existing approaches, the method has the advantages of: identifying all of the electrical parameters as well as the mechanical parameters, deriving identification algorithms that are guaranteed to converge, validating the analytical results with experimental data. The results obtained without position sen-sors are compared to those obtained with sensors following the approach from [START_REF] Blauch | High-speed parameter estimation of stepper motors[END_REF]. The theory is validated through experiments that were performed using a test bench available at the LAGIS laboratory at the École Centrale de Lille. The paper extends results presented at the 2012 American Control Conference [START_REF] Delpoux | Parameter estimation of permanent magnet stepper motors without position or velocity sensors[END_REF]. Compared to [START_REF] Delpoux | Parameter estimation of permanent magnet stepper motors without position or velocity sensors[END_REF], this paper provides a more extensive comparison of a method identifying the resistance based on the voltage equation to a method using the absorbed electrical power and the mechanical equation.

The article is divided into three parts. Section 2.1 presents the model of the PMSM in three different reference frames and the identification algorithms used in the paper. In Section 3, an identification procedure is developed for motors with position and velocity sensors, to be used as a basis for comparison. The last section 4 presents the new identification procedure and the results obtained experimentally.

Stepper Motor

d -q transform f -g transform Inverse d -q transform Inverse f -g transform v d v q θ θ v f v g θ r v a v b i a i b i d i q i f i g Ω Figure 1
: Global scheme of the PMSM with dq and fg transformations.

Preliminaries

PMSM Model

In this section, the model of the PMSM is presented in three different frames, including a new frame that is particularly useful for sensorless applications.

Model in the phase variables (a -b)

Equations ( 1) give the standard PMSM model in the phase (or winding) variables

                 L di a (t) dt = v a (t) -Ri a (t) + KΩ(t) sin(Nθ(t)), L di b (t) dt = v b (t) -Ri b (t) -KΩ(t) cos(Nθ(t)), J dΩ(t) dt = K (i b (t) cos(Nθ(t)) -i a (t) sin(Nθ(t))) -f v Ω(t) -C r sgn(Ω(t)).
(

) 1
where v a and v b are the voltages applied to the two phases of the PMSM, i a and i b are the two phase currents, L is the inductance of a phase winding, R is the resistance of a phase winding, K is the back-EMF constant (and also the torque constant), θ is the angular position of the rotor, Ω = dθ/dt is the angular velocity of the rotor, N is the number of pole pairs (or rotor teeth), J is the moment of inertia of the rotor (including the load), f v is the coefficient of viscous friction, and C r is the coefficient of Coulomb friction.

Model in the rotating frame (d -q)

The phase model can be transformed using Park's transformation [START_REF] Park | Two-reaction theory of synchronous machines generalized method of analysis-part I[END_REF]:

[i d , i q ] T = M p (θ) [i a , i b ] T , ( 2
) [v d , v q ] T = M p (θ) [v a , v b ] T , (3)
where

M p (θ) =   cos(Nθ) sin(Nθ) -sin(Nθ) cos(Nθ)   . (4)
Using this change of coordinates, the system (1) is transformed into the so

-called d -q model            L di d (t) dt = v d (t) -Ri d (t) + NLΩ(t)i q (t), L di q (t) dt = v q (t) -Ri q (t) -NLΩ(t)i d (t) -KΩ(t), J dΩ(t) dt = Ki q (t) -f v Ω(t) -C r sgn(Ω(t)).
(5)

The d -q transformation is commonly used for PMSM's (and synchronous motors in general), because it results in constant voltages and currents at constant speed (instead of the high-frequency phase variables). Also, the model highlights the role of the quadrature current i q in determining the torque. However, the d -q transformation is based on the position θ, which is not directly available in sensorless applications.

Model in the rotating reference frame (f -g)

The goal of this article is the identification of the parameters without the need for the position and the velocity. The d -q transformation is not suitable for such a purpose, given that the transformation M p (θ) uses the position θ. To overcome the problem, a solution is to use the model in the phase variables to identify the parameters. However, the high frequencies of the phase variables at high speeds pose problems for the identification, in particular by making it difficult to filter out measurement noise.

In this article, a different frame is proposed that uses a reference position instead of the real position. The model in these variables is obtained using the transformation (2) and (3) but the Park matrix is defined as follows:

M p (θ r ) =   cos(Nθ r ) sin(Nθ r ) -sin(Nθ r ) cos(Nθ r )   . (6)
The PMSM model in the transformed variables is

                 L di f (t) dt = v f (t) -Ri f (t) -KΩ(t) sin(N(θ r (t) -θ(t))) + LNΩ r (t)i g (t), L di g (t) dt = v g (t) -Ri g (t) -KΩ(t) cos(N(θ r (t) -θ(t))) -LNΩ r (t)i f (t), J dΩ(t) dt = K(i f (t) sin(N(θ r (t) -θ(t))) + i g (t) cos(N(θ r (t) -θ(t)))) -f v Ω(t) -C r sgn(Ω(t)), (7)
where Ω r = dθ r /dt.

The variables in the d -q and f -g frames are related through the transformation

  i f i g   =   cos(δ) sin(δ) -sin(δ) cos(δ)     i d i q   , (8)
where

δ = N(θ r (t) -θ(t)), (9)
is the difference between the angles of rotation of the two transformations.

The f -g frame is potentially useful in two ways. First, θ r may be defined as the reference position that the motor is supposed to track. For

example, an open-loop voltage control strategy consists in setting v f = V 0 , v g = 0 (or any other constant values). In this case, δ is the angle by which the rotor lags the reference position (scaled by N). θ r could also be defined as an estimate of θ, using a procedure to be determined. Then, the f -g model approximates the d -q model, with the advantage that it is valid and computable even if θ r is not exactly equal to θ. Note that when θ r = θ, the model in the f -g frame is the same as the model in the d -q frame. In this paper, we consider the first option, namely setting θ r as a reference position in an open-loop control strategy.

Constant speed operation with open-loop commands

For constant v f , v g , and Ω r , the phase voltages v a and v b are sinusoidal voltages with frequency NΩ r . The motor will reach a constant speed Ω = Ω r if constant values of i f , i g , and δ exist such that

         v f = Ri f + KΩ r sin(δ) -LNΩ r i g , v g = Ri g + KΩ r cos(δ) + LNΩ r i f , K(i f sin(δ) + i g cos(δ) = f v Ω r + C r sgn(Ω r ). (10)
We propose the following fact.

Fact: for an arbitrary Ω r , a steady-state equilibrium of the f -g model exists with Ω = Ω r provided that the voltage is sufficiently large.

Proof: let V , α, Z, β such that

v f = V cos(α), v g = V sin(α), (11)
R = Z cos(β), LNΩ r = Z sin(β). ( 12
)
Then, the solution of the first two equations of (10) gives

  i f i g   = 1 Z   V cos(α -β) -KΩ r sin(δ + β) V sin(α -β) -KΩ r cos(δ + β)   . ( 13
)
Substitution in the third equation of ( 10) then gives

sin(δ + α -β) = 1 V Z K (f v Ω r + C r sgn(Ω r )) + KΩ r cos(β) . (14)
Letting V be sufficiently large that the right-hand side has magnitude less than 1, the equation can be solved for δ, and the value of δ can be substituted in the previous equation to get i f and i g . Thus, an equilibrium state exists, which corresponds to a limit cycle in the orginal phase variables.

Comments: the equation for δ normally has two solutions. However, the second is typically unstable. Both can also be unstable. The local stability of the system in the f -g frame can be determined by linearizing the nonlinear f -g model around the equilibrium state. Such an analysis must typically be performed numerically (see [START_REF] Verghese | Analysis of Instability in Electrical Machines[END_REF] for a detailed analysis). Typically, analysis and experiments with stepper motors show that the equilibrium is stable at low speeds but becomes unstable above a certain threshold. Some motors become stable again at speeds above a sec-ond threshold, as is the case for the motor used in the experiments of this paper. Even for speeds corresponding to a stable equilibrium in the f -g frame, fluctuations in the friction torque and poorly damped modes result in oscillations of the position and speed around the reference. For this reason, it is useful to filter the data over some period of time. In practice, motors that do not reach a stable equilibrium at Ω = Ω r typically stall, so that lack of synchronism can be detected without position sensors (at least visually by the operator) . In the analysis and experiments of the paper, the assumption is made that, on the average, Ω = Ω r when using open-loop commands in the f -g frame.

Identification algorithms

Least-squares identification

The core of the identification procedures of this paper is the least-squares algorithm [START_REF] Blauch | High-speed parameter estimation of stepper motors[END_REF][START_REF] Söderström | System Identification[END_REF]. The least squares algorithm is applied when the equations are linear with respect to the parameters to be estimated. Specifically, assume that

y[n] = W T [n]p nom , (15)
where y[n] is the output vector, n is either an index or the time instant, W [n] is the regressor matrix, and p nom is the nominal (unknown) parameter vector. Given measurements of y and W, the objective is to determine p, an estimate of the nominal parameter vector p nom . The error equation is formed by subtracting the output from the estimated output. By definition, the residual error is equal to the sum of the norm squared of the error over an interval [N 0 , N 1 ]. The least-squares estimate minimizes the residual error.

In ideal conditions, the estimate is found by setting the derivative of the residual error with respect to p equal to zero, leading to

p = N 1 n=N 0 W[n]W T [n] -1 N 1 n=N 0 W[n]y[n] . ( 16
)

Elimination theory

In some cases, a linear parameterization can be obtained, but only with a set of parameters that are not independent from each other. Although one can identify the vector p using a linear least-squares algorithm, the problem is often poorly conditioned numerically, and the elements of the vector may not satisfy the constraints. When the parameters in the over-parameterized model are rationally related to the minimal set of parameters, the overparameterization can be handled by using elimination theory [START_REF] Chiasson | Elimination Theory for Nonlinear Parameter Estimation[END_REF]Wang et al., 2005a,b).

Experimental platform

Experiments were performed using a stepper motor test bench developed in the LAGIS at École Centrale de Lille (see Fig. 3). The parameters of the motor with coils in series furnished by the manufacturer are: L = 6.4mH,

R = 2.6Ω, K = 0.3N.m.A -1 , f v =

Parameter identification using a position sensor

The d -q frame is attractive to perform the identification. Most variables can be estimated using steady-state measurements at constant speed and, because the d -q variables can be averaged over batches of data, the effect of noise and unmodelled effects can be reduced. Such effects include the detent torque, amplifier biases in the voltages, measurement biases in the currents, and amplifier distortions.

Parameter estimation using steady-state measurements Identification scenario

Steady-state is reached when the velocity is constant. Using the d -q frame, this condition is easily obtained applying constant voltages v d and v q to the motor. After a short time, the variables of the motor are stabilized, the steady-state is reached, resulting in constant velocity and currents. In practice, there are actually small oscillations in the velocity and d-q currents, but these fluctuations can be averaged out. Thus, the necessary variables for the identification are recorded during a short period of time and averaged, yielding one data vector for the least-squares algorithm. At steady-state, the equations of the motor in the d -q frame are given by

         0 = v d -Ri d + NLΩi q , 0 = v q -Ri q -NLΩi d -KΩ, 0 = Ki q -f v Ω -C r sgn(Ω), (17)
where γ represents the average of a variable γ.

To obtain a sufficient data set, the procedure is repeated for different pairs of voltages. The set of equations in averaged steady-state is given by

         0 = v d [n] -Ri d [n] + NLΩ[n]i q [n], 0 = v q [n] -Ri q [n] -NLΩ[n]i d [n] -KΩ[n], 0 = Ki q [n] -f v [n]Ω[n] -C r sgn(Ω[n]), (18)
where n refers to the index of the experiment. Note that different pairs of voltages can yield the same velocity, and multiple values were deliberately applied to increase the richness of the data.

The paper considers two different approaches for the identification of the parameters in the presence of mechanical sensors.

First method

For this first method, it is proposed to identify the parameters directly from the system of equations ( 18). The parameters R, L, and K can be obtained from the two first equations. Then, using the estimate of K, the mechanical equation is used to identify the mechanical parameters. The first two equations of ( 18) can be written in the form of ( 15) with

y 1 [n] =   v d [n] v q [n]   , W T 1 [n] =   i d [n] -NΩ[n]i q [n] 0 i q [n] NΩ[n]i d [n] Ω[n]   , p nom,1 = R L K T . (19)
The estimate of K, denoted by K, allows one to write the mechanical equation as:

y 2 [n] = Ki q [n] W T 2 [n] = Ω[n] sgn(Ω[n]) p nom,2 = f v C r T . ( 20
)
By application of equation ( 16) to the two previous equations, one obtains estimates of the motor parameters. The results of identification using the two-step procedure are shown in Table 1.

In order to verify the estimated parameters, the set of data is reorganized Fig. 4 shows the output vector from (19) (which represents the averaged values of v d and v q of each experiment) as a function of Ω, as well as the fit from the least-squares algorithm W T 1 [n]p 1 . The jagged appearance is due to the fact that different pairs of voltages v d and v q applied to the motor can yield the same velocity, as mentioned before. The match of the least-squares fit with the data on the figure is found to be very good. The result of the fit of the mechanical equation ( 20) is shown on Fig. 5. In this case, the output variable is the torque. The figure shows that the fit is good, although some unmodeled effects are not matched by the estimated torque. However this uncertainty is not important. Typically, the integrator of a velocity or position control loop will compensate for such effects.

Second method

The parameters of the mechanical equation can be obtained using an alternative method that is more suited to the sensorless case. Specifically, the parameters R, f v and C r are identifiable based on the electric power

v d vd v q vq Ω[n] in rad.s -1 Ω[n] in rad.s -1 v d [n], vd [n] (in V ) v d [n], vd [n] (in V )
Figure 4: Representation of the identification of the electrical parameters in the dq frame.

converted to mechanical power

P em v d i d + v q i q -R(i 2 d + i 2 q ) = KΩi q . ( 21
)
In the steady-state Thus, the least-squares algorithm can be applied with

KΩi q = f v Ω 2 + C r |Ω|. ( 22
) 20 40 60 80 100 0 0.1 0.2 0.3 y 2 ŷ2 Ω[n] in rad.s -1 y 2 [n], ŷ2 [n] (in Nm)
y 3 [n] = v d [n]i d [n] + v q [n]i q [n], W 3 [n] =      i 2 d [n] + i 2 q [n] Ω 2 [n] |Ω[n]|      , p nom,3 =      R f v C r      . ( 23
)
Using the estimate of the parameter R, equation ( 19) is reformulated as the following least-squares fit

y 4 [n] =   v d [n] -Ri d [n] v q [n] -Ri q [n]   , W T 4 [n] =   -NΩ[n]i q [n] 0 NΩ[n]i d [n] Ω[n]   , p nom,4 =   L K   . ( 24
)
The results for this identification are reported in Table 1. The corresponding results are plotted on Fig. 6 significant.

Parameter estimation of inertia using transient measurements

The inertia J does not affect the dynamics in steady-state, and it is the only parameter that requires the measurement of a transient response. Using the estimates K, fv , and Ĉr obtained in steady-state, the mechanical equation in the d -q frame gives the equation

y 5 (t) = Ki q (t) -fv Ω(t) -Ĉr sgn(Ω(t)) = Ω(t) J . ( 25
)
which is linear in the parameter J to be estimated. In this case, t refers to the time instant. A large step in voltage is applied to produce a large acceleration dΩ/dt, which is useful to obtain a reliable estimate of the inertia.

The angular acceleration can be reconstructed from the velocity using the

Ω[k] = Ω[k] -Ω[k -1] T difference equation
y 2 ŷ2 Ω[n] in rad.s -1 Ω[n] in rad.s -1 y 4,1 [n], ŷ4,1 [n] (in V ) y 4,2 [n], ŷ4,2 [n] (in V )
Figure 7: Representation of the identification of the electrical parameters in the dq frame using the second method.

to be reconstructed from position using a similar procedure. However, the testbed used in the experiments of this paper included a tachometer as well as an encoder, so that it was not necessary to do so.

When step changes in voltages are applied, the velocity changes rapidly, yielding an acceleration profile such as shown in Fig. 8. The response exhibits spikes, which are well matched by the least-squares fit. The identification results are shown in Table 1. Using the d -q frame, all the parameters of the model have been identified. The results provide a reference to compare with in the sensorless identification section presented next.

Sensorless parameter identification

This section develops the results that form the contributions of the paper.

Parameter estimation using steady-state measurements Identification procedure at steady-state

Identification experiments are performed in the f -g frame similarly to the d -q frame, specifically with constant voltages v f and v g applied at constant speed. However, as opposed to the d-q frame, the PMSM cannot be operated for all speeds in this manner. There are regions of the speed range where the PMSM looses synchronism. The magnitude of the acceleration is also limited. Experiments were performed using the following proportional controller in the rotating reference frame

v f = max(k(I max -I), V max ), v g = 0, ( 26
)
where I 2 = i 2 a + i 2 b is the peak current in steady-state, I max and V max are the current and voltage limits, and k is an adjustable gain. This control law was used to maximize the currents and voltages while respecting the limits.

Because the method does not rely on position and velocity sensors or on estimates of these variables, the control algorithm will be referred to as the open-loop controller. It can then be assumed that Ω r = Ω on the average as long as the motor keeps synchronism. Indeed, equation ( 7 at Ω r = Ω for constant voltages v f and v g , and for some constant values i f , i g and θ r = θ. The equilibrium is not always stable, but lack of convergence can be verified without a sensor because loss of synchronism results in a stalled motor and/or considerable vibrations. Later in the section, Ω is replaced by Ω r in the equations.

In order to generate a sufficient data set for the identification, multiple values of Ω r were reached. For this purpose, a smooth reference velocity was applied from 0rad.s -1 to the desired velocity. Once the desired velocity was reached, the motor was allowed to settle at that speed and the measured variables were averaged over a small time interval, to get one data point for the identification. This procedure was repeated for different reference velocities.

To show that the averaged velocity reached the reference velocity, Fig. 9 shows the value of the averaged velocities together with the reference ones for each experiment. One can see that the averaged values are exactly equal to the desired ones, which shows that the assumption is experimentally verified.

The dashes represent a region between 14 and 28 rad/s with no data points.

In this region, the motor is able to reach the desired velocity, but then starts to oscillate and looses synchronism. However, speeds well above the unstable region could be reached as long as the speed reference was slowly increased without lingering at any speed in the region.

Almost all the parameters are estimated with steady-state experiments in the f -g frame, leaving the identification of the inertia to a well-designed transient experiment. The identifications at steady-state are obtained using the same data set. The data set, which comprises the reference speed, the currents and the voltages in the f -g frame, is obtained from the experiments, with the reference speed shown in Fig. 9. Note that different experiments are used later for the identification of inertia, as for the d -q experiments. In steady-state, the variables are related through

KΩ r sin(N(θ r -θ)) = v f -Ri f + LNΩ r i g , ( 27
) KΩ r cos(N(θ r -θ)) = v g -Ri g -LNΩ r i f , ( 28
) 0 = K(i f sin(N(θ r -θ)) + i g cos(N(θ r -θ))) (29) -f v Ω r -C r sign(Ω r ). ( 30
)
Although linear parameterizations are apparent, a different procedure must be used due to the unknown rotor position.

Consider (27) 2 + (28) 2 , which yields

K 2 Ω 2 r = v 2 f + v 2 g + R 2 (i 2 f + i 2 g ) -R(2v f i f + 2v g i g )+ L 2 (N 2 Ω 2 r (i 2 f + i 2 g )) + L(2NΩ r (v f i g -v f i g )).
(31) Substituting ( 27) and ( 28) in ( 30), one finds

v f i f + v g i g = R(i 2 f + i 2 g ) + f v Ω 2 r + C r |Ω r |. ( 32
)
Equations ( 31) and ( 32) are independent of the rotor position and provide linear parameterizations. However, the parameters R and L appear together with R 2 and L 2 .

First method

The parameters R, L, and K are identified using equation ( 31) and f v and C r using equation ( 32). ( 31) can be written for the set of averaged data points in the standard form (15) as

y 6 [n] = v 2 f [n] + v 2 g [n], (33)
W 6 [n] =            2(v f [n]i f [n] + v g [n]i g [n]
)

-(i 2 f [n] + i 2 g [n]) -2NΩ r [n](v f [n]i g [n] + v g [n]i f [n]) -N 2 Ω 2 r [n](i 2 f [n] + i 2 g [n]) Ω 2 r [n]            , (34)
p nom,6 = R R 2 L L 2 K 2 T = P 1 P 2 P 3 P 4 P 5 T . (35)
Although the representation is linear in the parameters, it is over-parameterized with P 2 = P 2 1 and P 4 = P 2 3 . These constraints can be handled using elimination theory (details given in Appendix A).

A standard least-squares algorithm is applicable for equation (32), which is linear with respect to the parameters. Using R, the estimate of R obtained previously, one finds that

y 7 [n] = v f [n]i f [n] + v g [n]i g [n] -R(i 2 f [n] + i 2 g [n]), (36)
W T 7 [n] = Ω 2 r [n] |Ω r [n]| , (37)
p nom,7 = f v C r T . (38)
The estimated parameters are shown in Table 2. We observe that the estimated parameters are close to the parameters obtained using mechanical sensors in the d -q frame. Moreover, Fig. 10 shows the output vector for the electrical equation as a function of Ω r , as well as W T 6 [n]p 6 as a function of Ω r . The estimate is close to the measured signal and the parameters are close to the parameters estimated using mechanical sensors. Fig. 11 shows the output signal for the mechanical equation, as well as its estimate.

Note that the signal is the converted electrical power and is negative at low frequencies. This is clearly not possible, and indicates that the resistance is over-estimated in the electrical equations (yielding an overestimation of the ohmic losses). This observation is the motivation for the second method, which was already used in the d -q frame. The second method estimates the resistance from the mechanical equation (expressed in terms of power), instead of the electrical equations.

Parameters

1 st Method 2 nd Method R(Ω) 3.0 2.86 L(mH)
10.4 10.4 K(Nm.A -1 ) 0.27 0.27 f v (Nms/rad) 5.35. 10 -4 2.69.10 -4 C r 0.0558 0.0742 J(kg.m 2 ) 3.13.10 -4 3.13.10 -4

Table 2: Results of the parameter estimation in the fg frame.

Second method

The parameters R, f v and C r are first identified using (32), while L and K are identified using (31). Specifically, equation ( 32) is linear in the parameters with Equation ( 31) can be written in the form of a linear equation using R, the estimate of R obtained in the first step, so that the following least-squares problem is obtained

y 8 [n] = v f [n]i f [n] + v g [n]i g [n], (39)
W T 8 [n] = i 2 f [n] + i 2 g [n] Ω 2 r [n] |Ω r [n]| , (40)
p nom,8 = R f v C r T . (41)
y 9 [n] = v 2 f [n] + v 2 g [n] -2 R(v f [n]i f [n] + v g [n]i g [n]) -R2 (i 2 f [n] + i 2 g [n]), (42)
W 9 [n] =      -2NΩ r [n](v f [n]i g [n] -v g [n]i f [n]) -N 2 Ω 2 r [n](i 2 f [n] + i 2 g [n]) Ω 2 r [n]      , (43)
p nom,9 = L L 2 K 2 T = P 1 P 2 P 3 T . ( 44
)
The system is overparameterized with P 2 = P 2 1 which can be resolved using elimination theory (details given in Appendix A).

The results for the identification are again shown in Table 2. The estimated parameters are close to the parameters estimated in the d -q frame using a position sensor. Fig. 12 shows the output variable, which is the power absorbed by the motor. The least-squares fit is excellent. The es- timated parameters are close to those obtained in the d -q frame using a position sensor. Moreover, Fig. 13 shows the output variable from (42) as a function of Ω r . The least-squares fit is very good. Overall, the second method gives better results than the first one. Further, the second method reduces the number of constraints in the nonlinear identification problem, yielding a simpler algorithm.

Identifiability of the parameters using the second method

Given the nonlinearity of the system, precise necessary and sufficient conditions for identifiability of the parameters are difficult to derive. However, partial yet useful answers can be obtained and show that it is generally easy to ensure identifiability. We limit the discussion to the second method, which was found to be preferable.

Since there are three parameters, a necessary condition for the identifiability of the mechanical parameters is that at least three data points are collected. Experiments must also include at least two separate velocities, in order to distinguish between viscous and Coulomb friction. A sufficient condition for identifiability is that three data points be collected with two experiments performed at the same velocity and different current levels, and a third experiment performed at a different velocity. Indeed, the matrix composed of the three regressor vectors is given by

     W T 8 [1] W T 8 [2] W T 8 [3]      =      I 2 [1] Ω 2 r [1] |Ω r [1]| I 2 [2] Ω 2 r [2] |Ω r [2]| I 2 [3] Ω 2 r [3] |Ω r [3]|      , (45)
where

I 2 [n] = i 2 f [n] + i 2 g [n] and, for Ω r [2] = Ω r [1], det      W T 8 [1] W T 8 [2] W T 8 [3]      = I 2 [1] -I 2 [2] |Ω r [1]||Ω r [3]| |Ω r [1]| -|Ω r [3]| , (46)
which is different from zero under the assumption. Therefore, the three parameters can be exactly estimated in theory.

For the identification of the electrical parameters in the second method, a necessary condition for identifiability is that at least two data points be collected. A sufficient condition for identifiability is that two data points be collected at the same speed but at different current magnitudes I 2 [n] (see subsection 2.2 for the open-loop command of the motor). Indeed, letting w 9,1 [n], w 9,2 [n], w 9,3 [n] be the elements of the vector W 9 [n], L can be obtained by solving the quadratic equation

y 9 [1] -y 9 [2] = (w 9,1 [1] -w 9,1 [2])L + (w 9,2 [1] -w 9,2 [2])L 2 , ( 47
)
which can be solved if

I 2 [1] = I 2 [2], since it implies that w 9,2 [1] -w 9,2 [2] = 0.
Having determined L, K can be obtained from

y 9 [1] -w 9,1 [1]L = w 9,2 [1]L 2 = w 9,3 [1]K 2 = Ω 2 r [1]K 2 . ( 48
)
As can be seen from this brief analysis, the conditions for identifiability are relatively weak. In general, one should collect significantly more data points than the minimum number required in order to reduce the effect of noise and unmodelled dynamics. However, the analysis suggests that, in addition to collecting data at speeds across the possible range, it may help to collect data at different current levels while at any given speed.

Parameter estimation of inertia using transient measurements

Without a sensor, the torque produced by the motor cannot be directly computed from the currents. In order to circumvent this problem, the technique discussed in this section is based on a computation of power rather than torque. Specifically, the electrical power converted to mechanical power is

P em v a i a + v b i b -Ri 2 a -Ri 2 b - L 2 d dt i 2 a + i 2 b . (49)
The transient inductive term is not necessarily zero as in the steady-state, but it is small in the experiments to come, and will be neglected. Using the mechanical equation of the model, one also has

P em = P kin + P f , (50)
where P kin is the kinetic power, i.e., the variation in kinetic energy P kin JΩ dΩ dt while P f is the power lost to friction

P f = f v Ω 2 + C r |Ω|.
The inertia can therefore be computed using a least-squares algorithm or, simply with

Ĵ = P em -P f Ω dΩ dt . (51)
Using (51) to estimate the inertia is difficult, because only low levels of acceleration can be achieved using the open-loop controller. Then, the converted electromechanical power is mostly equal to the power lost to friction (P em ≃ P f ) and dΩ/dt ≃ 0. However, with a carefully designed experiment, it is possible to identify the inertia even with the simple open-loop controller used for the other experiments.

Specifically, with constant v f and v g , the motor is operated at some arbitrary reference velocity Ω r = Ω 1 , and data is taken for a period [t 0 , t 1 ].

The friction power P f,1 is computed as the average of P em over [t 0 , t 1 ] (since dΩ/dt = 0). Then, the reference velocity is rapidly increased using Ω r (t) = Ω 2 1 + 2c(t -t 1 ). ( 52)

With Ω r dΩ r dt = c, where the constant c is adjusted experimentally to be as large as possible (while maintaining synchronism). The acceleration continues until some time t 2 where an arbitrary reference velocity Ω 2 is reached.

Then, a third batch of data is collected at constant speed Ω 2 for a period [t 2 , t 3 ]. The friction power P f,2 is computed as the average of P em over [t 2 , t 3 ].

The inertia estimate is then computed using

Ĵ = AV G P em -Pf c , (53)
where the average is taken over [t 1 , t 2 ] and Pf is an estimate of the friction.

Pf could be obtained using the estimated coefficients of friction obtained in the previous section. The experiments discussed here used the interpolation formula :

Pf = P f,1 + (P f,2 -P f,1 ) (t -t 1 ) (t 2 -t 1 ) . (54)
The condition for identifiability of the parameter J is simply that c = 0, i.e., that the motor actually accelerates.

The instantaneous value of the estimate of inertia and its average value for the experiment represented in Fig. 14 the detent torque are also visible. However, through averaging over a period of time, a suitable estimate of the inertia is obtained. The estimate of inertia was further refined by repeating the experiment with different velocities Ω 1 and Ω 2 , and with different values of c. The resulting estimate is shown in Table 2, and matches closely the value estimated with a position sensor.

Note that it is possible to include the inductance term that was neglected in (49) using the estimate of the inductance. However it was found that the impact of this term on the inertia estimate was no more than 0.3% of the inertia.

Conclusions

There has been considerable interest in developing sensorless control methods for synchronous motors and permanent magnet stepper motors in particular.

The objective is to replace position and velocity sensors with less costly and more reliable current sensors (which are often present anyway). Sensorless methods rely on the knowledge of a model and of its parameters, which raises the question of how these parameters can be determined without position sensors.

In this paper, different experimental approaches for the estimation of the parameters of PMSM's were presented, using position and velocity sensors as well as without any mechanical sensor. For each case, the identification was developed based on the model equations. Data was collected with the motor operated in open-loop, which is possible with PMSM's, although certain speed ranges may be prohibited and accelerations may be limited.

Although open-loop tracking of a speed reference may not be satisfactory from a control perspective, the paper showed that the information gathered was sufficient to determine all the parameters of the model. In the four cases under consideration, results were close to each other.

A special frame of reference was introduced based on the reference position. It presented the advantage of resulting in constant voltages and currents at constant speed, yet without the need for a position sensor. The reference frame may be useful for sensorless control as well. An interesting possibility is to have θ r be an estimate of the position θ, using an observer to be derived using the f -g model.

deg in P 1 deg in P 3 π 1 (P min ) 1 2 π 2 (P min ) 1 3

Table A.3: Degrees of π 1 and π 3 before elimination of P 1 .

deg in P 1 deg in P 3 r 1 (P 3 ) 0 3

Table A.4: Degree of the resultant.

Solving r 1 (P 3 ) = 0, one obtains the three solutions P 3 = {0.072, 0.032 ± 0.041j}. P 3 being a physical parameter, the solution must be real. The complex conjugate solution is removed. Substituting P 3 by its value into π 1 (P min ) (or into π 2 (P min )), and solving π 1 (P min ) = 0, gives the set of values for the parameter P 1 . Three real values are found which give three candidate solutions. The one that results in the smallest squared error is chosen. The three candidate solutions for the identification of section 4 as well as the squared error are given in the following table : The third solution is the one L K E 2 (P min ) Solution 1 0.07160 0.26758 1.7772e + 09 Solution 2 0.00246 0.26758 6.8341e + 05 Solution 3 0.01048 0.26758 1.6839e + 03 that results in the smallest squared error. The output vector and its estimate using the three candidate solutions are plotted on Fig. A.16. One observes that the third solution is indeed the one that best fits the measured data.

Fig. 1 Figure 2 :

 12 Fig.1shows the global scheme of the PMSM with the different variables used for identification. The different axes are represented in Fig.2. For the purpose of off-line parameter identification, one assumes that all the parameters are constant. On-line parameter estimation can then be used to correct for variations, if necessary, but is not considered in this paper.

  10 -3 and J = 5.10 -4 kg.m 2 . Note that these parameters are the parameters of the motor independently of the test bench, which alters some of the parameters. The number of pole pairs is N = 50. The input voltages v a and v b of each coil are delivered by two

Figure 3 :

 3 Figure 3: PMSM test-bench.

Figure 5 :

 5 Figure 5: Representation of the identification of the mechanical parameters in the dq frame.

Figure 6 :

 6 Figure6: Representation of the identification of the mechanical parameters in the dq frame using the second method.

  . The resulting signal is low-pass filtered to reduce noise (in the experiments, a third-order Butterworth with cut-off frequency at 500Hz was used). Using the Matlab function filtfilt, no delays are introduced by the filter. Note that, in general, velocity may have

Figure 8 :

 8 Figure 8: Inertia estimation in the dq frame.

Figure 9 :

 9 Figure 9: Speed tracking with open-loop controller.

)Figure 10 :Figure 11 :

 1011 Figure 10: Representation of the identification of the electrical parameters in the fg frame.

Figure 12 :

 12 Figure12: Representation of the identification of the mechanical parameters in the fg frame using the second method.

)Figure 13 :

 13 Figure13: Representation of the identification of the electrical parameters in the fg frame using the second method.

Figure 14 :

 14 Figure 14: Inertia identification experiment.

Figure 15 :

 15 Figure 15: Inertia identification result.

Table 1 :

 1 Results of the parameter estimation in the dq frame. according to increasing speed Ω[n]. Ideally, the output vector y[n] should be equal to W T [n]p, where p is the least-squares estimate of the parameters.

	Parameters	1 st Method 2 nd Method
	R(Ω)	2.86	2.88
	L(mH)	10.2	10.2
	K(Nm.A -1 )	0.26	0.26
	f v (Nms/rad) 2.37.10 -4	1.49.10 -4
	C r	0.0752	0.0805
	J(kg.m 2 )	3.18.10 -4	3.18.10 -4

Table A .

 A 5: Candidate solutions for the identification.

Acknowledgments

This work was supported by the Ministry of Higher Education and Research, Nord-Pas-de-Calais Regional Council and FEDER through the "Contrat de Projets Etat Region (CPER) 2007-2013". Marc Bodson acknowledges the support of the École Centrale de Lille through two visiting professorships.

Appendix A. Elimination Theory

More details are given in this appendix for the solution of the parameter identification of (42) using elimination theory. The problem is to minimize

subject to the constraints P 2 = P 2 1 and with p min = [

The extrema equations are π 1 (P min ) ∂R e (P min ) ∂P 1 = 0, π 2 (P min ) ∂R e (P min ) ∂P 3 = 0 , which are polynomials in the parameters P 1 and P 3 . The degrees of the polynomials are given in table A.3.

One computes the resultant to eliminate P 1 , r 1 (P 3 ) = Res(π 1 , π 2 , P 1 ) where r 1 (P 3 ) is a polynomial in P 3 . The degrees are reported in Table A.4.