行列転置の深淵：`TRANSPOSE`の甘い罠とキャッシュラインの物理的制約
スパコンの演算性能がペタからエクサへと向かう中、いまだに多くの研究者が「コードが動けばそれでいい」という幻想を抱いている。しかし、我々が扱う数千億自由度のシミュレーションにおいて、演算器（ALU）のスペックなど飾りだ。支配的なのは、常にメモリ階層との対話である。
今回は、最も初歩的に見えて、実はHPCの現場で最も頻繁にボトルネックとなる「行列転置」と、そこに潜む`TRANSPOSE`関数の物理挙動について深掘りする。
—
1. `TRANSPOSE`関数の本質と隠されたコスト
Fortranの組込み関数`TRANSPOSE(A)`。一見、非常にエレガントでモダンな記述だ。しかし、この関数を巨大な疎行列や密行列の計算ループ内で安易に呼ぶことは、キャッシュラインを自ら焼き払う行為に等しい。
`TRANSPOSE`は、関数の戻り値としてテンポラリ配列を生成する。この際、メモリ上のデータ配置は「列優先（Column-major）」の原則に従ってコピーされるわけだが、出力先のキャッシュ状況やアライメントによっては、メモリアクセスが完全に非連続（Strided access）になる。
なぜ「列優先」が神なのか
Fortranのメモリ配置は、第一添字が最も変化が速い。`A(i, j)`と`A(i+1, j)`は隣接しているが、`A(i, j)`と`A(i, j+1)`は、行列のサイズ（先頭次元長：Leading Dimension）だけ離れている。転置を行う際、行を列に変換しようとすると、CPUはキャッシュライン（通常64バイト）の半分以上を捨てながらデータをフェッチせざるを得ない。これが、VTuneで「L1/L2キャッシュミス」が真っ赤に表示される主原因だ。
—
2. キャッシュを味方につける：ブロッキング（Tiling）戦略
巨大行列の転置を最適化する際、単なるループ順序の変更だけでは不十分だ。我々がとるべきは「ブロッキング（タイリング）」である。行列をキャッシュサイズに収まる小さなタイルに分割し、そのタイル内で転置を完結させる。
以下のコードは、キャッシュヒット率を劇的に向上させるための、現場で叩き上げられた実装パターンである。
subroutine blocked_transpose(n, A, B, block_size)
implicit none
integer, intent(in) :: n, block_size
real(8), intent(in) :: A(n, n)
real(8), intent(out) :: B(n, n)
integer :: i, j, ii, jj
!$omp parallel do collapse(2) private(ii, jj, i, j)
do jj = 1, n, block_size
do ii = 1, n, block_size
! タイル内部の転置処理
! この範囲がL1/L2キャッシュに収まるようにblock_sizeを調整する
do j = jj, min(jj + block_size – 1, n)
do i = ii, min(ii + block_size – 1, n)
B(j, i) = A(i, j)
end do
end do
end do
end do
end subroutine blocked_transpose
極限のヒント: `block_size`は、ターゲットとするCPUのL2/L3キャッシュ容量を考慮し、`block_size block_size 8byte`がキャッシュの約半分に収まるようにチューニングせよ。AVX-512等のSIMD命令をコンパイラに自動生成させるには、ループの不変性を維持し、`!dir$ ivdep`などの最適化指示子を活用することが肝要だ。
—
3. スパコン環境でのリンク挙動と最適化フラグ
コンパイラに単に`-O3`を渡して満足しているなら、あなたはまだ素人だ。スパコンのプロファイラ（ScalascaやVTune）を使ってボトルネックを特定したあと、以下のフラグを検討せよ。
- `-xHost` / `-march=native`: 特定のCPUマイクロアーキテクチャの命令セット（AVX-512, FMA3等）をフル活用する。
- `-qopt-prefetch=4`: キャッシュミスを先読みで隠蔽する。ただし、メモリ帯域を飽和させるリスクがあるため、ストリーミング負荷が高い場合は調整が必要。
- `-fno-alias` / `-alias-args`: Fortranではポインタ（ポインタ変数ではなく実引数）のエイリアスを考慮する必要がある。コンパイラに「メモリの重なりはない」と断言させることで、レジスタへのデータ保持を強力に推し進める。
—
4. MPI/OpenMPハイブリッド並列化の深層
数万コア規模のHPC環境では、ノード間通信（MPI）のオーバーヘッドが最大の敵だ。転置処理において、タイル化したデータをMPIで送受信する場合、`MPI_Type_vector`や`MPI_Type_create_subarray`を使って、非連続なデータを一度にパックして送信せよ。
転置は、しばしば「All-to-All」の通信パターンを誘発する。このとき、単なる`MPI_Send/Recv`ではなく、ハードウェアトポロジーを意識したMPI_Alltoallvの利用が必須となる。ノード内部ではOpenMPでの共有メモリ並列を、ノード間では非同期通信（非ブロッキング通信）を組み合わせ、演算と通信のオーバーラップ（Overlap）を狙うのだ。
最後に：泥臭いデバッグの極意
「計算が合わない」のではない。「メモリが壊れている」のだ。
モダンFortranの`allocatable`な動的配列と、レガシーな固定配列の混在は、コンパイラのエイリアス解析を狂わせる。`ISO_FORTRAN_ENV`を使用して型を厳密に定義し、`implicit none`は強制すること。
究極の性能は、コードの美しさではなく、ハードウェアの物理的制約をいかに数学的にねじ伏せるかという、執念の積み重ねから生まれる。理論と実測値の乖離を埋めるその瞬間まで、プロファイラを回し続けろ。それが、真の数値計算アーキテクトの矜持である。
コメント