[論文レビュー] Rotation method for accelerating multiple-spherical Bessel function integrals against a numerical source function

Zachary Slepian, Yin Li+2 more|arXiv (Cornell University)|Nov 29, 2019

Radio Astronomy Observations and Technology参考文献 4被引用数 8

ひとこと要約

この論文は、三角関数の積-和恒等式を用いて、複数の球ベッセル関数（sBF）積分を再定式化することで、1次元FFTを用いた$N_{\rm g}\log N_{\rm g}$スケーリングを可能にする「回転法」という新規アルゴリズムを導入する。これにより、ナイーブな$N_{\rm g}^{p+1}$手法と比較して最大1000–3000倍の高速化が達成され、DESI や LSST のような大規模な宇宙線調査における相関関数および共分散行列の計算の効率が著しく向上する。

ABSTRACT

A common problem in cosmology is to integrate the product of two or more spherical Bessel functions (sBFs) with different configuration-space arguments against the power spectrum or its square, weighted by powers of wavenumber. Naively computing them scales as $N_{ m g}^{p+1}$ with $p$ the number of configuration space arguments and $N_{ m g}$ the grid size, and they cannot be done with Fast Fourier Transforms (FFTs). Here we show that by rewriting the sBFs as sums of products of sine and cosine and then using the product to sum identities, these integrals can then be performed using 1-D FFTs with $N_{ m g} \log N_{ m g}$ scaling. This "rotation" method has the potential to accelerate significantly a number of calculations in cosmology, such as perturbation theory predictions of loop integrals, higher order correlation functions, and analytic templates for correlation function covariance matrices. We implement this approach numerically both in a free-standing, publicly-available extsc{Python} code and within the larger, publicly-available package exttt{mcfit}. The rotation method evaluated with direct integrations already offers a factor of 6-10$ imes$ speed-up over the naive approach in our test cases. Using FFTs, which the rotation method enables, then further improves this to a speed-up of $\sim$$1000-3000 imes$ over the naive approach. The rotation method should be useful in light of upcoming large datasets such as DESI or LSST. In analysing these datasets recomputation of these integrals a substantial number of times, for instance to update perturbation theory predictions or covariance matrices as the input linear power spectrum is changed, will be one piece in a Monte Carlo Markov Chain cosmological parameter search: thus the overall savings from our method should be significant.

研究の動機と目的

• 宇宙論における数値的ソース関数との複数の球ベッセル関数積分を評価する際の計算ボトルネックを解消すること。
• p個の配置空間変数を持つ場合に生じるナイーブな積分手法の$N_{\rm g}^{p+1}$スケーリングを克服すること。
• 通常のFFT高速化と互換性が薄いsBF積分に対しても、高速フーリエ変換（FFT）を効果的に適用できるようにすること。
• 大規模構造調査における相関関数および共分散行列を高精度で数値的に安定して計算するための手法を提供すること。
• 繰り返し積分の再計算を要するMCMCパラメータ推定などの計算負荷の高いタスクを支援すること。

提案手法

• 三角関数の恒等式を用いて、球ベッセル関数を正弦および余弦の和に書き換える。
• 積-和恒等式を適用し、複数の配置空間変数に関する積分を、単一の差と和に関する積分に変換する。
• 配置空間平面における$(a,b)$から$(a-b, a+b)$への45°回転として問題を再定式化し、1次元FFTの適用を可能にする。
• 事前に生成された大きな照会テーブルを避けて、即座に1次元FFTを実行し、主記憶装置内での計算を可能にする。
• 回転変換を用いて、高次元のsBF積分を$N_{\rm g}\log N_{\rm g}$の複雑さを持つ一連の1次元FFTに変換する。
• 公開用のPythonパッケージおよびmcfitフレームワーク内に実装し、広範な利用可能性を実現する。

実験結果

リサーチクエスチョン

• RQ1 数値的ソース関数を伴う複数の球ベッセル関数積分は、ナイーブな$N_{\rm g}^{p+1}$スケーリングを上回る速度向上が可能か？
• RQ2 元来多次元的構造を持つsBF積分を、1次元FFTに適した形に再定式化することは可能か？
• RQ3 回転法は、直接数値積分法や従来の照会ベース手法と比較して、どの程度の高速化を達成できるか？
• RQ4 赤方偏移調査における共分散行列推定のような高精度な宇宙論的応用において、この手法はどの程度の性能を示すか？
• RQ5 数値的精度の損失を伴わずに、sBFの次数$\ell$をより高次に拡張することは可能か？

主な発見

• 回転法により、複数のsBF積分の計算スケーリングが$N_{\rm g}^{p+1}$から$N_{\rm g}\log N_{\rm g}$に低減され、FFTに基づく評価が可能になった。
• テストケースでは直接積分法と比較して6–10倍の高速化が達成され、FFTを組み合わせると最大1000–3000倍の高速化が実現された。
• ナイーブ法との数値的一致性が十分に保たれており、DESI解析のような応用において系統的誤差は無視できる水準にとどまる。
• 大きな照会テーブルを避けることで、効率的な主記憶内計算が可能になり、I/Oボトルネックが軽減された。
• 既に$\ell = 0$および$\ell = 2$のモーメントにおいて有効であり、2PCF共分散行列の主要な構成要素をカバーしている。
• 高次の$\ell$では被積分関数の動的範囲が大きくなるため課題が生じるが、現段階の調査規模の応用には依然として有効である。