傾き場と常微分方程式ソルバー

1階微分方程式を、小さな接線の場として読み取ろう。プロットをクリックして初期条件を追加すると、対応する解曲線が描かれます。

インタラクティブシミュレーションを読み込んでいます...

レッスン

理論 — 傾き場と常微分方程式ソルバー

1階微分方程式 dy/dx = f(x,y) は、y が何であるかを直接示すものではありません。各点において y がどのように変化しているか(変化の速さ)を示します。とは、すべての点における傾きがそこでの f の値と等しくなるような曲線のことであり — それぞれの始点を通る曲線からなる解の族が存在します。

各記号の意味

dy/dx
解曲線の傾き — 方程式が求める対象ではなく、あらかじめ規定している量。
f(x,y)
右辺: 任意の点における傾きを返す規則。ここでの5つのプリセットは yx*y-y/xy*(1-y)x^2-y です。
(x₀,y₀)
初期条件 — 解曲線が必ず通過しなければならない1点。解の族の中から1本の曲線を選び出すものであり、プロット上をクリックすることで指定されます。

表示の読み方

このパネルは描画されたものを説明するのではなく、カウントして表示します: デフォルトのグリッドにおける 21 x 21 = 441 個の傾きマーク、現在表示されている最大の有限傾き、プロットウィンドウの範囲、そして作成した解曲線の数です。クリックするたびに初期条件が追加され、RK4ソルバーが対応する曲線を x の正方向と負方向の両方に向かって追跡します。

前提
方程式が1階であり、すでに dy/dx について整理されていること、および f が各点に対してただ1つの傾きを返すことを前提としています。グリッドは傾き場そのものではなく、その標本(サンプル)です — 密度を上げると、1本の解曲線も変化させることなく、より多くのマークが描画されます。
成り立たない場合
-y/x を選択すると、カウンターは 421 finite, 20 vertical/singular に変化します。これらは y 軸上に位置するグリッド点であり、x = 0 であるため傾きが定義されません。その場所には場としての方向がまったく存在せず、いかなる解曲線もそこを通過できません。また、傾き場は、有限の x において解が無限大に発散する場所を示すこともできません — 解が急激に変化していても、画面の図はきわめて穏やかなままです。

この絵は、方程式を解かずに答えを示す 🖖

方向フィールド(傾き場)は、格子点上に傾き f(x,y) を持つ短い線分を描くことで、一階常微分方程式 dy/dx = f(x,y) を表します。これにより、方程式を解析的に解くことなく、潜在的な解曲線の流れを幾何学的に可視化できます。dy/dx = 0となる平衡解、安定・不安定なノード、および変曲点を示します。

解かずに矢印をたどる 🖖

傾き場は無数の小さな道しるべだと考えてください。各点で、解が次にどちらへ進むべきかを示しています。解曲線はこの向きにただ従うだけで、流れに漂う一枚の葉のようです。各点の向きはただ一つなので、出発点(初期条件)を選ぶと一本の曲線が定まります。どこかをクリックすれば、ツールがその唯一の経路を描きます。

RK4は密かにシンプソンの公式になる 🖖

このツールが使うRK4法は、各ステップで傾きを4回評価し、重み 1:2:2:1 で混ぜ合わせます。方程式が右辺に y を含まない dy/dx = f(x) の形なら、この混合は積分のシンプソンの公式とちょうど一致します:h/6·(f₀ + 4f_m + f₁)。つまり常微分方程式を数値的に解くことと定積分を計算することが、この場合は同じ計算になるのです。

練習

自分で確かめる

