[論文レビュー] A GEMM-based direct solver for finite-difference Poisson problems in non-uniform grids

Pedro Filipe da Silva Costa, Duarte Palancha+3 more|arXiv (Cornell University)|Mar 10, 2026

Parallel Computing and Optimization Techniques被引用数 0

ひとこと要約

この論文は、2方向で固有分解により演算子を対角化し、3番目の方向を直接解く3D非均一格子の直接ポアソン解法を提案する。前方/後方変換はCPU/GPU効率のためGEMMとして実装されている。

ABSTRACT

We present a direct Poisson solver for massively parallel simulations on three-dimensional Cartesian grids with non-uniform spacing. The method uses a tensor-based formulation in which the operator is diagonalized numerically along two directions through one-dimensional eigendecompositions, while the third direction is solved directly. The resulting dense transforms are evaluated efficiently as GEMMs (General Matrix--Matrix Multiplications), allowing many independent one-dimensional operations to be combined into matrix-matrix products that map well to modern CPU and GPU hardware. For uniform grids, the method reduces to the classical eigenfunction-expansion approach, and it naturally supports hybrid combinations of FFT-based and GEMM-based transforms depending on grid uniformity. After coupling the solver to an incompressible Navier-Stokes code, we assess its accuracy and performance against geometric multigrid and block cyclic reduction with FFT diagonalization. The results show that the proposed method is robust and consistently achieves the best time-to-solution. In strong scaling, the more compute-intensive GEMM-based variants attain higher parallel efficiency by better amortizing communication costs, while weak scaling highlights the expected trade-off between FFT-based and dense-transform formulations. Overall, the method enables efficient high-resolution stretched-mesh simulations on modern heterogeneous systems.

研究の動機と目的

• DNSの不可压縮性を保証する高速ポアソン解法の必要性を動機付ける。
• 非均一格子上で固有関数ベースの直接解法を拡張しつつ、GPU向けカーネルを維持する。
• 既存のCaNSインフラと統合し、ハイブリッド合成（均一方向はFFT、非均一方向はGEMM）をサポートするGEMMベースの実装を開発する。
• 最先端手法（幾何学的多重格子法およびFFTを用いたブロック巡回還元）と比較して性能を評価する。
• CPUとGPUで高解像度の伸長メッシュにおける強・弱スケーリングを Demonstrate する。

提案手法

• 非均一格子上の1Dポアソン演算子を三重対角系として定式化し、対角スケーリングDにより対称代替T̃ = D^{1/2} T D^{-1/2}を得て対称化する。
• 対称代替T̃ = Q̃ Λ Q̃^Tを固有分解し、Q = D^{-1/2} Q̃からQを回収し、Q^{-1} = Q̃^T D^{1/2}を得る。
• 3Dへ拡張するにはKronecker和を用い、T = T_x ⊗ I ⊗ I + I ⊗ T_y ⊗ I + I ⊗ I ⊗ T_zとする。xとy方向（T_x, T_y）を対角化し、z方向に沿って独立した1D三対角問題の集合を解く。
• 前方固有基底変換と後方逆変換はxおよびy方向でGEMM（均一方向はFFT）を用いて実行する。z方向では直接三対角解法（TDMA）を用いる。
• この手法はFFTベースの対角化と同じドメイン分割と通信パターンを維持し、均一格子にはFFT、非均一格子にはGEMMを組み合わせたハイブリッド合成を可能にする。
• ポアソン解法はナビエ・ストークス解法と組み合わせ、分数ステップのRK3スキームと、必要に応じて同じポアソン系/ヘルムホルツ系で解く暗黙拡散を統合する。

実験結果

リサーチクエスチョン

• RQ1 GEMMベースの二方向対角化が、FFTベースおよび多重格子法と比較して、非均一格子上の3Dポアソン問題の解時間に対して競争力または優位性を持つか。
• RQ2 伸長メッシュでのCPUおよびGPUでの強・弱スケーリングのパフォーマンスはどうか。
• RQ3 問題サイズが大きくなると、柔軟性（非均一格子）と漸近コストのトレードオフはどうなるか。
• RQ4 非圧縮性ナビエ–ストークスワークフローへの統合と、現代のハードウェア上での非均一格子機能のオーバーヘッドはどの程度か。

主な発見

• GEMMベースのアプローチは、試験した代替案の中でCPUコア1つあたりの最良の解時間を提供し、特に強度の高い勾配格子で有利である。
• CPUとGPUでの強いスケーリングは、GEMMベースの計算集約型変種が転置・通信のオーバーヘッドをFFT中心の変種よりも良く償却する。
• 弱いスケーリングではトレードオフが見られ、成長する方向でFFTを使用すると時間の増加が穏やになる一方、GEMMベースの変換は解像度が上がるとコストが増加する。
• GPUでは非均一格子機能を完全に有効にするとオーバーヘッドは中程度であり、多数GPUでの強いスケーリングを維持する。
• 総じて、この手法は伸長メッシュでの乱流流れシミュレーションの壁時計時間を均一格子実行と比較して大幅に削減でき、ヘテロジニアスなハードウェアでも堅牢な性能を示す。