【PySCF】Pythonで始める量子化学計算【一点計算】

【PySCF】Pythonで始める量子化学計算【一点計算】

この記事では、Pythonを使って手軽に始められる量子化学計算ツール「PySCF」について解説します。無料で利用できるPySCFを使いこなす方法と、主に一点計算によるエネルギー計算を行う手順について詳しく解説しています。機械学習や化学分野での応用に役立つため、PySCFの基本操作をマスターすることで、研究やプロジェクトに直接活用できる知識が得ることが出来ます。

動作検証済み環境

Google colab, pyscf 2.5.0

量子化学計算用のソフトウェアとして、PySCFは高い汎用性と優れた計算性能により、広く利用されています。ここでは、PySCFについて、その基本情報から特徴的な機能までを紹介します。

1. 量子化学計算用のソフトウェア: PySCFとは


PySCF(Python Simulations of Chemistry Framework)は、Pythonで書かれたオープンソースの量子化学計算フレームワークです。その主な目的は、量子化学計算のための柔軟で使いやすいツールを提供することにあります。PySCFは、分子システムや凝縮系の電子構造計算に広く使用されており、研究者や開発者による新しい計算手法の開発にも適しています[1]。

参考文献[1]から引用

2. PySCFの機能と特徴


PySCFの特徴は、その豊富な機能性と高い計算効率にあります。

2.1 機能の概要

PySCFは、ハートリー・フォック法、密度汎関数理論(DFT)、多配置自己無撞着場(MCSCF)、結合クラスター理論(CCSD)など、多様な計算手法を提供します。これにより、分子だけでなく、固体を対象とした平面波に基づく計算や有効核ポテンシャルを使用した計算など、表1のように、特定の研究ニーズに応じた柔軟な計算が可能です。

表1: PySCF+拡張パッケージで出来ること(ver.1.3リリース時)[3]

方法分子固体コメント
HF約5000 AO ※1
MP2約1500 MO ※1
CCSD約1500 MO ※1
EOM-CCSD約1500 MO ※1
CCSD(T)約1500 MO ※1
MCSCF〇※2約3000 AO ※1、30-50アクティブ軌道 ※3
MRPT〇※2約1500 MO^1、330-50アクティブ軌道 ※3
DFT約5000 AO ※1
TDDFT〇※2約5000 AO ※1
CISD〇※2約1500 MO ※1
FCI〇※2約(18e, 18o) ※4
LocalizerIAO, NAO, meta-Löwdin, Boys–Foster, Edmiston–Ruedenberg, Pipek–Mezey
Relativityすべての方法に対するECPおよびスカラー相対論的補正、2成分、4成分法
GradientsHF、DFT、CCSD、CCSD(T)、TDDFT
HessianHFおよびDFT
Properties非相対論的、4成分相対論的NMR
SymmetryD2hとサブグループ
AO, MO積分1電子積分、2電子積分
密度フィッティングHF、DFT、MP2

注釈:
※1 AO: 原子軌道、MO: 分子軌道
※2 サポートされているが、特定の機能に限定される
※3 アクティブ軌道:シミュレーションにおいて電子相関が考慮される軌道の範囲であり、pyscfではDMRG、SHCI、またはFCIQMCプログラムをアクティブ空間ソルバーとして使用します
※4 18e, 18o: 18電子系、18軌道系

CISD: シングルとダブルを持つ構成間相互作用、DFT: 密度汎関数理論、ECP: 有効コアポテンシャル、FCI: 完全構成相互作用、HF: Hartree-Fock、IAO: 固有原子軌道、MCSCF: 多配置自己無撞着場、MP2: Møller-Plesset二次摂動理論、NAO: 自然原子軌道、TDDFT: 時間依存DFT。

2.2 Gaussianとの比較

Gaussianは量子化学計算の分野で最も用いられているソフトウエアです。Gaussianは長年にわたって化学分野で盛んに利用されており、信頼性が高いソフトウエアである一方で、ライセンスは非常に高価であるという一面があります。

オープンソースソフトウェアのPySCFは、Gaussianとの比較で水分子のエネルギー計算において絶対値についても同等の結果を得るなど、Gaussianと比べても計算の精度面で高い信頼性を示すことも特徴です。

2.3 並列計算のサポート

PySCFは並列計算をフルサポートし、複数のCPUやGPUを利用して計算時間を短縮できます。これにより、効率的な大規模計算を可能にします。

2.4 開発言語としてのPythonとC

