# FPGA Realization of Sensorless PMSM Speed Controller Based on Extended Kalman Filter

Authors: Ying-Shieh Kung; Nguyen Vu Quynh; Nguyen Trung Hieu; Jin-Mu Lin
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

Based on extended Kalman filter (EKF), the design and FPGA implementation of a sensorless control intellectual property (IP)  for permanent magnet synchronous motor (PMSM) drive are presented in this paper. Firstly, the mathematical model for PMSM is derived and the vector control is built up. Secondly, the rotor flux angle (FA) and rotor speed which are estimated by using EKF are described. These estimated values are feedbacked to the current loop for vector control and to the speed loop for speed control. Thirdly, the very-high-speed IC hardware description language (VHDL) is adopted to describe the behavior of the sensorless speed control IP which includes the circuits of space vector pulse width modulation (SVPWM), coordinate transformation, EKF, and PI controller. Finally, to evaluate the effectiveness and correctness of the proposed of system, a cosimulation work performed by Simulink and ModelSim is firstly conducted. Then, an experimental system by FPGA chip and motor driving board is set up to further validate the performance of the proposed EKF-based sensorless speed control IP.


---

## Body

## 1. Introduction

PMSM has been increasingly used in many automation control fields as actuators, due to its advantages of superior power density, high-performance motion control with fast speed, and better accuracy. However, conventional motor control needs a speed sensor or an optical encoder to measure the rotor speed and feedback it to the controller for ensuring the precision speed control. Such sensor presents some disadvantages such as drive cost, machine size, reliability, and noise immunity. In recent years, a sensorless control without position and speed sensors for PMSM drive become a popular research topic in the literature [1–7]. Those sensorless control strategies have sliding mode observer (SMO), extended Kalman filter (EKF), reduced-order EKF, and so forth. The EKF is basically a full-order stochastic observer for the recursive optimum state estimation of a nonlinear dynamic system in real time by using signals that are in noisy environment [7]. Comparing with SMO, EKF can directly estimate the angular speed and it has high convergence rate which can give a more rapid speed response [5]. However, EKF requires heavy online 4×4 matrix computing; therefore, the complex computation becomes a challenge for a fixed-point processor system. In realization, a fixed-point processor using digital signal processor (DSP) and field programmable gate array (FPGA) both can provide a good solution in this issue. Particularly, FPGA with programmable hard-wired feature, fast computation ability, shorter design cycle, embedding processor, low power consumption, and higher density is better for the implementation of the digital system [8, 9] than DSP.Recently, a cosimulation work by Electronic Design Automation (EDA) Simulator Link has been gradually applied to verify the effectiveness of the Verilog and VHDL code in the motor drive system [10–12]. The EDA Simulator Link [13] provides a cosimulation interface between MATLAB or Simulink and HDL simulators-ModelSim [14]. Using it you can verify a VHDL, Verilog, or mixed-language implementation against your Simulink model or MATLAB algorithm [13]. Therefore, EDA Simulator Link lets you use MATLAB code and Simulink models as a test bench that generates stimulus for an HDL simulation and analyzes the simulation’s response [13]. In this paper, a cosimulation by EDA Simulator Link is applied to sensorless speed control for PMSM drive. The PMSM, inverter, and speed command are performed in Simulink and the EKF algorithm; current vector controller and speed PI controller described by VHDL code are executed in ModelSim. After successful verification in simulation, an FPGA-based experimental system is established with configuration as in Figure 1 for realizing the proposed IP code again, and experiment results will validate the effectiveness of the sensorless speed control system of PMSM drive.Figure 1
Architecture of FPGA-based sensorless PMSM drive system.


## 2. System Description of PMSM Drive and Sensorless Speed Control

The sensorless speed control block diagram for PMSM drive is shown in Figure1. The modeling of PMSM, the EKF-based rotor flux angle (FA) and rotor speed estimation, and the PI controller are introduced as follows.


### 2.1. Mathematical Model of PMSM

The mathematical model of a PMSM is described, in two-axisd-qsynchronous rotating reference frame, as follows:
(1)diddt=-rsLdid+ωeLqLdiq+1Ldvd,diqdt=-ωeLdLqid-rsLqiq-ωeλfLq+1Lqvq,
where vd,vq are the d and q axes voltages; id, iq, are thedandqaxes currents,rsis the phase winding resistance; Ld, Lq are the d and q axes inductance;ωeis the rotating speed of magnet flux; and λf is the permanent magnet flux linkage.The current loop control of PMSM drive in Figure1 is based on a vector control approach. That is, if the id is forced to 0, the PMSM will be decoupled, and controlling a PMSM is like controlling a DC motor. After decoupling, the motor’s torque is proportional to iq:
(2)Te=3NP4λfiq=ΔKtiq.
Considering the mechanical load, the overall dynamic equation of PMSM drive system is obtained by
(3)Jmddtωr+Bmωr=Te-TL,
where Te is the motor torque,NPis pole pairs,Ktis torque constant,Jm is the inertial value,Bm is damping ratio,TL is the external torque, and ωr is rotor speed.


