"計算装置を用いて物質の電子状態を計算する方法において、前記計算装置は、物質の電子状態に対する近似解を与える複数の演算モデルの夫々を構成要素として含む集合を設定し、複数の電子軌道上に存在する少なくとも2つ以上の電子を含む電子系に働く有効相互作用を含む有効ハミルトニアンを、各演算モデルを特定するために定め、前記有効ハミルトニアンの自己無撞着解を前記集合内の各演算モデルを用いて計算する過程で、計算した自己無撞着解同士が連続的に変化する方向を特定しながら、前記集合がなす空間内で距離が近接する複数の演算モデルのうちで最適な演算モデルを量子力学的変分法に基づいて定め、前記最適な演算モデルを逐次更新する際に、前記有効ハミルトニアンの自己無撞着解による電子系の変分エネルギを評価し、評価した変分エネルギが計算すべき厳密解のエネルギに近接してゆき、しかも変分エネルギが単調減少凸関数をなすように演算モデルを更新して、一又は複数の変分エネルギの系列から前記電子状態の厳密解を計算することを特徴とする電子状態計算方法。局所密度近似法、一般化勾配近似法、又はハートレーフォック法を基に、有効相互作用として揺らぎを示す非局所型演算子を定める手続きを通じて前記有効ハミルトニアンを決定することを特徴とする請求項1に記載の電子状態計算方法。前記揺らぎは、クーロン相互作用とハートレー平均場項のずれである密度揺らぎを含むことを特徴とする請求項2に記載の電子状態計算方法。前記厳密解に近接してゆく変分エネルギの系列を与える一又は複数の演算モデルの系列のうち、電子状態を示す相が転移せずに厳密解と同一視される相を与え、かつ最小の計算ステップで前記有効ハミルトニアンの自己無撞着解に至る演算モデルを含むものを定めることを特徴とする請求項1に記載の電子状態計算方法。各演算モデルによる自己無撞着解をLAPW法LinearizedAugmentedPlaneWave法、PAW法ProjectorAugmentedWave法、又は数値基底展開法を用いた並列計算により求めることを特徴とする請求項1から請求項4の何れか1つに記載の電子状態計算方法。物質の電子状態を計算する装置において、物質の電子状態に対する近似解を与える複数の演算モデルの夫々を構成要素として含む集合を設定する手段と、複数の電子軌道上に存在する少なくとも2つ以上の電子を含む電子系に働く有効相互作用を含む有効ハミルトニアンを、各演算モデルを特定するために定める手段と、前記有効ハミルトニアンの自己無撞着解を前記集合内の各演算モデルを用いて計算する過程で、計算した自己無撞着解同士が連続的に変化する方向を特定しながら、前記集合がなす空間内で距離が近接する複数の演算モデルのうちで最適な演算モデルを量子力学的変分法に基づいて定める手段と、前記最適な演算モデルを逐次更新する際に、前記有効ハミルトニアンの自己無撞着解による電子系の変分エネルギを評価する手段と、評価した変分エネルギが計算すべき厳密解のエネルギに近接してゆき、しかも変分エネルギが単調減少凸関数をなすように演算モデルを更新する手段と、一又は複数の変分エネルギの系列から前記電子状態の厳密解を計算する手段とを備えることを特徴とする電子状態計算装置。局所密度近似法、一般化勾配近似法、又はハートレーフォック法を基に、有効相互作用として揺らぎを示す非局所型演算子を定める手続きを通じて前記有効ハミルトニアンを決定するようにしてあることを特徴とする請求項6に記載の電子状態計算装置。前記揺らぎは、クーロン相互作用とハートレー平均場項のずれである密度揺らぎを含むことを特徴とする請求項7に記載の電子状態計算装置。前記厳密解に近接してゆく変分エネルギの系列を与える一又は複数の演算モデルの系列のうち、電子状態を示す相が転移せずに厳密解と同一視される相を与え、かつ最小の計算ステップで前記有効ハミルトニアンの自己無撞着解に至る演算モデルを含むものを定めることを特徴とする請求項6に記載の電子状態計算装置。各演算モデルによる自己無撞着解をLAPW法LinearizedAugmentedPlaneWave法、PAW法ProjectorAugmentedWave法、又は数値基底展開法を用いた並列計算により求めるようにしてあることを特徴とする請求項6から請求項9の何れか1つに記載の電子状態計算装置。コンピュータに、物質の電子状態を計算させるコンピュータプログラムにおいて、物質の電子状態に対する近似解を与える複数の演算モデルの夫々を構成要素として含む集合を設定させるステップと、複数の電子軌道上に存在する少なくとも2つ以上の電子を含む電子系に働く有効相互作用を含む有効ハミルトニアンを、各演算モデルを特定するために定めさせるステップと、前記有効ハミルトニアンの自己無撞着解を前記集合内の各演算モデルを用いて計算させる過程で、計算させた自己無撞着解同士が連続的に変化する方向を特定しながら、前記集合がなす空間内で距離が近接する複数の演算モデルのうちで最適な演算モデルを量子力学的変分法に基づいて定めさせるステップと、前記最適な演算モデルを逐次更新させる際に、前記有効ハミルトニアンの自己無撞着解による電子系の変分エネルギを評価させるステップと、評価させた変分エネルギが計算すべき厳密解のエネルギに近接してゆき、しかも変分エネルギが単調減少凸関数をなすように演算モデルを更新させ、一又は複数の変分エネルギの系列から前記電子状態の厳密解を計算させるステップとをコンピュータに実行させることを特徴とするコンピュータプログラム。前述した近似計算手法は、いずれも真の解厳密解が示す物性に有限回数の計算で常に到達する方法を与えていないという問題点を有していた。すなわち、LDAが与える高精度の一電子基底を活用しながら、常に真の解の物性を再現することが保障されている自己無撞着計算理論は、存在していなかった。一方、本願発明者の一名は、多配置参照密度汎関数法に基づく有効多電子計算として揺らぎ参照決定法を提案している例えば、特許文献1を参照。この方法により、一電子密度に加えて、局所密度揺らぎなどの正定値性を備えた揺らぎ変数相関関数を再現する拡張コーンシャム方程式を与えて、物理量の再現精度を任意の精度で高めることが可能である。この計算手法では、ハミルトニアンの最低エネルギ状態の決定を通して、物質の全エネルギ、一電子密度、特定した揺らぎ変数を同時に再現した。更に、拡張コーンシャム方程式が定めるモデルとその安定解を初期値として用いることで、全エネルギ、一電子密度、正準相関関数を再現する時間発展方程式を与えて物質の外場応答力学変形、電磁場印加に対する応答などを再現するための、高精度な第一原理電子状態計算の方法を与えた特許文献1を参照。しかしながら、前記の多配置参照密度汎関数法に基づく有効多電子計算は、量子モンテカルロ法、トランスコリレイティッド法、配置間相互作用法、摂動計算グリーン関数法、有効ポテンシャル法などの参照計算が不可避な計算手法であった。そのため、物理量の再現精度は参照計算の計算精度に依存することとなり、参照計算の計算精度を越えて真の解に到達する計算手法を与えるものではなかった。そこでさらに、本願発明者の一名は、密度汎関数を用いる量子力学的変分計算理論密度汎関数変分法及び厳密解を含むバナッハ空間となるモデル空間中で厳密解に到達する近似モデルの列を与える方法を与え、計算機資源の許す範囲内で安定した数値解を求めることで、厳密解の評価ができる電子状態計算方法、電子状態計算装置、コンピュータプログラム、及び記録媒体を提供した特許文献2を参照。しかしながら、前記モデル空間中で厳密解に到達する近似モデルの列により厳密解に到達する方法は、物質の物理的又は化学的性質のみを再現することに目的を限定すれば、より少ない有限回数の計算で厳密解の物性に到達することを用いてはいなかった。本発明は斯かる事情に鑑みてなされたものであり、局所密度近似法、一般化勾配近似法、ハートレーフォック法、配置間相互作用法、又は摂動MPMoellerPlesset法などの既存の近似モデルを出発点として、モデル空間中で与えられる密度汎関数変分法の変分原理に従って、複数の演算経路を探索しながら厳密解が示す物性に有限回の計算手続きで到達することができる電子状態計算方法、電子状態計算装置及びコンピュータプログラム提供することを目的とする。本発明に係る電子状態計算方法は、計算装置を用いて物質の電子状態を計算する方法において、前記計算装置は、物質の電子状態に対する近似解を与える複数の演算モデルの夫々を構成要素として含む集合を設定し、複数の電子軌道上に存在する少なくとも2つ以上の電子を含む電子系に働く有効相互作用を含む有効ハミルトニアンを、各演算モデルを特定するために定め、前記有効ハミルトニアンの自己無撞着解を前記集合内の各演算モデルを用いて計算する過程で、計算した自己無撞着解同士が連続的に変化する方向を特定しながら、前記集合がなす空間内で距離が近接する複数の演算モデルのうちで最適な演算モデルを量子力学的変分法に基づいて定め、前記最適な演算モデルを逐次更新する際に、前記有効ハミルトニアンの自己無撞着解による電子系の変分エネルギを評価し、評価した変分エネルギが計算すべき厳密解のエネルギに近接してゆき、しかも変分エネルギが単調減少凸関数をなすように演算モデルを更新して、一又は複数の変分エネルギの系列から前記電子状態の厳密解を計算することを特徴とする。本発明に係る電子状態計算方法は、局所密度近似法、一般化勾配近似法、又はハートレーフォック法を基に、有効相互作用として揺らぎを示す非局所型演算子を定める手続きを通じて前記有効ハミルトニアンを決定することを特徴とする。本発明に係る電子状態計算方法は、前記揺らぎが、クーロン相互作用とハートレー平均場項のずれである密度揺らぎを含むことを特徴とする。本発明に係る電子状態計算方法は、前記厳密解に近接してゆく変分エネルギの系列を与える一又は複数の演算モデルの系列のうち、電子状態を示す相が転移せずに厳密解と同一視される相を与え、かつ最小の計算ステップで前記有効ハミルトニアンの自己無撞着解に至る演算モデルを含むものを定めることを特徴とする。本発明に係る電子状態計算方法は、各演算モデルによる自己無撞着解をLAPW法LinearizedAugmentedPlaneWave法、PAW法ProjectorAugmentedWave法、又は数値基底展開法を用いた並列計算により求めることを特徴とする。本発明に係る電子状態計算装置は、物質の電子状態を計算する装置において、物質の電子状態に対する近似解を与える複数の演算モデルの夫々を構成要素として含む集合を設定する手段と、複数の電子軌道上に存在する少なくとも2つ以上の電子を含む電子系に働く有効相互作用を含む有効ハミルトニアンを、各演算モデルを特定するために定める手段と、前記有効ハミルトニアンの自己無撞着解を前記集合内の各演算モデルを用いて計算する過程で、計算した自己無撞着解同士が連続的に変化する方向を特定しながら、前記集合がなす空間内で距離が近接する複数の演算モデルのうちで最適な演算モデルを量子力学的変分法に基づいて定める手段と、前記最適な演算モデルを逐次更新する際に、前記有効ハミルトニアンの自己無撞着解による電子系の変分エネルギを評価する手段と、評価した変分エネルギが計算すべき厳密解のエネルギに近接してゆき、しかも変分エネルギが単調減少凸関数をなすように演算モデルを更新する手段と、一又は複数の変分エネルギの系列から前記電子状態の厳密解を計算する手段とを備えることを特徴とする。本発明に係る電子状態計算装置は、局所密度近似法、一般化勾配近似法、又はハートレーフォック法を基に、有効相互作用として揺らぎを示す非局所型演算子を定める手続きを通じて前記有効ハミルトニアンを決定するようにしてあることを特徴とする。本発明に係る電子状態計算装置は、前記揺らぎが、クーロン相互作用とハートレー平均場項のずれである密度揺らぎを含むことを特徴とする。本発明に係る電子状態計算装置は、前記厳密解に近接してゆく変分エネルギの系列を与える一又は複数の演算モデルの系列のうち、電子状態を示す相が転移せずに厳密解と同一視される相を与え、かつ最小の計算ステップで前記有効ハミルトニアンの自己無撞着解に至る演算モデルを含むものを定めることを特徴とする。本発明に係る電子状態計算装置は、各演算モデルによる自己無撞着解をLAPW法LinearizedAugmentedPlaneWave法、PAW法ProjectorAugmentedWave法、又は数値基底展開法を用いた並列計算により求めるようにしてあることを特徴とする。本発明に係るコンピュータプログラムは、コンピュータに、物質の電子状態を計算させるコンピュータプログラムにおいて、物質の電子状態に対する近似解を与える複数の演算モデルの夫々を構成要素として含む集合を設定させるステップと、複数の電子軌道上に存在する少なくとも2つ以上の電子を含む電子系に働く有効相互作用を含む有効ハミルトニアンを、各演算モデルを特定するために定めさせるステップと、前記有効ハミルトニアンの自己無撞着解を前記集合内の各演算モデルを用いて計算させる過程で、計算させた自己無撞着解同士が連続的に変化する方向を特定しながら、前記集合がなす空間内で距離が近接する複数の演算モデルのうちで最適な演算モデルを量子力学的変分法に基づいて定めさせるステップと、前記最適な演算モデルを逐次更新させる際に、前記有効ハミルトニアンの自己無撞着解による電子系の変分エネルギを評価させるステップと、評価させた変分エネルギが計算すべき厳密解のエネルギに近接してゆき、しかも変分エネルギが単調減少凸関数をなすように演算モデルを更新させ、一又は複数の変分エネルギの系列から前記電子状態の厳密解を計算させるステップとをコンピュータに実行させることを特徴とする。本発明にあっては、計算装置を用いて物質の電子状態を計算する方法、装置、コンピュータプログラムにおいて、物質の電子状態が示す物性を再現する有限回数の計算手続きを多配置参照密度汎関数理論に基づいて与える。本発明にあっては、多配置参照密度汎関数理論による演算モデルが与える決定方程式の自己無撞着解を決定することを通して、物質の厳密解が示す電子状態がもつ物性を、計算装置が自動的に決定する。本発明にあっては、厳密解が示す物性の再現に必要となる局所密度揺らぎなどの正定値性を備えた揺らぎ変数相関関数又は多重秩序変数が有限個であることが、密度汎関数理論における原理により保証される。本発明による場合は、物質の電子状態に対する近似解を与える演算モデルを出発点とし、電子間に働く有効相互作用を示す有効ハミルトニアンの自己無撞着解を求める過程で、モデル間の距離が接近していくように、計算過程で使用する演算モデルを決定する。これにより、新しい演算子空間の最適化手続きが与えられ、繰り込み思想に基づくダウンフォールディングの方法では不可能であった上位モデルからなる系列を形成することで、厳密解が与える物性の評価方法が有限回の計算手続きにより与えられる。本発明による場合、多種多様の演算モデルを含むモデル空間内で、厳密解に到達しうる複数の経路を与えることができるため、電子相関を取込む経路を適宜設定することができる。例えば、強相関効果をフェルミ準位から順次取り込み、エネルギスケール又は空間スケールに依存した全ての相関を考慮することでクーロン相互作用に漸近する演算子を導入することができる。本発明による場合は、自己無撞着解を求める演算を並列して実行することが可能であるため、並列計算を実現することにより計算速度が向上する。また、計算の各ステップで密度汎関数法を含む既存の第一原理計算法で実現されている並列化計算技法を用いることで、計算速度が向上する。本発明による場合は、局所密度近似法、一般化密度勾配近似法、ハートレーフォック法、配置間相互作用法、及び摂動MP法を含む既存の演算モデルを出発点として厳密解が示す物性に有限回の計算ステップで到達する手法を与えることができる。本発明による場合は、有効ハミルトニアンの表式が、フェルミ準位近傍の電子間のクーロン相互作用に漸近する有効相互作用を演算子として含む手法が与えられる。本発明による場合は、各演算モデルによる自己無撞着解をLAPW法LinearizedAugmentedPlaneWave法、PAW法ProjectorAugmentedWave法、又は数値基底展開の方法などを用いて、基底展開を行った高精度な一電子波動関数の表現方法を用いて求める手法が与えられる。本発明による場合は、モデル系とその安定解を初期値として用いることで、全エネルギ、一電子密度、正準相関関数を再現する時間発展方程式をさらに与え、物質の外場応答力学変形、電磁場印加に対する応答などを再現する高精度な第一原理電子状態計算の方法を与えることができる。本発明による場合は、与えられた計算機資源を利用する範囲内で効率良く厳密解が示す物性に到達し得る計算手続きを決定することができる。その結果、物質設計を計算機シミュレーションによる量子シミュレータを用いて実施する新しい量子デザイン手法において、スピンエレクロニクス、分子エレクトロニクスなどの量子素子形成要素の設計法の提供、エレクトロニクス応用によるエネルギ問題回避法の提供、低環境負荷又は戦略的元素選択による環境問題資源エネルギ問題回避型物質解の提供、新型センサー設計、生体親和材料設計、薬剤設計などの計算装置による物質設計が可能となる。以下、本発明をその実施の形態を示す図面に基づいて具体的に説明する。実施の形態1.本実施の形態に係る電子状態計算方法における計算原理は以下の通りである。本実施の形態では、電子状態の厳密解が示す物性を再現する最適化された演算モデルを決定する方法として、クーロン多体系が基底状態で発生する電子密度なる秩序変数を数値的に求めることを通して、物性の再現を行うことを目標とする。電子密度は実験で観測できる物理量であり、存在することが知られている。クーロン多体系における基底状態のエネルギE0は以下の表式を満たすため、数値的に評価を行うことができる。数式1エネルギ汎関数GバーXiεigiは、次の式により定義される。数式2また、ΔEバーXiεigiは、次の式により定義される。数式3ここで、GバーXiεigiは演算モデルが表現する電子の波動関数Ψに関するエネルギ汎関数を表し、ΔEバーXiεigiはクーロン系の変分エネルギを評価する際に発生するエネルギのずれを補正するエネルギ汎関数を表す。ここで、Ψは多粒子波動関数、Tは運動エネルギに係る演算子、VredXiは以下の数式4によって表される作用素である。なお、数式中では、演算子T及び作用素VredXiはハット付の文字で表記しているが、明細書中においてはそれぞれハットを付さずに表記する。また、数式2において、eは電子の電荷、nΨrは位置ベクトルrでのΨが与える電子密度、vextrは外部スカラーポテンシャルを表す。数式4ここで、Xiはモデルを定めるパラメータと演算子との組からなる集合を表す。Pは射影演算子である。また、::はノーマルオーダと呼ばれる演算子順序を指定するものであり、生成演算子を消滅演算子の前にフェルミオン演算子の交換規則に沿って順序交換してから用いることを指定するものである。また、Yn、Znは、一体波動関数の完全系が与える軌道l上での電子の消滅演算子clσから与えられる演算子である。fnlσ及びfnlσは複素定数である。αnは自己相互作用補正を制御するパラメータである。Ξnは、演算子Yn、Znに関する関数であり、期待値が有限な下界をもつことが分かっていることが望ましいが、特に多項式関数であってもよい。また、数式2におけるEεilocal[Ψ]及びEginonlocal[Ψ]は局所性をもつモデルエネルギ汎関数及び非局所性をもつモデルエネルギ汎関数であり、例えば、以下のように定義されるものを用いてよい。数式5ここで、εinは有界単調減少連続関数、gil1,l2,σは実数値係数であり、入力の際に与えられる。これらは計算プロセスの途中で定まるコーンシャム方程式の解を用いて最適化しなおすこともできるが、充分な計算規模を確保できたことにより数式4を通して取り込まれる揺らぎが充分であれば結果には影響しない。数式2,3,4及び5の有効相互作用によるモデルを導入していることが本願発明の特徴の1つとなる。演算子に関する関数Ξn及び有界単調減少連続関数εin及びgil1,l2,σを最適化することで、計算資源の範囲内で、高効率で高速かつ高精度な計算を可能とする最適化モデル系列を探索することができるようになる。数式4は、クーロン系の密度揺らぎをもとにチャンネル分解を行った有効多粒子相互作用を導入することによって、以下の式を具体的な表式として用いることが可能である。数式6数式1は、モデル空間において成立する変分原理を与えている。この原理は、Eginonlocal[Ψ]が特別な場合に限れば証明が公開されている例えば、KKusakabeJPhysSocJpn78114716を参照。そのため、不等式として数式1に類似の評価式が与えられる。等号を成立させる条件の存在は示されているが、有限の計算機資源による計算では厳密な等号は達成されない。しかしながら、変分エネルギに関する評価値の数列を得ることを通して、十分な計算結果の集積により基底状態のエネルギE0を数値的に評価することが可能となる特許文献2を参照。本出願では、この変分原理が広範なエネルギ汎関数に対して成立することを用いる。さらに、この計算過程を複数用意することによって、厳密解が示す物性を導出する有限回数の計算過程を複数与え、計算資源の範囲内で、高効率で高速かつ高精度な計算方法を見つけ出すことができる。数式1を与える原理は次のように厳密に示される。残差交換相関エネルギ汎関数ΔEXiεigi[Ψ]を、自己相互作用補正項を分離せずにクーロン揺らぎを与える表式を用いて、次式により導入する。数式7ここで、λに関する積分は、クーロン相互作用をλ倍に変化させて得られる制限付き探索法によるエネルギ汎関数数式8を参照のλ微分に関するルベーグ積分、及びFλΨのλディニ微分が有限個の飛びを生じている点での補正をもとに定義される量であり、存在が示されている。密度を変動させない範囲で、λディニ微分が飛びを生じさせる点において、残る相関による相転移が発生する。数式8このとき、数式1のΔEバーXiεigi[Ψ]を、数式7で定義されたΔEXiεigi[Ψ]に置き換えることによって、以下の表式を得ることが可能となる。数式9この数式9は、ΨXiεigiが変分原理を満たすことを通して得られるが、さらに次の不等式の成立を保証する。数式10図1は計算原理を説明する説明図である。数式1の表式に従い、まず、演算モデルを固定した上でエネルギ汎関数GバーXiεigiを最小化し、次に、最小化されたエネルギ汎関数GバーXiεigiと、この最小値を与える波動関数ΨXiεigiによって定められるΔEバーXiεigiとの和である変分エネルギの評価値を得る。続いて、演算モデルを変化させて同じ計算ステップを実行することで、変分エネルギの評価値を複数個得る。複数の演算モデルの夫々を構成要素として含む集合に、電子密度間の絶対値ノルムを用いて距離を導入することで得られる空間内で、演算モデルがもつ相互作用パラメータを連続的に変化させることを用いて、近接したモデルの系列を得る。その結果、相互作用パラメータに対して変分エネルギの変化が作る漸近曲線を得ることができる。この漸近曲線の収束値を数値的に評価することで、基底状態のエネルギE0の評価値を求める。図1は、この方法が複数個存在することを表している。演算モデルとは、GバーXiεigiと、この最小値を与える波動関数ΨXiεigiにより構成される。波動関数ΨXiεigiは、自己無撞着方程式を形成するGバーXiεigiの最小化問題を定める決定方程式に対して、自己無撞着性を備えるように決定することができる。後述するように、あるモデル例えば、図1のXiに関するモデル、及び他のモデル例えば、図1のXjに関するモデルの収束性について判定する。このとき、モデル系列{Xi}内でモデル列Xiが電子密度及び変分エネルギに関して収束しているか否かを判定する。また、モデル系列{Xj}内のモデル列Xjについても、電子密度及び変分エネルギに関して収束しているか否かを判定する。本実施の形態では、2種類の演算モデルを含むモデル列が共通の収束点として与える有効モデルを用いて物理量の算出を行う。演算モデルを構成するときには、物質の電子状態に対する近似解を与える既存のモデルである、局所密度近似法LDA、LDAU、一般化勾配近似法GGA、スピン依存GGA、メタGGA、ハートレーフォック法、配置間相互作用法、及び摂動MPMoellerPlesset法などを、初期条件として用いることができる。必要充分な揺らぎの寄与を数式4が表現する非局所演算子を用いてGバーXiεigiに取り込むことにより、計算が与えている近似解が、真の解厳密解の充分近傍にまで接近すると、ユニバーサルエネルギ汎関数Fλ1ΨからGバーXiεigiを定めるエネルギ汎関数数式11を参照まで、制限付き探索法によるエネルギ汎関数のパラメータ変化に関する微分が連続になり、演算モデルが与えるnΨが、クーロン系が与える密度と同一相に含まれるようになる。このとき、定めている数値的精度限界の範囲内で波動関数ΨXiεigiは厳密解が示す物性と同一の物性を与える。数式11次に、演算モデルが与えている近似解が、真の解厳密解の充分近傍にまで接近していることの確認を行う方法を開示する。複数の初期条件と、複数の揺らぎの導入順序を組み合わせて与えることにより、厳密解に漸近する演算モデル系列を複数与えることができる。このとき、系列ごとにモデル空間内でコーシー列を発生させる演算モデルの系列を常に構成することができる。それぞれの演算モデルにおいて明確に定義される一電子軌道が量子化を発生することを通して、上位モデル系列の発生手順が定義される。そこで、全ての系列が充分な揺らぎの取り込みがなされていることが系列ごとに定められる。同一系に対する複数系列が同じ物性を発生することは、電子密度が近接していくことを数値的に確認することでいつでも評価できる。よって、数値的に発生させられる全ての演算モデルが与えている近似解が、密度に関して同じ極限を与えるコーシー列を発生しているとき、厳密解に充分近傍にまで漸近していることが確認される。すなわち、本実施の形態では、多種多様の演算モデルを含むモデル空間内で、厳密解に到達しうる複数の経路を与えることができるため、電子相関を取込む経路を適宜設定することができる。例えば、強相関効果をフェルミ準位から順次取り込み、全スケールに依存した相関を考慮することでクーロン相互作用に近づけることができる。そこで、Xi,Xi1,Xi2,...のようにパラメータを変化させたときに表れる変分エネルギの点列に対して数値的収束を調べることにより、厳密解に漸近する方向を見出すことが可能となる。厳密解としての基底状態における物性が再現されていることは、電子密度がXiを変化させる全ての方向について変動しなくなることで定められる。このとき、変分エネルギの最小化が完了した時点の電子密度が得られている。量子力学系の定常状態に対する決定方程式を、連続座標空間を設定して、その上の微分方程式として得る際に、共変普遍性を満たす場の理論の表現から自由度の分離を行うことによる簡略化がなされた有効方程式を得る手続きが定義可能となるためには、時間的に不変な物理量期待値が与えるC数化された変量を与える自由度が見出されている必要がある。その例は、電子系のシュレディンガー方程式として採用される表現では、古典質点系として扱われている原子核配置であり、同時にマクスウェル方程式の解である古典場として扱われている静電場である。この課題設定によって電子系に対して付随して電子密度という秩序変数発生と、有効ハミルトニアンの決定が起こっている。このことは、現実の問題において自発的対称性の破れを伴う秩序変数が発生している状況では常に現れる。このようなC数化された物理変数発生を既知とし得ない純粋な多粒子系の量子状態決定問題に対しては、自由度の分離は不可能である。現実の物性を議論する多電子系では、低温での原子核の量子流体化が発生する場合を除くと、凝縮相において常に自由度の分離が可能であることが、電子密度の測定実験結果から確定されている。この不可能性を除いた量子力学系の定常状態に対する決定方程式を得る際にさらに別の時間的に不変な物理量期待値が与えるC数化されたある変量を与える自由度が見出されれば、密度だけでなくより一般化された秩序変数などC数化された物理変数を用いることで有効方程式を改善することが直ちに可能となる。そこで、本出願における密度汎関数変分法では、一般的に多体量子系に対して成立する方法論を与えることができる。密度汎関数の発生は外部スカラーポテンシャル中で多電子系の運動決定が行われることによる。エネルギ密度汎関数の評価をする際、従来の局所密度近似、局所スピン密度近似、一般化勾配近似、スピン依存一般化勾配近似、を含む汎用性が高いものの物質系によらずに関数形が汎関数として与えられている従来の方法を用いるのであれば、本願で開示する手法の適用が直ちに可能となる。本願発明者は、国際公開第2007/141942号パンフレットにて、多配置参照密度汎関数法に基づく有効多体電子計算として揺らぎ参照決定法を与えた。この決定法により、一電子密度に加えて、局所密度揺らぎなどの正定値性を備えた変数相関関数を再現するエネルギ汎関数を定義し数式12を参照、物理量の再現を行うことができる。数式12揺らぎ参照決定法を用いると、数式12の残差交換相関エネルギ汎関数に局所密度近似したモデルを用いることで、ハートレーフォック法、配置間相互作用法、及び摂動MP法などの近似的電子状態計算が与える局所揺らぎを再現する密度汎関数理論に基づいた演算モデルを直ちに構成することができる。したがって、ハートレーフォック法、配置間相互作用法、及び摂動MP法などの電子状態の近似解を与える演算モデルを初期条件として、演算モデルの系列を取り入れる揺らぎを定める演算子空間を拡張してゆくことにより与え、当該系列からコーシー列を見出すことで真の解の物性に到達する手続きを与えることができる。さらに、この決定プロセスにおいて、数値的精度限界を設定すると、有限の計算ステップで収束する計算プロセスを、点列{Xi}の発生を通して定めることができる。図2は本実施の形態に係る電子状態計算装置の内部構成を示すブロック図である。電子状態計算装置10は、CPU11、ROM13、RAM14、入力IF15、出力IF16、補助記憶装置17、及び記憶装置18を備えており、これらの各ハードウェアがバス12を介して互いに接続されている。ROM13には、ハードウェア各部を制御するための制御プログラムが格納されている。CPU11は、ROM13に格納されている制御プログラムをRAM14上にロードして実行することにより、前述したハードウェア各部の制御を行う。入力IF15には、マウス、キーボードなどの入力デバイス21が接続される。一方、出力IF16には、CRT、液晶ディスプレイなどの出力デバイス22が接続される。補助記憶装置17は、本実施の形態で説明する電子状態計算方法をコンピュータで実現するためのコンピュータプログラムを記録したFD、CD-ROMなどの記録媒体Mから前記コンピュータプログラムを読取るためのFDドライブ、CD-ROMドライブなどの読取装置を備える。補助記憶装置17によって読取られたコンピュータプログラムは記憶装置18に記憶される。演算手段としてのCPU11は、記憶装置18に記憶されたコンピュータプログラムを記憶手段としてのRAM14上にロードして実行することにより、装置全体を本発明に係る電子状態計算装置として機能させる。RAM14には、前記コンピュータプログムの他、入力デバイス21を通じて入力された各種情報、CPU11による演算の途中結果及び最終結果等が記憶される。なお、前記コンピュータプログラムを記録する記録媒体Mとしては、前述したFD及びCD-ROMの他に、MO、MD、DVD-ROM等の光ディスク、ハードディスク等の磁気記録媒体、ICカード、メモリカード、光カード等のカード型記録媒体、マスクROM、EPROMErasableProgrammableReadOnlyMemory、EEPROMElectricallyErasableProgrammableReadOnlyMemory、フラッシュROM等による半導体メモリの利用も可能である。また、インターネットを含む通信ネットワークを接続可能なシステムを構成することにより、通信ネットワークから前述のコンピュータプログラムをダウンロードするようにしてもよい。更に、ROM13の内部に前述のコンピュータプログラムを予め格納しておく構成であってもよい。記憶装置18は、例えばハードディスクドライブであり、補助記憶装置17によって読み取られたコンピュータプログラム、電子状態の計算に必要な初期データ、電子状態の計算により得られた途中結果、最終結果などが記憶される。また、記憶装置18が持つ記憶領域の一部は、物質の電子状態に対する近似解を与える既存のモデルを記憶するモデル記憶領域として利用される。この既存のモデルには、前述したように、LDA、LDAU、GGA、スピン依存GGA、メタGGA、ハートレーフォック法、配置間相互作用法、摂動MP法などが含まれる。記憶装置18に記憶されているコンピュータプログラム、初期データ、演算モデル等は、電子状態を計算する際に読み出され、RAM13に一時的に格納される。電子状態計算装置10のCPU11は、RAM13に格納されたコンピュータプログラムを実行することにより、電子状態の計算を行う。なお、本実施の形態では、記録媒体Mに記録された本発明に係るコンピュータプログラムを記憶装置18に記憶させて利用する構成としたが、予めROM13に組み込まれていてもよい。また、電子状態計算装置10が通信手段を備え、前記コンピュータプログラムを通信により取得する構成としてもよいことは勿論のことである。以下、具体的な計算手順について説明する。図3は電子状態計算装置10が実行する処理の手順を示すフローチャートである。電子状態計算装置10は、まず、初期値の設定を行うステップS11。初期値の設定は、計算対象となる物質の原子座標RIを設定し、有効相互作用の形式Xiと組み合わせるεi,giを定めることによって行う。初期値の設定を行うことにより、モデルが1つ定まる。原子座標RIの設定は、入力デバイス21を通じて原子座標RIのデータを受付け、RAM14に格納することによって行う。または、予め記憶装置18に原子座標RIのデータを記憶しておき、設定時に記憶装置18からそのデータを読み込み、RAM14に格納することによって行う。また、εi,giの関数形などはサブルーチンとして与えられているものとする。Xiは揺らぎ項の組を示す。揺らぎ項は、相互作用パラメータXnとαnとを含むように、例えば、数式6によって与えられる。次いで、電子状態計算装置10は、Xi=0i=1,2,3,...にセットするステップS12。ここで、Xi=0にセットするとは、数式6における相互作用パラメータXnをゼロにすることを意味する。次いで、電子状態計算装置10のCPU11は、設定された初期値を用いて、外部スカラーポテンシャルvextrを算出するステップS13。外部スカラーポテンシャルvextrの算出方法は公知であり、原子座標RIから一通りに決定される。次いで、電子状態計算装置10のCPU11は、有効系自己無撞着解を決定し、局在軌道の算出を行うステップS14。Xi=0における有効系自己無撞着解は、通常のLDAなどの計算ルーチンにより決定される。一体の展開基底がこの段階で初期値として与えられる。ユニタリ変換が原子座標RIに応じて決定され、局在軌道が算出される。次いで、電子状態計算装置10のCPU11は、ループカウンタiを1にセットした上でステップS15、作用素VredXiを算出するステップS16。このステップで、射影演算子Pとユニタリ変換の形式及び操作とを一連のものとして与える。すなわち、ループによって射影演算子とユニタリ変換の形式及び操作とが更新されてもよいこと意味する。更新によって、上位モデルに至る。次いで、電子状態計算装置10のCPU11は、Xiにおけるエネルギ汎関数の最適化を行うステップS17。自己無撞着解を持つモデルの列であるモデル列において、モデル同士が数値的に誤差範囲で同じ密度を与え、かつクーロン系の変分エネルギが最小化した値を持つことを確認したとき、収束したモデルXiが得られたと判定する。一体の展開基底がこの段階で再決定される。また、ユニタリ変換が原子座標RIに応じて決定され、局在軌道が算出される。次いで、電子状態計算装置10のCPU11は、電子密度nΨXiを算出するステップS18。電子密度nΨXiは、収束モデルに対してその波動関数ΨXiから算出される。次いで、電子状態計算装置10のCPU11は、Xi1におけるエネルギ汎関数の最適化を行うステップS19。本願では、射影演算子を変更してPとし、ユニタリ変換も同時に変更することにより、アップコンバージョンを行って、上位モデルXi1を作成する。この上位モデルの列において、収束したモデルXi1を得る。次いで、電子状態計算装置10のCPU11は、電子密度nΨXi1を算出するステップS20。電子密度nΨXi1は、収束モデルに対してその波動関数ΨXi1から算出される。次いで、電子状態計算装置10のCPU11は、モデルが収束したか否かを判定するステップS21。すでに下位モデルXiから上位モデルXi1までのモデル系列{Xi}を与えている。このモデル系列{Xi}内でモデル列Xiが電子密度及び変分エネルギに関して収束しているときに、収束したと判定する。モデル系列が収束していないと判定した場合S21:NO、電子状態計算装置10のCPU11は、ループカウンタiを1だけインクリメントした上でステップS22、処理をステップS16へ戻す。モデルが収束したと判定した場合S21YES、電子状態計算装置10のCPU11は、未計算のモデルが存在するか否かを判断するステップS23。未計算のモデル例えば、図1に示すXjに関するモデルが存在する場合にはS23:YES、処理をステップS11へ戻し、未計算のモデルについて初期値を設定する。そして、ステップS12~ステップS21の処理を新たに設定したモデルについて実行することにより、例えば、図1に示すXjに関するモデルについての収束性を判定する。未計算のモデルが存在しない場合S23:NO、すでに複数のモデル系列{Xi}、{Xj}が得られている。電子状態計算装置10のCPU11は、これらのモデル系列を含むモデル列が収束したか否かを判定するステップS24。この異なるモデル系列が収束先の上位モデルに関して、それらの与える電子密度が収束した値を得ているときに、モデル列が収束したと判定する。モデル列が収束していない場合S24:NO、まだ充分に揺らぎの取り込みに至っていないことから、モデル系列内の上位モデルを定める射影演算子とユニタリ変換とについて更新を行って、再びモデル系列をアップコンバージョンすべく、処理をステップS16へ戻す。モデル列が収束したと判定した場合S24:YES、得られた2種類の演算モデルを含むモデル列が共通の収束点として与える有効モデルを用いて物理量の算出を行うステップS25。以上のように、本実施の形態では、複数の一体有効モデルLDA,GGAなどから出発し、クーロン解に漸近するモデル系列を複数構成することができる。密度汎関数法が与えるエネルギ汎関数法が示す状態のレベル交差点が有限数である。また、一体有効ハミルトニアンが必ず量子化を示すことから、上位モデルを形成したときに有限回数でモデル系列が収束する。したがって、上記モデル系列は、有限回数の計算でクーロン解の示す物性に至る。また、本実施の形態では、モデル系列間の比較を電子密度に関する距離を用いて行うことができる点を利用しており、複数のモデル系列の中で最も計算効率が高い方法を有限の計算回数で決定することができる。なお、本実施の形態では、2種類の演算モデルを用いてモデル列の収束性を判定したが、3種類以上の演算モデルを用いてモデル列の収束性を判定するようにしてもよい。本実施の形態では、特定の演算モデル例えば、局所密度近似を出発点として、演算モデルを逐次更新しながら、演算を実行する構成としたが、各演算モデルによる自己無撞着解をLAPW法LinearizedAugmentedPlaneWaveMethod法を用いた並列計算、及び数値規定展開法により求める構成としてもよい。実施の形態2.実施の形態1で説明した電子状態計算方法に、アップコンバージョンのプロセスを導入することも可能である。実施の形態2では、実施の形態1の電子状態計算方法においてアップコンバージョンのプロセスを導入した形態について説明を行う。以下の数式13は、揺らぎ効果が重要になる一電子バンドや一電子軌道に対して射影演算子を導入し、その揺らぎ項が与える有効相互作用の強さを調整するパラメータとして遮蔽相互作用型の相互作用強度を導入したモデルによる、アップコンバージョン法の実施形態2に用いられるモデルの定義式である。数式13アップコンバージョンのプロセスは、初期モデルをLDAなどに取った場合には、2段階のプロセスとして纏めることができる。まず、揺らぎを導入する重要な一電子バンドや一電子軌道が、フェルミ準位近傍から順に取り入れられていく。ここで、遮蔽相互作用を定めるパラメータκを正の有限値に導入すると、中距離相関効果がまず取り込まれる。実施の形態1において説明したように、モデル列の収束点を示すモデルが得られたとき、このモデルを特徴付けるブロッホ軌道またはワニエ軌道を得ているため、これらの軌道により定められる射影演算子PAハットを導入する。ここで、λは、上記揺らぎを取り込んだモデルから最終的なクーロン系に至る単一連続パラメータによるアップコンバージョンの第二プロセスを示しているパラメータで、0から1までの数値をとる。具体的な実施の形態として、Sr2CuO3の電子状態計算を例にとる。GGAによるコーンシャム方程式の解を確認すると、フェルミ準位近傍にCu一原子当たり1つの1次元性バンドが生じていることが直ちに確認できる。このバンド上のワニエ軌道に対して射影演算子を導入する。遮蔽相互作用を取り込むことは、ワニエ軌道上の対角型相関項のみからアップコンバージョンの第一プロセスを始めることに相当する。そこで、相関パラメータとしては、遮蔽相互作用の積分が与えるUを単一パラメータとして採用することになる。相関パラメータを導入した多配置参照密度汎関数法計算は、自己無撞着に解を与えることが出来る。このパラメータUを変化させていくと、多配置参照密度汎関数法による多電子状態が、量子相関発生した状態ベクトルを自己無撞着解として、密度を伴って与えられる。密度汎関数変分法の変分原理に従って、ΔEバーを評価すると、変分エネルギがUの関数として得られる。図4は変分エネルギの評価結果の一例を示す図である。図4において、横軸は、相関パラメータU単位はエレクトロンボルト、縦軸は変分エネルギ単位はリュードベリを表す。ここで、モデルが与える変分エネルギは、相関パラメータUの連続関数となることが分かっている。そこで、この図から明らかなように、密度汎関数変分法の原理によって、変分エネルギを最小化する最適な相関パラメータを決定することができる。決定された相関パラメータのモデルをさらに初期条件として、アップコンバージョンの第二プロセスを実行することもできる。その場合には、λを変化させて残差相関を取り込んでいくことによって、モデル系列の変動を起こさせ、その変化を数値的に追跡することで、モデル空間での相転移発生の有無を検証することにより行うことができる。ここで、λを変化させたときのモデル系列が与えるエネルギの評価は、摂動展開グリーン関数法、モンテカルロ法などによって行うことができる。摂動展開グリーン関数法においては、相関パラメータUが有限のときに定まった局所相関効果を含むグリーン関数G0を求め、それをもとにして、長距離相互作用効果をダイソン方程式の解析を通して行うことができる。この場合の近似的評価法としては、GW近似における入力グリーン関数にG0を用いて、遮蔽効果を取り込むことによる、ワンショット計算が、密度を変化させないことを通して、アップコンバージョンプロセスの収束が確認される。密度変動が見えた場合には、相関項をより加えることを射影演算子PAハットの再定義を通して行い、再び多配置参照密度汎関数法計算を行うことで、アップコンバージョンの第一プロセスを再度行う。モンテカルロ法においては、相関パラメータUが有限のときに定まった局所相関効果を含む状態ベクトルを求め、それを基にした多体波動関数を拡散モンテカルロ法の節固定のための試行関数に採用することで、極短距離相関効果を取り込んだ計算を実施することができる。この場合の近似的評価法としては、節固定近似した拡散モンテカルロ法、変分モンテカルロ法等を採用することができる。実施の形態3.実施の形態1及び2で説明した電子状態計算方法を用いることにより、厳密解の示す物性を再現しながら、計算回数の最も少ない、最適化された演算モデルを定めることができる。実施の形態3では、最適化された演算モデルの選定手法について説明する。以下に示す数式14は、クーロン系と単一連続パラメータで凸性をもって繋がる連続モデル系列を与える、モデル汎関数を示したものである。Vハットチルダλは、λ微分から正定値演算子を発生する相関項として与えられたものである。λは0から1までの値をとり、λ=1において、もとのクーロン相互作用する系が再現されることになる。数式14λ=0においては、相関項に対して、凸性を保証するために設定された射影演算子を用いて定義されるモデル相互作用項が含まれることになる。この項は、数式13のλ=0での解を元にして構築することが出来る。数式14のλ=0でのモデルの最適解は、そのエネルギの評価を摂動展開グリーン関数法によって行うことができるものである。これが、射影演算子PAの再定義を用いるアップコンバージョンの第一プロセスである。このとき、密度を変動させずにエネルギ評価が可能な場合を探索することができる。次に、λを0から1まで変化させたときの効果を、変分計算によって評価する。この数式14が与えるアップコンバージョンの第二プロセスにおいては、nハットチルダrのチルダ付きの演算子が与えるモデル汎関数のλ微分が正定値であることを通して、相転移回数の有限性定理が成立する条件が与えられる。この条件は、アップコンバージョンプロセスにおける密度の不変性である。したがって、計算実施が可能な第二プロセスにおける密度の不変性確認によって、それ以上のアップコンバージョンプロセスが必要でないことが確認され、全プロセスの収束判定が、有限の数値計算プロセス内で定めた誤差範囲において可能となる。図5は最適モデル決定のプロセスを与える原理を説明する説明図である。既に、実施の形態1を通して、複数の収束モデル列が与えられているとする。このとき、λ変形されたモデルを数式14に従って構成することができる。このλ変形プロセスにおいて、既に充分に密度に関する収束が得られている。密度を変動させない範囲で残る相関による相転移点の発生は、高々有限回数である。なお、図5では明示されていないが、相転移点はパラメータλの変動に対して現れるグラフの屈曲点として生じ、この点においてλディニ微分が飛びを示す。このとき、厳密解の示す物性を再現するモデルにおいては、内部秩序変数の変化によるモデル空間における相転移発生が生じないことを用いている。収束した密度を再現し、充分な相関を取り込んだモデルの収束点から始めて、この相転移点が発生せずにλ=0まで至ることができるモデルを選択することによって、厳密解の示す物性を再現しながら、計算回数の最も少ない、最適化されたモデルを選び出すことができる。図5の原理からは、相転移発生を通した確認が出来ることになるが、これを数値計算で実現するための方法が、図6において示される、電子状態計算装置が最適モデル決定のために実行する処理の手順を示すフローチャートである。実施の形態1の図3により示したフローチャートのステップS11からステップS24により、収束モデル列が複数得られているとする。次に、電子状態計算装置は、収束密度を与える複数のモデル間でλ=0での相の差異を確認することにより、密度が収束していると判定されていながら、内部秩序変数の変動により相転移が発生しているか否かを判断するステップS31。λ変形において相転移発生があると確認した場合、例えば、図5の2つのモデル系列が得られている場合、S31:YES、それらのモデル列を棄却し、収束したモデル系列が同一の相を与えることを確認するまで、ステップS16に戻って揺らぎを更に取り込んだモデル系列の探索を行う。ステップS31で相転移発生がないと確認した場合S31:NO、同一相を与える複数のモデル系列を含むモデル列が得られる。このとき、λ変形に対する相転移発生は、アルゴリズム上、摂動展開グリーン関数法によるアップコンバージョンの第二プロセスが、密度変動を生じないことに加えて、各種応答の変動が無いことによって確認されることになる。このとき、複数のモデル系列が、λ変形の上位において同一の物性を示すことが系の真の解を得たことの確認を与えることになる。この条件が満たされたとき、λ変形プロセスの下位モデルまでにおいて、相転移点の発生が無いものを選択することで、最適モデル判定を行うステップS32。最後に、最適モデルを用いた物理量の算出を行いステップS33、計算規模に対して最も適切に選択されたモデル系による物性評価が可能となる。以上の実施の形態では、既知結晶構造に対する電子状態の最適化計算について説明したが、分子動力学の手法を用いて系のダイナミックスを計算することも可能である。系のダイナミックスを計算するためには、ある時間発展ステップにおける有効多体波動関数を求め、バルク秩序についてMD計算MDMolecularDynamicsを行う。このMD計算では、まず、拡張コーンシャム方程式により得られた電子密度、外部スカラーポテンシャル、及び局在軌道を用いて、原子間力及び応力を求め、得られた原子間力及び応力から物質構造変化率を算出し、微小時間後の物質の構造を順次決定する手法が採用される。また、本実施の形態、または分子動力学法の形態により得られる高精度に決定された基底状態を、引き続いて時間依存問題に対する初期条件として用いることができる。この方法によって、電子密度と電流密度を基本となる秩序変数として用いた時間依存電流密度汎関数法を実現することができる。"