計算のコア部分にはC言語を用いつつも、Pythonで書かれているため、PythonとC言語で書かれており、使いやすさと高速処理を両立しています。このバランスはPySCFを特に魅力的な選択肢にしています。

2.5 拡張性

高い拡張性を持つPySCFは、ユーザーが独自の計算手法を簡単に実装できるように設計されています。そのためPySCFは、計算のシンプルさ、汎用性、効率性を重視した開発が進められており、これらの特徴により、量子化学計算分野でのその使用範囲は広がり続けています。開発者や研究者が直面する多様な課題に対して柔軟に対応できる設計が、PySCFのユニークな特徴です。

2.6 商用利用が可能

オープンソースソフトウェアのPySCFはApache License 2.0というライセンスのもとで公開されています。このライセンスでは、ソフトウェアを商業目的を含むあらゆる目的で使用、再配布することが許可されています。ただし、原著作者のクレジットを表示し、ライセンスのコピーを提供するなどの条件を満たす必要があります。商用利用に際しても、追加の料金は発生しません。詳細については、Apache License 2.0をご覧ください:Apache License 2.0

3. PySCFの準備


では早速一点計算の準備を行っていきます。今回はGoogle Colabを使って実行します。

Google Colabは、ブラウザベースのPythonプログラミング環境であり、PySCFを使用するのに適したプラットフォームの一つです。Google Colabを使用することで、ローカル環境のセットアップなしにPySCFの計算をすぐに開始することができます。

Google colabのリンク:https://colab.research.google.com/?hl=ja

ColabノートブックでPySCFを使用するには、まず以下のコマンドをセルに入力して実行し、PySCFをインストールします。ローカル環境でインストールする際はを取って実行してください。

!pip install pyscf

インストールが完了したら、次のようにPySCFをインポートして、必要な計算を行うスクリプトを書き始めることができます。

この簡単な例では、水素とフッ素からなるHF分子のHartree-Fock計算を実行しています。

from pyscf import gto, scf

# 分子の定義
mol = gto.M(atom='H 0 0 0; F 0 0 1.1', basis='ccpvdz')

# Hartree-Fock計算の実行
mf = scf.RHF(mol)
energy = mf.kernel()
print('Energy:', energy)


実はこれだけで計算を行うことが出来ます。それでは入力ファイルの中身を詳しく見ていきます。

4. 入力ファイルの構造


PySCFを利用する際、Pythonスクリプトによって計算の全体像を定義することが出来ます。これには計算環境の設定、分子構造の指定、計算方法の設定が含まれます。密度汎関数理論(DFT)を例に、これらのステップを詳しく解説します。

4.1 計算環境の設定

計算を開始する前に、必要なPySCFのパッケージをインポートします。計算にはgtoモジュール(分子の構造を定義)とdftモジュール(DFT計算を実施)が必要です。

from pyscf import gto, dft

4.2 分子構造の指定

gto.M関数を用いて、PySCFで計算対象の分子構造を定義します。この関数は、分子の原子座標のセットアップ、使用する基底関数セットの選択、分子の全体的な電荷、そしてスピン多重度などの重要な計算パラメータの指定に利用されます。具体的には、以下のように水分子(H2O)の構造を設定することができます。

# 分子の定義
mol = gto.M(
    atom='''
        O 0.0 0.0 0.1173
        H 0.0 0.7572 -0.4692
        H 0.0 -0.7572 -0.4692
    ''',
    basis='cc-pVDZ',
    charge=0,
    spin=0
)

この設定例では、水分子が中性であることを示すために電荷を0に設定し、未対電子がない状態のシングレットを示すためにスピン多重度を0に設定しています。

  • 基底関数 基底関数セットはcc-pVDZという一般的なセットを使用していますが、計算の精度や目的に応じて、6-311+G(d)のような他の基底関数セットを使用することもできます。
  • スピン 注意点として、PySCFではスピンの指定にspin属性を使用しますが、これは分子の未対電子の総数(2S)であり、伝統的な量子化学の2S+1とは異なります。したがって、ダブレット状態の分子(1つの未対電子)ではspin=1、トリプレット状態(2つの未対電子)ではspin=2を使用します。
  • 電荷 計算対象の分子が陽イオンまたは陰イオンである場合は、適切にcharge属性を調整する必要があります。例えば、陽イオン(例:H3O+)であればcharge=1、陰イオン(例:OH-)であればcharge=-1と設定します。

4.3 計算方法の設定

PySCFでの計算方法の設定は、使用する理論や手法に基づいて異なります。以下に、密度汎関数理論(DFT)を用いる場合の設定方法を例として示します。