### 2.2. Extended Kalman Filter (EKF)

For a motor stochastic nonlinear equation, it can be written in the following form:(4)x˙(t)=f[x(t)]+Bu(t)+σ(t),(5)y(t)=h[x(t)]+μ(t),
wherex(t),u(t), and y(t) are system state, system input, and system output, respectively. Theσ(t)andμ(t)represent system noise and measurement noise which are zero-mean white Gaussian distribution with covarianceQ(t)andR(t), respectively. Once a nominal solution to a motor nonlinear equation can be found in (4)-(5), the linearized perturbation equation is
(6)δx˙(t)=F(x(t))δx(t)+Bδu(t)+σ(t),(7)δy(t)=H(x(t))δx+μ(t),
where the Jacobian and output matrices are defined as follows:
(8)F(x(t))=∂f∂x|x=x(t),H(x(t))=∂h∂x|x=x(t).
After discretization with sampling timeTc, (6) becomes
(9)x(tn)=Φ(tn,tn-1,x(tn-1))x(tn-1)+∫tn-1tnΦ(tn,tn-1,x(tn-1))Bdτ·u(tn-1)+ν(tn-1),
whereΦ(tn,tn-1,x(tn-1))is an exponential matrix and its computation can be simplified by
(10)Φ(tn,tn-1,x(tn-1))≅I+FTc.
Further,
(11)∫tn-1tnΦ(tn,tn-1,x(tn-1))Bdτ≅BTc.
Therefore, the discrete model of (6)-(7) becomes
(12)x(n)=(I+FTc)x(n-1)+BTc·u(n-1)+ν(n-1),(13)y(n)=H·x(n)+ξ(n-1),
wheretn=nTcandtn-1=(n-1)Tc.The EKF is an optimal estimator which searches the cost functionJ=∑n=1mE{x~2(n)}at the least square sense. Thex~(n) is defined byx~(n)=x^(n)-x(n), which is the difference of estimation of statex^(n)and system state x(n). The EKF algorithm is described by the following two-step recursive equations.Step 1 (prediction step).
From (4) and using a simple rectangular integration,
(14)x^n|n-1=x^n-1+(f[x^n-1]+B·un-1)Tc
or from (12)
(15)x^n|n-1=Φn-1x^n-1+BTc·un-1.
The covariance is updated by(16)Pn|n-1=Φn-1Pn-1Φn-1T+Qd.Step 2 (innovation step).
Consider(17)x^n=x^n|n-1+Kn(yn-Hx^n|n-1),(18)Pn=Pn|n-1-KnHPn|n-1.
The Kalman gain is calculated by(19)Kn=Pn|n-1HT[HPn|n-1HT+R]-1.
Figure 2 shows the procedure of EKF with two step above, after measurement update, it will be changed to time update.Figure 2
Demonstration of the EKF operation.


### 2.3. Design of the Sensorless PMSM Based on EKF

