用Matlab求解差分方程问题.ppt

《用Matlab求解差分方程问题.ppt》由会员分享，可在线阅读，更多相关《用Matlab求解差分方程问题.ppt（35页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、用Matlab求解差分方程问题 Still waters run deep.流静水深流静水深,人静心深人静心深 Where there is life,there is hope。有生命必有希望。有生命必有希望一、一阶线性常系数差分方程一、一阶线性常系数差分方程濒濒危物种的自然演危物种的自然演变变和人工孵化和人工孵化问题问题 Florida沙丘沙丘鹤鹤属于属于濒濒危物种，它在危物种，它在较较好好自然自然环环境下，年均增境下，年均增长长率率仅为仅为1.94%，而在中，而在中等和等和较较差差环环境下年均增境下年均增长长率分率分别为别为-3.24%和和 -3.82%，如果在某自然保，如果在某自然保护

2、护区内开始有区内开始有100只只鹤鹤，建立描述其数量建立描述其数量变变化化规规律的模型，并作律的模型，并作 数数值计值计算。算。模型建立模型建立记记第第k年沙丘年沙丘鹤鹤的数量的数量为为xk,年均增年均增长长率率为为r，则则第第k+1年年鹤鹤的数量的数量为为 xk+1=(1+r)xk k=0,1,2已知已知x0=100,在在较较好，中等和好，中等和较较差的自然差的自然环环境下境下 r=0.0194,-0.0324,和和-0.0382 我我们们利用利用Matlab编编程，程，递递推推20年后年后观观察沙丘察沙丘鹤鹤的的数量数量变变化情况化情况Matlab实现实现首先建立一个关于首先建立一个关于变

3、变量量n，r的函数的函数function x=sqh(n,r)a=1+r;x=100;for k=1:n x(k+1)=a*x(k);end在在command窗口里窗口里调调用用sqh函数函数 k=(0:20);y1=sqh(20,0.0194);y2=sqh(20,-0.0324);y3=sqh(20,-0.0382);round(k,y1,y2,y3)To Matlab（ff6）利用plot 绘图观察数量变化趋势可以用不同可以用不同线线型和型和颜颜色色绘图绘图r g b c m y k w 分分别别表示表示 红绿兰兰绿红绿兰兰绿洋洋红红黄黑白色黄黑白色：+o*.X s d 表示不同的表示不

4、同的线线型型 plot(k,y1,k,y2,k,y3)在同一坐在同一坐标标系下画系下画图图 plot(k,y2,:)plot(k,y2,-)plot(k,y2,r)plot(k,y2,y)plot(k,y2,y,k,y1,:)plot(k,y2,k,y1,:)plot(k,y2,oy,k,y1,:)用gtext(r=0.0194),gtext(r=-0.0324),gtext(r=-0.0382)在图上做标记。To Matlab（ff6）人工孵化是挽救人工孵化是挽救濒濒危物种的措施之一，如危物种的措施之一，如果每年孵化果每年孵化5只只鹤鹤放入保放入保护护区，区，观观察在中等察在中等自然条件下沙

5、丘自然条件下沙丘鹤鹤的数量如何的数量如何变变化化Xk+1=aXk+5 ,a=1+r如果我们想考察每年孵化多少只比较合适，可以令Xk+1=aXk+b ,a=1+rfunction x=fhsqh(n,r,b)a=1+r;x=100;for k=1:nx(k+1)=a*x(k)+b;endk=(0:20);%一个行向量一个行向量y1=(20,-0.0324,5);%也是一个行向量也是一个行向量round(k,y 1 )%对对k,y1四舍五入，但四舍五入，但%是是 不改不改变变变变量的量的值值 plot(k,y1)%k y1 是行向量列向量都可以是行向量列向量都可以也可以也可以观观察察20年的年的发