DFT計算を行うには、dft.RKSクラスを使用します。このクラスは、制限Kohn-Sham DFTを実行するためのものです。計算に使用する交換相関汎関数は、xc属性を通じて指定します。ここでは、広く使用されているB3LYP汎関数を例に設定方法を説明します。

from pyscf import gto, dft

# 分子の定義
mol = gto.M(atom='O 0 0 0; H 0 -0.757 0.587; H 0 0.757 0.587', basis='cc-pVDZ')

# DFT計算の設定
mf = dft.RKS(mol)  # RKSオブジェクトの生成
mf.xc = 'B3LYP'    # B3LYP汎関数を使用
energy = mf.kernel()  # 計算の実行

print('DFT Energy with B3LYP functional:', energy)

このスクリプトは、水分子のエネルギーをB3LYP汎関数を用いてDFTで計算します。mf.kernel()は計算を実行するメソッドで、計算から得られたエネルギーを返します。

5. 計算の実行と結果の確認


それでは最後に、量子化学計算ソフトウェアPySCFを使用して、分子のエネルギー計算を実行していきましょう。ここでは、PySCFを用いた密度汎関数理論(DFT)計算の一例として、ここまで解説した水分子(H2O)の計算を行い、計算結果を確認する方法を説明します。

計算の実行

PySCFを使用した計算の実行は、Pythonスクリプトを介して行われます。Google colab上での計算の実行はパッケージのインストールと計算のセットアップに分けて以下のように行いましょう。

!pip install pyscf
from pyscf import gto, dft

# 分子の定義
mol = gto.M(
    atom='''O 0.0 0.0 0.0; H 0.0 -0.757 0.587; H 0.0 0.757 0.587''',
    basis='cc-pVDZ'
)

# 計算方法の選択(DFT)と実行
mf = dft.RKS(mol)
mf.xc = 'B3LYP'
energy = mf.kernel()
print('Energy:', energy)


結果の確認

上記のDFT計算をPySCFで実行した後の出力は以下のようになります。

converged SCF energy = -76.4203783695937
Energy: -76.42037836959366

この出力は、計算が適切に収束したこと、そして得られた水分子の全エネルギーが-76.42037836959362 Hartreeであることを示しています。このエネルギー値は、分子の基底状態のエネルギーを反映しており、分子の安定性や反応性の初期評価に重要なデータを提供します。

エネルギーをより一般的な単位である電子ボルト(eV)に変換するには、以下のコードを使用します。

# ハートリーから電子ボルトへの換算係数
hartree_to_ev = 27.2114

# エネルギーをeVに変換
energy_ev = energy * hartree_to_ev
print('Energy:', energy_ev, 'eV')

実行結果は以下の通りです。

Energy: -2079.505483966361 eV

計算実行する度に微小な変動は見られますが、化学的な意義を考えると、これらの変動はeV単位で10^-10オーダーとなり、実用上の計算においては無視できる範囲です。したがって、このエネルギー値は分子の特性を理解する上で信頼性の高い参考値として機能します。

6. 計算レベルMP2, CCSD, MCSCFの使い方


ここでは、水分子に対して行えるより高精度な一点計算の種類と、それぞれの計算手法をどのように実行するのかについて解説します。

MP2計算

MP2(2次のMøller-Plesset摂動理論)は、基本的な電子相関を考慮した手法で、特にHartree-Fock手法で得られた結果を改善したい場合に適しています。単純な分子や中間体の反応エネルギーを計算する際に有効で、比較的計算コストが低いため、広範なシステムに適用可能です。

from pyscf import gto, scf, mp

# 分子を定義:指定されたジオメトリと基底セットを使用して水分子
mol = gto.M(atom='O 0 0 0; H 0 -0.757 0.587; H 0 0.757 0.587', basis='cc-pVDZ')

# ハートリーフォック計算を実行(MP2計算の前提条件)
mf = scf.RHF(mol)
mf.kernel()

# HF計算の結果を使用してMP2計算を設定し実行
mp2 = mp.MP2(mf)
mp2.kernel()

# MP2レベルで計算された全エネルギーを出力
print('MP2 Energy:', mp2.e_tot)


CCSD計算

CCSD(結合クラスターシングルズ・ダブルズ)は、より高度な電子相関を扱うことができる手法です。分子間相互作用や微妙なエネルギー差を求めたい場合に適しています。計算コストはMP2よりも高いですが、高精度な結果が得られるため、重要な化学反応の研究や複雑な分子系の分析に使用されます。

