[論文レビュー] A recurrence scheme for least-square optimized polynomials

C. Gebert, I. Montvay|ArXiv.org|Feb 28, 2003

Particle physics theoretical and experimental studies参考文献 2被引用数 3

ひとこと要約

この論文は、フェルミオン格子場理論シミュレーション、特に2段階マルチボソン（TSMB）アルゴリズムにおいて用いられる高次元最小二乗最適化多項式を計算する再帰的スキームを提示する。直交多項式展開と高精度算術を活用することで、軽いフェルミオン質量に相当するゼロ付近の逆累乗関数の正確な近似が可能となり、同時に単精度浮動小数点演算による効率的な評価も可能となる。C言語実装により、数千年にわたる多項式次数まで対応可能である。

ABSTRACT

A recurrence scheme is defined for the numerical determination of high degree polynomial approximations to functions as, for instance, inverse powers near zero. As an example, polynomials needed in the two-step multi-boson (TSMB) algorithm for fermion simulations are considered. For the polynomials needed in TSMB a code in C is provided which is easily applicable to polynomial degrees of several thousands.

研究の動機と目的

• ゼロ付近の逆累乗関数に対する高次元多項式近似を、数値的に安定かつ効率的に計算する手法を開発すること。これはフェルミオン格子場理論シミュレーションにおいて不可欠である。
• 区間 $[\epsilon, \lambda]$ において $\epsilon \ll \lambda$ となるとき、$\lambda/\epsilon$ が $10^6$ を超えるような場合に生じる、$x^{-\alpha}\bar{P}(x)$ のような関数に対する最小二乗最適化問題の悪条件性に対処すること。
• 多項式次数が数千年に達する場合でも、最小限の計算コストでこれらの多項式を計算可能な実用的なC言語実装を提供すること。
• 直交展開技術を用いて、実行時における高精度算術を必要とせず、標準的な32ビット浮動小数点演算で高精度に多項式を評価できるようにすること。

提案手法

• 直交多項式 $\Phi_\mu(x)$ は、係数 $\beta_\mu = -p_\mu / q_\mu$ および $\gamma_{\mu-1} = -q_\mu / q_{\mu-1}$ を用いた三項再帰関係で定義される。ここで $p_\mu$, $q_\mu$, $s_\mu$ は、重み関数 $w(x)$ を用いた区間 $[\epsilon, \lambda]$ 上での積分である。
• 展開係数 $d_{n\nu} = b_\nu / q_\nu$ は、内積 $b_\nu = \int_\epsilon^\lambda w(x)^2 f(x) \Phi_\nu(x) dx$ から計算され、最小二乗偏差ノルム $\delta_n^{\text{min}}$ を最小化する。
• 大規模な $\lambda/\varepsilon$ 比に起因する悪条件行列を避けるために、再帰計算段階で高精度算術が使用される。
• 多項式評価は直交展開 $P_n(x) = \sum_{\nu=0}^n d_{n\nu} \Phi_\nu(x)$ を通じて行われ、単精度浮動小数点演算でも正確な評価が可能である。
• TSMBアルゴリズムでは、3つの主要な多項式を計算する：$P^{(1)}$ は $x^{-\alpha}$ を近似し、$P^{(2)}$ は $P^{(1)}$ に多項式重みを乗じたものを近似し、$P^{(3)}$ は $P^{(2)}$ の逆平方根を近似する。これらはニュートン型反復を用いた反復的改善により得られる。
• Cコードは複数の多項式タイプと設定をサポートしており、Neubergerの公式または粗い初期近似からのニュートン反復を用いた $P^{(3)}$ の反復的改善が可能であり、調整可能な区間境界と重み関数を備えている。

実験結果

リサーチクエスチョン

• RQ1 \varepsilon \sim 10^{-6} かつ \lambda \sim 1 のとき、$x^{-\alpha}$ のような関数の条件数が非常に大きい場合、ゼロ付近の高次元最小二乗最適化多項式をどのように安定して計算できるか？
• RQ2 高精度算術を用いて再帰的計算を行うことで、標準的な32ビット浮動小数点演算で効率的に評価可能な、このような多項式の計算を可能にする再帰的スキームは何か？
• RQ3 TSMBアルゴリズムにおける補正不能な誤差を低減するために、区間下限部での $P^{(3)} \approx (P^{(2)})^{-1/2}$ の近似をどのように改善できるか？
• RQ4 大きな $K$ に対して、高価な Neubergerの公式に依存せずに、$P^{(3)}$ の多項式近似を効率的に改善する反復法は何か？
• RQ5 低次の初期推定値を用いたニュートン型反復により、計算コストを最小限に抑えつつ高精度な結果を得ることは可能か？

主な発見

• 再帰的スキームにより、条件数が $10^6$ を超える場合でも、数千年に達する多項式次数まで、高精度算術を再帰段階で使用することで、安定した計算が可能となった。
• 直交多項式展開の使用により、32ビット浮動小数点演算のみで多項式の正確な評価が可能となり、実行時間コストが顕著に削減された。
• 第3の多項式 $P^{(3)}$ については、粗い初期近似（例：$P^{(2)'}$）からのニュートン型反復が2〜3回の反復で高精度な結果を達成でき、高精度な状況では Neubergerの公式の使用が不要となった。
• C実装は複数の多項式タイプと設定をサポートしており、$P^{(1)}$、$P^{(2)}$、$P^{(3)}$ のそれぞれの再帰係数と多項式パラメータを別々に出力ファイルに保存している。
• 高精度計算に必要な桁数は多項式次数に線形に比例する：第1多項式では約 $40 + 1.6n_1$ 桁、第2および第4多項式では約 $70 + 2.2n_{2,4}$ 桁であり、内挿積分の複雑さに起因する。
• 近似結果は機械精度に達しており、$P^{(5)}$ のシーケンス（$P^{(3)}$ の反復的改善）は初期推定値に依存せず数回の反復で収束を示し、ニュートン反復アプローチの強固さが裏付けられた。