资源描述:
《一种对刚性常微分方程的有效算法》由会员上传分享,免费在线阅读,更多相关内容在行业资料-天天文库。
1、一种对刚性常微分方程的有效算法MarekM.Stabrowski华沙科技大学,IETiME,Koszykowa7.5,PLO0661华沙,波兰d收于1995年7月10日;1996年4月16日修订摘要Brayton-Gustavson-Hatchel法(BGH法)是向后微分公式法(BDF法)的一种,用于求解刚性常微分方程。本文介绍了BGH法的一些基础的细节,且描述了作者原创的一些对该算法的新改动。本文着重关注了减小计算量和改进误差控制两个方面。本文所举的两个例子,包括一个刚性和适度刚性的方程,体现出了BGH法相对于经典的吉尔法巨大的优越性。文中介绍的软件是可以升级的,它在PC/DOS
2、平台和两个UNIX系统下运行良好,没有出现任何问题。1997ElsevierScienceB.V.关键词:刚性微分方程,刚性微分系统一、介绍吉尔法是现在用于求解刚性常微分方程(ODES)的一种经典工具。它是一种相当稳定有效的办法(就求解速度而言)。它也是BDF法的一种。当然,其广泛使用不仅归功于在数学上的优势,也得益于一套完整fortran程序的发表。Byrne-Hindmarsh的吉尔法程序这样就成为了公共软件。但是,即使在同类的BDF方法里面,也有其他一些相当有竞争力的方法。似乎其中一些算法的潜在优势并没有被深入的研究和探索过,其中BGH法就理应得到更多的关注。其作者声称,BG
3、H法比吉尔法更稳定,尤其是在步长需要频繁变换的时候。并且,BGH法的计算量更小,误差控制也更灵活。在数学上,BGH法和吉尔法(甚至吉尔—诺德萨克法)地位并不相当。但是,吉尔—诺德萨克法按这样的方式进行一些改动之后,就和BGH法旗鼓相当了。所有这些表明,BGH法是吉尔—诺德萨克法的有力竞争者。可惜,目前发表的Fortran和C语言下的BGH法程序的特点是非常的刻板和特殊化。他们有些“特别”快速算法的特点,并没有真正库的风格。甚至,其多项式次数k被固定为两次,误差控制也被简化了,以至于对操作量的减少并没有完全的体现出来。本文将介绍作者自己改进的BGHstiff程序的一些细节。作者本着库
4、风格(就大小而言,完全可以升级)对BGH法进行了一些改进,增加了一些新特点。包括积分进程的有效流动,更加彻底和准确的误差控制,最后还有对减小计算量的完全探索。第一个版本的BGHstiff程序是在流行的PC/DOS平台下开发的。对于UNIX平台的兼容也正在考虑当中。二、Brayton-Gustavson-Hatchel法(BGH法)的基本原理似乎对BGH方法进行一个简要的介绍是明智的。我们设带求解的常微分方程如下:其中f是一个向量(方程组).t的取值空间[0,T]被分布不均的点t0,t1,t2……,tn-k,tn-k+1…,tN分为长度为h=Δtn=tn+1-tn的小区间.在(1)式
5、中的导数可以用下面的k次向后微分公式(修正)进行取代:这样,常微分方程系统(1)就被转化为了下一个时间点tn+1的非线性方程系统.开始时间t0,初始条件x0下k=1,(2)式可解,于是时间t1下的向量x1也可以求出来.在求解向量x2的时候,多项式的次数k可以取到2.在6步之后,微分方程的k可以在1到6之间变化,其值取决于精确性和计算速度的标准.方程式(1)的预告向量xPn+1也可以用BDF法计算,公式如下:其中,预告向量xPn+1最好用做在修正期间迭代求解下式的初始值:预告公式(3)同样适用于求解局部舍位误差.可以证明,局部舍位误差在第k步的方法下可用下式表示:这里是导数的真实值.
6、同时,xn-1=x(tn-1)+0(hk+1).预告向量xPn+1由(3)式定义,而(4)式可以求出xn+1.综上所述,BGH法可以归纳为以下五步:(BGH-1)建立新的时间点tn+1=tn+h(BGH-2)由预告公式(3)求出xPn+1的值.(BGH-3)用xPn+1作初值,再加上由(2)式求得的,解方程组(4),得到xn+1的值.(BGH-4)根据(5)式求得局部舍位误差,如果它不在规定的范围之内,那么需要改变积分步长或者所定k的值,回到步骤(BGH-1)(BGH-5)把时间n增加到n+1,进入下一步的(BGH-1)方程式(3)完全由k次多项式的预告系数γi来决定.这种情况下,
7、需要用到如下的k+1次的代数方程:其中并且我们规定(tj-tj)0=1同样的,校正系数αi也由时间点t0,t1,t2……,tn-k,tn-k+1…,tN来决定,用于求解(2)式.它也是由一个k+1词的代数方程来求解的:其中在ZimmermanandBaker的程序里面,(6)式和(11)式的求解都是用最简单的线性方程逼近法稍加改动而成的.三、实际算法对计算量的减少在Brayton,GustavsonandHatchel的原文中,他们提到了一种更为有效的对预告参数γi和