クレンショウ-カーティス(Clenshaw–Curtis)型の数値積分公式です。

[2]の論文”IS GAUSS QUADRATURE BETTER THAN CLENSHAW–CURTIS?”の結論だけかいつまみますと、

・クレンショウ-カーティスの求積法は、基本的にはガウス-ルジャンドル求積法と同じ性能
・クレンショウ-カーティスの求積法は実装が簡単

ということです。
本気で計算精度を求めようとするならばガウス-ルジャンドル求積法、簡単に、でもかなりの精度で求めたければクレンショウ-カーティスの求積法なのだと思います。
ガウス-ルジャンドル求積法は分点の位置と重みの計算が簡単ではないのです。

積分
\(
\displaystyle \int_{-1}^1 f(x) dx \sim \sum_{k=0}^{n} \omega_k f(x_k)
\)
Nが偶数の時、分点\(x_k\)と重み\(\omega_k\)は
\(
\displaystyle x_k=\cos\left(\frac{k\pi}{N}\right),\;\;\;k=0,1,\cdots,N
\)
\(
\begin{align}
\omega_0 &=\omega_N =\frac{1}{N^2 -1} \\
\omega_s &=\omega_{N-s} = \frac{4}{n^2}\sum^{n/2}_{j=0}{\prime\prime }\frac{1}{1-4j^2}\cos\left(\frac{2j\pi s}{n}\right),\;\;\;s=0,1,\cdots,\frac{N}{2}
\end{align}
\)

ここで、\(\sum{\prime\prime }\)は和の最初と最後の項に\(\frac{1}{2}\)を乗じることを意味します。

N点の分点に対し、誤差のオーダーは\(O(N+1)\)です。
\(x_k\)はチェビシェフ多項式\(T_n(x)\)の極値点になっています。

利点

・重みの計算に掛かる時間が\(O(N log N)\).
・2N点での積分値を求める際にN点で計算した評価点を用いることが出来る。

ガウス求積法と比べて

・ガウス型の重みの計算に掛かる時間は\(O(N^2)\).ただし、分点と重みは一度求めてしまえばずっと使いまわせるので、さほど重要な問題ではありません。
・ガウス-ルジャンドル求積法の分点を\(N\)点利用し、それに\(N\)点追加することでガウス-クロンロッド求積法が使えます。ガウス-ルジャンドルはN点の求積で\(O(x^{2N})\)、ガウス-クロンロッドは2N点の求積で\(O(x^{3N})\)の多項式まで正確に表すことが出来ます。

結論

本気で計算精度が欲しいならばガウス-ルジャンドル求積法。
簡単に、でもかなりの精度で求めたければクレンショウ-カーティスの求積法なのだと思います。

プログラム

最適化は特にしていませんが、fortran90ではこんな感じになります。

参考文献

[1]P. J. Davis, P. Rabinowitz著、森正武訳 『計算機による数値積分法』p. 131-132,初版 (1980)
[2]LLOYD N. TREFETHEN, IS GAUSS QUADRATURE BETTER THAN CLENSHAW–CURTIS?

Tweet