Improvement of Position Estimation of PMSMs Using an Iterative Vector Decoupling Algorithm

Abstract

This paper presents an improvement of sensorless techniques based on anisotropy for the estimation of the electrical angular position of synchronous machines by means of an iterative algorithm. The presented method reduces the effect of the fourth saliency harmonics on the measured signals avoiding the use of an observer or filter, thus, no additional dynamics are introduced on the system. Instead, a static algorithm based on iterative steps is proposed, minimizing the angular position error. The algorithm is presented and applied using the DFC (Direct Flux Control) technique but it is not limited to this choice. The advantages and limitations of this method are presented within this paper. The proof of the algorithm convergence is given. Simulations and experimental tests are performed in order to prove the effectiveness of the proposed algorithm.

1. Introduction

Permanent magnet synchronous machines (PMSMs) are widely used in various applications to realize drive systems with high power density and efficiency and highly dynamic behaviour. For efficient control, the electrical rotor position must be known, because it provides the orientation of the flux created by the rotor’s permanent magnets. Sensorless techniques have been developed in the literature that allow estimating the rotor position based on electrical quantities instead of measuring it with a mechanical position sensor, thus allowing reducing the cost and space or increase reliability of the overall drive system.

Originally, sensorless approaches extracted the position information from the back-electromotive force (back-EMF). Since the back-EMF is proportional to rotor speed, such estimation becomes increasingly difficult at low speeds and is impossible at standstill. To overcome this problem, anisotropy-based techniques have been developed, with the first one being the INFORM technique [1,2]. Anisotropy-based techniques exploit the position dependency of the machine inductances to estimate the position and can therefore be used even at very low speeds and standstill. Often, these approaches are also referred to as being based on machine saliency.

In fact, many different anisotropy-based techniques have been developed over the years that vary in the signals they inject, the measurements they use and how the position is estimated from these measurements. An overview can be found in References [3,4,5]. In many techniques, two-dimensional anisotropy vectors are obtained as intermediate results from which the position can be estimated. The latter can be achieved either via direct calculation using the arctangent function or by creating an error function that is fed into a controller for the position estimate. In an idealized representation of an anisotropy in PMSMs, the anisotropy vector rotates with two times the electrical rotor position around the origin. However, additional harmonic components can be present which in synchronous machines are either harmonics of the mechanical angular frequency or of the electrical angular frequency. In this case it is often referred to the machine with secondary or multiple saliencies [6,7].

Harmonics of the mechanical position are much harder to compensate because the absolute mechanical position is usually not known when using anisotropy-based approaches. Therefore, the compensation of harmonics in PMSMs concentrates usually only on harmonics with multiples of the electrical frequency, which will also be the case in this paper. Sources of such harmonics can be the motor geometry, the winding arrangement or the nonlinear behaviour of the soft-magnetic material, in particular saturation effects. These secondary saliencies, if they have a significant amplitude, lead to errors in the position estimation if not being considered in the estimation algorithm. To counteract, the secondary saliencies must first be identified either by FEM simulations or, preferably, via measurements on the real machine and then be compensated using suitable compensation structures.

Secondary saliencies at four times the electrical frequency have shown to be the most significant in PMSMs and are therefore usually the main focus of compensation efforts [8,9]. When the zero-sequence voltage is used as the source of measurements, fourth harmonics are even introduced as a systematic condition even if the self and mutual inductances of the machine phases are perfectly sinusoidal functions of two times the rotor position [10,11] and are therefore also of particular importance in these approaches.

In Reference [12], a work focusing on induction machines, three different approaches for compensation were suggested. The first approach used a scalar decoupling approach, subtracting an estimated error term from the estimated position. Such decoupling can therefore be used even in techniques that use only an error signal for making the position estimate converge, such as commonly used alternating injection techniques [8]. The second approach performed decoupling directly on the anisotropy vector by subtracting the estimated harmonics. The last approach did not subtract the secondary saliencies but considered them part of the expected anisotropy vector for calculation of the error signal. In all approaches, an error signal was created from the measured anisotropy and the expected anisotropy vector by a function that the authors refer to as a vector-cross-product. This error-signal was then used as feedback to an observer to estimate the position.

In Reference [8], the scalar and the vector decoupling approaches have been presented in a slightly different way, replacing the cross-product function by an arctangent function and a subtraction of the estimated position and simplifying the observer structure to a PI controller and an integrator. A new method was also presented that makes use not only of the estimated anisotropy angle but also of the magnitude and can therefore be used with machines where the secondary saliencies are so significant that the angle is no longer monotonous.

The aforementioned approaches rely on observer structures to make the position estimate converge to zero using controllers and integrators. Depending on the tuning of the controller, the bandwidth of these observer structures may limit the dynamics of the motor control if the gains are too low. For high gains, stability issues may arise, so that a compromise has to be found. In this paper, we propose a vector decoupling method that replaces the observer structure by a fixed number of iterations performed at each sample step of the time-discrete control and that we will refer to as Iterative Vector Decoupling (IVD). One focus of this paper is on the analysis of the initial error and the proof of convergence for the case of a secondary component of four times the electrical frequency. As the source of the anisotropy-information, the Direct Flux Control (DFC) technique is considered, which is based on zero-sequence voltage measurements and leads therefore to prominent fourth harmonic components in the anisotropy vector [13,14,15]. The focus on a single technique does not let the problem lose generality, since the approach can be applied to any anisotropy-based sensorless techniques that exhibit a fourth order anisotropy harmonic.

This article is organized in three parts. First, a brief introduction to PMSMs is given. PMSMs with accessible star-point are considered, since the investigated sensorless technique is based on voltage measurement at the neutral point. Then, the expression of the obtained two-dimensional saliency vector is presented. The second part introduces the IVD method applied to the DFC technique (IVD-DFC) and the analysis of the obtained signal expressions as well as the algorithm convergence proof. The analysis is presented considering saturation effects on the machine as well as variations on the parameters. Eventually, the results obtained from the experimental setup using the IVD-DFC are compared and discussed.

2. Estimation Principle of the DFC Technique

