龙格库塔公式的推导.pdf

《龙格库塔公式的推导.pdf》由会员上传分享，免费在线阅读，更多相关内容在行业资料-天天文库。

1、龙格库塔公式的推导为了避免计算各阶导数和偏导数，将式(3.5.12)写成ry(th)y(t)hbk(3.5.13)iii1其中r称为阶数，bi待定系数，k由下式i决定i1kiftch(i,yt()hakjj)(3.5.14)j1i1,2,3,rc011下面针对r的取值进行讨论：r=1，此时，有c110,kfty(,)因此有：yth()yt()bhfty1(,)（3－5－15）取b11即得一阶RK公式，它就是欧拉公式，因此可以说欧拉公式是RK公式的特例。2r＝2时此时有：k1f(t,y)（3－5－16）

2、k2f(tc2h,y(t)a1k1h)将f(tc2h,y(t)a1k1h)在点(t,y)展成泰勒级数ftchyt(2,()akh11)（3－5－17）fty(,)chfty(,)akhfty(,)211ty3将(3-5-17)式代入(3-5-16)式，再将式(3.5.16)代入式(3.5.13)，得yth()yt()hbkbk1122fty(,)fty(,)yt()hbfty(,)bfty(,)cbhabkh1222121ty22ffyt()(bbhfty)(,)bcha

3、bhfty(,)122212ty(3-5-18)4将式(3.5.18)与式(3.5.12)逐项比较，按照对应项系数相等比较可得：bb112(3-5-19)b2c21/2ab1/2125式(3.5.19)是一个不定方程，它有无穷多个解。取a11/2,b10,b21,c21/2，可得yn1ynhk2k1f(tn,yn)(3-5-20)hhkf(t,yk)2nn1226取a11,b1b21/2,c21，可得：hyy(kk)n1n122kf(t,y)(3-5-21)1nnk

4、f(th,yhk)2nn1显然式(3.5.21)正好是前面介绍的改进的欧拉公式。7r＝3时按照前面的推导方法可以得到常用的三阶RK公式：hyy(k3k)n1n134kf(t,y)1nn(3-5-22)hhkf(t,yk)2nn1332h2hk3f(tn,ynk2)338r=4时同样可以得到常用的四阶RK公式：hyy(k2k2kk)n1n12346kf(t,y)1nnhhk2f(tn,ynk1)(3-5-23)22hhkf(t,yk)3nn222kf(th,

5、yhk)4nn39最优步长控制策略由此作出判断：1、若en0，则本步积分成功，现确定下一步的最大步长hn1。假定hn1足够小，则(tnhn1)(tn)，下一步误差为：kk(t)h()thn＋1n1nn1en1yy11nn110e为使n10，即：(t)hknn10y1n则有1/k0(yn1)hn1(t)n11将式(3.5.37)代入上式得:1/k1/kkh(3-5-38)0n0hhn1eennn122、若en0，

6、则本步失败，按式(3.5.38)求出一个积分步长，它表示重新积分的本步步长，再算一遍，即：1/k0hhnn(3-5-39)en13由于假定了hn1足够小，因此(tn)基本不变，故必须限制步长的缩小与放大，一般限制的最大放缩系数为10，即要求：0.1hh10hnn1n14有关最优步长的控制，除此方法之外，还有吉尔(Gear)法等。采用最优步长控制后、f计算量有明显减少，但上述两种控制方法对于函数中含有间断特性的情况不适合。因为在间断点附近会出现步长频繁放大、缩小的振荡现象，由于最优步长控制法是以本步误差外推下一步步长，因此振

7、荡现象更为严重。15根匹配法关于最后一步，附加零点的说明：因为一般的传递函数模型有：mn因此在s平面上有n-m个零点在负无穷远处。不妨假设均在处，由此可见，T在z平面上尚有n-m个零点在e0处，即尚有n-m个零点在z平面的原点。16一个例子sG(s)设s22s1，试利用根匹配法，求与之相匹配的离散化模型，步骤如下：（1）G(s)s，求得2(s1)p1,p1,q0,n2,m1121（2）TTpe,pe,q11217K(z1)z（3）G(z)T2(ze)（4）对斜坡函数u(t)t,G(s)具有

8、非零且有限的稳态值：s1y()limsG(