第二章：微积分方法建模

overview:

涉及连续的变量，可以用微积分求解，求得解析式便于下一步分析。有些离散的变量也可以演变成连续变量进行求解。当我们描述实际对象的某些特性随时间（或空间）而演变的过程，分析它的变化规律，预测它的未来性态时，通常要建立对象的动态模型。建模时首先要根据建模目的和对问题的具体分析作出简化假设，然后按照对象内在的或可以类比的其它对象的规律列出微分方程，求出方程的解并将结果翻译回实际对象，就可以进行描述、分析或预测了。

2.1飞机降落曲线

飞机降落时，其曲线是一条三次抛物线，如下图：

水平速度为常数u，出于安全考虑，飞机的垂直加速度不得超过g/10，已知飞行高度h，要求在坐标原点降落，求开始下降点x₀允许的最小值。

求曲线解析式

设飞机降落曲线为：

y=ax3+bx2+cx+dy = ax^3+bx^2+cx+dy=ax3+bx2+cx+d

依据假设可知：
y(0)=0,y(x0)=h,y′(0)=0,y′(x0)=0y(0)=0,y(x_0)=h,y'(0)=0,y'(x_0)=0y(0)=0,y(x0)=h,y′(0)=0,y′(x0)=0，带入上式可得：

{y(0)=d=0,y′(0)=c=0,y(x0)=ax3+bx2+cx+d=h,y′(x0)=3ax2+2bx0+c=0x=0\small \begin{cases} y(0)=d=0,\\ y'(0)=c=0,\\ y(x_0)=ax^3+bx^2+cx+d=h,\\ y'(x_0)=3ax^2+2bx_0+c=0x=0\\ \end{cases} ⎩⎪⎪⎨⎪⎪⎧y(0)=d=0,y′(0)=c=0,y(x0)=ax3+bx2+cx+d=h,y′(x0)=3ax2+2bx0+c=0x=0
所以:
a=−2hx03,b=3hx02,c=0,d=0a = -\frac{2h}{{x_0}^3},b=\frac{3h}{{x_0}^2},c=0,d=0a=−x032h,b=x023h,c=0,d=0

其曲线为：
y=−hx02(2x0x3−3x2)y=-\frac{h}{{x_0}^2}(\frac{2}{x_0}x^3-3x^2)y=−x02h(x02x3−3x2)

求最佳着陆点

飞机垂直速度是y关于时间t的导数：
dydt=−hx02(6x0x2−6x)dxdt\frac{dy}{dt}=-\frac{h}{{x_0}^2}(\frac{6}{x_0}x^2-6x)\frac{dx}{dt}dtdy=−x02h(x06x2−6x)dtdx

其中，dxdt\frac{dx}{dt}dtdx是飞机水平速度，所以dxdt=u\frac{dx}{dt}=udtdx=u，则：
dydt=−6hux02(x2x0−x)\frac{dy}{dt}=-\frac{6hu}{{x_0}^2}(\frac{x^2}{x_0}-x)dtdy=−x026hu(x0x2−x)

垂直加速度为：d2ydt2=−6hux02(2xx0−1)dxdt=−6hu2x02(2xx0−1)\frac{d^2y}{dt^2}=-\frac{6hu}{{x_0}^2}(\frac{2x}{x_0}-1)\frac{dx}{dt}=-\frac{6hu^2}{x_0^2}(\frac{2x}{x_0}-1)dt2d2y=−x026hu(x02x−1)dtdx=−x026hu2(x02x−1)

设a(x)=d2ydt2a(x)=\frac{d^2y}{dt^2}a(x)=dt2d2y,则：

∣a(x)∣=6hu2x02∣2xx0−1∣,x∈[0,x0]|a(x)|=\frac{6hu^2}{{x_0}^2}|\frac{2x}{x_0}-1|,x\in[0,x_0]∣a(x)∣=x026hu2∣x02x−1∣,x∈[0,x0]

所以，垂直加速度最大值为：

max∣a(x)∣=6hu2x02,x∈[0,x0]max|a(x)|=\frac{6hu^2}{{x_0}^2},x\in[0,x_0]max∣a(x)∣=x026hu2,x∈[0,x0]

由于题目限制：max∣a(x)∣≤g10max|a(x)|\leq\frac{g}{10}max∣a(x)∣≤10g，所以：

x0≥u60hgx_0\geq u\sqrt{\frac{60h}{g}}x0≥ug60h，此为x0x_0x0的最小值。

2.2经济增长模型

发展经济、增加生产有两个重要因素，一是增加投资（扩大厂房、购买设备、技术革新等），二是增加劳动力。恰当调节投资增长和劳动力增长的关系，使增加的产量不致被劳动力的增长抵消，劳动生产率才能不断提高，从而真正起到发展经济的作用。为此，需要分析产量、劳动力和投资之间变化规律，从而保证经济正常发展。
假设：
Q(t)Q(t)Q(t)——地区、部门、企业在某时间t的产量
L(t)L(t)L(t)——地区、部门、企业在某时间t的劳动力
K(t)K(t)K(t)——地区、部门、企业在某时间t的资金
Z(t)Z(t)Z(t)——每个劳动力在某时间t所占有的产量（劳动力生产率）

道格拉斯（Douglas）生产函数

由于我们只关心产量、劳动力、投资的相对增长，而不是绝对量，所以定义：
产量指数：iQ(t)=Q(t)Q(0)i_Q(t) = \frac{Q(t)}{Q(0)}iQ(t)=Q(0)Q(t)，劳动力指数：iL(t)=L(t)L(0)i_L(t) = \frac{L(t)}{L(0)}iL(t)=L(0)L(t)，投资指数：iK(t)=K(t)K(0)i_K(t) = \frac{K(t)}{K(0)}iK(t)=K(0)K(t) （1）
令：ξ(t)=ln⁡iL(t)iK(t)\xi(t)=\ln\frac{i_L(t)}{i_K(t)}ξ(t)=lniK(t)iL(t),ψ(t)=ln⁡iQ(t)iK(t)\psi(t)=\ln\frac{i_Q(t)}{i_K(t)}ψ(t)=lniK(t)iQ(t) （2）
根据统计数据,散点(ξ,ψ)(\xi,\psi)(ξ,ψ)在直角坐标系下的图像大致为：

由图可知，大多数的点处在一条经过原点的直线附近，故ξ\xiξ和ψ\psiψ有如下关系：

ψ=γξ(0<γ<1)\psi=\gamma\xi (0<\gamma<1)ψ=γξ(0<γ<1) （3）
带入上式可得：
iQ(t)=iLγ(t)iK1−γ(t)i_Q(t)={i_L}^\gamma(t){i_K}^{1-\gamma}(t)iQ(t)=iLγ(t)iK1−γ(t) （4）
记a=Q(0)L−γ(0)K1−γ(0)a=Q(0)L^{-\gamma}(0)K^{1-\gamma}(0)a=Q(0)L−γ(0)K1−γ(0),则由（1）和（4）可知:
Q(t)=aLγ(t)K1−γ(t)Q(t)=aL^\gamma(t)K^{1-\gamma}(t)Q(t)=aLγ(t)K1−γ(t) （5）
式（5）就是经济学中著名的Douglas生产函数，它表明产量余劳动力、投资之间的关系，即：
Q˙Q=γL˙L+(1−γ)K˙K\frac{\dot{Q}}{Q}=\gamma\frac{\dot{L}}{L}+(1-\gamma)\frac{\dot{K}}{K}QQ˙=γLL˙+(1−γ)KK˙ （6）
其中Q˙\dot{Q}Q˙、L˙\dot{L}L˙、K˙\dot{K}K˙表示Q、L、K关于t的导数。
这个函数表明，相对增长量，Q˙Q\frac{\dot{Q}}{Q}QQ˙、L˙L\frac{\dot{L}}{L}LL˙、K˙K\frac{\dot{K}}{K}KK˙之间呈线性关系，当γ→1\gamma\rightarrow1γ→1时产量增长主要依靠劳动力增长，当γ→0\gamma\rightarrow0γ→0时产量增长主要依靠投资。亦称γ\gammaγ为产量对劳动力的弹性系数

劳动生产率增长条件

定义劳动生产率Z(t)=Q(t)L(t)Z(t)=\frac{Q(t)}{L(t)}Z(t)=L(t)Q(t)，则Z˙Z=Q˙Q−L˙L\frac{\dot{Z}}{Z}=\frac{\dot{Q}}{Q}-\frac{\dot{L}}{L}ZZ˙=QQ˙−LL˙,将（6）代入Z(t)Z(t)Z(t)可得：
Z˙Z=(1−γ)(K˙K−L˙L)\frac{\dot{Z}}{Z}=(1-\gamma)(\frac{\dot{K}}{K}-\frac{\dot{L}}{L})ZZ˙=(1−γ)(KK˙−LL˙) （7）
由（7）可见，只要K˙K>L˙L\frac{\dot{K}}{K}>\frac{\dot{L}}{L}KK˙>LL˙就能保证Z˙Z>0\frac{\dot{Z}}{Z}>0ZZ˙>0,即劳动生产率的提高需要由投资的相对增长大于劳动力的相对增长为前提条件。

2.3存贮模型

原料、商品的存贮问题，存的太多，存贮费用过高；存的太少，无法及时满足需求。目的：制定最有存贮策略，即：多长时间顶一次货，每次顶多少货，才能使总费用最小。

模型一

模型假设

每次订货费用为C1C_1C1，每天每吨货物贮存费用为C2C_2C2(已知)
每天货物的需求量r吨为已知
订货周期为T天，每次订货Q吨，当贮存量降到0时订货立即到达。

模型建立

订货周期、订货量余每天的需求量存在关系：
Q=rTQ=rTQ=rT
订货以后贮存量q(t)q(t)q(t)均匀地下降，即q(t)=Q−rtq(t)=Q-rtq(t)=Q−rt,如下图：

则一个订货周期的费用：
{订货费:C1,存贮费:∫0Tq(t)dt=12C2QT=12C2rT2\small \begin{cases} 订货费:C_1,\\ 存贮费:\int_{0}^{T}q(t)\,dt=\frac{1}{2}C_2QT=\frac{1}{2}C_2rT^2\\ \end{cases} {订货费:C1,存贮费:∫0Tq(t)dt=21C2QT=21C2rT2
即：C(T)=C1+12C2rT2C(T)=C_1+\frac{1}{2}C_2rT^2C(T)=C1+21C2rT2
则一个订货周期每天的费用为：C‾(T)=C1T+12C2rT\overline{C}(T)=\frac{C_1}{T}+\frac{1}{2}C_2rTC(T)=TC1+21C2rT

模型求解

令dC‾dT=0\frac{d\overline{C}}{dT}=0dTdC=0,可求得：T=2c1RC2T=\sqrt{\frac{2c_1}{RC_2}}T=RC22c1
进而：Q=2C1rC2Q=\sqrt{\frac{2C_1r}{C_2}}Q=C22C1r
此模型称为经济订货批量公式，简称EOQ公式

模型二允许缺货的存贮模型

模型假设

每次订货费用为C1C_1C1,每天每吨货物贮存费用C2C_2C2(已知)
每天的货物需求量r吨（已知）
订货周期为T天，订货量Q，允许缺货，每天每吨货物缺货费C3C_3C3（已知）

模型建立

缺货时的存储量q视为负值，则q(t)q(t)q(t)的图像如下，货物在t=T1t=T_1t=T1时用完，于时Q=rT1Q=rT_1Q=rT1

则一个订货周期内总费用为：
{订货费:C1,存贮费:∫0T1q(t)dt=12C2QT1=12C2Q2r,缺货费:∫T1T∣q(t)∣dt=C32r(T−T1)2=C32r(rT−Q)2\small \begin{cases} 订货费:C_1,\\ 存贮费:\int_{0}^{T_1}q(t)\,dt=\frac{1}{2}C_2QT_1=\frac{1}{2}C_2\frac{Q^2}{r},\\ 缺货费:\int_{T_1}^{T}|q(t)|\,dt=\frac{C_3}{2}r(T-T_1)^2=\frac{C_3}{2r}(rT-Q)^2 \end{cases} ⎩⎪⎨⎪⎧订货费:C1,存贮费:∫0T1q(t)dt=21C2QT1=21C2rQ2,缺货费:∫T1T∣q(t)∣dt=2C3r(T−T1)2=2rC3(rT−Q)2
即：C(T,Q)=C1+12C2Q21r+12rC3(rT−Q)2C(T,Q)=C_1+\frac{1}{2}C_2Q^2\frac{1}{r}+\frac{1}{2r}C_3(rT-Q)^2C(T,Q)=C1+21C2Q2r1+2r1C3(rT−Q)2
则平均每天的费用为：
C‾(T,Q)=C(T,Q)T=C1T+C2Q22rT+C3(rT−Q)22rT\overline{C}(T,Q)=\frac{C(T,Q)}{T}=\frac{C_1}{T}+\frac{C_2Q^2}{2rT}+\frac{C_3(rT-Q)^2}{2rT}C(T,Q)=TC(T,Q)=TC1+2rTC2Q2+2rTC3(rT−Q)2

模型求解

{∂C‾∂T=0∂C‾∂Q=0\small \begin{cases} \frac{\partial\overline{C}}{\partial T}=0\\ \frac{\partial\overline{C}}{\partial Q}=0 \end{cases} {∂T∂C=0∂Q∂C=0
求出T、Q的最优解，分别记为T′、Q′T'、Q'T′、Q′
T′=2C1rC2C2+C3C3,Q′=2C1rC2C3C2+C3T'=\sqrt{\frac{2C_1}{rC_2}\frac{C_2+C_3}{C_3}},Q'=\sqrt{\frac{2C_1r}{C_2}\frac{C_3}{C_2+C_3}}T′=rC22C1C3C2+C3,Q′=C22C1rC2+C3C3

分析

令μ=C2+C3C3\mu=\sqrt{\frac{C_2+C_3}{C_3}}μ=C3C2+C3,与模型一相比，有：
T′=μT,Q′=QμT'=\mu T,Q'=\frac{Q}{\mu}T′=μT,Q′=μQ
显然，T’ > T,Q’ < Q,即在允许缺货时应增大订货周期，减少订货批次；当缺货非C3C_3C3相对于存贮费C2C_2C2而言越大时，μ\muμ越小，T’和Q’越接近T和Q。

2.4城市人口统计模型

模型一（估算城市现有人口）

城市人口密度常用P(r)=br2+aP(r)=\frac{b}{r^2+a}P(r)=r2+ab或者P(r)=ae−br(a,b>0)P(r)=ae^{-br}(a,b>0)P(r)=ae−br(a,b>0)来近似表示，其中r是距城中心的距离。则计算距离市中心C区域内的人口数N可以这样：从城市中心画一条射线，把这条线上从0到C之间n等分，每小区间长度为Δr\Delta rΔr，每小区间确定一个圆环，第j个圆环面积为从城市中心为：
πrj2−πrj−12=πrj2−π(rj−Δr)2=2πrjΔr−π(Δr)2≈2πrjΔr\pi{r_j}^2-\pi{r_{j-1}}^2=\pi{r_j}^2-\pi(r_j-\Delta r)^2=2\pi r_j\Delta r-\pi(\Delta r)^2\approx2\pi r_j \Delta rπrj2−πrj−12=πrj2−π(rj−Δr)2=2πrjΔr−π(Δr)2≈2πrjΔr (Δr很小)(\Delta r很小)(Δr很小)

第j个圆环上的人口数近似为P(rj)2πrjΔrP(r_j)2\pi r_j \Delta rP(rj)2πrjΔr,所以：

N≈∑j=1nP(rj)2πrjΔrN\approx \sum_{j=1}^{n}P(r_j) 2\pi r_j \Delta rN≈∑j=1nP(rj)2πrjΔr

令n→∞n\rightarrow\inftyn→∞时，

N=∫0CP(r)2πrdrN=\int_{0}^{C}P(r)2\pi r\,drN=∫0CP(r)2πrdr

模型二（预测城市未来人口）

P(t)表示t时刻的城市人口数量，人口变化手下面规则的影响：

t时刻的净增长人口以每年r(t)的比率增加
在一段世界内，由于自然死亡和人口迁移，T1T_1T1时刻的人口数P(T1)P(T_1)P(T1)的一部分在T2T_2T2时刻仍然存在，用h(T2−T1)P(T1)h(T_2-T_1)P(T_1)h(T2−T1)P(T1)来表示，这里0<h(T2−T1)<10<h(T_2-T_1)<10<h(T2−T1)<1，T2−T1T_2-T_1T2−T1是这段时间的长度。

把(0,T]时间划分为n等分，每个小区间长度为Δt\Delta tΔt。
假设初始时刻人口为P(0)，到时刻T将只剩h(T)P(0)。当Δt\Delta tΔt很小的时候，从tj−1t_{j-1}tj−1到tjt_jtj，净增长的人口比率近似为常数r(tj)r(t_j)r(tj)。这段时间净增的人口数近似为r(tj)Δtr(t_j)\Delta tr(tj)Δt,在tjt_jtj时刻的人口到T时刻只剩h(T−tj)r(tj)Δth(T-t_j)r(t_j)\Delta th(T−tj)r(tj)Δt。所以在T时刻的总人口数近似为：

P(T)≈h(T)P(0)+∑j=1nh(T−tj)r(tj)ΔtP(T)\approx h(T)P(0)+\sum_{j=1}^{n}h(T-t_j)r(t_j)\Delta tP(T)≈h(T)P(0)+∑j=1nh(T−tj)r(tj)Δt

令n→∞n\rightarrow\inftyn→∞，得：

P(T)=h(T)P(0)+∫0Th(T−t)r(t)dtP(T)=h(T)P(0)+\int_{0}^{T}h(T-t)r(t)\,dtP(T)=h(T)P(0)+∫0Th(T−t)r(t)dt

2.7古生物年代确定

主要思路：C14C^{14}C14的半衰期为5568年，根据C14C^{14}C14的蜕变减少量的变化来判断生物的死亡时间。

模型假设

地球周围大气的C14C^{14}C14可以认为不变，现代生物和古代生物体内的C14C^{14}C14蜕变速度一致。
由原子物理学，C14C^{14}C14的蜕变速度与该时刻的C14C^{14}C14含量成正比。

模型建立

记，x(t)x(t)x(t)代表t时刻生物体内C14C^{14}C14的含量。由假设可知：dxdt=−kx,k>0（1）\frac{dx}{dt}=-kx,k>0 （1）dtdx=−kx,k>0（1）
（1）的通解为x(t)=Ce−ktx(t)=Ce^{-kt}x(t)=Ce−kt，设生物死亡时间为t0t_0t0，代入可知，C=x0C=x_0C=x0，于时：
x(t)=x0e−kt(2)x(t)=x_0e^{-kt} (2)x(t)=x0e−kt(2)
记C14C^{14}C14的半衰期为T，则有：
x(T)=x02（3）x(T)=\frac{x_0}{2} （3）x(T)=2x0（3）
将（3）代入（2）可得k=ln2Tk=\frac{ln2}{T}k=Tln2，故可解得：
t=Tln2lnx0x(t)（4）t=\frac{T}{ln2}ln\frac{x_0}{x(t)} （4）t=ln2Tlnx(t)x0（4）
由于x0和x(t)x_0和x(t)x0和x(t)难以测量，故使用另一种方法求t：
对（2）两边求导：
x′(t)=−x0ke−kt=−kx(t)x'(t)=-x_0ke^{-kt}=-kx(t)x′(t)=−x0ke−kt=−kx(t)
x′(0)=−kx(0)=−kx0x'(0)=-kx(0)=-kx_0x′(0)=−kx(0)=−kx0
两式相除，可得：x′(t)x′(0)=x0x(t)\frac{x'(t)}{x'(0)}=\frac{x_0}{x(t)}x′(0)x′(t)=x(t)x0,代入(4)可得：
t=Tln2lnx′(0)x′(t)t=\frac{T}{ln2}ln\frac{x'(0)}{x'(t)}t=ln2Tlnx′(t)x′(0)
所以只需要测出标本C14C^14C14的平均蜕变速度（单位：次/分钟），即x′(t)x'(t)x′(t)，和现在的C14C^14C14的蜕变速度x′(0)x'(0)x′(0)，即可求出t。

2.8预测人口的增长

模型一：Malthhus指数增长模型

模型假设

1.某国/地区再t时刻的人口数x(t)为连续可微函数。
2.人口的增长率r是常数，即，单位时间人口的增长量与当时的人口成正比

模型建立

记x0x_0x0为初始时刻的人口，即x(0)=x0x(0)=x_0x(0)=x0
则从ttt到t+Δtt+\Delta tt+Δt内人口的增长量为：
x(t+Δt)−x(t)=rx(t)Δtx(t+\Delta t)-x(t)=rx(t)\Delta tx(t+Δt)−x(t)=rx(t)Δt
可导出下面的微分方程：
{dxdt=rxx(0)=x0\small \begin{cases} \frac{dx}{dt}=rx\\ x(0)=x_0 \end{cases} {dtdx=rxx(0)=x0
解得：x(t)=x0ert,r>0x(t)=x_0e^{rt} ,r>0x(t)=x0ert,r>0

模型二：Logistic阻滞增长模型

模型假设

1.同模型一
2.当人口增加到一定数量后，增长率随着人口的继续增长而逐渐减少，且r(x)为x的线性函数r(x)=r−sxr(x)=r-sxr(x)=r−sx，其中r为x=0时的增长率，成为固有增长率。
3.自然资源和环境条件所能容纳的最大人口数量为xmx_mxm，称作最大人口容量

模型建立

当x=xmx=x_mx=xm时，增长率为0，即r(xm)=0r(x_m)=0r(xm)=0,进而：S=rxmS=\frac{r}{x_m}S=xmr,所以：r(x)=r(1−xxm)r(x)=r(1-\frac{x}{x_m})r(x)=r(1−xmx)
其中的r,xmr,x_mr,xm是根据人口统计数据确定的常数，xmx_mxm常由经验决定(模型缺点：xmx_mxm不易轻易地得到)。
模仿模型一可得：
{dxdt=r(1−xxm)xx(0)=x0\small \begin{cases} \frac{dx}{dt}=r(1-\frac{x}{x_m})x\\ x(0)=x_0 \end{cases} {dtdx=r(1−xmx)xx(0)=x0
解得：x(t)=xm1+(xmx0)e−rtx(t)=\frac{x_m}{1+(\frac{x_m}{x_0})e^{-rt}}x(t)=1+(x0xm)e−rtxm

2.9药物在体内地分布与排除

模型背景

药物进入机体后，在随血液输送到各器官和组织的过程中，不断地被吸收、分布、代谢，最终排出体外。药物在血液中的浓度(μg/mv\mu g/mvμg/mv)称血药浓度。血药浓度的大小直接影响到药物的疗效，浓度太低不能达到预期的治疗效果，浓度太高又可能导致中毒、副作用太强或造成浪费。
因此，研究药物在体内吸收、分布和排除的动态过程，对于新药研制时剂量的确定、给药方案设计等药理学和临床医学的发展具有重要的指导意义和实用价值。
为了研究目的，将一个机体划分成若干个房室，每个房室是机体的一部分，比如中心室和周边室。在一个房室内药物呈均匀分布，而在不同的房室之间按一定规律进行转移。

模型假设

1.药物进入机体后，全部进入中心室（血液较丰富的心、肺、肾等器官和组织），中心室的容积在给药过程中保持不变；
2.药物从中心室排出体外，与排除的数量相比，药物的吸收可以忽略；
3.药物排除的速率与中心室的血药浓度成正比。

模型建立

f0(t)：给药速度f_0(t)：给药速度f0(t)：给药速度
c(t)：中心室血药浓度c(t)：中心室血药浓度c(t)：中心室血药浓度
x(t)：中心室药量x(t)：中心室药量x(t)：中心室药量
V：中心室容积V：中心室容积V：中心室容积
k：排除速率系数k：排除速率系数k：排除速率系数

求各种给药方式下血药浓度变化情况

上述变量有如下关系:
x˙=f0(t)−kx\dot{x}=f_0(t)-kxx˙=f0(t)−kx
即x˙+kx=f0(x)\dot{x}+kx=f_0(x)x˙+kx=f0(x)
又x(t)=Vc(t)x(t)=Vc(t)x(t)=Vc(t)
可得：c˙(t)+kc(t)=f0(t)V（1）\dot{c}(t)+kc(t)=\frac{f_0(t)}{V} （1）c˙(t)+kc(t)=Vf0(t)（1）

1.快速静脉注射(指数模型)
给药量为D，则初始条件c(0)=DV,f0(t)=0c(0)=\frac{D}{V},f_0(t)=0c(0)=VD,f0(t)=0
解得：c(t)=DVe−kt（2）c(t)=\frac{D}{V}e^{-kt} （2）c(t)=VDe−kt（2）
2.恒速静脉注射
设持续时间为τ\tauτ,注射速率为k0k_0k0，则有：
当(x≤t≤τ)(x\le t\le\tau)(x≤t≤τ)时，f0(t)=k0f_0(t)=k_0f0(t)=k0,初始条件为c(0)=0,
当(t≥τ)(t\ge\tau)(t≥τ)时，f0(t)=0f_0(t)=0f0(t)=0,初始条件c(τ)=k0Vk(1−e−kτ)c(\tau)=\frac{k_0}{Vk}(1-e^{-k\tau})c(τ)=Vkk0(1−e−kτ),所以（1）的解为：
c(t)={k0Vk(1−e−kt),0≤t≤τk0Vk(1−e−kτ)e−k(t−τ),t≥τ（3）\small c(t)= \begin{cases} \frac{k_0}{Vk}(1-e^{-kt}) ,0\le t\le\tau \\ \frac{k_0}{Vk}(1-e^{-k\tau})e^{-k(t-\tau)},t\ge\tau （3） \end{cases} c(t)={Vkk0(1−e−kt),0≤t≤τVkk0(1−e−kτ)e−k(t−τ),t≥τ（3）
3.口服或肌肉注射
在药物输入中心室之前先有一个将药物吸入血液的过程，可以看作有一个吸收室，药物由吸收室进入中心室，药物由吸收室进入中心室额转移速率系数记为k1k_1k1，给药量D，吸收室药量x0(t)x_0(t)x0(t)。则有：

{x0˙=−k1x0x0(t)=D\small \begin{cases} \dot{x_0}=-k_1x_0\\ x_0(t)=D \end{cases} {x0˙=−k1x0x0(t)=D
上式可推出：x0(t)=De−k1tx_0(t)=De^{-k_1t}x0(t)=De−k1t
于时f0(t)=k1De−k1tf_0(t)=k_1De^{-k_1t}f0(t)=k1De−k1t初始条件c(0)=0,（1）的解为：
c(t)=k1DV(k1−k)(e−kt−e−k1t),k1>k（4）c(t)=\frac{k_1D}{V(k_1-k)}(e^{-kt}-e^{-k_1t}) , k_1>k （4）c(t)=V(k1−k)k1D(e−kt−e−k1t),k1>k（4）

2.10导弹跟踪

背景

在发射导弹时刻（t=0），导弹位于坐标原点(0,0)，飞机位于(a,b)，飞机研平行于x轴的方向以常速v0v_0v0飞行。导弹在时刻t的位置为点(x,y)，其速度为常值v1v_1v1，导弹在飞行过程中，按照制导系统时钟指向飞机。请确定导弹的飞行轨迹以及击中飞机所需的时间T。

模型建立与求解

首先建立导弹的运动方程。导弹飞行曲线在点M(x,y)处的切线方程为：
Y−y=dydx(X−x)=dydtdxdt(X−x)Y-y=\frac{dy}{dx}(X-x)=\frac{\frac{dy}{dt}}{\frac{dx}{dt}}(X-x)Y−y=dxdy(X−x)=dtdxdtdy(X−x)
其中(x,y)为切线上动点的坐标。由于点A(xA,b)A(x_A,b)A(xA,b)应位于切线上，且xA=a+v0tx_A=a+v_0txA=a+v0t,所以：
b−y=dydtdxdt(a+v0t−x)b-y=\frac{\frac{dy}{dt}}{\frac{dx}{dt}}(a+v_0t-x)b−y=dtdxdtdy(a+v0t−x)
从而，导弹的飞行轨迹为：
{dxdt(b−y)=dydt(a+v0−x)（1）(dxdt)2+(dydt)2=v12（2）\small \begin{cases} \frac{dx}{dt}(b-y)=\frac{dy}{dt}(a+v_0-x) （1）\\ (\frac{dx}{dt})^2+(\frac{dy}{dt})^2=v_1^2 （2） \end{cases} {dtdx(b−y)=dtdy(a+v0−x)（1）(dtdx)2+(dtdy)2=v12（2）
由（1）可得：dxdy(b−y)=a+v0−x\frac{dx}{dy}(b-y)=a+v_0-xdydx(b−y)=a+v0−x
两边对t求导，得：d2xd2ydydt(b−y)−dxdydydt=v0−dxdt\frac{d^2x}{d^2y}\frac{dy}{dt}(b-y)-\frac{dx}{dy}\frac{dy}{dt}=v_0-\frac{dx}{dt}d2yd2xdtdy(b−y)−dydxdtdy=v0−dtdx
即：d2xd2ydydt(b−y)=v0（3）\frac{d^2x}{d^2y}\frac{dy}{dt}(b-y)=v_0 （3）d2yd2xdtdy(b−y)=v0（3）
由（2）得：(dydt)2[1+(dxdtdydt)2]=v12(\frac{dy}{dt})^2[1+(\frac{\frac{dx}{dt}}{\frac{dy}{dt}})^2]=v_1^2(dtdy)2[1+(dtdydtdx)2]=v12
即：dydt=v1[1+(dxdy)2]12\frac{dy}{dt}=\frac{v_1}{[1+(\frac{dx}{dy})^2]^{\frac{1}{2}}}dtdy=[1+(dydx)2]21v1
代入（3）可得导弹得运动方程：
d2xd2y(b−y)=λ[1+(dxdy)2]12（4）\frac{d^2x}{d^2y}(b-y)=\lambda[1+(\frac{dx}{dy})^2]^{\frac{1}{2}} （4）d2yd2x(b−y)=λ[1+(dydx)2]21（4）
其中，λ=v0v1\lambda=\frac{v_0}{v_1}λ=v1v0，
又x(0)=0,x(b)=a+v0Tx(0)=0,x(b)=a+v_0Tx(0)=0,x(b)=a+v0T（在T时刻击中目标）（5）
接下来求（4）满足（5）的解：
设p=dxdyp=\frac{dx}{dy}p=dydx，则dpdy=d2xdy2\frac{dp}{dy}=\frac{d^2x}{dy^2}dydp=dy2d2x，
（4）可化为：dpdy(b−y)=λ(1+p2)12（6）\frac{dp}{dy}(b-y)=\lambda(1+p^2)^{\frac{1}{2}} （6）dydp(b−y)=λ(1+p2)21（6）
ln[p+(1+p2)12]=−λln(b−y)+c1（7）ln[p+(1+p^2)^{\frac{1}{2}}]=-\lambda ln(b-y)+c_1 （7）ln[p+(1+p2)21]=−λln(b−y)+c1（7）
（6）的初始条件为p(0)=abp(0)=\frac{a}{b}p(0)=ba，令c1=ln(kbλ)c_1=ln(kb^{\lambda})c1=ln(kbλ)，则：
ln(p+1+p2)=ln(b−y)−λ+ln(kbλ)=ln[kbλ(b−y)λ]ln(p+\sqrt{1+p^2})=ln(b-y)^{-\lambda}+ln(kb^\lambda)=ln[\frac{kb^\lambda}{(b-y)^\lambda}]ln(p+1+p2)=ln(b−y)−λ+ln(kbλ)=ln[(b−y)λkbλ]
p+1+p2=kbλ(b−y)λp+\sqrt{1+p^2}=\frac{kb^\lambda}{(b-y)^\lambda}p+1+p2=(b−y)λkbλ
于时可以得到降阶方程：
dxdy=12[kbλ(b−y)λ−(b−y)λkbλ]\frac{dx}{dy}=\frac{1}{2}[\frac{kb^\lambda}{(b-y)^\lambda}-\frac{(b-y)^\lambda}{kb^\lambda}]dydx=21[(b−y)λkbλ−kbλ(b−y)λ]
其通解为：
x=12[kbλ(λ−1)(b−y)(λ−1)+(b−y)λ+1(λ+1)kbλ]+c（8）x=\frac{1}{2}[\frac{kb^\lambda}{(\lambda-1)(b-y)^{(\lambda-1)}}+\frac{(b-y)^{\lambda+1}}{(\lambda+1)kb^\lambda}]+c （8）x=21[(λ−1)(b−y)(λ−1)kbλ+(λ+1)kbλ(b−y)λ+1]+c（8）
根据初始条件x(0)=0,可得：
c=b[(1+k2)λ+k2−1]/2k(1−λ2)c=b[(1+k^2)\lambda+k^2-1]/2k(1-\lambda^2)c=b[(1+k2)λ+k2−1]/2k(1−λ2)
所以导弹飞行轨迹方程为：x=12[kbλ(λ−1)(b−y)(λ−1)+(b−y)λ+1(λ+1)kbλ]+b[(1+k2)λ+k2−1]/2k(1−λ2)x=\frac{1}{2}[\frac{kb^\lambda}{(\lambda-1)(b-y)^{(\lambda-1)}}+\frac{(b-y)^{\lambda+1}}{(\lambda+1)kb^\lambda}]+b[(1+k^2)\lambda+k^2-1]/2k(1-\lambda^2)x=21[(λ−1)(b−y)(λ−1)kbλ+(λ+1)kbλ(b−y)λ+1]+b[(1+k2)λ+k2−1]/2k(1−λ2)
又由x(b)=a+v0Tx(b)=a+v_0Tx(b)=a+v0T得到导弹集中目标的时间为：
T=c−av0=a2+b2−aλv1(1−λ2)（10）T=\frac{c-a}{v_0}=\frac{\sqrt{a^2+b^2}-a\lambda}{v_1(1-\lambda^2)} （10）T=v0c−a=v1(1−λ2)a2+b2−aλ（10）

2.11食饵-捕食者系统

一个包含两个群体的系统，其中一个群体紧密地依赖于另一个群体，成为食饵-捕食者系统。假设：x(t):t时刻食饵的数量；y(t):t时刻捕食者的数量。
如果各自独立生活，则：
{dxdt=λxdydt=−μy(λ,μ>0)\small \begin{cases} \frac{dx}{dt}=\lambda x\\ \frac{dy}{dt}=-\mu y (\lambda,\mu>0) \end{cases} {dtdx=λxdtdy=−μy(λ,μ>0)
如今两者生活在一起，则有：
{dxdt=(λ−αy)x（1）dydt=−(μ−βx)y（2）(α,β>0)\small \begin{cases} \frac{dx}{dt}=(\lambda-\alpha y)x （1）\\ \frac{dy}{dt}=-(\mu-\beta x)y （2）(\alpha,\beta>0) \end{cases} {dtdx=(λ−αy)x（1）dtdy=−(μ−βx)y（2）(α,β>0)
上式称为Volterra-Lotka方程，初始条件为x(0)=x0,y(0)=y0x(0)=x_0,y(0)=y_0x(0)=x0,y(0)=y0(1)/(2)可得：
dydx=(βx−μ)y(λ−αy)x\frac{dy}{dx}=\frac{(\beta x-\mu)y}{(\lambda-\alpha y)x}dxdy=(λ−αy)x(βx−μ)y
可得通解：
−αy−βx+λlny+μlnx=lnc-\alpha y-\beta x+\lambda lny+\mu lnx=lnc−αy−βx+λlny+μlnx=lnc或yλeαyxμeβx=c\frac{y^\lambda}{e^{\alpha y}}\frac{x^\mu}{e^{\beta x}}=ceαyyλeβxxμ=c
将初始条件代入，可得特解，是xoy面上的一条闭轨线

当食饵较多时，捕食者增多因而食饵必定减少，使得捕食者也随之减少，从而食饵又会增多。两者的数量如此起伏，周而复始，维持着生态平衡。

2.12传染病模型

背景

建立传染病模型的目的是描述传染过程、分析受感染人数的变化规律、预报高潮来到的时间。
为了简单起见，假设传播期间内所观察地区人数N不变，不计生死迁移，时间以天为单位。

模型一 SI模型

模型假设

1.人群分为健康者和病人，在t时刻这两类人所占比例分别为s(t),i(t)，即s(t)+i(t)=1;
2.平均每个病人每天接触人数是常数λ\lambdaλ，即每个病人平均每天使得λs(t)\lambda s(t)λs(t)个健康者受感染变成病人，λ\lambdaλ称为日接触率。

模型建立

根据模型假设2，在T时刻，每个病人每天可以使得λs(t)\lambda s(t)λs(t)个健康者变成病人，病人人属为Ni(t)，故每天新增λNs(t)i(t)\lambda Ns(t)i(t)λNs(t)i(t)个患者，即：
Ndidt=λNsiN\frac{di}{dt}=\lambda NsiNdtdi=λNsi，假设t=0时患者比例i0i_0i0，可得模型：
{didt=λi(1−i)（1）i(0)=i0\small \begin{cases} \frac{di}{dt}=\lambda i(1-i) （1）\\ i(0)=i_0 \end{cases} {dtdi=λi(1−i)（1）i(0)=i0
式（1）的解为：i(t)=11+(1t0e−λt)（2）i(t)=\frac{1}{1+(\frac{1}{t_0}e^{-\lambda t})} （2）i(t)=1+(t01e−λt)1（2）

模型解释

1.当i=12i=\frac{1}{2}i=21时，didt\frac{di}{dt}dtdi达到最大值，此时t=tm=λ−1ln(1i0−1)t=t_m=\lambda^{-1}ln(\frac{1}{i_0}-1)t=tm=λ−1ln(i01−1)，也就是说，高潮到来时，λ\lambdaλ越大，则tmt_mtm越小。
2.当t→∞t\rightarrow\inftyt→∞时，i→1i\rightarrow 1i→1此时所有的人都被感染，因为SI模型没有考虑治愈病人。

模型二 SIS模型

在SI模型的基础上引入治愈，对SI模型进行修正。

模型假设

1.同SI模型假设1
2.同SI模型假设2
3.病人每天被治愈的占病人总数的比例为μ\muμ，称作日治愈率。

模型修正

SI模型可修正为，t时刻每天有NiμNi\muNiμ的病人转变为健康。
{didt=λi(1−i)−μi（3）i(0)=i0\small \begin{cases} \frac{di}{dt}=\lambda i(1-i)-\mu i （3）\\ i(0)=i_0 \end{cases} {dtdi=λi(1−i)−μi（3）i(0)=i0
（3）的解为：
i(t)={[λλ−μ+(1i0−λλ−μ)e−(λ−μ)t]−1,λ≠μ(λt+1i0)−1,λ=μ（4）\small i(t)= \begin{cases} [\frac{\lambda}{\lambda-\mu}+(\frac{1}{i_0}-\frac{\lambda}{\lambda-\mu})e^{-(\lambda-\mu)t}]^{-1},\lambda\not ={\mu}\\ (\lambda t+\frac{1}{i_0})^{-1},\lambda=\mu （4） \end{cases} i(t)={[λ−μλ+(i01−λ−μλ)e−(λ−μ)t]−1,λ=μ(λt+i01)−1,λ=μ（4）
由（3）可以计算出使得didt\frac{di}{dt}dtdi达到最大值的高潮时刻tmt_mtm（didt\frac{di}{dt}dtdi的最大值(didt)m(\frac{di}{dt})_m(dtdi)m在i=λ−μ2λi=\frac{\lambda-\mu}{2\lambda}i=2λλ−μ时达到）
记a=λμa=\frac{\lambda}{\mu}a=μλ，可知：
i(∞)={1−1a,a>10,a≤1\small i(\infty)= \begin{cases} 1-\frac{1}{a},a>1\\ 0,a\le 1 \end{cases} i(∞)={1−a1,a>10,a≤1

SIR模型

模型假设

1.人群分为健康者、病人、移出者（病愈免疫者），三类人在t时刻在总人数N中占比例分别为s(t)、i(t)、r(t)，即s(t)+i(t)+r(t)=1
2.病人日接触率为λ\lambdaλ，日治愈率为μ\muμ，传染期间接触数σ=λμ\sigma=\frac{\lambda}{\mu}σ=μλ

模型建立

i(t)随t的变化规律同模型二，对于r(t)：
Ndrdt=μNi,且dsdt+didt+drdt=0N\frac{dr}{dt}=\mu Ni,且\frac{ds}{dt}+\frac{di}{dt}+\frac{dr}{dt}=0Ndtdr=μNi,且dtds+dtdi+dtdr=0
于时可得模型：
{dsdt=−λsididt=λsi−μi（5）s(0)=s0,i(0)=i0\small \begin{cases} \frac{ds}{dt}=-\lambda si\\ \frac{di}{dt}=\lambda si-\mu i （5）\\ s(0)=s_0,i(0)=i_0 \end{cases} ⎩⎨⎧dtds=−λsidtdi=λsi−μi（5）s(0)=s0,i(0)=i0
从（5）中消去dt，结合σ\sigmaσ的实际意义，可得：
{dids=1σs−1（6）i∣s=s0=i0\small \begin{cases} \frac{di}{ds}=\frac{1}{\sigma s}-1 （6）\\ i|_{s=s_0}=i_0 \end{cases} {dsdi=σs1−1（6）i∣s=s0=i0
（6）的解为：
i=(s0+i0)−s+1σlnss0（7）i=(s_0+i_0)-s+\frac{1}{\sigma}ln\frac{s}{s_0} （7）i=(s0+i0)−s+σ1lns0s（7）
根据（5）（7）以及图像可分析s(t),i(t),r(t)的变化规律：

1.无论s0,i0s_0,i_0s0,i0为多少，i∞=0i_\infty=0i∞=0，即病人终将消失。
2.最终未被感染的健康者比例s∞s_\inftys∞时方程s0+i0−s∞+1σlns∞s0=0（8）s_0+i_0-s_\infty+\frac{1}{\sigma}ln\frac{s_\infty}{s_0}=0 （8）s0+i0−s∞+σ1lns0s∞=0（8）在(0,1σ)(0,\frac{1}{\sigma})(0,σ1)内的单根。
3.若s0>1σs_0>\frac{1}{\sigma}s0>σ1，则当s=1σs=\frac{1}{\sigma}s=σ1时，i(t)达到最大值im=s0+i0−1σ(1+lnσs0)i_m=s_0+i_0-\frac{1}{\sigma}(1+ln\sigma s_0)im=s0+i0−σ1(1+lnσs0)，i(t)先增后减至0.
4.若s0≤1σs_0\le\frac{1}{\sigma}s0≤σ1，则i(t)→0,s(t)→s∞i(t)\rightarrow0,s(t)\rightarrow s_\inftyi(t)→0,s(t)→s∞。

模型解释

1.1σ\frac{1}{\sigma}σ1是一个阈值，当s0>1σs_0>\frac{1}{\sigma}s0>σ1时传染病会蔓延，s0≤1σs_0\le\frac{1}{\sigma}s0≤σ1时就不会蔓延
2.σ=λμ\sigma=\frac{\lambda}{\mu}σ=μλ表示λ\lambdaλ越小，μ\muμ越大，σ\sigmaσ也越小，从而越有利。

数学建模笔记CH2：微积分方法建模相关推荐

数学建模笔记——插值拟合模型（二）
今天是8月21日,距离上次写文章好像将近一个月了--这段时间经历了建模校内选拔赛,考试周,以及与网络小说的斗智斗勇--好吧,其实也没干什么,除了考试就是荒废-- 我最近有在思考一个问题,就是我所关注的 ...
数学建模中常用的方法
数学建模中常用的方法:类比法.二分法.差分法.变分法.图论法.层次分析法.数据拟合法.回归分析法.数学规划(线性规划,非线性规划,整数规划,动态规划,目标规划).机理分析.排队方法.对策方法.决策方法 ...
数学建模的影响因素分析方法
数学建模的影响因素分析方法 PCA主成分分析灰色关联分析 AHP层次分析小结如果赶时间可直接看小结部分,再返回看正文作为萌新参加了数学建模,为解决影响煤炭价格的主要因素的问题, ...
[数学建模]学习笔记1：初等建模
初等模型: 1.研究对象的机理比较简单 2.用静态,线性,确定性模型即可达到建模的目的 3.可以利用初等数学方法来构造和求解模型注:尽量用简单的数学工具来建模 2.1 光盘的数据容量调查和分析经 ...
数学建模常见的一些方法【04拟合算法】
文章目录数学建模常见的一些方法 1. 拟合算法 1.1 插值和拟合的区别 1.2 求解最小二乘法 1.3 Matlab求解最小二乘 1.4 如何评价拟合的好坏 1.5 证明SST = SSE + S ...
数学建模常见的一些方法【03插值算法】
文章目录数学建模常见的一些方法 1. 插值算法 1.1 插值法的定义 1.2 插值法的分类 1.3 一般插值多项式原理 1.4 拉格朗日插值法 1.5 龙格现象(Runge phenomenon) ...
数学建模笔记之一起读论文2019年C题——机场的出租车问题
数学建模笔记之一起读论文--机场的出租车问题 2021-8-28 全国大学生数学建模竞赛 2019年C题 B站链接--国赛C题真题解析 1 赛题阅读与分析原题再现: 问题C 机场的出租车问题大多数 ...
【清风数学建模笔记】第七讲多元回归分析
回归分析是数据分析中最基础也是最重要的分析工具,绝大多数的数据分析问题,都可以使用回归的思想来解决.回归分析的任务就是,通过研究自变量X和因变量Y的相关关系,尝试去解释Y的形成机制,进而达到通过X去预 ...
建模笔记之maple学习
建模笔记之 maple 学习本笔记主要介绍基本的方程组求解,对于语法.画图功能不加以阐述.在数学建模中,大部分的编程工作还是由matlab或python来完成,而maple可以快速解决一些需要手算的 ...

数学建模笔记CH2：微积分方法建模

第二章：微积分方法建模

overview:

2.1飞机降落曲线

求曲线解析式

求最佳着陆点

2.2经济增长模型

道格拉斯（Douglas）生产函数

劳动生产率增长条件

2.3存贮模型

模型一

模型假设

模型建立

模型求解

模型二 允许缺货的存贮模型

模型假设

模型建立

模型求解

分析

2.4城市人口统计模型

模型一（估算城市现有人口）

模型二（预测城市未来人口）

2.7古生物年代确定

模型假设

模型建立

2.8预测人口的增长

模型一：Malthhus指数增长模型

模型假设

模型建立

模型二：Logistic阻滞增长模型

模型假设

模型建立

2.9药物在体内地分布与排除

模型背景

模型假设

模型建立

求各种给药方式下血药浓度变化情况

2.10导弹跟踪

背景

模型建立与求解

2.11食饵-捕食者系统

2.12传染病模型

背景

模型一 SI模型

模型假设

模型建立

模型解释

模型二 SIS模型

模型假设

模型修正

SIR模型

模型假设

模型建立

模型解释

数学建模笔记CH2：微积分方法建模相关推荐

最新文章

热门文章

模型二允许缺货的存贮模型