6、发展展趋势趋势，以及在，以及在较较差差条件下的条件下的发发展展趋势趋势，也可以考察每年孵，也可以考察每年孵化数量化数量变变化的影响。化的影响。To Matlab（ff7）高阶线性常系数差分方程高阶线性常系数差分方程 如果第如果第k+1时时段段变变量量Xk+1不不仅仅取决取决于第于第k时时段段变变量量Xk，而且与以前，而且与以前时时段段变变量有关，就要用高量有关，就要用高阶阶差分方程来描述差分方程来描述一年生植物的繁殖一年生植物的繁殖一年生植物春季一年生植物春季发发芽，夏天开花，秋季芽，夏天开花，秋季产产种，没有腐种，没有腐烂烂，风风干，干，没没被人被人为为掠取掠取的那些种子可以活的那些种子可以

7、活过过冬天，其中一部分冬天，其中一部分能在第能在第2年春季年春季发发芽，然后开花，芽，然后开花，产产种，种，其中的另一部分其中的另一部分虽虽未能未能发发芽，但如又能芽，但如又能活活过过一个冬天，一个冬天，则则其中一部分可在第三其中一部分可在第三年春季年春季发发芽，然后开花，芽，然后开花，产产种，如此种，如此继继续续，一年生植物只能活，一年生植物只能活1年，而近似的年，而近似的认认为为，种子最多可以活，种子最多可以活过过两个冬天，两个冬天，试试建建立数学模型研究立数学模型研究这这种植物数量种植物数量变变化的化的规规律，及它能一直繁殖下去的条件。律，及它能一直繁殖下去的条件。模型及其求解模型及其求

8、解记记一棵植物春季一棵植物春季产产种的平均数种的平均数为为c,种子能种子能活活过过一个冬天的一个冬天的(1岁岁种子种子)比例比例为为b,活活过过一个冬天没有一个冬天没有发发芽又活芽又活过过一个冬天的（一个冬天的（2岁岁种子）比例仍种子）比例仍为为b,1岁岁种子种子发发芽率芽率a1，2岁岁种子种子发发芽率芽率a2。设设c,a1,a2固定，固定，b是是变变量，考察能一直繁殖的条件量，考察能一直繁殖的条件记记第第k年植物数量年植物数量为为Xk，显显然然Xk与与Xk-1,Xk-2有关，由有关，由 Xk-1决定的部分是决定的部分是 a1bcXk-1,由由Xk-2决定的部分是决定的部分是 a2b(1-a1

9、)bcXk-2 Xk=a1bcXk-1 +a2b(1-a1)bcXk-2 Xk=a1bcXk-1 +a2b(1-a1)bcXk-2实际实际上，就是上，就是Xk=pXk-1+qXk-2 我我们们需要需要知道知道x0,a1,a2,c,考察考察b不同不同时时，种子繁殖，种子繁殖的情况。在的情况。在这这里假里假设设X0=100,a1=0.5,a2=0.25,c=10,b=0.180.20这样这样可以用可以用matlab计计算了算了Xk=a1bcXk-1 +a2b(1-a1)bcXk-2function x=zwfz(x0,n,b)c=10;a1=0.5;a2=0.25;p=a1*b*c;q=a2*b*

