第二章 经典的连续系统仿真建模方法学-8-3

第二章经典的连续系统仿真建模方法学本章教学目的及要求l掌握连续系统数值积分法的基本原理2掌握Euler法和Runge-Kutta法的概念和方法3掌握线性多步法的概念和方法4掌握数值积分法稳定性分析方法5熟悉数值积分法的选择方法第二章经典的连续系统仿真建模方法学2.1离散化原理及要求2.2Runge-Kutta积分法2.3线性多步法2.4数值积分法稳定性分析2.5数值积分法的选择和计算步距的确定2.1离散化原理及要求一离散化原理1问题求解的两种方法（1）仿真方法---通常是离散的数值解（仅能说明一些特定的状态的，如果要了解其他状态下的情况，就必须在其他给定条件下，重新进行仿真过程）（2）解析方法---得到的是更为一般的系统解的表示形式（可以表示系统所有可能的状态）（3）适用范围---对于复杂系统而言，求解析解是很复杂的，大多数情况下，求不出解析解，因此必须借助于数值解法，对连续系统进行仿真。2相似原理设系统模型为：，其中：u(t)为输入变量，y(t)为系统变量；令仿真时间间隔为h，离散化后的输入变量为系统变量为，tn=nh。若，即，（对所有n=0,1,2,…），则可认为两模型等价。),,(tuyfy)(ˆntu)(ˆnty)()(ˆnntutu)()(ˆnntyty0)()(ˆ)(nnnututute0)()(ˆ)(nnnytytyte原连续模型仿真模型u(t)hy(t)-+相似原理),,(tuyfy),ˆ,ˆ(ˆntuyfy0)(nyte)(ˆntu)(ˆnty基本要求（1）稳定性：不改变原系统的稳定性若原连续系统是稳定的，则离散化后得到的仿真模型也应是稳定的。若原连续系统是不稳定的，则离散化后得到的仿真模型也应是不稳定的。（2）准确性：有不同的准确性评价准则，最基本的准则是：绝对误差准则：相对误差准则：其中规定精度的误差量。)()(ˆ)(nnnytytyte)(ˆ)()(ˆ)(nnnnytytytyte3数值方法的基本数学描述将仿真系统表示成1阶微分方程组或状态方程的形式若1阶微分方程及其初值为：式中Y为n维状态变量，F(t,Y)n维向量函数。设方程的在tn处的形式上的连续解为：希望能够找到一个近似公式：来近似表示，其中:Yn为准确值Y(tn)的近似值；Qn是准确积分值的近似值。00)(),(YtYYtFY110),()(),()()(01nnnttttnndtYtFtYdtYtFtYtYnnQYY1),(nnttndtYtFQ4数值解法的含义（1）寻求初值问题真解在一系列离散点:t1t2···tn,···的近似解Y1,Y2,···,YN,···(即数值解)相邻两个时间离散点间隔:hn=tn+1-tn称为步长，通常hn=h（定值）。（2）数值积分法的主要问题就是对函数f(t,y)的数值积分问题（如何求出该函数定积分的近似解）。（3）连续系统数值积分方法（或称为数值解法）就是利用数值积分方法，对常微分方程（组）建立离散化的数学模型—差分方程，并求出其数值解。5主要步骤(1)将连续变量问题用数值积分方法,转化为离散的差分方程的初值问题；(2)根据已知的初始条件y0，逐步递推计算出以后各个离散点上的数值解yi(i=1,2,···)6主要分类（不同的递推算法）单步法，多步法和预估-校正法（不同的系统求解精度、速度和稳定性）二Euler法Euler法是最简单的一种数值积分方法。特点：导出简单几何意义明显，便于理解计算精度比较差1Taylor级数展开有如下标量微分方程：若y(t)=g(t,y)为其解析解，将y(t)展开成Taylor级数：从而：将上式写成差分方程：00)(),(ytyytfyhtytyhty)()()(),()()(ythftyhty),()(1nnnnythftyy2矩形近似解法对微分方程在区间上求积分，得：如果把积分间隔取得足够小，使得在[tn,tn+1]上f(t,y)可以近似地看成f(tn,yn)，就可以用矩形面积近似代替该区间上的曲线面积，有：即：将曲线f(t,y)看成是阶梯函数1),()()(1nnttnndtytftyty11),()(nnnnnyythfyty),2,1,0(1nhfyynnn3切线近似在tn的一个小邻域(tn-ε,tn+ε)内，曲线y(t)可以用tn处的切线来表示，如图，•y(t)在tn处的斜率和切线方程分别为：),(|nnttytfdtdyn))(,(nnnnttytfyy取tn+1=tn+∆t,使得点tn+1与tn非常靠近，利用上述切线方程可获得tn+1处y(t)的近似值：差分方程形式为：11),()(nnnnnyythfyty),3,2,1,0(1nhfyynnn4Euler公式使用方法：只要给定初始条件y0和步长h,就可以根据f(t0,y0)算出y1，再由y1算出y2，逐次递推，计算出y3,y4,•••。特点:(1)公式公式中任一个新的数值解yn+1都是基于前一个数值解及其导数f(tn,yn)求得的。(2)方法简单、计算量小(3)精度低),3,2,1,0(1nhfyynnn例对于线性微分方程试采用Euler法进行计算，取步长h=0.2。0)0(,1)0(,043YYYYY9904.19907.16576.16592.1376.1384.112.116.18.01,8.04.02.0)43(0)0(,431)0(,55443322111111zyzyzyzyzyyzzzyyyzhzzhzyyEulerzyzzyzynnnnnnnnnnnnn，，，，计算结果代入初值可得到法可得到：采用解：令三改进Euler法1提出目的：提高Euler公式的精度2改进的途径：使用梯形面积来代替每一个小区间的曲线面积（如图）3公式的推导（1）计算曲边梯形的面积（2）计算直边梯形的面积)()(),(111nntttytydtytfSnn)],(),([2)()()],(),([21111112nnnnnnnnnnytfytfhtytyytfytfhS（3）在步长h比较小情况下，用直边梯形面积代替曲边梯形面积，得到差分方程：（4）计算中的预估和校正原因：差分方程中右边含有待求量yn+1,公式无法起步计算。解决方法：A每次计算时先用Euler法计算出y(tn+1)的近似值ypn+1作为初值；B将ypn+1代入原微分方程，计算出函数fn+1的近似值fpn+1=f(tn+1,ypn+1)；C利用梯形公式，求出修正后的ycn+1。D重复B，C过程，直到ycn+1和yc-1n+1之间的差值在可接受范围内。)(21)],,(),([211111nnnnnnnnnnffyyytfytfhyy即4改进的Euler公式(预估-校正法)式中第一项为预估公式，第二项为校正公式特点：采用[tn,y(tn)]和[tn+1,y(tn+1)]两点斜率平均值的结果，利用了两点信息，从而提高了计算精度。)],(),([2),(1111pnnnnncnnnnpnytfytfhyyythfyy),3,2,1,0(1nhfyynnn)],(),([2),(1111pnnnnncnnnnpnytfytfhyyythfyy2.2Runge-Kutta积分法提出的原因：直接采用Taylor级数展开方法，需要计算函数y(t)在某一展开点的高阶导数，使用起来不方便。基本思想：用几个点上的y(t)的一阶导函数值的线性组合来近似代替y(t)在某一点的各阶导数，然后利用Taylor级数展开式确定线性组合中各加权系统----间接利用Taylor级数展开式目的：1避免计算高阶导数；2多取几项，提高数值积分的精度一阶微分方程：1假定y(t)是微分方程的解，将y(t)展开成Taylor级数式中：2若令：00)(),,(Ytyytfdtdy)(2)()()(2tyhtyhtyhty！yfytftftyyftfytfdtdtyytfty),()],([)(),()(*)(2),()()(),,(,2ytytfffhythftyhtyyffytfftff则：2.2.1基本原理3将y(t+h)写成线性组合形式：式中r称为阶数，bi待定系数，ki由下式决定：且定义a1=c1=0。4对y(t+h)的线性组合形式进行讨论（1）r=1为Euler公式（是一个特例）iriikbhtyhty1)()(),,2,1]()(,[11rikahtyhctfkjijjii),()()(1ythfbtyhty（2）r=2，由ki的取值公式可得到：将f[t+c2h,y(t)+a1k1h]在点(t,y)展开成Taylor级数：将k1和k2的取值代入线性组合公式，有：])(,[),(11221hkatyhctfkytfktyytfthkaytfthcytfhkatyhctf),(),(),(])(,[112112tyyfhbatfhcbfthfbbtykbhtyhtyiii2212222121),()()()()(将上式同*式逐项比较（对应系数相等），有：tyyfhbatfhcbfthfbbtykbhtyhtyiii2212222121),()()()()(21211212221bacbbb*)(2),()()(2ytfffhythftyhty得到一个不定方程组，有无穷多个解，若a1=1,b1=b2=0.5,c2=1，有：依照上述过程，可以得到r=3,r=4等情况下的对应公式。),(),()(2121211hkyhtfkytfkkkhyynnnnnn),()21,2()21,2(),()22(6342312143211hkyhtfkhkyhtfkhkyhtfkytfkkkkkhyynnnnnnnnnn)32,32()31,3(),()3(423121311hkyhtfkhkyhtfkytfkkkhyynnnnnnnn特点1单步法，可以自启动，存储量小在计算时只用到，而不直接用等项。2可变步长步长h在整个计算中并不要求固定，可以根据精度要求改变但是在一步中，为计算若干个系数ki，则必须用同一个步长h。3速度与精度四阶方法的h可以比二阶方法的h大10倍，每步计算量仅比二阶方法大一倍高于四阶的方法由于每步计算量将增加较多，而精度提高不快。4形式多样1ny21nnyy，1误差的两个准则绝对误差准则：相对误差准则：2RK法的误差估计和步长控制的基本思路每积分一步都设法估计出本步的计算误差，判断是否满足允许误差E，再选择相应的步长控制策略，调整步长，再作下一步积分运算。)()(ˆ)(kkkytytyte)(ˆ)()(ˆ)(kkkkytytytytek2.2.2误差估计与步长控制3误差估计（单步法）方法：找到另一个低阶（一般是低一阶）的RK公式，两个公式计算结果之差就可以看作是误差。常用的计算公式对：aRKM（Merson默森，1957）3-4法四阶公式与三阶公式对：两式中的k1,k3,k4相同，是一个4阶精度，3阶估计误差的公式，每次需要计算5次f，计算量较大。543111431154118926ˆ12936ˆ46KKKKhyyeKKK