The circuit equation of PMSM on thed-q rotating coordinate in (1) is re-formulated as
(20)[vdvq]=[rs+sLs-ωeLsωeLsrs+sLs][idiq]+[0ωeλf],
whereLs=ΔLd=Lq. Transforming (20), the circuit equation of PMSM on theα-βfixed coordinate can be derived by the following equation:
(21)[vαvβ]=[rs+sLs00rs+sLs][iαiβ]+ωeλf[-sinθecosθe],
where[vαvβ]Tis voltage on fixed coordinate;[iαiβ]Tis current on fixed coordinate;θeis angular position at magnet flux; and s is differential operator. If we choose the x(t)=[iα,iβ,ωe,θe]Tas the state variable of the PMSM mathematical model, then (21) can be expanded to the following formulation:
(22)ddt[iαiβωeθe]=[-rsiαLs+ωeλfLssinθe-rsiβLs-ωeλfLscosθe0ωe]+[1Ls001Ls0000][vαvβ]=Δf+B[vαvβ],(23)[iαiβ]=[10000100][iαiβωeθe]=Δh,
wherex(t)=[iαiβωeθe]T,y(t)=⌊iαiβ⌋, andu(t)=[vαvβ]T. From (8) and (22), the Jacobian and output matrices can be obtained:
(24)F(x(t))=∂f∂x|x=x(t)=[-rsLs0λfLssinθeωeλfLscosθe0-rsLs-λfLscosθeωeλfLssinθe00000010],(25)H(x(t))=∂h∂x|x=x(t)=[10000100].
Further, substituting (24) into (10), the exponential matrix is shown as follows:
(26)Φ(tn,tn-1,x(tn-1))≅I+FTc=[1-rsTcLs0λfTcLssinθeωeλfTcLscosθe01-rsTcLs-λfTcLscosθeωeλfTcLssinθe001000Tc1]=Δ[G0ϕ13ϕ140Gϕ23ϕ24001000Tc1].
As the results, the EKF algorithm in (14)–(19) can be carried out to estimate the state value. The initial values ofQd, R, and P0need to be chosen. Through the recursive calculation, the state value ofx^(n)=[i^α(n),i^β(n),ω^e(n),θ^e(n)]T is estimated at each sampling period; then the rotor speed can be derived by
(27)ω^r(n)=ω^e(n)NP.Finally, a summary for the estimation of rotor FA and rotor speed based on EKF is presented as follows.Step 1.
Set the initial values ofQd, R, P0, and n=1.Step 2.
Measure the values ofiα(n),iβ(n),vα(n), and vβ(n) from PMSM system.Step 3.
Estimate the temporary state variables from (14). In addition, refer to (22), and the scalar form of the prediction equation can be expressed as follows:
(28)i^α(n∣n-1)=(1-TcrsLs)i^α(n-1)+ω^e(n-1)Tcλfsin(θ^e(n-1))Ls+Tcvα(n-1)Ls,i^β(n∣n-1)=(1-TsrsLs)i^β(n-1)-ω^e(n-1)Tsλfcos(θ^e(n-1))Ls+Tcvβ(n-1)Ls,ω^e(nn-1)=ω^e(n-1),θ^e(nn-1)=θ^e(n-1)+ω^e(n-1)Tc.Step 4.
Calculate theΦn-1from (26).Step 5.
Obtain the temporary covariance matrixPn∣n-1 from (16). Because matrixPn∣n-1is symmetric matrix by pij=pji, it can be chosen by the following form:
(29)Pn∣n-1=[p11p12p13p14p21p22p23p24p31p32p33p34p41p42p43p44].Hence, from (16), individual elements ofPn|n-1 could be updated as follows:
(30)p11⟸G2p11+2Gϕ14p14+2Gϕ13p13+2ϕ14ϕ13p34+ϕ142p44+ϕ132p33+q11,p12⟸G2p12+Gϕ14p24+Gϕ13p23+ϕ23(Gp13+ϕ13p33+ϕ14p34)+ϕ24(Gp14+ϕ13p34+ϕ14p44),p13⟸Gp13+ϕ13p33+ϕ14p34,p14⟸(Gp13+ϕ13p33+ϕ14p34)Tc+(Gp14+ϕ13p34+ϕ14p44),p22⟸2G(ϕ23p23+ϕ24p24)+G2p22+ϕ232p33+ϕ242p44+2ϕ23ϕ24p34+q22,p23⟸Gp23+ϕ23p33+ϕ24p34,p24⟸(Gp23+ϕ23p33+ϕ24p34)Tc+Gp24+ϕ24p44+ϕ23p34,p33⟸p33+q33,p34⟸Tcp33+p34,p44⟸p44+Tc2p33+2Tcp34+q44.Step 6.
Calculate the Kalman gain from (19). Using (25) and (29), the formulation is further simplified as follows. Firstly,
(31)[HPn∣n-1HT+R]-1=[[10000100][p11p12p13p14p21p22p23p24p31p32p33p34p41p42p43p44][10010000]ddfdd+[R1100R22][p11p12p13p14p21p22p23p24p31p32p33p34p41p42p43p44][10010000]]-1=[[p11p12p21p22]+[R1100R22]]-1,
whereRrepresents the uncertainties of measurement variable which are very difficult to know. Hence, Ris often chosen as unity matrix for simplifying the task but does not affect the outcome (31). SinceP12=P21, (31) can be rewritten in the following form:
(32)[HPn∣n-1HT+R]-1=1p11p22+p11+p22-p122+1[p22+1-p12-p12p11+1].
To normalize the input value and to prevent the numerical overflow condition which occurred during computation, (32) is normalized and those values are all divided by 4. Therefore, (32) becomes
(33)[HPn∣n-1HT+R]-1=1M[p22+14-p124-p124p11+14]=Δ[t11t12t21t22],
whereM=(p11p22+p11+p22-p122+1)/4. Then, from (19), the Kalman gain can be calculated by
(34)Kn=Δ[k11k21k31k41k12k22k32k42]=[p11p21p31p41p12p22p32p42][t11t12t21t22]=[p11t11+p12t21p21t11+p22t21p31t11+p32t21p41t11+p42t21p11t12+p12t22p21t12+p22t22p31t12+p32t22p41t12+p42t22].Step 7.
Tune the present state variables from (17), and its scalar form is shown as follows:
(35)i^α(n)=i^α(nn-1)+k11i~α(n)+k12i~β(n),i^β(n)=i^β(nn-1)+k21i~α(n)+k22i~β(n),ω^e(n)=ω^e(nn-1)+k31i~α(n)+k32i~β(n),θ^e(n)=θ^e(n∣n-1)+k41i~α(n)+k42i~β(n),
wherekijis the element in Kalman gainKn·i~α(n)=iα(n)-i^α(n∣n-1)and i~β(n)=iβ(n)-i^β(n∣n-1).Step 8.
Update the present covariance matrixPnfrom (18), whereKnHPn∣n-1can be directly obtained as follows:(36)KnHPn∣n-1=[k11p11+k12p21k11p12+k12p22k11p13+k12p23k11p14+k12p24k21p11+k22p21k21p12+k22p22k21p13+k22p23k21p14+k22p24k31p11+k32p21k31p12+k32p22k31p13+k32p23k31p14+k32p24k41p11+k42p21k41p12+k42p22k41p13+k42p23k41p14+k42p24].Step 9.
Calculate the rotor speed from (27); then setn=n+1and go back to Step 2.


