遅延微分方程式の求解

遅延微分方程式 (DDE) は、現在の時間の解を過去の時間の解に関連付ける常微分方程式です。この遅延は定数、時間依存、状態依存、または微分依存のいずれかです。積分を開始するには、一般的に、ソルバーが最初の積分点より前の時点の解にアクセスできるように、解の履歴を提供しなければなりません。

定数遅延 DDE

定数の遅延を伴う微分方程式系は、以下の形式をもちます。

$$y^{\prime }\left(t\right)=f\left(t,y\left(t\right),y\left(t-{\tau }_1\right),\dots ,y\left(t-{\tau }_k\right)\right).$$

ここで、t は独立変数、y は従属変数の列ベクトル、y' は y の t に対する 1 次導関数を表します。遅延 (ラグ) τ₁,…,τ_k は、それぞれ正の定数です。

関数 dde23 は、t <t₀ に対する履歴 y(t) = S(t) を使用して、定数遅延を伴う DDE を解きます。

DDE の解は、一般に連続性がありますが、その導関数には不連続性が存在します。関数 dde23 は低次導関数の不連続部を追跡します。この関数は、同一の陽的ルンゲ・クッタ (2, 3) のペアと ode23 で使用する内挿を使って、微分方程式を積分します。ルンゲ・クッタ公式は、遅延よりもステップ サイズが大きい場合には陰的になります。ステップ サイズが大きくても問題がないほど y(t) が滑らかな場合、この陰的な公式は予測子修正子の反復を使って計算されます。

時間依存型および状態依存型 DDE

定数の時間遅延はより一般的な DDE 形式の特殊なケースです。

$$y^{\prime }\left(t\right)=f\left(t,y\left(t\right),y\left(d{y}_1\right),\mathrm{...},y\left(d{y}_p\right)\right).$$

時間依存型および状態依存型 DDE には遅延 dy₁、...、dy_k が含まれます。この遅延は、時間 t と状態 y の両方によって異なります。この遅延 dy_j(t, y) は、t₀ < t_f の場合の区間 [t₀, t_f] にて dy_j(t, y) ≤ t を満たさなければなりません。

関数 ddesd は、t < t₀ に対する履歴 y(t) = S(t) を使用して、時間依存型および状態依存型 DDE の解 y(t) を求めます。関数 ddesd は、従来の 4 段階法である 4 次の陽的ルンゲ・クッタ法を使用して積分し、自然内挿の残差のサイズを制御します。そのとき、遅れよりも大きいステップをとるために、反復を行います。

中立型 DDE

中立型遅延微分方程式は、y の遅延だけでなく、y′ の遅延も含みます。

$$y^{\prime }\left(t\right)=f\left(t,y\left(t\right),y\left(d{y}_1\right),\mathrm{...},y\left(d{y}_p\right),y^{\prime }\left(dy{p}_1\right),\mathrm{...},y^{\prime }\left(dy{p}_q\right)\right).$$

解の遅延は、dy_i(t,y) ≤ t を満たさなければなりません。1 次導関数の遅延は、y′ が方程式の両辺共に現れないように、dyp_j(t,y) < t を満たさなければなりません。

関数 ddensd は、時間依存および状態依存の遅延に対して与えられた形式の DDE を使用して中立型 DDE を近似することにより、その解を求めます。

$$y^{\prime }\left(t\right)=f\left(t,y\left(t\right),y\left(d{y}_1\right),\mathrm{...},y\left(d{y}_p\right)\right).$$

詳細については、Shampine の[1]を参照してください。

特定の点における解の計算

関数 deval と任意の DDE ソルバーの出力を使って、積分区間内の特定の点での解を計算できます。たとえば、y = deval(sol, 0.5*(sol.x(1) + sol.x(end))) は、積分区間の中点で解を計算します。

履歴と初期値

DDE を解く際、区間 [t₀,t_f] の解を t₀ < t_f で近似します。DDE は、y′(t) が t より前の時点でどのように解の値 (および場合によってはその導関数) に依存しているかを示します。たとえば、正の定数 τ_j の場合、定数の遅延を伴う y′(t₀) は、y(t₀ – τ₁),…,y(t₀ – τ_k) に依存します。このため、[t₀, t_k] での解は、t ≤ t₀ における値に依存します。これらの値は、t <t₀ に対する履歴関数 y(t) = S(t) を使って定義しなければなりません。

