GANSUで始める量子化学計算: 水素分子の結合距離を求めてみる

水素分子の結合距離を求めてみる

はじめに

前回の記事では、GANSUのインストールと基本的な使い方を紹介しました。
https://zenn.dev/comp_lab/articles/e809e6e01c4be4

今回は、実際に水素分子 ($H2\mathrm{H_2}$) の結合距離を求める計算を行います。

水素分子の結合距離

水素分子の結合距離は、量子化学計算によって求めることができます。具体的には、原子間距離を変えながらエネルギーを計算し、エネルギーが最小になる距離を見つけることで、最も安定な構造を決定します。これは、量子化学計算が「エネルギーを最小にする構造を探す」ことに本質的に関わっていることを直感的に理解する良い例です。

計算の流れ

1. 入力ファイルの準備
   • 水素分子($H2\mathrm{H_2}$)の原子間距離を変えた複数のXYZファイルを作成します。
2. GANSUでエネルギー計算を実行
   • 距離ごとにエネルギーを計算します。
3. エネルギー最小値を探す
   • 計算結果をプロットし、エネルギーが最小になる距離を特定します。

入力ファイルの準備

GANSUでは、XYZ形式のファイルを入力として受け付けます。XYZファイルは、分子を構成する原子の座標を記述したテキストファイルです。以下のような形式で記述されます。

H2.xyz

2
Hydrogen molecule (H2)
H 0.000 0.000 0.000
H 0.740 0.000 0.000

内容は以下の通りです。

• 1行目: 原子数
• 2行目: コメント（自由記述）
• 3行目以降: 原子の種類と3次元座標(x, y, z)をスペース区切りで記述。座標の単位はÅ(オングストローム)です。

この形式で水素分子の原子間距離を変えたファイルを作成し、それぞれのエネルギーを計算します。

GANSUを使ったエネルギー計算

準備したXYZファイルを使って、GANSUでエネルギーを計算します。

./HF_main -x H2.xyz -g ../basis/sto-3g.gbs

このコマンドは、上記XYZファイルで定義した水素分子のエネルギーを計算するものです。これを原子間距離を変えながら計算を繰り返します。

と言いたいところですが、XYZファイルを作成するのが面倒ですよね。そこで、Pythonスクリプトを使って自動的にXYZファイルを生成し、計算を行う方法を紹介します。

Pythonスクリプトによる自動計算

以下のPythonスクリプトを使って、水素分子の原子間距離を変えたXYZファイルを作成し、GANSUでエネルギーを計算し、結果をCSVファイルに保存します。

import subprocess
import csv
import os

def run_calculations(start=0.3, stop=2.0, step=0.05, output_csv="energy_results.csv"):
    results = []
    d = start
    while d <= stop:
        xyz_filename = "H2_temp.xyz"
        with open(xyz_filename, "w") as f:
            f.write("2\n")
            f.write(f"H2 molecule with distance {d} angstroms\n")
            f.write(f"H   0.000   0.000   0.000\n")
            f.write(f"H   {d:.3f}   0.000   0.000\n")

        command = ["./HF_main", "-x", xyz_filename, "-g", "../basis/sto-3g.gbs"]
        try:
            result = subprocess.run(command, capture_output=True, text=True, check=True)
            for line in result.stdout.split("\n"):
                if "Total Energy:" in line:
                    energy = float(line.split()[2])
                    results.append((d, energy))
                    print(f"Distance: {d:.2f} Å, Energy: {energy} hartree")
                    break
        except subprocess.CalledProcessError as e:
            print(f"Error running calculation for distance {d:.2f} Å: {e}")

        os.remove(xyz_filename)
        d += step

    with open(output_csv, "w", newline="") as f:
        writer = csv.writer(f)
        writer.writerow(["Distance (Å)", "Total Energy (hartree)"])
        writer.writerows(results)
    print(f"Results saved to {output_csv}")

# 実行
if __name__ == "__main__":
    run_calculations()

このスクリプトを実行すると、

• 0.3Å から 2.0Å まで 0.05Å 刻みで水素分子のXYZファイルを一時的に作成
• GANSU (HF_main) を実行し、エネルギーを取得
• XYZファイルを削除
• 結果を energy_results.csv に保存

上記スクリプトではHF_mainを繰り返し実行しており、実行ごとにGPUの初期化をするためその時間がかかります。このように繰り返し計算をする際はHF_mainコマンドを繰り返し実行せず、GANSUをライブラリとしてプログラムから直接呼び出す使い方が推奨されます。

結果の解析

上記で得たenergy_results.csvをエクセルなどを用いてプロットすると、エネルギーが最小になる原子間距離がわかります。この実行では0.7Åが最小エネルギーとなり、この値が水素分子の最安定構造、つまり、結合距離がわかります。

水素分子の原子間距離とエネルギーの関係

まとめ

本記事では、GANSUを使って水素分子の結合距離を求める方法を紹介しました。

• XYZ形式を用いて分子の座標を記述し、計算に使用する方法を学びました。
• Pythonスクリプトを利用して、異なる距離でのエネルギー計算を自動化しました。
• 計算結果を分析し、結合距離を決定する方法を理解しました。

次回は、計算結果の見方や可視化方法について紹介します。エネルギーの解釈や、Moldenを使った分子軌道の可視化を通じて、量子化学計算の結果をより直感的に理解できるようにしましょう。

https://zenn.dev/comp_lab/articles/29e73268f402b6