## 2.1. Mathematical Model of PMSM

The mathematical model of a PMSM is described, in two-axisd-qsynchronous rotating reference frame, as follows:
(1)diddt=-rsLdid+ωeLqLdiq+1Ldvd,diqdt=-ωeLdLqid-rsLqiq-ωeλfLq+1Lqvq,
where vd,vq are the d and q axes voltages; id, iq, are thedandqaxes currents,rsis the phase winding resistance; Ld, Lq are the d and q axes inductance;ωeis the rotating speed of magnet flux; and λf is the permanent magnet flux linkage.The current loop control of PMSM drive in Figure1 is based on a vector control approach. That is, if the id is forced to 0, the PMSM will be decoupled, and controlling a PMSM is like controlling a DC motor. After decoupling, the motor’s torque is proportional to iq:
(2)Te=3NP4λfiq=ΔKtiq.
Considering the mechanical load, the overall dynamic equation of PMSM drive system is obtained by
(3)Jmddtωr+Bmωr=Te-TL,
where Te is the motor torque,NPis pole pairs,Ktis torque constant,Jm is the inertial value,Bm is damping ratio,TL is the external torque, and ωr is rotor speed.


## 2.2. Extended Kalman Filter (EKF)

For a motor stochastic nonlinear equation, it can be written in the following form:(4)x˙(t)=f[x(t)]+Bu(t)+σ(t),(5)y(t)=h[x(t)]+μ(t),
wherex(t),u(t), and y(t) are system state, system input, and system output, respectively. Theσ(t)andμ(t)represent system noise and measurement noise which are zero-mean white Gaussian distribution with covarianceQ(t)andR(t), respectively. Once a nominal solution to a motor nonlinear equation can be found in (4)-(5), the linearized perturbation equation is
(6)δx˙(t)=F(x(t))δx(t)+Bδu(t)+σ(t),(7)δy(t)=H(x(t))δx+μ(t),
where the Jacobian and output matrices are defined as follows:
(8)F(x(t))=∂f∂x|x=x(t),H(x(t))=∂h∂x|x=x(t).
After discretization with sampling timeTc, (6) becomes
(9)x(tn)=Φ(tn,tn-1,x(tn-1))x(tn-1)+∫tn-1tnΦ(tn,tn-1,x(tn-1))Bdτ·u(tn-1)+ν(tn-1),
whereΦ(tn,tn-1,x(tn-1))is an exponential matrix and its computation can be simplified by
(10)Φ(tn,tn-1,x(tn-1))≅I+FTc.
Further,
(11)∫tn-1tnΦ(tn,tn-1,x(tn-1))Bdτ≅BTc.
Therefore, the discrete model of (6)-(7) becomes
(12)x(n)=(I+FTc)x(n-1)+BTc·u(n-1)+ν(n-1),(13)y(n)=H·x(n)+ξ(n-1),
wheretn=nTcandtn-1=(n-1)Tc.The EKF is an optimal estimator which searches the cost functionJ=∑n=1mE{x~2(n)}at the least square sense. Thex~(n) is defined byx~(n)=x^(n)-x(n), which is the difference of estimation of statex^(n)and system state x(n). The EKF algorithm is described by the following two-step recursive equations.Step 1 (prediction step).
From (4) and using a simple rectangular integration,
(14)x^n|n-1=x^n-1+(f[x^n-1]+B·un-1)Tc
or from (12)
(15)x^n|n-1=Φn-1x^n-1+BTc·un-1.
The covariance is updated by(16)Pn|n-1=Φn-1Pn-1Φn-1T+Qd.Step 2 (innovation step).
Consider(17)x^n=x^n|n-1+Kn(yn-Hx^n|n-1),(18)Pn=Pn|n-1-KnHPn|n-1.
The Kalman gain is calculated by(19)Kn=Pn|n-1HT[HPn|n-1HT+R]-1.
Figure 2 shows the procedure of EKF with two step above, after measurement update, it will be changed to time update.Figure 2
Demonstration of the EKF operation.


