Wolfram言語における数値積分

概要

Wolfram言語関数のNIntegrateは汎用数値積分器である.これは多様な一次元積分から多次元積分まで扱うことができる.

NIntegrate[f[x1,x2,,xn],{x1,a1,b1},{x2,a2,b2},,{xn,an,bn}]
範囲[a1,b1]×[a2,b2]××[an,bn]における関数 f の数値積分を求める

ある範囲における関数の数値積分を求める

一般にNIntegrateは積分領域上の被積分関数値のサンプリングを通して積分を推定する.さまざまな数値積分法によりサンプリングの初期ステップやサンプリングの展開方法が指定されている.

以下は積分 を数値的に計算する:
次のプロットは,被積分関数が評価される地点である の値の列をプロットする.サンプルは,被積分関数がより急速に変化する 付近でより密になる:

NIntegrateはユーザが指定した目標精度・確度を満足する積分推定を計算しようと試みる「積分ストラテジー」と呼ばれるアルゴリズムを使う.積分ストラテジーは重み付きの総和を頻繁に使用し,非積分関数値の集合から単独の積分推定を計算する「積分則」を使用する.

次は,積分法と積分則を指定して数値積分を計算する:

記号的前処理

NIntegrateは記号処理を使って特殊な構造の積分を簡約し,自動的に積分法を選択する.偶関数または奇関数である非積分関数や,区分関数を含む非積分関数では,積分領域が複数の別々の積分領域に変換されたり分割されたりする場合がある.

区間上で区分関数を積分する.積分領域は,特異点 で自動的に分割される:

高度に振動する被積分関数が見付けられ,特化された積分則が適用される.

区間上で高度に振動する関数を積分する.フーリエ正弦積分に特化されたメソッドが自動的に選ばれる:
以下は積分領域の上における上記の振動被積分関数のプロットである:

記号処理を行うことにより,不連続性や極度に速く変化する領域を含むさまざまな積分の自動計算が可能になる.

次は,1つの領域で高振動の区分関数を積分する.この積分領域は および で自動的に分割され,振動領域には特化されたフーリエ法が適用される:

求積法

NIntegrateには非常に古典的な一次元求積法が含まれている.

次は,さまざまな標準的な求積法で実行される一次元積分である:

古典的求積法は,抽出された被積分関数値の線形結合を形成することにより動作する.線形結合の重みのベクトルは,それぞれの求積法に対して決まっている.

次は,単位区間上の3点ニュートン・コーツ(Newton Cotes)公式に対するノードと重みである:

NIntegrateは多次元積分のために,疎な格子に基づく規則のクラスを含んでおり,また,規則を一次元規則の直積から形成することもできる.

デフォルト法の"MultidimensionalRule"は疎な格子則である:
直積則は,各次元につき1つの,一次元則のリストとして指定することができる:
下は上記の2つの多次元規則のそれぞれにおける最初のステップで使われるサンプル点を示している:

Booleはより複雑な多次元領域を指定するために使うことができる.この方法で指定された領域は,記号的前処理の間にさらに簡約化することもできる.

次は二次元関数である.底面が正方形×に内接する円錐である:
これが円錐関数の数値積分である:
これはNIntegrateにより使われるサンプル点である.サンプル点は積分領域の1つの象限にのみ存在し,すべてのサンプル点は領域 の内部にある:
次は記号的前処理なしでNIntegrateにより使われるサンプル点である.適応型アルゴリズムは積分領域全体をサンプルしなければならない.領域 の境界線付近は密にサンプリングされている.記号的前処理がない場合,NIntegrateはこの境界を実質的に数値的に見付けなければならない:

振動積分

NIntegrateには,一次元か多次元か,また有限領域上か無限領域上かを問わず,非常に広範囲の非積分関数に適用できる汎用の振動積分法が備わっている.さらにNIntegrateにはExpを含む特殊な形式の関数や,SinCos等の三角関数の他,BesselJ等の特殊関数の一次元非積分関数に特に適したメソッドもいくつか含まれている.

次は不規則に振動する関数の数値積分を計算する:
次は有限範囲上での上記の不規則な振動関数のプロットである:

特異点の自動処理

NIntegrateには,特異な被積分関数を扱う方法がいくつかある.決定的適応ストラテジー"GlobalAdaptive"および"LocalAdaptive"では,積分プロセスの収束を速めるために(変数変換に基づく)特異値処理法が使われる.ストラテジー"DoubleExponential"では,被積分関数に対して特殊な変数変換を行う台形積分法が使われる.この規則と変換の組合せにより,積分区間を含む複素平面の開集合における解析的な被積分関数に対する最適な収束が得られる.

特異値処理のある一次元積分である:
積分を100回繰り返すと,より正確な所要時間情報が得られる:
特異値処理がないと,上の積分計算はずっと遅くなる:

特殊ストラテジー

ストラテジー"DuffyCoordinates"は多次元積分におけるある種の特異点を簡約化あるいは排除する.ある球対称の積分は非常に速く収束することがある.

"DuffyCoordinates"を使った積分である:
同じ積分をNIntegrateのデフォルト設定を使って行う.正しい結果は得られるが,計算時間が格段に長くなる:

"Trapezoidal"ストラテジーでは,積分区間が厳密に1周期である場合に,解析的で周期的な被積分関数に対する最適収束が与えられる.

trapezoidalストラテジーで計算した積分である.結果を厳密値と比較する."Trapezoidal"で計算した結果はNIntegrateのデフォルト設定の結果よりも速く得られ,正確である:
以下は同じ積分をNIntegrateのデフォルトのMethod設定で計算したものであるが,上より時間がかかっている:

高次元積分の場合,あるいは大雑把な積分推定値だけが必要な場合,モンテカルロ法が便利である.NIntegrateには単純で適応型のモンテカルロ法と準モンテカルロ法とがある.

以下はモンテカルロ法を使って素早く計算された高次元積分である:

設計

NIntegrateフレームワークの主な機能は以下の通りである.

  • オブジェクト指向(メソッドプロパティ指定と通信)
  • メソッドの初期化段階とランタイム計算の分離
  • 階層的で再入可能な数値法
  • タイプおよび精度が動的なメソッド
  • プラグイン機能を使ったユーザ拡張性とプロトタイプ化
  • 特化されたデータ構造

ストラテジー,規則,プリプロセッサ

NIntegrate積分ストラテジーは,積分範囲のサンプリング方法,適用できる被積分関数のクラス,「規則ベース」のストラテジーか否かによって分類することができる.

"GlobalAdaptive"任意の被積分関数,適応型サンプリング,規則ベース
"LocalAdaptive"任意の被積分関数,適応型サンプリング,規則ベース
"DoubleExponential"任意の被積分関数,一様サンプリング
"Trapezoidal"任意の被積分関数,一様サンプリング
"MultiPeriodic"多次元被積分関数,一様サンプリング
"MonteCarlo"任意の被積分関数,一様ランダムサンプリング
"QuasiMonteCarlo"任意の被積分関数,一様擬似ランダムサンプリング
"AdaptiveMonteCarlo"任意の被積分関数,適応型ランダムサンプリング
"AdaptiveQuasiMonteCarlo"任意の被積分関数,適応型擬似ランダムサンプリング
"DoubleExponentialOscillatory"一次元無限区間振動型被積分関数
"ExtrapolatingOscillatory"一次元無限区間振動型被積分関数

NIntegrate積分ストラテジー

適応型サンプリングストラテジーは,誤差推定が大きい部分領域でより頻繁にサンプリングを行うことにより,積分推定を向上させようとするものである.通常,その部分領域を再分割して行われる.一様サンプリングストラテジーは,積分領域全体を通じて一様にサンプリングの密度を大きくすることにより,積分推定を向上させようとするものである.

規則ベースのストラテジーは指定された積分則を各部分領域に適用し,その領域に対する積分と誤差推定を得るものである.積分則はMethod->{"strategy",Method->"rule"}の設定で指定することができる.

次は,各部分領域にクレンショウ・カーチス(Clenshaw-Curtis)積分法を施す大域的適応型部分分割を使った積分の指定方法である:

NIntegrate積分則は,適用範囲が一次元領域か多次元領域かにより,また積分則のタイプにより分類される.

"BooleRule"一次元,重み付きの総和
"ClenshawCurtisRule"一次元,重み付きの総和
"GaussBerntsenEspelidRule"一次元,重み付きの総和
"GaussKronrodRule"一次元,重み付きの総和
"LobattoKronrodRule"一次元,重み付きの総和
"LobattoPeanoRule"一次元,重み付きの総和
"NewtonCotesRule"一次元,重み付きの総和
"TrapezoidalRule"一次元,重み付きの総和
"ClenshawCurtisOscillatoryRule"一次元,特化された振動則
"LevinRule"一次元および多次元,一般の振動則
"MultipanelRule"一次元,重み付きの総和,一次元則の組合せ
"CartesianRule"多次元,重み付きの総和,一次元則の積
"MultidimensionalRule"多次元,重み付きの総和

規則ベースのストラテジー"GlobalAdaptive"および"LocalAdaptive"とともに使える積分則

従来の「重み付きの総和」タイプの規則は,点の集合における関数値の決定済みの線形結合として積分を推定する.振動規則は被積分関数の特定の振動カーネルに依存する求積の重みを使って積分を推定する.

結合規則は1つ以上の部分規則から求積規則を構築する.これらは設定Method->{"rule",Method->{"subrule1",}}で指定する.

2つの一次元部分規則の直積多次元規則を指定する方法である:

すべてのストラテジーの機能は,被積分関数の記号的前処理で拡張される.前処理は最初に積分を変換,あるいは解析してから別のストラテジー(多くの場合別のプリプロセッサストラテジー)に積分を任せるプリプロセッサストラテジーにより制御される.

"SymbolicPreprocessing"全体的なプリプロセッサコントローラ
"EvenOddSubdivision"偶関数および奇関数の被積分関数を簡約化する
"InterpolationPointsSubdivision"関数の補間を含む被積分関数を部分分割する
"OscillatorySelection"振動被積分関数を検知し,適切なメソッドを選ぶ
"SymbolicPiecewiseSubdivision"区分関数を含む被積分関数を部分分割する
"UnitCubeRescaling"多次元被積分関数を単位立方体に再スケールする
"DuffyCoordinates"変換を削除する多次元特異点
"PrincipalValue"コーシー主値と等価の数値積分

NIntegrateのプリプロセッサストラテジー

プリプロセッサストラテジーは設定Method->{"preprocessor",Method->m}で指定される.ここで m は前処理が終了した後に積分が任せられるストラテジーあるいは規則である.

両方の変数について偶関数である被積分関数に,プリプロセッサストラテジーを明示的に適用する方法である.前処理の後,積分は"LocalAdaptive"ストラテジーに任される:

プロセッサストラテジーは,最終的な積分ストラテジーにより要求される仕事量を減らすことがよくある.

前の積分のサンプル点である.前処理がなければ,積分領域全体を覆い尽くすこの4倍のサンプル点が必要である:

ユーザ拡張性

組込みのメソッドは,特別な目的の積分器を効率的に構築するための構成要素として使うことができる.ユーザ定義の積分則,積分ストラテジー,プリプロセッサストラテジーを加えることもできる.