差分方程方法

1、第四章差分方程方法在实际中，许多问题所研究的变量都是离散的形式，所建立的数学模型也是离散的，譬如，像政治、经济和社会等领域中的实际问题。有些时候，即使所建立的数学模型是连续形 式，例如像常见的微分方程模型、积分方程模型等等，但是，往往都需要用计算机求数值解。这就需要将连续变量在一定条件下进行离散化，从而将连续型模型转化为离散型模型，因此，最后都归结为求解离散形式的差分方程解的问题。关于差分方程理论和求解方法在数学建模和解决实际问题的过程中起着重要作用。卜面就不同类型的差分方程进行讨论。所谓的差分方程是指：对于一个数列xn ,把数列中的前n 1项xi i 0,1,2, n关联起来所得到的方程。4

2、 . 1常系数线性差分方程4.1.1常系数线性齐次差分方程常系数线性齐次差分方程的一般形式为其中k为差分方程的阶数,XnaXn 1a2Xn 23kXn k 0（4.1）ai i 1,2, ,k为差分方程的系数，且 ak 0 k n。对应的代数方程kk 1k 2a1a2ak 0(4.2)称为差分方程的（4.1 ）的特征方程，其特征方程的根称为 特征根。常系数线性齐次差分方程的解主要是由相应的特征根的不同情况有不同的形式。下面分 别就特征根为单根、重根和复根的情况给出差分方程解的形式。1. 特征根为单根设差分方程（4.1）有k个单特征根1, 2, 3, k，则差分方程（4.1 ）的通解为1747

3、。7 ,XnC1 nC2Ck n，其中G,C2,70k为任意常数，且当给定初始条件Xi i i 1,2, ,k（4.3）时，可以唯一确定一个特解。2. 特征根为重根设差分方程（4.1 ）有l个相异的特征根1, 2, 3, I 1 k重数分别为Im17 m2,7 ml且mik则差分方程（4.1）的通解为xnmim2i 1 ni 1 nGE 1C2in2i 1i 1mii 1cHni 1同样的，由给定的初始条件（3.特征根为复根4.3 ）可以唯一确定一个特解。设差分方程（4.1 ）的特征根为一对共辆复根1,i和相异的k 2个单根Ckk，则差分方程的通解为XnC1 n cosn C2 nsinnc3

4、 nC4 4其中 22 ,arCtan.同样由给定的初始条件（4.3 ）可以惟一确定一个特解。另外，对于有多个共轲复根和相异实根，或共轲复根和重根的情况，都可以类似地给出差分方程解的形式。4. 1. 2常系数线性非齐次差分方程常系数线性非齐次差分方程的一般形式为Xn吨 1 a2Xn 2akXn kf n(4.4)其中k为差分方程的阶数, 为已知函数。ai i I2，k为差分方程的系数,ak 0 k n f (n)在差分方程（4.4）中，令f（n） 0,所得方程XnaXn1a2Xn 2akXn k 0(4.5)4.1 ）的形式相同。称为非齐次差分方程（4.4 ）对应的齐次差分方程，即与差分方程（

5、求解非齐次差分方程通解的一般方法为首先求对应的齐次差分万程（4.5 ）的通解Xn ,然后求非齐次差分方程（4.4 ）的一个特解Xn° ,则 参照常系数线性非齐次微分方程求特解的方法。Xn*Xn0Xn为非齐次差分方程（4.4 ）的特解。关于求Xn的方法同求差分方程（4.1 ）的方法相同。对于求非齐次方程（4.4）的特解Xn0的方法，可以用观察法确定，也可以根据f（n）的特性用待定系数法确定，具体方法可4.2 差分方程的平衡点及其稳定性一般来说，差分方程的求解是困难，实际中往往不需要求出差分方程的一般解，而只需要研究它的平衡点及其稳定性即可。4.2.1 一阶线性常系数差分方程一阶线性常系

