微分方程模型
在实际问题中,我们很难直接得出变量之间的数量关系,但是有时却很容易写出包括变量的导函数在内的一个方程,这就是微分方程,我们在一般的建模中常涉及常微分方程。
微分方程一般形式为:
0),...,'',',,()
(=n y
y y y x F 或),...,'',',,()
1()
(-=n n y
y y y x f y
。
若在某个范围内存在具有n 阶导数的函数)(x y ψ=使得
))(,...,')'(,)'(),(,()
(=x x x x x F n ψ
ψψψ,则称)(x y ψ=是微分方程的解。
微分方程所解决的问题通常可以分为两类:一类是用微分方程列出变量之间的关系式,求出位置函数的表达式,有时要借助软件进行数值分析;另一类是要了解未知变量或函数的某些性质即可,常需要根据微分方程的定性理论来研究,这两类建模问题我们将在后面进行讨论。
1. 微分方程简介
1.1. 简单的微分方程模型
一种比较简单的微分方程模型是变量的变化率与函数的即时值成正比,即
ky
y =',它的解就是kt e y t y 0)(=,这里0y 是初值,k 是待定常量。通常情况下,
如果0>k 称)(t y 指数递增;如果0 1.1.1. 放射性元素的自然衰变 放射性元素的自然衰变是化学上的一个基本事实,它常用于定碳测量,在考古学学上利用该方法测定古生物生存年代。存活于生物组织中占有确定比例的碳原子是放射性同位素14C ,一旦生物组织死亡,这种14C 不会增加,而会将一定比例的14C 衰变为12C ,并保持一定的速率(14C 的半衰期为5568年)按指数规律 下降。测定它现存的比例并与活的样品比较,就可以求得比例下降了多少,也就得到了被测样品的实际年代。 建立模型:假定)(t y 为t 时刻生物体内14C 原子的个数,经过相同的时间T ,y 的值减少为原值的1/2 (指数衰减)。002 1)(,)(y tT y e y t y kt ?= =,易得 T Ln k e T k 2 1,2 1= = ?。 如果给出一个具体时刻1t 的1y 值,由1011 )(y e y t y t k ==就可 以求出初值0y ,最后由kt e y t y 0)(=去指导实践。 1.1.2. 细菌数量的增长 细菌总数的增长也是与总数成正比的。如果培养体中的细菌总数在24小时之内由100(单位)增长导500(单位),问在10小时时细菌数量是多少? 建立模型:设t 时刻的细菌数目为)(t N ,细菌的增展率与总数成正比,马上得到如下微分方程模型: ????? ===500 )24(,100)0()()(N N t kN dt t dN 该方程的通解为kt e N t N 0)(=,k N 和0均为待定系数。 由500)24(,100)0(0===N N N 得 524 1,1000Ln k N ?= =, 即5 24 100)(Ln t e t N =,所以196)10(≈N ,可知10小时时细菌总数为196个单位。 1.1.3. 人口模型的演变 人口问题是当今世界上普遍关注的问题。人类在18世纪已经开始关注人口预报的工作,关于人口问题的建模很多,这里介绍两个较为简单的模型。 (1)Malthus 人口模型 马尔萨斯是英国人。1789年,他根据多年来的人口资料提出了一个人口模型,主要是假定单位时间内的人口增长数量与当时的人口数量称正比。 假定某地区t 时刻的人口数量为)(t N ,当考察一个国家或较大地域时)(t N 可 以假设为连续函数。记初始时刻的人口为0N ,人口增长率为r ,r 也可视为人口增长量与)(t N 的比例系数,通常假定r 为常数。 建立模型:在],[t t t ?+时间内,人口数量)(t N 的增量为: t t rN t N t t N ?=-?+)()()(, 整理得 ?????==0 )0() ()(N N t rN dt t dN (1-1) 解(1-1)式得rt e N t N 0)(=。该模型可以与19世纪前期欧洲一些地区的人口变化情况吻合的很好。但是随着时间的发展它和许多地区的实际人口增长之间差异越大,而且进入20世纪以来,人口的增长受到自然条件、环境条件、自身等问题的制约越来越显著,Malthus 人口模型中随时间发展,人口数量指数增长趋于无穷大也和实际不相符合。所以,有必要对Malthus 人口模型做相应的改进,这就是常见的Logistic 模型(阻滞增长模型)。 (2) Logistic 模型 根据人口增长的实际情况,人口增长率r 一般是随着当前人口数量的增加而减小的,做一个简单的假设)(t r 是)(t N 的线性函数:))(1()(m N t N r t r - =,其中r 是 固有的增长率,m N 是自然资源条件限制所能容纳的最大的人口数量;可以看到,当m N t N =)(时,人口的增长率0)(=t r 。从而得到如下阻滞增长模型: ?? ? ??=-=0)0()1() (N N N N N r dt t dN m (1-2) (2)式是非线性微分方程,利用分离变量或直接用Mathematica 软件可解得: rt m m e N N N t N --+= )1( 1)(0 根据(1-1)式,)('t N 实际上是)(t N 的二次函数,当2 m N N = 时,)('t N 取得 最大值。也即在m N N <<0上,)(t N 是递增的,当2 0m N N < <时,)(t N 增长速 度越来越大,当m m N N N <≤2 时增长速度越来越小;在m N N >时,人口是递减 的。 以美国的人口发展为例,1790年他美国的实际人口为3900=N 万,根据历史统计取21300,31.0==m N r 万,t e t N 31.0)1390 21300( 121300)(--+= ,m N 值是根据1800~ 1810美国的实际情况的一个预测值。 图 1 在图1中,)(t N 随着时间t 增长变化的趋势十分明显,函数的拐点出现在 10650 =N 附近。模型(1)一直到1930年左右仍然与实际统计数据吻合的很好, 但是以后就和实际有了较大差异。原因是m N 值的选取不尽合理,随着人类科学技术的进步,资源、环境等因素对人口自身的抑制作用相对减弱,所以一个国家的最大人口允许值不是一成不变的。 综上所述,我们给出微分方程模型建立的基本步骤: (1) 写出变量在(时间)单位内增量的变化值:这里需要把实际问题的具体信息转化为导数问题,变量可能是人口、速度、原子数、浓度、能力等等。一般问题都遵循“单位时间的净增量=单位时间的输入量-单位时间的输出量”,保证两端单位一致,然后再把等式转化为导数问题; (2) 边界条件:因为微分方程若有解,解常是一个函数族许多参量或系数都需要通过实际值回代才能得到。在实际问题中,初始条件和一些经验数据是非常重要的; (3)建立模型:审视(1)、(2),建立正确的表达式,写出 0),...,'',',,() (=n y y y y x F ; (4)求解和讨论:模型若有解析解,尽可能用准确巧妙的解法,若难以求解则尝试用软件进行求解或进行方程稳定性的讨论。 1.2. 微分方程稳定性的讨论 在利用微分方程进行数学建模活动中,有时建模的目的并非是为了寻求动态过程每一瞬时状态。我们往往需要了解系统某种稳定性的特征,特别是时间充分长以后动态过程的变化趋势,以及这种稳定状态是否容易收到破坏,这就需要进行微分方程平衡解与稳定性的讨论。 1.2.1. 简单微分方程的平衡点和稳定性 关于常微分方程平衡点和稳定性讨论,仅考察右端不显含有自变量t 的一阶微分方程 )(x f dt dx = (1-3) 称代数方程0)(=x f 的实根0x x =为原方程的平衡点或奇点,它也是(1-3)式的解。 在实际问题中,我们不仅要得到问题的解,还要得到∞→t 时解的变化趋势。如果从一定范围内的初始条件出发,(1-3)式的解)(t x 都满足)()(0∞→→t x t x ,则称平衡点0x 是稳定的。这种做法有时很难实现,我们经常利用如下判定0x 是否稳定的方法: 在0x 处将)(x f 作泰勒展开,只取线性部分得到(3)式的近似解: ))((00' x x x f dt dx -= (1-4) 显然0x 也是上述微分方程的平衡点,而(1-3)式的通解可写为 0)(0' )(x ce t x t x f +=, 可以看到: (1)当0)(0' 关于常微分方程组平衡点和稳定性讨论,仅考察右端不显含有自变量t 的微分方程组: ???? ?==),() ,(y x g dt dy y x f dt dx (1-5) 代数方程 ?? ?==0),(0 ),(y x g y x f 的实根 ?? ?==0 y y x x 称为上面微分方程组的平衡点,记 作),(000y x p ,也是原微分方程组的解。如果从一定范围内的初始条件出发,微分方程组的解))(),((t y t x 都满足00)(,)(y t y x t x t t →→∞ →∞ →,则称平衡点),(000y x p 是稳 定的。下面给出),(000y x p 是否稳定的判别准则。 设??? ?? ? ??????????????=??+??-=y p g x p g y p f x p f n y p g x p f m )() ()()(],)()([ 000000,则 (1)当00>>n m 且时,平衡点),(000y x p 是稳定的; (2)当00< 2. 传染病模型 本世纪初,瘟疫、AIDS 等传染病经常在某些地区流行。我们不可能做传染病传播的试验去获取数据,从医疗卫生部门那里取得的资料又是不充分、不完善的。但人们又十分关注这些问题:被传染的人数与哪些因素有关?如何预报传染病高潮的到来?为什么在同一地区一种传染病流行时,被传染的人数大致不变?我们 通过建立传染病传播的数学模型,可以较好的回答这些问题。 传染病的传播涉及因素很多,不可能通过一次简单的假设就能建立完善的数学模型。我们采取由简到繁,层层深入的方法,先做出简单假设,看结果是否与实际吻合,然后针对不合理不完善之处逐步修改和增加假设,得到最终满意的模型。 2.1. 模型1 基本假设: (1)每个病人在单位时间内传染的人数为常数; (2)某人得病之后经久不愈。 一些记号: )(t i :t 时刻的得病人数 i :刚开始时刻的得病人数 0k :一个病人在单位时间内能传染的人数(传染率) 建立模型:考察在单位时间t ?内,某疫区病人人数的变化情况,显然有 t t i k t i t t i ?=-?+)()()(0,整理得: ?????==0 0)0() () (i i t i k dt t di (2-1) 解(2-1)式得t k e i t i 0 0)(=。 模型分析:该结果表明,该疫区的病人人数将按照指数规律无限增加,即 ∞ →)(t i ,这显然和事实不相符合。实际上,一个地区的总人数大致可视为常数 (不考虑传染病传播期内出生和迁移的人数)。另外,在传染病的传播期,我们曾一个病人单位时间内能传染的人数0k ,0k 也不是一成不变的。在传染病发作的初期,0k 较大,随着病人数的增加,健康者逐渐减少,被传染的机会也减少了,于是0k 也就变小了,所以应对模型1加以修改,从而有模型2。 2.2. 模型2 基本假设: (1)疫区的总人数为一个定值; (2)每个病人在单位时间内传染的人数随着健康者人数的减少而减少; (3)某人得病之后经久不愈。 记号说明: )(t i :t 时刻的被传染的病人数,0i 为刚开始时刻的得病人数 )(t s :t 时刻的健康者人数 n :疫区的总人数,显然有)(t i +)(t s =n k :传染强度,假定单位时间内一个病人传染的人数与当时健康者人数 的成正比,比例系数记为常数k 建立模型:参考模型1中的(3。4)易得 ?? ???==0)0()()() (i i t i t ks dt t di (2-2) 把)()(t i n t s -=代入(7)式得下面微分方程, ?? ???=+=02)0() ()()(i i t kni t ki dt t di 解得nkt e i n n t i --+= )1( 1)(0 模型分析:首先利用)(t i 表达式,我们可以对传染病发作后的一段时间被传染者的数量变化做出预报,其次我们可以研究)(t i 的导函数, t dt t di ~)(被医学界 称为传染病曲线,它反映了病人数变化率与时间的关系,可以用来预报一些恶性传染病前期传染高峰的到来时间。 由nkt e i n n t i --+= )1( 1)(0 ,得 2 2 ] )1( 1[)1()(nkt nkt e i n e i n kn dt t di ---+-=. 取0003.0,100,100000===k i n ,做出t dt t di t t i ~)(~)(和 的图像如图2、图3。 图2 图3 令 ] )1( 1[)1(] )1( 1[)1(2)(2 3 2 3 220 3 2 2 =-+---+-= ----knt knt nkt nkt e i n e i n n k e i n e i n n k dt t di ,解得kn i n Ln t ) 1( 0*-=,* t 时刻附近就是传染病发作的最高峰期,同时2 )(*n t i = . 由此可见,当传染病的传染强度k 或n 增加时,*t 将变小,即传染病的高峰来得快些,这与实际情况是吻合的,这里的k 可由经验或统计数据求得。同时,我们也看到,随着∞→t ,∞→)(t i ,这意味着最终疫区内人人将被传染,这又与实际情况不相符合,这主要是因为疫区内人群的划分不尽合理,仍需改进。 2.3. 模型3 在模型2的基础上,我们将某一疫区所有的人分为三类: 第一类是传染类,记为)(t i ,它由能把疾病传染给他人的人组成; 第二类是受传染类,记为)(t s ,它由并非传染者但能够成为传染者的人组成; 第三类是疫区中排除了第一类和第二类后余下的人组成,记为)(t r ,它包括患病死去的人,病愈后具有长期免疫力的人以及在病愈且并未出现长期免疫能力之前被隔离的人等组成。 基本假设: (1)总人数不变,)(t i +)(t s +)(t r =n (2)传染强度:单位时间内一个病人能传染的人数与当时的健康者人数成正比,比例系数为k (3)恢复系数:单位时间内病愈且具有免疫力的人数与当时的病人人数成正比,比例系数为l 建立模型: 由基本假设(1)、(2)、(3),可将模型(7)式做如下改进, ?????===-=0 0)0(,0)0(,)0()()()() (s s r i i dt t dr t i t ks dt t di (2-3) 由假设(3)得 )()(t li dt t dr =,且)(t i =n -)(t s -)(t r 代入(2-3)式可得 ) ()() (1 )(t d t dr t s dt t ds ρ - =,其中k l =ρ, 直接用Mathematica 软件可解得 )](21)(1 1[)(2 2 0) (1 0t r t r s e s t s t r ρ ρ ρ + - ==- (2-4) 所以 )) ()(()()(t s t r n l t li dt t dr --?== )](2)()1()[(2 2 00 0t r s t r s s n l ρ ρ - -+-?= 解得: ) 2tan([21)(' c qt q b c t r +- ?-= (2-5) 其中)(0s n l a -=,)1( -=ρ s l b ,ac b q 42 +=,q b c = 'tan 模型分析:注意到在(2-5)式中,00)0(i s n i =-=, 故当0i 很小时,2 2 02 02) 1( 24-<=ρ ρ s l i s l ac , 所以b q ≈。 又因为 1)2 tan(lim ' -=+- ∞ →c t q t 所以,)1( 2)(21)(lim 0 2 -= ≈ += ∞ →ρ ρs s c b q b c t r t (2-6) (11)式表明,当0s 很小时(ρ≤0s )时,传染病根本传染不开,只有0s 很大时才会有传染病疫情。假定ερ+=0s (0≥ε,只有如此才能有传染问题)且ρε<<,则有 ε ρ εε ρρ 22)(lim 2 =? +≈ ∞ →t r t (2-7) 在这里我们看到,一个疫区最终免疫的人数,也即最终痊愈的人数和该传染病最终传染的人数是相等的,于是由(12)式得到: (1)对于同一地区,同一传染病所传染的人数大致是一个常数为2ε,这个与统计结果一致; (2)当恢复系数增大l ,传染强度k 减小时,k l =ρ变大,相应的ε变小, 从而传染病最终传染的总人数2ε也变小; (3)如果0s 相对于n 非常大,则0s 可视为总人数,则传染病最终传染的总人数为)(2ρ-n 。 在有的参考书中,把l k == ρ ρ1'叫做接触数,它表示一个传染期内每个病人 有效接触到的平均人数。为了降低传染的比例,通常需要我们不断提高卫生和医疗的水平,尽可能的降低'ρ,从而增大ρ,从而使得最终传染的总人数)(2ρ-n 不断减小。 另外一个降低传染病传染比例的方法是尽可能的降低0s ,这个可以通过预防、接种等方法提高群体的免疫力做到。但有实际操作的困难,因为免疫者一般不是均匀分布于疫区的,所以我们的做法存在有待改进的地方。本文的做法旨在让读者学会一种由简到繁,循序渐进、不断完善的数学建模的方法,把这种思想用到日后的建模实际中去。 3. 血液中的酒精浓度模型 酒后驾车的危害性,已经引起交通部门乃至全社会的高度关注。国家质量监督检查检疫局在2004年5月31日颁发了新的《车辆驾驶人员血液、呼气酒精含 量阀值与检验》国家标准,对酒后驾车人员血液中所允许的酒精浓度做了具体规定。 问题:对于一个驾驶员,他饮酒后多长时间不能开车?对于一定量的酒,在短时间饮完好还是在较长时间饮完比较好?本文在一定假设基础上建立数学建模,从量上给出了答案。 人饮完酒后,酒精通过胃肠的吸收扩散到人的体液(包括血液)中,同时体液中的酒精又通过汗液、尿液等排到体外。实际上,这种吸收、扩散、排除的过程相对复杂,干扰因素繁多。这里,我们把人体设想成一个含有两个房室的模型,只有胃肠道和体液两部分组成,并作如下假定。 模型假设: (1)酒精浓度: 一个基本的假设是酒精进入体液后迅速扩散到全身各个部位,任一时刻血液和体液中的酒精浓度相等,记为)(t p 。另外,可以认为人体体液的总体积是一个常数V ,且酒后人体的体液和酒精是正常排除的,没有意外呕吐等现象。 (2)吸收因子: 胃肠中的酒精被吸收到体液中的速率与胃肠道中的酒精的质量成正比,即若t 时刻胃肠中的酒精质量为)(t y ,那么此时的吸收速率为 )(1t y k ,其中1k 叫做吸收因子。一般情况下,1k 收到胃肠蠕动速度、体液 PH 值 等因素影响,但在一定时间段里1k 为常量。 (3)消除因子: 体液排出的速度受到诸多因素干扰,比如温度、运动量等等,但这里假定体液排处体外的速度是匀速的,并且单位时间内排出体外的体液体积为2k ,令V k k 23 ,称3k 为消除因子,它表示单位体积的体液在单位时间内 排除体外的量。 体液中酒精含量一方面通过胃肠吸收而得,另一方面又随体液排出体外。所以,体液中酒精得浓度与吸收因子1k 、消除因子3k 以及饮酒的量有密切联系,随着时间变化而发生变化。同时,我们看到,一定量的酒精进入体内有多种时间方式,譬如,可以在很短时间进入体内(短时饮酒模型)、在较长时间内进入(长时饮酒模型)或每隔一段时间分批量进入体内(间断饮酒模型)等。本节主要讨 论短时饮酒模型和长时饮酒模型。 3.1. 模型1 短时饮酒模型 在这个模型的实际问题中,饮酒者在很短的时间内(近似瞬间)M 毫升的酒,按照酒的度数计算,假设喝入的酒精质量为m (单位:mg ),相当于初始时刻体内酒精含量,那么计算t 时刻胃肠中的酒精质量)(t y 可以化为如下微分方程: ?? ???=-=m y t y k dt t dy )0() () (1 (3-1) 由分离变量解得 t k me t y 1)(-= (3-2) 建立酒精浓度模型: 在时间段],0[t 上考虑血液内酒精浓度的变化,由假设(2),在t ?内吸收到体液中的酒精质量约为t t y k ?)(1,因此在],0[t 时间内吸收到体液中的酒精质量为 ?t dt t y k 01)(。 假设在t 时刻体液中的酒精浓度为)(t p (单位:mg/100ml )。由假设(3),在 ],0[t 时间内排出体外的酒精质量为 ? t dt t p k 0 2)(100 ,从而t 时刻体液内酒精质量为 ? ?- t t dt t p k dt t y k 0 20 1)(100 )(, 故t 时刻体液中的酒精浓度为 100 )(100 )(00 21?- ?? V dt t p k dt t y k t t , 从而有 )(t p =V dt t p k dt t y k t t ?? -00 2 1)()(100 对上式两端同时求导可得: )(100)()(1 2' t y V k t p V k t p = + (3-3) 注意到V k k 23= ,令V E 100= ,从而(14)式化为 ? ? ?=+-00()()(113')=p me Ek t p k t p t k (3-4) 解(3-4)式得) 1 31 21()(t k t k e e k k Emk t p ----= (3-5) 模型分析: (1) 浓度分析 对(3-5)式求导可得在1 31 3 *k k Lnk Lnk t --= 处,0)(*='t p , 所以在区间], 0[131 3k k Lnk Lnk --上)(t p 单调递增, 在],[ 1 31 3∞--k k Lnk Lnk 上) (t p 单调递减。 这说明一开始血液酒精浓度以较快的速度增长,在1 31 3*k k Lnk Lnk t --= 时刻浓度 达到最大值,之后又逐渐降低,随着时间的推移,体液中酒精的浓度越练越低,直到完全消除为止。 (2) 数据拟合 为了检验我们建立模型的合理性,我们选取以下一组数据针对(3-5)中的1 k 和3k 进行数据拟和: 取)(5300mg m =,相当于两瓶酒精浓度为4。2g/100ml 的啤酒中酒精的含量;)(49000ml V =,相当于一个体重为70kg 的人的体液含量; t 与)(t p 的一组对应值表如下(选自 2004年大学生数学建模竞赛C 题数据)表1 t 0.25 0.5 0.75 1 1.5 2 2.5 3 3.5 4 4.5 5 P(t) 30 68 75 82 82 77 68 68 58 51 50 41 t 6 7 8 9 10 11 12 13 14 15 16 P(t) 38 35 28 25 18 15 12 10 7 7 4 表1 利用Matlab 软件下的非线性拟合命令nlinfit 对表中的数据进行拟和,得到对应于这组数据的吸收因子和消除因子199.0,98.131==k k ,对照这组数据的散点图,我们看到(3-3)、(3-4)式建立的模型基本符合人体血液中酒精浓度变化的规律,验证了我们所建模型的合理性。 但是,应该看到这里的)(t p 为光滑函数,说明我们建立的模型具有理想化的特点,图4是表1的散点图合)(t p 的图像,黑色的是利用软件拟合的曲线,灰色的是连线后散点图,可以看到两者吻合的很好,在快速饮酒的基础上我们再建立个人长时间的饮酒模型。 图4 两种类型的数据拟合曲线 3.2. 模型2 长时间饮酒模型 设某人在较长时间0T 内,摄入的酒精质量为m ,我们在这里假设这种摄入是匀速进行的,即在0T 时间内酒精以 T m 的速度进入胃肠道。 基本假设: 当0T t ≤时胃肠道中酒精质量为)(1t y ,体液中酒精浓度为)(1t p ,则)(1t y 的变化率为 )(110 t y k T m -,从而假设可以转化为如下初值问题, ?? ? ??=-=0)0()()(1110 1y t y k T m dt t dy (0T t ≤) (3-6) 解之得)1()(11 01t k e k T m t y --= (3-7) 建立模型: 当0T t ≤时, 我们首先来建立)(1t p 模型,类似于前面分析在],0[t 时间内排出体外的酒精质量为 ? t dt t p k 0 12)(100 , 从而t 时刻体液内酒精质量为? ?- t t dt t p k dt t y k 0 120 11)(100 )(, 故t 时刻体液中的酒精浓度为 100 )(100 )(0 1211?- ?? V dt t p k dt t y k t t , 从而有 )(1t p = V dt t p k dt t y k t t ?? -00 12 11)()(100, 整理得)1()()(10 13' 1t k e T Em t p k t p --= + (3-8) (3-8)式中初始条件为0)0(1=p ,解得(3-8)式的解为: 3 11 313 1 31113)11( )(k A e k k A e k k k A t p t k t k + -- --=--,其中0 1T Em A = 当0T t >时,设胃肠道中的酒精质量为)(2t y ,血液中酒精浓度为)(2t p ,则)(2t y 和)(2t p 的变化规律类似于(3-1)、(3-3)式,只是初值不同,直接写出如下方程: ???? ?===+-=) ()(),()()()()()()(0102010 22123' 221' 2T p T p T y T y t y Ek t p k t p t y k t y (3-9) 该一阶线性系统的解为) (1) (10120103 ))(()(T t k T t k e B e B T p t p ----+-= (3-10) 其中1 30111)(k k T y Ek B -= 。 这样在整个过程中,血液中的酒精浓度)(t p 的方程为 ???>≤=0 201),(),()(T t t p T t t p t p (3-11) 3.3. 模型1和2综合分析: (1) 浓度指标分析 我们选取如下一组参数进行实例分析:199.0,98.131==k k ,小时23=T , )(79500mg m =,)(49000ml V =,则181 .0,122.104)2(,576.79)2(111-==≈B p y , 利用软件计算出两中模型比较之下不能驾车的时间范围,时间长度,最高酒精浓度等指标见表2 最大浓度 最大浓度时 间 禁止时间范围 禁止时间长 度 短时饮酒 125.51 1.2901 [0.07,10.95] 10.88 长时饮酒 118.44 2.6145 [0.57,11.97] 11.40 表2 (2) 饮酒的个体差异 图5 短时和长时间饮酒的浓度比较 利用数学软件画出短时饮酒和长时饮酒的酒精浓度的函数图像如上图5。由图知,短时饮酒的酒精浓度最大值比长时饮酒的酒精浓度最大值要大。另外,对于短时饮酒来说在很短时间内,大约在1.3小时,血液中酒精浓度最大;对于长 时饮酒来说,这个时间要晚些,大约在2.6小时,血液中酒精浓度最大,之后又以较慢的速度降低,因此禁止驾驶的时间较长。 显然,喝酒的时间长短不同,其峰值浓度达到有明显的区别。对于酒量偏低的饮酒者更适于长时间的饮酒。对于需要在酒后能够快速的醒酒以便工作的饮酒者来说,那么不如在短时间内以较快的速度喝酒。 对于间断喝酒模型的讨论请参考文献。 4. 捕鱼业的收获模型 渔业资源是一种再生资源,再生资源应注意适度开发,不能完全为了的经济利益而“竭泽而渔”,渔业生产应该在持续稳产的前提下追求产量和经济效益。 实际问题 考察一个渔场,其中鱼量在天然环境下按一定的规律增长。如果捕捞量恰好等于增长量,那么鱼量将保持不,这个捕捞量就可以持续。本节要建立在捕捞情况下,渔场鱼量遵循的方程,分析鱼量稳定的条件,并且在稳定的前提下讨论如何控制捕捞使持续产量或经济效益达到最大,并研究过度捕捞的问题。 建模与求解 为建立模型,做如下假设和记号: ~ )(t x t 时刻渔场中的鱼量; ~ m x 渔场当前条件下所能允许养殖的最大鱼量; ~ )(x r 渔场鱼量的增长率,它通常是当前鱼量的函数,通常随着)(t x 的增大 而减少,这里假定) 1()(m x x r x r - =,r 为鱼量的固有增长率,m x t x =)(时, 0)(=x r 。 类似于人口增长的阻滞增长模型,得到自然增长条件下渔场鱼量变 化率的表达式, ) 1(m x x rx dt dx - = (4-1) 假定单位时间内的捕鱼量与渔场的鱼量)(t x 成正比,捕捞率为k ,则得到有 捕捞情况下渔场鱼量变化率的表达式, kx x x rx dt dx m -- =)1( (4-2) 根据微分方程平衡点的讨论方法,令0)1()(=--=kx x x rx x f m , 得到(4-2)式的平衡解 m x r k r x -= 0 (4-3) 则r k k x x r r x f x x m -=-?- ==0 2)(0' (4-4) 所以得到如下结论: (1)当r k <时,0)(0'<-=r k x f ,0x 是(4-2)式的稳定的平衡解,所以 r k <是渔业生产必须遵守的规则。 (2)当r k >时,0)(0'>-=r k x f ,0x 是(4-2)式的不稳定的平衡解,所以r k >时就必须改变捕捞方法,以保持鱼量稳定。 解的讨论 下面我们利用图解法讨论在保持鱼量稳定的前提下,如何选择捕捞系数k ,使捕捞量最大。考察rx x x r x x rx x g m m +- =- =2 )1()(和kx x h =)( 显然,)(x g 为抛物线在x 轴上方部分,只有)(x g 和)(x h 相交时,微分方程(4-2)式才会有平衡点0x 。 而r x x r x g m +-=2)(',显然r g =)0(',即)(x g 在原 点处切线斜率为r 。 只有当r k <<0时,)(x h 与)(x g 才会有交点。观察到,当 kx x h =)(过抛物 线顶点0P 时,渔场将获得最大的捕捞量m y ,此时稳定点 4 )1(,2 00m m m m rx x x rx y x x = - == 。 所以,当控制捕捞率在2 r k = 或控制捕捞率使渔场鱼量在 2 m x 时,就可以使渔 场鱼量稳定的条件下使捕捞量最大,最大捕捞量4 m m rx y = 效益分析 若我们希望在保持渔场鱼量稳定条件下得到利润函数L 最大,可以重新建模。设鱼单价为p ,且捕捞成本鱼捕捞率成正比,比率系数为c ,则在保持渔场鱼量稳定的条件下(r k <)单位时间内的利润为ck pkx L -=0。 由(4-3)式得)1(0m x x r k - =,代入 L 得到 )1)(()(000m x x c px r x L - -= (4-5) 令0)(0'=x L ,得p c x x m 22 0+ = 此时捕捞量0kx y =m m m x p rc rx x x x r 2 2 0044)1(- = - = 捕捞量) 1(2 )1(0m m px c r x x r k -=- = 在效益分析中,最优捕鱼率比最大捕鱼量少m x p rc 22 4,它与成本的平方成正比, 与鱼价的平方成反比。 捕捞过度 上面的分析是以独家捕捞为基础的,即捕捞行为是理性的。如果渔场开放,存在多家捕捞企业,每个经营者既不能控制价格,也不能控制捕捞量,这种情况下会导致捕捞过度。 实验07讲评、参考答案 讲评未按时交的同学 批改情况: 附参考答案: 实验07 微分方程模型(2学时) (第5章 微分方程模型) 1.(验证)传染病模型2(SI 模型)p136~138 传染病模型2(SI 模型): 0(1),(0)di k i i i i dt =-= 其中, i (t )是第t 天病人在总人数中所占的比例。 k 是每个病人每天有效接触的平均人数(日接触率)。 i 0是初始时刻(t =0)病人的比例。 1.1 画~di i dt 曲线图p136~138 取k =0.1,画出 i dt di ~的曲线图,求i 为何值时dt di 达到最大值,并在曲线图上标注。 %传染病模型2(SI 模型)的di/dt~i 曲线图 %文件名:p137fig2.m %λ=0.1 clear; clc; 《数学建模实验》 王平 提示:fplot, fminbnd, plot, text, title, xlabel 1)画曲线图 用fplot函数,调用格式如下: fplot(fun,lims) fun必须为一个M文件的函数名或对变量x的可执行字符串。 若lims取[xmin xmax],则x轴被限制在此区间上。 若lims取[xmin xmax ymin ymax],则y轴也被限制。 本题可用 fplot('0.1*x*(1-x)',[0 1.1 0 0.03]); 2)求最大值 用求解边界约束条件下的非线性最小化函数fminbnd,调用格式如下: x=fminbnd('fun',x1,x2) fun必须为一个M文件的函数名或对变量x的可执行字符串。 返回自变量x在区间x1 第五章微分方程建模案例 微分方程作为数学科学的中心学科,已经有三百多年的发展历史,其解法和理论已日臻完善,可以为分析和求得方程的解(或数值解)提供足够的方法,使得微分方程模型具有极大的普遍性、有效性和非常丰富的数学涵。微分方程建模包括常微分方程建模、偏微分方程建模、差分方程建模及其各种类型的方程组建模。微分方程建模对于许多实际问题的解决是一种极有效的数学手段,对于现实世界的变化,人们关注的往往是其变化速度、加速度以及所处位置随时间的发展规律,其规律一般可以用微分方程或方程组表示,微分方程建模适用的领域比较广,涉及到生活中的诸多行业,其中的连续模型适用于常微分方程和偏微分方程及其方程组建模,离散模型适用于差分方程及其方程组建模。本章主要介绍几个简单的用微分方程建立的模型,让读者一窥方程的应用。下面简要介绍利用方程知识建立数学模型的几种方法: 1.利用题目本身给出的或隐含的等量关系建立微分方程模型 这就需要我们仔细分析题目,明确题意,找出其中的等量关系,建立数学模 型。 例如在光学里面,旋转抛物面能将放在焦点处的光源经镜面反射后成为平行光线,为了证明具有这一性质的曲线只有抛物线,我们就是利用了题目中隐含的条件——入射角等于反射角来建立微分方程模型的。 2.从一些已知的基本定律或基本公式出发建立微分方程模型 我们要熟悉一些常用的基本定律、基本公式。例如从几何观点看,曲线 y y(x)上某点的切线斜率即函数y y(x)在该点的导数;力学中的牛顿第二运 动定律:F ma ,其中加速度a 就是位移对时间的二阶导数,也是速度对时间 的一阶导数等等。从这些知识出发我们可以建立相应的微分方程模型。 例如在动力学中,如何保证高空跳伞者的安全问题。对于高空下落的物体, 我们可以利用牛顿第二运动定律建立其微分方程模型, 设物体质量为m ,空气阻 力 系数为k ,在速度不太大的情况下,空气阻力近似与速度的平方成正比;设时 刻t 时物体的下落速度为v ,初始条件:v (o ) 0.由牛顿第二运动定律建立其微 分方程模型: 求解模型可得: 体在地面上的投影面积。根据极限速度求解式子,在m,, 一定时,要求落地速 度w 不是很大时,我们可以确定出s 来,从而设计出保证跳伞者安全的降落伞的 直径大小来 3?利用导数的定义建立微分方程模型 dv m 一 dt mg kv 2 ? k(exp[2t 由上式可知,当t 其中,阻力系数k 1) 时,物体具有极限速度: lim v t mg :k , s , 为与物体形状有关的常数, 为介质密度,s 为物 、mg(exp[2t 1) 第二章:动力学系统的微分方程模型 利用计算机进行仿真时,一般情况下要给出系统的数学模型,因此有必要掌握一定的建立数学模型的方法。在动力学系统中,大多数情况下可以使用微分方程来表示系统的动态特性,也可以通过微分方程可以将原来的系统简化为状态方程或者差分方程模型等。在这一章中,重点介绍建系统动态问题的微分方程的基本理论和方法。 在实际工程中,一般把系统分为两种类型,一是连续系统;其数学模型一般是高阶微分方程;另一种是离散系统,它的数学模型是差分方程。 §2.1 动力学系统统基本元件 任何机械系统都是由机械元件组成的,在机械系统中有3种类型的基本机械元件:惯性元件、弹性元件和阻尼元件。 1 惯性元件:惯性元件是指具有质量或转动惯量的元件,惯量可以定义为使加速度(或角加速度)产生单位变化所需要的力(或力矩)。 惯量(质量)= ) 加速度(力(2 /) s m N 惯量(转动惯量)= ) 角加速度(力矩(2/) s rad m N ? 2 弹性元件:它在外力或外力偶作用下可以产生变形的元件,这种元件可以通过外力做功来储存能量。按变形性质可以分为线性元件和非线性元件,通常等效成一弹簧来表示。 对于线性弹簧元件,弹簧中所受到的力与位移成正比,比例常数为弹簧刚度k 。 x k F ?= 这里k 称为弹簧刚度,x ?是弹簧相对于原长的变形量,弹性力的方向总是指向弹 簧的原长位移,出了弹簧和受力之间是线性关系以外,还有所谓硬弹簧和软弹簧,它们的受力和弹簧变形之间的关系是一非线性关系。 3 阻尼元件:这种元件是以吸收能量以其它形式消耗能量,而不储存能量,可以形象的表示为一个活塞在一个充满流体介质的油缸中运动。阻尼力通常表示为: α x c R = 阻尼力的方向总是速度方向相反。当1=α,为线性阻尼模型。否则为非线性阻 尼模型。应注意当α等于偶数情况时,要将阻尼力表示为: ||1--=αx x c R 这里的“-”表示与速度方向相反 第3章微分方程模型 3.1 微分方程模型的建模步骤 在自然科学以及工程、经济、医学、体育、生物、社会等学科中的许多系统,有时很难找到该系统有关变量之间的直接关系——函数表达式,但却容易找到这些变量和它们的微小增量或变化率之间的关系式,这时往往采用微分关系式来描述该系统——即建立微分方程模型。我们以一个例子来说明建立微分方程模型的基本步骤。 例1 某人的食量是10467(焦/天),其中5038(焦/天)用于基本的新陈代谢(即自动消耗)。在健身训练中,他所消耗的热量大约是69(焦/公斤?天)乘以他的体重(公斤)。假设以脂肪形式贮藏的热量100%地有效,而1公斤脂肪含热量41868(焦)。试研究此人的体重随时间变化的规律。 模型分析 在问题中并未出现“变化率”、“导数”这样的关键词,但要寻找的是体重(记为W )关于时间t 的 函数。如果我们把体重W 看作是时间t 的连续可微函数,我们就能找到一个含有的dt dW 微分方程。 模型假设 1.以)(t W 表示t 时刻某人的体重,并设一天开始时人的体重为0W 。 2.体重的变化是一个渐变的过程。因此可认为 )(t W 是关于t 连续而且充分光滑的。 3.体重的变化等于输入与输出之差,其中输入是指扣除了基本新陈代谢之后的净食量吸收;输出就是进行健身训练时的消耗。 模型建立 问题中所涉及的时间仅仅是“每天”,由此,对于“每天” 体重的变化=输入-输出。 由于考虑的是体重随时间的变化情况,因此,可得 体重的变化/天=输入/天—输出/天。 代入具体的数值,得 输入/天 = 10467(焦/天)—5038(焦/天)=5429(焦/天), 输出/天 = 69(焦/公斤?天)×W (公斤)= 69W (焦/天)。 体重的变化/天=t W ??(公斤/天)dt dW t =→?0 考虑单位的匹配,利用 “公斤/天=公斤焦天 焦/41868 /”, 可建立如下微分方程模型 第一章 绪论 §1、1 微分方程:某些物理过程的数学模型 §1、2 基本概念 习题1、2 1.指出下面微分方程的阶数,并回答方程就是否线性的: (1) y x dx dy -=24; (2)0122 2 2=+??? ??-xy dx dy dx y d ; (3)0322 =-+?? ? ??y dx dy x dx dy ; (4)x xy dx dy dx y d x sin 352 2=+-; (5) 02cos =++x y dx dy ; (6)x e dx y d y =+??? ? ??22sin . 解 (1)一阶线性微分方程; (2)二阶非线性微分方程; (3)一阶非线性微分方程; (4)二阶线性微分方程; (5)一阶非线性微分方程; (6)二阶非线性微分方程. 2.试验证下面函数均为方程022 2=+y dx y d ω的解,这里0>ω就是常数. (1)x y ωcos =; (2)11(cos C x C y ω=就是任意常数); (3)x y ωsin =; (4)22(sin C x C y ω=就是任意常数); (5)2121,(sin cos C C x C x C y ωω+=就是任意常数); (6)B A B x A y ,()sin(+=ω就是任意常数). 解 (1)y x dx y d x dx dy 2 222cos ,sin ωωωωω-=-=-=,所以022 2=+y dx y d ω,故 x y ωcos =为方程的解. (2)y x C y x C y 2 2 11cos , sin ωωωωω-=-=''-=',所以022 2=+y dx y d ω,故x C y ωcos 1=为方程的解. (3)y x dx y d x dx dy 2222sin ,cos ωωωωω-=-==,所以02 2 2=+y dx y d ω,故x y ωsin =为方程的解. (4)y x C y x C y 2 2 22sin , cos ωωωωω-=-=''=',所以022 2=+y dx y d ω,故x C y ωsin 2=为方程的解. (5)y x C x C y x C x C y 2222121sin cos , cos sin ωωωωωωωωω-=--=''+-=',所 以02 2 2=+y dx y d ω,故x C x C y ωωsin cos 21+=为方程的解. (6)y B x A y B x A y 2 2 )sin(, )cos(ωωωωω-=+-=''+=',故0222=+y dx y d ω,因 此)sin(B x A y +=ω为方程的解. 3.验证下列各函数就是相应微分方程的解: (1)x x y sin = ,x y y x cos =+'; (2)212x C y -+=,x xy y x 2)1(2 =+'-(C 就是任意常数); (3)x Ce y =,02=+'-''y y y (C 就是任意常数); (4)x e y =,x x x e ye y e y 2212-=-+'-; (5)x y sin =,0cos sin sin 22 2 =-+-+'x x x y y y ; (6)x y 1- =,12 22++='xy y x y x ; (7)12 +=x y ,x y x y y 2)1(2 2 ++-='; (8))()(x f x g y = ,) () ()()(2x f x g y x g x f y '-'='. 第一章 绪论 §1.1 微分方程:某些物理过程的数学模型 §1.2 基本概念 习题1.2 1.指出下面微分方程的阶数,并回答方程是否线性的: (1) y x dx dy -=24; (2)0122 2 2=+??? ??-xy dx dy dx y d ; (3)0322 =-+? ? ? ??y dx dy x dx dy ; (4)x xy dx dy dx y d x sin 352 2=+-; (5) 02cos =++x y dx dy ; (6)x e dx y d y =+??? ? ??22sin . 解 (1)一阶线性微分方程; (2)二阶非线性微分方程; (3)一阶非线性微分方程; (4)二阶线性微分方程; (5)一阶非线性微分方程; (6)二阶非线性微分方程. 2.试验证下面函数均为方程02 2 2=+y dx y d ω的解,这里0>ω是常数. (1)x y ωcos =; (2)11(cos C x C y ω=是任意常数); (3)x y ωsin =; (4)22(sin C x C y ω=是任意常数); (5)2121,(sin cos C C x C x C y ωω+=是任意常数); (6)B A B x A y ,()sin(+=ω是任意常数). 解 (1)y x dx y d x dx dy 2222cos ,sin ωωωωω-=-=-=,所以02 2 2=+y dx y d ω,故x y ωcos =为方程的解. (2)y x C y x C y 2 2 11cos , sin ωωωωω-=-=''-=',所以0222=+y dx y d ω,故 x C y ωcos 1=为方程的解. (3)y x dx y d x dx dy 2 222sin ,cos ωωωωω-=-==,所以022 2=+y dx y d ω,故x y ωsin =为方程的解. (4)y x C y x C y 2 2 22sin , cos ωωωωω-=-=''=',所以022 2=+y dx y d ω,故x C y ωsin 2=为方程的解. (5)y x C x C y x C x C y 2222121sin cos , cos sin ωωωωωωωωω-=--=''+-=', 所以022 2=+y dx y d ω,故x C x C y ωωsin cos 21+=为方程的解. (6)y B x A y B x A y 2 2 )sin(, )cos(ωωωωω-=+-=''+=',故02 2 2=+y dx y d ω,因此)sin(B x A y +=ω为方程的解. 3.验证下列各函数是相应微分方程的解: (1)x x y sin = ,x y y x cos =+'; (2)212x C y -+=,x xy y x 2)1(2 =+'-(C 是任意常数); (3)x Ce y =,02=+'-''y y y (C 是任意常数); (4)x e y =,x x x e ye y e y 2212-=-+'-; (5)x y sin =,0cos sin sin 22 2 =-+-+'x x x y y y ; (6)x y 1- =,12 22++='xy y x y x ; (7)12 +=x y ,x y x y y 2)1(2 2 ++-='; 微分方程模型建模实例 1.一个半球状雪堆,其体积融化的速率与半球面面积S成正比,比例系数k > 0。设融化中雪堆始终保持半球状,初始半径为R且3小时中融化了总体积的7/8,问雪堆全部融化还需要多长时间? 2.从致冰厂购买了一块立方体的冰块,在运输途中发现,第一小时大约融化了1/4 (1)求冰块全部融化要多长时间(设气温不变) (2)如运输时间需要2.5小时,问:运输途中冰块大约会融化掉多少? 3.一展开角为α的圆锥形漏斗内盛着高度为H的水,设漏斗底部的孔足够大(表面张力不计),试求漏斗中的水流光需要多少时间? 4.容器甲的温度为60度,将其内的温度计移入容器乙内,设十分钟后温度计读数为70度,又过十分钟后温度计读数为76度,试求容器乙内的温度。 5.一块加过热的金属块初始时比室温高70度,20分钟测得它比室温高60度,问:(1)2小时后金属块比室温高多少?(2)多少时间后,金属块比室温高10度? 6.设初始时容器里盛放着含净盐10千克的盐水100升,现对其以每分钟3升的速率注入清水,容器内装有搅拌器能将溶液迅时搅拌均匀,并同时以每分钟2升的速率放出盐水,求1小时后容器里的盐水中还含有多少净盐? 7.某伞降兵跳伞时的总质量为100公斤(含武器装备),降落伞张开前的空气阻力为0.5v,该伞降兵的初始下落速度为0,经8秒钟后降落伞打开,降落 伞打开后的空气阻力约为0.6 试球给伞降兵下落的速度v(t),并求其下落的极限速度。 8. 1988年8月5日英国人Mike McCarthy创建了一项最低开伞的跳伞纪录,它从比萨斜塔上跳下,到离地179英尺时才打开降落伞,试求他落地时的速度。 9.证明对数螺线r=A 上任一处的切线与极径的夹角的正切为一常 数,() 数学建模学习辅导 第三章 微分方程模型 本章重点: 车间空气清洁问题、减肥问题、单种群增长问题与多物种相互作用问题等数学模型的建立过程与所使用的方法 复习要求: 1.进一步理解建模基本方法与基本建模过程,掌握平衡原理与微元法在建模中的用法. 所谓平衡原理是指自然界的任何物质在其变化的过程中一定受到某种平衡关系的支配.注意发掘实际问题中的平衡原理是从物质运动机理的角度组建数学模型的一个关键问题.就象中学的数学应用题中等量关系的发现是建立方程的关键一样. 微元法是指在组建对象随着时间或空间连续变化的动态模型时,经常考虑它在时间或空间的微小单元变化情况,这是因为在这些微元上的平衡关系比较简单,而且容易使用微分学的手段进行处理.这类模型基本上是以微分方程的形式给出的. 例1 设警方对司机饮酒后驾车时血液中酒精含量的规定为不超过80%(mg/ml). 现有一起交通事故,在事故发生3个小时后,测得司机血液中酒精含量是56%(mg/ml), 又过两个小时后, 测得其酒精含量降为40%(mg/ml),试判断: 事故发生时,司机是否违反了酒精含量的规定? 解:模型建立 设)(t x 为时刻t 的血液中酒精的浓度, 则依平衡原理时间间隔],[t t t ?+内, 酒精浓度的改变量 t t x x ??∝?)(, 即 t t kx t x t t x ?-=-?+)()()( 其中k >0为比例常数, 式前负号表示浓度随时间的推移是递减的, 遍除以t ?, 并令0→?t , 则得到 ,d d kx t x -= 且满足40)5(,56)3(==x x 以及0)0(x x =. 模型求解 容易求得通解为kt c t x -=e )(, 代入0)0(x x =,得到 kt x t x -=e )(0 则)0(0x x =为所求. 又由,40)5(,56)3(==x x 代入0)0(x x =可得 17.04056e 40e 56e 25030=?=????==--k x x k k k 将17.0=k 代入得 25.93e 5656e 17.03017 .030≈?=?=??-x x >80 第二章 微 分 方 程 模 型 建立微分方程模型就是把物理、化学、生物科学、工程科学和社会科学中的规律和原理用含有待定函数的导数或微分的数学关系式表示出来。这一章我们由浅入深地介绍一些微分方程模型。 2.1 简单模型 例1 物体在空气中的下落与特技跳伞问题 假设质量为m 的物体在空气中下落,空气阻力与物体的速度平方成正比,阻尼系数为k (>0),求物体的运动规律。 解 所谓运动规律即下落距离与时间的关系,如图2.1.1, 建立坐标系。设x 为物体下落的距离,于是物体下落的速度为 dx v dt =, 加速度为 22d x a dt =, 根据牛顿第二定律F ma =,可以列出微分方程 2 22d x d x m k m g d t d t ?? =-+ ???, (2.1.1) 负号表示阻力方向与速度方向相反。 例2 单摆的自由振动问题。 如图2.1.2 为一个单摆,上端固定在O 点,M 为一质量为m 的质点,摆杆OM 之长为L (摆杆的质量忽略不计)。单摆的平衡位置为铅垂线'OO 。将质点M 拉开,使OM 与'OO 成一个角度0θ,然后放手任其自由运动,试求摆杆OM 和铅垂线'OO 的夹角θ与时间t 的关系。 解 将重力分解为径向力F 与切向力T ,T 的大小为sin mg θ,M 的切向加速 度为22d a L dt θ =,于是,由牛顿第二定律,列出微分方程 22s i n d m a m L m g dt θ θ== , 即 22s i n d g dt L θθ=-, (2.1.2) 设初始时刻0t =,摆杆的初始位置为0θ,初始角速度为0,则单摆的运动规律的研究就化为微分方程的初值问题 ()()22 00' 0s i n ,,0.t t d g dt L t t θθθθθ==?=-??? =??=??? (2.1.3) 图2.1.1 图2.1.2 例3 考古和地质学中文物和化石年代的测定问题。 考古、地质学等方面的专家常用14C (碳14)来估计文物或化石的年代。它们的依据是,宇宙射线不断轰击大气层,使之产生中子,中子与氧气作用生成具有放射性的14C 。这种放射性碳可以氧化成二氧化碳。二氧化碳被植物所吸收,而动物又以植物为食物,于是放射性碳就被带到各种动植物体内。由于14C 是放射性的,无论存在于空气中或生物体内它都在不断衰变,活着的生物通过新陈代谢不断地摄取14C ,使得生物体内的14C 与空气中的14C 有相同的百分含量。生物体死后它停止摄取14C ,因而尸体内的14C 由于不断衰变而不断减少。碳定年代法就是根据14C 的衰变减少量的变化情况来判定生物的死亡时间的。 基本假设 (1)现代生物体中14C 的衰变速度与古代生物体中14C 的衰变速度相同(依据是地球周围大气中14C 的百分含量可认为基本不变,即宇宙射线照射大气层的强度自古至今基本不变); (2)14C 的衰变速度与该时刻14C 的含量成正比(这条假设的根据来自于原子物理学理论)。 下面用微分方程建模。 设在时刻t (年)生物体中14C 的存量为()x t ,由假设(2)知 习题解答 1. 系统的微分方程为()4()2()x t x t u t '=-+,其中()u t 是幅度为1,角频率为1rad/s 的方波输入信号,试建立系统的Simulink 模型并进行仿真。 解:用积分器直接构造求解微分方程的模型 由原微分方程()4()2()x t x t u t '=-+可知 x '经积分模块作用就得x ,而x 经代数运算又产生x ',据此可以建立系统模型并仿真,实现建模与仿真步骤如下。 ⑴利用Simulink 模块库中的基本模块,不难建立系统模型,如题图1所示。 题图1 求解微分方程的模型 模型中各个模块说明如下。 ①()u t 输入模块:它的参数设置如题图1(a)所示,模块名称由原来的Pulse Generrator 改为()u t 。 题图1(a) ()u t 输入模块的参数设置 ②Gs 增益模块:增益参数Gain 设置为2。 ③求和模块:其图标形状Icon shape 选择rectangular ,符号列表Lisl of signs 设置为+-。 ④积分模块:参数不需改变。 ⑤G 1增益模块:增益参数设置为4,它的方向旋转可借助Format 菜单中的Rotate Block 命令实现。 ⑥Scope 示波器:在示波器参数设置窗口选择Data history 页,选中其中的Save data to workspace 复选框。这将使送入示波器的数据同时被保存在MA TLAB 工作空间的默认名为ScopeData 的结构矩阵或矩阵中。 ⑵设置系统仿真参数。单击模型编辑窗口Simulation 菜单中的Configuration Parameters 选项,打开仿真参数设置对话框,选择Solver 选项,把仿真的停止时间Sto ptime 设置为20。 ⑶仿真操作。双击示波器图标,打开示波器窗口。选择模型编辑窗口中Simulation 菜单中的Stan 命令,就可在示波器窗口中看到仿真结果的变化曲线,如题图1(b)所示。 题图1(b) 仿真结果 2. 建立使用阶跃信号为输入信号,经过传递函数为1 5.01 s 的一阶系统的Simulink 模型并进行仿真。要求:⑴查看其输出波形在示波器上的显示;⑵修改仿真参数Max step size 为2、Min step size 为1,在示波器上查看波形;⑶修改示波器Y 坐标轴范围为0~2,横坐标范围为0~15,查看波形。 解:⑴①利用Simulink 模块库中的基本模块,不难建立系统模型,如题图2所示。 题图2 一阶系统的Simulink 模型 模型中各个模块说明如下。 ()u t 输入模块:它的step time 被设置为0,模块名称由原来的step 改为()u t 。 Transfer Fon 传递函数模块:在Denominator coefficient 文本框中定义分母多项式系数向量为[0.5 1]。 实验07 微分方程模型(2学时) (第5章 微分方程模型) 1.(验证)传染病模型2(SI 模型)p136~138 传染病模型2(SI 模型): 0(1),(0)di k i i i i dt =-= 其中, i (t )是第t 天病人在总人数中所占的比例。 k 是每个病人每天有效接触的平均人数(日接触率)。 i 0是初始时刻(t =0)病人的比例。 1.1 画~di i dt 曲线图p136~138 取k =0.1,画出i dt di ~的曲线图,求i 为何值时dt di 达到最大值,并在曲线图上标注。 参考程序: 提示:fplot, fminbnd, plot, text, title, xlabel 1)画曲线图 用fplot函数,调用格式如下: fplot(fun,lims) fun必须为一个M文件的函数名或对变量x的可执行字符串。 若lims取[xmin xmax],则x轴被限制在此区间上。 若lims取[xmin xmax ymin ymax],则y轴也被限制。 本题可用 fplot('0.1*x*(1-x)',[0 1.1 0 0.03]); 2)求最大值 用求解边界约束条件下的非线性最小化函数fminbnd,调用格式如下:x=fminbnd('fun',x1,x2) fun必须为一个M文件的函数名或对变量x的可执行字符串。 返回自变量x在区间x1 (微分方程模型) .一个半球状雪堆,其体积融化地速率与半球面面积成正比,比例系数 > .设融化中雪堆始终保持半球状,初始半径为且小时中融化了总体积地,问雪堆全部融化还需要多长时间? .从致冰厂购买了一块立方体地冰块,在运输途中发现,第一小时大约融化了 ()求冰块全部融化要多长时间(设气温不变) ()如运输时间需要小时,问:运输途中冰块大约会融化掉多少? .一展开角为α地圆锥形漏斗内盛着高度为地水,设漏斗底部地孔足够大(表面张力不计),试求漏斗中地水流光需要多少时间? .容器甲地温度为度,将其内地温度计移入容器乙内,设十分钟后温度计读数为度,又过十分钟后温度计读数为度,试求容器乙内地温度. .一块加过热地金属块初始时比室温高度,分钟测得它比室温高度,问:()小时后金属块比室温高多少?()多少时间后,金属块比室温高度? .设初始时容器里盛放着含净盐千克地盐水升,现对其以每分钟升地速率注入清水,容器内装有搅拌器能将溶液迅时搅拌均匀,并同时以每分钟升地速率放出盐水,求小时后容器里地盐水中还含有多少净盐? .某伞降兵跳伞时地总质量为公斤(含武器装备),降落伞张开前地空气阻力为,该伞降兵地初始下落速度为,经秒钟后降落伞打开,降落伞打开后地空气阻力约为试球给伞降兵下落地速度(),并求其下落地极限速度. .年月日英国人创建了一项最低开伞地跳伞纪录,它从比萨斜塔上跳下,到离地英尺时才打开降落伞,试求他落地时地速度. .证明对数螺线上任一处地切线与极径地夹角地正切为一常数,().实验证明,当速度远低于音速时,空气阻力正比与速度,阻力系数大约为.现有一包裹从离地米高地飞机上落下,()求其落地时地速度()如果飞机高度更大些,结果会如何,包裹地速度会随高度而任意增大吗? .生态学家估计人地内禀增长率约为,已知年世界人口数为亿(×)而当时地人口增长率则为.试根据模型计算:()世界人口数地上限约为多少()何时将是世界人口增长最快地时候? .早期肿瘤地体积增长满足模型(λ,其中λ为常数),()求肿瘤地增倍时间 σ.根据统计资料,一般有σ()(单位为天),肺部恶性肿瘤地增倍时间大多大于天而小于天(发展太快与太慢一般都不是恶性肿瘤),故σ是确定肿瘤性质地重要参数之 习 题 2.1 什么是线性系统?其最重要的特性是什么?下列用微分方程表示的系统中,x o 表示系统输出,x i 表示系统输入,哪些是线性系统? (1) x x x x x i o o o o 222=++ (2) x tx x x i o o o 222=++ (3) x x x x i o 222o o =++ (4) x tx x x x i o o o 222o =++ 解: 凡是能用线性微分方程描述的系统就是线性系统。线性系统的一个最重要特性就是它满足叠加原理。该题中(2)和(3)是线性系统。 2.2 图(题2.2)中三同分别表示了三个机械系统。求出它们各自的微分方程,图中x i 表示输入位移,x o 表示输出位移,假设输出端无负载效应。 图(题2.2) 解: (1)对图(a)所示系统,由牛顿定律有 x m x c x x c i o o 2 o 1 )(=-- 即 x c x c c x m i 1 2 1 o o )(=++ (2)对图(b)所示系统,引入一中间变量x,并由牛顿定律有 )1()()(1 x x c k x x o i -=- )2()(2 x k x x c o o =- 消除中间变量有 x ck x k k x k k c i o 1 2 1 o 2 1 )(=-- (3)对图(c)所示系统,由牛顿定律有 x k x x k x x c o o i o i 2 1 )()(=-+- 即 x k x c x k k x c i i o o 1 2 1 )(+=++ 2.3求出图(题2.3)所示电系统的微分方程。 图(题2.3) 解:(1)对图(a)所示系统,设i 1为流过R 1的电流,i 为总电流,则有 ?+=idt C i R u o 12 2 i R u u o i 1 1=- 第2章 微分方程模型2(2课时) 教学目的 懂得如何根据实际情况建立微分方程。 教学内容 介绍应用微分方程方法建模。 模型Ⅱ 人口模型 1. 引言 在研究某些实际问题时,经常无法直接得到各变量之间的联系,问题的特性往往会给出关于变化率的一些关系。利用这些关系,我们可以建立相应的微分方程模型。在自然界以及工程技术领域中,微分方程模型是大量存在的。它甚至可以渗透到人口问题以及商业预测等领域中去,其影响是广泛的。 问题:人口问题是当今世界上最令人关注的问题之一,一些发展中国家的人口出生率过高,越来越威胁着人类的正常生活,有些发达国家的自然增长率趋于零,甚至变为负数,造成劳动力紧缺,也是不容忽视的问题。另外,在科学技术和生产力飞速发展的推动下,世界人口以空前的规模增长,统计数据显示: 年 1625 1830 1930 1960 1974 1987 1999 人口(亿) 5 10 20 30 40 50 60 可以看出,人口每增长十亿的时间,由一百年缩短为十二三年。我们赖以生存的地球,已经带着它的60亿子民踏入了21世纪。 长期以来,人类的繁衍一直在自发地进行着。只是由于人口数量的迅速膨胀和环境质量的急剧恶化,人们才猛然醒悟,开始研究人类和自然的关系,人口数量的变化规律,以及如何进行人口控制等问题。 我国是世界第一人口大国,地球上每九个人中就有一个中国人,在20世纪的一段时间内我国人口的增长速度过快,如下表: 年 1908 1933 1953 1964 1982 1990 2000 人口(亿) 3.0 4.7 6.0 7.2 10.3 11.3 12.95 有效地控制人口的增长,不仅是使我国全面进入小康社会、到21世纪中叶建 微分方程模型 1.一个半球状雪堆,其体积融化的速率与半球面面积S成正比,比例系数k > 0。设融化中雪堆始终保持半球状,初始半径为R且3小时中融化了总体积的7/8,问雪堆全部融化还需要多长时间? 2.从致冰厂购买了一块立方体的冰块,在运输途中发现,第一小时大约融化了1/4 (1)求冰块全部融化要多长时间(设气温不变) (2)如运输时间需要2.5小时,问:运输途中冰块大约会融化掉多少? 3.一展开角为α的圆锥形漏斗内盛着高度为H的水,设漏斗底部的孔足够大(表面张力不计),试求漏斗中的水流光需要多少时间? 4.容器甲的温度为60度,将其内的温度计移入容器乙内,设十分钟后温度计读数为70度,又过十分钟后温度计读数为76度,试求容器乙内的温度。 5.一块加过热的金属块初始时比室温高70度,20分钟测得它比室温高60度,问:(1)2小时后金属块比室温高多少?(2)多少时间后,金属块比室温高10度? 6.设初始时容器里盛放着含净盐10千克的盐水100升,现对其以每分钟3升的速率注入清水,容器内装有搅拌器能将溶液迅时搅拌均匀,并同时以每分钟2升的速率放出盐水,求1小时后容器里的盐水中还含有多少净盐? 7.某伞降兵跳伞时的总质量为100公斤(含武器装备),降落伞张开前的空气阻力为0.5v,该伞降兵的初始下落速度为0,经8秒钟后降落伞打开,降落伞打 开后的空气阻力约为0.6试球给伞降兵下落的速度v(t),并求其下落的极限速度。 8. 1988年8月5日英国人Mike McCarthy创建了一项最低开伞的跳伞纪录,它从比萨斜塔上跳下,到离地179英尺时才打开降落伞,试求他落地时的速度。 9.证明对数螺线r=A上任一处的切线与极径的夹角的正切为一常数, () 10.实验证明,当速度远低于音速时,空气阻力正比与速度,阻力系数大约为0.005。现有一包裹从离地150米高的飞机上落下,(1)求其落地时的速度(2)如果飞机高度更大些,结果会如何,包裹的速度会随高度而任意增大吗?11.生态学家估计人的内禀增长率约为0.029,已知1961年世界人口数为30.6 亿(3.06×)而当时的人口增长率则为0.02。试根据Logistic模型计算: 实验一SDH网元基本配置 一、实验目的: 通过本实验,了解SDH光传输的原理和系统组成,了解ZXMP S325设备的硬件构成和单板功能,学习ZXONM 300 网管软件的使用方法,掌握SDH 网元配置的基本操作。 二、实验器材: 1、SDH 设备:3 套ZXMP 325; 2、实验用维护终端。 三、实验原理 1、SDH 原理 同步数字体制(SDH)是为高速同步通信网络制定的一个国际标准,其基础在于直接同步复用。按照SDH 组建的网络是一个高度统一的、标准化的、智能化的网络,采用全球统一的接口以实现多环境的兼容,管理操作协调一致,组网与业务调度灵活方便,并且具有网络自愈功能,能够传输所有常见的支路信号,应用于多种领域(如光纤传输,微波和卫星传输等)。 SDH 具有以下特点: (1)接口:接口的规范化是设备互联的关键。SDH对网络节点接口(NNI)作了统一的规范,内容包括数字信号数率等级、帧结构、复接方法、线路接口、监控管理等。 电接口:STM-1是SDH的第一个等级,又叫基本同步传送模块,比特率为155.520Mb/s;STM-N是SDH第N个等级的同步传送模块,比特率是STM-1的N 倍(N=4n=1,4,16,- - -)。 光接口:采用国际统一标准规范。SDH仅对电信号扰码,光口信号码型是加扰的NRZ 码,信号数率与SDH 电口标准信号数率相一致。 (2)复用方式 a)低速SDH----高速SDH,字节间插; b) 低速PDH-----SDH,同步复用和灵活的映射。 (3)运行维护:用于运行维护(OAM)的开销多,OAM功能强——这也是线路编码不用加冗余的原因. (4)兼容性:SDH 具有很强的兼容性,可传送PDH 业务,异步转移模式信号(ATM)及其他体制的信号。 (5)SDH 复用映射示意图 第二章:动力学系统的微分方程模型 利用计算机进行仿真时,一般情况下要给出系统的数学模型,因此有必要掌 握一定的建立数学模型 的方法。在动力学系统中,大多数情况下可以使用微分方程 来表示系统的动态特性,也可以通过微分方程可以将原来的系统简化为状态方程或 者差分方程模型等。在这一章中,重点介绍建系统动态问题的微分方程的基本理论 和方法。 在实际工程中,一般把系统分为两种类型,一是连续系统;其数学模型一般 是高阶微分方程;另一 种是离散系统,它的数学模型是差分方程。 § 2.1 动力学系统统基本元件 任何机械系统都是由机械元件组成的,在机械系统中有 3种类型的基本机械元 件:惯性元件、弹性元件和阻尼元件。 1惯性元件:惯性元件是指具有质量或转动惯量的元件, 惯量可以定义为使加速度 (或角加速度)产生单位变化所需要的力(或力矩)。 2弹性元件:它在外力或外力偶作用下可以产生变形的元件, 力做功来储存能量。按变形性质可以分为线性元件和非线性元件,通常等效成一弹 簧来表示。 对于线性弹簧元件,弹簧中所受到的力与位移成正比, 比例常数为弹簧刚度 k 。 F Wx 这里k 称为弹簧刚度, 级是弹簧相对于原长的变形量,弹性力的方向总是指向弹 簧的原长位移,出了弹簧和受力之间是线性关系以外,还有所谓硬弹簧和软弹簧, 它们的受力和弹簧变形之间的关系是一非线性关系。 3阻尼元件:这种元件是以吸收能量以其它形式消耗能量, 而不储存能量,可以形 象的表示为一个活塞在一个充满流体介质的油缸中运动。阻尼力通常表示为: D a R = ex 阻尼力的方向总是速度方向相反。当 1,为线性阻尼模型。否则为非线性阻 尼模型。应注意当:等于偶数情况时,要将阻尼力表示为: R - -ex | x 4 | 这里的"-”表示与速度方向相反 惯量(质量) 力(N ) 加速度(m/ s 2 ) 惯量(转动惯量) 力矩(N m ) 角加速度(rad / s 2 ) 这种元件可以通过外 第二章微分方程管模型 第一节人口学模型 人口问题是当今世界人类面临的五大问题的首要问题。我国是世界上人口最多的国家,由于20世纪五六十年代人口政策方面的失误,不仅造成人口总数增长过快,而且年龄结构也不合理。人口的过份增长给我国经济发展造成沉重袍袱,严重地影响经济建设。能否有效地控制人口的增长,已成为本世纪中叶我国能否达到中等发达国家以至赶上发达国家的关键。 建立数学模型对人口发展过程进行描述、分析和预测,并进而研究控制人口增长的生育政策,已引起有关人口专家和官员的极大关注和兴趣。以下我们就如何建立人口数学模型作简要的介绍。 一. 马尔萨斯人口增长模型 对于如何预测人口的增长,早在8世纪人们就开始进行了。英国早期经济学家马尔萨斯(1766-1843)他在担任牧师期间,根据教堂100多年人口出生统计资料,他发现人口出生率是一个常数。于是在1798年发表的《人口原理》一书中提出了哄动于世的马尔萨斯人口模型:假设x(t)表示t时刻人口总数,r为人口增长率(常数),其他影响人口增长的因素均不考虑。则在t到t+△t这段时间内人口总数增长为 ()()() +?-=? x t t x t rx t t 两端同除以t ?,并令0t ?→,得 ()00dx rx dt x t x ?=???=? ……… (1) 其解为 ()() 00r t t x t x e -= (2) 方程(1)称为马尔萨斯模型,它的解是一个以r e 为公比的几何级数。马尔萨斯根据这个模型认为人口的增长是按几何级数增加,而物质的增长只能按算术级数增加。因此,人口必须加以控制。不幸的是马尔萨斯的这一忠告没有引起人们的足够重视。 马尔萨斯模型 (1) 和 (2) 与19世纪以前欧州一些地区的人口统计数据十分吻合。据估计1961年全世界人口总数为93.0610?,而在此之前的十来年间人口按年2%的速率增长。因此有 9001961, 3.0610,0.02t x r ==?= 于是 ()()0.02196193.0610t x t e -=?? (3) 根据统计资料,在1700—1961年间世界人口大约每35年增加一倍。由(3)式可以算出每34.6年人口增加一倍。 事实上,设在0T t t =-内地球上人口增加一倍,即当0t t =时, 90.0209 0 3.0610 3.0610x e ?=??=?,当0 T t t =-时,90.0202 3.0610,T x e =?所以 第5章 微分方程模型 一、讨论题 1. 某人每天由饮食获取10467焦热量,其中5038焦用于新陈代谢,此外每公斤体重需支付69焦热量作为运动消耗,其余热量则转化为脂肪,已知以脂肪形式贮存的热量利用率为100%,每公斤脂肪含热量41868焦,问此人的体重如何随时间而变化? 2. 根据罗瑟福的放射性衰变定律,放射性物质衰变的速度与现存的放射性物质的原子数成正比,比例系数成为衰变系数,试建立放射性物质衰变的数学模型。若已知某放射性物质经时间21T 放射物质的原子下降至原来的一半(21T 称为该物质的半衰期)试决定其衰变系数。 3. 建立耐用消费品市场销售量的模型。如果已知了过去若干时期销售量的情况,如何确定模型的参数。 4. 两棵不同类别的植物种在一起,按比例吸取养料,试建立它们的生长模型。 二、思考题 1、具有什么样特征的数学建模问题需要用微分方程方法建立模型? 2、用微分方程方法建立数学模型的基本步骤是什么?应注意哪些问题? 3、某种细菌的增长率不知道,但假设它是常数,试验开始时估计大约有11500个细菌,一时后有2000个,问四时后大约有多少细菌? 三、习题 1. 生活在阿拉斯加海滨的鲑鱼服从Malthus 增长模型 )(003.0) (t p dt t dp = 其中t 以分钟计。在0=t 时一群鲨鱼来到此水域定居,开始捕食鲑鱼。鲨鱼捕杀鲑鱼的速 率是)(001.02 t p ,其中)(t p 是t 时刻鲑鱼总数。此外,由于在它们周围出现意外情况,平 均每分钟有0.002条鲑鱼离开此水域。 (1)考虑到两种因素,试修正Malthus 模型。 (2)假设在0=t 是存在100万条鲑鱼,试求鲑鱼总数 )(t p ,并问∞→t 时会发生什么情况? 2. 用具有放射性的14 C 测量古生物年代的原理是:宇宙线轰击大气层产生中子,中子 与氮结合产生14C 。植物吸收二氧化碳时吸收了14C ,动物食用植物从植物中得到14 C 。在 活组织中14C 的吸收速率恰好与14 C 的衰变速率平衡。但一旦动植物死亡,它就停止吸收14C ,于是14C 的浓度随衰变而降低。由于宇宙线轰击大气层的速度可视为常数,既动物刚 死亡时14 C 的衰变速率与现在取的活组织样本(刚死亡)的衰变速率是相同的。若测得古生 物标本现在14C 的衰变速率,由于14 C 的衰变系数已知,即可决定古生物的死亡时间。试建 立用14C 测古生物年代的模型(14 C 的半衰期为5568年)。 3. 试用上题建立的数学模型,确定下述古迹的年代: (1)1950年从法国Lascaux 古洞中取出的碳测得放射性计数率为0.97计数(min ?g ),而活树木样本测得的计数为6.68计数(min ?g ),试确定该洞中绘画的年代; (2)1950年从某古巴比伦城市的屋梁中取得碳标本测得计数率为4.09计数(min ?g ),活数标本为6.68计数(min ?g ),试估计该建筑的年代。 4. 一容器用一薄膜分成容积为A V 和B V 的两部分,分别装入同一物质不同浓度的溶液。实验07讲评参考答案_微分方程模型(2学时)
微分方程建模案例
第二章动力学系统的微分方程模型
3.1 微分方程模型的建模步骤
常微分课后答案解析第二章
常微分课后答案解析第二章
微分方程模型建模实例
微分方程模型
第二章 微 分 方 程 模 型.
第2章习题解答
数学建模实验答案微分方程模型
微分方程模型习题
2机械控制工程基础第二章答案解析
微分方程模型2(2课时)
微分方程模型(二)
扩散问题的偏微分方程模型_数学建模
第二章:动力学系统的微分方程模型
[数学]第二章微分方程模型
第5章 微分方程模型