クイックスタート: 最初の計算

クイックスタート: 最初の計算#

この章では、空の Python プロンプトから数行で本物の量子化学の結果までたどり着きます。水分子のエネルギーを計算し、手法を DFT に切り替え、原子電荷の物性を読み出します — qc-rs の計算の全体像をつかむのに十分です。

注釈

始める前に qc-rs がインストールされ、インポートできる（import qc が通る）必要があります。まだなら先にインストールとビルドを行ってください。

最初の計算を実行する#

これを Python セッション（またはスクリプト、Jupyter ノートブックのセル）に入力します:

import qc

# 1. チェックポイントを作る:水、cc-pVDZ 基底、座標はオングストローム。
mychk = qc.chk.new(
    atom="O 0.0000 0.0000 0.1173; H 0.0000 0.7572 -0.4692; H 0.0000 -0.7572 -0.4692",
    ao="cc-pvdz",
    unit="angstrom",
)

# 2. 制限 Hartree–Fock（RHF）ステップを追加して実行。
mychk = mychk.scf(ref="r").run()

# 3. 結果を読む。
print(f"energy    = {mychk.scf.energy:.6f} hartree")
print(f"converged = {mychk.scf.converged}")

次のように表示されるはずです:

energy    = -76.026772 hartree
converged = True

これで Hartree–Fock 計算が完結です。おめでとうございます — たった今、水について電子シュレディンガー方程式を（近似的に！）解きました。

いま何が起きたのか？#

3 つの考え方、1 行につき 1 つ:

qc.chk.new(...) がチェックポイントを作る — 分子と、のちに全結果を保持するオブジェクト。ここでは構造（atom=）、基底関数（ao="cc-pvdz"）、座標の単位（unit="angstrom"）を渡しました。
.scf(ref="r") がステップを追加する — 自己無撞着場計算です。ref="r" は制限（閉殻）Hartree–Fock の意味。ステップの追加ではまだ重い計算は行いません。何をしたいかを記録するだけです。
.run() が実際に計算し、その後 mychk.scf.energy で結果を読み出します。

この「組み立ててから実行する」パターンが qc-rs の核心です。中核概念できちんと見ていきます。

Tip

実況を見る:run(log=...) 既定では .run() は沈黙します — 計算するだけです。run(log="stdout") を渡すと、代わりに実況の量子化学風トランスクリプトを流します:システム要約、SCF のサイクルごとのエネルギー表、収束チェック。早い段階で手を伸ばすオプション少数 — log_style="modern"/"orca" は見た目のスタイル、plot=True は SCF 収束曲線をインライン描画（ノートブックで %matplotlib inline）。実行後は mychk.log() で再計算なしにトランスクリプトを再表示できます。全体はエディタ設定とログと出力の章にあります。

Tip

速くする:nthread= 既定では計算は CPU コア 1 つで走ります。run(nthread=8) を渡すと、マシンの 8 コアに仕事を分散できます — ノート PC や単一ノードでほぼすべての実行を速くする、共有メモリのスレッド並列で、追加のセットアップは要りません。（複数のマシンにまたがる MPI プロセス並列 nmpi= はもっと高度で、後の並列計算と HPC で導入します。今は nthread= だけで十分です。）

結果はエネルギー以外もはるかに多く持ちます。例えば mychk.scf.ncycle は収束までの SCF 反復回数、 mychk.scf.energy_components は全エネルギーを物理的な内訳（運動・核引力・クーロン・交換 …）に分解します。必要になったところで登場します。

Tip

Python メモ f"...{mychk.scf.energy:.6f}..." は f-string です — Python が {...} 内の式の値をテキストに埋め込み、 :.6f は小数点以下 6 桁で整形します。

DFT に切り替える#

Hartree–Fock でなく密度汎関数理論を使いたい？ 1 か所だけ変えます — xc= で汎関数を指定:

mydft = qc.chk.new(
    atom="O 0.0000 0.0000 0.1173; H 0.0000 0.7572 -0.4692; H 0.0000 -0.7572 -0.4692",
    ao="cc-pvdz", unit="angstrom",
).scf(ref="r", xc="b3lyp").run()

print(f"{mydft.scf.energy:.6f} hartree")   # -76.420369

B3LYP のエネルギー（-76.420369）が Hartree–Fock より低いのは、電子相関を取り込んでいるからです — これは基礎で展開する概念です。

最初の物性#

エネルギーは始まりにすぎません。同じ収束計算から、数百の物性を求められます。ここでは Mulliken 原子電荷を:

q = qc.prop.chrg.mulliken(mychk)
print(q["charges"])       # [-0.3060, 0.1530, 0.1530]
print(q["atom_labels"])   # ['O', 'H', 'H']

酸素が部分負電荷、各水素が部分正電荷を帯びています — まさに水に期待される極性です。物性スイート全体は分子物性にあります。

開殻分子（ラジカル）#

すべての分子が閉殻ではありません。ラジカル — 電子数が奇数、または不対スピンがある — には非制限 Hartree–Fock （ref="u"）を使い、スピンを設定します。メチルラジカル（CH₃•、ダブレット）の例:

ch3 = qc.chk.new(
    atom="C 0 0 0; H 0 1.079 0; H 0.934 -0.539 0; H -0.934 -0.539 0",
    ao="cc-pvdz", unit="angstrom",
    spin=2,                       # スピン多重度 2S+1: 2 = ダブレット（不対電子 1 個）
).scf(ref="u").run()              # ref="u" = 非制限

print(f"{ch3.scf.energy:.6f} hartree")     # -39.563802
print(qc.prop.spin.s_squared(ch3))         # 0.7612   (理想ダブレット: 0.75)

Tip

スピンについて spin は多重度 2S+1:1 = シングレット（全対）、2 = ダブレット（不対 1）、3 = トリプレット（不対 2）。閉殻には ref="r"、開殻には ref="u"（非制限）または ref="ro"（制限開殻）を使います。

ここでの ⟨S²⟩ ≈ 0.76 は純粋なダブレットの理想値 0.75 に近い値です。わずかな超過はスピン汚染を表し、 SCF の章で説明します。

次はどこへ#

構造を緩和してから解析するには、チュートリアル: DFT 構造最適化 → 物性を進めてください。
チェックポイント/ステップ/実行のパターンがなぜこの形なのかを理解するには、中核概念を。
Hartree–Fock と DFT が実際に何を計算しているのかを理解するには、基礎を。