特許 7314097 シミュレーション方法、シミュレーション装置、及びプログラム

(19)【発行国】日本国特許庁(JP)

(12)【公報種別】特許公報(B2)

(11)【特許番号】

(24)【登録日】2023-07-14

(45)【発行日】2023-07-25

(54)【発明の名称】シミュレーション方法、シミュレーション装置、及びプログラム

(51)【国際特許分類】

   G16Z 99/00 20190101AFI20230718BHJP

【ＦＩ】

G16Z99/00

【請求項の数】 3

(21)【出願番号】P 2020101739

(22)【出願日】2020-06-11

(65)【公開番号】P2021196757

(43)【公開日】2021-12-27

【審査請求日】2022-07-13

(73)【特許権者】

【識別番号】000002107

【氏名又は名称】住友重機械工業株式会社

(74)【代理人】

【識別番号】100105887

【弁理士】

【氏名又は名称】来山 幹雄

(72)【発明者】

【氏名】松村 裕也

【審査官】宮地 匡人

(56)【参考文献】

【文献】特開２０１９－１６４６１０（ＪＰ，Ａ）

【文献】国際公開第１９９７／０４６９４９（ＷＯ，Ａ１）

【文献】特開平０９－１４６９２１（ＪＰ，Ａ）

【文献】特開２００８－０５２３０８（ＪＰ，Ａ）

【文献】橋本 保，氷に関する最近の分子動力学シミュレーション －氷の構造相転移－，日本雪氷学会誌雪氷，1996年09月，Vol.58 No.5，pp.422-424

【文献】上田 顕，分子シミュレーション，初版，株式会社裳華房，2003年10月25日，pp.38-49

(58)【調査した分野】(Int.Cl.，ＤＢ名)

Ｇ１６Ｚ ９９／００

ＪＳＴＰｌｕｓ（ＪＤｒｅａｍＩＩＩ）

ＪＳＴ７５８０（ＪＤｒｅａｍＩＩＩ）

(57)【特許請求の範囲】

【請求項1】

粒子間の相互作用ポテンシャルが異方性を持つ複数の粒子の挙動を、分子動力学法を用いて解析するシミュレーション方法であって、
粒子間の相互作用ポテンシャルとして、二体間相互作用ポテンシャルと三体間相互作用ポテンシャルとを合成したポテンシャルを用い、
前記二体間相互作用ポテンシャルとしてレナードジョーンズポテンシャルを用い、三体間相互作用ポテンシャルとして単原子水モデルの三体間相互作用ポテンシャルを用い、
シミュレーション対象の粒子が水分子であり、
前記二体間相互作用ポテンシャルＵ _２ （ｒ）として、

【数1】

（ｒは２つの粒子間の距離）
を用い、前記三体間相互作用ポテンシャルＵ _３ （ｒ _０１ ，ｒ _０２ ，θ）として、

【数2】

（ｒ _０１ 、ｒ _０２ は、着目する粒子と他の２つの粒子のそれぞれとの距離、θは、着目する粒子から他の２つの粒子に向かうベクトルのなす角度）
を用い、θ _０ ＝１０９．４７°、λ≦２３．１５×２、０．９≦γ≦１．５、１．７≦ａ≦１．９、１．０×１０ ^－２１ Ｊ≦ε≦１．０×１０ ^－１９ Ｊ、１．２５×１０ ^－１０ ｍ≦σ≦３．７５×１０ ^－１０ ｍであるシミュレーション方法。

【請求項2】

シミュレーション条件を取得するシミュレーション条件取得部と、
取得したシミュレーション条件に基づいて、粒子間の相互作用ポテンシャルとして、二体間相互作用ポテンシャルと三体間相互作用ポテンシャルとを合成したポテンシャルを用いて分子動力学法により、複数の粒子の挙動を解析する演算部と、
前記演算部による解析結果を出力する出力制御部と
を有し、
前記演算部は、前記二体間相互作用ポテンシャルとしてレナードジョーンズポテンシャルを用い、三体間相互作用ポテンシャルとして単原子水モデルの三体間相互作用ポテンシャルを用い、
シミュレーション対象の粒子が水分子であり、
前記演算部は、
前記二体間相互作用ポテンシャルＵ _２ （ｒ）として、

【数3】