In this section, the model of a PMSM and the dynamic behavior of its star point voltage are given. Moreover, the DFC sensorless technique is briefly presented since it is used for the purpose of this work. As will be discussed, the IVD can be applied to every sensorless techniques that obtains a two-dimensional anisotropy signal vector as an intermediate result.

2.1. Pmsm Model

The following equation describes the typical electrical dynamics of a general three-phase electric rotating machine.

$${\mathbf{v}}_{abc}\left(t\right)=\mathbf{R}{\mathbf{i}}_{abc}\left(t\right)+\frac{\mathrm{d}}{\mathrm{d}t}\left[{\mathbf{L}}_{abc}\left(t\right){\mathbf{i}}_{abc}\left(t\right)\right]+{\omega }_e\left(t\right)\mathbf{e}\left(t\right),$$

(1)

where $\mathbf{R}$ is the resistance matrix, ${\mathbf{L}}_{abc}$ is the inductance matrix, ${\omega }_e$ is the electrical rotor speed, ${\mathbf{v}}_{abc}={\left[\begin{array}{ccc}{v}_{AN}& v_{BN}& v_{CN}\end{array}\right]}^T$ and ${\mathbf{i}}_{abc}={\left[\begin{array}{ccc}{i}_a& i_b& i_c\end{array}\right]}^T$ are, respectively, voltages and current vectors of the phases $A,B,C$, and $\mathbf{e}$ is defined as follows:

$$\mathbf{e}\left(t\right)=-{\lambda }_{PM}\left[\begin{array}{c}sin\left({\theta }_e\left(t\right)\right)\\ sin({\theta }_e\left(t\right)-\frac{2}{3}\pi )\\ sin({\theta }_e\left(t\right)-\frac{4}{3}\pi )\end{array}\right]=\left[\begin{array}{c}{e}_a\\ e_b\\ e_c\end{array}\right],$$

(2)

where ${\theta }_e$ is the electrical rotor position and ${\lambda }_{PM}$ is the back-EMF constant. The schematic of the electrical machine is shown in Figure 1. Let us assume that ${\mathbf{L}}_{abc}$ is a function of ${\theta }_e$ and that it is invertible per each value of ${\theta }_e$. One can write the electrical equation as follows:

$${\mathbf{v}}_{abc}=\mathbf{R}{\mathbf{i}}_{abc}+{\mathbf{L}}_{abc}\frac{\mathrm{d}}{\mathrm{d}t}{\mathbf{i}}_{abc}+{\omega }_e\frac{\partial {\mathbf{L}}_{abc}}{\partial {\theta }_e}{\mathbf{i}}_{abc}+{\omega }_e\mathbf{e},$$

(3)

where the time dependency notation “$\left(t\right)$” has been neglected. The voltage vector ${\mathbf{v}}_{abc}$ is referred to the phase voltages ${\mathbf{v}}_{XO}$ and the star-point voltage $v_{NO}$ in this way: ${\mathbf{v}}_{abc}={\mathbf{v}}_{XO}-v_{NO}{\mathbf{T}}^T$, where:

$${\mathbf{v}}_{XO}=\left[\begin{array}{c}{v}_{AO}\\ v_{BO}\\ v_{CO}\end{array}\right],$$

(4)

and

$$\mathbf{T}=\left[\begin{array}{ccc}1& 1& 1\end{array}\right].$$

(5)

Hence, according to Reference [16], applying the Laplace transform to Equation (3) and after some mathematical manipulation we get:

$$V_{NO}\left(s\right)=\mathbf{G}\left(s\right)\left({\mathbf{V}}_{XO}-{\omega }_e\mathbf{e}\frac{1}{s}\right),$$

(6)

where $\mathbf{G}\left(s\right)$ is a transfer function vector that is dependent on the machine inductance, thus, it contains information about the electrical angular rotor position. Then, it is possible to extract this information directly from $v_{NO}$ by means of a modified phase voltage PWM pattern and a specific star-point voltage measurement procedure.

2.2. Position Extraction from the Star-Point Voltage Measurement

A modified PWM pattern, as shown in Figure 2, has to be generated in order to obtain the DFC signals. The voltage $v_{NV}$ has to be measured for every PWM period two times consecutively before and after $t_1$ inside the highlighted regions. Since this procedure has to be performed for every phase, three PWM period are needed in order to get the DFC signal vector ${\mathbf{\Gamma }}_{abc}$ defined as:

$${\mathbf{\Gamma }}_{abc}={\left[\begin{array}{ccc}{\Gamma }_a& {\Gamma }_b& {\Gamma }_c\end{array}\right]}^T={\left({\mathbf{L}}_{k_{\Sigma }}-\frac{1}{3}\mathbf{T}\right)}^T{v}_{DC},$$

(7)

where:

$${\Gamma }_x=v_{NV}\left(t_{1_x}^{+}\right)-v_{NV}\left(t_{1_x}^{-}\right),x\in \{a,b,c\},$$

(8)

$t_{1_x}$ is the time instant in which the PWM state excitation is switching referred to the phase $a,b\mathrm{or}c$, $v_{DC}$ is the power DC-line voltage, $v_{NV}$ is the star-point voltage referred to the virtual star-point and ${\mathbf{L}}_{k_{\Sigma }}$ is a vector that is expressed as follows:

$${\mathbf{L}}_{k_{\Sigma }}=\left[\begin{array}{ccc}\frac{L_{{\Sigma }_a}}{L_{{\Sigma }_a}+L_{{\Sigma }_b}+L_{{\Sigma }_c}}& \frac{L_{{\Sigma }_b}}{L_{{\Sigma }_a}+L_{{\Sigma }_b}+L_{{\Sigma }_c}}& \frac{L_{{\Sigma }_c}}{L_{{\Sigma }_a}+L_{{\Sigma }_b}+L_{{\Sigma }_c}}\end{array}\right],$$

(9)

where $L_{{\Sigma }_a}$, $L_{{\Sigma }_b}$, $L_{{\Sigma }_c}$ are respectively the sum of the elements along the first, second and third column of the adjugate matrix of ${\mathbf{L}}_{abc}$. To extract the position information, it is convenient to transform this vector into an orthogonal reference system such as $\alpha -\beta -\gamma$. Thus, one can obtain:

