24.1 常微分方程式
この関数は、lsode次のような形式の常微分方程式を解くのに使用できる。
dx -- = f (x, t) dt
Hindmarsh の ODE ソルバーLSODEを使用します。
: [x, istate, msg] = lsode (fcn, x_0, t) : [x, istate, msg] = lsode (fcn, x_0, t, t_crit) 常微分方程式 (ODE) ソルバー
解くべき微分方程式の集合は
dx -- = f (x, t) dt
と
x(t_0) = x_0
解は行列xで返され、各行はベクトルtの要素に対応します。 tの最初の要素はt_0で、システムの初期状態x_0に対応するため、出力の最初の行はx_0 になります。
最初の引数fcnは、方程式の右辺のベクトルを計算するために呼び出す 関数fを指定する文字列、インライン、または関数ハンドルです。関数は次の形式でなければなりません。
xdot = f (x, t)
ここで、xdotとx はベクトルであり、tはスカラーです。
fcnが2要素の文字列配列または文字列、インライン関数、または関数ハンドルの2要素のセル配列である場合、最初の要素は上記の関数f の名前になり、2番目の要素はfのヤコビアンを計算する関数の名前になります。ヤコビアン関数は次の形式でなければなりません。
jac = j (x, t)
ここでjacは偏微分行列である。
| df_1 df_1 df_1 |
| ---- ---- ... ---- |
| dx_1 dx_2 dx_N |
| |
| df_2 df_2 df_2 |
| ---- ---- ... ---- |
df_i | dx_1 dx_2 dx_N |
jac = ---- = | |
dx_j | . . . . |
| . . . . |
| . . . . |
| |
| df_M df_M df_M |
| ---- ---- ... ---- |
| dx_1 dx_2 dx_N |
2 番目の引数はシステムx_0の初期状態を指定します。3 番目の引数は、解を求める時間値を指定する ベクトルtです。
4 番目の引数はオプションであり、ODE ソルバーが積分しない過去の時間のセットを指定するために使用できます。これは、特異点や導関数に不連続がある点に関する問題を回避するのに役立ちます。
計算が成功すると、istateの値は 2 になります ( LSODEの Fortran バージョンと一致します)。
計算が成功しなかった場合、istate は2 以外になり、msg には追加情報が含まれます。
関数を使用してlsode_options、 のオプション パラメータを設定できますlsode。
の内部動作の詳細については、Alan C. Hindmarsh 著 「ODEPACK, A Systematized Collection of ODE Solvers」、Scientific Computing、RS Stepleman 編、(1983) またはhttps://computing.llnl.gov/projects/odepack をlsode参照してください。
例: ファンデルポール方程式を解く
fvdp = @(y,t) [y(2); (1 - y(1)^2) * y(2) - y(1)]; t = linspace (0, 20, 100); y = lsode (fvdp, [2; 0], t); See also: daspk, dassl, dasrt.
: lsode_options ()
: val = lsode_options (opt)
: lsode_options (opt, val)
lsode関数のオプションを照会または設定します。
引数なしで呼び出されると、使用可能なすべてのオプションの名前とその現在の値が表示されます。
引数を 1 つ指定すると、オプションoptの値を返します。
2 つの引数で呼び出されると、lsode_optionsオプション opt が値valに設定されます。
オプションには以下が含まれます
"absolute tolerance"
絶対許容値。ベクトルまたはスカラーのいずれかになります。ベクトルの場合は、状態ベクトルの次元と一致する必要があります。
"relative tolerance"
相対許容値パラメータ。絶対許容値とは異なり、このパラメータはスカラー値のみになります。
各積分ステップで適用されるローカルエラーテストは
abs (local error in x(i)) <= ...
rtol * abs (y(i)) + atol(i)
"integration method"
ODEシステムを解くために使用する積分法を指定する文字列。有効な値は次のとおりです。
"adams" "non-stiff"
ヤコビアンは使用されません (使用可能な場合でも)。
"bdf" "stiff"
スティフな後方微分公式 (BDF) 法を使用します。ヤコビアンを計算する関数が指定されていない場合は、lsodeヤコビアン行列の有限差分近似を計算します。
"initial step size"
最初のステップで試行されるステップ サイズ (デフォルトは自動的に決定されます)。
"maximum order"
解法の最大次数を制限します。Adams 法を使用する場合、このオプションは 1 から 12 までの範囲で指定する必要があります。それ以外の場合は、1 から 5 までの範囲で指定する必要があります。
"maximum step size"
最大ステップサイズを設定すると、非常に大きな領域を通過することを回避できます (デフォルトは指定されていません)。
"minimum step size"
許容される最小絶対ステップ サイズ (デフォルトは 0)。
"step limit"
許可される最大ステップ数 (デフォルトは 100000)。
"jacobian type"
スティフ後退微分公式（BDF）積分法で使用されるヤコビアンの種類を指定する文字列。有効な値は次のとおりです。
"full"
デフォルト。すべての偏導関数は、ユーザー指定のヤコビ関数から近似または使用されます。
"banded"
"lower jacobian subdiagonals"オプションとでそれぞれ指定された対角成分と下側および上側副対角成分の数のみが"upper jacobian subdiagonals"、ユーザー指定のヤコビ関数から近似または使用されます。ユーザー指定のヤコビ関数は、他のすべての偏導関数を任意の値に設定できます。
"diagonal"
ヤコビ関数がユーザーによって提供される場合、この設定は効果がありません。 によって近似されるヤコビアンはlsode対角に制限され、各偏導関数は状態のすべての要素に有限の変更をまとめて適用することによって計算されます。実際のヤコビアンが常に対角である場合、これは状態のそれぞれの要素にのみ有限の変更を適用するのと同じ効果がありますが、より効率的です。
"lower jacobian subdiagonals"
"jacobian type"オプションが に設定されている 場合に使用される下側サブ対角要素の数"banded"。デフォルトは 0 です。
"upper jacobian subdiagonals"
"jacobian type"オプションが に設定されている 場合に使用される上側サブ対角要素の数"banded"。デフォルトは 0 です。
以下は、3つの微分方程式をlsode関数を使って解く例です 。
## oregonator differential equation
function xdot = f (x, t)
xdot = zeros (3,1);
xdot(1) = 77.27 * (x(2) - x(1)*x(2) + x(1) ...
- 8.375e-06*x(1)^2);
xdot(2) = (x(3) - x(1)*x(2) - x(2)) / 77.27;
xdot(3) = 0.161*(x(1) - x(3));
endfunction
と初期条件x0 = [ 4; 1.1; 4 ]、方程式のセットはコマンドを使用して積分できる。
t = linspace (0, 500, 1000);
y = lsode ("f", x0, t);
これを試してみると、結果の値がt = 0と5の間、そしてt = 305付近で劇的に変化していることがわかります。より効率的な出力ポイントのセットは次のようになります。
t = 305. A more efficient set of output points might be
t = [0, logspace(-1, log10(303), 150), ...
logspace(log10(304), log10(500), 150)];
上記で使用した微分方程式のmファイルは、Octaveの配布物に含まれており、examplesディレクトリに次のオレゴネーター名前で含まれています。。