(57)【要約】【課題】分子の挙動を適切にかつ、効率的にシミュレーションすることができるシミュレーション方法及びその装置、プログラムを提供する。【解決手段】シミュレーションボックス中の分子集合状態の異方性を反映する異方性因子を変数とする、前記シミュレーションボックスの形状及び挙動を記述する式を設定する。設定された式を下に静水圧条件及びビリアルの定理を満足する演算式を生成する。生成された演算式を記憶媒体に記憶する。記憶された演算式を用いて、シミュレーションボックスの形状を制御し、該シミュレーションボックス及び分子の挙動を演算する。この演算によって、バルク状態の分子挙動をシミュレーションし、その演算結果を表示部に出力する。【選択図】  図３
2004192275.tif 000002
【特許請求の範囲】
【請求項１】
少なくとも１方向に周期境界条件を課したシミュレーションボックスを用いて、有限数の分子系でバルク状態の分子挙動をシミュレーションするシミュレーション方法であって、
シミュレーションボックス中の分子集合状態の異方性を反映する異方性因子を変数とする、前記シミュレーションボックスの形状及び挙動を記述する式を設定する設定工程と、
前記設定工程で設定された式を下に静水圧条件及びビリアルの定理を満足する演算式を生成する生成工程と、
前記生成工程で生成された演算式を記憶媒体に記憶する記憶工程と、
前記記憶工程で前記記憶媒体に記憶された演算式を用いて、前記シミュレーションボックスの形状を制御し、該シミュレーションボックス及び分子の挙動を演算する演算工程と、
前記演算工程による演算結果を表示部に出力する出力工程とを備え、
前記演算工程による演算によって、前記バルク状態の分子挙動をシミュレーションする
ことを特徴とするシミュレーション方法。
【請求項２】
前記設定工程は、前記シミュレーションボックスの挙動を、運動エネルギーで表現する表現式として設定する
ことを特徴とする請求項１に記載のシミュレーション方法。
【請求項３】
前記設定工程は、前記異方性因子及び該異方性因子との積で前記シミュレーションボックスの体積を示す等方体積因子を変数とする、前記シミュレーションボックスの形状を記述する式と該シミュレーションボックスの挙動を運動エネルギーで表現する表現式を設定する
ことを特徴とする請求項１に記載のシミュレーション方法。
【請求項４】
前記異方性因子をα、前記等方体積因子をＱ、
それぞれの速度を
2004192275.tif 000045
前記シミュレーションボックスの揺らぎを調整するためのパラメータをＭと定義した場合、
前記表現式とするＫｃｅｌｌは、以下のように表現される
2004192275.tif 000046
ことを特徴とする請求項３に記載のシミュレーション方法。
【請求項５】
圧力をＰｅｘ
前記シミュレーションボックス中の各分子の座標をｒｉ
前記シミュレーションボックス中の各分子のｘｙｚ空間におけるスケール座標をｓｉｘ，ｓｉｙ，ｓｉｚ
その速度を
2004192275.tif 000047
と定義した場合、
ｒｉ＝Ｑ１／３ｓｉｘ＋Ｑ１／３ｓｉｙ＋αＱ１／３ｓｉｚ
と表現され、かつ
前記シミュレーションボックスの揺らぎを調整するためのパラメータをＭ、
注目分子ｉに対する前記周期境界条件を満たす他分子ｊとの距離のスケール座標をｓｉｊ、
分子間相互作用ポテンシャルをφ、分子の質量をｍと定義した場合、
前記生成工程で生成する演算式であるラグランジアンは、以下のように表現される
2004192275.tif 000048
ことを特徴とする請求項４に記載のシミュレーション方法。
【請求項６】
前記異方性因子をα、
その速度を
2004192275.tif 000049
前記シミュレーションボックスの揺らぎを調整するためのパラメータをＭと定義した場合、
前記表現式とするＫｃｅｌｌは、以下のように表現される
2004192275.tif 000050
ことを特徴とする請求項２に記載のシミュレーション方法。
【請求項７】
前記シミュレーションボックスの体積をＶ、
前記シミュレーションボックス中の各分子の座標をｒｉ、
前記シミュレーションボックス中の各分子のｘｙｚ空間におけるスケール座標をｓｉｘ，ｓｉｙ，ｓｉｚ
その速度を
2004192275.tif 000051
と定義した場合、
2004192275.tif 000052
と表現され、かつ
前記シミュレーションボックスの揺らぎを調整するためのパラメータをＭ、
注目分子ｉに対する前記周期境界条件を満たす他分子ｊとの距離のスケール座標をｓｉｊ、
分子間相互作用ポテンシャルをφ、分子の質量をｍと定義した場合、
前記生成工程で生成する演算式であるラグランジアンは、以下のように表現される
2004192275.tif 000053
ことを特徴とする請求項６に記載のシミュレーション方法。
【請求項８】
前記異方性因子をα、
その速度を
2004192275.tif 000054
前記シミュレーションボックスの揺らぎを調整するためのパラメータをＭと定義した場合、
前記表現式とするＫｃｅｌｌは、以下のように表現される
2004192275.tif 000055
ことを特徴とする請求項２に記載のシミュレーション方法。
【請求項９】
前記シミュレーションボックスの体積をＶ、
前記シミュレーションボックス中の各分子の座標をｒｉ、
前記シミュレーションボックス中の各分子のｘｙｚ空間におけるスケール座標をｓｉｘ，ｓｉｙ，ｓｉｚ
その速度を
2004192275.tif 000056
と定義した場合、
2004192275.tif 000057
と表現され、かつ
前記シミュレーションボックスの揺らぎを調整するためのパラメータをＭ、
注目分子ｉに対する前記周期境界条件を満たす他分子ｊとの距離のスケール座標をｓｉｊ、
分子間相互作用ポテンシャルをφ、分子の質量をｍと定義した場合、
前記生成工程で生成する演算式であるラグランジアンは、以下のように表現される
2004192275.tif 000058
ことを特徴とする請求項８に記載のシミュレーション方法。
【請求項１０】
前記分子は、異方性液体状態を有する分子である
ことを特徴とする請求項１に記載のシミュレーション方法。
【請求項１１】
前記シミュレーションは、前記異方性液体状態間または前記異方性液体状態を含む相転移をシミュレーションする
ことを特徴とする請求項１０に記載のシミュレーション方法。
【請求項１２】
少なくとも１方向に周期境界条件を課したシミュレーションボックスを用いて、有限数の分子系でバルク状態の分子挙動をシミュレーションするシミュレーション装置であって、
シミュレーションボックス中の分子集合状態の異方性を反映する異方性因子を変数とする、前記シミュレーションボックスの形状及び挙動を記述する式を設定する設定手段と、
前記設定手段で設定された式を下に静水圧条件及びビリアルの定理を満足する演算式を生成する生成手段と、
前記生成手段で生成された演算式を記憶する記憶手段と、
前記記憶手段で記憶された演算式を用いて、前記シミュレーションボックスの形状を制御し、該シミュレーションボックス及び分子の挙動を演算する演算手段と、
前記演算手段による演算結果を表示部に出力する出力手段とを備え、
前記演算手段による演算によって、前記バルク状態の分子挙動をシミュレーションする
ことを特徴とするシミュレーション装置。
【請求項１３】
少なくとも１方向に周期境界条件を課したシミュレーションボックスを用いて、有限数の分子系でバルク状態の分子挙動をシミュレーションするシミュレーションを実現するプログラムであって、
シミュレーションボックス中の分子集合状態の異方性を反映する異方性因子を変数とする前記シミュレーションボックスの形状及び挙動を記述する式を下に生成した、静水圧条件及びビリアルの定理を満足する演算式をプログラミングするプログラミング工程のプログラムコードと、
前記生成工程で生成された演算式を記憶媒体に記憶する記憶工程のプログラムコードと、
前記記憶工程で前記記憶媒体に記憶された演算式を用いて、前記シミュレーションボックスの形状を制御し、該シミュレーションボックス及び分子の挙動を演算する演算工程のプログラムコードと、
前記演算工程による演算結果を表示部に出力する出力工程のプログラムコードとを備え、
前記演算工程による演算によって、前記バルク状態の分子挙動をシミュレーションする
ことを特徴とするプログラム。
【発明の詳細な説明】
【０００１】
【発明の属する技術分野】
本発明は、少なくとも１方向に周期境界条件を課したシミュレーションボックスを用いて、有限数の分子系でバルク状態の分子挙動をシミュレーションするシミュレーション方法及びその装置、プログラムに関するものである。
【０００２】
【従来の技術】
環境情報（圧力・温度など）をも含めた分子（あるいは粒子）集合の挙動・状態をシミュレーションする方法は、すでにいくつか開発されている。
【０００３】
その中でも、圧力が一定の環境下における分子（粒子）挙動をシミュレーションする方法は、主に以下の２つに分けられる。
【０００４】
（ａ）等方性液体のために開発された方法（Ａｎｄｅｒｓｅｎ法）
（ｂ）結晶のために開発された方法（Ｐａｒｒｉｎｅｌｌｏ  ＆  Ｒａｈｍａｎ法）
方法（ａ）は、液体が等方的である場合に適用され、方法（ｂ）は、弾性体である結晶または固体のシミュレーションに適する。
【０００５】
【発明が解決しようとする課題】
しかしながら、異方性液体は粘弾性など多くの物性に異方性を持つ液体であるため、別の目的のために開発された上記方法を用いて適切なシミュレーションを行うことは難しい。
【０００６】
例えば、上記の方法（ａ）を異方性液体に適用すると、系が静水圧にならない等、本来シミュレーションしたい環境ではなくなってしまう場合が多い。また、方法（ｂ）を異方性液体に適用すると、異方性液体の特徴である粘弾性の異方性により適切な平衡状態が得られない等の問題点が生じる。
【０００７】
そのため、異方性液体の適切な環境下における挙動を適切な数の分子（粒子）のシミュレーションで知ることは重要な課題となっている。
【０００８】
本発明は上記の課題に鑑みてなされたものであり、分子の挙動を適切にかつ、効率的にシミュレーションすることができるシミュレーション方法及びその装置、プログラムを提供することを目的とする。
【０００９】
【課題を解決するための手段】
上記の目的を達成するための本発明によるシミュレーション方法は以下の構成を備える。即ち、
少なくとも１方向に周期境界条件を課したシミュレーションボックスを用いて、有限数の分子系でバルク状態の分子挙動をシミュレーションするシミュレーション方法であって、
シミュレーションボックス中の分子集合状態の異方性を反映する異方性因子を変数とする、前記シミュレーションボックスの形状及び挙動を記述する式を設定する設定工程と、
前記設定工程で設定された式を下に静水圧条件及びビリアルの定理を満足する演算式を生成する生成工程と、
前記生成工程で生成された演算式を記憶媒体に記憶する記憶工程と、
前記記憶工程で前記記憶媒体に記憶された演算式を用いて、前記シミュレーションボックスの形状を制御し、該シミュレーションボックス及び分子の挙動を演算する演算工程と、
前記演算工程による演算結果を表示部に出力する出力工程とを備え、
前記演算工程による演算によって、前記バルク状態の分子挙動をシミュレーションする。
【００１０】
また、好ましくは、前記設定工程は、前記シミュレーションボックスの挙動を、運動エネルギーで表現する表現式として設定する。
【００１１】
また、好ましくは、前記設定工程は、前記異方性因子及び該異方性因子との積で前記シミュレーションボックスの体積を示す等方体積因子を変数とする、前記シミュレーションボックスの形状を記述する式と該シミュレーションボックスの挙動を運動エネルギーで表現する表現式を設定する。
【００１２】
また、好ましくは、前記異方性因子をα、前記等方体積因子をＱ、
それぞれの速度を
2004192275.tif 000003
前記シミュレーションボックスの揺らぎを調整するためのパラメータをＭと定義した場合、
前記表現式とするＫｃｅｌｌは、以下のように表現される
2004192275.tif 000004
ことを特徴とする。
【００１３】
また、好ましくは、圧力をＰｅｘ
前記シミュレーションボックス中の各分子の座標をｒｉ
前記シミュレーションボックス中の各分子のｘｙｚ空間におけるスケール座標をｓｉｘ，ｓｉｙ，ｓｉｚ
その速度を
2004192275.tif 000005
と定義した場合、
ｒｉ＝Ｑ１／３ｓｉｘ＋Ｑ１／３ｓｉｙ＋αＱ１／３ｓｉｚ
と表現され、かつ
前記シミュレーションボックスの揺らぎを調整するためのパラメータをＭ、
注目分子ｉに対する前記周期境界条件を満たす他分子ｊとの距離のスケール座標をｓｉｊ、
分子間相互作用ポテンシャルをφ、分子の質量をｍと定義した場合、
前記生成工程で生成する演算式であるラグランジアンは、以下のように表現される
2004192275.tif 000006
ことを特徴とする。
【００１４】
また、好ましくは、前記異方性因子をα、
その速度を
2004192275.tif 000007
前記シミュレーションボックスの揺らぎを調整するためのパラメータをＭと定義した場合、
前記表現式とするＫｃｅｌｌは、以下のように表現される
2004192275.tif 000008
ことを特徴とする。
【００１５】
また、好ましくは、前記シミュレーションボックスの体積をＶ、
前記シミュレーションボックス中の各分子の座標をｒｉ、
前記シミュレーションボックス中の各分子のｘｙｚ空間におけるスケール座標をｓｉｘ，ｓｉｙ，ｓｉｚ
その速度を
2004192275.tif 000009
と定義した場合、
2004192275.tif 000010
と表現され、かつ
前記シミュレーションボックスの揺らぎを調整するためのパラメータをＭ、
注目分子ｉに対する前記周期境界条件を満たす他分子ｊとの距離のスケール座標をｓｉｊ、
分子間相互作用ポテンシャルをφ、分子の質量をｍと定義した場合、
前記生成工程で生成する演算式であるラグランジアンは、以下のように表現される
2004192275.tif 000011
ことを特徴とする。
【００１６】
また、好ましくは、前記異方性因子をα、
その速度を
2004192275.tif 000012
前記シミュレーションボックスの揺らぎを調整するためのパラメータをＭと定義した場合、
前記表現式とするＫｃｅｌｌは、以下のように表現される
2004192275.tif 000013
ことを特徴とする。
【００１７】
また、好ましくは、前記シミュレーションボックスの体積をＶ、
前記シミュレーションボックス中の各分子の座標をｒｉ、
前記シミュレーションボックス中の各分子のｘｙｚ空間におけるスケール座標をｓｉｘ，ｓｉｙ，ｓｉｚ
その速度を
2004192275.tif 000014
と定義した場合、
2004192275.tif 000015
と表現され、かつ
前記シミュレーションボックスの揺らぎを調整するためのパラメータをＭ、
注目分子ｉに対する前記周期境界条件を満たす他分子ｊとの距離のスケール座標をｓｉｊ、
分子間相互作用ポテンシャルをφ、分子の質量をｍと定義した場合、
前記生成工程で生成する演算式であるラグランジアンは、以下のように表現される
2004192275.tif 000016
ことを特徴とする。
【００１８】
また、好ましくは、前記分子は、異方性液体状態を有する分子である。
【００１９】
また、好ましくは、前記シミュレーションは、前記異方性液体状態間または前記異方性液体状態を含む相転移をシミュレーションする。
【００２０】
上記の目的を達成するための本発明によるシミュレーション装置は以下の構成を備える。即ち、
少なくとも１方向に周期境界条件を課したシミュレーションボックスを用いて、有限数の分子系でバルク状態の分子挙動をシミュレーションするシミュレーション装置であって、
シミュレーションボックス中の分子集合状態の異方性を反映する異方性因子を変数とする、前記シミュレーションボックスの形状及び挙動を記述する式を設定する設定手段と、
前記設定手段で設定された式を下に静水圧条件及びビリアルの定理を満足する演算式を生成する生成手段と、
前記生成手段で生成された演算式を記憶する記憶手段と、
前記記憶手段で記憶された演算式を用いて、前記シミュレーションボックスの形状を制御し、該シミュレーションボックス及び分子の挙動を演算する演算手段と、
前記演算手段による演算結果を表示部に出力する出力手段とを備え、
前記演算手段による演算によって、前記バルク状態の分子挙動をシミュレーションする。
【００２１】
上記の目的を達成するための本発明によるプログラムは以下の構成を備える。即ち、
少なくとも１方向に周期境界条件を課したシミュレーションボックスを用いて、有限数の分子系でバルク状態の分子挙動をシミュレーションするシミュレーションを実現するプログラムであって、
シミュレーションボックス中の分子集合状態の異方性を反映する異方性因子を変数とする前記シミュレーションボックスの形状及び挙動を記述する式を下に生成した、静水圧条件及びビリアルの定理を満足する演算式をプログラミングするプログラミング工程のプログラムコードと、
前記生成工程で生成された演算式を記憶媒体に記憶する記憶工程のプログラムコードと、
前記記憶工程で前記記憶媒体に記憶された演算式を用いて、前記シミュレーションボックスの形状を制御し、該シミュレーションボックス及び分子の挙動を演算する演算工程のプログラムコードと、
前記演算工程による演算結果を表示部に出力する出力工程のプログラムコードとを備え、
前記演算工程による演算によって、前記バルク状態の分子挙動をシミュレーションする。
【００２２】
【発明の実施の形態】
以下、本発明の実施の形態について図面を用いて詳細に説明する。
【００２３】
本発明は、分子、特に、異方性液体（液晶、生体膜等）状態を有する分子の挙動を効率良く正確にシミュレーションすることを実現するものである。但し、本発明は、異方性液体に限定されず、様々な分子の挙動のシミュレーションを実現することができる。
【００２４】
実験で扱う分子数１０１６から１０２３程度）からなる物質の挙動を再現し、その物質の諸物性のシミュレーションをする場合、そのシミュレーションを実行するためのハードウエアの処理能力の限界から、実際の系より遙かに少ない分子（粒子）数しか扱えない。
【００２５】
そのため、本実施形態で、限られた分子数でなるべく現実の系の状態を精度良く知るために、分子（粒子）集合体（バルク状態の分子）の環境を表現するための周期境界条件を課した仮想的なシミュレーションボックスを用いてシミュレーションを行う。
【００２６】
この周期境界条件を課したシミュレーションボックスは、バルク状態の分子の一部分（シミュレーション対象部分）の周辺の環境情報を代表するものであり、そのシミュレーションボックスでシミュレーションすることによって、バルク状態の分子の一部分をシミュレーションしているにもかかわらずバルク状態の分子全体の性質・挙動を良い精度で再現することができる。
【００２７】
シミュレーションボックスには、シミュレーションしたい系の状況に合わせて周期境界条件を１から３方向に設定する。
【００２８】
例えば、バルク状態の分子の中央付近の振る舞いが知りたければ３方向、界面の振る舞いが知りたければ界面と平行方向の２方向、セル中などに閉じ込められた系なら開放されている方向のみに周期境界条件を設定する。
【００２９】
このような状況に合わせて、本実施形態では、図１に示すように、バルク状態の分子の分子数より少ない所定分子数からなるシミュレーションボックスを定義し、そのシミュレーションボックスと中の分子（粒子）集合を１つの系として、バルク状態の分子の挙動のシミュレーションを行う。
【００３０】
その際、系を記述する基本方程式であるラグランジアンは、シミュレーションボックスの動的挙動を含む方程式になっている。
【００３１】
以下、シミュレーションの対象となる異方性液体について説明した後、シミュレーションの具体的な方法を
（Ａ）圧力が一定の場合
（Ｂ）体積一定の場合
のそれぞれの場合について説明を行う。
【００３２】
異方性液体には、例えば、液晶が挙げられる。
【００３３】
液晶状態は、固体や液体のように物質の状態の一つである。また、液晶状態は、固体あるいは結晶の温度を上昇または濃度を下げる、あるいは液体の温度を下降または濃度を下げることによって、相転移（結晶（Ｃ）−スメクティック液晶（Ｓｍ）−ネマティック液晶（Ｎ）−等方性液体（Ｉ））を経て出現する。
【００３４】
そして、この液晶状態が出現する際の温度又は濃度、相転移の様子は重要な情報である。加えて、液晶状態になっている温度域・濃度域でどのような物性（配向状態、拡散係数、各種粘弾性など）を示すかを得ることも、例えば、液晶ディスプレイの設計等のデバイスの製造をする場合に重要な情報となる。
【００３５】
液晶相の中には大きく分けると層構造をもつスメクティック相（Ｓｍ相）（あるいはラメラ相）と層構造をもたないネマティック相（Ｎ相）がある。どちらの相も分子の配向が相関をもち、それによって様々な物性に異方性が現れる。
【００３６】
液晶状態になる物質には、低分子棒状分子、円盤状分子、高分子などの多岐にわたり、生体膜も液晶状態をとる。
【００３７】
従って、異方性液体の有効なシミュレーションは、液晶デバイスの開発にとって不可欠なだけではなく、高分子や生体現象の理解・応用にも欠かせない。
【００３８】
これら異方性液体の一部分（シミュレーション対象分子）を切り出してくると圧力一定の環境の下では、それらの異方性液体特有のシミュレーションボックスの形をとる。
【００３９】
また、少しの温度あるいは濃度変化でその形が変化するという共通の特徴がある。
【００４０】
体積一定の場合でも比較的少ない分子数（粒子数）で分子バルクの状態を実現するには、分子（粒子）状態の異方性に応じたシミュレーションボックスを用いる必要がある。
【００４１】
そこで、本実施形態において、方法（Ａ）と（Ｂ）に共通することは、シミュレーションボックスの動力学を記述するのに異方性因子αを導入している点である。
【００４２】
方法（Ａ）では、圧力を一定に保つようにシミュレーションボックスの形状および体積が揺らぐ。
【００４３】
シミュレーションボックスの形状は、立方体と比べてどの程度、異方的であるか（つまり、シミュレーションボックスをｘｙｚの３次元空間で定義した場合、そのｘｙｚ方向に対してどのように変形するか）を示す異方性因子αによって表現される。そして、本実施形態では、この異方性因子αと等方体積因子Ｑの積、即ち、Ｖ＝αＱがシミュレーションボックスの体積であり、かつその形状を定義する。
【００４４】
尚、この異方性因子αは、シミュレーションボックスあるいはその分子が存在する３次元空間を規定するｘｙｚ空間の形状を制御する時間を変数とする時間関数Ｆ（ｘ（ｔ），  ｙ（ｔ），  ｚ（ｔ））として表現されるものである。そのため、体積Ｖは、時間の経過に伴って変化する。これによって、時系列で体積と形状が変化するシミュレーションボックスを用いたシミュレーションを実行することができる。
【００４５】
シミュレーションボックスの運動エネルギーＫｃｅｌｌは、異方性因子α、等方体積因子Ｑそれぞれの変数の時間微分、即ち、速度を
2004192275.tif 000017
と記述すると、
【００４６】
【数１】
2004192275.tif 000018
【００４７】
となる表現式を生成する。
【００４８】
シミュレーションボックス中の各分子の座標をｒｉ、シミュレーションボックス中の各分子のｘｙｚ空間におけるスケール座標をｓｉｘ，ｓｉｙ，ｓｉｚとすると、各分子の座標ｒｉとスケール座標ｓｉは、
ｒｉ＝Ｑ１／３ｓｉｘ＋Ｑ１／３ｓｉｙ＋αＱ１／３ｓｉｚ
で関係づけられる。
【００４９】
また、スケール座標ｓｉｘ，ｓｉｙ，ｓｉｚの速度を、
2004192275.tif 000019
と記述する。
【００５０】
これにより、系の基本方程式であるラグランジアンは、
【００５１】
【数２】
2004192275.tif 000020
【００５２】
となる。
【００５３】
ここで、Ｓｉｊは、注目分子ｉに対する前記周期境界条件を満たす他分子ｊのスケール座標を示している。
【００５４】
ここで、φは分子（粒子）間相互作用ポテンシャルである。また、ｍは、分子または粒子の質量を示し、Ｐｅｘは、設定した圧力を示している。
【００５５】
更に、Ｎ、ｉ、ｊはそれぞれ自然数である。
【００５６】
これをラグランジュ方程式に代入して解くと、以下の体積Ｖの方程式を得る。
【００５７】
【数３】
2004192275.tif 000021
【００５８】
式（３）より、平衡状態
2004192275.tif 000022
では、ビリアルの定理である以下の
2004192275.tif 000023
と、静水圧である以下の
2004192275.tif 000024
が自動的に保証されることが解る（ここで、Ｐｘｙ＝（Ｐｘ＋Ｐｙ）／２）。
【００５９】
つまり、本実施形態のシミュレーションで使用する（２）式によれば、常に、静水圧条件およびビリアルの定理を満足することができる。これは、バルク状態の分子の分子数よりも少ない所定分子数のシミュレーションボックスに対するシミュレーションを実行するだけで、そのバルク状態の分子の挙動のシミュレーションを実現できることを意味する。
【００６０】
具体的に、コンピュータ等の計算機で解く方程式を得るには、例えば、次のハミルトニアンを用いる。
【００６１】
【数４】
2004192275.tif 000025
【００６２】
ここで、
2004192275.tif 000026
と、
2004192275.tif 000027
は、それぞれαとＱの共役運動量である。
【００６３】
ｓｉの共役運動量は、
2004192275.tif 000028
である。
【００６４】
最終的に解くべき式は、
【００６５】
【数５】
2004192275.tif 000029
【００６６】
【数６】
2004192275.tif 000030
【００６７】
【数７】
2004192275.tif 000031
【００６８】
【数８】
2004192275.tif 000032
【００６９】
である。
【００７０】
ここで、
2004192275.tif 000033
である。
【００７１】
一方、方法（Ｂ）では、設定した体積Ｖの下で異方性因子αが揺らぐ。シミュレーションボックスの運動エネルギーＫｃｅｌｌは、αの速度を
2004192275.tif 000034
と記述すると、
【００７２】
【数９】
2004192275.tif 000035
【００７３】
とおける。
【００７４】
シミュレーションボックス中の各粒子の座標ｒｉとスケール座標ｓｉは、
【００７５】
【数１０】
2004192275.tif 000036
【００７６】
で関係づけられる。
【００７７】
ラグランジアンは、
【００７８】
【数１１】
2004192275.tif 000037
【００７９】
となる。ここで、φは分子（粒子）間相互作用ポテンシャル、ｍは分子または粒子の質量である。そして、最終的に、コンピュータ等の計算機で解くべき式は、
【００８０】
【数１２】
2004192275.tif 000038
【００８１】
【数１３】
2004192275.tif 000039
【００８２】
【数１４】
2004192275.tif 000040
【００８３】
となる。ここで、式（１４）中のｆｉｊは、
2004192275.tif 000041
である。式（１４）を用いた場合、αはＭαの値に応じた周期でゆらぐ。以下のＫｃｅｌｌを用いると、Ｍの値のみ応じた周期でゆらぐ式を導き出せる。
【００８４】
尚、体積一定の場合は、変数がαだけなので式（１４）は、圧力一定の場合の式（３）に相当する。式（３）で平衡状態では、＜ｄＶ２／ｄ２ｔ＞＝０なので、Ｍの値によらずにビリアルの定理と静水圧は保証される。
【００８５】
同様に、式（１４）で、平衡状態で＜ｄα２／ｄ２ｔ＞＝０であり、Ｍαの値によらず、式（１４）中の項｛．．．｝＝０が保証される。ここで、式（１４）中の項｛．．．｝＝０がどのように揺らぐかを決めるのは、Ｍαである。
【００８６】
従って、揺らぎ方はＭの値だけでなく、その時々のαの値にも依存してしまう。すると、同じ現象をシミュレーションし、しかも同じＭの値を用いているのに、シミュレーションボックスの形状が異なると、異なるシミュレーション結果を得ることになる。例えば、Ｓｍ相を３層シミュレーションしたときと、６層シミュレーションしたときで同じ条件下でも異なるシミュレーション結果が得られてしまう。
【００８７】
従って、この点を解決したのが、式（１５）、（１６）を使用して導出した式（１７）である。
【００８８】
【数１５】
2004192275.tif 000042
【００８９】
ラグランジアンは、
【００９０】
【数１６】
2004192275.tif 000043
【００９１】
となる。よって、シミュレーションボックスの運動を記述する最終的な式は、
【００９２】
【数１７】
2004192275.tif 000044
【００９３】
となる。この方法でも分子の運動は式（１２）と式（１３）で記述される。
【００９４】
尚、式（１１）と式（１６）とでは、式（１０）を共通して使用していることに注意されたい。
【００９５】
次に、上記シミュレーションを実行するシミュレーション装置のハードウエア構成について、図２を用いて説明する。
【００９６】
図２は本発明の実施形態のシミュレーション装置のハードウエア構成を示す図である。
【００９７】
図２において、１０１はＣＰＵであり、シミュレーション装置の各種構成要素を制御する。１０２はＲＡＭであり、データの作業領域及び一時待避領域として機能する。１０３はＲＯＭであり、後述する図３のフローチャートを含む本実施形態のシミュレーションによる各種処理を実現するプログラムや、シミュレーションを実行するための各種テーブル、データ群、方程式群を記憶する。
【００９８】
１０４はキーボードやマウスからなる入力部であり、このキーボードやマウスを介して、シミュレーションを実行するための各種パラメータの入力や、シミュレーションの実行指示等の各種入力を行う。１０５はディスプレイであり、ＬＣＤやＣＲＴから構成され、本実施形態のシミュレーションに対する各種パラメータの入力／設定画面、シミュレーションの実行を指示するための指示画面、シミュレーション結果の出力画面、ステータス表示画面等の各種表示画面を表示する。
【００９９】
１０６は外部メモリであり、ハードディスクや、ＣＤ−ＲＯＭ、ＣＤ−Ｒ／ＲＷ、ＤＶＤ−ＲＯＭ、ＤＶＤ−ＲＡＭ、ＭＯ等から構成され、上記のプログラムやデータを記憶したり、本実施形態のシミュレーションプログラムを制御するためのオペレーティングシステムを記憶している。
【０１００】
１０７はネットワークインタフェースであり、ネットワーク（不図示）を介して外部機器とのデータの送受信を制御する。例えば、シミュレーションプログラムを搭載するシミュレーションサーバをネットワーク上に構成し、本シミュレーション装置から各種パラメータをシミュレーションサーバへ送信し、そのシミュレーションサーバによるシミュレーション結果をシミュレーション装置で受信するようにすることもできる。
【０１０１】
１０８は本シミュレーション装置の各種構成要素を相互に接続するバスである。
【０１０２】
尚、本実施形態のシミュレーションは、汎用コンピュータで動作可能なソフトウェアや、そのソフトウェアを記憶するファームウェアで実現されたり、専用ハードウェアとして実現されても良い。
【０１０３】
次に、本実施形態のシミュレーションの処理フローについて、図３を用いて説明する。
【０１０４】
図３は本発明の実施形態のシミュレーションの処理フローを示すフローチャートである。
【０１０５】
尚、図３では、汎用コンピュータで動作可能なソフトウェアを、ＣＰＵ１０１によって実行することによって、上記シミュレーションを実行するシミュレーション装置の処理フローとして説明する。また、図３では、（Ａ）圧力一定の場合でシミュレーションを実行する例について説明する。
【０１０６】
ステップＳ１０１で、シミュレーションボックス中の分子集合状態の異方性を反映する異方性因子及びその異方性因子との積でシミュレーションボックスの体積を示す等方体積因子を変数とする、シミュレーションボックスの形状及び挙動を記述する式を設定する。また、この設定された式を下に静水圧条件及びビリアルの定理を満足する演算式を生成する。そして、この演算式を、ＲＡＭ１０２に一時的に記憶する。
【０１０７】
次に、ステップＳ１０２で、入力部１０４より、この演算式を用いてシミュレーションを実行するために、シミュレーションボックス及び分子の初期状態（配置、速度、数、密度等）を設定する。
【０１０８】
次に、ステップＳ１０３で、ＣＰＵ１０１は、ＲＡＭ１０２に記憶される演算式を用いて、シミュレーションボックスの形状を制御し、シミュレーションボックス及び分子の挙動を演算する。
【０１０９】
ステップＳ１０１〜ステップＳ１０３の処理を具体的に説明すると、異方性因子αと等方体積因子Ｑを変数とする、シミュレーションボックス形状を記述する式（Ｖ＝αＱ）と、そのシミュレーションボックスの挙動を運動エネルギーで表現する表現式（式（１））を設定する。次に、この表現式をラグランジアンに導入することで、静水圧条件及びビリアルの定理を保証する演算式（２）を生成することができる。そして、この演算式（２）から導出される演算式群（式（４）〜式（８））を用いて、設定された初期状態からのシミュレーションボックス及びその中の分子の挙動（運動量）を演算する（バルク状態の分子の挙動を算出する（シミュレーションする））ことができる。
【０１１０】
そして、ステップＳ１０４で、ＣＰＵ１０１は、その演算結果をディスプレイ１０５に出力する。
【０１１１】
次に、本実施形態のシミュレーションによるシミュレーション結果の具体例について、以下に説明する。
【０１１２】
ここでは、軟らかい冠球円柱粒子系の結晶（Ｃ）−スメクティック液晶（Ｓｍ）−ネマティック液晶（Ｎ）−等方性液体（Ｉ）の相転移のシミュレーション結果について説明する。
【０１１３】
図４は、換算温度Ｔ＊に対する１粒子あたりの体積Ｖ／Ｎを示すグラグである。この図４より、各相転移温度（結晶（Ｃ）→スメクティック液晶（Ｓｍ）：６２℃、スメクティック液晶（Ｓｍ）→ネマティック液晶（Ｎ）：８０℃、ネマティック液晶（Ｎ）→等方性液体（Ｉ）：１１８℃）、各温度における密度、相転移の次数等の情報が得られる。
【０１１４】
次に、同じ系における換算温度Ｔ＊に対する配向秩序のパラメータ（Ｏｒｄｅｒ  ｐａｒａｍｅｔｅｒ）Ｓ２（○と●）と、σ（◇と◆）を図５に示す。  尚、図４及び図５の（○と●）と、σ（◇と●）は、それぞれ汎用的なアルゴリズムである
○、◇：ベルレ速度法  （Ｖｅｒｌｅｔ  ｖｅｌｏｃｉｔｙ  ａｌｇｏｒｉｔｈｍ）
●、◆：予測子修正子法（Ｇｅａｒ  ｐｒｅｄｉｃｔｏｒ−ｃｏｒｒｅｃｔｏｒ  ａｌｇｏｒｉｔｈｍ）
を用いて、上記シミュレーションを実行した場合のシミュレーション結果を示している。
【０１１５】
尚、ここでは、図４や図５には示していないが、従来のシミュレーション方法で上記アルゴリズムを使用したシミュレーションを行った場合は、その結果は、図４や図５に示すような結果は得られず、本シミュレーションが分子の挙動を知る上で有効であることが理解されるであろう。
【０１１６】
ここで、
Ｓ２＝＜（３ｃｏｓ２θｉｊ−１）／２＞
と、
σ＝＜ｃｏｓ（２πｚｉｊ／ｄ）（３ｃｏｓ２θｉｊ−１）／２＞
である。
【０１１７】
ここで、ｚｉｊは、層と垂直方向の距離を示している。
【０１１８】
これにより、各温度で分子の配向の度合いが分かる。
【０１１９】
図６は、Ｓｍ相（Ｔ＊＝７０℃）とＮ相（Ｔ＊＝１００℃）における分子配向方向（点線）および垂直方向（実線）の平均二乗変位（Ｍｅａｎ  Ｓｅｑｕａｒｅｄ  Ｄｉｓｐｌａｃｅｍｅｎｔ）を示している。これにより、分子の拡散の様子および拡散係数が解る。
【０１２０】
図７は、Ｎ相（Ｔ＊＝１００℃）におけるｘ軸（実線）、ｙ軸（点線）、ｚ軸（波線）方向の角度の変化の二乗である平均二乗角度変位（Ｍｅａｎ  Ｓｅｑｕａｒｅｄ  Ａｎｇｕｌａｒ  Ｄｉｓｐｌａｃｅｍｅｎｔ）を示している。ここでは、分子配向方向はｚ軸と平行である。これにより、回転粘度が解る。
【０１２１】
また、図６及び図７には、図４や図５と同様に、従来のシミュレーション方法によるシミュレーション結果は示していないが、これらの図６及び図７に示すような本シミュレーションによるシミュレーション結果と同様のシミュレーション結果が、従来のシミュレーション方法からは結晶相しか精度良く得られないことがシミュレーションにより明らかになっている。
【０１２２】
以上説明したように、本実施形態によれば、バルク状態の分子の環境を加味した上で、バルク状態の分子の分子数よりも少ない分子（粒子）数で静水圧条件およびビリアルの定理を満足する分子の挙動のシミュレーションを実行することができるので、精度良く分子バルクの挙動・物性をシミュレートすることができる。
【０１２３】
また、従来のシミュレーション方法に比べて、より少ない分子数でのシミュレーションでも、精度の高いシミュレーション結果を得ることができる。つまり、演算量を削減し、より効率的にシミュレーションを実行することができる。
【０１２４】
また、シミュレーションボックスに対して、異方性因子αを導入することで、シミュレーションボックスの形状や体積を可変に制御することができるので、従来の方法でシミュレーションが不可能だった系にも適応が可能である。例えば、異方性液体だけでなく結晶、等方性液体の挙動もシミュレーションでき、より詳しくは、結晶（Ｃ）−スメクティック液晶（Ｓｍ）−ネマティック液晶（Ｎ）−等方性液体（Ｉ）の一連の相転移もシミュレーションすることができる。
【０１２５】
尚、上記例は分子動力学法に適応した場合であるが、モンテカルロ法にも適応可能である。
【０１２６】
【発明の効果】
以上説明したように、本発明によれば、分子の挙動を適切にかつ、効率的にシミュレーションすることができるシミュレーション方法及びその装置、プログラムを提供できる。
【図面の簡単な説明】
【図１】本発明の実施形態のシミュレーションを説明するための図である。
【図２】本発明の実施形態のシミュレーション装置のハードウエア構成を示す図である。
【図３】本発明の実施形態のシミュレーションの処理フローを示すフローチャートである。
【図４】本発明の実施形態のシミュレーション結果の一例を示す図である。
【図５】本発明の実施形態のシミュレーション結果の一例を示す図である。
【図６】本発明の実施形態のシミュレーション結果の一例を示す図である。
【図７】本発明の実施形態のシミュレーション結果の一例を示す図である。
【符号の説明】
１０１  ＣＰＵ
１０２  ＲＡＭ
１０３  ＲＯＭ
１０４  入力部
１０５  ディスプレイ
１０６  外部メモリ
１０７  ネットワークインタフェース
１０８  バス
【図１】
2004192275.tif 000059
【図２】
2004192275.tif 000060
【図３】
2004192275.tif 000061
【図４】
2004192275.tif 000062
【図５】
2004192275.tif 000063
【図６】
2004192275.tif 000064
【図７】
2004192275.tif 000065