まず答えを予想し、それから上のコントロールで確かめてください。予想を決めてから答えを開くこと。それが練習になる条件です。

  1. 双曲線のプリセットは y′ = −y/x を解き、「有限な|傾き|の最大値」を 10 と表示します。格子密度を 30 に上げてこの行をもう一度読み、指数関数放物線のプリセットでも同じことをしてください。

    答えを表示
    双曲線では、密度 10・20・30 に対してこの行は 5 → 10 → 15 と動きます。指数関数では 3 のまま、放物線では 11 のままで、何をしても変わりません。特異な場では、この行が語っているのは方程式ではなくあなたの格子です。−y/x には最大の傾きなど存在せず、x → 0 で無限大へ走り去るので、ツールはたまたま置いた最も急な目盛りしか報告できません。軸に最も近い列は 8/密度 の位置にあり、密度 20 なら x = 0.4、そして 4/0.4 = 10 — 画面上の数値そのものです。垂直と印された 20 本は x = 0 の列まるごとで、そこでは傾きは数ですらありません。他の二つの場はその窓の中で本当に有界なので、格子は毎回同じ最大値を見つけます。
  2. ロジスティックのプリセットは y′ = y(1−y) を解きます。その成長は y = 1/2 で最も速く、そこでの変化率は 0.25 であることがよく知られています。ところがパネルは最大 |傾き| を 0.75 と表示します。三倍も急な傾きを、どこで拾っているのでしょうか。

    答えを表示
    個体群の外側です。このプリセットの窓は y = −0.5 から 1.5 までで、y = −0.5 では変化率が (−0.5)(1.5) = −0.75 になります。二本の平衡線の間にあるどれよりも急で、しかも減少方向です。行が報告するのは絶対値だからです。y min を 0、y max を 1 にすればこの行は 0.25 に落ち、生物学が置くとおり y = 0.5 で見つかります。そこまで厳密である必要すらありません。−0.2 から 1.2 でも 0.25 になります。(−0.2)(1.2) は 0.24 にしかならないからです。有名な数値が語るのはロジスティックな個体群であり、パネルが語るのはあなたが描かせた長方形です。

全プロセスの詳細解説

  1. y′ = −y/x の曲線群 6 ステップ

    y′ = −y/x(−4 ≤ x ≤ 4、−4 ≤ y ≤ 4、格子密度 20)。解曲線を 1 本も描く前に、そもそも傾きを定義できるマークがいくつあるかを数え、最も急な傾きを持つものを求め、さらにこの場がどのような曲線族から構成されているかを導き出しなさい。

    1. 格子は幅 20 区間であり、両端も数えるためサンプル線は 20 本ではなく 21 本存在する。したがって、幅 8 単位の表示範囲では 0.4 ごとにサンプリングが行われ、垂直方向も同じ間隔になる。

    2. 21 は奇数であるため、1 本のサンプル線が表示範囲のちょうど中央に位置し、この表示範囲の中央は x = 0 である。その列こそが、右辺の分母に 0 が現れる場所である。

    3. その列において y ≠ 0 であるすべての点は、値を持たない商を与える。すなわち、垂直に描かれる方向は存在するが、報告すべき数値は存在しない。これは 21 点のうち 20 点にあたる。最後の点は原点であり、そこでの比は 0/0 となるため、ツールはマークを水平に描くことで処理している。これは必然的な結果ではなく選択によるものであり、有限値を持つマークの個数が 420 ではなく 421 になるのはそのためである。

    4. 傾きの急さは |y|/|x| で表されるため、最も急なマークは、0 でない最も狭い列上の最高点、すなわち 4 を 0.4 で割った値となる。条件を満たすマークは特異線のすぐ隣に位置しており、特筆すべき点はそれだけである。

    5. 次に曲線そのものについて考える。方程式は変数分離形であり、両辺を積分すると対数の和が積へと変換される。すべての解曲線に沿って xy は一定となるため、この曲線族は直角双曲線群である。

    6. 点 (1, 2) で試してみよう。この点におけるマークの傾きは −2 である。この点を通る双曲線 xy = 2 のここでの傾きは −2 であり、適合する。同じ点を通る円 x² + y² = 5 の傾きは −x/y = −0.5 であり、適合しない。半径に垂直という条件に対応する方程式は x と y が逆になった y′ = −x/y である。

    解答

    有限なマーク 421 個、垂直なマーク 20 個、および最も急な傾き 10。これら 3 つの数値のうち、微分方程式そのものについて述べているのは 2 つだけである。最後の項目について格子の密度を記号 d のまま再計算すると、格子上で捉えられる最も急な傾きは 4 を 8/d で割ったもの、すなわち d/2 となる。この数値は格子密度の設定値の半分に過ぎない。方程式自体には最大の傾きなど存在せず、x → 0 において |y′| は無限に増大する。パネルが報告しているのは方向場がどこまで急であるかではなく、どれほど細かく観察したかである。密度 21 ではこの効果がより顕著になる。各辺 22 本のサンプル線が x = 0 の上に落ちるのではなく挟み込む形になるため、484 個すべてのマークが有限値となり、特異列は以前と全く同じ特異性を保ったままカウントから消滅する。

  2. y′ = y(1 − y) における y = 0.1 を通る解 6 ステップ

    y′ = y(1 − y)(−2 ≤ x ≤ 6、−0.5 ≤ y ≤ 1.5、格子密度 20)。方程式を解くことなく、場がどこで水平になるか、y = 0.1 を通る解が進む方向、およびこの表示範囲内でマークが最も急になる場所を答えよ。さらに、その上昇にどれだけの時間を要するかを数値で示せ。

    1. マークが水平になるのは右辺が 0 になるときだけであり、積が 0 になるのはいずれかの因数が 0 になるときに限られる。該当する値は 2 つあり、それぞれがそれ自体で解をなしている。定数関数 y ≡ 0 および y ≡ 1 は、すべての x に対してこの方程式を満たす。

    2. その両方がサンプル線上に位置する。21 個の y の値は −0.5 から始まり 0.1 刻みで変化するため、0 と 1 もその中に含まれており、2 行の水平な行は両側の傾きから示唆されるだけでなく実際に描画される。

    3. 水平な行の間では両方の因数が正となるため、すべてのマークは右上がりに傾く。その外側では一方の因数の符号が反転するため、マークは右下がりに傾く。0.1 から始まる解は上昇するが、1 に到達することは決してない。その線を越えるにはその上に存在する定数解と交差する必要があるが、この方程式の 2 つの解が交差することはあり得ないからである。

    4. 「上昇して有界であること」は「安定であること」と同義ではない。そのため、各平衡点をわずかにずらして符号の変化を観察する。0 のすぐ上では傾きはおよそ +ε であり、さらに遠ざかる方向へ押し出す。1 のすぐ上では傾きはおよそ −ε であり、押し戻す。1 つの方程式、2 本の水平線、相反する結論が得られる。

    5. 最も急なマークは、成長が起きている場所にあるわけではない。0 ≤ y ≤ 1 の内部では最も速い変化は y = 0.5 のときでその値は 0.25 であるが、表示範囲は −0.5 や 1.5 にまで達しており、そこでの積は −0.75 となる。これは 3 倍も急であり、このモデルが本来対象としていない個体数の値におけるものである。

    6. ある 1 つの量が曲線を直線へと変える。オッズの対数 y/(1 − y) を微分すると、合成関数の微分公式によって方程式が与える因数 y(1 − y) がちょうど現れ、その 2 つが相殺して、残った部分は一定の割合 1 で増加する。

    解答

    2 つの平衡点(0 は反発的、1 は吸引的)、および 0 ≤ y ≤ 1 の完全な外側に位置する最大 |傾き| 0.75。最後に示したステップこそ、覚えておくべきものである。対数オッズは x の 1 単位あたり正確に 1 ずつ増加するため、初期値が何であれ、10% から 90% への上昇には x の変化量として ln(9 ÷ 1/9) = ln 81 = 4.394 を要する。パネルにはこの数値は表示されないが、これこそがモデラーが実際に求める数値である。すなわち、曲線の高さではなく、成長率のみによって決まる遷移の長さである。この表示範囲の幅は 8 であるため、左端の 0 付近から入る解は、その表示範囲の半分強を費やして上昇の全体を行い、残りは y = 1 に沿って水平に推移する。