（ｒは２つの粒子間の距離）
を用い、前記三体間相互作用ポテンシャルＵ _３ （ｒ _０１ ，ｒ _０２ ，θ）として、

【数4】

（ｒ _０１ 、ｒ _０２ は、着目する粒子と他の２つの粒子のそれぞれとの距離、θは、着目する粒子から他の２つの粒子に向かうベクトルのなす角度）
を用い、θ _０ ＝１０９．４７°、λ≦２３．１５×２、０．９≦γ≦１．５、１．７≦ａ≦１．９、１．０×１０ ^－２１ Ｊ≦ε≦１．０×１０ ^－１９ Ｊ、１．２５×１０ ^－１０ ｍ≦σ≦３．７５×１０ ^－１０ ｍであるシミュレーション装置。

【請求項3】

シミュレーション条件を取得する機能と、
取得したシミュレーション条件に基づいて、粒子間の相互作用ポテンシャルとして、二体間相互作用ポテンシャルと三体間相互作用ポテンシャルとを合成したポテンシャルを用いて分子動力学法により、複数の粒子の挙動を解析する機能と、
解析結果を出力する機能と
をコンピュータに実現させるプログラムであって、
前記二体間相互作用ポテンシャルとしてレナードジョーンズポテンシャルを用い、三体間相互作用ポテンシャルとして単原子水モデルの三体間相互作用ポテンシャルを用い、
シミュレーション対象の粒子が水分子であり、
前記二体間相互作用ポテンシャルＵ _２ （ｒ）として、

【数5】

（ｒは２つの粒子間の距離）
を用い、前記三体間相互作用ポテンシャルＵ _３ （ｒ _０１ ，ｒ _０２ ，θ）として、

【数6】

（ｒ _０１ 、ｒ _０２ は、着目する粒子と他の２つの粒子のそれぞれとの距離、θは、着目する粒子から他の２つの粒子に向かうベクトルのなす角度）
を用い、θ _０ ＝１０９．４７°、λ≦２３．１５×２、０．９≦γ≦１．５、１．７≦ａ≦１．９、１．０×１０ ^－２１ Ｊ≦ε≦１．０×１０ ^－１９ Ｊ、１．２５×１０ ^－１０ ｍ≦σ≦３．７５×１０ ^－１０ ｍであるプログラム。

【発明の詳細な説明】

【技術分野】

【0001】

本発明は、シミュレーション方法、シミュレーション装置、及びプログラムに関する。

【背景技術】

【0002】

分子間（粒子間）の相互作用ポテンシャルとしてレナードジョーンズポテンシャルを用いた分子動力学法（ＭＤ法）によるシミュレーションは、不活性ガスの挙動を精度よく再現する。繰り込み群分子動力学法（ＲＭＤ法）も広い意味でＭＤ法の一種であり、本明細書において、ＭＤ法は、狭義のＭＤ法とＲＭＤ法との両方を含むものとする。

【0003】

レナードジョーンズポテンシャルを用いたＭＤ法は、不活性ガスの挙動を精度よく再現するが、分子間の相互作用ポテンシャルが異方性を持つ水分子等の挙動を精度よく再現できない場合がある。水分子を構成する酸素原子及び水素原子の挙動を解析することにより、水分子の挙動の解析精度を高めることができる。

【0004】

また、単原子水（ｍｏｎａｔｏｍｉｃ ｗａｔｅｒ）モデル（以下、ｍＷモデルという。）を用いた解析方法が知られている（例えば、下記の非特許文献１参照。）。ｍＷモデルでは、水分子間の相互作用ポテンシャルとして、二体間相互作用ポテンシャルと三体間相互作用ポテンシャルとを合成したポテンシャルが用いられる。なお、レナードジョーンズポテンシャルに、静電相互作用ポテンシャルを組み合わせたポテンシャルを作成し、このポテンシャルを高分子の分子運動シミュレーションに適用する技術が下記の特許文献１に開示されている。

【先行技術文献】

【特許文献】

【0005】

【文献】特開２０１１－２４８５４１号公報

【非特許文献】

【0006】

【文献】V. Molinero and E. B. Moore, "Water Modeled As an Intermediate Element between Carbon and Silicon", J. Phys. Chem. B 113, 4008-4016 (2009)

【発明の概要】

【発明が解決しようとする課題】

【0007】