6、数差分方程的一般形式为Xk 1 aXk b, k 0,1,2,. 、 . 、 . . . . . . . . . . . . . . . . . *其中a, b为常数，它的平衡点由代数方程x ax b求解得到，不妨记为 x .如果lim xkx ,则称平衡点x是稳定的，否则是 不稳定的。k k为了便于研究平衡点 x的稳定性问题，一般将其转化为求方程xk 1 axk 0的平衡点* x0的稳定性问题。事实上，由xk 1 axk 0可以解得kxkax。，一 .于是x0是稳定的平衡点的充要条件是：a 1.4.2.2 一阶线性常系数差分方程组一阶线性常系数齐次差分方程组的一般形式为x k 1 Ax k

7、B,k 0,1,2,其中x k为n维向量，A为n n阶常数矩阵。 = . . . . .它的平衡点x0是稳定的充要条件是A的所有特征根都有i 1 (i=1,2,n).对于一阶线性常系数非齐次差分方程组x(k 1) Ax(k) B , k 0,1,2,的情况同样给出4.2.3 二阶线性常系数差分方程二阶线性常系数齐次差分方程的一般形式为xk 2alxk 1 a2xk0,k 0,1,2,其中aha2为常数，其平衡点 x*0是稳定的充要条件是特征方程2aia20 ,的根1, 2满足11,21。对于一般的Xk 2 a1Xk1 a2Xk 0的平衡点的稳定性问题同样给出。类似地，也可直接推广到n阶线性差分

8、方程的情况。4.2.4 一阶非线性差分方程一阶非线性方程的一般形式为Xk1f Xk, k 0,1,2,其中f为已知函数，其平衡点定义为方程 X f,、. ， - . 事实上，将f Xk在X处作一阶的泰勒展开有(*Xk 1 f X Xk X则X也是一阶线性方程 Xk 1 f X Xk X的充要条件是| f(X )1。X的解X 。- * f X ,一*i. 一. i.*f X 的平衡点，故此，平稳衡点 X稳定4.3连续模型的差分方法4.3.1 微分的差分方程已知 f(X)在点 Xk处的函数值 f(Xk) k 0,1,.,n 1), 且. _ ax0X1. Xn 1 b ,试求函数的导致值f Xk