参考文献 (1)
  • Direction fields, initial-value problems and when a solution is guaranteed to exist: W. E. Boyce & R. C. DiPrima, Elementary Differential Equations and Boundary Value Problems, 10th ed. John Wiley & Sons, 2012. ISBN 978-0-470-45831-0.

例題

  • y' = y - 441個の線分があり、y = 0 に沿っては水平で、傾きが 3 を示す四隅で最も急になります。この水平な列自体が一つの解であり、どの曲線もそこを越えることはありません。そのため、最初にクリックした位置の符号によって曲線が上昇するか下降するかが完全に決まり、後から変わることはないのです。
  • ロジスティック - 二つの直線 y = 0 と y = 1 に沿って水平な線分が並んでおり、どちらも解となります。この二直線の間から始まる曲線はそこに閉じ込められ、必ず 1 へと向かいます。グリッド上で最も急な傾き 0.75 は、このS字曲線上には存在しません。下の解よりもさらに低い、y = −0.5 の位置にあるのです。
  • 双曲線 - カウンターの表示は、有限な傾きを示す421個の線分と、20個の垂直な線分に分かれます。この垂直な列は x = 0 に位置しており、ここでは −y/x は値を持ちません。それでは、密度を20から21へと少し動かしてみましょう。484個すべてが有限となり、特異点は一つも現れません。グリッドが座標軸と重ならなくなったのです。特異点が移動したわけではなく、サンプリングがそこをすり抜けただけなのです。
  • 放物線 - 水平な線分が放物線 y = x² を描き出し、どの解もそこを水平に横切ります。この放物線自体は解ではありません。方程式がその上で傾き 0 を要求するのに対し、曲線の傾きは 2x となるためです。最も急な線分は左下の隅にあり、その値は 11 です。