23.3 複数の変数の関数
Octave には、複数の変数の関数の積分を計算するための関数がいくつか含まれています。この手順は、通常、f をxに関して積分し、その関数をyに関して積分する関数を作成することによって実行できます。この手順は、関数を積分する次の例を使用して手動で実行できます。
f(x, y) = sin(pi*x*y) * sqrt(x*y)
xとyは0 から 1 の間です。
以下の例のように、を使用するとquadgk、二重積分を実行できます。( を除いて、任意の 1 次元求積関数をこの方法で使用できますが、quadこれは Fortran で記述されており、再帰的に呼び出すことはできないため、注意してください。)
function q = g(y)
q = ones (size (y));
for i = 1:length (y)
f = @(x) sin (pi*x.*y(i)) .* sqrt (x.*y(i));
q(i) = quadgk (f, 0, 1);
endfor
endfunction
I = quadgk ("g", 0, 1)
⇒ 0.30022
上記のアルゴリズムは、dblquad2 つの変数の積分に対して関数に実装されています。このプロセスの 3D 相当は、 triplequad3 つの変数の積分に対して関数に実装されています。たとえば、上記の結果は、dblquad次に示すように、 の呼び出しで再現できます。
I = dblquad (@(x, y) sin (pi*x.*y) .* sqrt (x.*y), 0, 1, 0, 1)
⇒ 0.30022
: q = dblquad (f, xa, xb, ya, yb)
: q = dblquad (f, xa, xb, ya, yb, tol)
: q = dblquad (f, xa, xb, ya, yb, tol, quadf)
: q = dblquad (f, xa, xb, ya, yb, tol, quadf, …)
fの二重積分を数値的に評価します。
f は、関数ハンドル、インライン関数、または評価する関数の名前を含む文字列です。関数fは、 z = f(x,y) の形式である必要があります 。ここで、xはベクトル、yはスカラーです。 xと同じ長さと方向のベクトルを返す必要があります。
xa、yaおよびxb、yb は、それぞれ x と y の積分の下限と上限です。基礎となる積分器によって、無限境界が受け入れられるかどうかが決定されます。
オプションの引数tol は、各サブ積分を積分するために使用される絶対許容値を定義します。デフォルト値は 1e-6 です。
オプションの引数quadf は、使用する基礎となる積分関数を指定します。 以外の任意の選択肢quadが使用可能で、デフォルトは ですquadcc。
追加の引数はfに直接渡されます。tolまたはquadfのデフォルト値を使用するには、または空の行列 ([]) を渡すことができます。':'
See also: integral2, integral3, triplequad, quad, quadv, quadl, quadgk, quadcc, trapz.
: q = triplequad (f, xa, xb, ya, yb, za, zb)
: q = triplequad (f, xa, xb, ya, yb, za, zb, tol)
: q = triplequad (f, xa, xb, ya, yb, za, zb, tol, quadf)
: q = triplequad (f, xa, xb, ya, yb, za, zb, tol, quadf, …)
fの三重積分を数値的に評価します。
fは、関数ハンドル、インライン関数、または評価する関数の名前を含む文字列です。関数fは、 w = f(x,y,z) の形式である必要があります 。ここで、xまたはyはベクトルであり、残りの入力はスカラーです。 xまたはyと同じ長さと方向のベクトルを返す必要があります。
xa、ya、zaおよびxb、yb、zb は、それぞれ x、y、z の積分の下限と上限です。基礎となる積分器によって、無限境界が受け入れられるかどうかが決定されます。
オプションの引数tol は、各サブ積分を積分するために使用される絶対許容値を定義します。デフォルト値は 1e-6 です。
オプションの引数quadf は、使用する基礎となる積分関数を指定します。 以外の任意の選択肢quadが使用可能で、デフォルトは ですquadcc。
追加の引数はfに直接渡されます。tolまたはquadfのデフォルト値を使用するには、または空の行列 ([]) を渡すことができます。':'
See also: integral3, integral2, dblquad, quad, quadv, quadl, quadgk, quadcc, trapz.
上で示した求積法の再帰アルゴリズムは と呼ばれます "iterated"。関数 には、別の 2-D 積分法が実装されていますquad2d。この関数は、"tiled"積分領域を長方形領域に分割し、それらの領域で個別の積分を実行することによって積分を実行します。領域は、必要な数値精度に達するために調整が必要な領域でさらに分割されます。特定の関数では、この方法は、上記の他の関数で使用される 2-D 反復よりも高速になる場合があります。
: q = quad2d (f, xa, xb, ya, yb)
: q = quad2d (f, xa, xb, ya, yb, prop, val, …)
: [q, err, iter] = quad2d (…)
タイル積分を使用して、 xa、xb、 ya、ybで定義される 2 次元領域上で適応求積法を使用してfの 2 次元積分を数値的に評価します。さらに、yaと yb はxのスカラー関数である可能性があり、非長方形領域での積分が可能になります。
fは、関数ハンドル、インライン関数、または評価する関数の名前を含む文字列です。関数f はz = f(x,y)の形式でなければならず 、すべての演算はベクトル化されて、x とy が配列入力を受け入れ、同じサイズの配列出力を返す必要があります。( xとyは同じサイズの配列か、一方がスカラーであると想定できます。) 基礎となる積分器は、積分点の配列をfに入力したり、内部ベクトル展開を使用して計算を高速化したりしますが、 f が要素ごとの演算に制限されていない場合は、予測できない結果が生じる可能性があります。これが避けられない被積分関数の場合、("Vectorized") option described below may produce より信頼性の高い結果が得られます。
追加のオプション パラメータはペアを使用して指定できます 。有効なプロパティは次のとおりです。 "property", value
AbsTol
求積法の絶対誤差許容値を定義します。デフォルト値は 1e-10 (単一の場合は 1e-5) です。
RelTol
求積法の相対誤差許容値を定義します。デフォルト値は 1e-6 (単一の場合は 1e-4) です。
MaxFunEvals
ベクトル化された関数fへの関数呼び出しの最大数。デフォルト値は 5000 です。
Singular
積分領域の端にある特異点を弱める変換を有効/無効にします。デフォルト値はtrueです。
Vectorized
ベクトル化された積分を有効または無効にします。 の値はfalse、積分関数を呼び出すときにスカラー入力のみを使用するように Octave に強制します。これにより、ベクトル化されていない、またはxまたはyのスカラー値のみを受け入れる積分関数f(x,y)が有効になります。デフォルト値は です。これは関数 で f(x,y)をラップすることによって実現されますが、これにより計算速度が大幅に低下する可能性がある ことに注意してください。truearrayfun
FailurePlot
MaxFunEvals に達する前に目的の誤差許容値に収束できない場合はquad2d、まだ調整が必要な領域のプロットが作成されます。デフォルト値はfalseです。
適応求積法は、次の条件が満たされるまで誤差の推定値を最小化するために使用されます。
error <= max (AbsTol, RelTol*|q|)
オプションの出力errは積分における誤差のおおよその境界です。ここで、I は積分の正確な値です。オプションの出力iter は、使用された関数fへのベクトル化された関数呼び出しの数です。 abs (q - I)
例1: xy平面上の長方形領域を積分する
f = @(x,y) 2*ones (size (x)); q = quad2d (f, 0, 1, 0, 1) ⇒ q = 2
結果は体積であり、この定数値積分関数の場合、ちょうど になります 。 Length * Width * Height
例2: xy平面上の三角形領域を積分する
f = @(x,y) 2*ones (size (x)); ymax = @(x) 1 - x; q = quad2d (f, 0, 1, 0, ymax) ⇒ q = 1
結果は体積であり、この定数値積分関数 f = 2の場合、三角形の面積 x 高さ、つまり になります 。 1/2 * Base * Width * Height
例3: 正方形領域上でベクトル化されていない関数を積分する
f = @(x,y) sinc (x) * sinc (y)); q = quad2d (f, -1, 1, -1, 1) ⇒ q = 12.328 (incorrect) q = quad2d (f, -1, 1, -1, 1, "Vectorized", false) ⇒ q = 1.390 (correct) f = @(x,y) sinc (x) .* sinc (y); q = quad2d (f, -1, 1, -1, 1) ⇒ q = 1.390 (correct)
最初の結果は正しくありません。これは、 fの sinc 関数間の非要素単位の演算子が、で使用される内部積分配列間で意図しない行列乗算を作成するためですquad2d。2 番目の結果では、"Vectorized"を false に設定するとquad2d、積分を計算するためにスカラー内部演算が強制的に実行されるため、計算時間が約 20 倍増加しますが、正しい数値結果が得られます。3 番目の結果では、要素単位の乗算演算子を使用して被積分関数fをベクトル化すると、 計算時間を増やすことなく正しい結果が得られます。
プログラミング ノート: 積分領域内に特異点がある場合は、積分を分割して特異点を境界上に配置するのが最適です。
既知のMATLAB の非互換性: 許容値が指定されず、積分限界が 型である場合single、Octave の積分関数は、上で指定されたデフォルトの絶対および相対誤差許容値を自動的に削減します。より厳しい許容値が必要な場合は、それを指定する必要があります。MATLAB はdouble、積分限界のクラスに関係なく、入力 に適したより厳しい許容値をそのまま残します。
参考文献: LF Shampine、 「 2D での求積法のためのMATLABプログラム」 、応用数学および計算、pp. 266–274、第 1 巻、2008 年。
See also: integral2, dblquad, integral, quad, quadgk, quadv, quadl, quadcc, trapz, integral3, triplequad.
最後に、関数integral2および は、integral3一般的な 2-D および 3-D 積分関数として提供されます。これらは、反復積分法とタイル積分法の間で自動的に選択され、dblquadおよび とは 異なりtriplequad、非長方形の積分領域で機能します。
: q = integral2 (f, xa, xb, ya, yb)
: q = integral2 (f, xa, xb, ya, yb, prop, val, …)
: [q, err] = integral2 (…)
xa、xb、 ya、yb(スカラーは有限または無限の可能性があります) で定義される 2 次元領域で適応求積法を使用してfの 2 次元積分を数値的に評価します。さらに、 yaとyb はxのスカラー関数である可能性があり、非長方形領域での積分が可能になります。
fは、関数ハンドル、インライン関数、または評価する関数の名前を含む文字列です。関数f はz = f(x,y)の形式でなければならず 、すべての演算はベクトル化されて、x とy が配列入力を受け入れ、同じサイズの配列出力を返す必要があります。( xとyは同じサイズの配列か、一方がスカラーであると想定できます。) 基礎となる積分器は、積分点の配列をfに入力したり、内部ベクトル展開を使用して計算を高速化したりしますが、 f が要素ごとの演算に制限されていない場合は、予測できない結果が生じる可能性があります。これが避けられない被積分関数の場合、("Vectorized") option described below may produce より信頼性の高い結果が得られます。
追加のオプション パラメータはペアを使用して指定できます 。有効なプロパティは次のとおりです。 "property", value
AbsTol
求積法の絶対誤差許容値を定義します。デフォルト値は 1e-10 (単一の場合は 1e-5) です。
RelTol
求積法の相対誤差許容値を定義します。デフォルト値は 1e-6 (単一の場合は 1e-4) です。
Method
使用する 2 次元積分法を指定します。有効なオプションは"auto"(デフォルト)、、"tiled"または です "iterated"。 を使用する場合、積分制限のいずれかが無限でない限り、 "auto"Octave は 方法を選択します 。"tiled"
Vectorized
ベクトル化された積分を有効または無効にします。 の値はfalse、積分関数を呼び出すときにスカラー入力のみを使用するように Octave に強制します。これにより、ベクトル化されていない、またはxまたはyのスカラー値のみを受け入れる積分関数f(x,y)が有効になります。デフォルト値は です。これは関数 で f(x,y)をラップすることによって実現されますが、これにより計算速度が大幅に低下する可能性がある ことに注意してください。truearrayfun
適応求積法は、次の条件が満たされるまで誤差の推定値を最小化するために使用されます。
error <= max (AbsTol, RelTol*|q|)
err は積分における誤差のおおよその境界であり 、ここでI は積分の正確な値です。 abs (q - I)
例1: xy平面上の長方形領域を積分する
f = @(x,y) 2*ones (size (x)); q = integral2 (f, 0, 1, 0, 1) ⇒ q = 2
結果は体積であり、この定数値積分関数の場合、ちょうど になります 。 Length * Width * Height
例2: xy平面上の三角形領域を積分する
f = @(x,y) 2*ones (size (x)); ymax = @(x) 1 - x; q = integral2 (f, 0, 1, 0, ymax) ⇒ q = 1
結果は体積であり、この定数値積分関数 f = 2の場合、三角形の面積 x 高さ、つまり になります 。 1/2 * Base * Width * Height
例3: 正方形領域上でベクトル化されていない関数を積分する
f = @(x,y) sinc (x) * sinc (y)); q = integral2 (f, -1, 1, -1, 1) ⇒ q = 12.328 (incorrect) q = integral2 (f, -1, 1, -1, 1, "Vectorized", false) ⇒ q = 1.390 (correct) f = @(x,y) sinc (x) .* sinc (y); q = integral2 (f, -1, 1, -1, 1) ⇒ q = 1.390 (correct)
最初の結果は正しくありません。これは、 fの sinc 関数間の非要素単位の演算子が、で使用される内部積分配列間で意図しない行列乗算を作成するためですintegral2。2 番目の結果では、"Vectorized"を false に設定するとintegral2、積分を計算するためにスカラー内部演算が強制的に実行されるため、計算時間が約 20 倍増加しますが、正しい数値結果が得られます。3 番目の結果では、要素単位の乗算演算子を使用して被積分関数fをベクトル化すると、 計算時間を増やすことなく正しい結果が得られます。
プログラミング ノート: 積分領域内に特異点がある場合は、積分を分割して特異点を境界上に配置するのが最適です。
既知のMATLAB の非互換性: 許容値が指定されず、積分限界が 型である場合single、Octave の積分関数は、上で指定されたデフォルトの絶対および相対誤差許容値を自動的に削減します。より厳しい許容値が必要な場合は、それを指定する必要があります。MATLAB はdouble、積分限界のクラスに関係なく、入力 に適したより厳しい許容値をそのまま残します。
参考文献: LF Shampine、 「 2D での求積法のためのMATLABプログラム」 、応用数学および計算、pp. 266–274、第 1 巻、2008 年。
See also: quad2d, dblquad, integral, quad, quadgk, quadv, quadl, quadcc, trapz, integral3, triplequad.
: q = integral3 (f, xa, xb, ya, yb, za, zb)
: q = integral3 (f, xa, xb, ya, yb, za, zb, prop, val, …)
xa、xb、ya、yb、za、zb (スカラーは有限または無限の可能性があります)で定義される 3 次元領域で適応求積法を使用して fの 3 次元積分を数値的に評価します。さらに、yaとyb はxとzaのスカラー関数である可能性があり、zb はxとyのスカラー関数である可能性があり、非長方形領域での積分が可能になります。
fは、関数ハンドル、インライン関数、または評価する関数の名前を含む文字列です。関数f はz = f(x,y,z)の形式でなければならず 、すべての演算はベクトル化されて、 x、y、z が配列入力を受け入れ、同じサイズの配列出力を返す必要があります。( x、y、z は同じサイズの配列またはスカラーのいずれかであると想定できます 。) 基礎となる積分器は、積分点の配列をfに入力したり、内部ベクトル拡張を使用して計算を高速化したりしますが、 f が要素ごとの演算に制限されていない場合は、予測できない結果が生じる可能性があります。これが避けられない被積分関数の場合、("Vectorized") option 以下に説明する方法の方が信頼性の高い結果をもたらす可能性があります。
追加のオプション パラメータはペアを使用して指定できます 。有効なプロパティは次のとおりです。 "property", value
AbsTol
求積法の絶対誤差許容値を定義します。デフォルト値は 1e-10 (単一の場合は 1e-5) です。
RelTol
求積法の相対誤差許容値を定義します。デフォルト値は 1e-6 (単一の場合は 1e-4) です。
Method
使用する 2 次元積分法を指定します。有効なオプションは"auto"(デフォルト)、、"tiled"または です "iterated"。 を使用する場合、積分制限のいずれかが無限でない限り、 "auto"Octave は 方法を選択します 。"tiled"
Vectorized
ベクトル化された積分を有効または無効にします。 の値はfalse、積分関数を呼び出すときにスカラー入力のみを使用するように Octave に強制します。これにより、ベクトル化されていない、またはx、y、またはzのスカラー値のみを受け入れる積分関数f(x,y,z) が有効になります。デフォルト値は です 。これは関数 でf(x,y,z)をラップすることによって実現されますが、これにより計算速度が大幅に低下する可能性があることに注意してください。 truearrayfun
適応求積法は、次の条件が満たされるまで誤差の推定値を最小化するために使用されます。
error <= max (AbsTol, RelTol*|q|)
err は積分における誤差のおおよその境界であり 、ここでI は積分の正確な値です。 abs (q - I)
例1: 長方形の体積を積分する
f = @(x,y,z) ones (size (x)); q = integral3 (f, 0, 1, 0, 1, 0, 1) ⇒ q = 1.00000
この定数値積分関数の場合、結果はちょうど の体積になります 。 Length * Width * Height
例2: 球状の体積を積分する
f = @(x,y) ones (size (x)); ymax = @(x) sqrt (1 - x.^2); zmax = @(x,y) sqrt (1 - x.^2 - y.^2); q = integral3 (f, 0, 1, 0, ymax, 0, zmax) ⇒ q = 0.52360
この定数値積分関数の場合、結果は単位球の 1/8 の体積、つまり になります1/8 * 4/3 * pi。
例3: 立方体上の非ベクトル化関数を積分する
f = @(x,y) sinc (x) * sinc (y), * sinc (z); q = integral3 (f, -1, 1, -1, 1, -1, 1) ⇒ q = 14.535 (incorrect) q = integral3 (f, -1, 1, -1, 1, -1, 1, "Vectorized", false) ⇒ q = 1.6388 (correct) f = @(x,y,z) sinc (x) .* sinc (y), .* sinc (z); q = integral3 (f, -1, 1, -1, 1, -1, 1) ⇒ q = 1.6388 (correct)
最初の結果は正しくありません。これは、 fの sinc 関数間の非要素単位の演算子が、で使用される内部積分配列間で意図しない行列乗算を作成するためですintegral3。2 番目の結果では、"Vectorized"を false に設定するとintegral3、積分を計算するためにスカラー内部演算が強制的に実行されるため、計算時間が約 30 倍増加しますが、正しい数値結果が得られます。3 番目の結果では、要素単位の乗算演算子を使用して被積分関数fをベクトル化すると、 計算時間を増やすことなく正しい結果が得られます。
プログラミング ノート: 積分領域内に特異点がある場合は、積分を分割して特異点を境界上に配置するのが最適です。
既知のMATLAB の非互換性: 許容値が指定されず、積分限界が 型である場合single、Octave の積分関数は、上で指定されたデフォルトの絶対および相対誤差許容値を自動的に削減します。より厳しい許容値が必要な場合は、それを指定する必要があります。MATLAB はdouble、積分限界のクラスに関係なく、入力 に適したより厳しい許容値をそのまま残します。
参考文献: LF Shampine、 「 2D での求積法のためのMATLABプログラム」 、応用数学および計算、pp. 266–274、第 1 巻、2008 年。
See also: triplequad, integral, quad, quadgk, quadv, quadl, quadcc, trapz, integral2, quad2d, dblquad.
上記の積分はかなり遅くなることがあり、その問題は積分の次元とともに指数関数的に増大します。2 次元積分のもう 1 つの可能な解決方法は、前のセクションで説明したように直交コロケーションを使用することです (直交コロケーションを参照)。関数f(x,y)のxとy の0 から 1 までの積分は、 n点を使用してi=1:nとj=1:nの和で近似できます。q(i)*q(j)*f(r(i),r(j))ここで、qとrは によって返される値です。3 つ以上の変数への一般化は簡単です。次のコードは、n=8点 colloc (n)を使用して、検討した積分を計算します。
f = @(x,y) sin (pi*x*y') .* sqrt (x*y');
n = 8;
[t, ~, ~, q] = colloc (n);
I = q'*f(t,t)*q;
⇒ 0.30022
点の数によって近似の質が決まることに注意してください。 0 と 1 の間ではなく、 aとbの間で積分を実行する必要がある場合は 、変数の変更が必要です。