ＭＤ法を適用するときの時間刻み幅は、一般的に粒子の質量の平方根に基づいて設定される。一般的に、粒子の質量が軽くなるにしたがって、時間刻み幅を短くする。水分子を酸素原子と水素原子で表して、酸素原子と水素原子の挙動を解析する方法では、ＭＤ法を適用するときの時間刻み幅が水素原子の質量に律速される。このため、水素原子の質量に基づいて時間刻み幅を短くしなければならない。また、クーロン力を考慮する必要があるため、カットオフ距離を長くしなければならない。例えば、原子間のカットオフ距離を、原子の直径の１０～２０倍程度にしなければならない。計算量はカットオフ距離の３乗にほぼ比例して増大するため、計算量が膨大になってしまう。

【0008】

従来のｍＷモデルで用いられる相互作用ポテンシャルは、カットオフ距離が気相の水分子間の距離に比べて短いため、気相の水分子の解析に適していない。このため、気相と液相との相転移を伴う条件下で高精度の解析を行うことができない。

【0009】

本発明の目的は、粒子間の相互作用ポテンシャルが異方性を持つ場合でも、計算量の増大を抑制し、かつ相転移を伴う解析を行うことが可能なシミュレーション方法、シミュレーション装置、及びプログラムを提供することである。

【課題を解決するための手段】

【0010】

本発明の一観点によると、
粒子間の相互作用ポテンシャルが異方性を持つ複数の粒子の挙動を、分子動力学法を用いて解析するシミュレーション方法であって、
粒子間の相互作用ポテンシャルとして、二体間相互作用ポテンシャルと三体間相互作用ポテンシャルとを合成したポテンシャルを用い、
前記二体間相互作用ポテンシャルとしてレナードジョーンズポテンシャルを用い、三体間相互作用ポテンシャルとして単原子水モデルの三体間相互作用ポテンシャルを用いるシミュレーション方法が提供される。

【0011】

本発明の他の観点によると
シミュレーション条件を取得するシミュレーション条件取得部と、
取得したシミュレーション条件に基づいて、粒子間の相互作用ポテンシャルとして、二体間相互作用ポテンシャルと三体間相互作用ポテンシャルとを合成したポテンシャルを用いて分子動力学法により、複数の粒子の挙動を解析する演算部と、
前記演算部による解析結果を出力する出力制御部と
を有し、
前記演算部は、前記二体間相互作用ポテンシャルとしてレナードジョーンズポテンシャルを用い、三体間相互作用ポテンシャルとして単原子水モデルの三体間相互作用ポテンシャルを用いるシミュレーション装置が提供される。

【0012】

本発明のさらに他の観点によると、
シミュレーション条件を取得する機能と、
取得したシミュレーション条件に基づいて、粒子間の相互作用ポテンシャルとして、二体間相互作用ポテンシャルと三体間相互作用ポテンシャルとを合成したポテンシャルを用いて分子動力学法により、複数の粒子の挙動を解析する機能と、
解析結果を出力する機能と
をコンピュータに実現させるプログラムであって、
前記二体間相互作用ポテンシャルとしてレナードジョーンズポテンシャルを用い、三体間相互作用ポテンシャルとして単原子水モデルの三体間相互作用ポテンシャルを用いるプログラムが提供される。

【発明の効果】

【0013】

粒子間の相互作用ポテンシャルが異方性を持つ場合でも、計算量の増大を抑制し、かつ相転移を伴う解析を行うことが可能である。

【図面の簡単な説明】

【0014】

【図1】図１Ａは、二体間相互作用ポテンシャルＵ_２が作用する２つの粒子Ｐ_ｉ，Ｐ_ｊを示す模式図であり、図１Ｂは、三体間相互作用ポテンシャルＵ_３が作用する３つの粒子Ｐ_ｉ、Ｐ_ｊ、Ｐ_ｋを示す模式図である。

【図2】図２は、二体間相互作用ポテンシャルＵ_２及び三体間相互作用ポテンシャルＵ_３を確定する手順のフローチャートである。

【図3】図３Ａは、ステップＳ２（図２）において分子動力学法により解析した結果の粒子の分布を示す斜視図であり、図３Ｂは、密度分布の一例を示すグラフである。

【図4】図４は、解析によって求めた飽和蒸気線を、解析対象である水の飽和蒸気線を比較して示すグラフである。

