经典的连续系统仿真建模方法学

《经典的连续系统仿真建模方法学》由会员分享，可在线阅读，更多相关《经典的连续系统仿真建模方法学（16页珍藏版）》请在金锄头文库上搜索。

1、32 第二章 经典的连续系统仿真建模方法学本章讨论经典的连续系统数字仿真的原理与方法，内容包括连续系统数字仿真的基本 概念、经典的数值积分法、经典的线性多步法等。在数字计算机上进行连续系统仿真，首 先要将连续模型离散化，因此，2.1 节首先讨论离散化原理及要求，这是连续系统仿真的基 础。然后，2.2 节对经典的数值积分法-龙格- 库塔法及其它典型的数值积分法仿真建模原 理进行详细分析，并通过实例说明其应用要点；而 2.3 节对经典的线性多步法进行了介 绍 2.1 离散化原理及要求 在数字计算机上对连续系统进行仿真时，首先遇到的问题是如何解决数字计算机在数 值及时间上的离散性与被仿真系统数值及时

2、间上的连续性这一基本问题。 从根本意义上讲，数字计算机所进行的数值计算仅仅是“数字”计算，它表示数值的 精度受限于字长，这将引入舍入误差；另一方面，这种计算是按指令一步一步进行的，因 而，还必须将时间离散化，这样就只能得到离散时间点上系统性能。用数字仿真的方法对 微分方程的数值积分是通过某种数值计算方法来实现的。任何一种计算方法都只能是原积 分的一种近似。因此，连续系统仿真，从本质上是对原连续系统从时间、数值两个方面对 原系统进行离散化，并选择合适的数值计算方法来近似积分运算，由此得到的离散模型来 近似原连续模型。如何保证离散模型的计算结果从原理上确能代表原系统的行为，这是连 续系统数字仿真首