$${\mathbf{\Gamma }}_{\alpha \beta \gamma }={\mathbf{T}}_C{\mathbf{\Gamma }}_{abc}={\left[\begin{array}{ccc}{\Gamma }_{\alpha }& {\Gamma }_{\beta }& {\Gamma }_{\gamma }\end{array}\right]}^T,$$

(10)

where ${\mathbf{T}}_C$ is the Clarke transformation matrix. If we consider, for example, the following inductance matrix:

$${\mathbf{L}}_{abc}=\left[\begin{array}{ccc}{L}_{aa}& M_{ab}& M_{bc}\\ M_{ba}& L_{bb}& M_{bc}\\ M_{ca}& M_{cb}& L_{cc}\end{array}\right],$$

(11)

where:

$$\begin{array}{cc}\hfill L_{aa}& =L_0-L_{2}cos\left(2{\theta }_e\right)\hfill \\ \hfill L_{bb}& =L_0-L_{2}cos\left(2\left({\theta }_e-\frac{2\pi }{3}\right)\right)\hfill \\ \hfill L_{cc}& =L_0-L_{2}cos\left(2\left({\theta }_e-\frac{4\pi }{3}\right)\right)\hfill \\ \hfill M_{ab}& =M_{ba}=M_0-M_{2}cos\left(2\left({\theta }_e-\frac{4\pi }{3}\right)\right)\hfill \\ \hfill M_{bc}& =M_{cb}=M_0-M_{2}cos\left(2{\theta }_e\right)\hfill \\ \hfill M_{ca}& =M_{ac}=M_0-M_{2}cos\left(2\left({\theta }_e-\frac{2\pi }{3}\right)\right).\hfill \end{array}$$

(12)

The expression of the vector signal ${\mathbf{\Gamma }}_{\alpha \beta \gamma }$ takes the following form:

$$\begin{array}{cc}\hfill {\Gamma }_{\alpha }& =acos\left(2{\theta }_e\right)+bcos\left(4{\theta }_e\right)\hfill \\ \hfill {\Gamma }_{\beta }& =-asin\left(2{\theta }_e\right)+bsin\left(4{\theta }_e\right)\hfill \\ \hfill {\Gamma }_{\gamma }& =0,\hfill \end{array}$$

(13)

where:

$$\begin{array}{cc}\hfill a& =\frac{\left(L_2-M_2\right)\left(L_0-M_0\right)}{3\left({\left(L_0-M_0\right)}^2-{\left(\frac{L_2}{2}+M_2\right)}^2\right)}{v}_{DC}\hfill \\ \hfill b& =\frac{\left(L_2-M_2\right)\left(\frac{L_2}{2}+M_2\right)}{3\left({\left(L_0-M_0\right)}^2-{\left(\frac{L_2}{2}+M_2\right)}^2\right)}{v}_{DC}\hfill \end{array}$$

(14)

Thus, substituting the term $2{\theta }_e$ with x we can express the double valued angular position estimation as follows:

$$\widehat{x}=-arctan\left({\displaystyle \frac{{\Gamma }_{\beta }}{{\Gamma }_{\alpha }}}\right)\approx x.$$

(15)

Starting from this point the term second harmonic will refer to the term $2x$. The position estimation procedure is summarized in Figure 3. The signals ${\Gamma }_{\alpha }$ and ${\Gamma }_{\beta }$ can be generalized to fit the results of any other anisotropy-based sensorless technique. One can notice that the variable $\widehat{x}$ is equal to x only if the parameter b is zero. However, if b is not zero, we can determine the estimation error $\Delta$ as follows:

$$\begin{array}{cc}\hfill \widehat{x}& =-arctan\left({\displaystyle \frac{-asin\left(x\right)+bsin\left(2x\right)}{acos\left(x\right)+bcos\left(2x\right)}}\right)\hfill \\ & =arctan\left({\displaystyle \frac{asin\left(x\right)-bsin\left(2x\right)}{acos\left(x\right)+bcos\left(2x\right)}}\right)\hfill \\ & =\angle \left(a{e}^{jx}+b{e}^{-j2x}\right)\hfill \\ & =\angle \left(a{e}^{jx}\left(1+{\displaystyle \frac{b}{a}}{e}^{-3jx}\right)\right)\hfill \\ & =x+\angle \left(1+{\displaystyle \frac{b}{a}}{e}^{-j3x}\right)\hfill \\ & =x-arctan\left({\displaystyle \frac{psin\left(3x\right)}{1+pcos\left(3x\right)}}\right)\hfill \\ & =x+\Delta \hfill \end{array}$$

(16)

where:

$$\begin{array}{cc}\hfill p& ={\displaystyle \frac{b}{a}}\hfill \\ \hfill \Delta & =-arctan\left({\displaystyle \frac{psin\left(3x\right)}{1+pcos\left(3x\right)}}\right).\hfill \end{array}$$

(17)

The estimation error term $\Delta$ is a periodic function of the electrical position. Even when the value of p is small enough allowing a correct drive of the electrical machine, the estimation of other variable such the electrical angular speed can be challenging since a third order harmonic is introduced on the estimated position. The aim is to eliminate or at least reduce as much as possible the error between the estimated and the real position (see Figure 4).

2.3. Inverse DFC Algorithm

A method for the elimination of the deviation was proposed in Reference [17]. The signals vector ${\mathbf{L}}_{k_{\Sigma }}$ is used instead of ${\mathbf{\Gamma }}_{abc}$ and its values are inverted:

$${\mathbf{\eta }}_{abc}=\left[\begin{array}{ccc}\frac{L_{{\Sigma }_a}+L_{{\Sigma }_b}+L_{{\Sigma }_c}}{L_{{\Sigma }_a}}& \frac{L_{{\Sigma }_a}+L_{{\Sigma }_b}+L_{{\Sigma }_c}}{L_{{\Sigma }_b}}& \frac{L_{{\Sigma }_a}+L_{{\Sigma }_b}+L_{{\Sigma }_c}}{L_{{\Sigma }_c}}\end{array}\right].$$

(18)

${\mathbf{\eta }}_{abc}$ is then transformed through the Clarke matrix into an $\alpha -\beta -\gamma$ frame. The angular position can be estimated as usual:

