中立型 DDE

遅延のコード化

最初に、関数を記述して方程式の遅延を定義します。遅延を伴う方程式の最初の項は $\mathit{y}\left(\frac{\mathit{t}}{2}\right)$ です。

function dy = dely(t,y) 
    dy = t/2;
end

遅延を伴う方程式の別の項は ${\mathit{y}}^{\prime }\left(\mathit{t}-\pi \right)$ です。

function dyp = delyp(t,y) 
    dyp = t-pi;
end

この例では、$\mathit{y}$ での遅延と ${\mathit{y}}^{\prime }$ での遅延がそれぞれ 1 つのみ存在します。より多くの遅延が存在する場合、それらの遅延を同じ関数ファイルに追加することで、関数がスカラーの代わりにベクトルを返すようにすることができます。

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

方程式のコード化

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

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

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

解の履歴のコード化

次に、解の履歴を定義する関数を作成します。解の履歴は時間 $\mathit{t}\le {\mathit{t}}_0$ での解です。

function y = history(t)
    y = cos(t);
end

方程式の求解

最後に、積分区間 $\left[{\mathit{t}}_0{\mathit{t}}_{\mathit{f}}\right]$ を定義し、ddensd ソルバーを使用して DDE を解きます。

tspan = [0 pi];
sol = ddensd(@ddefun, @dely, @delyp, @history, [0,pi]);

解のプロット

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

0 から pi までにある等間隔の 20 点で解を計算します。

tn = linspace(0,pi,20);
yn = deval(sol,tn);

計算された解と履歴を解析解に対してプロットします。

th = linspace(-pi,0);
yh = history(th);
ta = linspace(0,pi);
ya = cos(ta);

plot(th,yh,tn,yn,'o',ta,ya)
legend('History','Numerical','Analytical','Location','NorthWest')
xlabel('Time t')
ylabel('Solution y')
title('Example of Paul with 1 Equation and 2 Delay Functions')
axis([-3.5 3.5 -1.5 1.5])

ローカル関数

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

function yp = ddefun(t,y,ydel,ypdel) % equation being solved
    yp = 1 + y - 2*ydel^2 - ypdel;
end
%-------------------------------------------
function dy = dely(t,y) % delay for y
    dy = t/2;
end
%-------------------------------------------
function dyp = delyp(t,y) % delay for y'
    dyp = t-pi;
end
%-------------------------------------------
function y = history(t) % history function for t < 0
    y = cos(t);
end
%-------------------------------------------

参照

[1] Paul, C.A.H."A Test Set of Functional Differential Equations."Numerical Analysis Reports.No. 243.Manchester, UK:Math Department, University of Manchester, 1994.