来自论文"Predictive Path-Following Control for Fixed-Wing UAVs Using the qLMPC Framework in the Presence of Wind Disturbances"
控制架构
采用的是 ULTRA-Extra无人机,相关参数如下:
这里用于guidance law的无人机运动学模型为:
{
x
˙
p
=
V
a
cos
γ
cos
χ
+
V
w
cos
γ
w
cos
χ
w
y
˙
p
=
V
a
cos
γ
sin
χ
+
V
w
cos
γ
w
sin
χ
w
z
˙
p
=
V
a
sin
γ
+
V
w
sin
γ
w
χ
˙
=
g
tan
ϕ
/
V
a
γ
˙
=
g
(
n
z
cos
ϕ
−
cos
γ
)
/
V
a
\begin{cases} \dot{x}_p = V_a\cos\gamma\cos\chi + V_w\cos\gamma_w\cos\chi_w \\ \dot{y}_p = V_a\cos\gamma\sin\chi + V_w\cos\gamma_w\sin\chi_w \\ \dot{z}_p = V_a\sin\gamma + V_w\sin\gamma_w \\ \dot{\chi} = g\tan\phi/V_a \\ \dot{\gamma} = g(n_z\cos\phi-\cos\gamma)/V_a \end{cases}
⎩
⎨
⎧x˙p=Vacosγcosχ+Vwcosγwcosχwy˙p=Vacosγsinχ+Vwcosγwsinχwz˙p=Vasinγ+Vwsinγwχ˙=gtanϕ/Vaγ˙=g(nzcosϕ−cosγ)/Va
其中状态量为
(
x
p
,
y
p
,
z
p
,
γ
,
χ
)
(x_p,y_p,z_p,\gamma,\chi)
(xp,yp,zp,γ,χ),控制量为
(
V
a
,
n
z
,
ϕ
)
(V_a,n_z,\phi)
(Va,nz,ϕ)。在自动驾驶仪(Autopilot)中,采用 Successive-Loop-Closure (SLC)实现参考量
(
V
a
m
,
n
z
m
,
ϕ
m
)
(V_{a_m},n_{z_m},\phi_m)
(Vam,nzm,ϕm)的信号跟踪:
自动驾驶仪中依然采用横纵向通道的SLC实现控制,相应的控制逻辑如下:
Path Following 最优控制器
对运动学模型进行二阶求导可以得到:
(
x
˙
p
y
˙
p
z
˙
p
χ
˙
γ
˙
x
¨
p
y
¨
p
z
¨
p
χ
¨
γ
¨
)
=
(
O
5
×
5
I
5
−
V
a
cos
γ
sin
χ
−
V
a
sin
γ
cos
χ
V
a
cos
γ
cos
χ
−
V
a
sin
γ
sin
χ
O
5
×
5
O
5
×
3
0
V
a
cos
γ
0
0
0
g
sin
γ
V
a
)
(
x
p
y
p
z
p
χ
γ
x
˙
p
y
˙
p
z
˙
p
χ
˙
γ
˙
)
+
(
O
5
×
3
cos
γ
cos
χ
0
0
cos
γ
sin
χ
0
0
sin
γ
0
0
−
g
tan
ϕ
V
a
2
g
V
a
cos
2
ϕ
0
g
(
cos
γ
−
n
z
cos
ϕ
)
V
a
2
−
g
n
z
sin
ϕ
V
a
g
cos
ϕ
V
a
)
(
V
˙
a
ϕ
˙
n
˙
z
)
\left( \begin{matrix} {{{\dot{x}}}_{p}} \\ {{{\dot{y}}}_{p}} \\ {{{\dot{z}}}_{p}} \\ {\dot{\chi }} \\ {\dot{\gamma }} \\ {{{\ddot{x}}}_{p}} \\ {{{\ddot{y}}}_{p}} \\ {{{\ddot{z}}}_{p}} \\ {\ddot{\chi }} \\ {\ddot{\gamma }} \\ \end{matrix} \right)=\left( \begin{matrix} {{O}_{5\times 5}} & {} & {{I}_{5}} & {} \\ {} & {} & -{{V}_{a}}\cos \gamma \sin \chi & -{{V}_{a}}\sin \gamma \cos \chi \\ {} & {} & {{V}_{a}}\cos \gamma \cos \chi & -{{V}_{a}}\sin \gamma \sin \chi \\ {{O}_{5\times 5}} & {{O}_{5\times 3}} & 0 & {{V}_{a}}\cos \gamma \\ {} & {} & 0 & 0 \\ {} & {} & 0 & \frac{g\sin \gamma }{V_{a}^{{}}} \\ \end{matrix} \right)\left( \begin{matrix} {{x}_{p}} \\ {{y}_{p}} \\ {{z}_{p}} \\ \chi \\ \gamma \\ {{{\dot{x}}}_{p}} \\ {{{\dot{y}}}_{p}} \\ {{{\dot{z}}}_{p}} \\ {\dot{\chi }} \\ {\dot{\gamma }} \\ \end{matrix} \right)+\left( \begin{matrix} {} & {{O}_{5\times 3}} & {} \\ \cos \gamma \cos \chi & 0 & 0 \\ \cos \gamma \sin \chi & 0 & 0 \\ \sin \gamma & 0 & 0 \\ -\frac{g\tan \phi }{V_{a}^{2}} & \frac{g}{{{V}_{a}}{{\cos }^{2}}\phi } & 0 \\ \frac{g(\cos \gamma -{{n}_{z}}\cos \phi )}{V_{a}^{2}} & -\frac{g{{n}_{z}}\sin \phi }{V_{a}^{{}}} & \frac{g\cos \phi }{V_{a}^{{}}} \\ \end{matrix} \right)\left( \begin{align} & {{{\dot{V}}}_{a}} \\ & {\dot{\phi }} \\ & {{{\dot{n}}}_{z}} \\ \end{align} \right)
x˙py˙pz˙pχ˙γ˙x¨py¨pz¨pχ¨γ¨
=
O5×5O5×5O5×3I5−VacosγsinχVacosγcosχ000−Vasinγcosχ−VasinγsinχVacosγ0Vagsinγ
xpypzpχγx˙py˙pz˙pχ˙γ˙
+
cosγcosχcosγsinχsinγ−Va2gtanϕVa2g(cosγ−nzcosϕ)O5×3000Vacos2ϕg−Vagnzsinϕ0000Vagcosϕ
V˙aϕ˙n˙z
这里设
ρ
=
(
γ
,
χ
,
V
a
,
ϕ
,
n
z
)
T
\rho=(\gamma,\chi,V_a,\phi,n_z)^T
ρ=(γ,χ,Va,ϕ,nz)T,
x
=
(
x
p
,
y
p
,
z
p
,
χ
,
γ
,
x
˙
p
,
y
˙
p
,
z
˙
p
,
χ
˙
,
γ
˙
)
T
x=(x_p,y_p,z_p,\chi,\gamma,\dot{x}_p,\dot{y}_p,\dot{z}_p,\dot{\chi},\dot{\gamma})^T
x=(xp,yp,zp,χ,γ,x˙p,y˙p,z˙p,χ˙,γ˙)T,
u
=
(
V
˙
a
,
ϕ
˙
,
n
˙
z
)
T
u=(\dot{V}_a,\dot{\phi},\dot{n}_z)^T
u=(V˙a,ϕ˙,n˙z)T,得到:
x
˙
=
A
v
(
ρ
)
x
+
B
v
(
ρ
)
u
\dot{x}=A_v(\rho)x+B_v(\rho)u
x˙=Av(ρ)x+Bv(ρ)u
假设要跟踪的量为
r
=
(
x
r
,
y
r
,
z
r
)
T
r=(x_r,y_r,z_r)^T
r=(xr,yr,zr)T,构造跟踪向量
e
=
(
x
r
−
x
p
,
y
r
−
y
p
,
z
r
−
z
p
)
T
=
r
−
(
x
p
,
y
p
,
z
p
)
T
e=(x_r-x_p,y_r-y_p,z_r-z_p)^T=r-(x_p,y_p,z_p)^T
e=(xr−xp,yr−yp,zr−zp)T=r−(xp,yp,zp)T,有:
e
˙
=
r
˙
−
(
O
3
×
5
,
I
3
,
O
3
×
2
)
x
\dot{e}=\dot{r}-(O_{3\times 5},I_3,O_{3\times 2})x
e˙=r˙−(O3×5,I3,O3×2)x,得到:
(
x
˙
e
˙
)
=
(
A
v
(
ρ
)
O
10
×
3
−
(
O
3
×
5
∣
I
3
∣
O
3
×
2
)
O
3
)
(
x
e
)
+
(
B
v
(
ρ
)
O
3
×
3
)
u
+
(
O
10
×
1
r
˙
)
\begin{pmatrix} \dot{x} \\ \dot{e} \end{pmatrix} = \begin{pmatrix} A_v(\rho) & O_{10\times3} \\ -(O_{3\times 5}|I_3|O_{3\times 2}) & O_3 \end{pmatrix}\begin{pmatrix} x \\e \end{pmatrix} + \begin{pmatrix} B_v(\rho) \\O_{3\times 3} \end{pmatrix} u + \begin{pmatrix} O_{10\times 1} \\\dot{r} \end{pmatrix}
(x˙e˙)=(Av(ρ)−(O3×5∣I3∣O3×2)O10×3O3)(xe)+(Bv(ρ)O3×3)u+(O10×1r˙)
利用4阶Runge-Kutta法可以将上式可以离散化为一个LPV状态空间方程(linear parameter varying state-space representation):
x
e
,
k
+
1
=
A
e
(
ρ
k
)
x
e
,
k
+
B
e
(
ρ
k
)
u
e
,
k
+
c
r
,
k
x_{e,k+1} = A_e(\rho_k)x_{e,k}+B_e(\rho_k)u_{e,k}+c_{r,k}
xe,k+1=Ae(ρk)xe,k+Be(ρk)ue,k+cr,k
令
P
k
=
(
ρ
k
T
,
ρ
k
+
1
T
,
.
.
.
ρ
k
+
N
−
1
T
)
T
P_k=(\rho_{k}^T,\rho_{k+1}^T,...\rho_{k+N-1}^T)^T
Pk=(ρkT,ρk+1T,...ρk+N−1T)T,
X
k
=
(
x
e
,
k
+
1
T
,
x
e
,
k
+
2
T
,
.
.
.
,
x
e
,
k
+
N
T
)
T
X_k=(x_{e,k+1}^T,x_{e,k+2}^T,...,x_{e,k+N}^T)^T
Xk=(xe,k+1T,xe,k+2T,...,xe,k+NT)T,
U
k
=
(
u
e
,
k
+
1
T
,
u
e
,
k
+
2
T
,
.
.
.
,
u
e
,
k
+
N
T
)
T
U_k=(u_{e,k+1}^T,u_{e,k+2}^T,...,u_{e,k+N}^T)^T
Uk=(ue,k+1T,ue,k+2T,...,ue,k+NT)T,
c
k
=
(
c
r
,
k
+
1
T
,
c
r
,
k
+
2
T
,
.
.
.
c
r
,
k
+
N
T
)
T
c_k = (c_{r,k+1}^T,c_{r,k+2}^T,...c_{r,k+N}^T)^T
ck=(cr,k+1T,cr,k+2T,...cr,k+NT)T得到:
X
k
+
1
=
H
(
P
k
)
X
k
+
S
(
P
k
)
U
k
+
c
k
=
L
k
+
S
(
P
k
)
U
k
X_{k+1}=H(P_k)X_k + S(P_k)U_k+c_k \\ =L_k + S(P_k)U_k
Xk+1=H(Pk)Xk+S(Pk)Uk+ck=Lk+S(Pk)Uk
其中:
x
e
,
k
=
(
x
k
T
,
e
k
T
)
T
x_{e,k}=(x^T_k,e^T_k)^T
xe,k=(xkT,ekT)T,
H
(
P
k
)
=
d
i
a
g
(
[
A
e
(
ρ
k
)
,
A
e
(
ρ
k
+
1
)
.
.
.
A
e
(
ρ
k
+
N
−
1
)
]
)
H(P_k) =diag([A_e(\rho_{k}),A_e(\rho_{k+1})...A_e(\rho_{k+N-1})])
H(Pk)=diag([Ae(ρk),Ae(ρk+1)...Ae(ρk+N−1)]),
S
(
P
k
)
=
d
i
a
g
(
[
B
e
(
ρ
k
)
,
B
e
(
ρ
k
+
1
)
.
.
.
B
e
(
ρ
k
+
N
−
1
)
]
)
S(P_k)=diag([B_e(\rho_{k}),B_e(\rho_{k+1})...B_e(\rho_{k+N-1})])
S(Pk)=diag([Be(ρk),Be(ρk+1)...Be(ρk+N−1)])。采用MPC控制器进行设计时,
k
+
1
k+1
k+1时刻需要优化的目标函数:
J
k
+
1
=
∑
i
=
1
N
(
x
k
+
i
+
1
T
Q
x
k
+
i
+
1
+
u
k
+
i
T
R
u
k
+
i
+
e
k
+
i
+
1
T
T
e
k
+
i
+
1
)
=
X
k
+
1
T
H
X
X
k
+
1
+
U
k
T
H
U
U
k
=
[
L
k
+
S
(
P
k
)
U
k
]
T
H
X
[
L
k
+
S
(
P
k
)
U
k
]
+
U
k
T
H
U
U
k
=
U
k
T
(
S
(
P
k
)
T
H
X
S
(
P
k
)
+
H
U
)
U
k
+
2
L
k
T
H
X
S
(
P
k
)
U
k
+
L
k
T
H
X
L
k
J_{k+1}=\sum_{i=1}^N(x_{k+i+1}^TQx_{k+i+1} + u_{k+i}^TRu_{k+i} + e_{k+i+1}^TTe_{k+i+1}) \\ =X_{k+1}^TH_XX_{k+1} + U_{k}^TH_UU_{k}\\=[L_k + S(P_k)U_k]^TH_X[L_k + S(P_k)U_k] + U_{k}^TH_UU_{k}\\ =U_k^T(S(P_k)^TH_XS(P_k)+H_U)U_k + 2L_k^TH_XS(P_k)U_k + L_k^TH_XL_k
Jk+1=i=1∑N(xk+i+1TQxk+i+1+uk+iTRuk+i+ek+i+1TTek+i+1)=Xk+1THXXk+1+UkTHUUk=[Lk+S(Pk)Uk]THX[Lk+S(Pk)Uk]+UkTHUUk=UkT(S(Pk)THXS(Pk)+HU)Uk+2LkTHXS(Pk)Uk+LkTHXLk
其中:
Q
=
Q
T
>
0
,
P
=
P
T
>
0
,
R
=
R
T
>
0
Q=Q^T>0,P=P^T>0,R=R^T>0
Q=QT>0,P=PT>0,R=RT>0;
H
X
=
d
i
a
g
(
[
Q
,
Q
,
.
.
.
,
Q
]
)
H_X=diag([Q,Q,...,Q])
HX=diag([Q,Q,...,Q]),
H
U
=
d
i
a
g
(
[
R
,
R
,
.
.
.
,
R
]
)
H_U=diag([R,R,...,R])
HU=diag([R,R,...,R])。
而对于控制量
U
k
U_k
Uk和状态量
X
k
+
1
X_{k+1}
Xk+1有限幅,即:
U
min
≤
U
k
≤
U
max
U_{\min}\leq U_k\leq U_{\max}
Umin≤Uk≤Umax,
X
min
≤
X
k
+
1
≤
X
max
X_{\min} \leq X_{k+1} \leq X_{\max}
Xmin≤Xk+1≤Xmax,得到约束:
(
I
−
I
S
(
P
k
)
−
S
(
P
k
)
)
U
k
≤
(
U
max
−
U
min
X
k
+
1
−
L
k
−
X
k
+
1
+
L
k
)
\begin{pmatrix} I \\ -I\\ S(P_k)\\ -S(P_k) \end{pmatrix}U_k \leq \begin{pmatrix} U_{\max} \\ -U_{\min} \\ X_{k+1} - L_k \\ -X_{k+1} + L_k \end{pmatrix}
I−IS(Pk)−S(Pk)
Uk≤
Umax−UminXk+1−Lk−Xk+1+Lk
上面的假设是基于全状态反馈的,也是就是说对于 k + 1 k+1 k+1时刻的在线优化能获取 k k k时刻所有的状态信息和偏差信息。
若观测量为
Y
k
=
C
k
X
k
Y_k = C_kX_k
Yk=CkXk,
Y
min
≤
Y
k
≤
Y
max
Y_{\min}\leq Y_k\leq Y_{\max}
Ymin≤Yk≤Ymax,则上面的约束将被修正为:
(
I
−
I
C
k
S
(
P
k
)
−
C
k
S
(
P
k
)
)
U
k
≤
(
U
max
−
U
min
Y
max
−
C
k
L
k
−
Y
min
+
C
k
L
k
)
\begin{pmatrix} I \\ -I\\ C_kS(P_k)\\ -C_kS(P_k) \end{pmatrix}U_k \leq \begin{pmatrix} U_{\max} \\ -U_{\min} \\ Y_{\max} - C_kL_k \\ -Y_{\min} + C_kL_k \end{pmatrix}
I−ICkS(Pk)−CkS(Pk)
Uk≤
Umax−UminYmax−CkLk−Ymin+CkLk
无论如何,上述问题都可以被转化成QP问题,利用Matlab工具箱中的quadprog
函数进行求解,或者说是在线优化为以下问题:
min
U
k
1
2
U
k
T
F
k
U
k
+
G
k
U
k
A
k
U
k
≤
b
k
\min_{U_k}\frac{1}{2}U_k^TF_kU_k +G_kU_k \\ A_kU_k \leq b_k
Ukmin21UkTFkUk+GkUkAkUk≤bk
附带相应的伪代码如下图所示:
参考文献
@inproceedings{bib:Samir,
title={Predictive Path Following Control for Fixed Wing UAVs Using the qLMPC Framework in the Presence of Wind Disturbances},
author={Rezk, Ahmed S and Calder{\'o}n, Horacio M and Werner, Herbert and Herrmann, Benjamin and Thielecke, Frank},
booktitle={AIAA SCITECH 2024 Forum},
pages={1594},
year={2024}
}