数模微分方程的数值求解技术.docVIP

PAGE PAGE 15 常微分方程的数值求解技术 1 Euler方法 2 Runge-Kutta方法 3 线性多步方法 1 Euler方法 1.1 从一个人口增长（Malthus）模型谈起 英国人口统计学家马尔萨斯（Malthus,1766-1834年）调查了英国100多年人口出生统计资料，发现人口出生率稳定于一个常数。并提出了著名的Malthus人口增长模型。此模型的基本假设是：在考虑一个国家或地区的人口总数随时间变化的人口自然增长过程中，略去迁移、自然环境条件等因素对人口变化的影响，视净相对增长率（出生率与死亡率之差）是常数。即单位时间内人口的增长量与人口成正比，比例系数为。此假设就是一个常微分方程的问题。 记为时刻的人口，视为连续、可微函数。据Malthus的假设，在到时间内人口的增长量为，又设时的人口为，令时，可得满足方程 假设取1790年为初始状态，已知与1800年的实际人口，由最小二乘法拟合，可以定出。计算1790年～1880年间的人口（百万），。即求解微分方程 以下先介绍Euler方法，再给出此微分方程的数值解。 1）． Euler公式 （1.2） 2）.改进的Euler公式 (1.5) 此公式称为改进的Euler公式。它也可写成一个显式公式 (1.6) 或 (1.7) 3）．人口增长（Malthus）模型的题解 例1 取1790年为初始状态，已知与1800年的实际人口，由最小二乘法拟合，可以定出。用改进的Euler公式计算1790年～1880年间的人口（百万），即求解微分方程 其中，。 解 改进的Euler公式 具体的计算公式为 ，计算结果如下表6-1。 表6-1 用改进的Euler公式估计1790～1880年份的人口数据（106） 年份 1790 1800 1810 1820 1830 1840 1850 1860 1870 1880 ti 0 1 2 3 4 5 6 7 8 9 xi 3.9 5.2811 7.1512 9.6837 13.1129 17.7565 24.0445 32.5593 44.0893 59.7024 1.2 软件介绍 在本节里，将通过数值例子和Matlab程序介绍本节提出的Euler公式及改进的Euler公式的上机实现方法。 1.根据式（6.1.2）及算法6-1所表示的Euler公式，编写Matlab程序(函数名：Euler_f.m)。 function [x,y]=Euler_f(fun1,x0,y0,h,N) % fun1 为一阶微分方程的函数； % x0,y0为初始条件； % h为区间步长； % N为区间的个数； % x为Xn构成的向量； % y为Yn构成的向量. x=zeros(1,N+1); y=zeros(1,N+1); x(1)=x0; y(1)=y0; for n=1:N x(n+1)=x(n)+h; y(n+1)=y(n)+h*feval(fun1,x(n),y(n)); end 2.根据式（6.1.7）及算法6-2所表示的改进的Euler公式，编写Matlab程序(函数名：Euler_r.m)。 function [x,y]=Euler_r(fun1,x0,y0,h,N) % fun1 为一阶微分方程的函数； % x0,y0为初始条件； % h为区间步长； % N为区间的个数； % x为Xn构成的向量； % y为Yn构成的向量. x=zeros(1,N+1); y=zeros(1,N+1); x(1)=x0; y(1)=y0; for n=1:N x(n+1)=x(n)+h; ybar=y(n)+h*feval(fun1,x(n),y(n)); y(n+1)=y(n)+h/2*(feval(fun1,x(n),y(n))+feval(fun1,x(n+1),ybar)); end 例2 采用Euler公式（Euler_f.m）和改进的Euler公式（Euler_r.m）的Matlab程序上机实现人口指数增长（Malthus）模型所反映的微分方程的题解。其中，。 解 利用Matlsb编程过程如下： 求精确解 y=dsolve(‘Dy=0.307*y’,’y(0)=3.9’ 运行后显示精确解 y=39/10*exp(307/1000*x); （1）建立并保存名为fun11.m的M文件函数。 function f=fun1(x,y) f=0.307*y; （2）建立并保存名为Euler_f.m和Euler_r.m的M文件函