一、连续时间系统
1、系统增益指标
(1)系统性能指标定义
对于线性时不变的连续时间系统:
x
˙
(
t
)
=
A
x
(
t
)
+
B
w
(
t
)
z
(
t
)
=
C
x
(
t
)
+
D
w
(
t
)
\dot{x}(t)=Ax(t)+Bw(t)\\z(t)=Cx(t)+Dw(t)
x˙(t)=Ax(t)+Bw(t)z(t)=Cx(t)+Dw(t)
其中,
x
(
t
)
∈
R
n
x(t)\in R^n
x(t)∈Rn是系统的状态,
w
(
t
)
∈
R
q
w(t)\in R^q
w(t)∈Rq是外部扰动输入,
z
(
t
)
∈
R
r
z(t)\in R^r
z(t)∈Rr是系统被调输出。如果对于某一类外界扰动信号w(t),系统的被调输出z(t)能够保持小,则认为这样的系统具有好的性能,反映了系统抑制外部扰动的能力。
考虑系统的增益 Γ \Gamma Γ: Γ = s u p w ≠ 0 s i z e ( z ) ‾ s i z e ( w ) \Gamma=\mathop{{sup}}\limits_{w\ne0}\mathop{\underline{size(z)}}\limits_{size(w)} Γ=w=0supsize(w)size(z)
对于平方可积的信号,定义:
∣
∣
f
∣
∣
2
=
(
∫
0
∞
∣
∣
f
(
t
)
∣
∣
2
d
t
)
1
/
2
||f||_2=(\int_0^{\infty}||f(t)||^2dt)^{1/2}
∣∣f∣∣2=(∫0∞∣∣f(t)∣∣2dt)1/2
其中是向量的欧式范数,这样定义的
∣
∣
f
∣
∣
2
||f||_2
∣∣f∣∣2正好是信号f的能量。
对于幅值有界的信号,定义:
∣
∣
f
∣
∣
∞
=
s
u
p
t
≥
0
∣
∣
f
(
t
)
∣
∣
||f||_{\infty}=\mathop{{sup}}\limits_{t\ge0}||f(t)||
∣∣f∣∣∞=t≥0sup∣∣f(t)∣∣
其中是向量的欧式范数,当f是一个标量信号时,
∣
∣
f
∣
∣
∞
||f||_{\infty}
∣∣f∣∣∞是f的峰值。
因此,系统的一些性能指标可以定义为:
- IE增益: Γ i e = s u p w ( t ) = w 0 δ ( t ) , ∣ ∣ w 0 ∣ ∣ ≤ 1 ∣ ∣ Z ∣ ∣ 2 \Gamma_{ie}=\mathop{{sup}}\limits_{w(t)=w_0\delta(t),||w_0||\le1}||Z||_2 Γie=w(t)=w0δ(t),∣∣w0∣∣≤1sup∣∣Z∣∣2
- EP增益:
Γ e p = s u p ∣ ∣ w ∣ ∣ 2 ≤ 1 ∣ ∣ Z ∣ ∣ ∞ \Gamma_{ep}=\mathop{{sup}}\limits_{||w||_2\le1}||Z||_{\infty} Γep=∣∣w∣∣2≤1sup∣∣Z∣∣∞ - EE增益:
Γ e e = s u p ∣ ∣ w ∣ ∣ 2 ≤ 1 ∣ ∣ Z ∣ ∣ 2 \Gamma_{ee}=\mathop{{sup}}\limits_{||w||_2\le1}||Z||_2 Γee=∣∣w∣∣2≤1sup∣∣Z∣∣2 - PP增益:
Γ p p = s u p ∣ ∣ w ∣ ∣ ∞ ≤ 1 ∣ ∣ Z ∣ ∣ ∞ \Gamma_{pp}=\mathop{{sup}}\limits_{||w||_{\infty}\le1}||Z||_{\infty} Γpp=∣∣w∣∣∞≤1sup∣∣Z∣∣∞
(2)性能指标求解
那么,如何求这些性能指标呢?
-
求IE增益:
定理:对于原线性时不变连续时间系统,如果该系统是严格真的(D=0)和渐进稳定(矩阵A特征值均具有负实部)的,那么系统的IE增益是有限的,并且 Γ i e = ∣ ∣ B T Y B ∣ ∣ 1 / 2 \Gamma_{ie}=||B^TYB||^{1/2} Γie=∣∣BTYB∣∣1/2,其中矩阵的范数取为谱范数,即矩阵的最大奇异值,即A的谱范数为 A T A A^TA ATA的最大特征值的平方根。矩阵Y是以下李雅普诺夫方程的解: Y A + A T Y + C T C = 0 YA+A^TY+C^TC=0 YA+ATY+CTC=0
该定理可以转化为以下的优化问题: m i n γ s . t . P A + A T P + C T C < 0 B T P B ≤ γ I P > 0 min~~ \gamma\\ s.t.~~~~~~~PA+A^TP+C^TC<0\\B^TPB\le\gamma I\\P>0 min γs.t. PA+ATP+CTC<0BTPB≤γIP>0
该优化问题有一个最优值 γ ∗ \gamma^* γ∗,则 Γ i e = γ ∗ \Gamma_{ie}=\sqrt{\gamma^*} Γie=γ∗。该优化问题是一个具有线性不等式约束和线性目标函数的凸优化问题,可以应用LMI工具箱中的求解器mincx来求解。 -
求EP增益:
定理:对于原线性时不变连续时间系统,如果该系统是严格真的(D=0)和渐进稳定(矩阵A特征值均具有负实部)的,那么系统的EP增益是有限的,并且 Γ e p = ∣ ∣ C X C T ∣ ∣ 1 / 2 \Gamma_{ep}=||CXC^T||^{1/2} Γep=∣∣CXCT∣∣1/2,其中矩阵的范数取为谱范数,即矩阵的最大奇异值,即A的谱范数为 A T A A^TA ATA的最大特征值的平方根。矩阵X是以下李雅普诺夫方程的解: A X + X A T + B B T = 0 AX+XA^T+BB^T=0 AX+XAT+BBT=0
该定理可以转化为以下的优化问题: m i n γ s . t . A Q + Q A T + B B T < 0 C Q C T ≤ γ I Q > 0 min~~ \gamma\\ s.t.~~~~~~~AQ+QA^T+BB^T<0\\CQC^T\le\gamma I\\Q>0 min γs.t. AQ+QAT+BBT<0CQCT≤γIQ>0
该优化问题有一个最优值 γ ∗ \gamma^* γ∗,则 Γ e p = γ ∗ \Gamma_{ep}=\sqrt{\gamma^*} Γep=γ∗。该优化问题是一个具有线性不等式约束和线性目标函数的凸优化问题,可以应用LMI工具箱中的求解器mincx来求解。 -
求EE增益:
定理:对于原线性时不变连续时间系统,设 γ \gamma γ是一个给定的常数,则下列条件是等价的:
[1]系统渐进稳定,且 Γ e e < γ \Gamma_{ee}<\gamma Γee<γ
[2]存在一个对称矩阵P>0,使得: [ A T P + P A P B C T B T P − γ I D T C D − γ I ] < 0 \begin{bmatrix} A^TP+PA & PB & C^T\\\\ B^TP & -\gamma I & D^T\\\\ C & D & -\gamma I\\\\ \end{bmatrix}<0 ⎣⎢⎢⎢⎢⎢⎢⎡ATP+PABTPCPB−γIDCTDT−γI⎦⎥⎥⎥⎥⎥⎥⎤<0
条件[2]是一个线性矩阵不等式,因此可以应用LMI工具箱中的求解器feasp来判断系统增益 Γ e e \Gamma_{ee} Γee是否满足给定的约束条件。 -
求PP增益:
定理:对于原线性时不变连续时间系统,设 γ \gamma γ是一个给定的标量,如果存在对称矩阵R>0、标量 γ > 0 \gamma>0 γ>0和 μ > 0 \mu>0 μ>0,使得: [ A T R + R A + γ R R B B T R − μ I ] < 0 \begin{bmatrix} A^TR+RA+\gamma R & RB\\\\ B^TR & -\mu I\\\\ \end{bmatrix}<0 ⎣⎢⎢⎡ATR+RA+γRBTRRB−μI⎦⎥⎥⎤<0
和 [ λ R 0 C T 0 ( γ − μ ) I D T C D γ I ] > 0 \begin{bmatrix} \lambda R & 0 & C^T\\\\ 0 & (\gamma-\mu) I & D^T\\\\ C & D & \gamma I\\\\ \end{bmatrix}>0 ⎣⎢⎢⎢⎢⎢⎢⎡λR0C0(γ−μ)IDCTDTγI⎦⎥⎥⎥⎥⎥⎥⎤>0那么 Γ p p < γ \Gamma_{pp}<\gamma Γpp<γ。
注意到上述不等式并不是线性矩阵不等式。
2、 H 2 H_2 H2性能
传递函数矩阵
T
(
s
)
=
C
(
s
I
−
A
)
−
1
B
+
D
T(s)=C(sI-A)^{-1}B+D
T(s)=C(sI−A)−1B+D的
H
2
H_2
H2范数定义为:
∣
∣
T
(
s
)
∣
∣
2
=
T
r
a
c
e
(
∫
−
∞
+
∞
T
(
j
w
)
T
∗
(
j
w
)
d
w
/
2
π
)
1
/
2
||T(s)||_2=Trace(\int_{-\infty}^{+\infty}T(jw)T^*(jw)dw/2\pi)^{1/2}
∣∣T(s)∣∣2=Trace(∫−∞+∞T(jw)T∗(jw)dw/2π)1/2其中
T
∗
(
j
w
)
T^*(jw)
T∗(jw)表示矩阵
T
(
j
w
)
T(jw)
T(jw)的共轭转置,Trace代表矩阵的迹(n×n矩阵A的主对角线上各个元素的总和)。
T的
H
2
H_2
H2范数的平方等于系统脉冲响应的总的输出能量,还等于系统在白噪声输入信号激励下的稳态输出方差。
对于单输入单输出系统,
∣
∣
T
(
s
)
∣
∣
2
=
Γ
i
e
=
Γ
e
p
||T(s)||_2=\Gamma_{ie}=\Gamma_{ep}
∣∣T(s)∣∣2=Γie=Γep。
定理:对于原线性时不变连续时间系统,如果该系统是渐进稳定的,则:
[1]
∣
∣
T
∣
∣
2
<
∞
||T||_2<\infty
∣∣T∣∣2<∞当且仅当D=0;
[2]如果D=0,那么以下结论等价:
(i)
∣
∣
T
∣
∣
2
<
γ
||T||_2<\gamma
∣∣T∣∣2<γ
(ii)存在对称矩阵X>0,使得
A
X
+
X
A
T
+
B
B
T
<
0
,
T
r
a
c
e
(
C
X
C
T
)
<
γ
2
AX+XA^T+BB^T<0,Trace(CXC^T)<\gamma^2
AX+XAT+BBT<0,Trace(CXCT)<γ2
(iii)存在对称矩阵Y>0,使得
A
T
Y
+
Y
A
+
C
T
C
<
0
,
T
r
a
c
e
(
B
T
Y
B
)
<
γ
2
A^TY+YA+C^TC<0,Trace(B^TYB)<\gamma^2
ATY+YA+CTC<0,Trace(BTYB)<γ2
根据该定理,可以应用LMI工具箱的求解器feasp来验证系统是否满足给定的
H
2
H_2
H2范数约束条件。
3、 H ∞ H_\infty H∞性能
传递函数矩阵
T
(
s
)
=
C
(
s
I
−
A
)
−
1
B
+
D
T(s)=C(sI-A)^{-1}B+D
T(s)=C(sI−A)−1B+D的
H
∞
H_\infty
H∞范数定义为:
∣
∣
T
(
s
)
∣
∣
∞
=
s
u
p
w
σ
m
a
x
(
T
(
j
w
)
)
||T(s)||_\infty=\mathop{{sup}}\limits_w\sigma_{max}(T(jw))
∣∣T(s)∣∣∞=wsupσmax(T(jw))即系统频率响应的最大奇异值的峰值。
∣
∣
T
(
s
)
∣
∣
∞
=
Γ
e
e
||T(s)||_\infty=\Gamma_{ee}
∣∣T(s)∣∣∞=Γee。
因此,通过求解优化问题:
m
i
n
γ
s
.
t
.
[
A
T
P
+
P
A
P
B
C
T
B
T
P
−
γ
I
D
T
C
D
−
γ
I
]
<
0
P
>
0
min~~ \gamma\\ s.t.~~\begin{bmatrix} A^TP+PA & PB & C^T\\\\ B^TP & -\gamma I & D^T\\\\ C & D & -\gamma I\\\\ \end{bmatrix}<0\\P>0
min γs.t. ⎣⎢⎢⎢⎢⎢⎢⎡ATP+PABTPCPB−γIDCTDT−γI⎦⎥⎥⎥⎥⎥⎥⎤<0P>0可以得到系统最优
H
∞
H_\infty
H∞性能分析问题的解。由于该问题是一个具有线性矩阵不等式约束和线性目标函数的凸优化问题,因此可以应用LMI工具箱的求解器mincx来求解。
4、举例说明
某二阶系统:
G
(
s
)
=
k
s
2
+
2
ζ
ω
s
+
ω
2
G(s)=\frac{k}{s^2+2\zeta\omega s+\omega^2}
G(s)=s2+2ζωs+ω2k该系统的一个状态空间实现为:
[
A
B
C
D
]
=
[
−
2
ζ
ω
−
ω
2
1
1
0
0
0
k
0
]
\begin{bmatrix} A & B\\ C & D \end{bmatrix}=\begin{bmatrix} -2\zeta\omega & -\omega^2 & 1\\ 1 & 0 & 0\\ 0 & k & 0 \end{bmatrix}
[ACBD]=⎣⎡−2ζω10−ω20k100⎦⎤
系统1:
G
1
(
s
)
:
ζ
=
0.1
,
ω
=
1
,
k
=
1
G_1(s):\zeta=0.1,\omega=1,k=1
G1(s):ζ=0.1,ω=1,k=1;
系统2:
G
2
(
s
)
:
ζ
=
1
,
ω
=
3
,
k
=
20
G_2(s):\zeta=1,\omega=3,k=20
G2(s):ζ=1,ω=3,k=20;
(1)IE增益
求解该系统的IE增益代码实现为:
clc;clear all
%% 定义系统状态空间矩阵
A1=[-0.2,-1;1,0];B1=[1;0];C1=[0,1]; % 系统1
A2=[-6,-9;1,0];B2=[1;0];C2=[0,20]; % 系统2
%% 求解IE增益方法1:求解李雅普诺夫方程 YA+A'Y+C'C=0
% matlab X=lyap(A,C)可用于求解李雅普诺夫方程AX+XA(^T)=-C
CC1=C1'*C1;AA1=A1';
Y1=lyap(AA1,CC1);
ie_1=sqrt((norm(B1'*Y1*B1,2)))
CC2=C2'*C2;AA2=A2';
Y2=lyap(AA2,CC2);
ie_2=sqrt((norm(B2'*Y2*B2,2)))
%% 求解IE增益方法2:转化为优化问题用LMI工具箱求解
setlmis([]); % 定义一个矩阵不等式系统
X=lmivar(1,[2 1]); % P矩阵是一个2*2的对称矩阵
r=lmivar(1,[1 1]); % gamma
lmiterm([1 1 1 X],1,A1,'s'); % PA+A'P
lmiterm([1 1 1 0],CC1); % C'C
lmiterm([2 1 1 X],B1',B1); % B'PB
lmiterm([-2 1 1 r],1,1); % -gammaI
lmiterm([-3 1 1 X],1,1); % P
lmisys=getlmis;
n = decnbr(lmisys); % 系统决策变量个数
c = zeros(n,1); % 确定向量c的维数
for j=1:n
[rj]=defcx(lmisys,j,r);
c(j)=trace(rj);
end
[copt, xopt]=mincx(lmisys,c);
xx=dec2mat(lmisys,xopt,X);
gamma=dec2mat(lmisys,xopt,r);
ie_11=sqrt(gamma)
%2
setlmis([]); % 定义一个矩阵不等式系统
X=lmivar(1,[2 1]); % P矩阵是一个2*2的对称矩阵
r=lmivar(1,[1 1]); % gamma
lmiterm([1 1 1 X],1,A2,'s'); % PA+A'P
lmiterm([1 1 1 0],CC2); % C'C
lmiterm([2 1 1 X],B2',B2); % B'PB
lmiterm([-2 1 1 r],1,1); % -gammaI
lmiterm([-3 1 1 X],1,1); % P
lmisys=getlmis;
n = decnbr(lmisys); % 系统决策变量个数
c = zeros(n,1); % 确定向量c的维数
for j=1:n
[rj]=defcx(lmisys,j,r);
c(j)=trace(rj);
end
[copt, xopt]=mincx(lmisys,c);
xx=dec2mat(lmisys,xopt,X);
gamma=dec2mat(lmisys,xopt,r);
ie_22=sqrt(gamma)
输出为:
系统1的ie_1=1.5811(李雅普诺夫方程);ie_11=1.5818(LMI工具箱);
系统2的ie_2=1.9245(李雅普诺夫方程);ie_22=1.9261(LMI工具箱);
(2)EP增益
求解该系统的EP增益代码实现为:
clc;clear all
%% 定义系统状态空间矩阵
A1=[-0.2,-1;1,0];B1=[1;0];C1=[0,1]; % 系统1
A2=[-6,-9;1,0];B2=[1;0];C2=[0,20]; % 系统2
%% 求解EP增益方法1:求解李雅普诺夫方程 AX+XA'+BB'=0
% matlab X=lyap(A,C)可用于求解李雅普诺夫方程AX+XA(^T)=-C
BB1=B1*B1';AA1=A1;
X1=lyap(AA1,BB1);
ep_1=sqrt((norm(C1*X1*C1',2)))
BB2=B2*B2';AA2=A2;
X2=lyap(AA2,BB2);
ep_2=sqrt((norm(C2*X2*C2',2)))
%% 求解EP增益方法2:转化为优化问题用LMI工具箱求解
setlmis([]); % 定义一个矩阵不等式系统
X=lmivar(1,[2 1]); % Q矩阵是一个2*2的对称矩阵
r=lmivar(1,[1 1]); % gamma
lmiterm([1 1 1 X],A1,1,'s'); % AQ+QA'
lmiterm([1 1 1 0],BB1); % BB'
lmiterm([2 1 1 X],C1,C1'); % CQC'
lmiterm([-2 1 1 r],1,1); % -gammaI
lmiterm([-3 1 1 X],1,1); % Q
lmisys=getlmis;
n = decnbr(lmisys); % 系统决策变量个数
c = zeros(n,1); % 确定向量c的维数
for j=1:n
[rj]=defcx(lmisys,j,r);
c(j)=trace(rj);
end
[copt, xopt]=mincx(lmisys,c);
xx=dec2mat(lmisys,xopt,X);
gamma=dec2mat(lmisys,xopt,r);
ep_11=sqrt(gamma)
%2
setlmis([]); % 定义一个矩阵不等式系统
X=lmivar(1,[2 1]); % Q矩阵是一个2*2的对称矩阵
r=lmivar(1,[1 1]); % gamma
lmiterm([1 1 1 X],A2,1,'s'); % AQ+QA'
lmiterm([1 1 1 0],BB2); % BB'
lmiterm([2 1 1 X],C2,C2'); % CQC'
lmiterm([-2 1 1 r],1,1); % -gammaI
lmiterm([-3 1 1 X],1,1); % Q
lmisys=getlmis;
n = decnbr(lmisys); % 系统决策变量个数
c = zeros(n,1); % 确定向量c的维数
for j=1:n
[rj]=defcx(lmisys,j,r);
c(j)=trace(rj);
end
[copt, xopt]=mincx(lmisys,c);
xx=dec2mat(lmisys,xopt,X);
gamma=dec2mat(lmisys,xopt,r);
ep_22=sqrt(gamma)
输出为:
系统1的ep_1=1.5811(李雅普诺夫方程);ep_11=1.5822(LMI工具箱);
系统2的ep_2=1.9245(李雅普诺夫方程);ep_22=1.9260(LMI工具箱);
(3)EE增益
求解该系统的EE增益代码实现为:
clc;clear all
%% 定义系统状态空间矩阵
A1=[-0.2,-1;1,0];B1=[1;0];C1=[0,1];D1=0; % 系统1
A2=[-6,-9;1,0];B2=[1;0];C2=[0,20];D2=0; % 系统2
%% 求解EE增益方法:转化为优化问题用LMI工具箱求解
setlmis([]); % 定义一个矩阵不等式系统
X=lmivar(1,[2 1]); % P矩阵是一个2*2的对称矩阵
r=lmivar(1,[1 1]); % gamma
lmiterm([1 1 1 X],A1',1,'s'); % A'P+PA
lmiterm([1 2 1 X],B1',1); % B'P
lmiterm([1 2 2 r],-1,1); % -gammaI
lmiterm([1 3 1 0],C1); % C
lmiterm([1 3 2 0],D1); % D
lmiterm([1 3 3 r],-1,1); % -gammaI
lmiterm([-2 1 1 X],1,1); % P
lmisys=getlmis;
n = decnbr(lmisys); % 系统决策变量个数
c = zeros(n,1); % 确定向量c的维数
for j=1:n
[rj]=defcx(lmisys,j,r);
c(j)=trace(rj);
end
[copt, xopt]=mincx(lmisys,c);
xx=dec2mat(lmisys,xopt,X);
gamma=dec2mat(lmisys,xopt,r);
ee_11=gamma
%2
setlmis([]); % 定义一个矩阵不等式系统
X=lmivar(1,[2 1]); % P矩阵是一个2*2的对称矩阵
r=lmivar(1,[1 1]); % gamma
lmiterm([1 1 1 X],A2',1,'s'); % A'P+PA
lmiterm([1 2 1 X],B2',1); % B'P
lmiterm([1 2 2 r],-1,1); % -gammaI
lmiterm([1 3 1 0],C2); % C
lmiterm([1 3 2 0],D2); % D
lmiterm([1 3 3 r],-1,1); % -gammaI
lmiterm([-2 1 1 X],1,1); % P
lmisys=getlmis;
n = decnbr(lmisys); % 系统决策变量个数
c = zeros(n,1); % 确定向量c的维数
for j=1:n
[rj]=defcx(lmisys,j,r);
c(j)=trace(rj);
end
[copt, xopt]=mincx(lmisys,c);
xx=dec2mat(lmisys,xopt,X);
gamma=dec2mat(lmisys,xopt,r);
ee_22=gamma
输出为:
系统1的ee_11=5.035388(LMI工具箱);
系统2的ee_22=2.222583(LMI工具箱);
本文是作者在日常学习生活中所作,难免有遗漏或错误,遇到问题的读者请点击给我写信向我的邮箱反馈,不胜感激。