$$\widehat{x}=-arctan\left({\displaystyle \frac{{\eta }_{\beta }}{{\eta }_{\alpha }}}\right)\approx x.$$

(19)

Referring to Equation (14), it can be proven that Equation (19) gives the exact position when $a^2=b^2+b$. It easy to verify that this condition is fulfilled for $M_2=0$.

3. Fourth Harmonic Elimination: Ivd-Dfc Algorithm

Within this section the static iteration algorithm (IVD-DFC) is presented. First, the algorithm is presented when no saturation effects are considered on the model. In this case, the algorithm requires only the knowledge of the parameter b. It has to be remarked that in this case the algorithm not only requires the knowledge of the parameter b but also other parameters that will be introduced within this section.

3.1. Case of Study: No Saturation

As shown, the variable $\widehat{x}$ is deviated from x by the term $\Delta$ but we can consider it a good approximation of x if the value of p is not too large. If we suppose to know the parameter b one may try to reduce the effect of the second harmonic term of ${\Gamma }_{\beta }$ and ${\Gamma }_{\alpha }$ by algebraically subtracting the second harmonic on ${\Gamma }_{\beta }$ and ${\Gamma }_{\alpha }$ themselves. Let us suppose we are working at a time instant. First, we include in our notation a numbered subscript as index for the iterations. Thus, the value of the initial estimation is ${\widehat{x}}_0$ and it is calculated as in Equation (16). Then, the initial estimation obtained is used for the first iteration:

$${\widehat{x}}_1=-arctan\left({\displaystyle \frac{{\Gamma }_{\beta }-bsin\left(2{\widehat{x}}_0\right)}{{\Gamma }_{\alpha }-bcos\left(2{\widehat{x}}_0\right)}}\right).$$

(20)

The idea is to reduce the effect of the second harmonic using the estimated $\widehat{x}$ previously calculated. As it will be shown, the Equation (20) has another deviation term ${\Delta }_1$. If we find that $\parallel {\Delta }_1\parallel <\parallel {\Delta }_0\parallel$ for some value of p then we can construct the following algorithm:

$${\widehat{x}}_k=-arctan\left({\displaystyle \frac{{\Gamma }_{\beta }-bsin\left(2{\widehat{x}}_{k-1}\right)}{{\Gamma }_{\alpha }-bcos\left(2{\widehat{x}}_{k-1}\right)}}\right),$$

(21)

for $k=1,2,\dots ,n$ and we want to see if the following limit holds:

$$\underset{n\to \infty }{lim}{\widehat{x}}_n=x.$$

(22)

The IVD-DFC algorithm consists actually in Equation (21) and its schematic is shown in Figure 5. The first operation is to extract the raw estimated position ${\widehat{x}}_0$ using Equation (15). Once the first information is obtained, the iterations Equation (21) can be performed. A new estimation of the variable x is given after n iterations. Thus, a large number of iterations can be problematic for a low-cost microcontroller to be executed in real time. Anyway, as shown in Figure 6 and Figure 7, even a single iteration can already considerably improve the position estimation.

Let us now analyze the proposed algorithm by looking at the first iteration. The Equation (20) can be manipulated in the following manner:

$$\begin{array}{cc}\hfill {\widehat{x}}_1& =arctan\left({\displaystyle \frac{asin\left(x\right)-bsin\left(2x\right)+bsin\left(2{\widehat{x}}_0\right)}{acos\left(x\right)+bcos\left(2x\right)-bcos\left(2{\widehat{x}}_0\right)}}\right)\hfill \\ & =\angle \left(a{e}^{jx}+b{e}^{-j2x}-b{e}^{-j2{\widehat{x}}_0}\right)\hfill \\ & =x+\angle \left(1+p{e}^{-j3x}-p{e}^{-j\left(2{\widehat{x}}_0+x\right)}\right)\hfill \\ & =x+\angle \left(1+p{e}^{-j3x}-p{e}^{-j\left(3x+2{\Delta }_0\right)}\right)\hfill \\ & =x+\angle \left(1+p\left(1-e^{-j2{\Delta }_0}\right)e^{-j3x}\right).\hfill \end{array}$$

(23)

Let us consider the following proposition:

Proposition 1.

The sum of two unitary complex exponential numbers $\mathsf{\Phi }=e^{j{\gamma }_1}+e^{j{\gamma }_2}$ can be written as: $\mathsf{\Phi }=2cos\left({\displaystyle \frac{{\gamma }_2-{\gamma }_1}{2}}\right)e^{j\left({\displaystyle \frac{{\gamma }_2+{\gamma }_1}{2}}\right)}$.

From Equation (23) using the previous proposition and defining ${\gamma }_2=-2{\Delta }_0-\pi$ and ${\gamma }_1=0$ one gets:

$$\begin{array}{cc}\hfill {\widehat{x}}_1& =x+\angle \left(1+2pcos\left({\Delta }_0+{\displaystyle \frac{\pi }{2}}\right)e^{-j\left(3x+{\Delta }_0+{\displaystyle \frac{\pi }{2}}\right)}\right)\hfill \\ & =x+\angle \left(1-2psin\left({\Delta }_0\right)e^{-j\left(3x+{\Delta }_0+{\displaystyle \frac{\pi }{2}}\right)}\right)\hfill \\ & =x+arctan\left({\displaystyle \frac{2psin\left({\Delta }_0\right)sin\left(3x+{\Delta }_0+{\displaystyle \frac{\pi }{2}}\right)}{1-2psin\left({\Delta }_0\right)cos\left(3x+{\Delta }_0+{\displaystyle \frac{\pi }{2}}\right)}}\right)\hfill \\ & =x+arctan\left({\displaystyle \frac{2psin\left({\Delta }_0\right)cos\left(3x+{\Delta }_0\right)}{1+2psin\left({\Delta }_0\right)sin\left(3x+{\Delta }_0\right)}}\right)\hfill \\ & =x+{\Delta }_1.\hfill \end{array}$$

(24)

A relation between ${\Delta }_k$ and its previous value ${\Delta }_{k-1}$ can be found from Equation (24):

