一、加减速法
核心思想：围绕着电机的机械运动方程，通过测量已知量求解惯量
1、原理
机械运动方程：
J
d
ω
m
d
t
=
T
e
−
B
ω
m
−
T
L
(1-1)
J\frac{d\omega_m}{dt}=T_e-B\omega_m-T_L \tag{1-1}
Jdtdωm=Te−Bωm−TL(1-1)
上式中，可以通过测量得到的参数有电磁转矩 T e T_e Te以及电机机械角速度 ω m \omega_m ωm。负载转矩 T L T_L TL、黏滞摩擦系数 B B B、总转矩惯量 J J J是无法通过测量得到的，一般采用简化以及消除等方法排除掉未知项影响。
忽略系统摩擦，机械运动方程简化为：
J
d
ω
m
d
t
=
T
e
−
T
L
(1-2)
J\frac{d\omega_m}{dt}=T_e-T_L \tag{1-2}
Jdtdωm=Te−TL(1-2)
对上式两边取积分离散，得到：
J
(
ω
m
t
−
ω
m
0
)
+
∫
T
L
(
t
)
d
t
=
∑
n
=
1
k
T
e
(
n
)
t
c
(1-3)
J(\omega_{mt}-\omega_{m0})+\int{T_L(t)dt}=\sum_{n=1}^{k}T_{e}(n)t_c \tag{1-3}
J(ωmt−ωm0)+∫TL(t)dt=n=1∑kTe(n)tc(1-3)
在规定的时间
t
n
t_n
tn内，使电机以相同大小的加速度由0r/min做匀加速运动到
ω
\omega
ω，再匀减速到0r/min，在每个电流环控制周期内将电机的
T
e
T_e
Te累加保存。
加速段与减速段方程分别为：
{
J
ω
+
∫
T
L
(
t
)
d
t
=
∑
n
=
1
k
T
e
(
n
)
t
c
−
J
ω
+
∫
T
L
(
t
)
d
t
=
∑
n
=
1
k
T
e
′
(
n
)
t
c
(1-4)
\begin{cases} J\omega+\int{T_L(t)dt}=\sum_{n=1}^{k}T_e(n)t_c\\ -J\omega+\int{T_L(t)dt}=\sum_{n=1}^{k}T^{'}_e(n)t_c \end{cases} \tag{1-4}
{Jω+∫TL(t)dt=∑n=1kTe(n)tc−Jω+∫TL(t)dt=∑n=1kTe′(n)tc(1-4)
通过将两式作差，消掉负载转矩项，最后则可以得到转动惯量表达式：
J
=
∑
n
=
1
k
T
e
(
n
)
t
c
−
∑
n
=
1
k
T
e
′
(
n
)
t
c
2
ω
(1-5)
J=\frac{\sum_{n=1}^{k}T_e(n)t_c-\sum_{n=1}^{k}T^{'}_e(n)t_c}{2\omega} \tag{1-5}
J=2ω∑n=1kTe(n)tc−∑n=1kTe′(n)tc(1-5)
二、朗道离散递推模型参考自适应算法
1、原理
基于模型参考自适应的辨识方法是将实际系统作为参考模型，建立含有未知参数的可调模型。两模型具有相同物理意义的输入与输出量，比较两个模型的输出，通过某种自适应规律调整可调模型的参数，最终实现可调模型输出跟随参考模型输出。当可调模型的输出和参考模型的输出偏差不再改善时，可调模型的未知参数被作为待辨识参数的估计值，即辨识结果。
忽略摩擦阻力，将机械运动方程离散化：
ω m ( k ) = ω m ( k − 1 ) + T J [ T e ( k − 1 ) − T L ( k − 1 ) ] (2-1) \omega_m(k)=\omega_m(k-1)+\frac{T}{J}[T_e(k-1)-T_L(k-1)] \tag{2-1} ωm(k)=ωm(k−1)+JT[Te(k−1)−TL(k−1)](2-1)
ω m ( k − 1 ) = ω m ( k − 2 ) + T J [ T e ( k − 2 ) − T L ( k − 2 ) ] (2-2) \omega_m(k-1)=\omega_m(k-2)+\frac{T}{J}[T_e(k-2)-T_L(k-2)] \tag{2-2} ωm(k−1)=ωm(k−2)+JT[Te(k−2)−TL(k−2)](2-2)
由于采样时间短，电机所带负载转矩变化周期远大于惯量辨识控制周期，可以认为在k-1、k-2时刻负载转矩不变：
T L ( k − 1 ) = T L ( k − 2 ) (2-3) T_L(k-1)=T_L(k-2) \tag{2-3} TL(k−1)=TL(k−2)(2-3)
使用式（2-1）减去（2-2），得到：
ω r ( k ) = 2 ω m ( k − 1 ) − ω m ( k − 2 ) + b ( k ) Δ T e ( k − 1 ) (2-4) \omega_r(k)=2\omega_m(k-1)-\omega_m(k-2)+b(k)\Delta T_e(k-1) \tag{2-4} ωr(k)=2ωm(k−1)−ωm(k−2)+b(k)ΔTe(k−1)(2-4)
式（2-4）中，
b
(
k
)
=
T
J
b(k)=\frac{T}{J}
b(k)=JT,
Δ
T
e
(
k
−
1
)
=
T
e
(
k
−
1
)
−
T
e
(
k
−
2
)
\Delta T_e(k-1)=T_e(k-1)-T_e(k-2)
ΔTe(k−1)=Te(k−1)−Te(k−2)。
将式（2-4）作为参考模型，得到可调模型方程为：
ω g ( k ) = 2 ω m ( k − 1 ) − ω m ( k − 2 ) + b g ( k ) Δ T e ( k − 1 ) (2-5) \omega_g(k)=2\omega_m(k-1)-\omega_m(k-2)+b_g(k)\Delta T_e(k-1) \tag{2-5} ωg(k)=2ωm(k−1)−ωm(k−2)+bg(k)ΔTe(k−1)(2-5)
其中， ω g \omega_g ωg为估计速度， b g b_g bg为估计变量，参考模型和可调模型的输出误差为：
e ( k ) = ω m ( k ) − ω g ( k ) (2-6) e(k)=\omega_m(k)-\omega_g(k) \tag{2-6} e(k)=ωm(k)−ωg(k)(2-6)
Popov超稳定性理论利用函数来判断系统的全局和局部渐进的稳定性，无需求解系统的微分方程，可以准确得到系统的自适应律。由Popov超稳定性理论设计的转动惯量自适应律为：
b g ( k ) = b g ( k − 1 ) + β Δ T e ( k − 1 ) 1 + Δ T e ( k − 1 ) 2 e ( k ) (2-7) b_g(k)=b_g(k-1)+\beta\frac{\Delta T_e(k-1)}{1+\Delta T_e(k-1)^2}e(k) \tag{2-7} bg(k)=bg(k−1)+β1+ΔTe(k−1)2ΔTe(k−1)e(k)(2-7)
其中，
β
\beta
β为自适应增益系数。
β
\beta
β越大，辨识过程收敛速度越快，但辨识结果波动大；
β
\beta
β越小，辨识精度越高，但收敛时间变长。
对于参数
β
\beta
β的选取，可以考虑变增益和PI等方式。
-
变增益
当系统转动惯量发生变化时，增益系数 β \beta β跟随变化。转动惯量突变， β \beta β迅速增大，当跟踪到真实值后，再逐步减小 β \beta β。
β ( k ) = β ( k − 1 ) − β k − 1 2 Δ T e ( k − 1 ) 2 λ + β ( k − 1 ) Δ T e ( k − 1 ) 2 (2-8) \beta(k)=\beta(k-1)-\frac{\beta{k-1}^2\Delta T_e(k-1)^2}{\lambda+\beta(k-1)\Delta T_e(k-1)^2} \tag{2-8} β(k)=β(k−1)−λ+β(k−1)ΔTe(k−1)2βk−12ΔTe(k−1)2(2-8)
其中， λ \lambda λ为调整系数， λ \lambda λ越小 β \beta β衰减越快。
-
PI
使用PI控制器进行 β \beta β的动态调整。
β = β m i n + K ∑ k = 1 n e ( k ) (2-9) \beta=\beta_{min}+K\sum_{k=1}^{n}e(k) \tag{2-9} β=βmin+Kk=1∑ne(k)(2-9)
其中， K K K为积分系数。
三、最小二乘递推
1、概述
最小二乘辨识算法（RLS）可以解决线性定常系统、线性时变系统、含有色噪声的线性系统等参数辨识问题。在利用最小二乘法对电机进行参数辨识过程中，主要将电机的非线性模型线性化，得到与电机参数有着直接关系的线性化模型，可以实时辨识出电机的所有主要参数。
2、原理
对于一个典型的离散单输入单输出系统，其传递函数表示如下：
H ( Z ) = Y ( Z ) X ( Z ) = b 1 z − 1 + b 2 z − 2 + ⋯ + b n z − n 1 + a 1 z − 1 + a 2 z − 2 + ⋯ + a m z − m (3-1) H(Z)=\frac{Y(Z)}{X(Z)}=\frac{b_1z^{-1}+b_2z^{-2}+\cdots +b_nz^{-n}}{1+a_1z^{-1}+a_2z^{-2}+\cdots +a_mz^{-m}} \tag{3-1} H(Z)=X(Z)Y(Z)=1+a1z−1+a2z−2+⋯+amz−mb1z−1+b2z−2+⋯+bnz−n(3-1)
其中，
n
≤
m
;
z
=
e
s
T
,
T
n\le m;z=e^{sT},T
n≤m;z=esT,T为采样周期。n表示输出量的采样次数，
b
1
、
b
2
、
⋯
b
n
b_1、b_2、\cdots b_n
b1、b2、⋯bn表示输出量在此次采样时刻的因数，m表示输入量的采样次数，
a
1
、
a
2
、
⋯
a
m
a_1、a_2、\cdots a_m
a1、a2、⋯am表示输入量在此次采样时刻的因数。
将式（3-1）写成差分方程形式：
{ y ( k ) + ∑ 1 n a i y ( k − i ) = ∑ 1 n b i x ( k − i ) y ( k ) = y ( k T ) x ( k ) = x ( k T ) (3-2) \begin{cases} y(k)+\sum_{1}^{n}a_iy(k-i)=\sum_{1}^{n}b_ix(k-i)\\ y(k)=y(kT)\\ x(k)=x(kT) \tag{3-2} \end{cases} ⎩⎪⎨⎪⎧y(k)+∑1naiy(k−i)=∑1nbix(k−i)y(k)=y(kT)x(k)=x(kT)(3-2)
在实际系统中，模拟电路的精度不足、系统模型误差以及噪声干扰等均可以归为RLS理论中的误差e(k)。对于n阶差分方程，需要前n个时刻的输入输出量，则式（3-2）可以整理成如下：
y ( k ) = − ∑ 1 n a i y ( k − i ) + ∑ 1 n b i x ( k − i ) + e ( k ) = ∑ 1 n Φ k T Θ ^ + e ( k ) (3-3) y(k)=-\sum_{1}^{n}a_iy(k-i)+\sum_{1}^{n}b_ix(k-i)+e(k)\\ =\sum_{1}^n\Phi_k^T\hat\Theta+e(k) \tag{3-3} y(k)=−1∑naiy(k−i)+1∑nbix(k−i)+e(k)=1∑nΦkTΘ^+e(k)(3-3)
其中，
Φ
k
T
=
(
−
y
(
k
−
1
)
,
⋯
,
−
y
(
k
−
n
)
,
x
(
k
−
1
)
,
⋯
,
x
(
k
−
n
)
)
\Phi_k^T=(-y(k-1),\cdots ,-y(k-n),x(k-1),\cdots ,x(k-n))
ΦkT=(−y(k−1),⋯,−y(k−n),x(k−1),⋯,x(k−n))，为输出量矩阵。
Θ
^
k
=
(
a
1
,
⋯
,
a
n
,
b
1
,
⋯
,
b
n
)
T
\hat\Theta_k=(a_1,\cdots ,a_n,b_1,\cdots ,b_n)^T
Θ^k=(a1,⋯,an,b1,⋯,bn)T，为因数矩阵。
最小二乘的指标是使：
C
=
∑
1
n
e
(
k
)
2
=
∑
1
n
[
y
(
k
)
−
Φ
k
T
Θ
^
k
]
2
C=\sum_1^ne(k)^2=\sum_1^n[y(k)-\Phi_k^T\hat\Theta_k]^2
C=∑1ne(k)2=∑1n[y(k)−ΦkTΘ^k]2最小。
令
E
=
[
e
(
1
)
,
e
(
2
)
,
⋯
,
e
(
n
)
]
T
E=[e(1),e(2),\cdots ,e(n)]^T
E=[e(1),e(2),⋯,e(n)]T，那么：
C = E T E = ( Y k − Φ k T Θ ^ k ) T ( Y k − Φ k T Θ ^ k ) (3-4) C=E^TE=(Y_k-\Phi_k^T\hat\Theta_k)^T(Y_k-\Phi_k^T\hat\Theta_k) \tag{3-4} C=ETE=(Yk−ΦkTΘ^k)T(Yk−ΦkTΘ^k)(3-4)
为了使误差和最小，对C求导： ∂ C ∂ Θ ^ k = − 2 Φ T Y k + 2 Φ T Φ Θ ^ k = 0 \frac{\partial C}{\partial \hat\Theta_k}=-2\Phi^TY_k+2\Phi^T\Phi\hat\Theta_k=0 ∂Θ^k∂C=−2ΦTYk+2ΦTΦΘ^k=0，则：
Θ ^ k = ( Φ k T Φ k ) − 1 Φ k T Y k (3-5) \hat\Theta_k=(\Phi_k^T\Phi_k)^{-1}\Phi_k^TY_k \tag{3-5} Θ^k=(ΦkTΦk)−1ΦkTYk(3-5)
递推最小二乘在最小二乘的基础上进行推广，其基本思想可以表示为:
Θ ^ k = Θ ^ k − 1 + 修 正 项 (3-6) \hat\Theta_k=\hat\Theta_{k-1}+修正项 \tag{3-6} Θ^k=Θ^k−1+修正项(3-6)
记 P k − 1 = Φ k T Φ k P_k^{-1}=\Phi_k^T\Phi_k Pk−1=ΦkTΦk，经过推导，最终得到递推方程为：
{ Θ ^ k = Θ ^ k − 1 + K k ( y ( k ) − Φ k T Θ ^ k − 1 ) P k = ( I − K k Φ k T ) P k − 1 K k = P k − 1 Φ k 1 + Φ k T P k − 1 Φ k (3-7) \begin{cases} \hat\Theta_k=\hat\Theta_{k-1}+K_k(y(k)-\Phi_k^T\hat\Theta_{k-1})\\ P_k=(I-K_k\Phi_k^T)P_{k-1}\\ K_k=\frac{P_{k-1}\Phi_k}{1+\Phi_k^TP_{k-1}\Phi_k} \tag{3-7} \end{cases} ⎩⎪⎨⎪⎧Θ^k=Θ^k−1+Kk(y(k)−ΦkTΘ^k−1)Pk=(I−KkΦkT)Pk−1Kk=1+ΦkTPk−1ΦkPk−1Φk(3-7)
整个递推最小二乘法的计算步骤可以归纳如下：
（1）设置初始值
Θ
^
0
,
P
0
\hat\Theta_0,P_0
Θ^0,P0，一般采用如下初值：
{ P 0 = α I Θ ^ 0 = z e r o (3-8) \begin{cases} P_0=\alpha I\\ \hat\Theta_0=zero \tag{3-8} \end{cases} {P0=αIΘ^0=zero(3-8)
其中，
α
\alpha
α为充分大的正实数，一般取
1
0
4
∼
1
0
6
10^4\sim 10^6
104∼106,
I
I
I为单位矩阵，
z
e
r
o
zero
zero为零向量。
（2）确定输出矩阵，组建
Φ
k
\Phi_k
Φk。
（3）对数据进行采集。
（4）由公式（3-7）进行递推运算，计算参数值。
随着运算次数的增多， P k P_k Pk逐渐趋近于0，从而失去其修正能力，可以考虑引入遗忘因子 λ \lambda λ。 λ \lambda λ介于0和1之间，其实质是在迭代过程中引入权值概念，在迭代过程中增加此时刻的数据权值，以此来减缓 P k 、 K k P_k、K_k Pk、Kk趋近于0的速度。引入遗忘因子后， P k 、 K k P_k、K_k Pk、Kk计算公式发生变化：
{ P k = 1 λ ( I − K k Φ k T ) P k − 1 K k = P k − 1 Φ k λ + Φ k T P k − 1 Φ k (3-9) \begin{cases} P_k=\frac{1}{\lambda}(I-K_k\Phi_k^T)P_{k-1}\\ K_k=\frac{P_{k-1}\Phi_k}{\lambda+\Phi_k^TP_{k-1}\Phi_k} \tag{3-9} \end{cases} {Pk=λ1(I−KkΦkT)Pk−1Kk=λ+ΦkTPk−1ΦkPk−1Φk(3-9)
3、运动方程离散
机械运动方程（1-1）拉式变换为：
T e ( s ) − T L ( s ) = B ω m ( s ) + J s ω m ( s ) (3-10) T_e(s)-T_L(s)=B\omega_m(s)+Js\omega_m(s) \tag{3-10} Te(s)−TL(s)=Bωm(s)+Jsωm(s)(3-10)
令 y ( s ) = ω m ( s ) , u ( s ) = T e ( s ) − T L ( s ) y(s)=\omega_m(s),u(s)=T_e(s)-T_L(s) y(s)=ωm(s),u(s)=Te(s)−TL(s)，则运动方程传递函数为：
H ( s ) = y ( s ) u ( s ) = 1 J s + B J (3-11) H(s)=\frac{y(s)}{u(s)}=\frac{\frac{1}{J}}{s+\frac{B}{J}} \tag{3-11} H(s)=u(s)y(s)=s+JBJ1(3-11)
已知零阶保持器传递函数：
G 0 ( s ) = 1 − e − T e s s (3-12) G_0(s)=\frac{1-e^{-T_es}}{s} \tag{3-12} G0(s)=s1−e−Tes(3-12)
带有零阶保持器的开环传递函数为：
H
(
Z
)
=
(
1
−
z
−
1
)
[
1
B
∗
z
(
1
−
e
−
B
T
c
/
J
)
(
z
−
1
)
(
z
−
e
−
B
T
c
/
J
)
]
)
=
z
−
1
(
1
−
e
−
B
T
c
/
J
)
B
(
1
−
Z
−
1
e
−
B
T
c
/
J
)
=
ω
m
(
z
)
T
e
(
z
)
−
T
L
(
z
)
(3-13)
H(Z)=(1-z^{-1})[\frac{1}{B}*\frac{z(1-e^{-BT_c/J})}{(z-1)(z-e^{-BT_c/J})}])\\ =\frac{z^{-1}(1-e^{-BT_c/J})}{B(1-Z^{-1}e^{-BT_c/J})} \\ =\frac{\omega_m(z)}{T_e(z)-T_L(z)} \tag{3-13}
H(Z)=(1−z−1)[B1∗(z−1)(z−e−BTc/J)z(1−e−BTc/J)])=B(1−Z−1e−BTc/J)z−1(1−e−BTc/J)=Te(z)−TL(z)ωm(z)(3-13)
其中 T c T_c Tc为采样周期。对上式进行变形：
ω m ( k ) − e − B T c / J ω m ( k − 1 ) = 1 − e − B T c / J B [ T e ( k − 1 ) − T L ( k − 1 ) ] (3-14) \omega_m(k)-e^{-BT_c/J}\omega_m(k-1)=\frac{1-e^{-BT_c/J}}{B}[T_e(k-1)-T_L(k-1)] \tag{3-14} ωm(k)−e−BTc/Jωm(k−1)=B1−e−BTc/J[Te(k−1)−TL(k−1)](3-14)
为了简化表示，令：
{ b = − e B T c / J a = ( 1 + b ) / B c = a T L ( k − 1 ) (3-15) \begin{cases} b=-e^{BT_c/J}\\ a=(1+b)/B\\ c=aT_L(k-1) \end{cases} \tag{3-15} ⎩⎪⎨⎪⎧b=−eBTc/Ja=(1+b)/Bc=aTL(k−1)(3-15)
那么式（3-14）可以表示为：
ω m ( k ) = a T e ( k − 1 ) − b ω m ( k − 1 ) + c ∗ ( − 1 ) (3-16) \omega_m(k)=aT_e(k-1)-b\omega_m(k-1)+c*(-1) \tag{3-16} ωm(k)=aTe(k−1)−bωm(k−1)+c∗(−1)(3-16)
用最小二乘的方式表示：
{ y ( k ) = Φ k T Θ y ( k ) = ω m ( k ) Φ k T = [ T e ( k − 1 ) , − ω m ( k − 1 ) , − 1 ] Θ = [ a , b , c ] T (3-17) \begin{cases} y(k)=\Phi_k^T\Theta\\ y(k)=\omega_m(k)\\ \Phi_k^T=[T_e(k-1),-\omega_m(k-1),-1]\\ \Theta=[a,b,c]^T \end{cases} \tag{3-17} ⎩⎪⎪⎪⎨⎪⎪⎪⎧y(k)=ΦkTΘy(k)=ωm(k)ΦkT=[Te(k−1),−ωm(k−1),−1]Θ=[a,b,c]T(3-17)
由于 T c T_c Tc足够小，所以可以有以下近似：
{ a = 1 + b B = 1 − e B T c / J B ≈ T c J b = − e B T c / J = lim T c → 0 ( − e B T c / J ) = − 1 c = a T L ( k − 1 ) ≈ T c J T L ( k − 1 ) (3-18) \begin{cases} a=\frac{1+b}{B}=\frac{1-e^{BT_c/J}}{B}\approx\frac{T_c}{J}\\ b=-e^{BT_c/J}=\lim_{T_c\to 0}(-e^{BT_c/J})=-1\\ c=aT_L(k-1)\approx\frac{T_c}{J}T_L(k-1) \end{cases} \tag{3-18} ⎩⎪⎨⎪⎧a=B1+b=B1−eBTc/J≈JTcb=−eBTc/J=limTc→0(−eBTc/J)=−1c=aTL(k−1)≈JTcTL(k−1)(3-18)
根据最小二乘法的辨识结果 Θ ^ = [ a ^ ( k ) , b ^ ( k ) , c ^ ( k ) ] T \hat\Theta=[\hat a(k),\hat b(k),\hat c(k)]^T Θ^=[a^(k),b^(k),c^(k)]T，可以计算出转动惯量 J ^ \hat J J^、摩擦系数 B ^ \hat B B^、负载转矩 T ^ L \hat T_L T^L。
4043
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