## 2.3. Design of the Sensorless PMSM Based on EKF

The circuit equation of PMSM on thed-q rotating coordinate in (1) is re-formulated as
(20)[vdvq]=[rs+sLs-ωeLsωeLsrs+sLs][idiq]+[0ωeλf],
whereLs=ΔLd=Lq. Transforming (20), the circuit equation of PMSM on theα-βfixed coordinate can be derived by the following equation:
(21)[vαvβ]=[rs+sLs00rs+sLs][iαiβ]+ωeλf[-sinθecosθe],
where[vαvβ]Tis voltage on fixed coordinate;[iαiβ]Tis current on fixed coordinate;θeis angular position at magnet flux; and s is differential operator. If we choose the x(t)=[iα,iβ,ωe,θe]Tas the state variable of the PMSM mathematical model, then (21) can be expanded to the following formulation:
(22)ddt[iαiβωeθe]=[-rsiαLs+ωeλfLssinθe-rsiβLs-ωeλfLscosθe0ωe]+[1Ls001Ls0000][vαvβ]=Δf+B[vαvβ],(23)[iαiβ]=[10000100][iαiβωeθe]=Δh,
wherex(t)=[iαiβωeθe]T,y(t)=⌊iαiβ⌋, andu(t)=[vαvβ]T. From (8) and (22), the Jacobian and output matrices can be obtained:
(24)F(x(t))=∂f∂x|x=x(t)=[-rsLs0λfLssinθeωeλfLscosθe0-rsLs-λfLscosθeωeλfLssinθe00000010],(25)H(x(t))=∂h∂x|x=x(t)=[10000100].
Further, substituting (24) into (10), the exponential matrix is shown as follows:
(26)Φ(tn,tn-1,x(tn-1))≅I+FTc=[1-rsTcLs0λfTcLssinθeωeλfTcLscosθe01-rsTcLs-λfTcLscosθeωeλfTcLssinθe001000Tc1]=Δ[G0ϕ13ϕ140Gϕ23ϕ24001000Tc1].
As the results, the EKF algorithm in (14)–(19) can be carried out to estimate the state value. The initial values ofQd, R, and P0need to be chosen. Through the recursive calculation, the state value ofx^(n)=[i^α(n),i^β(n),ω^e(n),θ^e(n)]T is estimated at each sampling period; then the rotor speed can be derived by
(27)ω^r(n)=ω^e(n)NP.Finally, a summary for the estimation of rotor FA and rotor speed based on EKF is presented as follows.Step 1.
Set the initial values ofQd, R, P0, and n=1.Step 2.
Measure the values ofiα(n),iβ(n),vα(n), and vβ(n) from PMSM system.Step 3.
Estimate the temporary state variables from (14). In addition, refer to (22), and the scalar form of the prediction equation can be expressed as follows:
(28)i^α(n∣n-1)=(1-TcrsLs)i^α(n-1)+ω^e(n-1)Tcλfsin(θ^e(n-1))Ls+Tcvα(n-1)Ls,i^β(n∣n-1)=(1-TsrsLs)i^β(n-1)-ω^e(n-1)Tsλfcos(θ^e(n-1))Ls+Tcvβ(n-1)Ls,ω^e(nn-1)=ω^e(n-1),θ^e(nn-1)=θ^e(n-1)+ω^e(n-1)Tc.Step 4.
Calculate theΦn-1from (26).Step 5.
Obtain the temporary covariance matrixPn∣n-1 from (16). Because matrixPn∣n-1is symmetric matrix by pij=pji, it can be chosen by the following form:
(29)Pn∣n-1=[p11p12p13p14p21p22p23p24p31p32p33p34p41p42p43p44].Hence, from (16), individual elements ofPn|n-1 could be updated as follows:
(30)p11⟸G2p11+2Gϕ14p14+2Gϕ13p13+2ϕ14ϕ13p34+ϕ142p44+ϕ132p33+q11,p12⟸G2p12+Gϕ14p24+Gϕ13p23+ϕ23(Gp13+ϕ13p33+ϕ14p34)+ϕ24(Gp14+ϕ13p34+ϕ14p44),p13⟸Gp13+ϕ13p33+ϕ14p34,p14⟸(Gp13+ϕ13p33+ϕ14p34)Tc+(Gp14+ϕ13p34+ϕ14p44),p22⟸2G(ϕ23p23+ϕ24p24)+G2p22+ϕ232p33+ϕ242p44+2ϕ23ϕ24p34+q22,p23⟸Gp23+ϕ23p33+ϕ24p34,p24⟸(Gp23+ϕ23p33+ϕ24p34)Tc+Gp24+ϕ24p44+ϕ23p34,p33⟸p33+q33,p34⟸Tcp33+p34,p44⟸p44+Tc2p33+2Tcp34+q44.Step 6.
Calculate the Kalman gain from (19). Using (25) and (29), the formulation is further simplified as follows. Firstly,
(31)[HPn∣n-1HT+R]-1=[[10000100][p11p12p13p14p21p22p23p24p31p32p33p34p41p42p43p44][10010000]ddfdd+[R1100R22][p11p12p13p14p21p22p23p24p31p32p33p34p41p42p43p44][10010000]]-1=[[p11p12p21p22]+[R1100R22]]-1,
whereRrepresents the uncertainties of measurement variable which are very difficult to know. Hence, Ris often chosen as unity matrix for simplifying the task but does not affect the outcome (31). SinceP12=P21, (31) can be rewritten in the following form:
(32)[HPn∣n-1HT+R]-1=1p11p22+p11+p22-p122+1[p22+1-p12-p12p11+1].
To normalize the input value and to prevent the numerical overflow condition which occurred during computation, (32) is normalized and those values are all divided by 4. Therefore, (32) becomes
(33)[HPn∣n-1HT+R]-1=1M[p22+14-p124-p124p11+14]=Δ[t11t12t21t22],
whereM=(p11p22+p11+p22-p122+1)/4. Then, from (19), the Kalman gain can be calculated by
(34)Kn=Δ[k11k21k31k41k12k22k32k42]=[p11p21p31p41p12p22p32p42][t11t12t21t22]=[p11t11+p12t21p21t11+p22t21p31t11+p32t21p41t11+p42t21p11t12+p12t22p21t12+p22t22p31t12+p32t22p41t12+p42t22].Step 7.
Tune the present state variables from (17), and its scalar form is shown as follows:
(35)i^α(n)=i^α(nn-1)+k11i~α(n)+k12i~β(n),i^β(n)=i^β(nn-1)+k21i~α(n)+k22i~β(n),ω^e(n)=ω^e(nn-1)+k31i~α(n)+k32i~β(n),θ^e(n)=θ^e(n∣n-1)+k41i~α(n)+k42i~β(n),
wherekijis the element in Kalman gainKn·i~α(n)=iα(n)-i^α(n∣n-1)and i~β(n)=iβ(n)-i^β(n∣n-1).Step 8.
Update the present covariance matrixPnfrom (18), whereKnHPn∣n-1can be directly obtained as follows:(36)KnHPn∣n-1=[k11p11+k12p21k11p12+k12p22k11p13+k12p23k11p14+k12p24k21p11+k22p21k21p12+k22p22k21p13+k22p23k21p14+k22p24k31p11+k32p21k31p12+k32p22k31p13+k32p23k31p14+k32p24k41p11+k42p21k41p12+k42p22k41p13+k42p23k41p14+k42p24].Step 9.
Calculate the rotor speed from (27); then setn=n+1and go back to Step 2.


