北京工业大学-数学建模2 -201310

微分方程实验作业一、 基本实验1. 微分方程稳定性分析绘出下列自治系统相应的轨线，并标出随t 增加的运动方向，确定平衡点， 并按稳定的、渐近稳定的、或不稳定的进行分类：,,,,(1)(2)(3)(4);2;2;2.dx dx dx dxx x y x dt dt dt dtdy dy dy dy y y x y dt dt dt dt ⎧⎧⎧⎧==-==-⎪⎪⎪⎪⎪⎪⎪⎪⎨⎨⎨⎨⎪⎪⎪⎪===-=-⎪⎪⎪⎪⎩⎩⎩⎩解：（1）易知12(),(,)f x y x f x y y ==，，则上述代数方程组的实根为x=0,y=0,可知微分方程的平衡点为（0,0）. 线性微分方程组的系数矩阵为1001A ⎡⎤=⎢⎥⎣⎦，求得特征值λ1=1>0，λ2=1>0；p=-(λ1+λ2)=-2<0,q=(λ1)* (λ2)=1>0,所以根据平衡点与稳定性的各种情况表可知，平衡点（0,0）是不稳定的。

由于2121df d x dx x x dt dx dt =⋅=⋅=，2221df d y dy y y dt dy dt=⋅=⋅=，可知''''(),()x t y t 在x=0和y=0处变号，且随着t 的增加远离远点。

根据x,y 一二阶导数的符号，图形（手画）如下1-1：图1-1（2）易知12(),(,)2f x y x f x y y =-=，，则上述代数方程组的实根为x=0,y=0,可知微分方程的平衡点为（0,0）. 线性微分方程组的系数矩阵为1002A -⎡⎤=⎢⎥⎣⎦，求得特征值λ1=1>0，λ2=-2<0；p=-(λ1+λ2)=-1<0,q=(λ1)* (λ2)=-2<0,所以根据平衡点与稳定性的各种情况表可知，平衡点（0,0）是不稳定的。

2121()df d x dx x x dt dx dt =⋅=-⋅-=，222224df d y dyy y dt dy dt=⋅=⋅=，根据x,y 一二阶导数的符号，可得如下图形1-2：图1-2（3）易知12(),(,)2f x y y f x y x ==-，，则上述代数方程组的实根为x=0,y=0,可知微分方程的平衡点为（0,0）. 线性微分方程组的系数矩阵为0120A ⎡⎤=⎢⎥-⎣⎦， p=-(λ1+λ2)=-（a 11+a 22）=0,q=(λ1)* (λ2)= a 11a 22- a 21a 12=2>0,所以根据平衡点与稳定性的各种情况表可知，平衡点（0,0）是不稳定的。

2122df d x dyx dt dt dt ===-,22222df d y dxy dt dt dt==-=-.根据x,y 一二阶导数的符号，可得如下图形1-3：图1-3(4) 易知12()1,(,)2f x y x f x y y =-+=-，，则上述代数方程组的实根为x=1,y=0,可知微分方程的平衡点为（1,0）. 线性微分方程组的系数矩阵为1002A -⎡⎤=⎢⎥-⎣⎦， p=-(λ1+λ2)=-（a 11+a 22）=3>0,q=(λ1)* (λ2)= a 11a 22- a 21a 12=2>0,所以根据平衡点与稳定性的各种情况表可知，平衡点（1,0）是稳定的,随着t 的增加，x,y 趋向稳定点（1,0）。

2121(1)1df d x dxx x dt dx dt =⋅=-⋅-+=-,222224df d y dyy y dt dy dt=⋅=-⋅-=.根据x,y 一二阶导数的符号，x <1时，'0,''0x x ><,x>1时，'0,''0x x <>；y>0时，'0,''0y y <>，y<0时，'0,''0y y ><，可得如下图形1-4：图1-42. 种群增长模型一个片子上的一群病菌趋向于繁殖成一个圆菌落.设病菌的数目为N ，单位成员的增长率为r1，则由Malthus 生长律有1dNr N dt=⋅，但是，处于周界表面的那些病菌由于寒冷而受到损伤，它们死亡的数量与N 1/2成比例，其比例系数为r2，求N 满足的微分方程.不用求解，图示其解族.方程是否有平衡解，如果有，是否为稳定的?解：由题意可得病菌的数目N 满足的微分方程为：1212()dNr N r N f N dt=-=，r 1为单位成员的增长率，r 2为比例系数。

令()f N =0，可得方程的两个平衡点为N1=0,N2=2221/r r ；对N 进行二阶求导得到：111222212122()1()()2d N df N r r N r N r N dt dt -==-⋅-⋅。

令220d Ndt=，可求得N 1=2221/4r r ，N 2=2221/r r ，N 3=0；则N=0，N=2221/4r r ，N=2221/r r 把第一象限（因为数量N 和时间t 都为正值）划为三部分，且从下到上三部分中分别有：0<N<2221/4r r 时,0,dN dt >220d N dt>，递增下凸函数； 2221/4r r <N<2221/r r 时,0,dN dt >220d N dt <，递增上凸函数； 2221/r r <N 时，0,dN dt <220d N dt>，为递减下凸函数。

则可以得到N （t ）的图形2-1：图2-1对于方程的两个平衡点，其中N 1=0是不稳定的；N 2=2221/r r 是稳定的。

3. 单种群开发模型考虑单种群开发方程： (1)().dx xr x Ex E r dt N=--< （1） 在不求解的情况下，绘出其解族曲线；（2） 用数学表达式证明：在稳定状态下，最优捕捞率*2r E =解：（1）令(1)=()'()0().dx xr x Ex f x x t E r dt N=--==< 可以求得120,()0Nx x r E r==->；且0()'()0;r x r E x t N <<->时，()'()0;rx r E x t N>-<时，求x 的二阶导数：22212()()''()().d x df dx r E r E r x x x x t E r dt dx dt r N r N--=⋅=--=<有直线0,()()2N Nx x r E x r E r r ==-=-和将第一象限划分为三部分。

且：0()''()0;2r x r E x t N <<-<时，()()''()0;2r r r E x r E x t N N-<<-<时， ()''()0;rx r E x t N >->时， '()''()x t x t ，的符号已经确定，可以得到解族曲线图3-1：图3-1由图像可以看出120,()Nx x r E r==-不是稳定点是稳定点。

（E<r ） （2）(1)=(),()xr x g x Ex f x N-=记 可知抛物线=()(1),xg x r x N =-y 两根为120,x x N ==，抛物线开口向下，在2Nx =处取得最大值记为m h 。

直线()f x Ex =的斜率为E,且抛物线=()(1),xg x r x N=-y 在x=0处的切线为y rx =，只有E<r ，存在()()f x g x 与的交点，即有稳定的平衡点。

稳定状态下：()-()=0=r(1-)x-Ex xf xg x N即：r(1-)x=Ex xN，()()f x g x 与相交（如图3-2）。

在2Nx =时达到最大捕捞量，此时： 1r(1-)-E =0222N N N N ⋅，求解得到：2r E =，最优捕捞率即为2r 。

图3-24. Compertz 模型设渔场鱼量自然增长服从Compertz 模型ln dx N rx dt x= 其中，r 为固有增长率，N 为最大的种群数量。

若单位时间内捕捞量为h Ex =，试讨论鱼量的平衡点及其稳定性，求最大持续产量m h ，此时捕捞强度m E 和鱼量水平*0x 。

解：（1）平衡点及其稳定性： ln =()defdx N rx g x dt x=, r 为固有增长率，N 为最大的种群数量。

鱼量x （t ）变化规律为：()()=f(x)defdxg x h x dt=-，令f(x)=0,可知： (ln )0N x r E x -=，平衡点即为12/0,E r N x x e ==。

/0E r N x e <<时，'()0x t <；/E r Nx e>时，'()0x t >。

求x 的二阶导数：''()(ln )(ln )N Nx t r E x r E x x =+⋅⋅-123/0,,E E r N x x x Ne e ===，由于x<=N 。

可知/0,,E r Nx x x N e===将第一象限分割成两部分。

/0E r N x e <<时，''()0x t >；/E r Nx N e<<时，''()0x t <。

根据'()''()x t x t ，的符号，可得到图4-1图4-1当初始条件02/E rNx x e<=，随着时间推移很快到达0，表明渔场中的鱼很快会被捕光。

当初始条件02/E r Nx x e>=时，渔场中的鱼量可以自动调节而趋于N 。

对于方程来说0是稳定的平衡点，2/E r Nx e=不是稳定的平衡点。

（2）最大持续产量、捕捞强度和鱼量水平：分析ln =()defN rx g x x，由于21'()ln ()ln N N Ng x r rx r r x x x x =+⋅⋅-=-'()0g x =，求解得到Nx e=。

当0Nx e<<，'()0g x >；N x N e <<时，'()0g x <；Nx e=时，'()0g x =； 所以，()g x 在N x e =时，取得极大值为ln N N rNr N e e e=。

(),()g x h x 交于Nx e =处便可得到平衡点，此时()=m rN Nh x h E e e==，求得E=r 。