4次ルンゲクッタ法で2階常微分方程式を解く計算機

使い方

1. 右辺Fをx、y、pの安全な式で入力し、モデルの符号と実数定義域を確認します。
2. 初期点x0、初期状態、x0と異なる終点xnを入力します。2階ではp0=y'(x0)も指定します。
3. 分割数nを1〜200で指定し、固定刻みh=(xn−x0)/nとして計算します。
4. SVGの形、代替数値表、TSVを確認し、nを増やした結果や解析解・別solverと比較します。

計算式・考え方

(y,p)系の始点・中点2回・終点で組のstageを作り、yとpをそれぞれ(k1+2k2+2k3+k4)/6で同時更新します。

注意事項

• このページは古典的Runge–Kutta 4次法で2階常微分方程式を扱います。方式名と階数を固定し、途中で別の更新式へ自動変更しません。
• 古典的RK4は始点、中点2回、終点の4stageをk1、k2、k3、k4として評価し、(k1+2k2+2k3+k4)/6を加えます。滑らかな問題で大域誤差は一般にhの4次で減り、step半減で誤差がおよそ1/16になるのが目安です。4次という名称は無条件の正確さや安定性を意味せず、非滑らかなF、硬い問題、丸め誤差支配では理論次数から外れます。
• 2階ページではp=y'と置き、y'=p、p'=F(x,y,p)の連立1階系として同時に進めます。既定例F=−4p−4y、x0=0、y0=0、p0=1の解析解はy=x exp(−2x)、p=(1−2x)exp(−2x)です。yだけでなくpも解析式と比較すると、位置に相当する状態と導関数状態の両方を検算できます。
• 連立系の同時更新は重要です。Euler法でも先にpを新しくしてからその新しいpで同じstepのyを更新してはいけません。中点RK2と古典RK4も、各stageのy増分とp増分を同じstage状態から組にして計算します。逐次的に片方だけ早く更新すると別のアルゴリズムになり、期待する収束次数や解析解との一致を失います。
• この2階古典RK4ページの既定n=50ではx=1でyは約0.135335274181835、pは約−0.135335259156936です。4組のstageでy用増分とp用増分を同じ時刻・同じstage状態から作ります。片方だけ新しいstageへ進めると古典RK4の連立系公式ではなくなるため、1分割テストで各組を再構成すると実装検算に役立ちます。
• この計算機は初期値問題を有限個の格子点で近似する学習・検算用ツールです。入力するのは右辺F、初期点x0、初期状態、終点xn、分割数nです。刻み幅はh=(xn−x0)/nで全区間を通して一定とし、初期点を含むn+1点を表示します。グラフだけで形を判断せず、展開できる数値表とTSVで各格子点の値も確認してください。
• 実装は固定刻み・非適応です。途中で解が急変しても刻み幅を自動で細かくせず、局所誤差推定による再試行も行いません。分割数を増やした結果が安定するか、nと2nで終点値や系列を比較することが最低限の確認になります。小さいhほど常に正確とは限らず、打切り誤差が減る一方でbinary64の丸め誤差や式評価の相殺が目立つ場合があります。
• この実装はstiffness非対応です。減衰時間尺度が大きく異なる硬い問題では、陽的Euler法や陽的Runge–Kutta法は真の解が滑らかでも不安定になり、刻みを非常に小さくしないと発散や振動が生じます。硬さの疑いがある問題、保存則が重要な長時間積分、業務上の重要判断では、stiff solverや誤差制御を備えた検証済み数値計算ライブラリを使用してください。
• 表示値には誤差保証なしという制約があります。計算機は局所誤差、累積した大域誤差、条件数、解の存在・一意性を証明しません。複数の解法、異なる分割数、既知の解析解、保存量や単調性など問題固有の性質と照合してください。入力したFが数式モデルと一致しているか、符号、係数、初期条件、変数の単位を別に点検することも必要です。
• 右辺式は変数xと状態変数、定数pi・π・e、四則演算、べき乗、丸括弧、abs、sqrt、exp、ln、log、log10、三角・逆三角・双曲線関数、floor、ceilだけを安全なパーサーで処理します。evalやnew Function、代入、配列、プロパティ参照、セミコロン、未知名は使いません。式は200文字、128トークン、入れ子32段に制限します。
• x0とxnの絶対値および区間幅は1,000,000以下、初期状態の絶対値は1,000,000以下、分割数は1〜200です。各stageのF、増分、途中状態、更新後状態は有限値かつ絶対値1e100以下でなければなりません。平方根・対数の定義域外、0除算、指数オーバーフロー、状態の暴走が起きた場合は、無理に丸めたり上限へ貼り付けたりせず計算を停止します。
• 終点xnをx0より小さくすれば逆向き区間も計算できます。その場合hは負になり、同じ更新式を逆方向へ適用します。逆積分が元の初期状態へ正確に戻る保証はなく、散逸系や不安定な逆時間問題では誤差が急増します。逆向き機能は式と実装の検算に便利ですが、物理的に逆向きの解が意味を持つかはモデルごとに判断してください。
• グラフは有限な格子値を直線で結んだ表示で、格子間の厳密な解や高次補間ではありません。急峻なピーク、振動、符号変化、eventは格子点の間に隠れることがあります。この計算機にはdense output、event location、零点検出、最大値探索がありません。必要なら区間を狭めて分割数を増やし、専用solverの連続出力と比較してください。
• 1分割は更新式を手計算で再構成する境界テストに向きますが、一般に精度を期待する設定ではありません。200分割はこのページの安全上限であり、十分な精度を意味しません。上限はSVG、代替表、TSVをスマートフォンでも扱える201点に抑え、ユーザー端末のCPU・メモリ使用量と操作性を予測可能にするためのものです。より細かい格子が必要なら専門ライブラリへ切り替えてください。
• すべての構文解析、最大4stage×200stepの反復、SVG描画、数値表生成、TSVコピーはclient-onlyで利用者のブラウザ内だけで実行します。右辺式、初期値、計算系列をサーバーや外部APIへ送信・保存せず、実行時fetchや外部CDNも使いません。ページを閉じれば入力と結果はサイト側に残らず、サーバー負荷を増やさずに端末の有限リソースを利用します。
• 更新式と固定刻み法の位置付けは2026年7月13日にNIST Digital Library of Mathematical Functionsの数値微分方程式解法とSciPy solve_ivpの公開説明を照合しました。実装は既知の解析解、定数右辺、平衡解、逆向き区間、分割数1・200、step半減時の収束次数、式注入拒否、PythonとJavaScriptの一致で確認します。公開説明との照合は個別入力の精度保証ではありません。
• 固定刻みbinary64による教育用近似で、適応刻み・stiffness検出・大域誤差保証はありません。

よくある質問