プログラムの高速化・並列化サービスの事例

＊ここに掲載するのは、STFCラザフォード・アップルトン研究所のDominik B. Jochym博士によるHECToRレポート「Hybrid Time-Dependent Density Functional Theory in CASTEP: a dCSE Project Part 2: Atomic Forces in Excited Electronic States, Dominik B. Jochym, STFC Rutherford Appleton Laboratory, February 13, 2012」を要約したものです。

概要

時間依存密度汎関数法$TDDFT$は、分子系の励起状態の特性の解析に有用な手法で、様々な量子化学コードに実装されています。CASTEPへのTDDFT実装[2]は、英国電子構造コミュニティに対して最先端の問題に対処する機会を提供します。これは例えば、無機あるいは有機物質の光起電性や、表面上の触媒反応、光学ディスプレイのための発光ポリマー物質、フェムト秒レーザー化学などです。CASTEPは特に、HPC領域で優秀な利用実績があり、HECToRシステムで最も利用頻度の高いコードです。ハイブリッド汎関数の実装は、以前にTDDFTで困難であった問題に対処することを可能にします[3,1,5]。
この報告は2年間のプロジェクトの後半部で、Tamm-Dancoff近似の線形応答TDDFTで決定される励起状態のフォース計算の実装に焦点を当てています。繰り返し法ソルバーを用いた励起状態の計算については、レポート：http://www.hector.ac.uk/cse/distributedcse/reports/castep02/を参照してください。

ハイブリッドTDDFT