## 3. Design of Sensorless Speed Controller IP

Figure3 illustrates the internal architecture of the proposed FPGA implementation of a sensorless speed control IP for PMSM drive system. The internal circuit comprises a Nios II embedded processor IP and an application IP. The Nios II processor is depicted to generate the speed command, collect the response data, and communicate with external device. All programs in Nios II processor are coded in the C language. The application IP includes mainly a circuit of an EKF-based rotor FA and rotor speed estimation, a speed PI controller, a circuit for current controllers and coordinate transformation (CCCT), a SVPWM circuit, and an ADC interface circuit. The sampling frequencies of the speed control loop and current control loop are designed with 2 kHz and 16 kHz, respectively. The frequency divider generates 50 Mhz (Clk) and 12.5 Mhz (Clk-step) clock to supply all circuits in Figure 3. The internal circuit of CCCT performs the function of two PI controllers, table look-up for sin/cos function, and the coordinate transformation for Clark, Park, inverse Park, and modified inverse Clarke. The detailed circuit design of CCCT, speed PI controller, ADC interface, and SVPWM generation in Figure 3 refers to [8]. Herein, only the digital hardware circuit design of the EKF-based rotor FA and rotor speed estimation is described. To reduce the FPGA resource usage, a finite-state machine (FSM) is employed to model the EKF-based algorithm whose design procedure and mathematical formulation are presented in the previous section. The design result is shown in Figure 4 which uses one adder, one multiplier, one divider, two sin/cos look-up table, shifters, registers, and so forth and manipulates 131 steps machine to carry out the overall computation. The data type adopts 16-bit length with Q15 format and 2’s complement operation. The multiplier, adder, and divider apply Altera LPM (Library Parameterized Modules) standard. In Figure 4, the steps s0~s12 execute the calculation of Jacobian matrix and predict temporary state variables; s13~s68 perform the computation of temporary covariance matrix; steps s69~s70 are for state error calculation; steps s71~s120 update of the present covariance matrix Pn and calculate the Kalman gain Kn; s121~s130 describe the present state tuning as well as executing the computation of rotor FA and rotor speed. According to (30), the computation of temporary covariance matrixPn∣n-1in Figure 4 is described in Figure 5. Further, according to (18), (33)-(34), and (36), the update of the present covariance matrix Pn and calculation the Kalman gain Kn in Figure 4 are illustrated in Figure 6. The operation of each step in Figure 4 is 80 ns (12.5 MHz) in FPGA; therefore a total of 131 steps need only 10.48 μs operation times. The FPGA (Altera) resource usages of EKF-based rotor FA and rotor speed estimation in Figure 4 is 4,158 LEs (logic elements) and 2,890 RAM bits. Additionally, the FPGA resource usages of the circuits of CCCT, speed PI controller, ADC interface, and SVPWM in Figure 3 are 864 LEs and 24,576 RAM, 2,043 LEs and 0RAM bits, 136 LEs and 0RAM bits, and 1,221 LEs and 0RAM bits, respectively.Figure 3
Internal architecture of a sensorless speed control IP in FPGA.Figure 4
State diagram of FSM for describing the EKF-based algorithm.Figure 5
State diagram of FSM for prediction of covariance matrixPn∣n-1.State diagram of FSM for (a) the calculation of matrixKn and (b) the update of covariance matrixPn.
(a)
(b)


