中立型初期値 DDE

この例では、ddensd を使用して、時間依存の遅延を伴う初期値 DDE (遅延微分方程式) 系を解く方法を説明します。この例は、元々 Jackiewicz [1] により提示されました。

この方程式は次のとおりです。

$${\mathit{y}}^{\prime }\left(\mathit{t}\right)=2\cos \left(2\mathit{t}\right){\mathit{y}\left(\frac{\mathit{t}}{2}\right)}^{2\cos \left(\mathit{t}\right)}+\log \left({\mathit{y}}^{\prime }\left(\frac{\mathit{t}}{2}\right)\right)-\log \left(2\cos \left(\mathit{t}\right)\right)-\sin \left(\mathit{t}\right).$$

この方程式は、${\mathit{t}}_0$ における時間遅延が 0 のため、"初期値" DDE になります。したがって、解の履歴は解の計算には必要ありません。必要なのは次の初期値のみです。

$$\mathit{y}\left(0\right)=1,$$

$${\mathit{y}}^{\prime }\left(0\right)=\mathit{s}.$$

$\mathit{s}$ は $2+\log \left(\mathit{s}\right)-\log \left(2\right)=0$ の解です。この方程式を満たす $\mathit{s}$ の値は ${\mathit{s}}_1=2$ と ${\mathit{s}}_2=0.4063757399599599$ です。

方程式の時間遅延は ${\mathit{y}}^{\prime }$ の項に存在しているため、この方程式は "中立型 DDE" と呼ばれます。

MATLAB® でこの方程式を解くには、遅延微分方程式ソルバー ddensd (中立型方程式のソルバー) を呼び出す前に、方程式と遅延をコード化する必要があります。(ここで行ったように) 必要な関数をファイルの最後にローカル関数として含めることも、あるいは個別のファイルとして MATLAB パスのディレクトリに保存することもできます。

遅延のコード化

最初に、方程式の遅延を定義する無名関数を作成します。$\mathit{y}$ と ${\mathit{y}}^{\prime }$ のどちらにも形式 $\frac{\mathit{t}}{2}$ の遅延があるため、必要な関数定義は 1 つのみです。この遅延関数は後でソルバーに 2 回渡されます。1 回は $\mathit{y}$ の遅延を示すため、もう 1 回は ${\mathit{y}}^{\prime }$ の遅延を示すためです。

delay = @(t,y) t/2; 

方程式のコード化

ここで、方程式をコード化する関数を作成します。この関数にはシグネチャ yp = ddefun(t,y,ydel,ypdel) がなければなりません。ここでは以下のようになります。

• t は時間 (独立変数) です。

• y は解 (従属変数) です。

• ydel には y の遅延が含まれます。

• ypdel には ${\mathit{y}}^{\prime }=\frac{\mathrm{dy}}{\mathrm{dt}}$ の遅延が含まれます。

これらの入力はソルバーによって自動的に関数に渡されますが、変数名によって方程式のコード化方法が決まります。この場合、

• ydel$\to \mathit{y}\left(\frac{\mathit{t}}{2}\right)$

• ypdel $\to {\mathit{y}}^{\prime }\left(\frac{\mathit{t}}{2}\right)$

function yp = ddefun(t,y,ydel,ypdel) 
   yp = 2*cos(2*t)*ydel^(2*cos(t)) + log(ypdel) - log(2*cos(t)) - sin(t);
end

メモ: 関数はすべて例の最後にローカル関数として含まれます。

方程式の求解

最後に、積分区間 $\left[{\mathit{t}}_0{\mathit{t}}_{\mathit{f}}\right]$ と初期値を定義してから、ddensd ソルバーを使用して DDE を解きます。4 番目の入力引数の cell 配列内に指定することで、初期値をソルバーに渡します。

tspan = [0 0.1];
y0 = 1;
s1 = 2;
sol1 = ddensd(@ddefun, delay, delay, {y0,s1}, tspan);

方程式を再度解きます。今回は初期条件に対して $\mathit{s}$ の代替値を使用します。

s2 = 0.4063757399599599;
sol2 = ddensd(@ddefun, delay, delay, {y0,s2}, tspan);

解のプロット

解の構造体 sol1 および sol2 には x フィールドと y フィールドがあり、ソルバーが受け入れる内部タイム ステップと、その時の対応する解が含まれます。ただし、deval を使用して特定の点における解を計算することができます。

結果を比較するために 2 つの解をプロットします。

plot(sol1.x,sol1.y,sol2.x,sol2.y);
legend('y''(0) = 2','y''(0) = .40637..','Location','NorthWest');
xlabel('Time t');
ylabel('Solution y');
title('Two Solutions of Jackiewicz''s Initial-Value NDDE');

ローカル関数

ここでは、DDE ソルバー ddensd が解を計算するために呼び出すローカル補助関数を紹介しています。あるいは、これらの関数を独自のファイルとして MATLAB パスのディレクトリに保存することもできます。

function yp = ddefun(t,y,ydel,ypdel) 
   yp = 2*cos(2*t)*ydel^(2*cos(t)) + log(ypdel) - log(2*cos(t)) - sin(t);
end

参照

[1] Jackiewicz, Z."One step Methods of any Order for Neutral Functional Differential Equations."SIAM Journal on Numerical Analysis.Vol. 21, Number 3. 1984. pp. 486–511.

参考

ddensd | ddesd | dde23 | deval

トピック

• 遅延微分方程式の求解
• 中立型 DDE