最初の作業目標は、ハイブリッド汎関数に対するXC応答カーネルを実装する事でした。これについては前回のレポート(http://www.hector.ac.uk/cse/distributedcse/reports/castep02/)を参照してください。要約すれば、励起状態に対するハートレーフォック交換項、およびTDDFTの断熱交換相関カーネルのためのハイブリッド交換相関汎関数が実装されました。

I/Oボトルネック

前回のレポートで示したように、固有値ソルバーがHECToRでテストされましたが、比較的小さいサイズ$50原子程度$の系でも光吸収スペクトル計算が困難な問題に直面していました。特に、大きな状態数$100以上$のリスタートにおいて、100GBオーダーのデータ読み込みが発生します。各状態は、次回の状態との直交性を確保するために、全波動関数と全ての低エネルギー状態を必要とします。この問題は、もう一つのプロジェクト「Boosting the scaling performance of CASTEP: enabling next generation HPC for next generation science」において解決しました。要約すると、800MBのファイルを100個読み込むのに、以前は2.5時間かかっていましたが、改善後は約10倍の15分程度に減少しました。

励起状態フォース

次の目標は、励起状態フォース計算の実装です。原子に働く力の計算は、全エネルギーの原子核座標による微分計算で得ることが出来ます。Hutterの論文[4]に従えば、フォースの計算では、コーン-シャム軌道の変化もまた軌道に影響を与えて微分を複雑にすることを考慮しなくてはなりません。そのラグランジアンから、フォースへは次の3種からの寄与があります：基底状態軌道、線形応答TDDFT軌道、所謂Zベクトル。Zベクトルは、「コーン-シャム軌道の定常性に関するラグランジュ未定乗数の行列」です[4]。このベクトルはレスポンス波動関数と同族であり、Sternheimer様方程式の解が要求されます。基底状態からの寄与は既存のCASTEP内に実装されています。
よって、フォースへのTDDFTからの寄与はTDDFTレスポンス軌道とZベクトルからのフォースの和と成ります。

フォースの計算値の正確性は、原子移動距離に関する励起状態エネルギーの数値微分を用いて検証されました。

· 励起状態の座標最適化

TDDFTフォースを既存の座標最適化に組み込む際には、モジュール依存性に注意が必要でした。CASTEP開発グループとの議論の結果、全エネルギーやフォース計算が必要となった場合は、モジュール「ｔｄｄｆｔ」を「ｅｌｅｃｔｒｏｎｉｃ」よりも低レベルの「ａｌｇｏｒ」から呼び出すこととしました。
Hutterの論文に合わせて、ホルムアルデヒド分子でテストを行いました。Hutterの論文は、静電ポテンシャルの周期イメージ分割にMartyna-Tuckerman法を用いる[6]ため、正確な結果の一致は期待されません。異なる擬ポテンシャルも用いています。こうした相違点もありますが、CASTEPが生成した座標は参照座標とよく一致しました。

性能の観点では、励起状態の最適化は基底状態に比べてよりコスト高です。明らかな理由の一つは、各座標においてTDDFT固有値問題を解かなければならないためです。もう一つの理由は、線形応答理論と連携する際には基底状態の波動関数を決めるためのエネルギー条件がより厳しいことです。これは、応答波動関数の精度を保証するために、基底状態計算の繰り返し数が増加することが要因です。TDDFTフォース計算もまた固有値問題に比べ無視できません。
Zベクトルを得るには、追加の一方向の密度汎関数摂動理論$DFPT$ E-fieldソルバー[7]による計算が必要です。この計算時間はほぼ基底状態計算と同等です。しかしながら、最も時間が掛かるのはTDDFT固有値/固有ベクトルの計算です。この計算は前回の座標を初期値とすることで削減することが可能です。無論、規格直交性に注意が必要ですが、これを実装したところ、TDDFTソルバーの繰り返し数は少なくとも半分に削減されました。前回の固有ベクトルに小さな乱数$0.001$を追加することで安定性を確保できましたが、TDDFTの繰り返し数には影響はありませんでした。ホルムアルデヒドのテストケースでは、第二励起状態の構造最適化で7回のBFGS繰り返し数を要しますが、この実装はランダムな初期試行固有ベクトルに比べ40%削減されました。

補足

CASTEPはよく整備された構造を持っています。元来は専用用途向けに作成されたコードでしたが、少ない修正で一般化できるものでした。また、サブルーチンヘッダとコード全体に豊富なコメントがあるため、コードの修正には極めて役立ちました。
ここで開発されたコードは、CASTEPの将来のメジャーリリースバージョンに含まれる予定です。

謝辞

このプロジェクトは、nAG Ltd.が運営するHECToRの分散計算科学および工学$CSE$サービスの基に実行されました。英国の国立スーパーコンピューティング・サービスである、HECToR：英国リサーチ・カウンシル・ハイエンド計算サービスは、リサーチ・カウンシルを代行するEPSRCが管理しています。そのミッションは英国学術界の科学および工学の研究支援です。HECToRスーパーコンピューターは、UoE HPCx Ltd.およびnAG Ltd.のCSEサポートサービスにより管理運営されています。

文献

[1]	L. Bernasconi, M. Sprik, and J. Hutter, Hartree-Fock exchange in time dependent density functional theory: application to charge transfer excitations in solvated molecular systems, Chemical Physics Letters 394 $2004$, no. 1-3, 141-146.
[2]	S.J. Clark, MD Segall, CJ Pickard, PJ Hasnip, MJ Probert, K. Refson, and MC Payne, First principles methods using CASTEP, Z. Kristallogr 220 $2005$, no. 5-6, 567-570.
[3]	A. Dreuw, J.L. Weisman, and M. Head-Gordon, Long-range charge-transfer excited states in time-dependent density functional theory require non-local exchange, The Journal of Chemical Physics 119 $2003$, 2943.
[4]	J. Hutter, Excited state nuclear forces from the Tamm-Dancoff approximation to time-dependent density functional theory within the plane wave basis set framework, The Journal of Chemical Physics 118 $2003$, 3928.
[5]	A.F. Izmaylov and G.E. Scuseria, Why are time-dependent density functional theory excitations in solids equal to band structure energy gaps for semilocal functionals, and how does nonlocal Hartree-Fock-type exchange introduce excitonic effects?, The Journal of Chemical Physics 129 $2008$, 034101.
[6]	G.J. Martyna and M.E. Tuckerman, A reciprocal space based method for treating long range interactions in ab initio and force-field-based calculations in clusters, The Journal of Chemical Physics 110 $1999$, 2810.
[7]	K. Refson, P.R. Tulip, and S.J. Clark, Variational density-functional perturbation theory for dielectrics and lattice dynamics, Physical Review B 73 $2006$, no. 15, 155114.