9、, (k 1,2,.n)。根据导数的定义，用差商代替微商，则有下面的差分公式。向前差：f (Xk)f(Xk1) f(Xk)(k 1,2,.,n),Xk 1Xk向后差:f (Xk)f(Xk) f(Xk)Xk Xk 1中心差:f (Xk)f(Xk1) f(Xk)xk 1 xk 1(k 1,2,.,n),4.3.2 定积分的差分方法1,2,., n),且在a,b上可积，试求函已知函数f(x)在点Xk处的函数值f (Xk) , (kb数在a,b上的积分值a f (x)dxo根据定积分的定义，则有一般的求积公式bnf(x)dxAk f(Xk)ak 0其中Ak为求积系数，它与Xk的选取方法有关。取不同的求

10、积系数，可以得不同的求积公式。对于等距节点xkakh (k0,1,.,n),其中步长hba为很小的数，则有如 nki(2)复化梯形矩阵;ba f(x)dxh n 1f(Xk)2 k 0f (Xk1)I / f (a) 21f (Xk) f (b)1(3)复化辛普森(Simpson)矩阵公式;bh n 1a f(x)dx 6 f (Xk) 4f(xk1)2f(Xk1)(4)n 1f(a) 4 f(xk 1)k 0 k 2f(Xk) f(b)1.一(Xk Xk 1) 为子区间Xk,Xk 1 (k2复化柯特斯(Cotes)公式;0,1,., n 1)的中点。下的求积公式。(1)复化矩阵公式;ba f

11、(x)dxn 112 f(Xk) 7f(b) k 1n 1n 1f (x)dx 7f (a) 32 f(x 1) 12 f(x 1) 32 f (x 3)90k 0 k4 k 0 k2 k 0 k 4其中x 1,x 1 ,x 3为子空间xk, Xk i(k0,1,.,n 1)中的四等分点。k k k ： 2444.3.3 常微分方程的差分方法1. 一阶常微分方程的差分方法设一阶常微分方程的定解问题为(y f(x,y),y(x°) y°,其中函数f(x, y)关于y满足李普希兹条件，即保证问题解的存在唯一性。现在的问题是求方程在一系列节点x1x2 . xn.处的近似数值解yi

12、,y2,., yn,.不妨假设步长为h xn i xn为常数。在此，我们根据微分的差分方法，即用差商来近似代替微商，再利用“步进式”方法，可以给出求解问题(4-6)的差分方法。(1)单步欧拉(Euler)公式用差商y(xn 1)一以立2近似代替y'(xn)f (xn,y(xn)中的导数，则可以得差分公式hyn 1yn hf(xn,yn),n 0,1,2,.2 .其精度为O(h )阶的。(2)两步欧拉公式用差商y(xn1) y(xn1)近似代替y'(xn)f(xn,y(xn)中的导数，则可得差分公式2hyn 1yn 1 2hf (xn , yn), n 0,1,2,.两步法需要用

13、到前两步的方信息，一般不能自行起步，需先用单步方法求出y1 ,其精度是O(h2)阶的。(3)梯形公式对于方程y' f(x, y)的两边在xn,xn 1上求积分得xn 1 ry(xn 1) y(xn) y(xn)f(x,y(x)dxxn利用积分的差分方法中梯形公式求解积分xn 1hf(x,y(x)dxf (xn,y(xn) f(xn 1,y(xn 1)xn2h .y(Xn 1)y(Xn) -f(Xn, y(Xn) f (Xn 1 , y(Xn 1)离散化即可得到微分方程的梯形差分公式hyn 1Yn 21f(Xn,yN)f (Xn 1,yN 1), N 0,1,2,.这是一个隐式格式，计算

14、量大，一般不单独使用。其镜的也是O(h3)阶的。(4)改进的欧拉公式由于单步欧拉公式色精度低，但计算量小；矩形公式精度高,但是计算量大，为此,我们综合运用这两种方法舅老爷得到改进的欧拉公式，其精度为一 3 一 .O h 阶的。预报:ynynhfxn, yn , n 0,1,校正:ynynf Xn,yn fxn 1, yn 1,n 0,1,或写成平均化形式;ypYnhfXn, Yn ,ycYnhfXn 17yp ,3ypyc), n0,1,2,(5)龙格一库塔法龙格一库塔方法的基本思想:对于微分方程的定解问题(4.6 ),考虑差商y(Xn 1)h也),根据阿格朗日微分中值定理可得y(Xn 1)

15、y Xnhyy Xn hf，y ， Xn,y*,称为 Xn, Xn 1上的平均变化率，则 y(Xn 1) y XnhY。现在一、 一 . . 一、 一 .一 *的问题只要找到寻找一种好的计算Y的近似方法。t . Mr-t T-TT* 一*_.如果取Y f Xn , yn丫1 ,则就是欧拉公式。,*1如果取Y - f Xn, ynf Xn 1,yn 1 Y2,则相应的就是改进的欧拉公式。现在，我们取m个点Xn ih, yn ihXn,Xn1i 1,2, , m，用f在这m个点的函数值的加权平均作为 Y的近似值，即-*Yi f Xnih,yn ihi 1其中i为权系数。则有myn 1 yn h i

16、 f Xn Efn ih i 1(4.7)其中 一i, i为待定系数。实际上，适当选择i , i , i ,使得公式有更高的精度，这是龙格-库塔方法的思想。二阶龙格-库塔公式:在Xn, Xn i内取中点X i n _ 21Xn -h,则可取 10, 2 1 ,21一代人（4.7 ）2式得到二阶龙格-库塔公式，其精度为 O h3阶。yn 1 yn hY2,Y f Xn , yn ,工 h h.y2 f Xn -,yn -Y ，n 0,1,2,三阶龙格-库塔公式：在Xn, Xn 1内任取点 Xn p Xnph, Xn q Xn qh, 0 p q 1类似的方法可得到三阶的龙格-库塔公式yn 1yn

17、 7 Yi 4丫2 丫3 ,6Y f Xn,yn , h hY>f Xn -,yn -x ,n 0,1,2,K f Xn h,yn h Y1 2Y2其精度是Oh4阶的常用的是三阶的情况。四阶龙格-库塔公式:类似的方法可以得到四阶龙格 -库塔公式，其精度是 O h5阶的yni yn -h Yl 2Y2 2Y3 Y4 ,Y f Xn,yn ,Y2 f Xn h,yn hY ,n 0,1.2,22f Xnh2，yn2Y2Y4f(Xni,yn hK)2. 一阶常微方程组的差分方法将前面的单个方程中的变量和函数视为向量, 一阶方程组的情形。对于二个方程的方程组相应的差分方法即可用于由多个方程组的f

18、 X.y,z ,y Xoyo,(4.8)g X, y, z ,z XoZo,设以yn,zn表示函数在节点Xnx0 nh,n1,2,上的近似解，则有改进的欧拉公式:预报:VnynhfXn, yn, zn校正:四阶龙格-库塔公式znhfynynf Xn, ynZf Xn 1,n1 yn 1,zn 1znzng Xn,ynZg xn 1,yn 1,zn 1Ynznh-Y1 2Y2 2Y3hZ1 2Z2 2Z36123丫4 ,n 0,1,2, Z4其中YiZi丫2Z2丫3Z3f Xn, yn,Zn ；Xn,yn,Zn ；xnXnXnXnh2,ynh2，ynh2,ynh2,ynh丫1, Zn22Yi,z

19、n2丫2,Zn2丫2,Zn2Z2丫4f Xn i, yn hY3,Zn hZ3Z4 g Xn i，yn hY3,Zn hZs其他的公式也都可以类似得到，即相当于同时求解多个一阶方程， 从方法上没有本质的差别。3.高阶常微分方程的差分方法对于某些高阶方法的定解问题，原则上可以转化为一阶方程组来求解。譬如，对于如下的二阶微分方程的定解问题'''y f X, y, y''y(Xo) yo,y x°y若令Z y ,则可以化为一阶方程组的定解问题 (Z f X, y,Zy Z,(4.9)''y(Xo)yo,y Xov。实际上，(4.9)式

20、可以视为(4.8)式的特例，类似地可以得到相应的求解差分公式。4.4最优捕鱼问题4.3.i问题的提出假设艇鱼可分为 4个年龄组：称i、2、3、4龄鱼。各年龄组每条鱼的平均重量分别为 5.o7 , ii.55 , i7.86 , 22.99 ( g );各年龄组鱼的自然死亡率均为o.8 (i/年)；这种鱼为季节性集中产卵繁殖，产卵孵化期为每年的最后4个月，平均每条4龄鱼的产卵量为5i.io9* io (个)，3龄鱼的广卵量为这个数的一半，2龄和i龄鱼不广卵。卵孵化并成活为i龄鱼，成活率(i龄鱼条数与产卵量n之比)为 iii.22 io一 ii.22 io n渔业部门规定，每年只允许在产卵孵化期前

21、的8个月内进行捕捞作业。 如果每年投入的捕捞能力固定不变，即固定努力量捕捞，这时单位时间捕捞量将与各年龄组鱼群条数成正比， 比例系数称为捕捞强度系数。通常使用13mm网眼的拉网，这种网只能捕捞3、4龄鱼，其两个捕捞系数之比为 o.42 : i.要解决的问题是：建立数学模型，分析如何实现可持续性捕捞（即每年开始捕捞时渔场中各年龄组鱼群 条数不变），并且在此前提下得到最高的年收获量（总质量）。4.3.2 模型的假设与符号说明1 .模型的假设（1）只考虑鱼的繁殖和捕捞的变化，不考虑鱼群迁入与迁出；（2）各龄鱼在一年的任何时间都会发生自然死亡；（3）所有鱼都在每年最后四个月内（后 1/3年）完成产卵孵

22、化的过程，成活的幼鱼在 下一年初成为1龄鱼；（4）产卵发生于后四个月之初，产卵鱼的自然死亡发生与产卵之后；（5）相邻两个年龄组的鱼群在相邻两年之间变化是连续的，即第 k年底i龄鱼的条数 等于第k 1年初i 1龄鱼的条数；（6） 4龄以上的鱼全部死亡；（7）采用固定努力量捕捞意味着捕捞的速率正比于捕捞时各龄鱼群的条数，比例系数 为捕捞强度系数。2 .符号的说明用x"t）表示t时刻（年）i龄鱼的条数；r表示鱼的平均自然死亡率，即 r 0.8; %A .表小龄i鱼的广卵数，即（f1, f2,f3, f4） （0,0,-, A） , A 1.109 105; Wi表示龄i鱼群的捕捞强度系数，

23、即（Wi，W2，W3，W4）（5.07,11.55,17.86,22.99）； qi表示i龄鱼群的捕捞2强度系数，即 &72q。） （Q0,0.42e,E） , E为捕捞努力量；t 表示产卵开始的 3Y月份；Y表示i龄鱼的捕捞量；Ci表示i龄鱼的捕捞率，即 Ci oXi4.4.3模型的建立与求解1.无捕捞时鱼群的自然增长模型dxi tdt由假设（1）和（2）得rxi t , i 1,2,3,4 ; k t k 1 , k 0,1,2,又由假设（3）和（4）得_11x1 k 11.22*10,Xo k_111.22*10x0 k txo k tx3 k tAx4 k t由假设（5）和（6

24、）得x, 0x1 , x 1k 1x, k 1limi (t)(i 1,2,3;k 0,1,)2。固定努力量捕捞鱼群的增长和捕捞模型由假设知，捕捞期为k t k t,则有则有dxi tdtrxi t qi E xi t , k t k t,(4.10)dxi tdtrxi t , k t t k 1(4.11)xi k 1 , i 1,2,3,(4.12)x1 k 1_ 111.22*10,x° k_ 111.22*10x0 k t(4.13)A x0 k tx3 k tAx4 k t , k 0,1,2,2(1)鱼群的增长规律求解方程(4.10 )和(4.11 ),并利用连续条件(

25、4.12 )式可得xi 1 k 1 sli E Xi k , i 1,2,3,(4.15)x(k 1)b x0 kx0(4.16)x0 k 1A ttsEx3 k t As I4 E x4 k t 2(4.17)其中 s e r 0.4993, I1EI2 E1, I3 Ee 0.42tE, A 1.109 105, b 1.22 1011。(2)捕捞量单位时间第i龄鱼的捕捞量（条数）为yi tqi E Xi t , k t k t ,第k年全年（8个月）第i龄鱼的捕捞量（条数）为tt（1） qi EtYi knyit dtQiEXitdtqi匚（1 stli E）Xi k00r qi E于是

26、，第k年总捕捞量（质量）为W（k）W3Y3*） W4Y4*）（3）可持续性捕捞模型可持续捕捞，即意味着由于自然死亡和捕捞使鱼群减少，而通过产卵繁殖补充，使得鱼群能够在每年初开始捕捞时保持平衡不变，这样的捕捞策略就可以年复一年地一直持续下去。因此，可持续捕捞的鱼群数应是（4.15 ）、（4.16 ）、（4.17）式的平衡解，即模型不依赖于时间 t 的解 xi （i 0,1,2,3,4）。求解（4.15）、（4.16）、（4.17）得*一*/ C Cxi 1sli E xi , i1,2,3,* *x2slx，x3sx2.*bx0xi r,b bx0s xi , x4S13（E）x3 s313 E

27、 xi*_ t _ *_ t _ *x0 0.5Astl3（E）x3 Astl4（E）x4将（4.18）式代入（4.20）式得8*o*x。0.5 s14（E）As313（E）xi代入（4.19 ）式有8*b0.5 s14（E）As313（E）x*x18，二 *b 0.5 sl4（E）As3l3（E）x1求解可得*, / Ixb 1 -B E代入（4.18）式得到*x2sb 1*X3sb 1*3X4s I3 E b 18其中 B E 0.5 sl4(E)As3l3(E).当B(E) 1时，xi 0即意味着捕捞过度，致使鱼群火绝。当B(E) 1时,Eo 31.4称之为过度捕捞力量，因此，可以在E Eo的范围内找最优捕捞策略。在可持续性捕捞白条件下，第i龄鱼的年捕捞量(条数)为t*qi E 1 s li E XiY , i 3,4。r qi E整个鱼群的年捕捞量