华中科技大学电气学院matlab选修课大作业PI控制器

13/13华中科技大学电气学院matlab选修课大作业PI控制器华中科技大学电气学院matlab选修课大作业PI控制器

2018级

《MATAB语言与操纵系统仿真》课程

大作业

姓名赖智鹏

学号u201811806

所在院系电气与电子工程学院

班级电气0809

日期1月16日

同组人员

作业评分评阅人

设计报告评分表

ＰＩ操纵器设计与讨论

1.引言

本文讨论的对象是智能交通/高速公路系统IVHS(习题7.3)，系统通过电子技术为驾驶者和操纵系统提供实时路面情形，该系统还提供更方便服务，驾驶者可不用自己操纵驾驶，而系统自动操纵车辆，且保持车相对速度，以实现交通的有序进行。

本文前半部分讨论了PI操纵器参数的设计咨询题，依照性能指标推导操纵器参数的约束条件，通过相关算法找到咨询题的解，后半部分讨论系统参数变化时对系统性能的阻碍，通过比较得出综合性能较优的操纵器参数。

2.给定系统的操纵器设计

2.1性能分析：

（1）阶跃响应零稳态误差。由原系统为零型系统，故必须通过操纵器在原点加入至少一个开环极点，可选择的由PI操纵器、PID操纵器等，而单独使用比例操纵器已达不到目标。（2）要求系统对爬坡响应的跟踪误差小于25%，故Kv>4，从而有：

001

4

16ssKvsGcGsGc->->==>∑∑。（3）阶跃响应的超调量小于5%，因此假如校正后系统近似为2阶系统，要求阻尼比ξ>0.7。

（4）调剂时刻Ts=4/nξω4/1.5=2.67。

就目前分析，单独使用比例操纵器不能满足要求，而在考虑选择PID操纵器之前我们选择相对简单的PI操纵器：

i

PK

GcKS=+。

图1原系统

做出根轨迹图并作出满足性能指标的区域：

rlocus(G);

holdon;

plot([-2.66-2.66],[-2020]);％指定性能指标在根轨迹图中所在区域zeta=0.7;

plot([0-20*zeta],[020*sqrt(1-zeta^2)],[0-20*zeta],[0-20*sqrt(1-zeta^2)]);％指定性能指标在根轨迹图中所在区域。

-20-15

-10

-5

510

15

20

RootLocus

RealAxis

ImaginaryAxis

图2原系统根轨迹图（黑线表示所要达到的指标边界）

系统波特图及其近似画法：figure;

sys2=tf([1/8],[10]);%第二段直线近似sys3=tf([1],[100]);%第三段直线近似bodemag(sys1,{0.1,2})holdon;

bodemag(sys2,{2,8})bodemag(sys3,{8,100})bode(num,den,{0.1,100})grid;

-100-80

-60

-40

-20

Magnitude(dB)10

10

10

10

-180

-135-90-450Phase(deg)

BodeDiagram

Frequency(rad/sec)

图3原系统波特图

由波特图可看出原系统开环波特式增益小于1，低频区水平，转折频率为2rad/s和8rad/s，是

一个稳固的二阶系统。

下图为其开环阶跃响应和斜坡响应：

subplot(2,1,1);

step(feedback(G,1))

subplot(2,1,2);

t=[0:0.01:5];

lsim(feedback(G,1),t,t)

00.511.522.533.544.55

2

4

6

System:untitled1

Time(sec):3.79

Amplitude:0.183

Input:In(1)

Time(sec):3.79

Amplitude:3.79

Input:In(1)

Time(sec):4.6

Amplitude:4.6

System:untitled1

Time(sec):4.6

Amplitude:0.228

00.20.40.60.811.21.41.61.8

0.02

0.04

0.06

StepResponse

Time(sec)

A

m

p

l

i

t

u

d

e

LinearSimulationResults

Time(sec)

A

m

p

l

i

t

u

d

e

图4原系统响应

阶跃响应稳态误差接近1，斜坡响应不能跟踪输入信号变化，可通过PI操纵器增加开环系统类型以排除这两个缺陷。

图5校正后系统框图

2.2．参数KP、KI约束条件：

开环传函：

(2)(8)

KpsKi

GH

sss

?+

=

++