## 4. Simulink/ModelSim Cosimulation and Its Simulation Results

The Simulink/ModelSim cosimulation architecture for sensorless PMSM speed control system is shown in Figure7. The SimPowerSystem-Blockset in the Simulink executes the PMSM and the IGBT-based inverter. The EDA simulator link for ModelSim executes the cosimulation using VHDL code running in ModelSim program. It executes the function of the sensorless speed controller by three works. The work-1 to work-3 of ModelSim in Figure 7, respectively, perform the function of speed loop PI controller, the function of current controller and coordinate transformation (CCCT) and SVPWM, and the function of EKF-based rotor FA and rotor speed estimation. The sampling frequency in current control loop and EKF-based estimation algorithm is designed with 16 kHz, but in speed control loop is 2 kHz. The clocks of 50 MHz and 12.5 MHz will supply all works of ModelSim. The FPGA (Altera) resource usages of work-1 to work-3 are 2,043 LEs and 0 RAM bits, 2,085 LEs and 24,576 RAM bits, and 4,158 LEs and 2,890 RAM bits, respectively.Figure 7
The Simulink/ModelSim cosimulation architecture for sensorless speed control of PMSM drive.In simulation, the designed PMSM parameters used in Figure7 are as follows: pole pairs are 4, stator phase resistance is 1.3 Ω, stator inductance is 6.3 mH, inertia is J=0.000108 kg*m2, and friction factor is F=0.0013 N*m*s. PI gains in speed control loop are designed with Kp=0.2442 and Ki=0.00458. To evaluate the proposed controller performance, PMSM running at low speed (200 rpm) and inverse speed (from 200 rpm to −200 rpm or vice versa) condition is firstly considered, and its simulation results regarding as the speed step response, the current response as well as the actual rotor flux angle (FA) and the estimated rotor flux angle are shown in Figure 8. It shows that the motor speed gives a good dynamic response performance with a little overshoot and 25 ms rising time. It also shows that the estimated rotor speed can track actual rotor speed well. Figure 8(b) shows a successful vector control becauseid is controlled to zero. Figure 8(c) shows that the estimated rotor FA can fast track the actual rotor FA even in the inverse speed command condition. Further, another simulation case while the PMSM running at wide-range high speed (from 0 rpm → 1000 rpm → 2000 rpm → 1000 rpm → −1000 pm → −2000 rpm) and inverse speed (from 1000 rpm to −1000 rpm) condition is tested, and its simulation results are shown in Figure 9. It presents that not only at speed tracking but also at rotor FA tracking, the estimated value adopted by EKF method can give a good follow to the real rotor value. However, except in the initial condition, the actual rotor speed gives a fast step response by 18 ms rising time, near 0 mm steady-state value, and maximum 17.5% overshoot.PMSM running at low speed (200 rpm) and inverse speed (±200 rpm) condition. (a) Speed step response, (b) current response, and (c) actual rotor FA and estimated rotor FA.
(a)
(b)
(c)PMSM running at wide-range high speed (2000 rpm) and inverse speed (±1000 rpm) condition. (a) Speed step response, (b) current response, AND (c) actual rotor FA and estimated rotor FA.
(a)
(b)
(c)