10、(1-a1)*b*c;x(1)=x0;x(2)=p*(x(1);for k=3:nx(k)=p*(x(k-1)+q*(x(k-2);endk=(0:20);y1=zwfz(100,21,0.18);y2=zwfz(100,21,0.19);y3=zwfz(100,21,0.20);round(k,y1,y2,y3);plot(k,y1,k,y2,:,k,y3,o),gtext(b=0.18);gtext(b=0.19);gtext(b=0.20)To Matlab（ff8）结果分析：Xk=pXk-1+qXk-2 (1)x1+px0=0 (2)差分方程的特征方程差分方程的特征方程差分方程的特征根

11、：差分方程的特征根：方程方程(1)的解可以表的解可以表为为C1,c2 由初始条件由初始条件x0,x1确定。确定。本例中，用待定系数的方法可以求出本例中，用待定系数的方法可以求出b=0.18时时,c1=95.64,c2=4.36 ，这样这样实际实际上，上，植物能一直繁殖下去的条件是植物能一直繁殖下去的条件是b0.191线性常系数差分方程组线性常系数差分方程组汽汽车车租租赁赁公司的运公司的运营营一家汽一家汽车车租租赁赁公司在公司在3个相个相邻邻的城市运的城市运营营，为为方便方便顾顾客起客起见见公司承公司承诺诺，在一个城市租，在一个城市租赁赁的汽的汽车车可以在任意一个城市可以在任意一个城市归还归还。

12、根据。根据经验经验估估计计和市和市场调查场调查，一个租一个租赁赁期内在期内在A市租市租赁赁的汽的汽车车在在A,B,C市市归还归还的比例分的比例分别为别为0.6,0.3,0.1;在在B市租市租赁赁的汽的汽车归还车归还比例比例0.2,0.7,0.1;C市租市租赁赁的的归还归还比比例分例分别为别为0.1,0.3,0.6。若公司开。若公司开业时业时将将600辆辆汽汽车车平均分配到平均分配到3个城市，建立运个城市，建立运营过营过程中汽程中汽车车数量在数量在3个城市个城市间转间转移的模型，并移的模型，并讨论时间讨论时间充充分分长长以后的以后的变变化化趋势趋势。0.60.3A B CA B CA B C假假

13、设设在在每个租每个租赁赁期开期开始能把始能把汽汽车车都都租出去，租出去，并都在并都在租租赁赁期期末末归还归还0.10.70.20.10.60.30.1模型及其求解模型及其求解记记第第k个租个租赁赁期末公司在期末公司在ABC市的汽市的汽车车数量数量分分别为别为x1(k),x2(k),x3(k)（也是第（也是第k+1个租个租赁赁期开始各个城市租出去的汽期开始各个城市租出去的汽车车数量），很数量），很容易写出第容易写出第k+1个租个租赁赁期末公司在期末公司在ABC市的市的汽汽车车数量数量为为(k=0,1,2,3)用矩用矩阵阵表示表示用用matlab编编程，程，计计算算x(k),观观察察n年以后的年以

14、后的3个城市的个城市的汽汽车车数量数量变变化情况化情况function x=czqc(n)A=0.6,0.2,0.1;0.3,0.7,0.3;0.1,0.1,0.6;x(:,1)=200,200,200;for k=1:n x(:,k+1)=A*x(:,k);end如果直接看如果直接看10年或者年或者20年年发发展展趋势趋势，可以直接在命令窗，可以直接在命令窗口（口（commond window）作，而不是必）作，而不是必须编须编一个函数一个函数A=0.6,0.2,0.1;0.3,0.7,0.3;0.1,0.1,0.6;n=10;for k=1:nx(:,1)=200,200,200;x(:,

15、k+1)=A*x(:,k);end round(x)作图观察数量变化趋势 k=0:10;plot(k,x)，gridgtext(x1(k),gtext(x2(k),gtext(x3(k)可以看到可以看到时间时间充分充分长长以后以后3个城市汽个城市汽车车数量数量趋趋于于180，300，120可以考察可以考察这这个个结结果与初始条件是否有关果与初始条件是否有关若最开始若最开始600辆辆汽汽车车都在都在A市，可以看到市，可以看到变变化化时间时间充分充分长长以后，各城市汽以后，各城市汽车车数量数量趋趋于于稳稳定，与初始定，与初始值值无关无关直接输入x（:,1）的值即可x(:,1)=600,0,0;ro

16、und(x);plot(k,x),grid6.6 按年龄分组的人口模型按年龄分组的人口模型 不同年龄组的繁殖率和死亡率不同不同年龄组的繁殖率和死亡率不同.建立差分方程模型，讨论稳定状况下种群的增长规律建立差分方程模型，讨论稳定状况下种群的增长规律.假设与建模假设与建模 种群按年龄大小等分为种群按年龄大小等分为n个年龄组，记个年龄组，记i=1,2,n 时间离散为时段，长度与年龄组区间相等，记时间离散为时段，长度与年龄组区间相等，记k=1,2,以雌性个体数量为对象以雌性个体数量为对象.第第i 年龄组年龄组1雌性个体在雌性个体在1时段内的时段内的繁殖率繁殖率为为bi 第第i 年龄组在年龄组在1时段内

17、的死亡率为时段内的死亡率为di,存活率存活率为为si=1-di假设假设与与建模建模xi(k)时段时段k第第i 年龄组的种群数量年龄组的种群数量按年龄组的分布向量按年龄组的分布向量预测任意时段种群预测任意时段种群按年龄组的分布按年龄组的分布Leslie矩阵矩阵(L矩阵矩阵)(设至少设至少1个个bi0)按年龄分组的种群增长野生或野生或饲饲养的养的动动物因繁殖而增加，因自然死亡物因繁殖而增加，因自然死亡和人和人为为屠屠杀杀而减少，不同年而减少，不同年龄动龄动物的繁殖率，物的繁殖率，死亡率有死亡率有较较大差大差别别，因此在研究某一种群数量，因此在研究某一种群数量的的变变化化时时，需要考，需要考虑虑年年

18、龄龄分分组组的种群增的种群增长长。将种群按年将种群按年龄龄等等间间隔的分成若干个年隔的分成若干个年龄组龄组，时时间间也离散化也离散化为时为时段，段，给给定各年定各年龄组龄组种群的繁殖种群的繁殖率和死亡率，建立按年率和死亡率，建立按年龄龄分分组组的种群增的种群增长长模型，模型，预测预测未来各年未来各年龄组龄组的种群数量，并的种群数量，并讨论时间讨论时间充充分分长长以后的以后的变变化化趋势趋势。模型及其求解模型及其求解设设种群按年种群按年龄龄等等间间隔的分成隔的分成n个年个年龄组龄组，记记i=1,2,，n,时时段段记记作作k=0,1,2,且年且年龄组龄组区区间间与与时时段段长长度相等度相等(若若5

19、岁为岁为一个年一个年龄组龄组，则则5年年为为一一个个时时段段)。以雌性个体。以雌性个体为为研究研究对对象象记记在在时时段段k第第i年年龄组龄组的数量的数量为为xi(k);第第i年年龄组龄组的的繁殖率繁殖率为为bi，表示每个个体在一个，表示每个个体在一个时时段内繁殖段内繁殖的数量；第的数量；第i年年龄组龄组死亡率死亡率为为di，表示一个，表示一个时时段段内死亡数与内死亡数与总总数的比，数的比，si=1-di是存活率。是存活率。记记在在时时段段k种群各年种群各年龄组龄组的数量的数量为为X(k)=x1(k),x2(k),xn(k)这样这样，有，有x(k+1)=Lx(k),k=0,1,给给定在定在0时

20、时段，各年段，各年龄组龄组的初始数量的初始数量x(0)就可以就可以预测预测任意任意时时段段k,各年各年龄组龄组的数量的数量设设一种群分成一种群分成5个年个年龄组龄组，繁殖率繁殖率b1=0,b2=0.2,b3=1.8,b4=0.8,b5=0.2存活率存活率s1=0.5,s2=0.8,s3=0.8,s4=0.1各年各年龄组现龄组现有数量都是有数量都是100只，只，用用matlab计计算算x(k)b=0,0.2,1.8,0.8,0.2;s=diag(0.5,0.8,0.8,0.1);L=b;s,zeros(4,1);x(:,1)=100*ones(5,1);n=30;for k=1:nx(:,k+1

21、)=L*x(:,k);end round(x)k=0:30;subplot(1,2,1),plot(k,x),gridTo Matlab（ff9）将将x(k)归归一化后的向量一化后的向量记记做做x(k),称称为为种群种群按年按年龄组龄组的分布向量的分布向量，即各年，即各年龄组龄组在在k时时段在数量上占段在数量上占总总数的百分比。数的百分比。y=diag(1./sum(x);%sum(x)对对列求和列求和Z=x*y Subplot(1,2,2),plot(k,z),grid结结果分析：果分析：时间时间充分充分长长以后，种群按年以后，种群按年龄组龄组的分布的分布x(k)趋趋向向稳稳定。定。To Matlab（ff9）