ドキュメンテーション

最新のリリースでは、このページがまだ翻訳されていません。 このページの最新版は英語でご覧になれます。

捕食者と被食者の方程式の求解

この例では、ode23ode45 の両方を使用して捕食者/被食者モデルを表す微分方程式を解く方法を示します。これらはルンゲ・クッタ積分手法の可変ステップ サイズを使用する常微分方程式の数値解に対する関数です。ode23 は、中規模な精度に対して単純に 2 次、3 次の公式の組み合わせを使用し、ode45 は、高い精度に対して 4 次、5 次の組み合わせを使用します。

Lotka-Volterra の方程式または捕食者と被食者のモデルとして知られる 1 次の常微分方程式の組み合わせについて考えます。

変数 および は、被捕食者と捕食者のそれぞれの個体群のサイズを表します。2 次のクロス項は、種の間の相互関係を説明します。捕食者が存在しない場合に被捕食者の個体群は増加し、被捕食者が少なくなると捕食者の個体群は減少します。

方程式のコード作成

方程式系のシミュレーションを実行するために、状態と時間の値を与えて、状態導関数の列ベクトルを返す関数を作成します。2 つの変数 および は、MATLAB でベクトル y の最初の 2 つの値として表すことができます。同様に、導関数はベクトル yp の最初の 2 つの値です。この関数は t および y の値を受け入れて、yp の方程式によって生成された値を返さなければなりません。

yp(1) = (1 - alpha*y(2))*y(1)

yp(2) = (-1 + beta*y(1))*y(2)

この例では、方程式は lotka.m というファイルに含まれています。このファイルは、 および のパラメーター値を使用します。

type lotka
function yp = lotka(t,y)
%LOTKA  Lotka-Volterra predator-prey model.

%   Copyright 1984-2014 The MathWorks, Inc.

yp = diag([1 - .01*y(2), -1 + .02*y(1)])*y;

方程式系のシミュレーション

ode23 を使用して、区間 lotka で定義された微分方程式を解きます。捕食者と被捕食者の個体群が等しくなるように の初期条件を使用します。

t0 = 0;
tfinal = 15;
y0 = [20; 20];   
[t,y] = ode23(@lotka,[t0 tfinal],y0);

結果のプロット

時間に対して結果の個体群をプロットします。

plot(t,y)
title('Predator/Prey Populations Over Time')
xlabel('t')
ylabel('Population')
legend('Prey','Predators','Location','North')

次に、相互に個体群をプロットします。結果の位相面プロットにより、個体群間の巡回的な関係が非常に明らかになります。

plot(y(:,1),y(:,2))
title('Phase Plane Plot')
xlabel('Prey Population')
ylabel('Predator Population')

さまざまなソルバーの結果の比較

2 度目は ode23 ではなく ode45 を使用して方程式系を解きます。ode45 ソルバーでは、各手順において時間がかかりますが、手順もより大きくなります。それにもかかわらず、既定の設定で、ソルバーは連続的な拡張式を利用して、使用した各ステップの期間において 4 つの等間隔な時間点で出力を生成するため、ode45 の出力は滑らかです ('Refine' オプションを使用して点の数を調整できます)。比較のために両方の解をプロットします。

[T,Y] = ode45(@lotka,[t0 tfinal],y0);

plot(y(:,1),y(:,2),'-',Y(:,1),Y(:,2),'-');
title('Phase Plane Plot')
legend('ode23','ode45')

結果は、異なる数値的手法を使用して微分方程式を解くと、若干異なる解が生成されることを示します。

参考

|

関連するトピック