$${\Delta }_k=arctan\left({\displaystyle \frac{2psin\left({\Delta }_{k-1}\right)cos\left(3x+{\Delta }_{k-1}\right)}{1+2psin\left({\Delta }_{k-1}\right)sin\left(3x+{\Delta }_{k-1}\right)}}\right).$$

(25)

A graphical visualization of the iterations effect is shown in Figure 8. Next step is to prove that the deviation $\left|{\Delta }_k\right|$ is convergent to zero, namely:

$$\underset{k\to \infty }{lim}\left|{\Delta }_k\right|=0.$$

(26)

Let us introduce the following definition:

Definition 1.

Let be given a real value sequence $X_k$. The sequence is said to be linearly convergent (at least) to zero if exist a value $\sigma \in \left[0,1\right)$ such that $\left|X_k\right|\le \sigma \left|X_{k-1}\right|$$\forall k\in \left[1,2,\dots ,n\right]$. Then, ${lim}_{k\to \infty }\left|X_k\right|=0.$

To prove (26) the following lemma is presented.

Lemma 1.

Let us be given a real value sequence $X_k$ defined on the interval $\left[-{\displaystyle \frac{\pi }{2}},{\displaystyle \frac{\pi }{2}}\right]$. If $tan\left(X_k\right)$ is at least linearly convergent to zero, then also the original sequence $X_k$ is at least linearly convergent to the same value.

The previous lemma will be useful for the proof of the following theorem.

Theorem 1.

Let be the following initial function: ${\Delta }_0\left(x,p\right)=-arctan\left({\displaystyle \frac{psin\left(3x\right)}{1+pcos\left(3x\right)}}\right)$ parameterized by $p\in \mathbb{R}$ and defined $\forall x\in \mathbb{R}$. If the following sequence:

${\Delta }_k\left(x,p\right)=arctan\left({\displaystyle \frac{2psin\left({\Delta }_{k-1}\left(x,p\right)\right)cos\left(3x+{\Delta }_{k-1}\left(x,p\right)\right)}{1+2psin\left({\Delta }_{k-1}\left(x,p\right)\right)sin\left(3x+{\Delta }_{k-1}\left(x,p\right)\right)}}\right)$ is defined $\forall k\in {\mathbb{N}}^0$. Then the sequence ${\Delta }_k$ is linearly convergent if and only if $\left|p\right|<{\displaystyle \frac{1}{2}}$.

Proof of Theorem 1.

Let us proof first the following statement:

$$\left|tan\left({\Delta }_k\left(x,p\right)\right)\right|\le \gamma \left(p\right)\left|tan\left({\Delta }_{k-1}\left(x,p\right)\right)\right|,$$

(27)

where:

$$\gamma \left(p\right)=2\left|p\right|.$$

(28)

It is trivial to see that the previous sequence is linearly convergent only if $\gamma \left(p\right)<1$ or $\left|p\right|<{\displaystyle \frac{1}{2}}$. (necessary condition)

Let us elaborate the equations considering $\left|p\right|<{\displaystyle \frac{1}{2}}$. (sufficient condition)

First, let us develop the left side of the inequality (27), for simplicity we put ${\Delta }_k\left(x,p\right)={\Delta }_k$:

$$\begin{array}{cc}\hfill \left|tan\left({\Delta }_k\right)\right|& =\left|tan\left(arctan\left({\displaystyle \frac{2psin\left({\Delta }_{k-1}\right)cos\left(3x+{\Delta }_{k-1}\right)}{1+2psin\left({\Delta }_{k-1}\right)sin\left(3x+{\Delta }_{k-1}\right)}}\right)\right)\right|\hfill \\ & =\left|{\displaystyle \frac{2psin\left({\Delta }_{k-1}\right)cos\left(3x+{\Delta }_{k-1}\right)}{1+2psin\left({\Delta }_{k-1}\right)sin\left(3x+{\Delta }_{k-1}\right)}}\right|\hfill \\ & ={\displaystyle \frac{\left|2psin\left({\Delta }_{k-1}\right)cos\left(3x+{\Delta }_{k-1}\right)\right|}{\left|1+2psin\left({\Delta }_{k-1}\right)sin\left(3x+{\Delta }_{k-1}\right)\right|}}\hfill \\ & ={\displaystyle \frac{2\left|p\right|\left|sin\left({\Delta }_{k-1}\right)\right|\left|cos\left(3x+{\Delta }_{k-1}\right)\right|}{\left|1+2psin\left({\Delta }_{k-1}\right)sin\left(3x+{\Delta }_{k-1}\right)\right|}}.\hfill \end{array}$$

(29)

Then, substituting Equation (29) in Equation (27), we obtain:

$${\displaystyle \frac{2\left|p\right|\left|sin\left({\Delta }_{k-1}\right)\right|\left|cos\left(3x+{\Delta }_{k-1}\right)\right|}{\left|1+2psin\left({\Delta }_{k-1}\right)sin\left(3x+{\Delta }_{k-1}\right)\right|}}\le 2\left|p\right|{\displaystyle \frac{\left|sin\left({\Delta }_{k-1}\right)\right|}{\left|cos\left({\Delta }_{k-1}\right)\right|}}$$

(30)

The previous inequality will be used in order to find the values of p that satisfies the inequality. The terms can be elaborated as follows:

$$\begin{array}{cc}\hfill & \left|cos\left({\Delta }_{k-1}\right)\right|{\displaystyle \frac{\bcancel{2\left|p\right|}\bcancel{\left|sin\left({\Delta }_{k-1}\right)\right|}\left|cos\left(3x+{\Delta }_{k-1}\right)\right|}{\left|1+2psin\left({\Delta }_{k-1}\right)sin\left(3x+{\Delta }_{k-1}\right)\right|}}\le \bcancel{2\left|p\right|}\bcancel{\left|sin\left({\Delta }_{k-1}\right)\right|}\hfill \\ & {\displaystyle \frac{\left|cos\left({\Delta }_{k-1}\right)\right|\left|cos\left(3x+{\Delta }_{k-1}\right)\right|}{\left|1+2psin\left({\Delta }_{k-1}\right)sin\left(3x+{\Delta }_{k-1}\right)\right|}}\le 1\hfill \\ & |\stackrel{{\chi }_1\left(x,{\Delta }_{k-1}\right)}{\overbrace{cos\left({\Delta }_{k-1}\right)cos\left(3x+{\Delta }_{k-1}\right)}}|\le |\stackrel{{\chi }_2\left(x,{\Delta }_{k-1}\right)}{\overbrace{1+2psin\left({\Delta }_{k-1}\right)sin\left(3x+{\Delta }_{k-1}\right)}}|.\hfill \end{array}$$