【図5】図５は、実施例によるシミュレーション装置のブロック図である。

【図6】図６は、上記実施例及び比較例によるシミュレーション方法を用いて求めた水の温度と圧力との関係を示すグラフである。

【図7】図７は、実施例及び比較例によるシミュレーション方法を用いて求めた水の飽和蒸気線、及び水の実際の物性値としての飽和蒸気線を示すグラフである。

【図8】図８は、実施例及び比較例によるシミュレーション方法を用いて求めた水の表面張力、及び水の実際の物性値としての表面張力を示すグラフである。

【発明を実施するための形態】

【0015】

以下に説明する実施例では、例えば水分子のように、複数の原子で構成され、相互作用ポテンシャルが異方性を持つ粒子を一つの粒子に粗視化し、相互作用ポテンシャルのパラメータを調整することにより、粒子の挙動の高精度の解析を行う。図１Ａ～図４Ｂを参照して、解析に用いる相互作用ポテンシャルのパラメータの決定方法について説明する。

【0016】

本実施例では、粒子間の相互作用ポテンシャルとして、二体間相互作用ポテンシャルＵ_２と、三体間相互作用ポテンシャルＵ_３とを合成したポテンシャルを用いる。以下、図１Ａ及び図１Ｂを参照して、二体間相互作用ポテンシャルＵ_２及び三体間相互作用ポテンシャルＵ_３について説明する。

【0017】

図１Ａは、二体間相互作用ポテンシャルＵ_２が作用する２つの粒子Ｐ_ｉ，Ｐ_ｊを示す模式図である。２つの粒子Ｐ_ｉ、Ｐ_ｊの間の距離をｒ_ｉｊと表記する。２つの粒子Ｐ_ｉ、Ｐ_ｊの間に作用する二体間相互作用ポテンシャルＵ_２として、下記の一般形で表されるレナードジョーンズポテンシャルを用いる。

【数1】

εは、エネルギの次元［Ｊ］を持つパラメータであり、σは距離の次元［ｍ］を持つパラメータである。一例として、斥力項の次数ｐとしてｐ＝１２またはｐ＝８が用いられ、引力項の次数ｑとして、ｑ＝６が用いられる。ｐ＝１２、ｑ＝６としたポテンシャルは、レナードジョーンズポテンシャルの代表例であり、ｐ＝１２、ｑ＝６のポテンシャルを狭義のレナードジョーンズポテンシャルという場合がある。本明細書では、ｐ＝１２、ｑ＝６に限定しない一般形をレナードジョーンズポテンシャルという。

【0018】

図１Ｂは、三体間相互作用ポテンシャルＵ_３が作用する３つの粒子Ｐ_ｉ、Ｐ_ｊ、Ｐ_ｋを示す模式図である。着目する粒子Ｐ_ｉと他の２つの粒子Ｐ_ｊ、Ｐ_ｋとの間の距離を、それぞれｒ_ｉｊ、ｒ_ｉｋと表記する。着目する粒子Ｐ_ｉから他の２つの粒子Ｐ_ｊ、Ｐ_ｋのそれぞれに向かう２つのベクトルのなす角度をθ_ｉｊｋと表記する。このとき、２つの粒子Ｐ_ｊ、Ｐ_ｋによって、着目する粒子Ｐ_ｉに作用する三体間相互作用ポテンシャルＵ_３として、単原子水モデルで用いられる三体間相互作用ポテンシャルを用いる。三体間相互作用ポテンシャルＵ_３は下記の式で表される。

【数2】

ここで、ε、σは、式（１）のパラメータε、σと同一である。λ、γ、ａは、無次元のパラメータである。

【0019】

ＭＤ法を用いて粒子の挙動を解析するときには、二体間相互作用ポテンシャルＵ_２及び三体間相互作用ポテンシャルＵ_３を合成した合成相互作用ポテンシャルＵを用いる。着目する粒子を粒子Ｐ_ｉと表記し、粒子Ｐ_ｉに作用する他の粒子を粒子Ｐ_ｊ、粒子Ｐ_ｋと表記する。位置ベクトルｒ_ｉの位置の粒子Ｐ_ｉに作用する相互作用ポテンシャルＵ（ｒ_ｉ）は、以下の式に示すように、二体間相互作用ポテンシャルＵ_２と三体間相互作用ポテンシャルＵ_３との和で表される。

