《微分方程模型 (2)幻灯片.ppt》由会员分享,可在线阅读,更多相关《微分方程模型 (2)幻灯片.ppt(116页珍藏版)》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、微分方程模型课件第1页,共116页,编辑于2022年,星期六第一部分微分方程模型简介第2页,共116页,编辑于2022年,星期六 在研究实际问题时,常常会联系到某些变量的变化率或导数,这样所得到变量之间的关系式就是微分方模型。微分方程模型反映的是变量之间的间接关系,因此,要得到直接关系,就得求微分方程。求解微分方程有三种方法:1)求精确解;2)求数值解(近似解);3)定性理论方法。第3页,共116页,编辑于2022年,星期六 微分方程模型是一类应用十分广泛而且最常见的数学模型,其建模方法在数学模型竞赛中占有重要的地位。1.微分方程模型是机理分析建模方法的最佳体现:根据对象的特征和目的,对问题进
2、行必要的简化和假设,用精确的语言做出适合的数学语言表达的假设,然后根据所作的假设来分析对象的因果关系,利用对象的内在规律建立 各个量之间的等式关系或其他数学结构。2.微分方程模型是物理、生物进化等定律精确的定量描述。第4页,共116页,编辑于2022年,星期六(1)根据规律列方程 利用数学、力学、物理、化学等学科中的定理或经过实验检验的规律等来建立微分方程模型。如牛顿第二定律、放射性物质的放射性规律等.我们常利用这些规律对某些实际问题列出微分方程.(2)微元分析法 利用已知的定理与规律寻找微元之间的关系式,与第一种方法不同的是对微元而不是直接对函数及其导数应用规律。建立微分方程模型的方法第5页
3、,共116页,编辑于2022年,星期六(3)模拟近似法 在生物、经济等学科的实际问题中,许多现象的规律性不很清楚,即使有所了解也是极其复杂的,建模时在不同的假设下去模拟实际的现象,建立能近似反映问题的微分方程,然后从数学上求解或分析所建方程及其解的性质,再去同实际情况对比,检验此模型能否刻画、模拟某些实际现象。第6页,共116页,编辑于2022年,星期六微分模型的建模原理在建立微分方程的时候,所要求的其实是微分方程在建立微分方程的时候,所要求的其实是微分方程的一条解曲线,通过它来反映某些我们所要寻求的的一条解曲线,通过它来反映某些我们所要寻求的规律。微分方程曲线思想是,如果知道曲线上每一规律。
4、微分方程曲线思想是,如果知道曲线上每一点处的导数以及它的起始点,那么就能构造这条曲点处的导数以及它的起始点,那么就能构造这条曲 线。具体步骤如下:线。具体步骤如下:1、转化、转化:实际问题中,有许多表示实际问题中,有许多表示“导数导数”的常用词,如的常用词,如“速率速率”、”增长增长“(在生物学以及人口问题研究中)、(在生物学以及人口问题研究中)、”衰衰变变“(在放射性问题中)以及(在放射性问题中)以及”边际的边际的“(在经济学中)(在经济学中)等。这些词就是信号,这个时候要注意是哪些研究对象等。这些词就是信号,这个时候要注意是哪些研究对象在变化,对这些规律的表示微分方程也许就能用得上。在变化
5、,对这些规律的表示微分方程也许就能用得上。第7页,共116页,编辑于2022年,星期六2、准确性和总体特征、准确性和总体特征微分方程式一个在任何时刻都必须正确的微分方程式一个在任何时刻都必须正确的瞬时表达式,这是问题的核心。建立微分瞬时表达式,这是问题的核心。建立微分方程模型,首先要把注意力放在方程文字方程模型,首先要把注意力放在方程文字形式的总关系上:形式的总关系上:净变化率净变化率=输入率输入率输出率输出率或者:或者:变化率(微商)变化率(微商)=单位增加量单位增加量-单位减少量单位减少量等式通常是利用已有的原则或定律。等式通常是利用已有的原则或定律。第8页,共116页,编辑于2022年,
6、星期六3、单位、单位一旦确定了哪些子项应该列入微分方程中,一旦确定了哪些子项应该列入微分方程中,就要确保每一项都采用同样的物理单位。就要确保每一项都采用同样的物理单位。这是在建立微分方程过程中容易疏忽的问题。这是在建立微分方程过程中容易疏忽的问题。第9页,共116页,编辑于2022年,星期六4、约束条件、约束条件约束条件是关于所研究对象在某一特定时刻约束条件是关于所研究对象在某一特定时刻的信息(比如初始时刻),它们独立于微分的信息(比如初始时刻),它们独立于微分方程而存在。在建立微分方程模型后,利用方程而存在。在建立微分方程模型后,利用它们来确定模型中有关的常数,这些常数包它们来确定模型中有关
7、的常数,这些常数包括比例系数、原微分方程的其他参数和解中括比例系数、原微分方程的其他参数和解中的积分常数。为了完整,充分地给出问题的的积分常数。为了完整,充分地给出问题的数学陈述,建模过程中应该将这些约束条件数学陈述,建模过程中应该将这些约束条件和微分方程一起写出。和微分方程一起写出。第10页,共116页,编辑于2022年,星期六5、概念框架、概念框架前面阐述的都是使用微分方程建模的关键问题。当面临前面阐述的都是使用微分方程建模的关键问题。当面临一个典型问题是,首先必须有一个明确的概念框架一个典型问题是,首先必须有一个明确的概念框架(建立其他模型也是如此),这个概念框架就是关键步骤。(建立其他
8、模型也是如此),这个概念框架就是关键步骤。具体如下:具体如下:(1)把用语言描述的情况转化为文字方程。)把用语言描述的情况转化为文字方程。(2)陈述出所涉及的原则或物理定律。)陈述出所涉及的原则或物理定律。(3)建立微分方程,配备方程各子项的单位。)建立微分方程,配备方程各子项的单位。(4)给定约束条件,包括初始条件或其他条件。)给定约束条件,包括初始条件或其他条件。(5)给出微分方程的解。)给出微分方程的解。(6)求出微分方程的常数。)求出微分方程的常数。(7)给出问题答案。)给出问题答案。(8)检验答案是否满足问题的要求。)检验答案是否满足问题的要求。在建模过程中,明确了概念框架,然后就是
9、依次完成在建模过程中,明确了概念框架,然后就是依次完成框架中每一步所要做的事情。框架中每一步所要做的事情。第11页,共116页,编辑于2022年,星期六建立微分方程模型的一般准则(1)转化翻译:有许多表示导数的常用词,如速率、增长、衰变、边际、弹性等。改变、变化、增加、减少这些词可能是一种暗示信号,只需弄清楚什么在变,随什么而变,这时也许导数就用得上。(2)机理分析:将所研究的问题看成一个封闭和系统,思考研究的问题是否遵循什么原理或物理定律,是应该用已知的定律还是去推导问题的合适结果?在不知道问题的机理时,合理的想象和类比是很重要的!对平衡式:变化率=输入率-输出率。能理解它,并且能使用正确的
10、物理量纲,或许就得到了需要的微分方程。第12页,共116页,编辑于2022年,星期六(3)微分方程模型:微分方程是在任何时刻必须正确的瞬时表达式。如看到了表示导数的关键词,就要寻找y(t)与y(t),t的关系。首先将注意力集中在文字形式的关系式变化率=输入率-输出率上,然后准确填好式中的所有项。(4)单位:一旦确信了哪些项应该列入微分方程中,就要确保第一项都彩同样的物理单位,保证式子的平衡。(5)定解条件:系统在某一特定时刻的信息,独立于微分方程而成立。利用它们来确定有关的常数,包括比例系数、原微分方程的其它参数以及解中的积分系数。第13页,共116页,编辑于2022年,星期六小案例-求解析解
11、1 国民生产总值2 环境污染问题3 刑事侦察中死亡时间的鉴定第14页,共116页,编辑于2022年,星期六 1999年我国的国民生产总值(年我国的国民生产总值(GDP)为)为80 423亿元,如果亿元,如果我国能保持每年我国能保持每年8%的相对增长率,问到的相对增长率,问到2010年我国的年我国的GDP是多少?是多少?国民生产总值国民生产总值 第15页,共116页,编辑于2022年,星期六解解 (1)(1)建立微分方程建立微分方程记记t=0代表代表1999年,并设第年,并设第t年我国的年我国的GDP为为P(t)由题意由题意知,从知,从1999年起,年起,P(t)的的相对增长率为相对增长率为8%
12、,即,即得微分方程得微分方程第16页,共116页,编辑于2022年,星期六(2 2)求通解)求通解分离变量得分离变量得方程两边同时积分,得方程两边同时积分,得第17页,共116页,编辑于2022年,星期六(3 3)求特解)求特解将将p(0)=80423代入通解,得代入通解,得C=80423,所以从,所以从1999年起年起第第t年年我国的我国的GDP为为将将t=2010-1999=11代入上式,得代入上式,得2010年我国的年我国的GDP的的预测值为预测值为 第18页,共116页,编辑于2022年,星期六 环境污染问题环境污染问题 某水塘原有某水塘原有50000t清水清水(不含有害杂质不含有害杂
13、质),从时间,从时间t=0开始,开始,含有有害杂质含有有害杂质5%的浊水流入该水塘流入的速度为的浊水流入该水塘流入的速度为2t/min,在塘中充分混合,在塘中充分混合(不考虑沉淀不考虑沉淀)后又以后又以2t/min的速度流出水的速度流出水塘问经过多长时间后塘中有害物质的浓度达到塘问经过多长时间后塘中有害物质的浓度达到4%?第19页,共116页,编辑于2022年,星期六解解 (1 1)建立微分方程)建立微分方程 设在时刻设在时刻t塘中有害物质的含量为塘中有害物质的含量为 ,此时塘中,此时塘中单位时间内有害物质的变化量单位时间内有害物质的变化量=(=(单位时间内流进塘内有害物质的量单位时间内流进塘
14、内有害物质的量)-(-(单位时间内流出塘的有害物质的量单位时间内流出塘的有害物质的量)有害物质的浓度为有害物质的浓度为 ,于是有,于是有第20页,共116页,编辑于2022年,星期六即即 ,(1 1)初始条件为初始条件为Q(0)=0.(2 2)求通解)求通解 式(式(1 1)是可分离变量方程,分离变量得)是可分离变量方程,分离变量得第21页,共116页,编辑于2022年,星期六积分,得积分,得即即 第22页,共116页,编辑于2022年,星期六(3 3)求特解)求特解由初始条件由初始条件t=0,Q=0得得C=-2500,故,故当塘中有害物质浓度达到当塘中有害物质浓度达到4%时,应有时,应有由此
15、解得由此解得(min)即经过即经过670.6min后,塘中有害物质浓度达到后,塘中有害物质浓度达到4%,塘中有害物质的最终浓度为,塘中有害物质的最终浓度为由于由于第23页,共116页,编辑于2022年,星期六 刑事侦察中死亡时间的鉴定刑事侦察中死亡时间的鉴定 当一次谋杀发生后,尸体的温度从原来的当一次谋杀发生后,尸体的温度从原来的3737按照牛顿冷却定按照牛顿冷却定律开始下降,如果两个小时后尸体温度变为律开始下降,如果两个小时后尸体温度变为3535,并且假定,并且假定周围空气的温度保持周围空气的温度保持2020不变,试求出尸体温度不变,试求出尸体温度H H随时间随时间t的变的变化规律又如果尸体
16、发现时的温度是化规律又如果尸体发现时的温度是3030,时间是下午,时间是下午4 4点点整,那么谋杀是何时发生的?整,那么谋杀是何时发生的?第24页,共116页,编辑于2022年,星期六注注:牛顿冷却定律指出牛顿冷却定律指出,物体在空气中冷却的速度与物体在空气中冷却的速度与物体温度和空气温度之差成正比,现将牛顿冷却物体温度和空气温度之差成正比,现将牛顿冷却 定律应用于刑事侦察中死亡时间的鉴定定律应用于刑事侦察中死亡时间的鉴定第25页,共116页,编辑于2022年,星期六解解(1 1)建立微分方程)建立微分方程的冷却速度的冷却速度其中其中k0是常数,初始条件为是常数,初始条件为H(0)=37 设尸
17、体的温度为设尸体的温度为H(t)(t从谋杀后计从谋杀后计),根据题意,尸体,根据题意,尸体正比即正比即与尸体温度与尸体温度H和空气温度和空气温度20之差成之差成第26页,共116页,编辑于2022年,星期六分离变量得分离变量得 (2 2)求通解)求通解积分得积分得 第27页,共116页,编辑于2022年,星期六 把把初初值值条条件件H(0)=37代代入入通通解解,求求得得C=17于于是是该该初初值值问问题题的解为的解为 为求出为求出k值,根据两小时后尸体温度为值,根据两小时后尸体温度为3535这一条件,这一条件,有有(3 3)求特解)求特解第28页,共116页,编辑于2022年,星期六求得求得
18、 ,于是温度函数为,于是温度函数为将将H=30代入式代入式(1)有有 ,即得,即得 (h)。于是,于是,(1)(1)可以判定谋杀发生在下午可以判定谋杀发生在下午4点尸体被发现前的点尸体被发现前的8.4h,即,即8小时小时24分钟,所以谋杀是在上午分钟,所以谋杀是在上午7点点36分发生的分发生的第29页,共116页,编辑于2022年,星期六微分方程的解析解的MATLAB实现 求微分方程(组)的解析解命令求微分方程(组)的解析解命令:dsolve(方程方程1,方程方程2,方程方程n,初始条件初始条件,自变量自变量)结结 果:果:u=tan(t+C1)第30页,共116页,编辑于2022年,星期六
19、解解 输入命令输入命令:y=dsolve(D2y+4*Dy+29*y=0,y(0)=0,Dy(0)=15,x)结结 果果 为为:y=3*exp(-2*x)*sin(5*x)第31页,共116页,编辑于2022年,星期六解解 输入命令输入命令:x,y,z=dsolve(Dx=2*x-3*y+3*z,Dy=4*x-5*y+3*z,Dz=4*x-4*y+2*z,t);结结 果果 为:为:x=C2*exp(2*t)+C3*exp(-t)y=C2*exp(2*t)+C3*exp(-t)+exp(-2*t)*C1z=C2*exp(2*t)+exp(-2*t)*C1第32页,共116页,编辑于2022年,星
20、期六微分方程的数值解(一)常微分方程数值解的定义 在生产和科研中所处理的微分方程往往很复杂且大多得不出一般解。而在实际上对初值问题,一般是要求得到解在若干个点上满足规定精确度的近似值,或者得到一个满足精确度要求的便于计算的表达式。因此,研究常微分方程的数值解法是十分必要的。第33页,共116页,编辑于2022年,星期六(二)建立数值解法的一些途径(二)建立数值解法的一些途径1、用差商代替导数、用差商代替导数 若步长若步长h较小,则有较小,则有故有公式:故有公式:此即欧拉法。此即欧拉法。第34页,共116页,编辑于2022年,星期六对于欧拉公式(1)我们看到,当 时公式右端的 都是近似的,所以用
21、它计算的 会有累积误差,分析累积误差比较复杂,这里先讨论比较简单的所谓局部截断误差。假定用(1)式时右端的 没有误差,即 ,那么由此算出 (2)局部截断误差指的是,按(2)式计算由 到 这一步的计算值 与精确值 之差 为了估计它,由Taylor展开得到的精确值 是 (3)(2)-(3),注意到 ,得 (4)第35页,共116页,编辑于2022年,星期六即局部截断误差是 阶的,而数值算法的精度定义为:若一种算法的局部截断误差为 ,则称该算法具有p阶精度。显然p 越大,方法的精度越高。式(4)说明,欧拉方法是一阶方法,因此它的精度不高。第36页,共116页,编辑于2022年,星期六欧拉方法的Mat
22、lab 实现程序如下:function x,y=euler(fun,x0,xfinal,y0,n);if nargin5,n=50;endh=(xfinal-x0)/n;x(1)=x0;y(1)=y0;for i=1:nx(i+1)=x(i)+h;y(i+1)=y(i)+h*feval(fun,x(i),y(i);end第37页,共116页,编辑于2022年,星期六2、使用数值积分、使用数值积分对方程对方程 两边由两边由 到到 积分,并利用梯形公式,有:积分,并利用梯形公式,有:故有公式:故有公式:(5)实际应用时,与欧拉公式结合使用:实际应用时,与欧拉公式结合使用:此即改进的欧拉法。此即改进
23、的欧拉法。第38页,共116页,编辑于2022年,星期六由(5)式可知,所以改进的欧拉方法是二阶的,精度较欧拉方法要高,实用性更加广泛。第39页,共116页,编辑于2022年,星期六为了编程方便,常把(3)式改写为第40页,共116页,编辑于2022年,星期六则改进欧拉方法的matlab实现程序为:function x,y=eulerpro(fun,x0,xfinal,y0,n);if nargin5,n=50;endh=(xfinal-x0)/n;x(1)=x0;y(1)=y0;for i=1:nx(i+1)=x(i)+h;y1=y(i)+h*feval(fun,x(i),y(i);y2=y
24、(i)+h*feval(fun,x(i+1),y1);y(i+1)=(y1+y2)/2;end第41页,共116页,编辑于2022年,星期六龙格-库塔方法MATLAB运算程序为 x,y=ode23(f,ts,y0,tol)ode23为2、3阶龙格-库塔公式 x,y=ode45(f,ts,y0,tol)ode45为4、5阶龙格-库塔公式 x,y=ode15s(f,ts,y0,tol)ode15s同ode45 f为常微分方程(组)函数,ts表示自变量的求解区间或范围,y0为函数的初值,tol为误差限 调出结束后输出变量x及其函数值y第42页,共116页,编辑于2022年,星期六基本思想:利用基本思
25、想:利用 在某些特殊点上的函数值的在某些特殊点上的函数值的线性组合来构造高阶单步法的平均斜率。线性组合来构造高阶单步法的平均斜率。什么叫平均斜率?什么叫平均斜率?对差商对差商 应用微分中值定理,有,应用微分中值定理,有,利用微分方程利用微分方程 ,有,有这里的这里的 称为平均斜率。称为平均斜率。第43页,共116页,编辑于2022年,星期六可将改进的欧拉格式改写成可将改进的欧拉格式改写成的算术平均值作为平均斜率。的算术平均值作为平均斜率。该公式可以看作是用该公式可以看作是用 和和 两个点处的斜率两个点处的斜率 和和由改进型欧拉公式我们可以猜想,如果在由改进型欧拉公式我们可以猜想,如果在内多预测
26、几个点的斜率,再对他们进行加权平均,内多预测几个点的斜率,再对他们进行加权平均,可能得到精度更好的平均斜率!可能得到精度更好的平均斜率!第44页,共116页,编辑于2022年,星期六下面以下面以2阶龙格阶龙格-库塔方法为例来阐述这种思想库塔方法为例来阐述这种思想考察区间考察区间 上的一点上的一点 ,用用 和和 的斜率的斜率 和和 的加权平均作为平均的加权平均作为平均斜率斜率 的近似值:的近似值:即取即取其中其中 和和 是待定常数。若取是待定常数。若取 ,则,则问题在于如何确定问题在于如何确定 处的斜率处的斜率 和常数和常数 和和 。第45页,共116页,编辑于2022年,星期六仿照改进的欧拉方
27、法,用欧拉方法预测仿照改进的欧拉方法,用欧拉方法预测 的值,的值,并用它来估计斜率并用它来估计斜率 :于是得到如下形式的算法:于是得到如下形式的算法:通过适当选取参数通过适当选取参数 和和 的值,使得公式具有的值,使得公式具有2阶精度!阶精度!第46页,共116页,编辑于2022年,星期六由泰勒公式展开,要使公式具有由泰勒公式展开,要使公式具有 2 阶精度,只需阶精度,只需方程组有无穷多解:方程组有无穷多解:二级方法有无穷多种二级方法有无穷多种常见的常见的3种二级方法:种二级方法:中点法(修正的中点法(修正的Euler法法)取取 二阶二阶龙格库塔方法龙格库塔方法取取第47页,共116页,编辑于
28、2022年,星期六三级方法:三级方法:N=3类似于类似于N=2的推导方法,可得到的推导方法,可得到常见的常见的2种三阶方法:种三阶方法:库塔库塔三阶方法三阶方法第48页,共116页,编辑于2022年,星期六 四级方法:四级方法:N=4局部截断误差局部截断误差常见的常见的2种四阶方法:种四阶方法:经典经典龙格龙格-库塔库塔方法方法第49页,共116页,编辑于2022年,星期六 1、在解、在解n个未知函数的方程组时,个未知函数的方程组时,x0和和x均为均为n维向量,维向量,m-文件中的待解方程组应以文件中的待解方程组应以x的分量形式写成的分量形式写成.2、使用、使用Matlab软件求数值解时,高阶
29、微分方程必须等价软件求数值解时,高阶微分方程必须等价地变换成一阶微分方程组地变换成一阶微分方程组.注意注意:第50页,共116页,编辑于2022年,星期六解解:令令 y1=x,y2=y11、建立、建立m-文件文件vdp1000.m如下:如下:function dy=vdp1000(t,y)dy=zeros(2,1);dy(1)=y(2);dy(2)=1000*(1-y(1)2)*y(2)-y(1);2、取、取t0=0,tf=3000,输入命令:,输入命令:T,Y=ode15s(vdp1000,0 3000,2 0);plot(T,Y(:,1),-)3、结果如图、结果如图第51页,共116页,编
30、辑于2022年,星期六解解 1、建立、建立m-文件文件rigid.m如下:如下:function dy=rigid(t,y)dy=zeros(3,1);dy(1)=y(2)*y(3);dy(2)=-y(1)*y(3);dy(3)=-0.51*y(1)*y(2);2、取取t0=0,tf=12,输入命令:,输入命令:T,Y=ode45(rigid,0 12,0 1 1);plot(T,Y(:,1),-,T,Y(:,2),*,T,Y(:,3),+)3、结果如图、结果如图图中,图中,y1的图形为实线,的图形为实线,y2的图形为的图形为“*”线,线,y3的图形为的图形为“+”线线.第52页,共116页,
31、编辑于2022年,星期六第二部分微分方程模型理论介绍第53页,共116页,编辑于2022年,星期六1 微分方程的一般理论微分方程的一般理论微分方程的一般形式微分方程的一般形式一阶微分方程 (1)其中f(t,x)是t和x的已知函数,x(t0)=为初始条件,又称定解条件。一阶微分方程组 (2)又称为一阶正规方程组。如果引入向量 则方程组(2)可以写为简单的形式 第54页,共116页,编辑于2022年,星期六 (3)即与方程(1)的形式相同,当n=1时为方程(1)对于任一高阶的微分方程 如果记 则方程为 即可化为一阶方程组的形式。因此,下面主要对正规方程组(3)进行讨论。第55页,共116页,编辑于
32、2022年,星期六微分方程解的存在唯一性微分方程解的存在唯一性 定理定理1(Cauchy-Peano)如果函数 在R:上连续,则方程组(3)在 上有解 满足初值条件 此处定理定理2 如果函数 在 R:上连续,且满足利普希茨(Lipschitz)条件(即存在正常数L使得 其中 则方程组(3)满足初值条件 的解是唯一的。第56页,共116页,编辑于2022年,星期六2 微分方程的平衡点及稳定性微分方程的平衡点及稳定性(1)微分方程的平衡点微分方程的平衡点 设有微分方程组(3),对于 ,在某个区域内连续,且满足解的存在唯一性条件。如果存在某个常数 使得 则称点 为方程组(3)的平衡点(或奇点),且称
33、 为方程组的平凡解(或奇解)。第57页,共116页,编辑于2022年,星期六 如果对所有可能初值条件,方程组(3)的解 都满足 ,则称平衡点 是稳定的(渐近稳定);否则是不稳定的。实际中,判断平衡点的稳定性有两种方法:间接方法和直接方法。间接方法:首先求出方程的解 然后利用定义 来判断。直接方法:不用求方程的解直接地来研究其稳定性。第58页,共116页,编辑于2022年,星期六(2)一阶方程的平衡点及稳定性)一阶方程的平衡点及稳定性 设有微分方程 ,其相应的平衡点为代数方程 的实根 .其稳定性可以用间接方法判断,下面说明直接方法。首先,将函数 在点 作一阶泰勒(Taylor)展开,即方程可以近
34、似地表示为 显然,也是该方程的一个平衡点,其稳定性主要取决于 符号,即有下面结论:若 则平衡点是稳定的;若 则平衡点是不稳定的。第59页,共116页,编辑于2022年,星期六(3)平面方程的平衡点及稳定性)平面方程的平衡点及稳定性 设平面方程组的一般形式为 (6)则称代数方程组 的实根 为平面方程组(6)的平衡点,记为 .如果对所有可能的初值条件方程的解为 满足 则称平衡点 是稳定的;否则是不稳定的。第60页,共116页,编辑于2022年,星期六也可以用直接方法讨论:将方程组(6)的右边的函数作一阶泰勒展开,即可表示为近似的线性方程组 (7)记系数矩阵为 ,且假设其行列式 ,则方程组(7)的特
35、征方程为 ,即第61页,共116页,编辑于2022年,星期六其中 为特征根。不妨设特征根分别为 ,即 根据特征根 和系数 的取值情况可以确定平衡点 的稳定性。事实上,当 时平衡点是稳定的;当 或 时平衡点是不稳定的。对于一般微分方程的平衡点和稳定性问题可以类似的讨论。第62页,共116页,编辑于2022年,星期六第三部分微分方程建模实例第63页,共116页,编辑于2022年,星期六第64页,共116页,编辑于2022年,星期六第65页,共116页,编辑于2022年,星期六第66页,共116页,编辑于2022年,星期六第67页,共116页,编辑于2022年,星期六第68页,共116页,编辑于20
36、22年,星期六第69页,共116页,编辑于2022年,星期六第70页,共116页,编辑于2022年,星期六例3 饿狼追兔问题现有一只兔子,一只狼,兔子位于狼的正西现有一只兔子,一只狼,兔子位于狼的正西100100米处。假设兔子与狼同时发现对方并一米处。假设兔子与狼同时发现对方并一起起跑,兔子往正北起起跑,兔子往正北6060米处的巢穴跑,而米处的巢穴跑,而狼在追兔子,已知兔子、狼是匀速跑且狼狼在追兔子,已知兔子、狼是匀速跑且狼的速度是兔子的两倍。问题是兔子能否安的速度是兔子的两倍。问题是兔子能否安全回到巢穴?全回到巢穴?yxhy=f(x)B-60A(100,0)C(x,y)O第71页,共116页
37、,编辑于2022年,星期六解解 首先建立坐标系,兔子在首先建立坐标系,兔子在O O处,狼在处,狼在A A处。处。由于狼要盯着兔子追,所以狼行走的是一条曲由于狼要盯着兔子追,所以狼行走的是一条曲线,且在同一时刻,曲线上狼的位置与兔子的线,且在同一时刻,曲线上狼的位置与兔子的位置的连线为曲线上该点处的切线。设狼的行位置的连线为曲线上该点处的切线。设狼的行走轨迹是走轨迹是y=f(x)y=f(x),则有则有 ,又因狼的速度是兔子的两倍,所以在相同时间又因狼的速度是兔子的两倍,所以在相同时间内狼走的距离为兔子走的距离的两倍。假设在内狼走的距离为兔子走的距离的两倍。假设在某一时刻,兔子跑到某一时刻,兔子跑
38、到(0,h)(0,h)处,而狼在处,而狼在(x,y)(x,y)处,处,则有则有第72页,共116页,编辑于2022年,星期六整理得到下述模型整理得到下述模型这属于可降阶的二阶微分方程,解得狼的行走轨迹这属于可降阶的二阶微分方程,解得狼的行走轨迹因因 ,所以狼追不上兔子。,所以狼追不上兔子。第73页,共116页,编辑于2022年,星期六例4 传染病模型传染病模型 随随着着卫卫生生设设施施的的改改善善、医医疗疗水水平平的的提提高高以以及及人人类类文文明明的的不不断断发发展展,诸诸如如霍霍乱乱、天天花花等等曾曾经经肆肆虐虐全全球球的的传传染染性性疾疾病病已已经经得得到到有有效效的的控控制制。但但是是
39、一一些些新新的的、不不断断变变异异着着的的传传染染病病毒毒却却悄悄悄悄向向人人类类袭袭来来。20世世纪纪80年年代代十十分分险险恶恶的的爱爱滋滋病病毒毒开开始始肆肆虐虐全全球球,至至今今带带来来极极大大的的危危害害。长长期期以以来来,建建立立制制止止传传染染病病蔓蔓延延的的手手段段等等,一一直直是是各各国国有有关关专专家家和和官员关注的课题。官员关注的课题。第74页,共116页,编辑于2022年,星期六 不不同同类类型型传传染染病病的的传传播播过过程程有有其其各各自自不不同同的的特特点点,弄弄清清这这些些特特点点需需要要相相当当多多的的病病理理知知识识,这这里里不不可可能能从从医医学学的的角角
40、度度一一一一分分析析各各种种传传染染病病的的传传播播,而而只只是是按按照照一一般般的的传传播播模模型机理建立几种模型。型机理建立几种模型。模模型型1 在在这这个个最最简简单单的的模模型型中中,设设时时刻刻t的病人人数的病人人数x(t)是连续、可微函数,是连续、可微函数,第75页,共116页,编辑于2022年,星期六方程(方程(1)的解为)的解为 结结果果表表明明,随随着着t的的增增加加,病病人人人人数数x(t)无无限限增长,这显然是不符合实际的。增长,这显然是不符合实际的。建建模模失失败败的的原原因因在在于于:在在病病人人有有效效接接触触的的人人群群中中,有有健健康康人人也也有有病病人人,而而
41、其其中中只只有有健健康康人人才才可可以以被被传传染染为为病病人人,所所以以在在改改进进的的模型中必须区别这两种人。模型中必须区别这两种人。第76页,共116页,编辑于2022年,星期六 模型模型2(SI模型)模型)假设条件为假设条件为 1.在在疾疾病病传传播播期期内内所所考考察察地地区区的的总总人人数数N不不变变,即即不不考考虑虑生生死死,也也不不考考虑虑迁迁移移。人人群群分分为为易易感感染染者者(Susceptible)和和已已感感染染者者(Infective)两两类类(取取两两个个词词的的第第一一个个字字母母,称称之之为为SI模模型型),以以下下简简称称健健康康者者和和病病人人。时时刻刻t
42、这这两两类类人人在在总总人人数数中中所所占占比比例例分别记作分别记作s(t)和和i(t)。2.每每个个病病人人每每天天有有效效接接触触的的平平均均人人数数是是常常数数,称称为为日日接接触触率率。当当病病人人与与健健康康者者接接触触时,使健康者受感染变为病人。时,使健康者受感染变为病人。第77页,共116页,编辑于2022年,星期六又因为又因为 第78页,共116页,编辑于2022年,星期六方程(方程(5)是)是Logistic模型。它的解为模型。它的解为 第79页,共116页,编辑于2022年,星期六第80页,共116页,编辑于2022年,星期六这这时时病病人人增增加加的的最最快快,可可以以认
43、认为为是是医医院院的的门门诊诊量量最最大大的的一一天天,预预示示着着传传染染病病高高潮潮的的到到来来,是是医疗卫生部门关注的时刻。医疗卫生部门关注的时刻。其其原原因因是是模模型型中中没没有有考考虑虑到到病病人人可可以以治治愈愈,人人群群中中的的健健康康者者只只能能变变成成病病人人,病病人人不不会会再再变成健康者。变成健康者。第81页,共116页,编辑于2022年,星期六 模模型型3(SIS模模型型)有有些些传传染染病病如如伤伤风风、痢痢疾疾等等愈愈后后免免疫疫力力很很低低,可可以以假假定定无无免免疫疫性性,于于是是病病人人被被治治愈愈后后变变成成健健康康者者,健健康康者者还还可可以以被被感感染
44、染再再变变成成病病人人,所所以以这这个个模模型型称称SIS模型。模型。SIS模模型型的的假假设设条条件件1,2与与SI模模型型相相同同,增加的条件为增加的条件为 3每每天天被被治治愈愈的的病病人人数数占占病病人人总总数数的的比比例例为为常常数数,称称为为日日治治愈愈率率,病病人人治治愈愈后后成成为为仍仍可可被被感感染染的的健健康康者者。显显然然1/是是这这种种传传染病的平均传染期。染病的平均传染期。第82页,共116页,编辑于2022年,星期六 不不难难看看出出,考考虑虑到到假假设设3,SI模模型型的的(3)式应修正为式应修正为(4)式不变,于是()式不变,于是(5)式应改为)式应改为 我我们
45、们不不去去求求解解方方程程(9)(虽虽然然它它的的解解可可以以解解析析地地表表出出),而而是是通通过过图图形形分分析析i(t)的的变化规律。定义变化规律。定义 第83页,共116页,编辑于2022年,星期六注注意意到到和和1/的的含含义义,可可知知 是是整整个个传传染染期期内内每每个个病病人人有有效效接接触触的的平平均均人人数数,称称为为接接触数。触数。利用利用,方程(,方程(9)可以改写作)可以改写作 第84页,共116页,编辑于2022年,星期六第85页,共116页,编辑于2022年,星期六不难看出,接触数不难看出,接触数=1是一个阈值。是一个阈值。SI模模型型可可视视为为本本模模型型的的
46、特特例例,请请读读者者考考虑虑它它相当于本模型中相当于本模型中 或或 取何值的情况。取何值的情况。第86页,共116页,编辑于2022年,星期六 模模型型4(SIR模模型型)大大多多数数传传染染病病如如天天花花、流流感感、肝肝炎炎、麻麻疹疹等等治治愈愈后后均均有有很很强强的的免免疫疫力力,所所以以病病愈愈的的人人即即非非健健康康者者(易易感感染染者者),也也非非病病人人(已已感感染染者者),他他们们已已经经退退出出传传染染系系统统。这这种种情情况况比比较较复复杂杂,下下面面将将详详细细分分析析建建模模过程。过程。模型假设模型假设 1.总总人人数数N不不变变。人人群群分分为为健健康康者者、病病人
47、人和和病病愈愈免免疫疫的的移移出出者者(Removed)三三类类,称称SIR模模型型。三三类类人人在在总总数数N中中占占的的比比例例分分别别记记作作s(t),i(t)和和r(t)。第87页,共116页,编辑于2022年,星期六病病人人的的日日接接触触率率为为,日日治治愈愈率率为为(与与SI模模型型相同),传染期接触为相同),传染期接触为 =/。模型构成模型构成由假设由假设1显然有显然有 s(t)+i(t)+r(t)=1 (12)根根据据条条件件2方方程程(8)仍仍然然成成立立。对对于于病病愈愈免免疫疫的移出者而言有的移出者而言有第88页,共116页,编辑于2022年,星期六方方程程(14)无无
48、法法求求出出s(t)和和i(t)的的解解析析解解,我我们们先作数值计算。先作数值计算。第89页,共116页,编辑于2022年,星期六例例5 人口模型人口模型第90页,共116页,编辑于2022年,星期六1 问题的提出问题的提出 人口问题是当今世界上最令人关注的问题之一,一些人口问题是当今世界上最令人关注的问题之一,一些发展中国家的人口出生率过高,越来越威胁着人类的正常生发展中国家的人口出生率过高,越来越威胁着人类的正常生活,有些发达国家的自然增长率趋于零,甚至变为负数,造活,有些发达国家的自然增长率趋于零,甚至变为负数,造成劳动力紧缺,也是不容忽视的问题。另外,在科学技术和成劳动力紧缺,也是不
49、容忽视的问题。另外,在科学技术和生产力飞速发展的推动下,世界人口以空前的规模增长,统生产力飞速发展的推动下,世界人口以空前的规模增长,统计数据显示:计数据显示:年 1625 1830 1930 1960 1974 1987 1999 人口(亿)5 10 20 30 40 50 60 例例5 人口模型人口模型第91页,共116页,编辑于2022年,星期六 可以看出,人口每增长十亿的时间,由一百年缩可以看出,人口每增长十亿的时间,由一百年缩短为十二三年。我们赖以生存的地球,已经带着它的短为十二三年。我们赖以生存的地球,已经带着它的60亿子民踏入了亿子民踏入了21世纪。世纪。长期以来,人类的繁衍一直
50、在自发地进行着。只是由于长期以来,人类的繁衍一直在自发地进行着。只是由于人口数量的迅速膨胀和环境质量的急剧恶化,人们才猛然人口数量的迅速膨胀和环境质量的急剧恶化,人们才猛然醒悟,开始研究人类和自然的关系,人口数量的变化规律,醒悟,开始研究人类和自然的关系,人口数量的变化规律,以及如何进行人口控制等问题。以及如何进行人口控制等问题。第92页,共116页,编辑于2022年,星期六 我国是世界第一人口大国,地球上每九个人中就有我国是世界第一人口大国,地球上每九个人中就有二个中国人,在二个中国人,在20世纪的一段时间内我国人口的增长速世纪的一段时间内我国人口的增长速度过快,如下表:度过快,如下表:年