系统存在3个确定的极点在原点（两个自带的极点，一个PI操纵器引入的极点），一个变化的零点，增益也是可变的，由性能指标要求即可确定主导极点所在区域。

（1）根轨迹渐近线与实轴的交点小于-2.66，即

280(/)5/22.662

pzKpKiKpKimnα+--=

==-+4，即

1

416sKvsGcGKi->==>∑

得出Ki64>

（3）特点方程为32

10(16)

0ssKpsKi++++=闭环传函

3210(16)KpsKi

ssKpsKi?+++++

由劳斯判据，有

3210

116010

0160100010

0sKpsKi

KpKi

ssKi

++-

得

160100KpKi+->0Ki>

综上所述，关于Kp，Ki有4个约束条件（实为3个），做出取值区间如下图：

symskpki

f1='ki-64+0*kp';f2='kp-ki/10+16';f3='ki-14/3*kp';figure

ezplot(f1,[0,160,-20,100]);holdon;ezplot(f2,[0,160,-20,100]);holdon;ezplot(f3,[0,160,-20,100]);holdon;Ki和Kp取值区域：

204060

80100120140160

-200

20

40

60

80

100

ki

kp

ki-14/3kp=0

图6Kp，Ki取值区域

图中箭头方向为解区域方向。

2.3考察边界值KI=64,KP=1

3.8

%通过linmod()函数将方框图转换成状态空间形式，再有状态空间得到传递函数[a,b,c,d]=linmod('untitled');sys=tf(ss(a,b,c,d))得到传函：

13.8s+64

s^3+4s^2+17.8s+64

BodeDiagramGm=Inf,Pm=Inf

Frequency(rad/sec)

-150-100-50

050

100

Magnitude(dB)

10

10

10

10

10

Phase(deg)

图7系统各环节的比较和对性能的阻碍

由图7波特图知加入PI校正器增加了系统的类型。

-25-20

-15-10-5

0510

152025

RootLocus

RealAxis

ImaginaryAxis

图8校正后根轨迹

由图8根轨迹看出，现在主导极点不能同时满足阻尼比和阻尼频率的条件。

阶跃和斜坡响应：

Time(sec)

Amplitude

00.511.5

22.5

33.5

0.2

0.4

0.6

0.8

1

1.2

1.4

图9校正后阶跃响应（Ki=64,Kp=13.8）

由图7知超调量20.2%大于5%，稳固时刻为2.68s也大于1.5s，故性能不满足要求。

0.511.522.53

3.54

4.55Time(sec)

Amplitude

图10校正后斜坡响应（Ki=64,Kp=13.8）

由上图8知现在已达到斜坡响应的稳态误差要求要求由代码得到性能指标数值PO、Ts和斜坡响应稳态误差Kp=13.8;Ki=64;s=tf('s');PI=Kp+Ki/s;

G_PI=series(G,PI);

G_PI_closed=feedback(G_PI,1);figure

step(G_PI_closed);

[y_s,t_s]=step(G_PI_closed);yss=dcgain(Gc);

po=(max(y_s)-yss)*100%超调量POm=find(abs(y_s-yss)>0.02);ts=t_s(length(m))%调剂时刻t_r=[0:0.001:5];u=t_r;%斜坡信号figure;

lsim(G_PI_closed,u,t_r);%斜坡响应y_r=lsim(G_PI_closed,u,t_r);

err=u(length(u))-y_r(length(y_r))由代码得po=20.2016ts=2.3831err=0.2500而由图得出的po=20.2Ts=2.68err=0.25

误差分析：二者并不完全一致，其中Ts差不较大0.3/2=11%左右。注意到m数组如下，其中25——28缺省2个数，46——52缺省了6个数，因此将导致此方法运算的Ts有差异：

ts=t_s(length(m)+3+6)>>ts=2.6852

现在就与图中值吻合了

关于一样情形，可用以下代码实现：m=find(abs(y_s-1)>0.02);