【数3】

【0020】

粒子Ｐ_ｉの運動方程式は、以下の式で表される。

【数4】

ここで、ｍは粒子の質量である。

【0021】

次に、図２を参照してパラメータε、σ、λ、γ、ａを決定し、二体間相互作用ポテンシャルＵ_２及び三体間相互作用ポテンシャルＵ_３を確定する手順について説明する。

【0022】

図２は、二体間相互作用ポテンシャルＵ_２及び三体間相互作用ポテンシャルＵ_３を確定する手順のフローチャートである。まず、二体間相互作用ポテンシャルＵ_２（式（１））及び三体間相互作用ポテンシャルＵ_３（式（２））、及び運動方程式（式（４））をパラメータε、σ、及びｍで無次元化し、残りのパラメータλ、γ、及びａに適当な値を設定する（ステップＳ１）。例えば、無次元化することにより得られるポテンシャルＵ’、長さｒ’、時間ｔ’、温度Ｔ’、質量ｍ’を、以下の式で表す。

【数5】

【0023】

複数の温度条件で、無次元化された運動方程式を用いてＭＤ法により粒子の挙動を解析し、解析結果から解析対象の無次元化された物性値を求める（ステップＳ２）。本実施例では、解析対象として水を採用し、水の物性値として飽和蒸気線を採用する。無次元化に用いたパラメータε、σに適当な値を設定して、無次元化された物性値に次元を付与する。次元が付与された物性値が、解析対象の実際の物性値に近づくように、パラメータε、σをフィッティングする（ステップＳ３）。パラメータε、σのフィッティング処理を行うには、例えば、ユーザが、パラメータの上限値と下限値、及び刻み幅を指定し、演算部３２が、与えられた範囲内で、与えられた刻み幅でパラメータの値を変化させて飽和蒸気線を求めるようにするとよい。

【0024】

パラメータε、σを最適化して得られた次元が付与された物性値と、解析対象の実際の物性値とを比較し、解析結果から得られた物性値が適切か否か判定する（ステップＳ４）。解析結果から得られた物性値が適切ではない場合、ステップＳ１からステップＳ３までを再度実行する。ステップＳ１を再度実行する際に、パラメータλ、γ、ａの値を調整する。

【0025】

解析結果から得られた物性値が適切である場合は、パラメータε、σ、λ、γ、ａを、適切な物性値を導出した時の値に決定し、二体間相互作用ポテンシャルＵ_２及び三体間相互作用ポテンシャルＵ_３を確定する（ステップＳ５）。

【0026】

次に、図３Ａ～図４を参照して、ステップＳ２、Ｓ３（図２）の処理について説明する。
図３Ａは、ステップＳ２（図２）においてＭＤ法により解析した結果の粒子の分布を示す斜視図である。解析領域１０を、高さ方向の寸法が最も大きな直方体とした。解析領域１０の高さ方向の中央部の約１／３の領域に粒子を配置し、その上下の領域には粒子を配置しない状態を初期状態とした。解析領域１０の６個の面に、それぞれ周期境界条件を適用した。なお、重力は考慮していない。定常状態において、解析領域１０の高さ方向のほぼ中央部に、粒子密度の高い液相領域１１が現れ、その上下に、粒子密度の低い気相領域１２が現れている。

【0027】

高さ方向の密度分布は、高さ方向の座標をｘとして以下の式にフィットさせることができる。

【数6】

ここで、ρ_ｇ、ρ_ｌ、ｘ_ｃ、ｗは、それぞれ気相領域１２の密度、液相領域１１の密度、境界の位置、及び境界の厚さである。

【0028】

図３Ｂは、密度分布の一例を示すグラフである。横軸は、高さ方向の無次元化した位置を表し、縦軸は、無次元化密度を表す。グラフ中の黒丸記号は、解析結果から求まる密度を示し、曲線は式（６）で近似した密度を示す。複数の温度条件で解析を行って得られた密度分布から、飽和蒸気線の温度及び密度を得ることができる。さらに、温度ごとに、圧力テンソルから圧力及び表面張力を求めることができる。

【0029】

