(57)【特許請求の範囲】
【請求項１】
コンピュータにより、ＴＤＤＦＴ法によって物質の電子状態を計算する電子状態計算方法であって、
前記物質の電子状態の初期状態を示す情報を含む前記ＴＤＤＦＴ法によって前記物質の電子状態を計算するためのパラメータを入力して、前記ＴＤＤＦＴ法における時間ステップｊ毎の、前記物質の電子状態のハミルトニアンＨ（ｊ）におけるハートレー項Ｖｈ（ｊ）を含む当該物質の電子状態を示す情報を記憶する記憶手段に記憶させる入力ステップと、
時間ステップｎにおいて、時間ステップｎ−１から時間ステップｎ−Ｉ−１までの前記ハートレー項Ｖｈ（ｎ−１）〜Ｖｈ（ｎ−Ｉ−１）を前記記憶手段から読み出して、予め設定されるＩ個の係数α１〜αＩを用いて、前記時間ステップｎにおけるハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）を予め記憶した式
【数１】
0005098480.tif 000034
により演算をして算出する初期値算出ステップと、
前記初期値算出ステップにおいて算出された前記ハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）から、予め記憶したハートレー項算出ルールに基づき差分法により当該ハートレー項Ｖｈ（ｎ）を算出すると共に、算出したハートレー項Ｖｈ（ｎ）を前記記憶手段に記憶させるハートレー項算出ステップと、
前記ハートレー項算出ステップにおいて算出された前記物質の前記ハートレー項Ｖｈ（ｎ）及び前記記憶手段により記憶された情報から、予め記憶した電子状態算出ルールに基づき前記物質の電子状態を算出すると共に、算出した電子状態を示す情報を前記記憶手段に記憶させる電子状態算出ステップと、
を含む電子状態計算方法（但し、Ｉは予め定められる１以上の整数であり、ｎ，ｊは時間ステップを示すインデックスである）。
【請求項２】
前記電子状態算出ステップにおいて算出された前記物質の電子状態を示す情報から、予め記憶した物性量算出ルールに基づき当該物質の物性量を算出する物性量算出ステップと、
前記物性量算出ステップにおいて算出された前記物質の物性量を示す情報を出力する出力ステップと、
を更に備える請求項１に記載の電子状態計算方法。
【請求項３】
前記入力ステップにおいて、前記物質の電子状態を計算するためのパラメータとして、前記物質の電子状態の初期状態を示す情報に与えられる外部刺激を示す情報を入力することを特徴とする請求項１又は２に記載の電子状態計算方法。
【請求項４】
前記初期値算出において、前記時間ステップｎにおけるハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）を予め記憶した式
【数２】
0005098480.tif 000035
により演算をして算出する、
ことを特徴とする請求項１〜３の何れか一項に記載の電子状態計算方法。
【請求項５】
前記記憶手段から前記ハートレー項Ｖｈ（ｊ）を含む前記物質の電子状態を示す情報を読み出して、当該ハートレー項Ｖｈ（ｊ）を含む当該物質の電子状態を示す情報から予め記憶した係数値算出ルールに基づき、係数αｊを算出して当該記憶手段に記憶させる係数算出ステップを更に含み、
前記初期値算出ステップにおいて、係数αｊを前記記憶手段から読み出して、前記初期値Ｖｈｉｎｉｔ（ｎ）の算出に用いるものとして設定する、
ことを特徴とする請求項１〜４の何れか一項に記載の電子状態計算方法。
【請求項６】
コンピュータに、ＴＤＤＦＴ法によって物質の電子状態を計算させる電子状態計算プログラムであって、
前記ＴＤＤＦＴ法における時間ステップｊ毎の、前記物質の電子状態のハミルトニアンＨ（ｊ）におけるハートレー項Ｖｈ（ｊ）を含む当該物質の電子状態を示す情報を記憶する記憶手段と、
前記物質の電子状態の初期状態を示す情報を含む前記ＴＤＤＦＴ法によって前記物質の電子状態を計算するためのパラメータを入力して、前記記憶手段に記憶させる入力手段と、
時間ステップｎにおいて、時間ステップｎ−１から時間ステップｎ−Ｉ−１までの前記ハートレー項Ｖｈ（ｎ−１）〜Ｖｈ（ｎ−Ｉ−１）を前記記憶手段から読み出して、予め設定されるＩ個の係数α１〜αＩを用いて、前記時間ステップｎにおけるハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）を予め記憶した式
【数３】
0005098480.tif 000036
により演算をして算出する初期値算出手段と、
前記初期値算出手段により算出された前記ハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）から、予め記憶したハートレー項算出ルールに基づき差分法により当該ハートレー項Ｖｈ（ｎ）を算出すると共に、算出したハートレー項Ｖｈ（ｎ）を前記記憶手段に記憶させるハートレー項算出手段と、
前記ハートレー項算出手段により算出された前記物質の前記ハートレー項Ｖｈ（ｎ）及び前記記憶手段により記憶された情報から、予め記憶した電子状態算出ルールに基づき前記物質の電子状態を算出すると共に、算出した電子状態を示す情報を前記記憶手段に記憶させる電子状態算出手段と、
を前記コンピュータに機能させるための電子状態計算プログラム（但し、Ｉは予め定められる１以上の整数であり、ｎ，ｊは時間ステップを示すインデックスである）。
【請求項７】
ＴＤＤＦＴ法によって物質の電子状態を計算する電子状態計算装置であって、
前記ＴＤＤＦＴ法における時間ステップｊ毎の、前記物質の電子状態のハミルトニアンＨ（ｊ）におけるハートレー項Ｖｈ（ｊ）を含む当該物質の電子状態を示す情報を記憶する記憶手段と、
前記物質の電子状態の初期状態を示す情報を含む前記ＴＤＤＦＴ法によって前記物質の電子状態を計算するためのパラメータを入力して、前記記憶手段に記憶させる入力手段と、
時間ステップｎにおいて、時間ステップｎ−１から時間ステップｎ−Ｉ−１までの前記ハートレー項Ｖｈ（ｎ−１）〜Ｖｈ（ｎ−Ｉ−１）を前記記憶手段から読み出して、予め設定されるＩ個の係数α１〜αＩを用いて、前記時間ステップｎにおけるハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）を予め記憶した式
【数４】
0005098480.tif 000037
により演算をして算出する初期値算出手段と、
前記初期値算出手段により算出された前記ハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）から、予め記憶したハートレー項算出ルールに基づき差分法により当該ハートレー項Ｖｈ（ｎ）を算出すると共に、算出したハートレー項Ｖｈ（ｎ）を前記記憶手段に記憶させるハートレー項算出手段と、
前記ハートレー項算出手段により算出された前記物質の前記ハートレー項Ｖｈ（ｎ）及び前記記憶手段により記憶された情報から、予め記憶した電子状態算出ルールに基づき前記物質の電子状態を算出すると共に、算出した電子状態を示す情報を前記記憶手段に記憶させる電子状態算出手段と、
を備える電子状態計算装置（但し、Ｉは予め定められる１以上の整数であり、ｎ，ｊは時間ステップを示すインデックスである）。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、ＴＤＤＦＴ法によって物質の電子状態を計算する電子状態計算方法、電子状態計算プログラム及び電子状態計算装置に関する。
【背景技術】
【０００２】
電子状態によって定まる時間依存の物質の物性量を算出する方法として、その物性量の要因となっている電子状態を記述する基本方程式を用い、外部刺激に対する応答の時間発展を求める方法が広く用いられている。この種の算出方法では、シミュレーションを行なう空間内に外部刺激を与え、その電子状態の揺らぎを記述する方法として、例えば、時間依存密度汎関数法（以下、ＴＤＤＦＴ（Time-Dependent Density Functional Theory）法と呼ぶ）がある（例えば、非特許文献１参照）。
【０００３】
この計算手法の大きな特徴として、励起状態の情報を提供できる特徴をもっていることが知られている。このため、物質の物性量の中でも励起情報が必要な発光・吸収に関するスペクトルあるいは光学定数の算出には、重要な手法として使われている。
【０００４】
この種の方法では電子状態は一般に波動関数で表され、Kohn-Sham方程式を用いセルフコンシステントな状態を求めておき、これにある時間において微小な摂動として外部刺激を与え、時間依存のKohn-Sham方程式を用い、その応答をその解である波動関数の変化として、時間を追って追跡する。一般に、これらの時間変化を記述する方法として、空間的、時間的に差分化して計算する方法が採られている。
【非特許文献１】K. Yabana and G.F. Bertsch,“Time-dependent local-density approximation in real time”, Phys. Rev B54,4484(1996)
【発明の開示】
【発明が解決しようとする課題】
【０００５】
この時間発展の各ステップで波動関数を求めるために、その時間ステップより以前のステップでの電子状態を初期値として用い安定な状態となったところで、次の時間ステップへ進めることにより、時間発展計算が行なわれる。この安定な状態とするために、特にハートレー項などのハミルトニアンの各項を、数値的に繰り返し計算することが行われている。計算対象となる物質によっては、この繰り返し計算の回数が多くなるため、計算時間が非常に長くなるという問題があった。繰り返し計算の回数は、物質の分子形状に大きく影響され、特に高分子では鎖長が長くなればなるほど、繰り返し回数が増えていた。
【０００６】
本発明は、以上の問題点を解決するためになされたものであり、ＴＤＤＦＴ法による物質の電子状態の計算を効率的に行うことができる電子状態計算方法、電子状態計算プログラム及び電子状態計算装置を提供することを目的とする。
【課題を解決するための手段】
【０００７】
上記の目的を達成するために、本発明に係る電子状態計算方法は、コンピュータにより、ＴＤＤＦＴ法によって物質の電子状態を計算する電子状態計算方法であって、物質の電子状態の初期状態を示す情報を含むＴＤＤＦＴ法によって物質の電子状態を計算するためのパラメータを入力して、ＴＤＤＦＴ法における時間ステップｊ毎の、物質の電子状態のハミルトニアンＨ（ｊ）におけるハートレー項Ｖｈ（ｊ）を含む当該物質の電子状態を示す情報を記憶する記憶手段に記憶させる入力ステップと、時間ステップｎにおいて、時間ステップｎ−１から時間ステップｎ−Ｉ−１までのハートレー項Ｖｈ（ｎ−１）〜Ｖｈ（ｎ−Ｉ−１）を記憶手段から読み出して、予め設定されるＩ個の係数α１〜αＩを用いて、時間ステップｎにおけるハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）を予め記憶した式
【数１】
0005098480.tif 000002
により演算をして算出する初期値算出ステップと、初期値算出ステップにおいて算出されたハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）から、予め記憶したハートレー項算出ルールに基づき差分法により当該ハートレー項Ｖｈ（ｎ）を算出すると共に、算出したハートレー項Ｖｈ（ｎ）を記憶手段に記憶させるハートレー項算出ステップと、ハートレー項算出ステップにおいて算出された物質のハートレー項Ｖｈ（ｎ）及び記憶手段により記憶された情報から、予め記憶した電子状態算出ルールに基づき物質の電子状態を算出すると共に、算出した電子状態を示す情報を記憶手段に記憶させる電子状態算出ステップと、を含むことを特徴とする（但し、Ｉは予め定められる１以上の整数であり、ｎ，ｊは時間ステップを示すインデックスである）。
【０００８】
本発明に係る電子状態計算方法では、時間ステップｎにおけるハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）が、当該時間ステップの１つ前の時間ステップｎ−１から時間ステップｎ−Ｉ−１までのハートレー項Ｖｈ（ｎ−１）〜Ｖｈ（ｎ−Ｉ−１）に基づいて算出される。このように算出された初期値Ｖｈｉｎｉｔ（ｎ）は算出するハートレー項Ｖｈ（ｎ）の値に近くなるため、ハートレー項Ｖｈ（ｎ）の算出を従来と比べてより少ない繰り返し計算で行うことが可能になる。これにより、本発明に係る電子状態計算方法によれば、ＴＤＤＦＴ法による物質の電子状態の計算を効率的に行うことができる。
【０００９】
電子状態計算方法は、電子状態算出ステップにおいて算出された物質の電子状態を示す情報から、予め記憶した物性量算出ルールに基づき当該物質の物性量を算出する物性量算出ステップと、物性量算出ステップにおいて算出された物質の物性量を示す情報を出力する出力ステップと、を更に備えることが望ましい。この構成によれば、効率的に物質の物性量を算出することができる。
【００１０】
入力ステップにおいて、物質の電子状態を計算するためのパラメータとして、物質の電子状態の初期状態を示す情報に与えられる外部刺激を示す情報を入力することが望ましい。この構成によれば、時間変化に応じた物質の電子状態あるいは物性量を適切に算出することができる。
【００１１】
初期値算出において、時間ステップｎにおけるハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）を予め記憶した式
【数２】
0005098480.tif 000003
により演算をして算出する、ことが望ましい。この構成によれば、簡易かつ確実に物質の電子状態の計算を効率的に行うことができる。
【００１２】
電子状態計算方法は、記憶手段からハートレー項Ｖｈ（ｊ）を含む物質の電子状態を示す情報を読み出して、当該ハートレー項Ｖｈ（ｊ）を含む当該物質の電子状態を示す情報から予め記憶した係数値算出ルールに基づき、係数αｊを算出して当該記憶手段に記憶させる係数算出ステップを更に含み、初期値算出ステップにおいて、係数αｊを記憶手段から読み出して、初期値Ｖｈｉｎｉｔ（ｎ）の算出に用いるものとして設定する、ことが望ましい。この構成によれば、ハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）を、より算出するハートレー項Ｖｈ（ｎ）の値に近づけることができ、ハートレー項Ｖｈ（ｎ）の算出を従来と比べて更に少ない繰り返し計算で行うことが可能になる。
【００１３】
ところで、本発明は、上記のように電子状態計算方法の発明として記述できる他に、以下のように電子状態計算プログラム及び電子状態計算装置の発明としても記述することができる。これはカテゴリが異なるだけで、実質的に同一の発明であり、同様の作用及び効果を奏する。
【００１４】
即ち、本発明に係る電子状態計算プログラムは、コンピュータに、ＴＤＤＦＴ法によって物質の電子状態を計算させる電子状態計算プログラムであって、ＴＤＤＦＴ法における時間ステップｊ毎の、物質の電子状態のハミルトニアンＨ（ｊ）におけるハートレー項Ｖｈ（ｊ）を含む当該物質の電子状態を示す情報を記憶する記憶手段と、物質の電子状態の初期状態を示す情報を含むＴＤＤＦＴ法によって物質の電子状態を計算するためのパラメータを入力して、記憶手段に記憶させる入力手段と、時間ステップｎにおいて、時間ステップｎ−１から時間ステップｎ−Ｉ−１までのハートレー項Ｖｈ（ｎ−１）〜Ｖｈ（ｎ−Ｉ−１）を記憶手段から読み出して、予め設定されるＩ個の係数α１〜αＩを用いて、時間ステップｎにおけるハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）を予め記憶した式
【数３】
0005098480.tif 000004
により演算をして算出する初期値算出手段と、初期値算出手段により算出されたハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）から、予め記憶したハートレー項算出ルールに基づき差分法により当該ハートレー項Ｖｈ（ｎ）を算出すると共に、算出したハートレー項Ｖｈ（ｎ）を記憶手段に記憶させるハートレー項算出手段と、ハートレー項算出手段により算出された物質のハートレー項Ｖｈ（ｎ）及び記憶手段により記憶された情報から、予め記憶した電子状態算出ルールに基づき物質の電子状態を算出すると共に、算出した電子状態を示す情報を記憶手段に記憶させる電子状態算出手段と、をコンピュータに機能させることを特徴とする。
【００１５】
本発明に係る電子状態計算装置は、ＴＤＤＦＴ法によって物質の電子状態を計算する電子状態計算装置であって、ＴＤＤＦＴ法における時間ステップｊ毎の、物質の電子状態のハミルトニアンＨ（ｊ）におけるハートレー項Ｖｈ（ｊ）を含む当該物質の電子状態を示す情報を記憶する記憶手段と、物質の電子状態の初期状態を示す情報を含むＴＤＤＦＴ法によって物質の電子状態を計算するためのパラメータを入力して、記憶手段に記憶させる入力手段と、時間ステップｎにおいて、時間ステップｎ−１から時間ステップｎ−Ｉ−１までのハートレー項Ｖｈ（ｎ−１）〜Ｖｈ（ｎ−Ｉ−１）を記憶手段から読み出して、予め設定されるＩ個の係数α１〜αＩを用いて、時間ステップｎにおけるハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）を予め記憶した式
【数４】
0005098480.tif 000005
により演算をして算出する初期値算出手段と、初期値算出手段により算出されたハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）から、予め記憶したハートレー項算出ルールに基づき差分法により当該ハートレー項Ｖｈ（ｎ）を算出すると共に、算出したハートレー項Ｖｈ（ｎ）を記憶手段に記憶させるハートレー項算出手段と、ハートレー項算出手段により算出された物質のハートレー項Ｖｈ（ｎ）及び記憶手段により記憶された情報から、予め記憶した電子状態算出ルールに基づき物質の電子状態を算出すると共に、算出した電子状態を示す情報を記憶手段に記憶させる電子状態算出手段と、を備えることを特徴とする。
【発明の効果】
【００１６】
本発明によれば、ハートレー項Ｖｈ（ｎ）の算出を従来と比べてより少ない繰り返し計算で行うことが可能になるため、ＴＤＤＦＴ法による物質の電子状態の計算を効率的に行うことができる。
【発明を実施するための最良の形態】
【００１７】
以下、図面と共に本発明に係る電子状態計算方法、電子状態計算プログラム及び電子状態計算装置の好適な実施形態について詳細に説明する。なお、図面の説明においては同一要素には同一符号を付し、重複する説明を省略する。
【００１８】
本実施形態に係る電子状態計算方法は、ＴＤＤＦＴ法（ＴＤＤＦＴシミュレーション）によって物質の電子状態を計算するものである。即ち、電子状態計算方法は、時間ステップを設定して、当該時間ステップ毎の物質の電子状態（電子状態の時間発展）を求めるものである。計算された物質の電子状態は、後述するように当該物質の物理的あるいは化学的物性量の算出に用いられる。
【００１９】
本実施形態において電子状態を示す波動関数ψ（ｎ）は、以下の式のように算出対象の時間ステップｎ（ｎ番目の時間ステップ）の前の時間ステップｎ−１の波動関数ψ（ｎ−１）に対する時間発展演算子Ｔ（ｎ−１）により算出される（なお、上記では、時間発展演算子を便宜上Ｔと記載しているが、実際には以下の式に記載されるようにハット記号（＾）が付けられたＴである。以降についても同様）。
【数５】
0005098480.tif 000006
なお、ｎは整数であり時間ステップを示すインデックスである。ここで時間発展演算子Ｔ（ｎ）は、以下の式により求められる。
【数６】
0005098480.tif 000007
ここで、ｉは虚数を、Ｈ（ｎ）は時間ステップｎのハミルトニアンを、Δｔは時間ステップの間隔をそれぞれ示している。上記の時間発展は、Δｔに関して４次以上の展開を行なっても良いし、Suzuki-Trotter形式で演算を行なっても良い（参考文献：Masuo Suzuki, “General NonsymmetricHigher-Order Decomposition of Exponential Operators and Symplectic Integrator”、J．Phys.Soc. Jpn. 61, p3015(1992)）。
【００２０】
時間ステップｎにおけるハミルトニアンＨ（ｎ）は、以下の式のように時間ステップｎにおける電子状態（例えば波動関数ψ（ｎ））に依存している。
【数７】
0005098480.tif 000008
ここで、ψｉ（ｎ）は、時間ステップｎの物質におけるｉ番目の電子状態を、｛ψｉ（ｎ）｝は、その扱う電子状態一式をそれぞれ示している。例えば、擬ポテンシャルを使った密度汎関数法では、ハミルトニアンＨ（ｎ）は、
【数８】
0005098480.tif 000009
と表される。第１項は電子の運動エネルギー項であり、第２項はイオンと電子の相互作用ポテンシャルである。また、第３項がハートレー項であり、第４項は交換相関相互作用を表す。この場合は、時間ステップｎのハミルトニアンＨ（ｎ）は時間ステップｎの電子密度（波動関数から算出可能）に依存している。
【００２１】
時間発展演算の中で一番変化の大きいハートレー項Ｖｈ（ｎ）は、実空間解析では次の式
【数９】
0005098480.tif 000010
を空間差分により差分法を用いて連立方程式に変換し、共役勾配法等の繰り返し法を用いて求める。従来の方法では、この際、初期値Ｖｈｉｎｉｔ（ｎ）として、以下の式
【数１０】
0005098480.tif 000011
に示すように、１つ前の時間ステップｎ−１のハートレー項の値が用いられて、上述の繰り返し計算を行いＶｈ（ｎ）が求められていた。
【００２２】
引き続いて、本実施形態に係る電子状態計算方法が実行される電子状態計算装置１０を説明する。図１に電子状態計算装置１０を示す。電子状態計算装置１０は、具体的には、ワークステーションやＰＣ（Personal Computer）等のコンピュータである。電子状態計算装置１０は、例えばＣＰＵ（Central ProcessingUnit）やメモリ等のハードウェアにより構成されており、これらの構成要素が動作することにより後述する電子状態計算装置１０としての機能が発揮される。なお、本実施形態に係る電子状態計算方法をコンピュータに実行させるプログラムが電子状態計算装置１０において実行されることにより本方法が行われてもよい。
【００２３】
図１に示すように電子状態計算装置１０は、記憶部１１と、入力部１２と、初期値算出部１３と、ハートレー項算出部１４と、電子状態算出部１５と、物性量算出部１６と、出力部１７と、係数算出部１８と、計算制御部１９とを備えて構成される。また、電子状態計算装置１０は、外部装置２０と接続されており、外部装置２０から情報が入力される。
【００２４】
記憶部１１は、物質の電子状態を計算するための情報を記憶する記憶手段である。記憶部１１は、例えば、電子状態計算装置１０が備えるメモリにより実現される。具体的には、物質を示す情報である原子配置やＴＤＤＦＴ法の計算に用いるためのパラメータ（例えば、時間発展のステップ数、ステップ間隔、汎関数の選択、計算空間の大きさ（差分格子間隔と格子点数））、及び物質の電子状態を示す情報等である。
【００２５】
物質の電子状態を示す情報としては、例えば、波動関数及び電子密度等である。記憶部１１に記憶される物質の電子状態を示す情報として、本実施形態に係る電子状態計算方法による算出の初期状態と、本実施形態に係る電子状態計算方法により算出される時間ステップ毎の電子状態を示す情報とを含む。時間ステップｊ毎の物質の電子状態を示す情報には、当該物質の電子状態のハミルトニアンＨ（ｊ）におけるハートレー項Ｖｈ（ｊ）を含む。特に、記憶部１１は、ＴＤＤＦＴ法における時間ステップｊ毎の、物質の電子状態のハミルトニアンＨ（ｊ）におけるハートレー項Ｖｈ（ｊ）を含む当該物質の電子状態を示す情報を記憶する。ここで、ｊは１以上の整数であり、時間ステップを示すインデックスである。なお、記憶部１１は、電子状態を示す情報を格納する際には、何れの時間ステップにおける電子状態であるかが識別できるように、時間ステップを示す情報に対応付けておく。
【００２６】
入力部１２は、物質の電子状態の初期状態を示す情報を含むＴＤＤＦＴ法によって物質の電子状態を計算するためのパラメータを入力して、記憶部１１に記憶させる入力手段である。入力部１２に入力される物質の電子状態の初期状態としては、物質の電子状態を記述する基本方程式を満たすセルフコンシステントなもの（波動関数）とするのが望ましい。
【００２７】
また、入力部１２は、物質の電子状態を計算するためのパラメータとして、物質の電子状態の初期状態を示す情報に与えられる外部刺激を示す情報を入力することとしてもよい。外部刺激を示す情報が入力されると、外部刺激が与えられた物質の電子状態が計算される。ここで外部刺激は、計算対象の物質を構成する原子あるいは分子の位置を時間変化させるものである。入力される外部刺激は、摂動の範囲の刺激であり、物質を構成する原子あるいは分子の位置の時間変化が電子の運動に比べて十分に小さいと見なせるものとするのがよい。セルフコンシステントな初期状態に対して、上記の外部刺激を与えると、物質の時間変化に応じた電子状態、及び物性量が適切に計算できるからである。
【００２８】
外部刺激は、セルフコンシステントな電子状態を示す波動関数ψに対して、時間ステップｊ＝０にて、
【数１１】
0005098480.tif 000012
とすることで与えられる。ここで、ψ´は外部刺激が与えられた波動関数を、Ｅは十分小さい電場を、ｒは刺激を与える方向をそれぞれ示している。なお、参考文献K. Yabana and G.F. Bertsch, “Time-dependent local-densityapproximation in real time”, Phys. Rev B54, 4484(1996)では、解析系の大きさの逆数を指標として、これよりも十分小さな波数とすることが必要であると記述されている。
【数１２】
0005098480.tif 000013
これは全ての電子に均一な速度場を与えることに等しく、これによって生じる電場を扱うといっても良いので上記標識を用いている。
【００２９】
より具体的には、原子単位系を用いて、ｒはｘ方向、ｙ方向、ｚ方向の何れかの座標である。ここで、原子単位系はａｕで表され、ｈ／２π＝ｍ＝ｅ＝１（但し、ここではｈはプランク定数を表す）で求められる単位系である。ボーア半径を１とするので、１ａｕ=0.529177Å、エネルギーは１ハートレーを１とするため１ａｕ＝27.12ｅＶとなる。表式が簡便となるため量子力学計算では良く用いる単位系である。
【００３０】
上記の情報の入力は、例えば、ユーザによる操作によって外部装置２０から行われる。また、所定の演算を行って上記の情報（例えば、セルフコンシステントな電子状態を表す情報等）が生成されてもよい。入力部１２は、入力した各情報を記憶部１１に記憶させる。また、外部装置２０から入力されたパラメータを予め電子状態計算装置１０の記憶部１１に記憶させておき、電子状態計算方法の実行時に記憶された情報が用いられてもよい。入力部１２は、各情報を入力するとその旨を計算制御部１９に通知する。
【００３１】
初期値算出部１３は、時間ステップｎにおいて、ハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）を算出する初期値算出手段である。また、初期値算出部１３は、算出した初期値Ｖｈｉｎｉｔ（ｎ）を記憶部１１に格納する。具体的には、初期値算出部１３は、時間ステップｎ−１から時間ステップｎ−Ｉ−１までのハートレー項Ｖｈ（ｎ−１）〜Ｖｈ（ｎ−Ｉ−１）を記憶部１１から読み出して、予め設定されるＩ個の係数α１〜αＩを用いて、初期値Ｖｈｉｎｉｔ（ｎ）を予め記憶した式
【数１３】
0005098480.tif 000014
により演算をして算出する。ここで、Ｉは予め定められる１以上の整数である。即ち、初期値算出部１３は、初期値Ｖｈｉｎｉｔ（ｎ）の算出に、１つ前の時間ステップｎ−１でのハートレー項Ｖｈ（ｎ−１）の値だけでなく、２つ以上前の時間ステップでのハートレー項の値を用いる。ハートレー項を用いる時間ステップの数は、Ｉ＋１である。
【００３２】
特に、初期値算出部１３は、初期値Ｖｈｉｎｉｔ（ｎ）の算出に、２つ前までの時間ステップでのハートレー項を用いることとするのが好ましい。即ち、初期値算出部１３は、初期値Ｖｈｉｎｉｔ（ｎ）を予め記憶した式
【数１４】
0005098480.tif 000015
により演算をして算出する（上記の一般化した式（１１）におけるＩ＝１）。なお、上記の式（１２）における係数αは、式（１１）における係数α１に対応する。これは、計算の容易さ、記憶容量、及び確実に本発明を実施する等の観点によるものである。
【００３３】
上記のように、初期値Ｖｈｉｎｉｔ（ｎ）の算出に、２つ前までの過去の時間ステップでのハートレー項の値を用いることとすると、初期値Ｖｈｉｎｉｔ（ｎ）は、算出するハートレー項Ｖｈ（ｎ）の値に近くなる。これは、図２に示すように、２つ前までの過去の時間ステップでのハートレー項の値を用いることによって、ハートレー項Ｖｈ（ｎ）の変動割合（図２においてはα［Ｖｈ（ｎ−１）−Ｖｈ（ｎ−２）］が相当する）が考慮された初期値Ｖｈｉｎｉｔ（ｎ）となるためである。
【００３４】
また、各［Ｖｈ（ｎ−ｊ）−Ｖｈ（ｎ−ｊ−１）］の係数αｊは、予め設定されておいてもよいし、後述するように係数算出部１８により算出されてもよい。係数αｊは記憶部１１に記憶されており、初期値算出部１３は、係数αｊを記憶部１１から読み出して、初期値Ｖｈｉｎｉｔ（ｎ）の算出に用いるものとして設定する。係数αｊは、計算に用いるハートレー項Ｖｈ（ｎ）の数に応じて、適切に設定されることが望ましい。上記のように２つ前までの過去の時間ステップでのハートレー項の値を用いる場合には、例えば、係数α１を１とするのが望ましい。また、計算機実験等によって、経験的に適切な値を設定するのが望ましい。
【００３５】
ハートレー項算出部１４は、初期値算出部１３により算出されたハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）から、予め記憶したハートレー項算出ルールに基づき差分法により当該ハートレー項Ｖｈ（ｎ）を算出するハートレー項算出手段である。また、ハートレー項算出部１４は、算出したハートレー項Ｖｈ（ｎ）を記憶部１１に記憶させる。ハートレー項算出部１４は、上記の算出のための演算に応じたパラメータ等を記憶部１１から読み出して、演算に用いる。
具体的には、ハートレー項算出部１４は、記憶部１１から算出するハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）を読み出して、上述したように、ハートレー項Ｖｈ（ｎ）を、次の式
【数１５】
0005098480.tif 000016
を空間差分により差分法を用いて連立方程式に変換し、共役勾配法等の繰り返し法を用いて求める。ハートレー項算出部１４は、算出したハートレー項Ｖｈ（ｎ）の値を記憶部１１に記憶させる。
【００３６】
電子状態算出部１５は、算出対象の１つ前の時間ステップでの物質の電子状態を示す情報を記憶部１１から読み出して、当該物質の電子状態を示す情報とハートレー項算出部１４によって算出されたハートレー項Ｖｈ（ｎ−１）とから、予め記憶した電子状態算出ルールに基づき物質の電子状態を算出する電子状態算出手段である。また、電子状態算出部１５は、算出した電子状態を示す情報を記憶部１１に記憶させる。また、電子状態算出部１５は、上記の算出のための演算に応じたパラメータ等を記憶部１１から読み出して、演算に用いる。
【００３７】
具体的には、電子状態算出部１５は、算出対象の１つ前の時間ステップｎ−１において、電子状態を示す波動関数ψ（ｎ−１）からハミルトニアンＨ（ｎ−１）を算出する。続いて、電子状態算出部１５は、ハミルトニアンＨ（ｎ−１）から時間発展演算子Ｔ（ｎ−１）を算出する。続いて、電子状態算出部１５は、時間ステップｎにおける波動関数ψ（ｎ）を、以下の式により算出する。
【数１６】
0005098480.tif 000017
【００３８】
また、電子状態算出部１５は、電子状態として波動関数ψ（ｎ）だけでなく、例えば、以下の式により電子密度ρ（ｎ）を算出することとしてもよい。
【数１７】
0005098480.tif 000018
【００３９】
ここで、ｗｉは、重み係数を示す値であり、予め記憶部１１に記憶されており、電子密度ρ（ｎ）の算出の際に記憶部１１から読み出されて算出に用いられる。
【００４０】
物性量算出部１６は、電子状態算出部１５により算出された物質の電子状態を示す情報から、予め記憶した物性量算出ルールに基づき当該物質の物性量を算出する物性量算出手段である。物性量算出部１６は、算出した物性量を示す情報を記憶部１１に記憶させる。
【００４１】
具体的には例えば、物性量算出部１６は、物質の物性量として、発光・吸収スペクトルを算出するための双極子モーメントμ（ｎ）を、物性量算出部１６に算出された電子密度ρ（ｎ）を記憶部１１から読み出して、以下の式
【数１８】
0005098480.tif 000019
により算出する。ここで、ｒは、空間上の位置を示すベクトルである。物性量算出部１６は、双極子モーメントμ（ｎ）を各時間ステップ毎に算出する。
【００４２】
また、物性量算出部１６は、物質の物性量として、上記以外にも、ベクトルポテンシャルＡ（ｎ）を用いて以下の式により表面電荷σ（ｎ）を算出することとしてもよい。
【数１９】
0005098480.tif 000020
【００４３】
出力部１７は、物性量算出部１６によって算出された物質の物性量を示す情報を、記憶部１１から読み出して、出力する出力手段である。出力部１７は、例えば、電子状態計算装置１０が備えるディスプレイ等の表示装置に物性量の情報を表示することによって出力を行う。また、電子状態計算装置１０から他の装置に出力することとしてもよい。また、出力部１７は、電子状態算出部１５によって算出された物質の電子状態を示す情報を、記憶部１１から読み出して、出力することとしてもよい。
【００４４】
係数算出部１８は、初期値Ｖｈｉｎｉｔ（ｎ）の算出に用いるハートレー項を含む物質の電子状態を示す情報を記憶部１１から読み出して、当該ハートレー項を含む物質の電子状態を示す情報から予め記憶した係数値算出ルールに基づき、Ｉ個の係数αｊを算出して記憶部１１に記憶させる係数算出手段である。ここで算出されるαｊは適切なαｊ、即ち、算出するハートレー項Ｖｈ（ｎ）により近い初期値Ｖｈｉｎｉｔ（ｎ）を算出させるものである。係数算出部１８は、算出した係数αｊを記憶部１１に記憶させる。
【００４５】
具体的には、係数算出部１８は、以下のような最適化解を求める処理により、係数αｊを求める。以下では、初期値Ｖｈｉｎｉｔ（ｎ）の算出に、２つ前までの時間ステップでのハートレー項を用いる場合の例を示す。まず、空間をｍ個に差分化し表現すると、ハートレー項Ｖｈ（ｎ）はｍ次元のベクトルで表される。差分化に伴う定数をｍ×ｍの正定値対称行列Ａとし、電子密度ρ（ｎ）に関連した既知量をｍ次元のベクトルｂ（ｎ）とすると、式（７）は、
【数２０】
0005098480.tif 000021
という連立一次方程式の形に表すことができる。このとき、共役勾配法と同様に重みつき誤差（Ｖｈｉｎｉｔ（ｎ）−Ｖｈ（ｎ），Ａ（Ｖｈｉｎｉｔ（ｎ）−Ｖｈ（ｎ）））を最小とする最適な係数αｊは、関数
【数２１】
0005098480.tif 000022
の最小化を
【数２２】
0005098480.tif 000023
により計算することで、求めることができる。
【００４６】
具体的には、表式を簡単にするため時間ステップｊにおけるハートレー項Ｖｈ（ｊ）は、ｘｊ、また、算出対象の初期値Ｖｈｉｎｉｔ（ｎ）をｘｎ＊と書き表すと、２つ前までの時間ステップを用いる場合ｘｎ＊＝ｘｎ−１＋α（ｘｎ−１−ｘｎ−２）（ここで、α＝α１である）であるので、この関数ｆ（ｘｎ＊）は、
【数２３】
0005098480.tif 000024
と書き表せる。この関数ｆ（ｘｎ＊）を最小化する係数αは式（１８）から、
【数２４】
0005098480.tif 000025
と求めることができる。ここで、ハートレー項算出部１４により時間ステップｎより前の時間ステップにおいて算出されたハートレー項Ｖｎ（ｊ）は、連立一次方程式（１６）の近似解であってもよいとし、以下の関係式
【数２５】
0005098480.tif 000026
が成り立つものとした。
【００４７】
また、３つ前までの時間ステップを用いる場合その係数をα（＝α１）及びβ（＝α２）とすると、ｘｎ＊＝ｘｎ−１＋α（ｘｎ−１−ｘｎ−２）＋β（ｘｎ−２−ｘｎ−３）であるので、この関数ｆ（ｘｎ＊）は、
【数２６】
0005098480.tif 000027
となる。係数α及びβについて、この関数ｆ（ｘｎ＊）の停留点は、
【数２７】
0005098480.tif 000028
及び
【数２８】
0005098480.tif 000029
から、
【数２９】
0005098480.tif 000030
と求めることができる。また、上記算出された係数α及びβにおいて、関数ｆ（ｘｎ＊）が最小値を取る条件は、
【数３０】
0005098480.tif 000031
を満たす場合である。
【００４８】
更に、最適な係数αｊを求める別の方法として、残差ベクトルｂ（ｎ）−ＡＶｈｉｎｉｔ（ｎ）の大きさを最小化しても良い。即ち、関数
【数３１】
0005098480.tif 000032
を用い、
【数３２】
0005098480.tif 000033
の計算から、残差を最小とする係数αｊを求めることも可能である。
【００４９】
計算制御部１９は、本実施形態に係る電子状態計算方法のＴＤＤＦＴ法による処理を全体的に統括する手段である。具体的には、上記の各処理を行う構成要素１２〜１９からの通知を受けて、各処理を行う構成要素１２〜１９に対して次に行うべき処理を指示する。計算制御部１９は、予め設定されたステップ数の計算が完了したか否かを判断して、計算を終了させる等を行う。ここで、計算制御部１９は、設定されたステップ数の情報を記憶部１１から読み出して、上記の判断に用いる。
【００５０】
上記の各処理を行う構成要素１２〜１９は、例えば、電子状態計算装置１０が備えるＣＰＵにより実現される。以上が電子状態計算装置１０の構成である。
【００５１】
引き続いて、図３のフローチャートを用いて、電子状態計算装置１０により実行される本実施形態に係る電子状態計算方法を説明する。
【００５２】
電子状態計算装置１０では、まず、入力部１２によって、物質の電子状態の算出に必要な、上述したパラメータ等の情報が入力される（Ｓ０１、入力ステップ）。入力された各情報は、記憶部１１に記憶される。入力される情報には、物質の電子状態を示す波動関数の初期値ψ（０）及び外部刺激が含まれる。ここで、ＴＤＤＦＴ法における時間ステップｎは０（初期値）とされる。
【００５３】
続いて、時間ステップｎは、それまでの時間ステップｎに１が加えられたものとされる。続いて、電子状態算出部１５によって、時間ステップｎにおける波動関数ψ（ｎ）が算出される（Ｓ０２、電子状態算出ステップ）。具体的には、前の時間ステップｎ−１において算出されたハミルトニアンＨ（ｎ−１）から時間発展演算子Ｔ（ｎ−１）が算出される。続いて、前の時間ステップｎ−１において算出された波動関数ψ（ｎ−１）と時間発展演算子Ｔ（ｎ−１）とが上述したように演算されることにより、時間ステップｎにおける波動関数ψ（ｎ）が算出される。なお、時間ステップｎが１の場合は、時間ステップｎ−１（即ち時間ステップが０の場合）における各値は、算出された値でなく、記憶部１１に初期値として記憶されているものが用いられる。算出された波動関数ψ（ｎ）は、記憶部１１に記憶される。
【００５４】
続いて、初期値算出部１３によって、時間ステップｎにおける、ハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）が算出される（Ｓ０３、初期値算出ステップ）。初期値Ｖｈｉｎｉｔ（ｎ）は、上述したように、時間ステップｎの１つ前の時間ステップｎ−１からＩ＋１前の時間ステップｎ−Ｉ−１までのハートレー項Ｖｈ（ｎ−１）〜Ｖｈ（ｎ−Ｉ−１）から算出される。なお、時間ステップｎがＩ以下である場合は、各ハートレー項Ｖｈ（ｎ−１）〜Ｖｈ（ｎ−Ｉ−１）が取得できないため、上記の方法による初期値Ｖｈｉｎｉｔ（ｎ）の算出は行わない。その場合、従来と同様に、前の時間ステップｎ−１のハートレー項Ｖｈ（ｎ−１）を初期値Ｖｈｉｎｉｔ（ｎ）とすることとしてもよい。算出された初期値Ｖｈｉｎｉｔ（ｎ）は、記憶部１１に記憶される。
【００５５】
ここで、初期値Ｖｈｉｎｉｔ（ｎ）の算出の前に、係数算出部１８によって、初期値Ｖｈｉｎｉｔ（ｎ）の算出に用いられる係数αｊの算出が行われることとしてもよい（Ｓ０４、係数算出ステップ）。係数αｊの算出が行われた場合、係数αｊは記憶部１１に記憶されて、初期値算出部１３に初期値Ｖｈｉｎｉｔ（ｎ）の算出の際に読み出されて用いられる。なお、係数αｊの算出が行われない場合、予め係数αｊが記憶部１１に記憶されている。
【００５６】
続いて、ハートレー項算出部１４によって、初期値算出部１３により算出されたハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）から、予め記憶したハートレー項算出ルールに基づき差分法により当該ハートレー項Ｖｈ（ｎ）が算出される（Ｓ０５、ハートレー項算出ステップ）。算出されたハートレー項Ｖｈ（ｎ）は、記憶部１１に記憶される。
【００５７】
続いて、電子状態算出部１５によって、算出された波動関数ψ（ｎ）からハミルトニアンＨ（ｎ）が算出される（Ｓ０６、電子状態算出ステップ）。算出されたハミルトニアンＨ（ｎ）は、記憶部１１に記憶される。
【００５８】
続いて、電子状態算出部１５によって、算出された波動関数ψ（ｎ）から電子密度ρ（ｎ）が算出される（Ｓ０７、電子状態算出ステップ）。算出された電子密度ρ（ｎ）は、記憶部１１に記憶される。
【００５９】
続いて、物性量算出部１６によって、算出された電子密度ρ（ｎ）から、物質の物性量である双極子モーメントμ（ｎ）が、予め記憶した物性量算出ルールに基づき算出される（Ｓ０８、物性量算出ステップ）。なお、算出される物質の物性量は、必ずしも双極子モーメントである必要は無く、物質の電子状態から算出されるものであればよい。算出された双極子モーメントμ（ｎ）は、記憶部１１に記憶される。
【００６０】
続いて、出力部１７によって、物性量算出部１６によって算出された物質の物性量を示す情報が出力される（Ｓ０９、出力ステップ）。また、併せて、算出された波動関数ψ（ｎ）等の物質の電子状態を示す情報が出力されてもよい。なお、この出力は、各時間ステップｎにおいて行われる必要は必ずしも無く、本処理の最後にまとめて行われてもよい。
【００６１】
続いて、計算制御部１９によって、現在の時間ステップｎが、処理を終了させるために予め設定された時間ステップ数であるか否かが判断される（Ｓ１０）。現在の時間ステップｎが、設定された時間ステップ数であった場合は、ここで処理は終了される。一方、現在の時間ステップｎが、設定された時間ステップ数でなかった場合は、上述したように、それまでの時間ステップｎに１が加えられたものして、即ち、時間ステップを更新して、Ｓ０２〜Ｓ１０の処理を繰り返し行う。以上が、本実施形態に係る電子状態計算方法である。
【００６２】
上述したように、本実施形態では、時間ステップｎにおけるハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）が、当該時間ステップの前の複数の時間ステップにおけるハートレー項に基づいて算出される。このように算出された初期値Ｖｈｉｎｉｔ（ｎ）は算出するハートレー項Ｖｈ（ｎ）の値に近くなるため、ハートレー項Ｖｈ（ｎ）の算出を従来と比べてより少ない繰り返し計算で行うことが可能になる。これにより、本発明に係る電子状態計算方法によれば、ＴＤＤＦＴ法による物質の電子状態の計算を効率的に行うことができる。
【００６３】
また、本実施形態のように、物質の物性値を算出することとすれば、効率的に物質の物性量を算出することができる。
【００６４】
また、本実施形態のように、物質のセルフコンシステントな電子状態に対して外部刺激を与えることにより、当該電子状態の外部刺激に対する応答を計算することとすれば、時間変化に応じた物質の電子状態あるいは物性量を適切に算出することができる。
【００６５】
また、本実施形態のように、ハートレー項Ｖｈ（ｎ）を算出する時間ステップｎの２つ前までのハートレー項を用いて、初期値Ｖｈｉｎｉｔ（ｎ）を算出することとすれば、簡易かつ確実に物質の電子状態の計算を効率的に行うことができる。
【００６６】
また、本実施形態のように初期値Ｖｈｉｎｉｔ（ｎ）を算出するための式に含まれる係数αｊを算出することとすれば、ハートレー項Ｖｈ（ｎ）の初期値Ｖｈｉｎｉｔ（ｎ）を、より算出するハートレー項Ｖｈ（ｎ）の値に近づけることができ、ハートレー項Ｖｈ（ｎ）の算出を従来と比べて更に少ない繰り返し計算で行うことが可能になる。しかしながら、この場合、係数αｊを算出する処理が増えることとなるので、係数αｊの算出によるハートレー項Ｖｈ（ｎ）の算出のための繰り返し計算を減少させる効果が顕著な場合に用いることが望ましい。
【００６７】
引き続いて、上述した一連の物質の電子状態の計算を行う処理をコンピュータに実行させるための電子状態計算プログラムを説明する。図４に示すように、電子状態計算プログラム３１は、コンピュータが備える記録媒体３０に形成されたプログラム格納領域３０ａ内に格納されている。
【００６８】
電子状態計算プログラム３１は、記憶モジュール３１ａと、入力モジュール３１ｂと、初期値算出モジュール３１ｃと、ハートレー項算出モジュール３１ｄと、電子状態算出モジュール３１ｅと、物性量算出モジュール３１ｆと、出力モジュール３１ｇと、係数算出モジュール３１ｈと、計算制御モジュール３１ｉとを備えて構成される。記憶モジュール３１ａ、入力モジュール３１ｂ、初期値算出モジュール３１ｃ、ハートレー項算出モジュール３１ｄ、電子状態算出モジュール３１ｅ、物性量算出モジュール３１ｆ、出力モジュール３１ｇ、係数算出モジュール３１ｈ及び計算制御モジュール３１ｉを実行させることにより実現される機能は、上述した電子状態計算装置１０の記憶部１１、入力部１２、初期値算出部１３、ハートレー項算出部１４、電子状態算出部１５、物性量算出部１６、出力部１７、係数算出部１８及び計算制御部１９の機能とそれぞれ同様である。
【００６９】
なお、電子状態計算プログラム３１は、その一部若しくは全部が、通信回線等の伝送媒体を介して伝送され、他の機器により受信されて記録（インストールを含む）される構成としてもよい。
【実施例】
【００７０】
次に本発明を実施例によって更に詳しく解説するが、本発明は実施例に限定されるものではない。
【００７１】
［実施例１］
フルオレン１分子を本発明における電子状態の算出対象とした。本実施例では、密度汎関数法を用いた実空間実時間ＴＤＤＦＴ法において、初期値Ｖｈｉｎｉｔの算出に、式（１２）を用いた。また、初期値Ｖｈｉｎｉｔとして式（８）を用いて算出した場合を比較例とする。各場合について、1,000ステップまでの各時間ステップにおけるハートレー項Ｖｈの算出に要する繰り返し回数、及び10,000ステップまでの平均繰り返し回数を計測した。
【００７２】
本実施例で用いる入力等として、密度汎関数法としては交換相関相互作用にP.J.Perdew and Alex Zunger, ”Self-interaction correction to density-functionalapproximation for many-electron systems”, Phys. Rev. B23, p5048(1981)（以後、Perdew-Zungerとして引用）により与えられるパラメータを用い、また擬ポテンシャルはN.Troullior and Jose Luis Martins, ”Efficient pseudopotentials for plane-wavecalculations”, Phys. Rev. B43, p1993(1991)（以後、Troullior-Martinsとして引用）により示された手法により予め生成しておいたものを用いた。原子配置は最適化された位置に固定し、空間の差分格子間隔は0.3Å、格子点数は、65×53×47である。これによりセルフコンシステントになるまでKohn-Sham方程式に従って電子状態を予め生成しておいた。これに外部刺激として0.01Ｖ／Åのδ関数型の電場をｔ＝０において印加した。
【００７３】
また、波動関数の時間発展計算は、時間間隔２πΔｔ／ｈ（但し、ここではｈはプランク定数を表す）を0.002ｅＶ−１として、時間発展演算子をΔｔについて４次まで展開することにより行われた。また、上述したように電子状態として電子密度ρ（ｎ）の算出を行った。また、それに基づき、物質の物性量として、物理量である光の吸収発光に関る双極子モーメントμ（ｎ）を算出した。また、これらの計算は時間ステップが10、000ステップとなることを判定条件として計算を終了させることとした。
【００７４】
ハミルトニアンＨ（ｎ）の更新は式（６）により行なわれるが、ここでのハートレー項Ｖｈ（ｎ）を求めるにあたり式（７）を実空間での差分法により連立方程式に置き換え、共役勾配法を用いた。このときの初期値として、本発明における予測手法を用い、最適な初期値Ｖｈｉｎｉｔ（ｎ）を提供することとした。初期値Ｖｈｉｎｉｔ（ｎ）の予測には式（１２）を用い、α＝１とした。
【００７５】
この結果、10,000ステップまでのハートレー項Ｖｈの算出に要する平均繰り返し回数は、比較例では１７回であるのに対し、本発明では５回であった。本発明は、従来の方法に比べて繰り返し回数が非常に小さくて済み、計算時間の短縮に大きく寄与していることが認識される。
【００７６】
また初期1,000ステップまでの各時間ステップにおいてハートレー項Ｖｈの算出に要する繰り返し回数の時間発展の様子を図５に示す。図５において横軸はＴＤＤＦＴ法による時間ステップを示す値であり、縦軸は各時間ステップにおいてハートレー項Ｖｈの算出に要する繰り返し回数である。本発明による方法が、初期の時間ステップではより効果があることが認識される。
【００７７】
［実施例２］
実施例１と同様に、本発明による計算と、比較例による計算とを行い、10,000ステップまでの平均繰り返し回数を計測した。本実施例では、フルオレン１量体、２量体、３量体、４量体、６量体及び８量体に対して、それぞれ計算を行った。
【００７８】
本実施例で用いる入力及び処理内容等は、実施例１とほぼ同様である。交換相関相互作用にPerdew-Zungerにより与えられるパラメータを用い、また擬ポテンシャルはTroullior-Martinsにより示された手法により予め生成しておいたものを用いた。原子配置はいずれも最適化された位置に固定し、空間の差分格子間隔は0.3Åで、格子点数はフルオレン１量体の場合は65×53×47、２量体の場合は91×55×51、３量体の場合は119×55×53、４量体の場合は147×55×55、６量体の場合は203×55×55、８量体の場合は257×55×55である。これによりセルフコンシステントになるまでKohn-Sham方程式に従って電子状態を予め生成しておいた。
【００７９】
また、実施例１と同じ外部刺激を加え、時間間隔２πΔｔ／ｈ（但し、ここではｈはプランク定数を表す）を0.002ｅＶ−１とし、物理量として双極子モーメントμ（ｎ）を算出した。
【００８０】
特に本発明で重要なハミルトニアンＨ（ｎ）の更新は、実施例１と同様に式（６）により行なわれるが、ここでのハートレー項Ｖｈを求めるにあたり式（７）を実空間での差分法により連立方程式に置き換え、共役勾配法を用いた。このときの初期値として、本発明における予測手法を用い、本発明による初期値Ｖｈｉｎｉｔ（ｎ）を提供することとした。初期値Ｖｈｉｎｉｔ（ｎ）の予測には式（１２）を用い、α＝１とした。
【００８１】
この結果、10,000ステップまでのハートレー項Ｖｈの算出に要する平均繰り返し回数と量数との関係について、従来法と本発明とによる予測手法を用いることにより得られた結果を、図６に示す。
【００８２】
従来法では分子形状（例えば高分子鎖長）が長くなるにつれて、計算にかかるコストが非常に大きくなるが、本発明を用いると共役勾配法の平均繰り返し回数は鎖長に依らずほぼ一定の値、この場合５回程度となっていることがわかる。このように、本発明は、分子形状による影響を最小限にとどめる効果があることもわかる。
【図面の簡単な説明】
【００８３】
【図１】本発明の実施形態に係る電子状態計算装置の構成図である。
【図２】本発明の実施形態において算出されるハートレー項の初期値を概念的に示した図である。
【図３】本発明の実施形態に係る電子状態計算方法を示すフローチャートである。
【図４】本発明の実施形態に係る電子状態計算プログラムを示す図である。
【図５】本発明と比較例とによる電子状態の算出の際の、各時間ステップにおけるハートレー項Ｖｈの算出に要する共役勾配法の繰り返し回数を示したグラフである。
【図６】本発明と比較例とによる電子状態の算出の際の、フルオレンの量体数に応じた、ハートレー項Ｖｈの算出に要する共役勾配法の平均繰り返し回数を示したグラフである。
【符号の説明】
【００８４】
１０…電子状態計算装置、１１…記憶部、１２…入力部、１３…初期値算出部、１４…ハートレー項算出部、１５…電子状態算出部、１６…物性量算出部、１７…出力部、１８…係数算出部、１９…計算制御部、２０…外部装置、３０…記録媒体、３０ａ…プログラム格納領域、３１…電子状態計算プログラム、３１ａ…記憶モジュール、３１ｂ…入力モジュール、３１ｃ…初期値算出モジュール、３１ｄ…ハートレー項算出モジュール、３１ｅ…電子状態算出モジュール、３１ｆ…物性量算出モジュール、３１ｇ…出力モジュール、３１ｈ…係数算出モジュール、３１ｉ…計算制御モジュール。
【図１】
0005098480.tif 000038
【図２】
0005098480.tif 000039
【図３】
0005098480.tif 000040
【図４】
0005098480.tif 000041
【図５】
0005098480.tif 000042
【図６】
0005098480.tif 000043