DDE の不連続性

問題が不連続性をもつ場合は、オプション構造体を使用して不連続性をソルバーに伝達するのが最良の方法です。これを行うには、関数 ddeset を使用して、問題の不連続性を含んだ options 構造体を作成します。

options 構造体には、不連続性を指定するために使用できる 3 つのプロパティ InitialY、Jumps、Events があります。不連続性の位置と特性によって選択するプロパティが決まります。

不連続性の特性	プロパティ	コメント
初期値は t = t0	InitialY	一般的に初期値 y(t0) は、履歴関数によって返された値 S(t0) であり、つまり解が初期点で連続的であることを意味しています。そうでない場合には、InitialY プロパティを使って異なる初期値を設定します。
履歴すなわち、t <t0 の場合の解、または t >t0 の場合の方程式の係数	Jumps	不連続の既知の位置 t を Jumps プロパティの値としてベクトルで与えます。dde23 に対してのみ適用します。
状態に依存	Events	dde23、ddesd および ddensd では、これらの不連続の位置を検出するために指定するイベント関数を使います。ソルバーがそのような不連続を検出すると、積分を再スタートして計算を続けます。次の積分用の履歴として、現在の積分に関する解の構造体を指定します。ソルバーは各再スタートの後で解の構造体の個々の要素を追加していき、最終の構造体が全区間の積分に対する解を与えます。次の積分の解に変化がある場合、InitialY プロパティを使って、次の積分の初期値を指定します。

不連続の伝播

一般に、解の 1 次導関数は、初期点で不連続になります。これは、履歴関数 S(t) の 1 次導関数が、一般にこの点で DDE を満たさないためです。遅延が定数の場合、y(t) の導関数の不連続性は、τ₁,…, τ_k の間隔で未来に伝播します。遅延が定数でない場合、不連続性の伝播はさらに複雑になります。定数遅延 DDEおよび時間依存型および状態依存型 DDEの形式の中立型 DDE の場合、これが伝播するたびに次に高い次数の導関数に不連続性が出現します。この意味において、解は積分が進むにつれて滑らかになります。中立型 DDEで与えられた形式の中立型 DDE の解は、定性的に異なります。解の不連続性は、より高い次数の導関数に伝播しません。特に、t₀ での y′(t) における典型的な不連続性は、[t₀, t_f] の全体を通じて y′(t) の不連続性として伝播します。

DDE の例とファイル

DDE のほとんどの問題に対して開始点として効果的に使用できるサンプル ファイルがいくつかあります。例を簡単に参照したり実行したりするには、微分方程式の例アプリを使用します。このアプリを実行するには、次のように入力します。

odeexamples

個々のサンプル ファイルを編集用に開くには、次のように入力します。

edit exampleFileName.m

例を実行するには、次のように入力します。

exampleFileName

次の表に、利用可能な DDE のサンプル ファイルと、それらで使用されているソルバーおよびオプションのリストを示します。

サンプル ファイル	使用されているソルバー	オプション指定	説明	リンクの例
ddex1	dde23	—	定数の履歴をもつ DDE	定数の遅延をもつ DDE
ddex2	dde23	'Jumps'	不連続性をもつ DDE	不連続をもつ心臓血管モデル DDE
ddex3	ddesd	—	状態依存の遅延をもつ DDE	状態依存の遅延をもつ DDE
ddex4	ddensd	—	2 つの遅延をもつ中立型 DDE	中立型 DDE
ddex5	ddensd	—	初期値をもつ中立型 DDE	中立型初期値 DDE

参照

[1] Shampine, L.F. “Dissipative Approximations to Neutral DDEs.” Applied Mathematics & Computation, Vol. 203, 2008, pp. 641–648.

参考

dde23 | ddesd | ddensd | ddeset

関連するトピック

• 定数の遅延をもつ DDE
• 状態依存の遅延をもつ DDE