第十一章 常微分方程边值问题的数值解法工程技术与科学实验中提出的大量问题是常微分方程边值问题.本章将研究常微分方程边值问题的数值求解方法.主要介绍三种边界条件下的定解问题和两大类求解边值问题的数值方法,打靶法算法和有限差分方法.11.1 引言在很多实际问题中都会遇到求解常微分方程边值问题. 考虑如下形式的二阶常微分方程),,(y y x f y '='', b x a <<, (11.1.1)在如下三种边界条件下的定解问题: 第一种边界条件:α=)(a y , β=)(b y (11.1.2)第二种边界条件:α=')(a y , β=')(b y (11.1.2)第三种边界条件:⎩⎨⎧=-'=-'1010)()()()(b b y b y a a y a y βα, (11.1.13) 其中0 0, ,00000>+≥≥b a b a .常微分方程边值问题有很多不同解法, 本书仅介绍打靶方法和有限差分方法.11.2 打靶法对于二阶非线性边值问题()()().,,βα==≤≤'=''b y a y b x a y y x f y ,,, (11.2.1)打靶法近似于使用初值求解的情况. 我们需要利用一个如下形式问题初值解的序列:()()v a w a w b x a w w x f w ='=≤≤'='')(,,,,,α, (11.2.2)引进参数v 以近似原边界值问题的解.选择参数k v v =,以使:()()β==∞→b y v b w k k ,lim , (11.2.3)其中),(k v x w 定义为初值问题(11.2.2)在k v v =时的解,同时()x y 定义为边值问题(11.2.1)的解.首先定义参数0v ,沿着如下初值问题解的曲线,可以求出点),(αa 对应的初始正视图()()v a w a w b x a w w x f w ='=≤≤'='')(,,,,,α. (11.2.4)如果),(0v b w 不严格收敛于β,那么我们选择1v 等值以修正近似值,直到),(0v b w 严格逼近β.为了取得合适的参数k v ,现在假定边值问题(11.2.1)有唯一解,如果),(v x w 定义为初始问题(11.2.2)的解,那么v 可由下式确定:0),(=-βv b w . (11.2.5)由于这是一个非线性方程,我们可以利用Newton 法求解.首先选择初始值0v ,然后由下式生成序列),)(()),((111-----=k k k k v b dvdw v b w v v β,此处),(),)((11--=k k v b dv dwv b dv dw , (11.2.6)同时要求求得),)((1-k v b dvdw,因为),(v b w 的表达式未知,所以求解这个有一点难度;我们只能得到这么一系列的值。
,,,),(),(),(),(1210-⋯⋯k v b w v b w v b w v b w假如我们如下改写初值问题(11.2.2),使其强调解对x 和v 的依赖性()()v v a w v a w b x a v x w v x w x f w ='=≤≤'=''),(,),(),,(,,,,α,(11.2.7)保留初始记号以显式与x 的微分相关.既然要求当k v v =时),)((v b dvdw的值,那么我们需要求出表达式(11.2.7)关于v 的偏导数.过程如下:)),(),,(,(),(v x w v x w x v fv x v w '∂∂=∂''∂),()),(),,(,()),(),,(,(v x v w v x w v x w x w f v x v x w v x w x x f ∂∂'∂∂+∂∂'∂∂=),()),(),,(,(v x v w v x w v x w x w f ∂'∂''∂∂+又因为x 跟v 相互独立,所以当b x a ≤≤上式如下;),()),(),,(,(),()),(),,(,(),(v x vw v x w v x w x w f v x v w v x w v x w x w f v x v w ∂'∂''∂∂+∂∂'∂∂=∂''∂(11.2.8)初始条件为:1),(0),(=∂'∂=∂∂v a v w v a v w , 如果简单地用),(v x z 定义),)((v x vw∂∂,并且假定x 和v 的微分阶翻转,(11.2.8)转化为初值问题:1)(0)(),,(),,(....==≤≤''∂∂+'∂∂=a z a z b x a z w w x w fz w w x w f z ,,,,(11.2.9) 因此,对于(11.2.2)和(11.2.9)式的每次迭代需要求解两个初值问题.那么从(11.2.6)式可得:),(),(111-----=k k k k v b z v b w v v β, (11.2.10)事实上,这些初值问题很难精确求解,而将这些解近似为一个初值问题的解.同样,我们可以按以上步骤考虑对于三阶非线性边值问题的打靶法算法. 对于三阶非线性边值问题()()().)(,,,,γβα='='=≤≤'''='''b y a y a y b x a y y y x f y ,,,(11.2.11)转变形式:()()v a w a w a w b x a w w w x f w =''='=≤≤'''=''')()(,,,,,,,βα,(11.2.12)选择参数k v v =,以使:()()β=='∞→b y v b w k k ,lim , (11.2.13)其中),(k v x w '定义为初值问题(11.2.12)在k v v =时的解,同时()x y 定义为边值问题(11.2.11)的解.定义参数0v ,沿着如下初值问题解的曲线,可以求出点),(αa 对应的初始正视图()()v a w a w a w b x a w w w x f w =''='=≤≤'''=''')()(,,,,,,,βα.(11.2.14)如果),(0v b w '不严格收敛于β,那么我们选择1v 等值以修正近似值,直到),(0v b w '严格逼近β.为了取得合适的参数k v ,现在假定边值问题(11.2.11)有唯一解,如果),(v x w 定义为初始问题(11.2.12)的解,那么v 可由下式确定:0),(=-'βv b w . (11.2.15)利用Newton 法求解.首先选择初始值0v ,然后由下式生成序列),)(()),((111---'-'-=k k k k v b dvw d v b w v v β,此处),(),)((11--'='k k v b dv w d v b dv w d , (11.2.16)同时要求求得),)((1-'k v b dvw d ,因为),(v b w '的表达式未知,所以只能得到这么一系列的值。
,,,),(),(),(),(1210-'⋯⋯'''k v b w v b w v b w v b w改写初值问题(11.2.12),使其强调解对x 和v 的依赖性()(),(,),(,),(,),(,)(,)w f x w x v w x v w x v a x b w a v w a v w a v vαβ''''''='''≤≤===,,,, (11.2.17)要求当k v v =时),)((v b dvw d '的值,那么我们需要求出表达式(11.2.17)关于v 的偏导数.所以当b x a ≤≤上式如下;vw w f v w w f v w w f v x v w ∂''∂''∂∂+∂'∂'∂∂+∂∂∂∂=∂'''∂),( (11.2.18) 初始条件为:1),(0),(0),(=∂''∂=∂'∂=∂∂v a vw v a v w v a v w ,, 对于第三种边界条件, 注意到 01)(a a y'(a)a y -=, 取 s a y =)('得到 ⎪⎩⎪⎨⎧=-==s a , y'a a a y'a y x,y,y'f y )()()()(01, b x a <<, (11.2.19) 设它的解为)(x,s y , 代入第二个边界条件得到s βs F s b y βb,s y b y βb y ==-=-')() ,()(')()(00即s βs F =)(当1s ββ= 时, 即为所求的解. 这样,同第一种边界条件一样, 可以利用打靶方法求解.例11.1 利用打靶法求解两点边值问题2()0(0)0,(0.4) 1.8y y x y y '''⎧-+=⎨==⎩ 解:(1) 为应用打靶方法,需要假定初值(0)y '来先求解初值问题,取初始参数0 1.8120.40t -==-, 解2()(0)0,(0)2y y x y y '''⎧=-⎨'==⎩令1y u =,2y u '=,则上述初值问题变为1212222,(0)1,(0)2u u u u u x u ⎧'==⎪⎨'=-=⎪⎩, 由标准的R-K 方法,取0.2h =,可得0 2.540β=。
(2)再令11t =, 解1212222,(0)1,(0)1u u u u u x u ⎧'==⎪⎨'=-=⎪⎩, 取0.2h =。
仍由标准的R-K 方法,可得1 1.497 1.8β=<。
(3)令2121(1.8 1.497) 1.291.497 2.540t -=+-=-再解1212222,(0)1,(0) 1.29u u u u u x u ⎧'==⎪⎨'=-=⎪⎩, 取0.11h =,得2 1.7094β=。