前言

今天我们要说的就是利用MATLAB对常微分方程进行数值求解的实例应用-Lotka-Volterra模型及其改进模型。本文是科学计算与MATLAB语言专题六第6小节的学习笔记，如果大家有时间的话，可以去听听课，没有的话，可以看看我的笔记，还是很不错的。

1 Lotka−VolterraLotka-VolterraLotka−Volterra模型

20世纪40年代，Lotka（1925）和Volterra（1926）奠定了种间竞争关系的理论基础，他们提出的种间竞争方程对现代生态学理论的发展有着重大影响。
可以用一阶非线性微分方程组表示：
{drdt=2r−λrf，r(0)=r0dfdt=−f+λrf，f(0)=f0\left\{ \begin{aligned} \frac{dr}{dt}&=2r-\lambda rf，r(0)=r_0\\ \frac{df}{dt}&=-f+\lambda rf，f(0)=f_0 \end{aligned} \right. ⎩⎪⎪⎨⎪⎪⎧dtdrdtdf=2r−λrf，r(0)=r0=−f+λrf，f(0)=f0

其中，ttt为时间，r（t）r（t）r（t）为兔子数量，f（t）f（t）f（t）为狐狸数量，λ\lambdaλ为一个正常数。该系统解具有周期性，其周期取决于初始条件。也就是说，对于任意数量的r(0)r(0)r(0)和f(0)f(0)f(0)，总存在一个时间t=tpt=t_pt=tp，使这两个种群的数量回归初始值。
①在r0=300、f0=150、λ=0.01r_0=300、f_0=150、\lambda =0.01r0=300、f0=150、λ=0.01的条件下，求该系统的解，结果会发现tpt_ptp接近5。进一步绘制r(t)r(t)r(t)和f(t)f(t)f(t)函数的曲线，以及以r和f为坐标轴画相平面图。
相平面相关概念：
设一个二阶系统可以用下列常微分方程来描述：
d2xdt2+a1(x,dxdt)dxdt+a0(x,dxdt)x=0，x¨=f(x,x˙)\frac{d^2x} {dt^2}+a_1(x,\frac{dx}{dt})\frac{dx}{dt}+a_0(x,\frac{dx}{dt})x=0 ，\ddot{x}=f(x,\dot{x}) dt2d2x+a1(x,dtdx)dtdx+a0(x,dtdx)x=0，x¨=f(x,x˙)
令x=x1，dx/dt=x2x=x1，dx/dt=x2x=x1，dx/dt=x2
dx2dx1=f(x1,x2)x2\frac{dx_2}{dx_1}=\frac{f(x_1,x_2)}{x_2}dx1dx2=x2f(x1,x2)
{dx1dt=x2dx2dt=f(x1,x2)\left\{ \begin{aligned} \frac{dx_1}{dt}&=x_2\\ \frac{dx_2}{dt}&=f(x_1,x_2)\\ \end{aligned} \right. ⎩⎪⎪⎨⎪⎪⎧dtdx1dtdx2=x2=f(x1,x2)
以x1x_1x1为自变量，以x2x_2x2为因变量的一阶微分方程。二阶系统常微分方程方程的解既可用x与t的关系来表示，也可用x2x_2x2与x1x_1x1的关系来表示。实际上，看作一个质点的运动方程，则x1(t)代表质点的位置，x2(t)x_2(t)x2(t)代表质点的速度。以x1x_1x1为自变量，以x2x_2x2为因变量的一阶微分方程。二阶系统常微分方程方程的解既可用x与t的关系来表示，也可用x，与x1的关系来表示。实际上，看作一个质点的运动方程，则x1(t)x_1(t)x1(t)代表质点的位置，x2(t)x_2(t)x2(t)代表质点的速度。
用x1、x2x_1、x_2x1、x2，描述二阶系统常微分方程方程的解，也就是用质点的状态来表示该质点的运动。在物理学中，状态又称为相。把由x1，x2x_1，x_2x1，x2，所组成的平面坐标系称为相平面，系统的一个状态则对应于相平面上的一个点。
当t变化时，系统状态在相平面上移动的轨迹称为相轨迹。
而与不同初始状态对应的一簇相轨迹所组成的图叫做相平面图。

②在r0=15、f0=22、λ=0.01r_0=15、f_0=22、\lambda=0.01r0=15、f0=22、λ=0.01时,求解并绘制图形，发现会发现tpt_ptp接近6.62。
③在r0=102、f0=198、λ=0.01r_0=102、f_0=198、\lambda=0.01r0=102、f0=198、λ=0.01时求解并绘制图形，并确定周期tpt_ptp。
④点（r0，f0）=（1/λ，2/λ）（r_0，f_0）=（1/\lambda，2/\lambda）（r0，f0）=（1/λ，2/λ）是一个稳定平衡点。若以此为初值，那么种群数量不发生变化。若初值接近此平衡点，那么数量不会发生大的变化。试作图验证。
模型分析
x1(t)x_1(t)x1(t):兔子在t时刻的数量
x2(t）x_2(t）x2(t）一狐狸在t时刻的数量
r1r_1r1一兔子独立生存时的增长率
r2r_2r2一狐狸独自存在时的死亡率
λ1\lambda_1λ1，一狐狸掠取兔子的能力
λ2\lambda_2λ2一兔子对狐狸的供养能力
这里假设r1=2，r2=1，λ1=λ2=0.01r1=2，r2=1，\lambda_1=\lambda_2=0.01r1=2，r2=1，λ1=λ2=0.01，即系统模型如下：
{dx1dt=x1(2−0.01x2)dx2dt=x2(−1+0.01x1)\left\{ \begin{aligned} \frac{dx_1}{dt}&=x_1(2-0.01x_2)\\ \frac{dx_2}{dt}&=x_2(-1+0.01x_1) \end{aligned} \right. ⎩⎪⎪⎨⎪⎪⎧dtdx1dtdx2=x1(2−0.01x2)=x2(−1+0.01x1)

第①问：兔子数量初始值x1（0）=300x_1（0）=300x1（0）=300，狐狸数量初始值x2(0)=150x_2(0)=150x2(0)=150。

rabbitFox=@(t,x)[x(1)*(2-0.01*x(2));x(2)*(-1+0.01*x(1))];
[t,x]=ode45(rabbitFox,[0,30],[300,150])
subplot(1,2,1);
plot(t,x(:,1),'-',t,x(:,2),'-*');
legend('xl（t）','x2（t）');
xlabel('时间');
ylabel('物种数量');
grid on
subplot(1,2,2);
plot(x(:,1),x(:,2))
grid on

我们可以看到，兔子和狐狸两者相互制约，当狐狸逐渐增多时，兔子逐渐减少，当狐狸增加到一定数量时，由于种群内部竞争激烈，又导致其数量减少，从而兔子天敌减少，进而导致兔子数量增加，当兔子增加到一定数量时，种间竞争加剧，从而又导致兔子数量减少，如此循环，相互制约发展。
第②问：兔子初始值x1（0）=15x_1（0）=15x1（0）=15，狐狸初始值x2（0）=22x_2（0）=22x2（0）=22。

rabbitFox=@(t,x)[x(1)*(2-0.01*x(2));x(2)*(-1+0.01*x(1))];
[t,x]=ode45(rabbitFox,[0,30],[15,22])
subplot(1,2,1);
plot(t,x(:,1),'-',t,x(:,2),'-*');
legend('xl（t）','x2（t）');
xlabel('时间');
ylabel('物种数量');
grid on
subplot(1,2,2);
plot(x(:,1),x(:,2))
grid on

第③问：兔子初始值x1（0）=102x_1（0）=102x1（0）=102，狐狸初始值x2（0）=198x_2（0）=198x2（0）=198。

rabbitFox=@(t,x)[x(1)*(2-0.01*x(2));x(2)*(-1+0.01*x(1))];
[t,x]=ode45(rabbitFox,[0,30],[102,198])
subplot(1,2,1);
plot(t,x(:,1),'-',t,x(:,2),'-*');
legend('xl（t）','x2（t）');
xlabel('时间');
ylabel('物种数量');
grid on
subplot(1,2,2);
plot(x(:,1),x(:,2))
grid on

第④问：验证（1/λ，2/λ）（1/\lambda，2/\lambda）（1/λ，2/λ）是稳定平衡点。
取λ\lambdaλ=0.01，所以稳定平衡点（1/λ,2/λ）（1/\lambda,2/\lambda）（1/λ,2/λ）即是(100,200)(100,200)(100,200)，以此点作为初值时，画出其图像。

rabbitFox=@(t,x)[x(1)*(2-0.01*x(2));x(2)*(-1+0.01*x(1))];
[t,x]=ode45(rabbitFox,[0,30],[100,200])
plot(t,x(:,1),'-o',t,x(:,2),'-*');
legend('x1（t）-兔子','x2（t）-狐狸');
xlabel('时间');
ylabel('物种数量');

当将初始值变为(98,195)(98,195)(98,195)时，即向下十分接近平衡点，画出其图像。

rabbitFox=@(t,x)[x(1)*(2-0.01*x(2));x(2)*(-1+0.01*x(1))];
[t,x]=ode45(rabbitFox,[0,30],[98,195])
plot(t,x(:,1),'-o',t,x(:,2),'-*');
legend('x1（t）-兔子','x2（t）-狐狸');
xlabel('时间');
ylabel('物种数量');

两者图像没有交点，兔子和狐狸的数量分别以各自的平衡点进行周期性波动，但是两者的数量变化不是很大。
当将初始值变为(70,150)(70,150)(70,150)时（向下偏离平衡点比较远时），画出其图像。

可以发现兔子和狐狸的数量变化比较剧烈，但是依然成周期性变化。
当将初始值变为(900,1600)(900,1600)(900,1600)时（向上偏离平衡点十分远时），画出其图像。

rabbitFox=@(t,x)[x(1)*(2-0.01*x(2));x(2)*(-1+0.01*x(1))];
[t,x]=ode45(rabbitFox,[0,500],[900,1560])
plot(t,x(:,1),t,x(:,2));
legend('x1（t）-兔子','x2（t）-狐狸');
xlabel('时间');
ylabel('物种数量');

可以发现，此时兔子和狐狸的数量变化相当的剧烈，但是依然存在周期。

2 Lotka−VolterraLotka-VolterraLotka−Volterra改进模型

修改后的方程：
{drdt=2（1−rR)−λrf，r(0)=r0dfdt=−f+λrf，f(0)=f0\left\{ \begin{aligned} \frac{dr}{dt}&= 2（1-\frac{r}{R})-\lambda rf，r(0)=r0 \\ \frac{df}{dt}&= -f+\lambda rf，f(0)=f0 \end{aligned} \right. ⎩⎪⎪⎨⎪⎪⎧dtdrdtdf=2（1−Rr)−λrf，r(0)=r0=−f+λrf，f(0)=f0
其中，t为时间；r（t）为兔子数量；f（t）为狐狸数量；λ\lambdaλ为一个正常数，R也是一个正常数。方程在r0=300、f0=150r0=300、f0=150r0=300、f0=150初始条件下的50个时间单位上求解并绘制图形。
①原模型下，狐狸数量和兔子数量的时间函数曲线
②改进模型下，狐狸数量和兔子数量的时间函数曲线。
③原模型下，狐狸数量相对兔子数量的关系曲线。
④改进模型下，狐狸数量相对兔子数量的关系曲线。
第①问：在原模型下，狐狸数量和兔子数量的时间函数曲线

rabbitFox=@(t,x)[x(1)*(2-0.01*x(2));x(2)*(-1+0.01*x(1))];
[t,x]=ode45(rabbitFox,[0,50],[300,150])
plot(t,x(:,1),t,x(:,2));
legend('x1（t）-兔子','x2（t）-狐狸');
xlabel('时间');
ylabel('物种数量');
title('原模型下，狐狸和兔子数量的函数曲线');

第②问：在改进模型下，狐狸数量和兔子数量的时间函数曲线。

rabbitFox=@(t,x) [2*x(1)*(1-x(1)/400-0.005*x(2));...x(2)*(-1+0.01*x(1))];%...代表换行
[t,x]=ode45(rabbitFox,[0,50],[300,150]);
plot(t,x(:,1),t,x(:,2));
legend('x1（t）-兔子','x2（t）-狐狸');
xlabel('时间');
ylabel('物种数量');
title('改进模型下，狐狸和兔子数量的函数曲线');

第③问：在原模型下，狐狸数量相对兔子数量的关系曲线。

rabbitFox=@(t,x)[x(1)*(2-0.01*x(2));x(2)*(-1+0.01*x(1))];
[t,x]=ode45(rabbitFox,[0,30],[300,150])
plot(x(:,1),x(:,2));
xlabel('兔子数量');
ylabel('狐狸数量');
title('原模型下，狐狸数量相对于兔子数量的关系曲线');

第④问：在改进模型下，狐狸数量相对于兔子数量的关系曲线。

rabbitFox=@(t,x) [2*x(1)*(1-x(1)/400-0.005*x(2));...x(2)*(-1+0.01*x(1))];%...代表换行
[t,x]=ode45(rabbitFox,[0,50],[300,150]);
plot(x(:,1),x(:,2));
xlabel('时间');
ylabel('物种数量');
title('改进模型下，狐狸数量相对于兔子数量的关系曲线');

模型的进一步改进：
考虑狐狸的环境容纳量
考虑到兔子其他天敌和狐狸天敌的影响
考虑自然灾害、人为捕捉等因素

结语

欢迎大家点赞

如何利用MATLAB建立Lotka-Volterra模型及其改进模型相关推荐

MATLAB/Simulink——利用S-Function建立高超声速飞行器的纵向模型
参考文献:朱平. 高超声速飞行器容错控制算法研究[D].南京航空航天大学,2020. 一.高超声速飞行器的纵向模型气动力以及力矩表达式为: ps:其中具体参数见论文or下面的S函数二.利用S-Fu ...
利用matlab建立含间隙关节的曲柄滑块动力学分析
好的. 建立含间隙关节的曲柄滑块动力学分析可以使用 Matlab 的机械设计工具箱 (Mechanical Toolbox) 进行.你需要了解曲柄滑块机构的几何特征和动力学模型,然后用 Matlab ...
怎样用matlab建立igbt的仿真分析模型,基于MATLAB/Simulink的IGBT导通模型研究
<电气自动化)2o15年第37卷第6期电力系统及其自动· Power System & Automation 基于 MATLAB/Simu¨nk的 IGBT导通模型研究朱永超, ...
matlab里的仿真模型块,搭建simulink模型（如何利用MATLAB/SIMULINK搭建简单的仿真模型）...
如何利用MATLAB/SIMULINK搭建简单的仿真模型安装完MATLAB软件后,在电脑桌面点击MATLAB快捷方式打开MATLAB后,点击Simulink Library按钮之后会进入Simu ...
MATLAB建立ar模型,matlab关于ar模型
基于参数建摸的功率谱估计是现代功率谱估计的重要内容,其目的就是为了改善功率谱估计的频率分辨率,它主要包括 AR 模型.MA 模型.ARMA 模型,其中基于 AR 模型...... MATLAB 仿 ...
利用 AssemblyAI 在 PyTorch 中建立端到端的语音识别模型
作者 | Comet 译者 | 天道酬勤,责编 | Carol 出品 | AI 科技大本营(ID:rgznai100) 这篇文章是由AssemblyAI的机器学习研究工程师Michael Nguyen ...
matlab建立的发动机的模型,基于MATLAB∕Simulink的摩托车发动机仿真模型建立.pdf
·信息技术· 王磊,等·基于MATLAB/simulink的摩托车发动机仿真模型建立 DOI:10．19344/j．cnki．issnl671-5276．2017．00．037 仿真模型建立王磊,申 ...
利用 Matlab Simulink 平台搭建双馈风力发电机在电网中的模型
利用 Matlab Simulink 平台搭建双馈风力发电机在电网中的模型,双馈风力发电机在风速变化的影响下转矩.电流.电压等参数波形变化. 适用于风电并网时对风电场影响的研究. 详情请见文档. ID ...
matlab 线性回归参数显著性,matlab建立多元线性回归模型并进行显著性检验及预测问题...
matlab建立多元线性回归模型并进行显著性检验及预测问题例子; x=[143 145 146 147 149 150 153 154 155 156 157 158 159 160 162 164 ...
matlab2018中变压器模块,利用MATLAB中Sim+Power+Systems模库时变压器模型的参数计算及其仿真结果比较...
[实例简介] 变压器模型 matlab 仿真参数计算第21卷第1期向秋风,等:利用 MATLAB中 Sim Power System模库时变压器模型的参数计算及其仿真结果比较 17 其标幺值:R= ...

如何利用MATLAB建立Lotka-Volterra模型及其改进模型

文章目录

前言

1 Lotka−VolterraLotka-VolterraLotka−Volterra模型

2 Lotka−VolterraLotka-VolterraLotka−Volterra改进模型

结语

如何利用MATLAB建立Lotka-Volterra模型及其改进模型相关推荐

最新文章

热门文章