スレッド・BLAS・LAPACK

スレッド・BLAS・LAPACK#

最初で最も易しい高速化は、1 台のマシンの複数コア — 共有メモリのスレッド（入門の第 1 モデル）を使います。本章は必要な唯一のノブ run(nthread=N) を示し、qc-rs が実際に何をスレッド化するかを説明し、その下で重い線形代数を担う BLAS ライブラリを解き明かします。

唯一のノブ:`run(nthread=N)`#

計算が使ってよい CPU コア数を、.run() の nthread= で設定します:

import qc
water = "O 0 0 0.117; H 0 0.757 -0.469; H 0 -0.757 -0.469"

m = qc.chk.new(atom=water, ao="cc-pvdz", unit="angstrom").scf(ref="r").run(nthread=8)

nthread は**（MPI）ランクごとのコア数**です:qc-rs カーネルのワーカースレッドと BLAS/LAPACK のスレッド数を、環境変数の調整なしに両方まとめて設定します。数は nthread 引数 > QC_THREADS 環境変数 > 利用可能な全コア の順で解決します。QC_THREADS は nthread を渡せないとき（例:自分で制御できないスクリプト）の外部フォールバックです:

export QC_THREADS=8      # run(nthread=...) が与えられないときに使われる

重要なのは — index で約束したとおり — スレッドが増えても変わるのは速さだけで、答えではないこと:

e1 = qc.chk.new(atom=water, ao="cc-pvdz", unit="angstrom").scf(ref="r").run(nthread=1).scf.energy
e4 = qc.chk.new(atom=water, ao="cc-pvdz", unit="angstrom").scf(ref="r").run(nthread=4).scf.energy
# e1 == e4 == -76.026794   （ビット単位で同一）

qc-rs が実際にスレッド化するもの#

すべてのステップが等しく恩恵を受けるわけではありません（Amdahl を思い出して）。現在の qc-rs では:

スレッド化:KS-DFT 交換相関グリッド — グリッドブロックごとのループ（LDA・GGA・meta-GGA、RKS・UKS とも）が std::thread ワーカーで走り、各々が自分のスクラッチを持ち、共有の V_xc に簡約します。主要なスレッド並列領域です。
BLAS 経由でスレッド化:密な線形代数 — Fock の対角化、大きな行列積 — がマルチスレッド BLAS/LAPACK ライブラリ（下記）で走ります。
まだ逐次:1 電子・2 電子の積分組み立て。小さな分子ではこれが支配しうるので、小さな水のジョブがコアを増やして大きく速くなるとは期待しないでください — そのボトルネックは（まだ）スレッド化されていません。

1 つの大きなスクラッチバッファ（Workspace）が並列ループの前にスレッドごとのレーンに切り分けられるので、スレッドは 1 つの物理確保を反復ごとのメモリ churn なしに共有します — qc-rs のスレッド化が安価な理由です。

下層:BLAS と LAPACK#

ほぼすべての量子化学プログラムは、行列積と固有値解法に BLAS（Basic Linear Algebra Subprograms）と LAPACK に頼ります — 数十年チューニングされた、密な代数の真のエンジンです。qc-rs は次の 2 つのいずれかをリンクします:

ライブラリ	いつ
Intel MKL	Intel/AMD CPU、特に HPC クラスタ（そこでは通常最速）
OpenBLAS	オープンソースの既定;移植性が高く、Intel ツールチェーン不要

ビルド時に 1 つを選びました（intel-mkl-system 対 openblas-* 機能）。どちらもマルチスレッドで、run(nthread=) がそのスレッド数を設定します。

重要

Linux/HPC での MKL スレッド層 MKL ビルドはスレッド層を固定する必要があり、さもないとインポートエラー（undefined symbol: omp_get_num_procs）や、静かに誤ったスレッド化になります。環境で MKL_THREADING_LAYER=GNU（と LP64 インターフェース MKL_INTERFACE_LAYER=LP64）を設定してください — make setup がこれらを .vscode/my.env に書き、シェルの env.sh がエクスポートします。HPC セットアップで最もよくある落とし穴です。

オーバーサブスクリプション:唯一の罠#

オーバーサブスクリプションは、コア数より多くのスレッドを要求したときに起こります — 例えば \(W\) 個の qc-rs ワーカースレッドが各々 \(M\) 個の BLAS スレッドを使うルーチンを呼ぶと、\(W \times M\) 個のスレッドがコアを奪い合います。CPU が thrash して、速くなるどころか遅くなります。qc-rs は設計上これを避けます（スレッド化された XC 領域は単一スレッド GEMM を使うので、オーバーサブスクライブする入れ子のマルチスレッド BLAS がない）し、 run(nthread=) が 1 つの整合した数を設定します。あなたへの規則:MPI と組み合わせるとき（次章）、 ノードあたりランク数 × nthread ≤ ノードあたりコア数 を保つこと。

総合例:コアは答えを変えない#

import qc
water = "O 0 0 0.117; H 0 0.757 -0.469; H 0 -0.757 -0.469"

for n in (1, 2, 4):
    e = qc.chk.new(atom=water, ao="cc-pvdz", unit="angstrom").scf(ref="r").run(nthread=n).scf.energy
    print(f"nthread={n}: E = {e:.6f}")
# nthread=1: E = -76.026794
# nthread=2: E = -76.026794
# nthread=4: E = -76.026794

エネルギーは同一です;十分大きな DFT ジョブなら実時間は nthread とともに下がります（測って — 比だけでなく絶対秒数を報告してください）。

練習 17

水/cc-pVDZ の RHF に run(nthread=16) を設定したのに nthread=2 と変わらない速さでした。qc-rs は壊れていますか。考えられる理由は。
64 コアのノードで 8 MPI ランクを走らせ、オーバーサブスクライブせずにノード全体を使いたいです。どの nthread を渡しますか。
MKL ビルドで import qc が undefined symbol: omp_get_num_procs で失敗します。対処は。

解答練習 17

壊れていません。水/cc-pVDZ は小さく、そのボトルネックはスレッド化された XC/BLAS 領域ではなく逐次の積分組み立て（Amdahl）です — 16 コアで恩恵を受けるほどの並列仕事がありません。大きな分子 / DFT グリッドはスケールします。
nthread=8 — 8 ランク × 8 スレッド = 64 = ノードあたりコア数、オーバーサブスクライブなしでノードを完全に使います。
環境で MKL_THREADING_LAYER=GNU（と MKL_INTERFACE_LAYER=LP64）を設定;make setup がこれらを .vscode/my.env / env.sh に入れます。

スレッドは 1 台のマシン内でスケールさせます。それを超えて — 多数のマシンにわたって — 行くには MPI と相互接続が必要です。