[論文レビュー] Practical Implementation of High-Order Multiple Precision Fully Implicit Runge-Kutta Methods with Step Size Control Using Embedded Formula

Tomonori Kouya|arXiv (Cornell University)|Jun 11, 2013

Numerical methods for differential equations参考文献 6被引用数 3

ひとこと要約

本稿では、埋め込み公式とステップサイズ制御を備えた完全陰的ルンゲ・クッタ（IRK）法を用いた高次精度ODEソルバを提案する。混合精度反復改善とW変換に基づく行列低減を活用し、内側の反復を高速化する。実装により、1Dブルセラトルのような剛性・大規模問題において、最大29桁の精度を達成し、従来のDPオンリーのソルバーよりも精度と効率性に優れている。

ABSTRACT

We propose a practical implementation of high-order fully implicit Runge-Kutta(IRK) methods in a multiple precision floating-point environment. Although implementations based on IRK methods in an IEEE754 double precision environment have been reported as RADAU5 developed by Hairer and SPARK3 developed by Jay, they support only 3-stage IRK families. More stages and higher-order IRK formulas must be adopted in order to decrease truncation errors, which become relatively larger than round-off errors in a multiple precision environment. We show that SPARK3 type reduction based on the so-called W-transformation is more effective than the RADAU5 type one for reduction in computational time of inner iteration of a high-order IRK process, and that the mixed precision iterative refinement method is very efficient in a multiple precision floating-point environment. Finally, we show that our implementation based on high-order IRK methods with embedded formulas can derive precise numerical solutions of some ordinary differential equations.

研究の動機と目的

• IEEE754倍精度の限界を超えて、剛性および悪条件問題を扱える実用的な高次精度複数精度ODEソルバの開発。
• 高次IRK法におけるRADAU5型の行列低減の非効率性を解消するため、SPARK3型のW変換を採用し、内側反復の高速化を図る。
• 任意のIRK法に対して埋め込み公式を導出することで、複数精度における適応的ステップサイズ制御を可能にする。
• 特に大規模問題におけるスパースヤコビ行列に対して、混合精度反復改善を用いて内側反復における線形システムの解法を高速化する。
• 1Dブルセラトル（1000個のDAEを含む）を含むベンチマーク剛性ODEにおいて、ソルバの高精度と性能を実証する。

提案手法

• ソルバは、複数精度算術における局所誤差推定とステップサイズ制御のための埋め込み公式を備えた高次完全陰的ルンゲ・クッタ法を用いる。
• 内側反復における線形システムの係数行列を、RADAU5の複素数対角化よりも効率的な実数非対称三重対角形に低減するためにW変換を適用する。
• 混合精度反復改善を採用し、線形システムを倍精度で解き、複数精度で補正することで、解法時間を著しく短縮する。
• バンド行列格納と左前処理付きDP-MP BiCGSTABソルバを用いることで、スパースヤコビアン行列に対応し、大規模問題に適応する。
• CPU最適化キーナルを備えたMPFRとGMPを基盤とする複数精度数値ライブラリBNCpackに基づく実装を行う。
• 主なIRK公式に付随する埋め込み公式により、ステップサイズの適応が可能となり、A安定かつ高次精度の積分が実現される。

実験結果

リサーチクエスチョン

• RQ1 埋め込み公式を備えた高次IRK法を複数精度環境に効果的に実装することで、従来の倍精度ソルバーよりも高い精度を達成できるか？
• RQ2 複数精度算術下で高次IRK法の内側反復を高速化する観点から、SPARK3型のW変換はRADAU5型の低減法よりも効果的か？
• RQ3 混合精度反復改善は、高次IRK法の内側反復における大規模線形システムの解法コストを著しく低減できるか？
• RQ4 1Dブルセラトルのような剛性および大規模問題において、提案されたMP ODEソルバの性能と精度は、既存のソルバーよりも優れているか？
• RQ5 MP ODEソルバにおける主なボトル neck は何か、特に内側反復における線形ソルバの収束性に起因するものであり、それらはどのように緩和できるか？

主な発見

• 1Dブルセラトル問題において20段階・20次精度の解法で、最大相対誤差が$1.3 \times 10^{-29}$に達し、25桁を超える正確な有効数字を達成した。
• N=500の1Dブルセラトルにおいて、20段階法は計算時間を3642秒にまで短縮した（非前処理BiCGSTABでは64573秒）。左前処理の恩恵が明確に示された。
• 混合精度反復改善により、特にスパースヤコビアンのバンド行列格納と組み合わせた場合、線形解法の収束が著しく高速化された。
• 複数精度下での高次IRK法の内側反復時間短縮において、SPARK3型W変換はRADAU5型低減法よりも優れた効果を示した。
• 100桁のMP結果との比較により、20桁を超える真の小数桁数の精度が保証された。
• 埋め込み公式のおかげでステップサイズの適応が安定かつ効率的であり、左前処理付きDP-MP BiCGSTABの使用により収束性が向上し、計算時間が短縮された。