dm=find(abs(y_s(1:(length(m)-1))-1)5%，不满足性能要求，因此刚才界定的区域还不足以完全确定满足要求的解。下图得出现在闭环系统的零极点图：

RealAxis

ImaginaryAxis

图12KP=40，KI=100时，零极点图

由图12知系统主导极点在规定区域，但零点与极点的位置专门接近，而零点有增大超调量的作用，下面为去掉该零点之后系统的响应：

Time(sec)

Amplitude

LinearSimulationResults

Time(sec)

Amplitude

图13系统去除零点后响应

由图13知现在超调量为0%4），超调

量要求不大于10%（ξ>0.6），仿上述设计流程，得出：

Ki/Kp2

160100KpKi+->

做出该区间：symskpki

f1='ki-64+0*kp';f2='kp-ki/10+16';f3='ki-2*kp';holdon

ezplot(f1,[0,160,-20,100]);holdon;ezplot(f2,[0,160,-20,100]);holdon;ezplot(f3,[0,160,-20,100]);holdon;

0204060

80100120140160

ki

kp

ki-2kp=0

图21

同样的方法能够得出一系列满足条件离散点，如下图

020406080100120

140160180200ki

kpfigure

ki-2kpfigure=0

图22

性能指标见下表：

E

rr0.25

0.249

985

0.249

923

0.231

884

0.231

873

0.231

83

0.216

216

0.216

208

0.216

177

0.202

532

0.202

525

K

I

79848484848489898994K

P

49.5424752576749.554.559.552

PO6.202

628

9.218

092

7.833

353

7.851

369

8.035

544

8.782

046

8.934

857

8.672

828

9.273

72

9.935

864

Ts0.950

096

0.918

025

0.773

279

0.817

075

0.924

344

0.996

271

0.762

198

0.918

862

0.915

286

0.738

163

Err0.202

503

0.190

476

0.190

471

0.190

454

0.190

415

0.190

231

0.179

771

0.179

758

0.179

728

0.170

21

2．7系统参数变化时，系统性能的研究

图23参数ab可变系统

关于系统

1

()()

G

sasb

=

++，其参数a,b可能发生改变，现在我们假定PI操纵器参数是恒定

的，给定Kp=33，Ki=69。其位于KI-KP坐标系中位置如下图，由图可知Kp=33，Ki=69时（Kp，Ki）（红点标记）与a=2,b=8时满足性能指标的解空间的相对位置，即处于解空间的临界区。

0204060

80100120140160

-20

20

40

60

80

100

ki

kp

ki-14/3kp=0

图24给定PI参数位置

图25

校验现在系统性能，下图为阶跃响应和斜坡响应，由图可知阶跃响应超调量略大于5%，其余指标满足要求。

Time(sec)

Amplitude

0.2

0.4

0.6

0.8

1

1.2

1.4

图26阶跃响应a=2b=8

LinearSimulationResults

Time(sec)

Amplitude

0.511.522.53

3.54

4.55

图27斜坡响应

类似于前面关于Kp、Ki取值范畴的分析，分析KPKI不变时，为满足性能指标要求的a、b取值区间。

开环传函

,33,69

()()KpsKi

GHKpKissasb?+===++

系统存在1个确定的零点（-Ki/Kp）和在原点的极点（0），两个变化的极点（-a和-b），由性能指标要求可确定主导极点所在区域，且

0(/)0.24()/22.662pzabKpKiabmnα+--=

==-+

由Kv>4，得出ab

0ab+>

综上所述，关于a,b有四个约束条件，做出a,b取值区间如下图：figure;symsabkp=33;ki=69;

f1='a^2*b+a*b^2+46*(a+b)-69';f2='a+b-5.8';f3='a*b-17.25'

ezplot(f1,[-100,100,-100,100]);holdon;ezplot(f2,[-100,100,-100,100]);holdon;ezplot(f3,[-100,100,-100,100]);holdon;

a

b

图29约束条件界定区域

a

b

ab-17.25=0

图30放大后

考虑到a,b之间是等价的，取区间的一半，如下图，箭头方向标注满足要求区域：

b

a

图31放大后再次采纳穷举法求a，b：

a=2;%初值

b=8;%初值

s=tf('s');%

kp=33;

ki=69;

PI=kp+ki/s;%PI

AB=zeros(100,2);%储备a,b

PO=zeros(100,1);%

Ts=zeros(100,1);%

Err=zeros(100,1);%

i=1;

whilea0.02);

dm=find(abs(y_s(1:(length(m)-1))-1)0.02);

ts=t_s(length(m))%?