(31)

We can see that for $\left|p\right|<{\displaystyle \frac{1}{2}}$ the statement ${\chi }_2\left(x,{\Delta }_{k-1}\right)>0$ is always true.

In the first case: ${\chi }_1\left(x,{\Delta }_{k-1}\right)\ge 0$:

$$\begin{array}{cc}\hfill cos\left({\Delta }_{k-1}\right)cos\left(3x+{\Delta }_{k-1}\right)& \le 1+2psin\left({\Delta }_{k-1}\right)sin\left(3x+{\Delta }_{k-1}\right)\hfill \\ \hfill cos\left({\Delta }_{k-1}\right)cos\left(3x+{\Delta }_{k-1}\right)-2psin\left({\Delta }_{k-1}\right)sin\left(3x+{\Delta }_{k-1}\right)& \le 1\hfill \\ \hfill \left(1-2p\right)cos\left(-3x\right)+\left(1+2p\right)cos\left(3x+2{\Delta }_{k-1}\right)& \le 2.\hfill \end{array}$$

(32)

We can rewrite the last inequality as follows:

$$c_{1}cos\left(y_1\right)+c_{2}cos\left(y_2\right)\le 2,$$

(33)

with $c_1=1-2p$, $c_2=1+2p$, $y_2=-3x$ and $y_1=3x+2{\Delta }_{k-1}$.

One can state that:

$$c_{1}cos\left(y_1\right)+c_{2}cos\left(y_2\right)\le \left|c_{1}cos\left(y_1\right)+c_{2}cos\left(y_2\right)\right|.$$

(34)

Then, applying the triangular inequality statement, we can write:

$$\left|c_{1}cos\left(y_1\right)+c_{2}cos\left(y_2\right)\right|\le \left|c_{1}cos\left(y_1\right)\right|+\left|c_{2}cos\left(y_2\right)\right|\le \left|c_1\right|+\left|c_2\right|.$$

(35)

Thus, putting together Equations (34) and (35):

$$c_{1}cos\left(y_1\right)+c_{2}cos\left(y_2\right)\le \left|c_1\right|+\left|c_2\right|.$$

(36)

Thus, if we find a range of values of p that satisfies $\left|c_1\right|+\left|c_2\right|\le 2$ the Equation (33) is also satisfied. Let us define the range for the p values considering $\left|p\right|\le {\displaystyle \frac{1}{2}}$:

$$\left|c_1\right|+\left|c_2\right|\le 2$$

(37)

Since $c_1$ and $c_2$ are positive for $\left|p\right|\le {\displaystyle \frac{1}{2}}$, we can avoid the absolute value notation:

$$\begin{array}{cc}\hfill c_1+c_2& \le 2\hfill \\ \hfill 1-2p+1+2p& \le 2\hfill \\ \hfill 2& \le 2.\hfill \end{array}$$

(38)

That is feasible. Then, $\left|p\right|\le {\displaystyle \frac{1}{2}}$ is a possible range in order for Equation (30) to be satisfied. Thus, $\left|p\right|\le {\displaystyle \frac{1}{2}}$ is a sufficient condition for Equation (30).

The case ${\chi }_1\left(x,{\Delta }_{k-1}\right)\le 0$ can be manipulated as before and the result is the same. Then, $\left|p\right|\le {\displaystyle \frac{1}{2}}$ is a possible range for Equation (30) in order to be satisfied.

Then, if we put together the necessary and sufficient conditions found, one can say that:

$$\underset{k\to \infty }{lim}\left|tan\left({\Delta }_k\right)\right|=0,\mathrm{iff}\left|p\right|<{\displaystyle \frac{1}{2}}.$$

(39)

Thus, using the Lemma 1:

$$\underset{k\to \infty }{lim}\left|{\Delta }_k\right|=0,\mathrm{iff}\left|p\right|<{\displaystyle \frac{1}{2}}.$$

(40)

□

3.2. Case of Study: With Saturation

If saturation effects are taken into account the equations in (13) can be modified as follow:

$$\begin{array}{cc}\hfill {\Gamma }_{\alpha }& =acos\left(x+{\phi }_a\right)+bcos\left(2x+{\phi }_b\right)\hfill \\ \hfill {\Gamma }_{\beta }& =-asin\left(x+{\phi }_a\right)+bsin\left(2x+{\phi }_b\right),\hfill \end{array}$$

(41)

where the value of ${\phi }_a$ and ${\phi }_b$ depend on the currents that are flowing through the windings. The parameter a and b also vary according to the current amplitudes. An estimation of the variable x can be given by:

$${\widehat{x}}_0=-arctan\left({\displaystyle \frac{{\Gamma }_{\beta }}{{\Gamma }_{\alpha }}}\right)-{\phi }_a.$$

(42)

The following mathematical elaboration can be done in order to find the deviation ${\Delta }_0$, the process is similar to Equation (16):

$$\begin{array}{cc}\hfill {\widehat{x}}_0& =-arctan\left({\displaystyle \frac{-asin\left(x+{\phi }_a\right)+bsin\left(2x+{\phi }_b\right)}{acos\left(x+{\phi }_a\right)+bcos\left(2x+{\phi }_b\right)}}\right)-{\phi }_a\hfill \\ & =arctan\left({\displaystyle \frac{asin\left(x+{\phi }_a\right)-bsin\left(2x+{\phi }_b\right)}{acos\left(x+{\phi }_a\right)+bcos\left(2x+{\phi }_b\right)}}\right)-{\phi }_a\hfill \\ & =\angle \left(a{e}^{j\left(x+{\phi }_a\right)}+b{e}^{-j\left(2x+{\phi }_b\right)}\right)-{\phi }_a\hfill \\ & =x+{\phi }_a+\angle \left(1+{\displaystyle \frac{b}{a}}{e}^{-j\left(3x+{\phi }_a+{\phi }_b\right)}\right)-{\phi }_a\hfill \\ & =x-arctan\left({\displaystyle \frac{psin\left(3x+{\phi }_a+{\phi }_b\right)}{1+pcos\left(3x+{\phi }_a+{\phi }_b\right)}}\right)\hfill \\ & =x+{\Delta }_0.\hfill \end{array}$$

