2011车辆动力学(7)- 数值求解

第四章微分方程数值积分算法一、高阶常微分方程转化为一阶常微分方程组二、数值积分算法的类型三、数值积分算法选用原则1.求系统最高阶的导数，即将动态方程写成如下形式：()(1)(,,,,)nnxftxxx2.作如下变换12(,,,,)nnxftxxx(1)123,,,,nnxxxxxxxx将动态方程转换为一阶常微分方程形式：一、高阶常微分方程转化为一阶常微分方程组3.写出如下每一个变量的一阶常微分方程：12231,,,nnxxxxxx例1将高阶常微分方程转化为一阶常微分方程组。()()()()tctktMt解：12,xx一、高阶常微分方程转化为一阶常微分方程组()()()()tctktMt12221()xxxcxkxMt112212{}{}TTXX122M1M1J2Jck11112122212122()()()()JMckJckM例2一、高阶常微分方程转化为一阶常微分方程组111112222200JMcckkJMcckk二、数值积分算法的类型显式算法：计算yn+1时，没有用到tn时刻以后的状态或输入。隐式算法:计算yn+1时，要用tn时刻以后的状态或输入。单步法：计算yn+1时只与tn时刻的状态或输入有关。多步法：计算yn+1时要用到tn时刻以前的状态或输入。1、欧拉(Euler)算法计算公式：1,,1,0),,(1Mnythfyynnnn()yt()yt1nyny1ntntht单步显式一阶算法前向1,,1,0),,(111Mnythfyynnnn后向()yt()ytny1ny1ntntht单步隐式一阶算法112121()2(,)(,)nnnnnnhyykkkftykfthyhk2、改进的欧拉公式（预测-校正法）单步显式二阶算法间接利用泰勒展开式。用在若干个点上函数值f(t,y)的线性组合来代替高阶导数项，既可以避免计算高阶导数，又可以提高数值计算精度。3、龙格-库塔法（Runge-Kutta）1）算法思路2）经典显式四阶Runge-Kutta法(RK4)),()2,2()2,2(),()22(6342312143211hkyhtfkkhyhtfkkhyhtfkytfkkkkkhyynnnnnnnnnn单步显式四阶算法例：设系统方程为2()()0,ytyt1)0(y试用显式四阶Runge—Kutta法求其数值解0.1,01)ht（取仿真步距解：原方程为2()(),ytyt)(),(2tyytf递推公式可写为1123421221123222433(22)6(,)(,)()222(,)()222(,)()nnnnnnnnnnnnnnhyykkkkkftyyhhhkftykykhhhkftykykkfthyhkyhk3、龙格-库塔法（Runge-Kutta）对于线性定常系统)()()(tttBuAxx)()()(tttBuAxfRK4的4个k向量可表示为)(),()2(2)2,2()2(2)2,2()(),(3342231121hthhxhththhxhththhxhttxtmmmmmmmmmmmmmmmmBukxAkfkBukxAkfkBukxAkfkBuAxfk)22(643211kkkkxxhmm3）矩阵分析法（RK4解状态方程）3、龙格-库塔法（Runge-Kutta）1）高阶隐式积分算法三阶隐式积分算法：1,,1,0)],(5),(8),([1211111Mnytfytfytfhyynnnnnnnn4、线性多步法多步隐式算法2）亚当姆斯（Adams）积分算法•二阶亚当姆斯积分算法1,,1,0)],,(3),([2111Mnytfytfhyynnnnnn多步显式算法•计算步长h；•算法阶数；•仿真对象的动态特性（主要是时间常数）。解决好如下矛盾：•减小步长，可以减小截断误差；•减小步长，会导致计算步数增加，累积次数增多，舍入误差增大。合理步长步长舍入误差综合误差局部截断误差误差O1、积分误差三、数值积分算法选用原则2、数值解法稳定性1）稳定性概念：如果误差积累能够得到抑制，不随计算时间的增加而无限增大，则认为相应的积分算法是数值稳定的。例：()()ytyt0，的数值解法稳定性讨论。前向欧拉公式110(1)nnnnyyhyhy要保证稳定，必须：11h即2h要求：积分步长必须小于系统时间常数的2倍。前向欧拉法稳定区域后向欧拉公式1(1)10(1)(1)nnnyhyhy11h要保证稳定，必须：后向欧拉法是恒稳定的。•计算步长h；•仿真对象的动态特性（主要是时间常数）。2）数值稳定性的主要影响因素2、数值解法稳定性minmin51201ThT40nthminT为系统的最小时间常数积分步长经验公式为系统在阶跃函数作用下的过渡过程时间nt•当状态变化缓慢时，加大步长以减少计算时间；•当状态变化很快时，减小步长以提高积分精度。3）变步长积分算法2、数值解法稳定性•分段变步长:将过渡过程分成几段，每段使用不同的步长；•根据每步积分的误差，自动调整下一步的积分步长；•最优步长法，即使每一步积分步长在保证精度的前提下取最大的步长。变步长积分算法类型：不适合于实时仿真。如果一个系统中包括快变和慢变行为，系统模型的时间常数反映了这种特性。一个系统的刚度可以通过系统的最长时间常数与最短时间常数之比来反映。当此比率超过三个数量级时，该系统称为刚性系统，否则称非刚性系统。1x01s+1object10.001s+1controlleryToWorkspace4）刚性系统仿真2、数值解法稳定性非刚性系统时，精度决定于积分步长。刚性系统，数值解的稳定性和精度决定于积分算法和积分步长。1）初始条件•每一个一阶微分方程和差分方程都有与之对应的初始条件；•所有初始条件必须在求解开始前设定。3、数值求解的其他问题2）驱动信号•驱动系统运行的输入信号；•信号的选择依赖于目标。例如，油门开度、方向盘、刹车输入。3）终止条件•设定的时间到；•用户的“终止”操作；•设定的目标达到。