図４は、解析によって求めた飽和蒸気線を、解析対象である水の物性値としての飽和蒸気線と比較して示すグラフである。横軸は、密度を単位「ｋｇ／ｍ^３」で表し、縦軸は、温度を単位「Ｋ」で表す。グラフ中の丸記号は、温度条件ごとに行った解析結果から求めた密度を示し、破線は、水の物性値としての飽和蒸気線を示す。グラフの左側の破線より左側の領域では気体のみが存在し、右側の破線より右側の領域では液体のみが存在する。解析結果から得られる飽和蒸気線が、水の実際の物性値としての飽和蒸気線にほぼ一致するように、パラメータε、σ、λ、γ、ａに適切な値が設定されている。解析対象が水の場合、解析結果から求まる飽和蒸気線と、実際の飽和蒸気線とがほぼ一致するように最適化されたパラメータの値は以下のとおりである。
λ＝２３．１５
γ＝０．９５
ａ＝１．８
ε＝２．３６３×１０^－２０Ｊ
σ＝２．３５４×１０^－１０ｍ
なお、式（１）の（ｐ、ｑ）＝（１２，６）とし、水分子の質量ｍ＝２．９８８×１０^－２６ｋｇとし、式（２）のθ_０＝１０９．４７°とした。θ_０の値は、正四面体の中心と２つの頂点とを結ぶ線分のなす角度である。

【0030】

次に、図５を参照して、実施例によるシミュレーション装置について説明する。
図５は、実施例によるシミュレーション装置のブロック図である。実施例によるシミュレーション装置は、処理装置３０、入力装置３８、及び出力装置３９を有する。処理装置３０は、シミュレーション条件取得部３１、演算部３２、及び出力制御部３３を含む。

【0031】

図５に示す各ブロックは、ハードウェア的には、コンピュータの中央処理ユニット（ＣＰＵ）をはじめとする素子や機械装置で実現することができ、ソフトウェア的にはコンピュータプログラム等によって実現することができる。図５では、ハードウェア及びソフトウェアの連携によって実現される機能ブロックが示されている。従って、これらの機能ブロックは、ハードウェア及びソフトウェアの組み合わせによって、種々の態様で実現することが可能である。

【0032】

処理装置３０は入力装置３８及び出力装置３９と接続される。入力装置３８は、処理装置３０で実行される処理に関係するユーザからのコマンド及びデータの入力を受ける。入力装置３８として、例えばユーザが操作を行うことにより入力を行うキーボード、ディスプレイとポインティングデバイス、インターネット等のネットワークを介して入力を行う通信装置、種々のリムーバブルメディアから入力を行う読取装置等を用いることができる。

【0033】

シミュレーション条件取得部３１は、ユーザが入力装置３８に入力したシミュレーション条件を取得する。シミュレーション条件には、シミュレーションに必要な種々の情報が含まれる。例えば、シミュレーション対象の粒子（分子）に関する物理量、シミュレーションの初期条件、境界条件等が含まれる。

【0034】

演算部３２に、図２に示す方法で決定された二体間相互作用ポテンシャルＵ_２及び三体間相互作用ポテンシャルＵ_３を定義するパラメータが記憶されている。演算部３２は、パラメータが決定されている二体間相互作用ポテンシャルＵ_２及び三体間相互作用ポテンシャルＵ_３を用い、入力されたシミュレーション条件に基づいてＭＤ法により粒子の挙動をシミュレーションする。

【0035】

出力制御部３３は、シミュレーション結果を出力装置３９に出力する。例えば、粒子の位置の変化を、出力装置３９の表示画面に図形で表示する。

【0036】

次に、図６～図８を参照して、上記実施例の優れた効果について説明する。
図６は、上記実施例及び比較例によるシミュレーション方法を用いて求めた水の温度と圧力との関係を示すグラフである。横軸は、温度を単位「Ｋ」で表し、縦軸は圧力を単位「ＭＰａ」で表す。グラフ中の丸記号及び四角記号は、それぞれ上記実施例及び比較例によるシミュレーション方法を用いて求めた水の温度と圧力との関係を示し、破線は、水の物性値としての温度と圧力との関係を示す。

【0037】