(43)

Let us suppose that we know the parameter b, ${\phi }_a$ and ${\phi }_b$, as before one can build the algorithm:

$${\widehat{x}}_k=-arctan\left({\displaystyle \frac{{\Gamma }_{\beta }-bsin\left(2{\widehat{x}}_{k-1}+{\phi }_b\right)}{{\Gamma }_{\alpha }-bcos\left(2{\widehat{x}}_{k-1}+{\phi }_b\right)}}\right)-{\phi }_a.$$

(44)

That is, for the first iteration:

$$\begin{array}{cc}\hfill {\widehat{x}}_1& =arctan\left({\displaystyle \frac{asin\left(x+{\phi }_a\right)-bsin\left(2x+{\phi }_b\right)+bsin\left(2{\widehat{x}}_0+{\phi }_b\right)}{acos\left(x+{\phi }_a\right)+bcos\left(2x+{\phi }_b\right)-bcos\left(2{\widehat{x}}_0+{\phi }_b\right)}}\right)-{\phi }_a\hfill \\ & =\angle \left(a{e}^{j\left(x+{\phi }_a\right)}+b{e}^{-j\left(2x+{\phi }_b\right)}-b{e}^{-j\left(2{\widehat{x}}_0+{\phi }_b\right)}\right)-{\phi }_a\hfill \\ & =x+{\phi }_a+\angle \left(1+p{e}^{-j\left(3x+{\phi }_a+{\phi }_b\right)}-p{e}^{-j\left(2{\widehat{x}}_0+x+{\phi }_a+{\phi }_b\right)}\right)-{\phi }_a\hfill \\ & =x+\angle \left(1+p{e}^{-j\left(3x+{\phi }_a+{\phi }_b\right)}-p{e}^{-j\left(3x+2{\Delta }_0+{\phi }_a+{\phi }_b\right)}\right)\hfill \\ & =x+\angle \left(1-2psin\left({\Delta }_0\right)e^{-j\left(3x+{\Delta }_0+{\phi }_a+{\phi }_b+\pi /2\right)}\right)\hfill \\ & =x+{\Delta }_1.\hfill \end{array}$$

(45)

Thus:

$${\Delta }_k=arctan\left({\displaystyle \frac{2psin\left({\Delta }_{k-1}\right)cos\left(3x+{\Delta }_{k-1}+{\phi }_a+{\phi }_b\right)}{1+2psin\left({\Delta }_{k-1}\right)sin\left(3x+{\Delta }_{k-1}+{\phi }_a+{\phi }_b\right)}}\right).$$

(46)

The mathematics follows the same procedure as in the previous subsection. Using the previous theorem we can prove the following limit:

$$\underset{k\to \infty }{lim}\left|{\Delta }_k\right|=0,$$

(47)

and:

$$\underset{k\to \infty }{lim}{\widehat{x}}_k=x.$$

(48)

On Figure 9 the flowchart related to the IVD-DFC operations in case of saturation effect is presented.

3.3. Sensitivity Analysis

If the parameters b, ${\phi }_a$ and ${\phi }_b$ are not precisely known the following algorithm can be proposed:

$${\widehat{x}}_k=-arctan\left({\displaystyle \frac{{\Gamma }_{\beta }-\widehat{b}sin\left(2{\widehat{x}}_{k-1}+{\widehat{\phi }}_b\right)}{{\Gamma }_{\alpha }-\widehat{b}cos\left(2{\widehat{x}}_{k-1}+{\widehat{\phi }}_b\right)}}\right)-{\widehat{\phi }}_a.$$

(49)

The previous equation can be easily related to (44). The position information x can be extracted from Equation (49) as presented in the previous subsections. After some mathematical elaborations we obtain:

$${\widehat{x}}_k=x+{\tilde{\phi }}_a+arctan\left({\displaystyle \frac{2psin\left({\Delta }_{k-1}+{\tilde{\phi }}_b\right)cos\left(3x+{\kappa }_1\right)+\tilde{p}sin\left(3x+{\kappa }_2\right)}{1-2psin\left({\Delta }_{k-1}+{\tilde{\phi }}_b\right)sin\left(3x+{\kappa }_1\right)+\tilde{p}cos\left(3x+{\kappa }_2\right)}}\right),$$

(50)

where:

$$\begin{array}{cc}\hfill {\kappa }_1& ={\phi }_a+\frac{{\phi }_b+{\widehat{\phi }}_b}{2}+{\Delta }_{k-1}\hfill \\ \hfill {\kappa }_2& ={\phi }_a+{\widehat{\phi }}_b+2{\Delta }_{k-1}\hfill \\ \hfill \tilde{p}& ={\displaystyle \frac{b-\widehat{b}}{a}}\hfill \\ \hfill {\tilde{\phi }}_a& ={\phi }_a-{\widehat{\phi }}_a\hfill \\ \hfill {\tilde{\phi }}_b& ={\phi }_b-{\widehat{\phi }}_b,\hfill \end{array}$$

(51)

where $\widehat{b}$, ${\widehat{\phi }}_a$ and ${\widehat{\phi }}_b$ are respectively the estimated value of b, ${\phi }_a$ and ${\phi }_b$.

From Equation (50) we cannot state that for increasing number of iteration the deviation term goes to zero. In fact, considering the results of Figure 10, a greater number of iteration does not mean a better convergence if the p value is too small. Figure 11 shows the deterioration of the performance due to the phase estimation error.

4. Experimental Results