3、先必须解决的问题。 设系统模型为： ，其中 u(t)为输入变量，y(t) 为系统变量；令仿真时间 ) , , ( t u y f y & 间隔为 h，离散化后的输入变量为 ，系统变量为 ，其中 表示 t=kh 。如果 ) ( k t u ) ( k t y k t ， ，即 ， ) ( ) ( k k t u t u ) ( ) ( k k t y t y 0 ) ( ) ( ) ( k k k u t u t u t e （对所有 k=0,1,2, ） ，则可认为两模型等价，这称为相似原 0 ) ( ) ( ) ( k k k y t y t y t e 理（参见图 2.1） 。 实际上，要

4、完全保证 是很困难的。进一步分析离散化引入的误差， 0 ) ( , 0 ) ( k y k u t e t e 随着计算机技术的发展，由计算机字长引入的舍入误差可以忽略，关键是数值积分算法， 也称为仿真建模方法。相似原理用于仿真时，对仿真建模方法有三个基本要求： （1）稳定性：若原连续系统是稳定的，则离散化后得到的仿真模型也应是稳定的。关于稳 定性的详细讨论将在 2.4 节中进行。 原连续模型 ) , , ( t u y f y & 仿真模型 ) , , ( k t u y f y & u(t) h y(t) ) ( k t y) ( k t u - + 0 ) ( k y t e 图 2.1

5、 相 似 原 理33 f(t,y) f(t 0 ,y o ) t t t 0 t 1 图 2.2 数值积分法原理 （2）准确性：有不同的准确性评价准则，最基本的准则是：绝对误差准则： ) ( ) ( ) ( k k k y t y t y t e相对误差准则： ) ( ) ( ) ( ) ( k k k k y t y t y t y t e 其中 规定精度的误差量。 （3）快速性：如前所述，数字仿真是一步一步推进的，即由某一初始值 出发，逐步 ) ( 0 t y 计算，得到 ，每一步计算所需时间决定了仿真速度。若第 k ) ( , ), ( ), ( 2 1 k t y t y t y L

6、步计算对应的系统时间间隔为 计算机由 计算 需要的时 , 1 k k k t t h ) ( k t y ) ( 1 k t y 间为 T k ，则，若 T k =h k 称为实时仿真，T k h k 称为超实时仿真，而大多数情况是 T k h k，对应于离线仿真。 连续系统数字仿真中离散化最基本的算法是数值积分算法。对于形如 的系统，已知系统变量 的初始条件 ，现在要求 随时间变化的 ) , , ( t u y f y & y y t y ( ) 0 0 y 过程 计算过程可以这样考虑（参见图 2.2）：首先求出初始点 的 y t ( ) 。 y t y ( ) 0 0 ，微分方程可以写作：

7、 f t y ( ) 0 0 ，（2.1） y t y f t y dt t t ( ) ( , ) 0 0 图 2.2 所示曲线下的面积就是 ，由于 y t ( ) 难以得到 f(y,u,t)积分的数值表达式，人们对数 值积分方法进行了长期探索，其中欧拉法是最 经典的近似方法。欧拉法用矩形面积近似表示积分结果，也 就是当 t=t 1 时， 的近似值为 ： y t ( ) 1 y 1(2.2) y y t y t f t y 1 1 0 0 0 ( ) ( ) ， 重复上述作法，当 时 t t 2y y t y t t f t y 2 2 1 2 1 1 1 ( ) ( ) ( ) ， 所以，

8、对任意时刻 t k+1 ，有：（2.3） y y t y t t f t y 1 1 1 ( ) ( ) ( ) ， 令 称为第 步的计算步距。若积分过程中步距不变 ，可以证明，欧拉 t t h 1 h h 法的截断误差正比于 。 h 2 为进一步提高计算精度，人们提出了“梯形法” 。 梯形法近似积分形式如式(2.4)所示， 令： 已知： 时 的近似值 ,那么： t t h h 1 t t y t ( ) y (2.4) y y t y h f t y f t y 1 1 1 1 1 2 ( ) ( , ) ( , ) 可见，梯形法是隐函数形式。采用这种积分方法最简单的预报校正方法是用欧拉法估

9、计 初值，用梯形法校正，即：(2.5) y y h f t y f t y i i 1 1 1 1 1 2 ( ) ( ) ( , ) ( , )34(2.6) y y h f t y i 1 ( ) ( , )式(2.6)称作预报公式，采用欧拉法，式(2.5)为校正公式，采用梯形法。用欧拉法估计一 次 的值，代入校正公式得到 的校正值 。设是规定的足够小正整数，称作允 y i 1 ( ) y 1 y i 1 1 ( ) 许误差，若 i=0, i+1=1 称作第一次校正； 称作第二次校正；通过反复迭代， i i 1 1 2 ， 直到满足 ，这时 是满足误差要求的校正值。 y y i i 1 1

10、 1 y i 1 1 ( )上述方法是针对(2.3) 式所示的微分方程在已知初值情况下进行求解，因此也称为微分方 程初值问题数值计算法，为统一起见，本书中称为数值积分法。 连续系统数字仿真的离散化方法有两类，它们是数值积分方法和离散相似方法，本章 讨论数值积分法。 数值积分方法采用递推方式进行运算，而采用不同的积分方法会引进不同的计算误差， 为了提高计算精度，往往会增加运算量。就同一种积分算法而言，为提高计算精度，减小 积分步距 ,计算量增大，影响系统运算速度。因此，计算精度和速度是连续系统仿真中常 h 迂到的一对矛盾，也是数字仿真中要求解决的问题之一。也就是说，选择合适的算法、合 适的软、硬

11、件环境，在保证计算精度的前提下，考虑怎样提高仿真的速度。 经典的数值积分法可分为两类：单步法与多步法，下面我们将分别来介绍这两类方法 中的最常用的算法。 2.2 龙格库塔法 2.2.1 龙格- 库塔法基本原理在上一节中我们已经讲过，在连续系统的仿真中，主要的数值计算工作是对 的一阶微分方程进行求解。因为 dy dt f t y ( , )y t y t f t y dt t t ( ) ( ) ( , ) 1 1 若令： , y y t ( ) Q f t y dt t t k ( , ) 1 则有 (2.7) y t y y ( ) 1 1 Q 因此主要问题就是如何对 进行数值求解，即如何对

12、 进行积分。通常称作 Q f t y ( , ) “右端函数”计算问题。已知 ，假设我们从 跨出一步， ， 时刻 y t y ( ) 0 0 t 0 t t h 1 0 t 1 为 ，可以在 附近展开台劳级数，只保留 项，则有： y y t h 1 0 ( ) t 0 h 2(2.8) 2 0 0 0 1 0 ) ( 2 1 ) , ( h t f dt dy y f h y t f y y t 在式(2.8)中括号后的下标 0 表示括号中的函数将用 均同。我们假设这 t t y y 0 0 ， 个解可以写成如下形式：(2.9) y y a k a k h 1 0 1 1 2 2 ( ) 其中

13、， k f t y 1 0 0 ( , )k f t b h y b k h 2 0 1 0 2 1 ( ) ，35 对 式右端的函数在 处展成台劳级数，保留 项，可得： k 2 t t y y 0 0 ， hh k b t f b y t f k t 0 ) ( ) , ( 1 2 1 0 0 2 将 代入(2.9) 式，则有： k k 1 2 ， ) ( ) , ( ) , ( 0 1 2 1 0 0 2 0 0 1 0 1 h y f k b t f b y t f h a y t hf a y y t 将上式与(2.8)式进行比较，可得： 可见有四个未知数 但只有三个方程，因此有无穷多

14、个解，若限定 a a b b 1 2 1 2 ， ， ， ， ，则可得其中一个解： a a 1 2 a a b b 1 2 1 2 1 2 1 ， 将它们代入(2.9) 式可得一组计算公式：(2.10) y y h 1 0 1 2 2 ( ) 其中， 1 0 0 2 0 0 1 f t y f t h y h ( , ) ( ) ， ， 若写成一般递推形式，即为：(2.11) y t y y h k k k ( ) ( ) 1 1 1 2 2 其中， 1 2 1 f t y f t h y h k k k k ( , ) ( , ) ，由于(2.8)式只取了 两项，而将 以上的高阶项略去了，所以这种递推公式的截 h h ， 2 h 3 断误差正比于 的三次方，又由于计算时只取了 及 项，故这种方法被称为二阶龙格 h h h 2 库塔法（简称 RK-2） 。根据上述原理，若在展成台劳级数时保留 及 项，那么可得一套截断误差 h h h ， ， 2 3 h 4 正比于 的四阶龙格- 库塔法（简称 RK-4）公式： h 5(2.12) y t y y k k k ( ) ( )