オプションに numeric と method=dverk78 を指定した dsolve コマンドは、 7-8 次 Runge-Kutta 法を用いて数値解を求めます。

以下のオプションは、 dverk78 形式をベースとする解法で使用することができます。 :

'output'

dsolve の出力形式を指定します。さらに、known オプションはユーザが定義した既知の関数を指定します。詳細は dsolve[numeric] をご参照ください。

'number', 'procedure', 'start', 'initial', と 'procvars'

これらのオプションは、関数を用いて初期条件 (IVP) を指定する場合に使用します。詳細は dsolve[numeric,IVP] をご参照ください。

'startinit','implicit', と 'optimize'

こららのオプションは、解法と計算上の振る舞いを制御します。'startinit' と 'implicit' の詳細は dsolve[numeric,IVP] を、'optimize' の詳細は dsolve[numeric] をご参照ください。

'abserr', 'relerr', 'initstep', 'minstep', と 'maxstep'

解の精度 (accuracy) とステップの幅を指定します。詳細は dsolve[Error_Control] をご参照ください。デフォルトの値は abserr=Float(1,-8), relerr=Float(1,-8), minstep=0, maxstep=2.0, そして initstep=maxstep*(relerr)^(1/6) となっています。

dverk78 による解法は、Digits の設定に基づいた任意の精度 (precision) で計算することができ、連立常微分方程式の高精度な解を得ることができます。しかしながら注意として、不正確な誤差評価によって与えられた丸め誤差を抑えるために、要求された誤差の許容範囲より大きな Digits の数を扱うことが、しばしば必要となります。Digits の設定として誤差の許容が厳しすぎる場合、アルゴリズムによって検出され、誤差が出力されます。dverk78 は Digits の設定として tol = O(Float(1,-Digits)^(9/8)) を指定しています。

注意 : 他の全ての数値を扱う dsolve の解法では、全ての Digits の設定はハードウェアの精度と等しいか、またはそれ以下で実行されています。つまり、アルゴリズムは trunc(evalhf(Digits)) で Digits < evalhf(Digits) として実行されています ( 詳細は evalhf をご参照ください。) 。

'errorest'

誤差評価を得るために、値の補間 (interpolant) や組 (pair) を持つことができ、使用する解法を制御することを許します。補間 (interpolant) オプションはデフォルトであり、補間による誤差評価を実行するように dverk78 に指示します。組 (pair) オプションは、計算による解と、より低い階数の解による誤差評価の違いを用いて、誤差評価を実行するように dverk78 に指示します ( つまり '組を用いる' という意味になります。) 。通常、補間 (interpolant) オプションは組 (pair) オプションよりも高精度な誤差評価を与えますが、より多くの計算を要することになります。このことは例えば、補間 (interpolant) オプションを伴う dverk78 が 1 ステップあたり 20 の評価関数を実行するのに対し、ペア (pair) オプションを伴う dverk78 は 13 の評価関数しか必要としないこともあり、高級なシステムに対してはより良いオプションであるかもしれません。

'maxfun'

一回の呼び出しで、常微分方程式および連立常微分方程式の評価数の最大値を指定する整数です。この計算量の制限は、全ての個々の呼び出しにおいて実行されます ( 詳細は dsolve[maxfun] をご参照ください。 ) 。このオプションは、手続き型の出力にのみ影響を与えます。デフォルトでは 0 に設定されています。

結果は plots パッケージの odeplot 関数を用いてプロットすることができます。

dsol1 := dsolve({D(y)(x) = y(x), y(0)=1}, numeric,
                method=dverk78, abserr=1e-8, relerr=1e-8);

dsol1(0);

dsol1(0.4);

dsol1(1.2);

dsol1(1.5);

dsys2 := {diff(x(t),t)=-y(t), diff(y(t),t)=x(t)};

init2 := {x(0)=1, y(0)=0};

dsol2 := dsolve(dsys2 union init2, numeric, method=dverk78,
                output=array([0,.7,1.1,1.5,2.2,2.5]));

dsol2[1,1][1]=dsol2[2,1][6,1], dsol2[1,1][2]=dsol2[2,1][6,2];

dsys3 := {diff(x(t),t$3) = x(t)+y(t), diff(y(t),t)=-x(t)+y(t)};

init3 := {x(0)=1, D(x)(0)=-1, (D@@2)(x)(0)=.5, y(0)=1.3};

dsol3 := dsolve(dsys3 union init3, numeric, method=dverk78,
                output=listprocedure, maxstep=0.01):
dsol3x := subs(dsol3,x(t)): dsol3y := subs(dsol3,y(t)):
dsol3x(0), dsol3y(0);

dsol3x(1), dsol3y(1);

dsol3x(1.5), dsol3y(1.5);

dsol4a := dsolve({diff(y(x),x,x) = -y(x), y(0)=0, D(y)(0)=1},
                numeric, method=dverk78, tolerance=1e-9,
                errorest=interpolant);

tt := time():
dsol4a(1000);

time()-tt;

dsol4b := dsolve({diff(y(x),x,x) = -y(x), y(0)=0, D(y)(0)=1},
                numeric, method=dverk78, tolerance=1e-9,
                errorest=pair);

tt := time():
dsol4b(1000);

time()-tt;

Enright, W.H. "The Relative Efficiency of Alternative Defect Control Schemes for High Order Continuous Runge-Kutta Formulas." Technical Report 252/91, Dept. of Computer Science, University of Toronto, June, 1991.

Verner, J.H. "Explicit Runge-Kutta Methods with Estimates of the Local Truncation Error." SIAM Journal of Numerical Analysis, Aug. 1978.

Much appreciation and thanks go to Dr. Jim Verner who provided us with the coefficients for the dverk78 method, and the related interpolant. Without these coefficients the method would be restricted to hardware precision.