Within this section, the IVD-DFC algorithm will be validated by means of experimental tests. It has been shown that the new algorithm does not need any additional dynamical system to be implemented. The parameters to be set for the correct functioning of the algorithm are: the amplitude of the fourth harmonic (b), the phase shifting due to the saturation effects (${\phi }_a$ and ${\phi }_b$) and the number of iterations. The first two parameters can be normally identified using either online or offline identification methods, the third one should be chosen considering a trade-off between accuracy (more iterations) and computational effort (less iterations). In fact, a large number of iterations can result in a considerable computing effort for the microcontroller. Nevertheless, it has to be remarked that several approaches can be adopted to optimize the calculation of trigonometric functions either in software or in hardware. Moreover, as it is shown in this work, usually a relatively low number of terations is required in order to allow the algorithm to converge to an acceptable value.

4.1. Test Setup

For the experimental validation, a test bench composed of a custom PMSM coupled to a servo motor and a Baumer GBA2H 18-bit encoder has been considered (see Figure 12, parameters of the PMSM listed in

Table 1. Custom PMSM parameter list.

Parameter	Value
Phase Resistance	1.1 $\Omega$
$L_d$ inductance	394 µH
$L_q$ inductance	475 µH
magnetic flux linkage ${\lambda }_{PM}$	9.89 mVs
Number of pole pairs	8
Nominal voltage	24 V

).

A specific electronic board has been used for the purpose. That board features a 32-bit microcontroller, a three-phase inverter bridge, a dedicated electronic for the star-point measurement and a USB communication port (see Figure 13). The parameters needed for the IVD-DFC have been identified offline for the complete range of the q-axis current using the MATLAB System Identification Toolbox. Different tests were performed in order to prove the correct algorithm functioning during typical operations and under stress conditions.

4.2. DFC and IVD-DFC Signals Comparison

Within this subsection the measured DFC signals are compared to the IVD-DFC ones when the motor is operating at nominal speed without external load torque applied on the shaft. The saliency effects of the examined motor are not considerably large and a single iteration is enough to eliminate almost entirely their component. As shown in Figure 14, the DFC signals ${\Gamma }_{\alpha }$ and ${\Gamma }_{\beta }$ are forming a circular 2-D plot where the saliency component can easily be seen compared to the reference circle in black, that represent the ideal case with no multiple harmonics.

Two IVD-DFC are compared using respectively one and two iterations. One can notice that the two IVD-DFC circles are more coherent with the reference circle than the DFC one. The second order harmonic seems to be almost completely suppressed. The thickness of the signals does not depend on higher order saliencies but on the unideal magnetic characteristic of the motor. A single iteration seems to present already a good response in terms of second saliency harmonic elimination. The FFT amplitude response calculated for both DFC signals and IVD-DFC signals with one iteration is presented in Figure 15.

As shown, two amplitude peaks are prevalent among the whole frequency spectrum. They represent respectively the first and the second saliency harmonics. The peak of the second harmonic is reduced using the IVD-DFC about more than 80% of the original value. The FFT shows small amplitude harmonics before the first peaks. We suppose that these subharmonics are responsible for the imperfect estimation of the angular position (see the thickness of the signals in Figure 14).

4.3. Position and Speed Estimation

In the previous subsection, we focused on the obtained $\alpha$ and $\beta$ signals and we proved the functioning of the IVD-DFC algorithm. Within this subsection, we present the comparison between the encoder information and the estimated angular position and speed using DFC and IVD-DFC.

The whole presented results are referring to one step iteration for the IVD-DFC. In order to prove the correct functioning of the IVD-DFC algorithm with saturation effects, a considerable amount of current into the q-axis has been generated. As already said, the currents change the value of the parameters a, b, ${\phi }_a$ and ${\phi }_b$, generating signals in the presented form in Equation (41). A previous identification of these parameters has been performed and look-up tables are used in this work to feed the algorithm with the correct parameter values. In Figure 16 the DFC estimated position error is compared with the IVD-DFC one as increasing amplitude steps of current on the q-axis are applied. The saturation effect introduces in the estimated position an offset of ${\phi }_a$ and, at the same time, modifies ${\Delta }_k$ since this depends both on ${\phi }_a$ and ${\phi }_b$. Actually, the shifted position can be adjusted subtracting ${\phi }_a$ from the estimation. Anyhow, the ripple around the mean value of the position error can be reduced only if both parameter ${\phi }_a$ and ${\phi }_b$ are known.

In Figure 17, the estimated DFC and IVD-DFC angular speed are compared with the encoder speed when opposite directions of speed reference are given for the speed control of the machine. The test is performed using the encoder as angular position information for the drive of the machine. The standard deviation of the DFC speed signal in respect to the encoder signal is more than the double of the IVD-DFC one. Therefore, the use of the IVD-DFC algorithm reduces the speed signal noise but not the bandwidth of the estimation. Thus, the low-pass filter used for the DFC speed signal can be tuned differently in order to have a wider bandwidth granting a better speed control performance. In Figure 18, the encoder angular position is compared to the DFC and IVD-DFC estimation. It can be easily seen that the IVD-DFC line follows a straight line over the encoder position.

5. Conclusions

This paper presented a new algorithm for the elimination of the second order saliency component from anisotropy-based sensorless techniques. The new method avoids the use of an additional dynamic system for the estimation of angular speed and position of a PMSM. A thoroughly analysis of the algorithm was presented, the convergence condition has been proved and experimental results have been shown in order to validate the concept. In our specific case, the IVD algorithm was proposed for the improvement of the DFC signals. Anyway, every anisotropy-based sensorless technique using two-dimensional anisotropy information could be theoretically used for the same purpose, since the method is based on the measured signals and not on the specific technique used to extract them. The functioning of the algorithm was also proven under stress conditions, such as high-current operation, in order to bring the machine to saturation and switching of the speed control reference between opposite value of nominal speed in order to test the robustness of the algorithm. The IVD-DFC method is able to eliminate almost completely the second order saliency effect on the considered electrical machine. Moreover, the noise over the estimated speed signal has been reduced more than 80% using the IVD-DFC compared to DFC. The IVD-DFC approaches the problem of the position shifting in a different way than the Lorentz method [7]. The main limitation of the presented algorithm is the inability to correct signals whose second order saliency harmonics are greater than the half of the first harmonic (that is $p>0.5$). Anyway, solutions to this problem have been already presented and we suppose that any of these solutions could be applied to the examined case. Furthermore, the research on online parameter estimation for b, ${\phi }_a$ and ${\phi }_b$ to be combined with IVD-DFC is ongoing.