数值计算9章.docVIP

第9章　常微分方法的数值解法 科学技术与工程问题常常需要建立微分方程形式的数学模型,下面是这类问题的例子. 设N(t)为某物种的数量,α为该物种的出生率与死亡率之差,β为生物的食物供给及它们所占空间的限制,描述该物种增长率的数学模型是 设Q是电容器上的带电量,C为电容,R为电阻,E为电源的电动势,描述该电容器充电过程的数学模型是 以上两个例子是常微分方程的初值问题,下面是一个两点边值问题的例子. 设一根长为L的矩形截面的梁,两端固定.E是弹性模量,S是端点作用力,I(x)是惯性矩,q是均匀载荷强度,梁的桡度y(x)满足如下方程 针对实际问题建立的数学模型,要找出模型解的解析表达式往往是困难的,甚至是不可能的.因此,需要研究和掌握微分方程的数值解法,即计算域内离散点上解的近似值的方法.本章讨论常微分方程数值解的基本方法和理论. 9.1　Euler方法 9.1.1　Euler方法及其有关的方法 考虑一阶常微分方程的初值问题 (9.1) 设f(x,y)是连续函数,对y满足Lipschitz条件,即存在正数L,使得对于任意两点(x,y)与(x,),有|f(x,y)-f(x,)|≤L|y-|,这样初值问题的解是存在唯一的,而且连续依赖于初始条件. 为了求得离散点上的函数值,将微分方程的连续问题(9.1)式进行离散化.一般是引入点列{xn},这里xn=xn-1+hn,n=1,2,….称hn为步长,经常考虑定长的情形,即hn=h,xn=x0+nh, n=0,1,…. 记y(xn)为初始问题(9.1)的准确解y(x)在xn处的值,用均差近似代替(9.1)式中的导数得 令yn为y(xn)的近似值,将上面两个近似式写成等式,整理后得 yn+1=yn+hf(xn,yn),n=0,1,…, (9.2) yn+1=yn+hf(xn+1,yn+1),n=0,1,…. (9.3) 从x0处的初值y0开始,按(9.2)式可逐步计算以后各点上的近似值.称(9.2)式为显式Euler公式.由于(9.3)式的右端隐含有待求函数值yn+1,不能逐步显式计算,称(9.3)式为隐式Euler公式或后退Euler公式. 如果将(9.2)式和(9.3)式作算术平均,就得梯形公式 yn+1=yn+[f(xn,yn)+f(xn+1,yn+1)],n=0,1,…. (9.4) 梯形公式也是隐式公式.以上公式都是由yn去计算yn+1,故称它们为单步法. 例9.1　取h=0.1,用Euler方法?隐式Euler方法和梯形方法解 解　本题有f(x,y)=x-y+1,y0=1.如果用Euler方法,由(9.2)式并代入h=0.1得 yn+1=0.1xn+0.9yn+0.1. 同理,用隐式Euler方法有 yn+1=(0.1xn+1+yn+0.1). 用梯形公式有 yn+1=(0.1xn+0.95yn+0.105). 3种方法及准确解y(x)=x+e-x的数值结果如表9-1所示.从表中可看到,在xn=0.5处,Euler方法和隐式Euler方法的误差|y(xn)-yn|分别是1.6×10-2和1.4×10-2,而梯形方法的误差却是2.5×10-4. 表9-1 xn Euler方法 隐式Euler方法 梯形法 准确解 0 1 1 1 1 0.1 1.000000 1.009091 1.004762 1.004837 0.2 1.010000 1.026446 1.018594 1.018731 0.3 1.029000 1.051315 1.040633 1.040818 0.4 1.056100 1.083013 1.070096 1.070320 0.5 1.090490 1.120921 1.106278 1.106531 在例9.1中,由于f(x,y)对y是线性的,所以对隐式公式也可方便地计算yn+1.但是,当f(x,y)是y的非线性函数时,如y′=5x+,其隐式Euler公式为yn+1=yn+h(5xn+1+).显然,它是yn+1的非线性方程,可以选择非线性方程求根的迭代法求解yn+1.以梯形公式为例,可用显式Euler公式提供迭代初值,用迭代公式 反复迭代,直到, 其中,步长h成为迭代参数,它需要满足一定的条件迭代公式才能收敛.若将(9.4)式减去该迭代公式,得yn+1-=(f(xn+1,yn+1)-f(xn+1,)). 假设f(x,y)关于y满足Lipschitz条件,则有 yn+1-≤yn+1-, 这里,L是Lipschitz常数.由上式可见,当1即h时,迭代序列{}收敛于yn+1. 对于隐式公式,通常采用预估－校正技术,即先用显式公式计算,得到预估值,然后以预估值为隐式公式的迭代初值,用隐式公式迭代一次得到校正值,称为