CN104993765A - Method for estimating rotating speed of brushless direct current motor - Google Patents

Abstract

The invention discloses a method for estimating the rotating speed of a brushless direct current motor. In allusion to a problem that a Kalman filter system noise matrix and a measurement noise matrix are difficult to be acquired in estimation of a Kalman filter for the rotating speed of the brushless direct current motor, the invention provides a method in which the system noise matrix and the measurement noise matrix are optimized at first by using an ant colony algorithm and a particle swarm algorithm, and then estimation for the rotating speed is carried out by using an extended Kalman filter. The algorithm integrates advantages of the ant colony algorithm and the particle swarm algorithm, and the particle swarm algorithm is introduced into the ant colony algorithm in allusion to defects that a prematurity phenomenon easily occurs in the ant colony algorithm and that the particle swarm algorithm is poor in local searching ability in a later period, thereby enabling ant colony algorithm to have characteristics of the particle swarm algorithm. The optimized Kalman filter is applied to estimation of the rotating speed of the brushless direct current motor, thereby being capable of improving the precision in estimation for the rotating speed of the brushless direct current motor.

Description

A kind of method for estimating rotating speed of brshless DC motor

Technical field

The invention belongs to motor control method field, relate to a kind of brshless DC motor method for estimating rotating speed.

Background technology

Motor speed detects mainly through physical sensors and software algorithm two kinds of methods obtain.Utilize physical sensors detection can increase fault point and physics cost, and carry out Speed Identification by software algorithm, the deficiencies such as fault is many, volume is large, difficult in maintenance can be overcome in hardware approach.Therefore the research of the virtual speed probe utilizing software algorithm to test the speed is become to the important research direction of modern transmission control technology, Kalman filter is considered to the important method that virtual speed probe medium speed is estimated.

Using extended Kalman filter to carry out one of subject matter of speed estimate existence is exactly choosing of noise matrix, the precision of speed estimate and the degree of convergence depend on the priori knowledge of user to noise to a great extent, so can improve the precision of speed estimate for the accurate planning of system noise matrix and measurement noises matrix.Utilize particle cluster algorithm to be optimized noise matrix, although can improve the precision of speed estimate result, particle cluster algorithm is poor in the local search ability in later stage, and feedback information utilizes insufficient.So, ant group algorithm is combined with particle cluster algorithm, a kind of integration algorithm is proposed, utilize quick, the global convergence of particle cluster algorithm and the pheromones positive feedback mechanism of ant group algorithm, reach mutual supplement with each other's advantages, noise matrix is optimized, thus improves precision and the overall performance of extended Kalman filter, make the estimated value of the extended Kalman filter after optimization to brshless DC motor rotating speed more accurate.

Summary of the invention

Goal of the invention: compared to physical sensors, utilizes virtual speed probe to test the speed can effectively reduce costs motor, improves reliability.In addition, carry out in the process of speed estimate at employing extended Kalman filter, initial priori value can not represent noise during real work, after utilizing particle cluster algorithm and ant group algorithm to be optimized system noise matrix and calculation matrix, effectively can improve the precision of speed estimate value and the overall performance of Kalman filter.

Technical scheme: when carrying out brshless DC motor speed estimate with extended Kalman filter, first utilize particle cluster algorithm and ant group algorithm to be optimized system noise matrix and measurement noises matrix, the correction value obtained carries out speed estimate as the input parameter of extended Kalman filter.The method for estimating rotating speed of brshless DC motor of the present invention mainly comprises the following steps:

Step 1: be based upon the square wave brshless DC motor model under static abc coordinate system, winding is wye connection.By threephase stator phase current i _a, i _b, i _cwith rotor speed ω and rotor position as state variable, the Mathematical Modeling (1) of structure motor:

$\left\{\begin{array}{c}\frac{{\mathrm{di}}_a}{dt}=\frac{1}{L_a}\[\frac{u_{ab}-u_{ac}}{3}-R_a{i}_a-e_a\]\\ \frac{{\mathrm{di}}_b}{dt}=\frac{1}{L_b}\[\frac{u_{bc}-u_{ab}}{3}-R_b{i}_b-e_b\]\\ \frac{{\mathrm{di}}_c}{dt}=\frac{1}{L_c}\[\frac{u_{ca}-u_{bc}}{3}-R_c{i}_c-e_c\]\\ \frac{d\mathrm{\ω}}{dt}=\frac{{\mathrm{Pk}}_e}{J}(T_e-T_L)\\ \frac{d\mathrm{\θ}}{dt}=\mathrm{\ω}\end{array}\right.---\left(1\right)$

U in formula _ab, u _bc, u _cabe respectively the threephase stator line voltage of brshless DC motor; R _a, R _b, R _cfor threephase stator resistance, and R _a=R _b=R _c=R; L _a, L _b, L _cfor threephase stator equivalent inductance, and L _i=L _is-L _im, i=a, b, c, wherein L _isand L _imbe respectively the self-induction and mutual inductance that represent i phase stator winding, L simultaneously _is=L _s, L _im=L _mif, L=L _s-L _m, then threephase stator equivalent inductance has L _a=L _b=L _c=L; P is number of pole-pairs; J is rotor inertia matrix, k _cfor back EMF coefficient; e _a, e _b, e _cbeing each phase back-emf, is the trapezoidal wave with 120 ° of electrical degree platforms, and differs 120 ° of electrical degrees successively.Because electromagnetic torque T _e=(e _ai _a+ e _bi _b+ e _ci _c)/ω, so formula (1) can be written as:

$\stackrel{\·}{x}=f\left(x\right(t),u(t),t)---\left(2\right)$

Wherein, x (t)=[i _ai _bi _cω θ] ^t; Formula (2) discretization is obtained:

x _k+1＝F _k(x _k)x _k+G _ku _k(3)

So motor Random Discrete model is:

$\left\{\begin{array}{c}{x}_{k+1}=F_k\left(x_k\right)x_k+G_k{u}_k+w_k\\ y_k={\mathrm{Hx}}_k+v_k\end{array}\right.---\left(4\right)$

Wherein: $x_k={[i_a,i_b,i_c,\mathrm{\ω},\mathrm{\θ}]}_k^T,u_k={[u_{\mathrm{ab}},u_{\mathrm{bc}},u_{\mathrm{ca}}]}_k^T,y_k={[i_a,i_b,i_c]}_k^T,$

$F_k\left(x_k\right)={\left[\begin{array}{ccccc}1-\frac{RT}{L}& 0& 0& \frac{k_{e}T(a_{11}+a_{12\mathrm{\θ}})}{L}& 0\\ 0& 1-\frac{RT}{L}& 0& \frac{k_{e}T(a_{21}+a_{22\mathrm{\θ}})}{L}& 0\\ 0& 0& 1-\frac{RT}{L}& \frac{k_{e}T(a_{31}+a_{32\mathrm{\θ}})}{L}& 0\\ 0& 0& 0& 1& 0\\ 0& 0& 0& T& 1\end{array}\right]}_k,$

$G_k\left(x_k\right)={\left[\begin{array}{ccccc}\frac{T}{3L}& -\frac{T}{3L}& 0& 0& 0\\ 0& \frac{T}{3L}& -\frac{T}{3L}& 0& 0\\ -\frac{T}{3L}& 0& \frac{T}{3L}& 0& 0\end{array}\right]}_k^T,H=\left[\begin{array}{ccccc}1& 0& 0& 0& 0\\ 0& 1& 0& 0& 0\\ 0& 0& 1& 0& 0\end{array}\right]$

Wherein, T is the sampling time, w _kand v _kbe respectively system noise vector sum measurement noises vector, be zero mean Gaussian white noise, have covariance matrix to be respectively Q and R.

Step 2: in particle cluster algorithm process, the corresponding particle of often group parameter of noise matrix Q and R of Kalman filter.Need to carry out initialization to noise matrix Q and R in initial, and need position vector, the velocity vector of random initializtion particle at first at particle cluster algorithm.Owing to having five state variables and three output variables in motor model, so get respectively:

Q=diag (Q ₁₁q ₂₂q ₃₃q ₄₄q ₅₅), R=diag (R ₁₁r ₂₂r ₃₃), here conveniently by calculating, make f (Q, R)=[Q ₁₁q ₂₂q ₃₃q ₄₄q ₅₅r ₁₁r ₂₂r ₃₃] ^t

Then judge whether particle is optimum fitness function and is:

$fitness(Q,R)=\frac{1}{k}\underset{i=0}{\overset{k}{\Σ}}{({\mathrm{\ω}}_{ri}^{*}-{\mathrm{\ω}}_{ri})}^2---\left(5\right)$

Wherein: and ω _ribe respectively the actual speed of moment i and estimate rotating speed, k equals the ratio in simulation time and sampling period;

Step 3: upgrade particle rapidity and position according to formula (6) and (7):

v _id(t+1)＝wv _id(t)+c ₁r ₁(p _id(t)-x _id(t))+c ₂r ₂(p _gd(t)-x _id(t)) (6)

x _id(t+1)＝x _id(t)+v _id(t+1) (7)

Wherein, p _idfor the optimal solution that each particle searches oneself; p _gdfor the optimal solution that whole population searches; v _idrepresent the speed of i-th particle in d dimension; x _idrepresent i-th position of particle in d dimension; W is inertia weight; r ₁and r ₂for being evenly distributed on the random number between [0,1]; c ₁and c ₂for acceleration limits the factor.The process of renewal speed and position is: the history adaptive optimal control angle value of more each particle and global optimum's fitness value, if the fitness value of certain particle current location is better than history value, then the history optimum position of this particle and fitness function value are replaced; If the history optimum position of certain particle is better than global optimum's adaptive value, Ze Zhi global optimum fitness value is the history adaptive optimal control value of this particle, records this global optimum's particle position.Until it is enough little to reach error, now, some optimization solutions are obtained.

Step 4: the some optimization solutions utilizing particle cluster algorithm to obtain generate the pheromones initial distribution of ant group algorithm.Each particle history optimal location P (i) that the corresponding particle cluster algorithm in position X (i) of ant i solves, that is:

X(i)＝P(i) (8)

Ant pheromones initial distribution utilizes the evaluation function value of each ant position, and formula is:

$T\left(i\right)={\mathrm{ka}}^{-f\left(X_i\right)}---\left(9\right)$

Wherein, k be greater than 0 constant, 0 < a < 1, f (x _i) be the target function value of ant group algorithm.Transition probability P (k, j) is:

$P(k,j)=\frac{T^{\mathrm{\&alpha;}}\left(j\right)Z^{\mathrm{\&beta;}}\left(j\right)}{\underset{i}{\&Sigma;}{T}^{\mathrm{\&alpha;}}\left(i\right)Z^{\mathrm{\&beta;}}\left(i\right)}---\left(10\right)$

Wherein, Δ F _ij=F _j-F _i, F _ifor the functional value of ant i position, T (j) is the pheromones quantity of ant j present position, and a < 1, α and β is non-negative parameter; In ant group algorithm, the transition rule of ant is: move by transition probability and random search in the field of radius r, is exactly carry out global search, if the pheromones increment of ant j adds S pheromones, then by probability transfer

ΔT _j＝ΔT _j+S (11)

In field, search is exactly carry out Local Search, and random search in the field of radius r, if the target function value of reposition is larger than initial value, then getting reposition is ant current location, otherwise do not change, radius r often terminates once to circulate and just reduces, the rear lastest imformation that once circulated element:

T(k)＝dT(k)+ΔT(k) (12)

Wherein: T (k) represents track persistence, 0≤d < 1 for pheromones track intensity, d.

Step 5: within each sampling period, extended Kalman filter carries out iterative process, and the correction value Q obtained and R carries out speed estimate as the input parameter of extended Kalman filter.Status predication formula is:

${\hat{x}}_{k+1|k}=F_k\left({\hat{x}}_{k|k}\right)\&CenterDot;{\hat{x}}_{k|k}{\mathrm{+G}}_k\left({\hat{x}}_{k|k}\right)\&CenterDot;u_k---\left(13\right)$

Error covariance matrix estimation formulas is:

$P_{k+1|k}=F_k^{\&prime;}\left({\hat{x}}_{k|k}\right)\&CenterDot;P_{k|k}\&CenterDot;{\left(F_k^{\&prime;}\right({\hat{x}}_{k|k}\left)\right)}^T+Q_k---\left(14\right)$

Wherein, $F_k^{\&prime;}=\frac{\&part;\&lsqb;F_k\left({\hat{x}}_{k|k}\right)\&CenterDot;{\hat{x}}_{k|k}\&rsqb;+G_k\left({\hat{x}}_{k|k}\right)\&CenterDot;u_k}{\&part;x_k}{|}_{x_k=x_{k|k}}---\left(15\right)$

Kalman gain formula is:

$K_k=P_{k|k-1}{H}^T{\&lsqb;H\&CenterDot;P_{k|k-1}\&CenterDot;H^T+R_k\&rsqb;}^{-1}---\left(16\right)$

Renewal state equation is:

${\hat{x}}_{k|k}={\hat{x}}_{k|k-1}+K_k\&lsqb;y_k-H\&CenterDot;{\hat{x}}_{k|k-1}\&rsqb;---\left(17\right)$

Upgrade error co-variance matrix:

P _k|k＝[1-K _k·H]P _k|k-1(18)

Just can carry out Kalman's interative computation according to formula (13) to formula (18) and recursion initial value, thus estimate motor speed.

Beneficial effect: the present invention is described in detail the technical problem involved by the method for estimating rotating speed of brshless DC motor, when carrying out speed estimate for employing extended Kalman filter, be difficult to the problem obtaining the noise matrix affecting extended Kalman filter speed estimate performance, add the algorithm that noise matrix is optimized, combine by ant group algorithm and particle cluster algorithm, be incorporated in extended Kalman filter.Quick, the global convergence of this algorithm synthesis particle cluster algorithm and the pheromones positive feedback mechanism of ant group algorithm, reach the object of mutual supplement with each other's advantages, effectively can improve the precision of extended Kalman filter, thus the precision improved based on the brshless DC motor method for estimating rotating speed of extended Kalman filter and overall performance.

Accompanying drawing explanation

Fig. 1 is the structure chart of the extended Kalman filter based on ant group particle group optimizing;

Fig. 2 is the integration algorithm flow chart of ant group algorithm and particle cluster algorithm.

Embodiment

Below with reference to accompanying drawing and embodiment, the present invention is described in further detail, and this enforcement does not form restriction to the present invention.

Brshless DC motor method for estimating rotating speed mainly comprises and being made up of extended Kalman filter 1, target function 2 and ant group particle cluster algorithm 3.Its structure is as shown in Figure 1: first utilize extended Kalman filter 1 to carry out speed estimate, obtain posteriority speed estimate value, according to target function 2 and filtered parameter value, ant group's particle cluster algorithm 3 couples of Q and R are used to carry out optimizing, the correction value obtained is carried out parameter Estimation, until the speed estimate value be optimized as the input parameter of extended Kalman filter 1.

Be based upon the square wave brshless DC motor model under static abc coordinate system, winding is wye connection.By threephase stator phase current i _a, i _b, i _cwith rotor speed ω and rotor position as state variable, the Mathematical Modeling (1) of structure motor:

$\left\{\begin{array}{c}\frac{{\mathrm{di}}_a}{dt}=\frac{1}{L_a}\&lsqb;\frac{u_{ab}-u_{ac}}{3}-R_a{i}_a-e_a\&rsqb;\\ \frac{{\mathrm{di}}_b}{dt}=\frac{1}{L_b}\&lsqb;\frac{u_{bc}-u_{ab}}{3}-R_b{i}_b-e_b\&rsqb;\\ \frac{{\mathrm{di}}_c}{dt}=\frac{1}{L_c}\&lsqb;\frac{u_{ca}-u_{bc}}{3}-R_c{i}_c-e_c\&rsqb;\\ \frac{d\mathrm{\&omega;}}{dt}=\frac{{\mathrm{Pk}}_e}{J}(T_e-T_L)\\ \frac{d\mathrm{\&theta;}}{dt}=\mathrm{\&omega;}\end{array}\right.---\left(1\right)$

U in formula _ab, u _bc, u _cabe respectively the threephase stator line voltage of brshless DC motor; R _a, R _b, R _cfor threephase stator resistance, and R _a=R _b=R _c=R; L _a, L _b, L _cfor threephase stator equivalent inductance, and L _i=L _is-L _im, i=a, b, c, wherein L _isand L _imbe respectively the self-induction and mutual inductance that represent i phase stator winding, L simultaneously _is=L _s, L _im=L _mif, L=L _s-L _m, then threephase stator equivalent inductance has L _a=L _b=L _c=L; P is number of pole-pairs; J is rotor inertia matrix, k _cfor back EMF coefficient; e _a, e _b, e _cbeing each phase back-emf, is the trapezoidal wave with 120 ° of electrical degree platforms, and differs 120 ° of electrical degrees successively, then define e _a, e _b, e _cfor:

$\left\{\begin{array}{c}{e}_a=k_e\mathrm{\&omega;}(a_{11}+a_{21}(\mathrm{\&theta;}-{\mathrm{\&theta;}}_{\mathrm{int}}))\\ e_b=k_e\mathrm{\&omega;}(a_{21}+a_{22}(\mathrm{\&theta;}-{\mathrm{\&theta;}}_{\mathrm{int}}))\\ e_c=k_e\mathrm{\&omega;}(a_{31}+a_{32}(\mathrm{\&theta;}-{\mathrm{\&theta;}}_{\mathrm{int}}))\end{array}\right.---\left(2\right)$

A in formula ₁₁to a ₃₂definition as shown in table 1, θ _int=2 π fix (θ/2 π), fix is bracket function.

Table 1 a ₁₁to a ₃₂definition value

Because electromagnetic torque T _e=(e _ai _a+ e _bi _b+ e _ci _c)/ω, so formula (1) can be written as:

$\stackrel{\&CenterDot;}{x}=f\left(x\right(t),u(t),t)---\left(3\right)$

Wherein, x (t)=[i _ai _bi _cω θ] ^t; Formula (3) discretization is obtained:

x _k+1＝F _k(x _k)x _k+G _ku _k(4)

So motor Random Discrete model is:

$\left\{\begin{array}{c}{x}_{k+1}=F_k\left(x_k\right)x_k+G_k{u}_k+w_k\\ y_k={\mathrm{Hx}}_k+v_k\end{array}\right.---\left(5\right)$

Wherein: $x_k={[i_a,i_b,i_c,\mathrm{\&omega;},\mathrm{\&theta;}]}_k^T,u_k={[u_{\mathrm{ab}},u_{\mathrm{bc}},u_{\mathrm{ca}}]}_k^T,y_k={[i_a,i_b,i_c]}_k^T,$

$F_k\left(x_k\right)={\left[\begin{array}{ccccc}1-\frac{RT}{L}& 0& 0& \frac{k_{e}T(a_{11}+a_{12\mathrm{\&theta;}})}{L}& 0\\ 0& 1-\frac{RT}{L}& 0& \frac{k_{e}T(a_{21}+a_{22\mathrm{\&theta;}})}{L}& 0\\ 0& 0& 1-\frac{RT}{L}& \frac{k_{e}T(a_{31}+a_{32\mathrm{\&theta;}})}{L}& 0\\ 0& 0& 0& 1& 0\\ 0& 0& 0& T& 1\end{array}\right]}_k,$

$G_k\left(x_k\right)={\left[\begin{array}{ccccc}\frac{T}{3L}& -\frac{T}{3L}& 0& 0& 0\\ 0& \frac{T}{3L}& -\frac{T}{3L}& 0& 0\\ -\frac{T}{3L}& 0& \frac{T}{3L}& 0& 0\end{array}\right]}_k^T,H=\left[\begin{array}{ccccc}1& 0& 0& 0& 0\\ 0& 1& 0& 0& 0\\ 0& 0& 1& 0& 0\end{array}\right]\&CenterDot;$

Wherein, T is the sampling time, w _kand v _kbe respectively system noise vector sum measurement noises vector, be zero mean Gaussian white noise, have covariance matrix to be respectively Q and R.

The corresponding particle of often group parameter of Q and R, needs to carry out initialization to noise matrix Q and R, and needs position vector, the velocity vector of random initializtion particle at first at particle cluster algorithm in initial.Owing to having five state variables and three output variables in motor model, so get respectively:

Q=diag (Q ₁₁q ₂₂q ₃₃q ₄₄q ₅₅), R=diag (R ₁₁r ₂₂r ₃₃), here conveniently by calculating, make f (Q, R)=[Q ₁₁q ₂₂q ₃₃q ₄₄q ₅₅r ₁₁r ₂₂r ₃₃] ^t

Then judge whether particle is optimum fitness function and is:

$fitness(Q,R)=\frac{1}{k}\underset{i=0}{\overset{k}{\&Sigma;}}{({\mathrm{\&omega;}}_{ri}^{*}-{\mathrm{\&omega;}}_{ri})}^2---\left(6\right)$

Wherein: and ω _ribe respectively the actual speed of moment i and estimate rotating speed, k equals the ratio in simulation time and sampling period;

Particle rapidity and position is upgraded according to formula (1.10) and (1.11):

v _id(t+1)＝wv _id(t)+c ₁r ₁(p _id(t)-x _id(t))+c ₂r ₂(p _gd(t)-x _id(t)) (7)

x _id(t+1)＝x _id(t)+v _id(t+1) (8)

Wherein, p _idfor the optimal solution that each particle searches oneself; p _gdfor the optimal solution that whole population searches; v _idrepresent the speed of i-th particle in d dimension; x _idrepresent i-th position of particle in d dimension; W is inertia weight; r ₁and r ₂for being evenly distributed on the random number between [0,1]; c ₁and c ₂for acceleration limits the factor.The process of renewal speed and position is: the history adaptive optimal control angle value of more each particle and global optimum's fitness value.If the fitness value of certain particle current location is better than history value, then the history optimum position of this particle and fitness function value are replaced; If the history optimum position of certain particle is better than global optimum's adaptive value, Ze Zhi global optimum fitness value is this history adaptive optimal control value, records this global optimum's particle position.Until it is enough little to reach error, now, some optimization solutions are obtained.

The some optimization solutions utilizing particle cluster algorithm to obtain generate the pheromones initial distribution of ant group algorithm.Each particle history optimal location P (i) that the corresponding particle cluster algorithm in position X (i) of ant i solves, that is:

X(i)＝P(i) (9)

Ant pheromones initial distribution utilizes the evaluation function value of each ant position, and formula is:

$T\left(i\right)={\mathrm{ka}}^{-f\left(X_i\right)}---\left(10\right)$

Wherein, k be greater than 0 constant, 0 < a < 1, f (x _i) be the target function value of ant group algorithm.Transition probability P (k, j) is:

$P(k,j)=\frac{T^{\mathrm{\&alpha;}}\left(j\right)Z^{\mathrm{\&beta;}}\left(j\right)}{\underset{i}{\&Sigma;}{T}^{\mathrm{\&alpha;}}\left(i\right)Z^{\mathrm{\&beta;}}\left(i\right)}---\left(11\right)$

Wherein, Δ F _ij=F _j-F _i, F _ifor the functional value of ant i position, T (j) is the pheromones quantity of ant j present position, and a < 1, α and β is non-negative parameter; In ant group algorithm, the transition rule of ant is: move by transition probability and random search in the field of radius r, is exactly carry out global search, if the pheromones increment of ant j adds the pheromones of S unit, then by probability transfer

ΔT _j＝ΔT _j+S (12)

In field, search is exactly carry out Local Search, and random search in the field of radius r, if the target function value of reposition is larger than initial value, then getting reposition is ant current location, otherwise do not change, radius r often terminates once to circulate and just reduces, the rear lastest imformation that once circulated element:

T(k)＝dT(k)+ΔT(k) (13)

Wherein: T (k) represents track persistence, 0≤d < 1 for pheromones track intensity, d.

Within each sampling period, extended Kalman filter carries out iterative process, and the correction value Q obtained and R carries out speed estimate as the input parameter of extended Kalman filter.Status predication formula is:

${\hat{x}}_{k+1|k}=F_k\left({\hat{x}}_{k|k}\right)\&CenterDot;{\hat{x}}_{k|k}{\mathrm{+G}}_k\left({\hat{x}}_{k|k}\right)\&CenterDot;u_k---\left(14\right)$

Error covariance matrix estimation formulas is:

$P_{k+1|k}=F_k^{\&prime;}\left({\hat{x}}_{k|k}\right)\&CenterDot;P_{k|k}\&CenterDot;{\left(F_k^{\&prime;}\right({\hat{x}}_{k|k}\left)\right)}^T+Q_k---\left(15\right)$

Wherein, $F_k^{\&prime;}=\frac{\&part;\&lsqb;F_k\left({\hat{x}}_{k|k}\right)\&CenterDot;{\hat{x}}_{k|k}\&rsqb;+G_k\left({\hat{x}}_{k|k}\right)\&CenterDot;u_k}{\&part;x_k}{|}_{x_k=x_{k|k}}---\left(16\right)$

Kalman gain formula is:

$K_k=P_{k|k-1}{H}^T{\&lsqb;H\&CenterDot;P_{k|k-1}\&CenterDot;H^T+R_k\&rsqb;}^{-1}---\left(17\right)$

Renewal state equation is:

${\hat{x}}_{k|k}={\hat{x}}_{k|k-1}+K_k\&lsqb;y_k-H\&CenterDot;{\hat{x}}_{k|k-1}\&rsqb;---\left(18\right)$

Upgrade error co-variance matrix:

P _k|k＝[1-K _k·H]P _k|k-1(19)

Just can carry out Kalman's interative computation according to formula (14) to formula (19) and recursion initial value, thus estimate the rotating speed of brshless DC motor.

Claims (2)

1. a method for estimating rotating speed for brshless DC motor, is characterized in that, specifically comprises the following steps:

Step 1: set up the square wave brshless DC motor model under static abc coordinate system, by threephase stator phase current i _a, i _b, i _cwith rotor speed ω and rotor position as state variable, the Mathematical Modeling of structure motor, obtains motor discrete model by equation discretization;

Step 2: init state error co-variance matrix P (0), state and the corresponding particle of the often group parameter of noise matrix Q and R, Q and R, the position vector of random initializtion particle, velocity vector, calculate fitness function value;

Step 3: utilize particle cluster algorithm, more the position of new particle and speed, the individual pole figure of merit and position, the overall pole figure of merit and position thereof, carry out the recursive iteration of particle algorithm, generate some groups of optimization solutions;

Step 4: according to optimization solution information generated element initial distribution, utilize ant group algorithm, calculating moves to next node probability and lastest imformation is plain, carries out ant group algorithm recursive iteration, obtain the optimal value of Q and R;

Step 5: according to expanded Kalman filtration algorithm, using noise matrix Q and R after optimization as the input parameter of extended Kalman filter, through Kalman filtering iteration, estimates motor speed value.

2., in the brshless DC motor method for estimating rotating speed based on ant group particle group optimizing extended Kalman filter according to claim 1, it is characterized in that:

In step 1, the discrete model of brshless DC motor:

$\left\{\begin{array}{c}{x}_{k+1}=F_k\left(x_k\right)x_k+G_k{u}_k+w_k\\ y_k={\mathrm{Hx}}_k+v_k\end{array}\right.---\left(1\right)$

Wherein: $x_k={\&lsqb;i_a,i_b,i_c,\mathrm{\&omega;},\mathrm{\&theta;}\&rsqb;}_k^T,u_k={\&lsqb;u_{ab},u_{bc},u_{ca}\&rsqb;}_k^T,y_k={\&lsqb;i_a,i_b,i_c\&rsqb;}_k^T;$ W _kand v _kbe respectively system noise vector sum measurement noises vector, be zero mean Gaussian white noise, have covariance matrix to be respectively Q and R; U in formula _ab, u _bc, u _cabe respectively the threephase stator line voltage of brshless DC motor, i _a, i _b, i _cfor threephase stator phase current, ω is rotor speed, and θ is rotor-position;

In described step 2, owing to having five state variables and three output variables in motor status equation, so get respectively: Q=diag (Q ₁₁q ₂₂q ₃₃q ₄₄q ₅₅), R=diag (R ₁₁r ₂₂r ₃₃), then fitness function is:

$fitness(Q,R)=\frac{1}{k}\underset{i=0}{\overset{k}{\&Sigma;}}{({\mathrm{\&omega;}}_{ri}^{*}-{\mathrm{\&omega;}}_{ri})}^2---\left(2\right)$

Wherein: and ω _ribe respectively the actual speed of moment i and estimate rotating speed, k equals the ratio in simulation time and sampling period;

In step 3, target function is fitness function fitness (Q, R), and the parameter of corresponding one group of Q and R of each particle, upgrades particle rapidity and position according to formula (1.10) and (1.11):

v _id(t+1)＝wv _id(t)+c ₁r ₁(p _id(t)-x _id(t))+c ₂r ₂(p _gd(t)-x _id(t)) (3)

x _id(t+1)＝x _id(t)+v _id(t+1) (4)

Wherein, p _idfor the optimal solution that each particle searches oneself; p _gdfor the optimal solution that whole population searches; v _idrepresent the speed of i-th particle in d dimension; x _idrepresent i-th position of particle in d dimension; W is inertia weight; r ₁and r ₂for being evenly distributed on the random number between [0,1]; c ₁and c ₂for acceleration limits the factor; The process of renewal speed and position is: the fitness value of each particle respectively with individual extreme value p _idwith global extremum p _gdrelatively, if better, replace individual extreme value or global extremum with fitness value, until it is enough little to reach error, now, obtain some optimization solutions;

In step 4, according to each particle history optimal location P (i) that the process of the optimization solution information generated element initial distribution corresponding particle cluster algorithm in position X (i) that is: ant i solves, that is:

X(i)＝P(i) (5)

Ant pheromones initial distribution utilizes the evaluation function value of each ant position, and formula is:

$T\left(i\right)={\mathrm{ka}}^{-f\left(X_i\right)}---\left(6\right)$

Wherein, k be greater than 0 constant, 0 < a < 1, f (x _i) be the target function value of ant group algorithm, X (i) is the position of ant i, and transition probability P (k, j) is:

$P(k,j)=\frac{T^{\mathrm{\&alpha;}}\left(j\right)Z^{\mathrm{\&beta;}}\left(j\right)}{\underset{i}{\&Sigma;}{T}^{\mathrm{\&alpha;}}\left(i\right)Z^{\mathrm{\&beta;}}\left(i\right)}---\left(7\right)$

Wherein, Δ F _ij=F _j-F _i, F _ifor the functional value of ant i position, T (j) is the pheromones quantity of ant j present position, and a < 1, α and β is non-negative parameter; In ant group algorithm, the transition rule of ant is: move by transition probability and random search in the field of radius r, is exactly carry out global search, if the pheromones increment of ant j adds S pheromones, then by probability transfer

ΔT _j＝ΔT _j+S (8)

In field, search is exactly carry out Local Search, and random search in the field of radius r, if the target function value of reposition is larger than initial value, then getting reposition is ant current location, otherwise do not change, radius r often terminates once to circulate and just reduces, the rear lastest imformation that once circulated element:

T(k)＝dT(k)+ΔT(k) (9)

Wherein: T (k) represents track persistence, 0≤d < 1 for pheromones track intensity, d,

In described step 5, within each sampling period, extended Kalman filter carries out iterative process, and status predication formula is:

${\hat{x}}_{k+1|k}=F_k\left({\hat{x}}_{k|k}\right)\&CenterDot;{\hat{x}}_{k|k}{\mathrm{+G}}_k\left({\hat{x}}_{k|k}\right)\&CenterDot;u_k---\left(10\right)$

Error covariance matrix estimation formulas is:

$P_{k+1|k}=F_k^{\&prime;}\left({\hat{x}}_{k|k}\right)\&CenterDot;P_{k|k}\&CenterDot;{\left(F_k^{\&prime;}\right({\hat{x}}_{k|k}\left)\right)}^T+Q_k---\left(11\right)$

Wherein, $F_k^{\&prime;}=\frac{\&part;\&lsqb;F_k\left({\hat{x}}_{k|k}\right)\&CenterDot;{\hat{x}}_{k|k}\&rsqb;+G_k\left({\hat{x}}_{k|k}\right)\&CenterDot;u_k}{\&part;x_k}{|}_{x_k=x_{k|k}}---\left(12\right)$

Kalman gain formula is:

K _k＝P _k|k-1H ^T[H·P _k|k-1·H ^T+R _k] ^-1(13)

Renewal state equation is:

${\hat{x}}_{k|k}={\hat{x}}_{k|k-1}+K_k\&lsqb;y_k-H\&CenterDot;{\hat{x}}_{k|k-1}\&rsqb;---\left(14\right)$

Upgrade error co-variance matrix:

P _k|k＝[1-K _k·H]P _k|k-1(15)

Just can carry out Kalman's interative computation according to formula (14) to formula (19) and recursion initial value, thus estimate brshless DC motor rotating speed.