from pyscf import gto, scf, cc

mol = gto.M(atom='O 0 0 0; H 0 -0.757 0.587; H 0 0.757 0.587', basis='cc-pVDZ')
mf = scf.RHF(mol).run()
ccsd = cc.CCSD(mf).run()
print('CCSD Energy:', ccsd.e_tot)

MCSCF計算

MCSCF(多配置自己無撞着場)は、電子構造の変化が大きい化学反応や励起状態の研究に適しています。特に電子の再配置が伴う反応メカニズムの解析や、電子相関が顕著に影響する系の計算に使用されます。計算コストは高いですが、反応経路や励起状態の詳細な解析が可能です。

from pyscf import gto, scf, mcscf

mol = gto.M(atom='O 0 0 0; H 0 -0.757 0.587; H 0 0.757 0.587', basis='cc-pVDZ')
mf = scf.RHF(mol).run()
mc = mcscf.CASSCF(mf, 4, 4).run()
print('MCSCF Energy:', mc.e_tot)

これらの計算手法は、それぞれ異なる状況や目的に適しています。基本的な相互作用や反応エネルギーの計算にはMP2、より精密なエネルギー差が求められる場合にはCCSD、電子の再配置が重要な場合にはMCSCFを選択すると良いでしょう。計算の目的とコストを考慮しつつ、最適な計算手法を選択してください。

最後に


いかがでしたでしょうか!PySCFとの最初の出会いは、その使い勝手の良さに本当に驚かされました。私たちが何か新しいことを始めるとき、その第一歩はしばしば難しいものです。しかし、PySCFの直感的なインターフェースと柔軟性は、量子化学の複雑さを大きく低減させてくれると思います。是非気軽に試してみてくださ

参考文献


[1] PySCF 2.4 公式サイト, https://pyscf.org/, (参照2024-04-10)

[2] Q. Sun, X. Zhang, S. Banerjee, P. Bao, M. Barbry, N. S. Blunt, N. A. Bogdanov, G. H. Booth, J. Chen, Z.-H. Cui, J. J. Eriksen, Y. Gao, S. Guo, J. Hermann, M. R. Hermes, K. Koh, P. Koval, S. Lehtola, Z. Li, J. Liu, N. Mardirossian, J. D. McClain, M. Motta, B. Mussard, H. Q. Pham, A. Pulkin, W. Purwanto, P. J. Robinson, E. Ronca, E. R. Sayfutyarova, M. Scheurer, H. F. Schurkus, J. E. T. Smith, C. Sun, S.-N. Sun, S. Upadhyay, L. K. Wagner, X. Wang, A. White, J. Daniel Whitfield, M. J. Williamson, S. Wouters, J. Yang, J. M. Yu, T. Zhu, T. C. Berkelbach, S. Sharma, A. Yu. Sokolov, and G. K.-L. Chan, Recent developments in the PySCF program package, J. Chem. Phys. 153, 024109 (2020)

[3] Q. Sun, T. C. Berkelbach, N. S. Blunt, G. H. Booth, S. Guo, Z. Li, J. Liu, J. McClain, S. Sharma, S. Wouters, and G. K.-L. Chan, PySCF: the Python-based simulations of chemistry framework, WIREs Comput. Mol. Sci. 8, e1340 (2018)

[4]Ulf Ekström, Lucas Visscher, Radovan Bast, Andreas J. Thorvaldsen, and Kenneth Ruud, Arbitrary-order density functional response theory from automatic differentiation, J. Chem. Theory Comput. 6, 1971 (2010).

[5] Susi Lehtola, Conrad Steigemann, Micael J.T. Oliveira, and Miguel A.L. Marques, Recent developments in libxc — A comprehensive library of functionals for density functional theory, SoftwareX 7, 1 (2018).

2 COMMENTS

ORDE

MCSCFとMRPTの30-50アクティブ軌道というのはDMRG-CASSCFやDMRG-NEVPT2などですよね?
表1には拡張パッケージで出来ること,とは書いてありますが,DMRGと記載していないと初心者に誤解を与えかねない気がしました。

返信する
cliuses

ORDE様

ご連絡ありがとうございます!おっしゃる通り、表記にDMRGを使用したアクティブ軌道であることを追加すべきでした。ご意見を基に記事の変更を行いました。混乱を招く可能性がある点をご指摘いただき、誠にありがとうございました。

返信する

コメントを残す

メールアドレスが公開されることはありません。 が付いている欄は必須項目です