仿真方法的介绍与对比

现代工程控制理论实验报告

学生姓名：任课老师：

学号：班级：

实验六：仿真方法的介绍与对比

一、实验内容及目的

实验内容：选用五种仿真方法对以下系统进行仿真，比较仿真精度和仿真时间的不同。

前四种仿真方法为：欧拉法、局部离散法、梯形法、龙格—库塔法，之后利用整体离散法作为衡量标准。

实验目的：1.分析比较各种仿真方法的特点以便将来能够选择最符合要求的仿真手段。

2. 练习使用多种仿真方法，提高仿真程序的书写能力。

二、实验内容及仿真结果

方法一：欧拉法

（1）仿真原理及程序：

○1首先将系统方框图化成积分图，求出系统的空间状态方程。

积分图如下：

1ki ⎰

2ki ⎰

⎰

⎰

⎰

⎰

⎰

1x 2x 3

x 4x 5

x 6

x x 8

x

从积分图中可以得到系统的状态空间方程。

=ｘＡｘ＋ＢＲ，其中

之后利用欧拉公式：

11(*)k k

k k x x x dt

x I dt A x BdtR ++-=

=++ 1k k x x

x dt

+-=，得1*k k x x x dt +=+。 将状态空间方程代入可得：1(*)k k x I dt A x BdtR +=++

通过循环、迭代可以得到每次的X k+1，对应的程序如下：

而系统最终的输出y=x(8)，因此记录下每次循环得到的x(8)，作为系统的输出。

（2）仿真结果：

系统的输出曲线如下：

方法二：局部离散法

（1）仿真原理及部分语句

局部离散法的思想是利用积分图将高阶系统转化为一阶系统，

之后对每一个一阶环节利用离散法求其递推关系。

该系统的一阶环节有两类，积分控制器中的积分环节和一阶惯

性环节。

对于积分环节应用离散相似法的公式(1)()()m x k x k u k φφ+=+，可

以得到递推关系式为：(1)()**()x k x k ki dt u k +=+。

对于一阶惯性环节

1k

Ts

+，通过离散相似法得到的递推关系式为：(1)()(1)()dt

dt T

T

x k e x k k e x k --+=+- 。

局部离散相似法仿真程序的主要语句如下：

（2）系统的仿真曲线：

方法三：梯形法

（1）仿真原理及部分语句:

首先观察欧拉法的主要公式(1)()()

+=+，其中的误差主

x k x k x k dt

要来源在于此时的()

x k，误差在于将上一时刻的变化率当作了这一

时刻的增长率，称此时的误差为截断误差。

如左图，s2的大小即为截断误差的大小。设E （k ）表示为从k 到k+1时刻的平均增长率。则提高仿真精度的实质，即为提高E(k)的精度，使其更接近实际值。

而梯形法的实质便是： 令()(1)()2x k x k E k ++=

，则()(1)

(1)()2

x k x k x k x k T +++=+ 而()()()x k Ax k Bu k =+，可以求得准确值。

但(1)x k +无法得到准确值，只能借助欧拉公式求其估计值：

(1)()()(1)(1)(1)

x k x k Tx k x k Ax k Bu k +=++=+++

而整理以上几条公式，可以得到直接的递推关系式：

222(1)()()()()(1)2222

T T T T

x k I TA A x k A Bu k Bu k +=++++++

与整体离散法相比较：

()

x k

(1)()()m x k x k u k φφ+=+

2221320

1

...

2!111

(...)2!3!!AT T

AT k k m e I AT T A e Bdt T T A T A T A B

k φφ==+++==+++⎰ 可以发现梯形法实际为离散相似法中的φ、m φ取前三项和前两项。 对应的主要语句如下：

（2）系统的仿真曲线：

方法四：龙格—库塔法

（1）仿真原理及部分语句：

龙格—库塔法实际叫做四阶龙格—库塔法。为了提高公式

(1)()E()x k x k k dt +=+

中E(k)的精度，这个方法通过四个不同的变化率（E ）加权平均作为最终的E(k)。

设k 时刻系统的状态变量x 的一阶导数为E1，则由状态空间方程可知：E1()()Ax k Bu k =+。

之后为了寻找12

k +时刻系统的变化率，需要用到12

k +时刻的状态变量1()2

x k +。但1()2x k +在此时是未知的，需要预测。这个方法的关键就在于1

()2x k +的求取。

既然k 时刻x 的变化率E1已知，可以通过类似欧拉法的方式

11()()122

T

x k x k E +=+，

将1

()2x k +求出。 之后借助11()2x k +可以求得12

k +时刻系统的变化率E2。

111

E1()()22

k Ax k Bu +=++

而E2已知的同时又可以通过

21()()222

T x k x k E +=+

求得另一个1

()2x k +，由此又可以得到12

k +时刻大小不同的变化率E3。

最后借助E3可以求出k+1时刻系统的状态变量x(k+1)，

(1)()*3x k x k T E +=+ 4(1)(1)E Ax k Bu k =+++

这样又可以得到k+1时刻系统变化率的估计值E4(k+1)。

这样集齐了4个变化率E1、E2、E3、E4，为了求得从k 到1

()

2

k k +时系统的平均变化率的近似值，需要将这四个变化率求平均。而E2、E3为中间时刻求到的值，更为逼近k 到1

()2

k k +时系统的平均变化率，因此权重要大一点。所以：