比較例では、水分子を従来のＴＩＰ４Ｐモデルで表す。ＴＩＰ４Ｐモデルでは、水分子を、１個の酸素原子、２個の水素原子、及び１個のダミー原子で表す。これらの複数の粒子について運動方程式を解くことにより、水分子の挙動を解析する。計算対象となる粒子の個数が水分子の個数に比べて多くなるため、実施例による方法と比べて計算量が多くなり、計算時間が長大化する。例えば、比較例による方法でシミュレーションを行う場合の計算時間は、実施例による方法でシミュレーションを行う場合の計算時間の約３００倍になる。

【0038】

図６に示すように、実施例による方法、及び比較例による方法のいずれも，実際の水の物性値から求まる圧力を精度よく再現していることがわかる。

【0039】

図７は、実施例及び比較例によるシミュレーション方法を用いて求めた水の飽和蒸気線、及び水の実際の物性値としての飽和蒸気線を示すグラフである。横軸は、密度を単位「ｋｇ／ｍ^３」で表し、縦軸は、温度を単位「Ｋ」で表す。グラフ中の丸記号及び四角記号は、それぞれ実施例及び比較例によるシミュレーション方法を用いて求めた水の飽和蒸気線を示し、破線は、水の物性値としての飽和蒸気線を示す。なお、実施例による方法で求めた飽和蒸気線、及び水の物性値としての飽和蒸気線は、図４に示したグラフと同一である。

【0040】

図７に示すように、実施例による方法で求めた飽和蒸気線は、従来のＴＩＰ４Ｐモデルを用いた場合と比べて、同等またはそれ以上の精度で水の飽和蒸気線を再現していることがわかる。

【0041】

図８は、実施例及び比較例によるシミュレーション方法を用いて求めた水の表面張力、及び水の実際の物性値としての表面張力を示すグラフである。横軸は、温度を単位「Ｋ」で表し、縦軸は、表面張力を単位「Ｎ／ｍ」で表す。グラフ中の丸記号及び四角記号は、それぞれ実施例及び比較例によるシミュレーション方法を用いて求めた水の表面張力を示し、破線は、水の物性値としての表面張力を示す。

【0042】

実施例による方法で求めた表面張力は、従来のＴＩＰ４Ｐモデルを用いた場合と比べて、同等またはそれ以上の精度で水の表面張力を再現していることがわかる。

【0043】

従来の単原子水モデルで用いられる二体間相互作用ポテンシャル及び三体間相互作用ポテンシャルは、粒子間距離が粒子の直径程度である状態（液相状態）に適用できるが、粒子間距離が遠い状態（気相状態）には適用が困難である。これに対して上記実施例では、二体間相互作用ポテンシャルＵ_２としてレナードジョーンズポテンシャルを用いているため、液相と気相とが混在し、相転移が生じ得る系を精度良く再現することができる。このように、相転移を伴う系の解析を行うことが可能である。

【0044】

また、水分子を１つの粒子として解析するため、１つの水分子を複数の粒子で表す従来の水モデルを用いたシミュレーションと比べて、計算量を削減することができる。また、図６～図８に示したように、計算量を削減しても、従来の水モデルを用いる方法と同等の解析精度を維持することができる。

【0045】

次に、上記実施例の変形例について説明する。
上記実施例では、解析対象を水にしているが、二体間相互作用ポテンシャル及び三体間相互作用ポテンシャルのパラメータの値を変更することにより、分子間に異方性ポテンシャルが作用する水以外の分子からなる対象物のシミュレーションを行うことも可能である。

【0046】

水の解析を行う場合のパラメータとして、図４に示したデータを得るときに用いたパラメータの値を含むある範囲内から選択してもよい。例えば、パラメータの値を、
λ≦２３．１５×２、
０．９≦γ≦１．５、
１．７≦ａ≦１．９、
１．０×１０^－２１Ｊ≦ε≦１．０×１０^－１９Ｊ、
１．２５×１０^－１０ｍ≦σ≦３．７５×１０^－１０ｍ
を満たす範囲から選択してもよい。

【0047】

上記実施例は例示であり、本発明は上記実施例に制限されるものではない。例えば、種々の変更、改良、組み合わせ等が可能なことは当業者に自明であろう。

【符号の説明】

【0048】

１０ 解析領域
１１ 液相領域
１２ 気相領域
３０ 処理装置
３１ シミュレーション条件取得部
３２ 演算部
３３ 出力制御部
３８ 入力装置
３９ 出力装置

【図1】

【図2】

【図3】

【図4】

【図5】

【図6】

【図7】

【図8】