## 5. Experimental System and Results

After confirming the correctness of the proposed EKF-based sensorless control IP by simulation, the VHDL codes are directly applied to the experimental FPGA-based sensorless PMSM drive system. The block diagram and real experimental system are depicted in Figure10. The main devices include a PMSM, a DE2 board with Altera CycloneII FPGA, a motor driver, and a power supplier. The parameters of the PMSM are rs=0.63 Ω, L=2.77 mH, and 4 pole pairs. The input voltage, continuous current, rating torque, rating speed, and continuous power of the PMSM are 220 V, 12 A, 2.3 N-m, 3000 rpm, and 750 W, respectively. The Altera CycloneII EP2C35 chip adopted in the design possesses 33,216 LEs, maximum 475 available I/O pins, 483,840 RAM, and 35 embedded multipliers. The chip can be embedded with a Nios II multicore processor that is equipped with several 32-bit CPU, flexibility of core size, 1 to 16 Mbytes of flash memory in the available memory chip, 1 Mbyte SRAM, 16 M byte SDRAM, and 4 Gbytes memory outside of the chip. In implementation, except the VHDL codes of CCCT, EKF, and PI controller, the VHDL code of ADC interface circuit should be added and integrated and then downloaded into FPGA. The PI gains of speed control loop in experiment are designed with Kp=0.512 and Ki=0.021.Experimental system (a) with block diagram and (b) in real system.
(a)
(b)In the experimental case, PMSM running at low speed (200 rpm) and inverse speed (from 200 rpm to −200 rpm or vice versa) condition are firstly considered, and its experiment results regarding the speed step response as well as the actual rotor FA and the estimated rotor FA are shown in Figure11. Figure 11(a) shows that the motor speed gives a good dynamic response performance with a little oscillation and overshoot and it has 15 ms rising time. Figure 11(b) shows that the estimated FA can fast track the actual rotor FA even in the inverse speed command condition. Further, another experimental case while the PMSM running at medium speed (from 0 rpm → 300 rpm → 600 rpm → 900 rpm → 600 rpm → 900 rpm) condition is evaluated, and its experiment results are shown in Figure 12. The Figure 12(a) presents a fast speed tracking with 15 ms rising time and no occurrence of overshoot condition. Figure 12(b) shows a successful vector control due to theid being controlled to zero. Therefore, from the simulation results in Figures 8 and 9 and the experimental results in Figures 11 and 12, it is demonstrated that the EKF-based estimation algorithm used in sensorless PMSM drive and the behavior description by using VHDL are effective and correct.(a) Step speed response of actual rotor speed and estimated rotor speed and (b) actual rotor FA and estimated rotor FA.
(a)
(b)(a) Step speed response of actual rotor speed and estimated rotor speed and (b) current response.
(a)
(b)


## 6. Conclusions

This study has presented a sensorless PMSM drive based on EKF and successfully demonstrated its performance through cosimulation by using Simulink/ModelSim and implementation by using FPGA. In realization aspect, the VHDL is used to describe the behavior of EKF algorithm, and FSM is used to reduce the FPGA resource usage; therefore, it only needs 4,158 LEs and 2,890 RAM bits. In computational power aspect, the operation time to complete the computation of EKF algorithm is only 10.48μs, which is less than the 62.5 μs (16 KHz) sampling time in current control loop. In simulated and experimental results, it shows that the use of EKF in sensorless PMSM drive can accurately estimate the rotor FA and rotor speed, and it can give a good step response performance in case of low speed control, inverse speed control, and high speed control as well.



---

*Source: 919318-2013-12-04.xml*