(57)【要約】      （修正有）【課題】安定した計算が可能な分子モデル作成の方法をする。【解決手段】分子動力学法とＯＮＩＯＭ法とを組み合わせて、分子鎖の切断を伴うシミュレーションを行う際に、高分子鎖の全原子モデルからその一部を取り出して、分子鎖が切断した後のＱＭ計算用の分子モデルを作成する。全原子モデルは、複数の粒子モデルと粒子モデル間を結合するボンドモデルとを含んでいる。複数の粒子モデルから選択された一対の切断粒子モデル間で切断して、全原子モデルから互いに異なる第１分子モデル及び第２分子モデルを取り出す工程と、第１分子モデルの切断粒子モデルに第１水素粒子モデルを結合する工程と、第２分子モデルの切断粒子モデルに第２水素粒子モデルを結合する工程と、第１水素粒子モデルと第２水素粒子モデルとの間の距離が予め定められた基準距離よりも大きくなるように、第１水素粒子モデルと第２水素粒子モデルとを互いに離間させる工程と、を有する。【選択図】図６
2020086773.tif 000002
【特許請求の範囲】
【請求項１】
分子動力学法（ＭＤ法）とＯＮＩＯＭ法とを組み合わせて、分子鎖の切断を伴うシミュレーションを行う際に、高分子鎖の全原子モデルからその一部を取り出して、前記分子鎖を切断した後のＱＭ計算用の分子モデルの構造を、コンピュータを用いて作成するための方法であって、
前記全原子モデルは、複数の粒子モデルと、前記粒子モデル間を結合するボンドモデルとを含み、
前記方法は、
前記コンピュータが、複数の前記粒子モデルから選択された一対の切断粒子モデル間で切断して、前記全原子モデルから互いに異なる少なくとも第１分子モデル及び第２分子モデルを取り出す工程と、
前記コンピュータが、前記第１分子モデルの前記切断粒子モデルに、第１水素粒子モデルを結合する工程と、
前記コンピュータが、前記第２分子モデルの前記切断粒子モデルに、第２水素粒子モデルを結合する工程と、
前記コンピュータが、前記第１水素粒子モデルと、前記第２水素粒子モデルとの間の距離が、予め定められた基準距離よりも大きくなるように、前記第１水素粒子モデルと、前記第２水素粒子モデルとを互いに離間させる工程とを含む、
分子モデルの作成方法。
【請求項２】
前記離間させる工程は、前記切断粒子モデルの少なくとも一方と、その切断粒子モデルと隣接する前記粒子モデルとを結合する前記ボンドモデルの長手方向に沿った軸線回りで、
前記切断粒子モデルの少なくとも一方と、その切断粒子モデルに結合する前記第１水素粒子モデル又は前記第２水素粒子モデルとを回転させる工程を含む、請求項１記載の分子モデルの作成方法。
【請求項３】
前記基準距離は、水素原子のファンデルワールス半径に基づいて設定され、
前記第１水素粒子モデルと、前記第２水素粒子モデルとの間には、前記ファンデルワールス半径に基づく斥力を定義するための第１ポテンシャルが設定される、請求項１又は２記載の分子モデルの作成方法。
【請求項４】
前記離間させる工程は、前記第１ポテンシャルの前記ファンデルワールス半径を小さくする工程と、
前記第１ポテンシャルの前記ファンデルワールス半径を徐々に大きくしながら、前記第１水素粒子モデル、及び、前記第２水素粒子モデルを対象とするＭＤ計算を行う工程とを含む、請求項３記載の分子モデルの作成方法。
【請求項５】
前記離間させる工程は、前記第１水素粒子モデルと、前記第２水素粒子モデルとの間に、前記第１ポテンシャルよりも小さい斥力を定義するための第２ポテンシャルを設定する工程と、
前記第２ポテンシャルに基づいて、前記第１水素粒子モデル、及び、前記第２水素粒子モデルを対象とするＭＤ計算を行う工程とを含む、請求項３又は４記載の分子モデルの作成方法。
【請求項６】
前記取り出す工程は、ＭＤ計算後の前記全原子モデルについて、前記ボンドモデルを介して隣り合う前記粒子モデル間の距離が、予め定められた第１距離よりも大きい部分で切断して、前記第１分子モデル及び前記第２分子モデルをＱＭ計算対象として取り出す工程を含む、請求項１ないし５のいずれかに記載の分子モデルの作成方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、分子動力学法（ＭＤ法）とＯＮＩＯＭ法とを組み合わせて、分子鎖の切断を伴うシミュレーションに関し、詳しくは、高分子鎖の全原子モデルの一部を取り出して、分子鎖を切断した後のＱＭ計算用の分子モデルの構造を作成するための方法に関する。
【背景技術】
【０００２】
近年、ゴム等の高分子材料の開発のために、高分子材料の物理現象等を、コンピュータを用いて評価するためのシミュレーション方法（数値計算）が種々提案されている。この種のシミュレーション方法では、複数の粒子モデルと、粒子モデル間を結合するボンドモデルとを含む全原子モデルが用いられる。隣接する粒子モデル間には、例えば、引力又は斥力が作用する相互作用ポテンシャルが定義される。そして、上記シミュレーション方法では、計算コストの低いＭＭ計算（分子力学計算）や、計算精度の高いＱＭ計算（量子力学計算）が行われている。
【０００３】
下記非特許文献１は、ＭＭ計算とＱＭ計算とを組み合わせて最適解を得るＯＮＩＯＭ法を提案している。ＯＮＩＯＭ法では、全原子モデルの全体を対象とするＭＭ計算の後に、全原子モデルから計算の重要度の高い部分を取り出した分子モデルを作成し、その分子モデルを対象とするＱＭ計算が行われている。分子モデルの末端に配される粒子モデルには、例えば、原子の原子価や混成軌道に応じて、水素粒子モデルが結合される。
【先行技術文献】
【特許文献】
【０００４】
【非特許文献１】Stefan Dapprich、外４名、"A new ONIOM implementation in Gaussian98. Part I. The calculation of energies, gradients, vibrational frequencies and electric field derivatives"、Journal of Molecular Structure（Theochem）、461-462、1999、p1-21
【発明の概要】
【発明が解決しようとする課題】
【０００５】
しかしながら、ＱＭ計算が行われた分子モデルが、全原子モデルの切断後の構造としてそのまま用いられると、ＭＤ計算（分子動力学計算）が不安定になり、途中終了する等、目的とする計算結果を得ることができない場合があった。
【０００６】
発明者らは、鋭意研究を重ねた結果、分子モデルの末端に配される粒子モデルに、水素粒子モデルを結合する際に、水素粒子モデル間の距離が、現実ではあり得ない距離（例えば、水素原子のファンデルワールス半径よりも小さい距離）まで短くなる場合があり、それらの水素粒子モデル間に大きな反発力が生じ、計算が不安定になることを知見した。
【０００７】
本発明は、以上のような実状に鑑み案出されたもので、安定した計算が可能な分子モデルを作成することができる方法を提供することを主たる目的としている。
【課題を解決するための手段】
【０００８】
本発明は、分子動力学法（ＭＤ法）とＯＮＩＯＭ法とを組み合わせて、分子鎖の切断を伴うシミュレーションを行う際に、高分子鎖の全原子モデルからその一部を取り出して、前記分子鎖を切断した後のＱＭ計算用の分子モデルの構造を、コンピュータを用いて作成するための方法であって、前記全原子モデルは、複数の粒子モデルと、前記粒子モデル間を結合するボンドモデルとを含み、前記方法は、前記コンピュータが、複数の前記粒子モデルから選択された一対の切断粒子モデル間で切断して、前記全原子モデルから互いに異なる少なくとも第１分子モデル及び第２分子モデルを取り出す工程と、前記コンピュータが、前記第１分子モデルの前記切断粒子モデルに、第１水素粒子モデルを結合する工程と、前記コンピュータが、前記第２分子モデルの前記切断粒子モデルに、第２水素粒子モデルを結合する工程と、前記コンピュータが、前記第１水素粒子モデルと、前記第２水素粒子モデルとの間の距離が、予め定められた基準距離よりも大きくなるように、前記第１水素粒子モデルと、前記第２水素粒子モデルとを互いに離間させる工程とを含むことを特徴とする。
【０００９】
本発明に係る前記分子モデルの作成方法において、前記離間させる工程は、前記切断粒子モデルの少なくとも一方と、その切断粒子モデルと隣接する前記粒子モデルとを結合する前記ボンドモデルの長手方向に沿った軸線回りで、前記切断粒子モデルの少なくとも一方と、その切断粒子モデルに結合する前記第１水素粒子モデル又は前記第２水素粒子モデルとを回転させる工程を含んでもよい。
【００１０】
本発明に係る前記分子モデルの作成方法において、前記基準距離は、水素原子のファンデルワールス半径に基づいて設定され、前記第１水素粒子モデルと、前記第２水素粒子モデルとの間には、前記ファンデルワールス半径に基づく斥力を定義するための第１ポテンシャルが設定されてもよい。
【００１１】
本発明に係る前記分子モデルの作成方法において、前記離間させる工程は、前記第１ポテンシャルの前記ファンデルワールス半径を小さくする工程と、前記第１ポテンシャルの前記ファンデルワールス半径を徐々に大きくしながら、前記第１水素粒子モデル、及び、前記第２水素粒子モデルを対象とするＭＤ計算を行う工程とを含んでもよい。
【００１２】
本発明に係る前記分子モデルの作成方法において、前記離間させる工程は、前記第１水素粒子モデルと、前記第２水素粒子モデルとの間に、前記第１ポテンシャルよりも小さい斥力を定義するための第２ポテンシャルを設定する工程と、前記第２ポテンシャルに基づいて、前記第１水素粒子モデル、及び、前記第２水素粒子モデルを対象とするＭＤ計算を行う工程とを含んでもよい。
【００１３】
本発明に係る前記分子モデルの作成方法において、前記取り出す工程は、ＭＤ計算後の前記全原子モデルについて、前記ボンドモデルを介して隣り合う前記粒子モデル間の距離が、予め定められた第１距離よりも大きい部分で切断して、前記第１分子モデル及び前記第２分子モデルをＱＭ計算対象として取り出す工程を含んでもよい。
【発明の効果】
【００１４】
本発明の分子モデルの作成方法は、コンピュータが、複数の粒子モデルから選択された一対の切断粒子モデル間で切断して、全原子モデルから互いに異なる少なくとも第１分子モデル及び第２分子モデルを取り出す工程と、前記第１分子モデルの前記切断粒子モデルに、第１水素粒子モデルを結合する工程と、前記第２分子モデルの前記切断粒子モデルに、第２水素粒子モデルを結合する工程と、前記第１水素粒子モデルと、前記第２水素粒子モデルとの間の距離が、予め定められた基準距離よりも大きくなるように、前記第１水素粒子モデルと、前記第２水素粒子モデルとを互いに離間させる工程とを含んでいる。
【００１５】
本発明の作成方法は、前記第１水素粒子モデルと、前記第２水素粒子モデルとの間の距離が前記基準距離以下になるような、現実ではあり得ない距離まで接近するように配置されたとしても、前記距離が前記基準距離よりも大きくなるように、前記第１水素粒子モデルと、前記第２水素粒子モデルとを互いに離間させることができる。これにより、本発明の作成方法は、前記第１水素粒子モデルと前記第２水素粒子モデルとの間に、大きな反発力が作用するのを防ぐことができる。したがって、本発明の作成方法は、安定した計算が可能な分子モデルを作成することができる。
【図面の簡単な説明】
【００１６】
【図１】本発明の分子モデルの作成方法を実行するコンピュータの一例を示す斜視図である。
【図２】ポリイソプレンの構造式である。
【図３】セルに配置された全原子モデルの一例を示す概念図である。
【図４】全原子モデルの一例を示す概念図である。
【図５】モース型関数のポテンシャル曲線の一例を示すグラフである。
【図６】分子モデルの作成方法の処理手順の一例を示すフローチャートである。
【図７】ＭＤ計算後の各全原子モデルの一例を示す部分拡大図である。
【図８】ＱＭ計算対象として取り出された第１分子モデル及び第２分子モデルの一例を示す図である。
【図９】第１水素粒子モデルが結合された第１分子モデル、及び、第２水素粒子モデルが結合された第２分子モデルの一例を示す図である。
【図１０】大きな反発力の影響を受けた第１水素粒子モデルの一例を示す図である。
【図１１】離間工程の処理手順の一例を示すフローチャートである。
【図１２】第１水素粒子モデルと、第２水素粒子モデルとを互いに離間させる工程の一例を説明する図である。
【図１３】本発明の他の実施形態の離間工程の処理手順の一例を示すフローチャートである。
【図１４】本発明のさらに他の実施形態の離間工程の処理手順の一例を示すフローチャートである。
【図１５】第１水素粒子モデルと、第２水素粒子モデルとの間に定義された第２ポテンシャルを示す図である。
【図１６】（ａ）は、第１ポテンシャルの斥力と距離との関係を示すグラフ、（ｂ）は、第２ポテンシャルの斥力と距離との関係を示すグラフである。
【発明を実施するための形態】
【００１７】
以下、本発明の実施の一形態が図面に基づき説明される。
本実施形態の分子モデルの作成方法（以下、単に「作成方法」ということがある。）は、分子動力学法（ＭＤ法）とＯＮＩＯＭ法とを組み合わせて、分子鎖の切断を伴うシミュレーションを行う際に、全原子モデルからその一部を取り出して、分子鎖を切断した後のＱＭ計算用の分子モデルの構造を、コンピュータを用いて作成するための方法である。
【００１８】
本実施形態では、ＭＤ計算（分子動力学計算）で求めた全原子モデルから、計算の重要度の高い部分（ＱＭ計算（量子力学計算）対象）を取り出す。このＱＭ計算において、全原子モデル（分子鎖）が切断したと判断された場合には、ＯＮＩＯＭ法に基づいて、切断された部分を水素等でキャップした分子モデルが作成される。
【００１９】
図１は、本発明の分子モデルの作成方法を実行するコンピュータの一例を示す斜視図である。コンピュータ１は、本体１ａ、キーボード１ｂ、マウス１ｃ及びディスプレイ装置１ｄを含んで構成されている。この本体１ａには、例えば、演算処理装置（ＣＰＵ）、ＲＯＭ、作業用メモリ、磁気ディスクなどの記憶装置、及び、ディスクドライブ装置１ａ１、１ａ２が設けられている。また、記憶装置には、本実施形態の作成方法を実行するための処理手順（プログラム）が予め記憶されている。
【００２０】
高分子鎖を有する高分子材料としては、例えば、ゴム、樹脂又はエラストマー等が含まれる。本実施形態の高分子材料は、cis-１，４ポリイソプレン（以下、単に「ポリイソプレン」ということがある。）である場合が例示される。図２は、ポリイソプレンの構造式である。
【００２１】
ポリイソプレンを構成する高分子鎖２は、メチン基等（例えば、−ＣＨ＝、＞Ｃ＝）、メチレン基（−ＣＨ２−）、及び、メチル基（−ＣＨ３）によって構成されるイソプレンのモノマー（イソプレン分子）３が、重合度ｎで連結されて構成されている。なお、高分子材料には、ポリイソプレン以外の高分子材料が用いられてもよい。
【００２２】
本実施形態の作成方法では、例えば、高分子材料に基づいて設定された数値計算用の高分子材料モデルが用いられる。高分子材料モデルは、高分子材料の一部に対応する仮想空間であるセルに、高分子鎖２（図２に示す）をモデル化した全原子モデルが配置されている。図３は、セル４に配置された全原子モデル５の一例を示す概念図である。図４は、全原子モデル５の一例を示す概念図である。
【００２３】
図３に示されるように、セル４は、少なくとも互いに向き合う一対の面８、８、本実施形態では、互いに向き合う三対の面８、８を有しており、直方体又は立方体（本実施形態では、立方体）として定義されている。各面８、８には、周期境界条件が定義されている。これにより、セル４は、一方側の面８ａと、他方側の面８ｂとが連続している（繋がっている）ものとして取り扱うことができる。セル４の一辺の各長さＬａは、例えば、全原子モデル５の拡がりを示す量である慣性半径（図示省略）の３倍以上に設定されるのが望ましい。また、セル４の大きさは、例えば１気圧で安定な体積に設定される。これにより、セル４は、解析対象の高分子材料の少なくとも一部の体積を定義することができる。
【００２４】
全原子モデル５は、セル４の内部において、少なくとも一つ（例えば、１０個〜１，０００，０００個）配置される。全原子モデル５は、一つの硫黄原子が架橋されたポリイソプレンとして設定されている。全原子モデル５は、複数の粒子モデル６と、粒子モデル６、６間を結合するボンドモデル７とを含んで構成されている。
【００２５】
図４に示されるように、粒子モデル６及びボンドモデル７は、図２に示した高分子鎖２のモノマー３をなす単位構造に基づいて連結される。これにより、モノマーモデル９が設定される。このモノマーモデル９は、分子量（重合度）Ｍnに基づいて連結される。これにより、全原子モデル５が設定される。
【００２６】
粒子モデル６は、ＭＤ計算、ＯＮＩＯＭ法に基づくＭＭ計算（分子力学計算）、及び、ＱＭ計算（量子力学計算）において必要となるパラメータが定義される。パラメータの一例としては、原子名、質量、直径、電荷、又は、初期座標などが含まれる。本実施形態の粒子モデル６は、図２に示した高分子鎖２の炭素原子をモデル化した炭素粒子モデル６Ｃ、水素原子をモデル化した水素粒子モデル６Ｈ、及び、硫黄原子をモデル化した硫黄粒子モデル６Ｓを含んでいる。
【００２７】
ボンドモデル７は、粒子モデル６、６間を拘束するものである。本実施形態のボンドモデル７は、主鎖７ａと側鎖７ｂとを含んでいる。主鎖７ａは、炭素粒子モデル６Ｃ、６Ｃを連結するものである。側鎖７ｂは、炭素粒子モデル６Ｃと水素粒子モデル６Ｈとの間、炭素粒子モデル６Ｃと硫黄粒子モデル６Ｓとの間、及び、硫黄粒子モデル６Ｓと水素粒子モデル６Ｈとの間を連結するものである。さらに、ボンドモデル７には、一対の原子間の一重結合を定義した一重結合ボンドモデル７Ａと、一対の原子間の二重結合を定義した二重結合ボンドモデル７Ｂとを含んでいる。
【００２８】
ボンドモデル７は、隣り合う２つの粒子モデル６、６間に、相互作用（斥力及び引力を含む）が生じさせるボンドポテンシャルＢＰによって定義される。ボンドポテンシャルＢＰは、適宜設定することができる。本実施形態のボンドポテンシャルＢＰは、例えば、下記式（１）で定義されるモース型関数を用いて設定される。図５は、モース型関数のポテンシャル曲線１０の一例を示すグラフである。
【００２９】
【数１】
2020086773.tif 000003
ここで、各定数及び変数は、次のとおりである。
Ｅ：ポテンシャルエネルギー
ForceC、Dlim：フィッティングパラメータ
Ｒ：２つの粒子モデルの間の距離
ＲM：２つの粒子モデルの間の平衡距離
なお、上記式（２）の距離Ｒ、及び、平衡距離ＲMは、粒子モデル６、６の中心６ｃ、６ｃ（図４に示す）間の距離として定義される。
【００３０】
図５のグラフにおいて、モース型関数のポテンシャルエネルギーＥは、２つの粒子モデル６、６（図４に示す）の間の距離Ｒが、平衡距離ＲMよりも小さくなるほど、粒子モデル６、６間の斥力に起因して大きくなる。また、ポテンシャルエネルギーＥは、距離Ｒが平衡距離ＲMになるときに最小となり、粒子モデル６、６間に斥力及び引力は作用しない。さらに、ポテンシャルエネルギーＥは、距離Ｒが平衡距離ＲMよりも大きくなると、粒子モデル６、６間の引力に起因して大きくなる。このように、モース型関数のポテンシャル曲線１０は、距離Ｒに応じて、粒子モデル６、６間の斥力及び引力を定義することができる。本明細書において、「ポテンシャルエネルギーＥ」は、平衡距離ＲMでの全エネルギーをゼロとした相対エネルギーである。
【００３１】
モース型関数のポテンシャル曲線１０は、粒子モデル６、６間の距離Ｒが平衡距離ＲMよりも大きい領域において、ポテンシャルエネルギーの増分が小さくなる変曲点１３を有している。変曲点１３は、変曲点１３よりも距離Ｒが大きいと予想される領域Ｔａ、及び、変曲点１３よりも距離Ｒが小さいと予想される領域Ｔｂにおいて、ポテンシャル曲線１０に近似する一対の直線１４ａ、１４ｂの交点で特定することができる。
【００３２】
変曲点１３よりも距離Ｒが小さい領域Ｔｂでは、変曲点１３よりも距離Ｒが大きい領域Ｔａに比べて、モース型関数のポテンシャルエネルギーの増分が大きい。これは、変曲点１３よりも距離Ｒが小さい領域Ｔｂにおいて、図４に示した粒子モデル６、６間の結合が強固に維持されることを示している。一方、変曲点１３よりも距離Ｒが大きい領域Ｔａでは、粒子モデル６、６間で切断、又は、結合次数の減少が生じており、ポテンシャルエネルギーの増分が小さくなっている。従って、変曲点１３は、粒子モデル６、６間の結合と解離との境界点と考えることができる。
【００３３】
上記式（１）のフィッティングパラメータForceCは、モース型関数のポテンシャル曲線１０において、平衡距離ＲM近傍の曲率に相当する。また、上記式（１）のフィッティングパラメータDlimは、モース型関数のポテンシャル曲線１０において、粒子モデル６、６間の結合と解離との境界でのポテンシャルエネルギー（即ち、結合解離エネルギー）に相当する。これらのフィッティングパラメータForceC、Dlimが特定されることにより、ボンドポテンシャルＢＰを定義することができる。
【００３４】
本実施形態では、図４に示した炭素粒子モデル６Ｃ、６Ｃ間、炭素粒子モデル６Ｃと水素粒子モデル６Ｈとの間、炭素粒子モデル６Ｃと硫黄粒子モデル６Ｓとの間、及び、硫黄粒子モデル６Ｓと水素粒子モデル６Ｈとの間に、それぞれ異なるボンドポテンシャルＢＰ（即ち、それぞれ異なるフィッティングパラメータForceC、Dlim）がそれぞれ定義される。なお、フィッティングパラメータForceC、Dlimは、従来の方法に基づいて、適宜設定することができる。
【００３５】
図３に示されるように、隣接する全原子モデル５において、ボンドモデル７を介さずに隣り合う粒子モデル６、６間には、第１ポテンシャルＰ１が定義される。本実施形態の第１ポテンシャルＰ１は、下記式（２）で定義されるＬＪポテンシャルＵLJ（ｒij）である。
【００３６】
【数２】
2020086773.tif 000004
ここで、各定数及び変数は、Lennard-Jones ポテンシャルのパラメータであり、次のとおりである。
ｒij：粒子モデル間の距離
εij：粒子モデル間に定義されるＬＪポテンシャルの強度
σ0：粒子モデルの直径に相当
ｒc：カットオフ距離（Lennard-Jonesポテンシャルの計算を打ち切るための距離）
Ｕcutoff：ｒ＝ｒcの場合に、ＵLJ（ｒij）がゼロとなるようにポテンシャルをシフトさせる値
なお、距離ｒij及びカットオフ距離ｒcは、各粒子モデル６、６の中心６ｃ、６ｃ（図４に示す）間の距離として定義される。
【００３７】
本実施形態では、炭素粒子モデル６Ｃ、６Ｃ（図４に示す）間、硫黄粒子モデル６Ｓ、６Ｓ（図４に示す）間、炭素粒子モデル６Ｃと硫黄粒子モデル６Ｓとの間、炭素粒子モデル６Ｃと水素粒子モデル６Ｈとの間、及び、硫黄粒子モデル６Ｓと水素粒子モデル６Ｈとの間に、それぞれ異なる第１ポテンシャルＰ１が設定されている。各第１ポテンシャルＰ１は、上記式（２）の定数がそれぞれ異なっている。なお、各定数は、例えば、論文（S. L. Mayo, B. D. Olafson & W. A. Goddard III 著、「DREIDING: A Generic Force Field for Molecular Simulations」、J. Phys. Chem. 1990, 94, 8897）に基づいて、適宜設定することができる。
【００３８】
本実施形態のパラメータσ0は、ファンデルワールス半径に基づいて設定される。ファンデルワールス半径は、原子や分子等がお互いに入り込むことができない半径である。これにより、第１ポテンシャルＰ１は、粒子モデル６、６間の距離が、パラメータ２1/6σ0（即ち、ＬＪポテンシャルＵLJ（ｒij）が極小値を与える各粒子モデル６の中心間の距離２1/6σ0）未満になった場合に、ファンデルワールス半径に基づく斥力（反発力）を定義することができる。
【００３９】
このような全原子モデル５が配置されたセル４について、ＭＭ計算やＭＤ計算に基づく構造緩和が計算されることにより、高分子材料モデル１６（図３に示す）が設定される。ＭＤ計算では、例えば、セル４について所定の時間、全原子モデル５が古典力学に従うものとして、ニュートンの運動方程式が適用される。そして、各時刻での粒子モデル６の動きが、単位時間ステップ毎に追跡される。このような構造緩和の計算は、例えばソフトマテリアル総合シミュレーター（ＯＣＴＡ）に含まれるＣＯＧＮＡＣを用いて処理することができる。
【００４０】
ＭＤ計算では、例えば、セル４において、圧力及び温度が一定、又は、体積及び温度が一定に保たれる。また、ＭＤ計算では、全原子モデル５の人為的な初期配置が十分に排除されたとみなせるまで行われる。これにより、実際の高分子材料の分子運動に近似させた全原子モデル５のシミュレーションを行うことができる。
【００４１】
本実施形態のシミュレーション（ＭＤ法とＯＮＩＯＭ法とを組み合わせたシミュレーション）では、先ず、高分子材料モデル１６を一軸方向（例えば、Ｚ軸方向）に伸長させて、全原子モデル５の切断を計算する伸長シミュレーションが実施される。
【００４２】
伸長シミュレーションでは、高分子材料モデル１６を構成する全ての粒子モデル６を対象とするＭＤ計算によって、シミュレーションの単位時間ステップ毎に、伸長した高分子材料モデル１６に作用する力及びエネルギーが求められる。
【００４３】
ＭＤ計算では、粒子モデル６に作用する力を、分子力学計算（ＭＭ計算）で求めた後に、時間差分計算で次の時間ステップにおける座標と速度を計算している。通常のＭＭ計算では、高分子材料モデル１６の全体領域を対象として、上述のボンドポテンシャルＢＰ（図４に示す）、アングルポテンシャルＡＰ、トーションポテンシャルＴＰ、及び、第１ポテンシャルＰ１（図３に示す）を用いたＭＭ計算（Ｒｅａｌ：Ｌｏｗ）によって、粒子モデル６に作用する力が計算される。ここで、「ＭＭ計算（Ｒｅａｌ：Ｌｏｗ）」とは、ＯＮＩＯＭ法において、全領域を対象とするＭＭ計算を意味している。
【００４４】
本実施形態のシミュレーション（ＭＤ計算）では、粒子モデル６、６間で切断しそうな部分（即ち、計算の重要度が高い部分）が存在する場合に、ＭＭ計算（Ｒｅａｌ：Ｌｏｗ）に加え、その切断しそうな部分を含むように全原子モデル５から切り出したモデルに対しても、ＭＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）とＱＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）とが実施される。粒子モデル６、６間で切断しそうな部分としては、例えば、ボンドモデル７を介して隣り合う粒子モデル６、６間の距離Ｒ（図４に示す）が、予め定められた第１距離Ｒａ（図示省略）よりも大きい部分として定義される。「ＭＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）」とは、ＯＮＩＯＭ法において、切り出したモデルを対象とするＭＭ計算を意味している。また、「ＱＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）」とは、ＯＮＩＯＭ法において、切り出したモデルを対象とするＱＭ計算を意味している。粒子モデル６に作用する力は、ＭＭ計算（Ｒｅａｌ：Ｌｏｗ）、ＭＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）及びＱＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）を用いて、ＯＮＩＯＭ法に基づき計算される。
【００４５】
第１距離Ｒａは、図５に示したモース型関数のポテンシャル曲線１０において、平衡距離ＲMよりも大、かつ、変曲点１３での距離Ｒｐよりも小に設定されている。上述したように、変曲点１３は、粒子モデル６、６間の結合と解離との境界点と考えることができる。従って、第１距離Ｒａよりも大きく離間した粒子モデル６、６間において、切断する可能性が高いと考えることができる。
【００４６】
ＱＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）では、取り出された部分を水素原子等でキャップした後に、全原子モデル５（分子モデル）が、原子核と電子とに基づく相互作用に従うものとして、粒子モデル６の電子状態が計算される。ＱＭ計算では、全原子モデル５に、量子力学計算に使用する理論及び基底関数が定義される。本実施形態で定義される理論としては、例えば、密度汎関数法（汎関数：Ｂ３ＬＹＰ）が定義される。また、本実施形態の基底関数としては、例えば、基底関数（６−３１Ｇ（ｄ））が定義される。ＱＭ計算は、Ｇａｕｓｓｉａｎ社製の量子化学計算プログラムＧａｕｓｓｉａｎ０３、又は、Ｇａｕｓｓｉａｎ０９を用いて処理することができる。
【００４７】
ＱＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）及びＭＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）の結果から、ＯＮＩＯＭ法に基づいて、粒子モデル６に作用する力が計算される。この粒子モデル６に作用する力を用いて、粒子モデル６に作用する加速度が求められる。これにより、単位時間ステップ後の粒子モデル６の座標を計算することができる（ＭＤ計算）。次に、この計算結果に、伸長による移動量を加算した値が、全原子モデル５に反映される。そして、ボンドモデル７を介して隣り合う粒子モデル６、６間の距離Ｒ（図４に示す）が、予め定められた第２距離Ｒｂ（図示省略）よりも大きい部分で切断が定義される。切断後の構造としては、ＯＮＩＯＭ法に基づくＱＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）やＭＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）と同様に、切断した部分を水素原子等でキャップした構造が用いられる。第２距離Ｒｂは、図５に示したモース型関数のポテンシャル曲線１０において、第１距離Ｒａよりも大、かつ、変曲点（粒子モデル６、６間の結合と解離との境界点）１３での距離Ｒｐ以上に設定されている。このようなＱＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）、ＭＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）、及び、ＭＤ計算が、単位時間ステップ毎に繰り返し行われることで、予め定められた長さ（上限値）まで伸長した高分子材料モデル１６が計算される。
【００４８】
このように、分子動力学法とＯＮＩＯＭ法とを組み合わせたシミュレーションでは、ＱＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）が、粒子モデル６、６間で切断しそうな部分（即ち、計算の重要度が高い部分）のみを対象に実施される。このため、例えば、高分子材料モデル１６の全体領域を対象に量子力学計算を実施する場合に比べて、計算時間を大幅に短縮することができる。
【００４９】
本実施形態の作成方法は、ＭＤ計算後の全原子モデル５から、分子鎖が切断しそうな部分（計算の重要度の高い部分）を取り出して、互いに異なる少なくとも第１分子モデル及び第２分子モデルが作成される。作成された第１分子モデル及び第２分子モデルは、ＯＮＩＯＭ法に基づくＱＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）及びＭＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）に用いられる。図６は、分子モデルの作成方法の処理手順の一例を示すフローチャートである。
【００５０】
本実施形態の作成方法では、先ず、コンピュータ１が、全原子モデル５から互いに異なる少なくとも第１分子モデル２１Ａ及び第２分子モデル２１Ｂを取り出す（工程Ｓ１）。図７は、ＭＤ計算後の各全原子モデル５の一例を示す部分拡大図である。図８は、ＱＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）対象として取り出された第１分子モデル２１Ａ及び第２分子モデル２１Ｂの一例を示す図である。なお、図８では、第１分子モデル２１Ａ及び第２分子モデル２１Ｂのみが示されているが、計算の重要度の高い部分が他にも存在する場合には、第１分子モデル２１Ａ及び第２分子モデル２１Ｂ以外に、複数の分子モデル２１（図示省略）が取り出される。
【００５１】
図７に示されるように、本実施形態の工程Ｓ１では、ＭＤ計算後の各全原子モデル５において、複数の粒子モデル６から選択された一対の切断粒子モデル１１、１１間（即ち、図７において、一点鎖線とボンドモデル７とが交差する部分）で切断する。これにより、少なくとも第１分子モデル２１Ａ及び第２分子モデル２１Ｂ（図８に示す）がＱＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）対象として取り出される。
【００５２】
図７に示されるように、切断粒子モデル１１、１１は、ボンドモデル７を介して隣り合う切断しそうな（即ち、距離Ｒが第１距離Ｒａ（図５に示す）よりも大きい）一対の粒子モデル（図７において、一点鎖線で囲まれた粒子モデル）６と、その粒子モデル６に隣接して結合されている粒子モデル（図７において、一点鎖線で囲まれていない粒子モデル）６とで構成されている。なお、切断粒子モデル１１、１１は、シミュレーションに応じて適宜設定することができる。第１分子モデル２１Ａ及び第２分子モデル２１Ｂを含む取り出された分子モデル２１は、コンピュータ１に記憶される。
【００５３】
次に、本実施形態の作成方法では、コンピュータ１が、第１分子モデル２１Ａの切断粒子モデル１１に、第１水素粒子モデル２２Ａを結合する（工程Ｓ２）。図９は、第１水素粒子モデル２２Ａが結合された第１分子モデル２１Ａ、及び、第２水素粒子モデル２２Ｂが結合された第２分子モデル２１Ｂの一例を示す図である。
【００５４】
工程Ｓ２では、図８及び図９に示されるように、第１分子モデル２１Ａの切断粒子モデル１１を構成する原子の原子価や混成軌道に応じて、切断粒子モデル１１に、少なくとも一つの第１水素粒子モデル２２Ａが結合される。第１水素粒子モデル２２Ａと切断粒子モデル１１との間には、ボンドモデル７が配置されている。ボンドモデル７は、上記式（１）のボンドポテンシャルＢＰで定義される。
【００５５】
第１水素粒子モデル２２Ａの結合には、セル４（図３に示す）での第１水素粒子モデル２２Ａの位置情報（位置ベクトル）が予め特定されている必要がある。第１水素粒子モデル２２Ａの位置情報は、例えば、上記非特許文献１に基づいて定義される。
【００５６】
第１水素粒子モデル２２Ａは、工程Ｓ２で新たに追加されたものであるため、第１分子モデル２１Ａを構成する他の粒子モデル６とは異なり、上記式（２）で定義される第１ポテンシャルＰ１（図３に示す）が設定されていない。このため、図９に示されるように、第１水素粒子モデル２２Ａと、カットオフ距離ｒc以内に存在する第１分子モデル２１Ａの粒子モデル６との間に、第１ポテンシャルＰ１が設定される。本実施形態の第１ポテンシャルＰ１は、第１水素粒子モデル２２Ａから予め定められた個数（個数は、力場ポテンシャルの種類に依存（一例として２〜３個））のボンドモデル７を介して隣接する粒子モデル６を除いて定義される。第１水素粒子モデル２２Ａが結合された第１分子モデル２１Ａは、コンピュータ１に記憶される。
【００５７】
次に、本実施形態の作成方法では、図８及び図９に示されるように、コンピュータ１が、第２分子モデル２１Ｂの切断粒子モデル１１に、第２水素粒子モデル２２Ｂを結合する（工程Ｓ３）。工程Ｓ３では、工程Ｓ３と同様に、第２分子モデル２１Ｂの切断粒子モデル１１を構成する原子の原子価や混成軌道に応じて、切断粒子モデル１１に、ボンドモデル７を介して、少なくとも一つの第２水素粒子モデル２２Ｂが結合される。
【００５８】
第２水素粒子モデル２２Ｂは、工程Ｓ３で新たに追加されたものであるため、第２分子モデル２１Ｂを構成する他の粒子モデル６とは異なり、上記式（２）で定義される第１ポテンシャルＰ１（図３に示す）が設定されていない。このため、第２水素粒子モデル２２Ｂと、カットオフ距離ｒc以内に存在する第２分子モデル２１Ｂの粒子モデル６との間に、第１ポテンシャルＰ１が設定される。本実施形態の第１ポテンシャルＰ１は、第２水素粒子モデル２２Ｂから予め定められた個数（個数は、力場ポテンシャルの種類に依存（一例として２〜３個））のボンドモデル７を介して隣接する粒子モデル６を除いて定義される。
第２水素粒子モデル２２Ｂが結合された第２分子モデル２１Ｂは、コンピュータ１に記憶される。
【００５９】
ところで、切断後の構造として、図９に示した第１分子モデル２１Ａ及び第２分子モデル２１Ｂを、ＭＤ計算にそのまま用いた従来のシミュレーションでは、そのＭＤ計算が不安定になり、途中終了する等、目的とする計算結果を得ることができない場合があった。
【００６０】
発明者らは、鋭意研究を重ねた結果、ＭＤ計算とＯＮＩＯＭ法とを組み合わせたシミュレーションにおいて、切断粒子モデル１１に隣接する第１水素粒子モデル２２Ａ、及び、第２水素粒子モデル２２Ｂが、現実ではあり得ない距離（例えば、水素原子のファンデルワールス半径よりも小さい距離）まで接近するように配置されると、それらの間に大きな反発力が生じ、計算が不安定になることを知見した。図１０は、大きな反発力の影響を受けた第１水素粒子モデル２２Ａの一例を示す図である。大きな反発力の影響を受けた第１水素粒子モデル２２Ａは、その結合距離が異常に長くなり、第１水素粒子モデル２２Ａを結合するボンドモデル７が伸縮して、大きな振動が生じる場合がある。
【００６１】
本実施形態の作成方法では、上記のような不安定な計算を防ぐために、コンピュータ１が、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとを互いに離間させている（離間工程Ｓ４）。離間工程Ｓ４では、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとの間の距離Ｌ１が、予め定められた基準距離よりも大きくなるように、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとを互いに離間させている。なお、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとを互いに離間させる際に、第１水素粒子モデル２２Ａ及び第２水素粒子モデル２２Ｂ以外の粒子モデル６を固定するのが望ましい。これにより、第１分子モデル２１Ａ及び第２分子モデル２１Ｂの構造が変化するのを防ぐことができる。
【００６２】
基準距離については、適宜設定することができる。上述したように、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間には、第１ポテンシャルＰ１（図９に示す）が設定されている。第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間の距離Ｌ１が、パラメータσ0の２1/6倍未満になった場合に、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間に、斥力（即ち、反発力）が生じる。したがって、本実施形態の基準距離は、水素原子のＬＪポテンシャルＵLJ（ｒij）のパラメータσ0（ファンデルワールス半径に相当）に基づいて設定される。図１１は、離間工程Ｓ４の処理手順の一例を示すフローチャートである。
【００６３】
本実施形態の離間工程Ｓ４では、先ず、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間の距離Ｌ１（図９に示す）が、基準距離未満か否かが判断される（工程Ｓ４１）。第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間の距離Ｌ１が、基準距離未満であると判断された場合（工程Ｓ４１において、「Ｙ」）、第１水素粒子モデル２２Ａ及び第２水素粒子モデル２２Ｂの間に大きな反発力が生じる。このため、次の工程Ｓ４２において、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとを互いに離間させる。
【００６４】
他方、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間の距離Ｌ１（図９に示す）が、基準距離以上であると判断された場合（工程Ｓ４１において、「Ｎ」）、第１水素粒子モデル２２Ａ及び第２水素粒子モデル２２Ｂの間に大きな反発力は生じない。このため、離間工程Ｓ４の一連の処理が終了する。
【００６５】
次に、本実施形態の離間工程Ｓ４では、切断粒子モデル１１の少なくとも一方と、その切断粒子モデル１１に結合する第１水素粒子モデル２２Ａ又は第２水素粒子モデル２２Ｂとを回転させる（工程Ｓ４２）。図１２は、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとを互いに離間させる工程の一例を説明する図である。
【００６６】
図１２に示されるように、工程Ｓ４２では、第１分子モデル２１Ａにおいて、切断粒子モデル１１、１１の少なくとも一方（第１水素粒子モデル２２Ａが結合されている切断粒子モデル１１）と、その切断粒子モデル１１と隣接する粒子モデル６とを結合するボンドモデル７の長手方向に沿った軸線回りで、切断粒子モデル１１の少なくとも一方と、その切断粒子モデル１１に結合する第１水素粒子モデル２２Ａ又は第２水素粒子モデル２２Ｂとを回転させている。図１２では、第１分子モデル２１Ａの切断粒子モデル１１と、第１水素粒子モデル２２Ａとを回転させているが、第２分子モデル２１Ｂの切断粒子モデル１１と、第２水素粒子モデル２２Ｂとを回転させてもよい。
【００６７】
本実施形態の工程Ｓ４２では、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとの間の距離Ｌ１が基準距離よりも大きくなるまで、切断粒子モデル１１及び第１水素粒子モデル２２Ａを、ボンドモデル７の軸線Ｓの回りで回転させている。これにより、工程Ｓ４２では、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとを簡単に離間させることができる。
【００６８】
なお、切断粒子モデル１１と隣接する粒子モデル６とを結合するボンドモデル７が二重結合ボンドモデル７Ｂ（図４に示す）であると、容易に回転させることができない。このため、上記のような回転は、切断粒子モデル１１と隣接する粒子モデル６とを結合するボンドモデル７が、一重結合ボンドモデル７Ａである場合のみに実施されるのが望ましい。
【００６９】
このように、本実施形態の作成方法は、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとの間の距離Ｌ１が基準距離以下になるような、現実ではあり得ない距離まで接近するように配置されたとしても（図１２において、二点鎖線で示す）、距離Ｌ１が基準距離よりも大きくなるように、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとを離間させることができる。これにより、本実施形態の作成方法は、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間に、大きな反発力が作用するのを防ぐことができる。したがって、本実施形態の作成方法は、安定した計算が可能な分子モデル２１を作成することができる。
【００７０】
離間工程Ｓ４では、工程Ｓ１において、第１分子モデル２１Ａ及び第２分子モデル２１Ｂとは異なる他の分子モデル（図示省略）が取り出され、かつ、第１水素粒子モデル２２Ａ及び第２水素粒子モデル２２Ｂ以外の水素粒子モデル（図示省略）が結合された場合、これらの水素粒子モデルについても、近接する水素粒子モデルから離間させるのが望ましい。
【００７１】
上述したように、本実施形態の作成方法で作成された分子モデル２１は、ＱＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）及びＭＭ計算（ＭＯＤＥＬ：Ｈｉｇｈ）に用いられ、これらの計算結果に基づくＭＤ計算が行われた後に、ボンドモデル７を介して隣り合う粒子モデル６、６間の距離Ｒ（即ち、ボンド長）が、第２距離Ｒｂよりも大きい場合に、分子鎖が切断した後の構造として、全原子モデル５（図７に示す）に反映される。そして、このようなＯＮＩＯＭ法と組み合せたＭＤ計算が所定の回数繰り返されたシミュレーション結果に基づいて、高分子材料の物理現象（本実施形態では、全原子モデル５の切断）が評価される。その評価結果が良好であれば、全原子モデル５（分子鎖）の構造に基づいて、高分子材料が製造される。他方、評価結果が良好でない場合には、高分子鎖の構造等を変更し、ＯＮＩＯＭ法に基づくシミュレーションが再度実施される。これにより、所望の性能を有する高分子材料を効率よくかつ確実に製造することができる。
【００７２】
本実施形態の離間工程Ｓ４では、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとを離間させるために、第１水素粒子モデル２２Ａ又は第２水素粒子モデル２２Ｂを回転させる工程Ｓ４２が例示されたが、このような態様に限定されない。図１３は、本発明の他の実施形態の離間工程Ｓ４の処理手順の一例を示すフローチャートである。なお、この実施形態において、これまでの実施形態と同一の構成については、同一の符号を付し、説明を省略することがある。
【００７３】
この実施形態の離間工程Ｓ４では、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間の距離Ｌ１（図９に示す）が、基準距離未満であると判断された場合（工程Ｓ４１において、「Ｙ」）、図９に示した第１ポテンシャルＰ１（ＬＪポテンシャルＵLJ（ｒij））のパラメータσ0（ファンデルワールス半径に相当）を小さくする（工程Ｓ４３）。  
【００７４】
工程Ｓ４３で定義されるパラメータσ0の値については、適宜設定することができる。この実施形態のパラメータσ0の値は、工程Ｓ４１で基準距離未満と判断された距離Ｌ１よりも小さく（例えば、距離Ｌ１の半分より小さく）設定されるのが望ましい。これにより、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間に、斥力（反発力）が生じるのを防ぐことができる。
【００７５】
次に、この実施形態では、図９に示した第１ポテンシャルＰ１（ＬＪポテンシャルＵLJ（ｒij））のパラメータσ0（ファンデルワールス半径に相当）を徐々に大きくし（元に戻し）ながら、第１水素粒子モデル２２Ａ及び第２水素粒子モデル２２Ｂを対象とするＭＤ計算を行う（工程Ｓ４４）。工程Ｓ４４では、予め定められた時間間隔毎に、ＬＪポテンシャルＵLJ（ｒij）のパラメータσ0（ファンデルワールス半径に相当）を一定の増分値で大きくしながら、第１水素粒子モデル２２Ａ及び第２水素粒子モデル２２Ｂを対象とするＭＤ計算が行われる。これにより、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間に大きな反発力が生じるのを防ぎつつ、第１ポテンシャルＰ１に基づいて、第１水素粒子モデル２２Ａ及び第２水素粒子モデル２２Ｂの配置が緩和する方向へと離間させることができる。
【００７６】
本実施形態の工程Ｓ４４では、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間の距離Ｌ１（図９に示す）が、基準距離よりも大きくなるまで、ＭＤ計算が実施される。これにより、工程Ｓ４４では、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとを離間させることができる。
【００７７】
ＬＪポテンシャルＵLJ（ｒij）のパラメータσ0（ファンデルワールス半径に相当）を大きくする時間間隔については、適宜設定することができる。なお、時間間隔が小さいと、第１水素粒子モデル２２Ａ及び第２水素粒子モデル２２Ｂが十分に離間する前に、ファンデルワールス半径が大きくなってしまい、第１水素粒子モデル２２Ａ及び第２水素粒子モデル２２Ｂを離間できないおそれがある。他方、時間間隔が大きいと、第１水素粒子モデル２２Ａ及び第２水素粒子モデル２２Ｂを離間させるまでに、大きな時間を要するおそれがある。このため、時間間隔は、０．１〜５０fsが望ましい。
【００７８】
ＬＪポテンシャルＵLJ（ｒij）のパラメータσ0（ファンデルワールス半径に相当）の増分値についても、適宜設定することができる。なお、増分値が大きいと、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとが十分に離間する前に、ファンデルワールス半径が増分値に基づいて大きくなると、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間に非常に強い斥力が生じるおそれがある。他方、増分値が小さいと、第１水素粒子モデル２２Ａ及び第２水素粒子モデル２２Ｂを離間させるまでに、大きな時間を要するおそれがある。このため、増分値は、０．０１〜０．０５nmに設定されるのが望ましい。
【００７９】
このように、この実施形態の作成方法でも、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとを離間させることができるため、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間に、大きな反発力が作用するのを防ぐことができる。したがって、この実施形態の作成方法は、安定した計算が可能な分子モデル２１を作成することができる。
【００８０】
さらに、この実施形態の作成方法では、切断粒子モデル１１と隣接する粒子モデル６とを結合するボンドモデル７が、その長手方向に沿った軸線回りで回転し難い二重結合ボンドモデル７Ｂ（図４に示す）であったとしても、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとを容易に離間させることができる。したがって、この実施形態の作成方法は、安定した計算が可能な分子モデル２１を、確実に作成することができる。
【００８１】
これまでの実施形態の離間工程では、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとを離間させるために、第１水素粒子モデル２２Ａ又は第２水素粒子モデル２２Ｂを回転させる工程Ｓ４２（図１１に示す）や、ＬＪポテンシャルＵLJ（ｒij）のパラメータσ0（ファンデルワールス半径に相当）を一定の増分値で大きくする工程Ｓ４３及び工程Ｓ４４（図１３に示す）が例示されたが、これらの態様に限定されない。図１４は、本発明のさらに他の実施形態の離間工程Ｓ４の処理手順の一例を示すフローチャートである。なお、この実施形態において、これまでの実施形態と同一の構成については、同一の符号を付し、説明を省略することがある。
【００８２】
この実施形態の離間工程Ｓ４では、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間の距離Ｌ１（図９に示す）が、基準距離未満であると判断された場合（工程Ｓ４１において、「Ｙ」）、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとの間に、第２ポテンシャルＰ２が設定される（工程Ｓ４５）。第２ポテンシャルＰ２は、距離Ｌ１が小さい場合において、第１ポテンシャルＰ１（図９に示す）よりも小さい斥力を定義するためのものである。図１５は、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとの間に定義された第２ポテンシャルＰ２を示す図である。
【００８３】
第２ポテンシャルＰ２としては、適宜設定することができる。本実施形態の第２ポテンシャルＰ２は、とりわけ距離Ｌ１が小さい場合において、第１ポテンシャルＰ１として定義されるＬＪポテンシャルＵLJ（ｒij）（上記式（２））よりも小さい斥力を定義可能なソフトコアポテンシャルＵsoftcore（ｒij）である。
【００８４】
【数３】
2020086773.tif 000005
ここで、各定数及び変数は、 ソフトコアポテンシャルのパラメータであり、次のとおりである。
ｒij：粒子モデル間の距離
ａij：粒子モデル間に定義される斥力の強度に対応する定数
ｒc：カットオフ距離
なお、距離ｒij及びカットオフ距離ｒcは、各粒子モデル６、６の中心６ｃ、６ｃ（図４に示す）間の距離として定義される。
【００８５】
ソフトコアポテンシャルＵsoftcore（ｒij）は、上記式（３）に限定されるわけではなく、ＬＪポテンシャルのようにｒijが０となった時に発散しない関数であれば適宜用いることができる。ソフトコアポテンシャルＵsoftcore（ｒij）の他の例として、下記式（４）に示す。
【００８６】
【数４】
2020086773.tif 000006
ここで、各定数及び変数は、ソフトコアポテンシャルのパラメータであり、変数ｒを除いて、上記式（３）と同一である。なお、変数ｒは、粒子モデルi, j間の距離を示している。
【００８７】
定数ａijの値としては、例えば、１０［ε］を中心に、５〜３０［ε］の範囲で適宜設定することができる。図１６（ａ）は、第１ポテンシャルＰ１（ＬＪポテンシャルＵLJ（ｒij））の斥力と距離との関係を示すグラフである。図１６（ｂ）は、第２ポテンシャルＰ２（ソフトコアポテンシャルＵsoftcore（ｒij））の斥力と距離との関係を示すグラフである。
【００８８】
図１６（ｂ）に示した第２ポテンシャルＰ２（ソフトコアポテンシャルＵsoftcore（ｒij））は、図１６（ａ）に示した第１ポテンシャルＰ１（ＬＪポテンシャルＵLJ（ｒij））に比べて、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間の距離Ｌ１が小さくなっても（０に近づいても）、斥力の上昇が穏やかである。このため、第２ポテンシャルＰ２は、次のＭＤ計算を行う工程Ｓ４６において、図１５に示した第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間に大きな反発力が生じるのを防ぎつつ、第１水素粒子モデル２２Ａ及び第２水素粒子モデル２２Ｂの配置が緩和する方向へと離間させることができる。
【００８９】
次に、この実施形態の離間工程Ｓ４では、第２ポテンシャルＰ２（図５に示す）に基づいて、第１水素粒子モデル２２Ａ、及び、第２水素粒子モデル２２Ｂを対象とするＭＤ計算が行われる（工程Ｓ４６）。本実施形態の工程Ｓ４６では、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間の距離Ｌ１（図１５に示す）が、基準距離よりも大きくなるまで、ＭＤ計算が実施される。これにより、工程Ｓ４６では、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとの間の距離Ｌ１が基準距離よりも大きくなるように、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとを離間させることができる。
【００９０】
ＭＤ計算が実施される単位時間ステップ数（ＭＤステップ）については、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間の距離Ｌ１を、基準距離よりも大きくすることができれば、適宜設定することができる。なお、単位時間ステップ数が小さいと、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとが十分に離間しないおそれがある。逆に、単位時間ステップ数が大きいと、第１水素粒子モデル２２Ａ及び第２水素粒子モデル２２Ｂ以外の構造が大きく変化したり、計算時間が増大したりするおそれがある。このような観点より、単位時間ステップ数は、５００〜５００００ＭＤステップが望ましい。
【００９１】
このように、この実施形態の作成方法は、これまでの実施形態の作成方法と同様に、第１水素粒子モデル２２Ａと、第２水素粒子モデル２２Ｂとを離間させることができるため、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとの間に、大きな反発力が作用するのを防ぐことができる。したがって、この実施形態の作成方法は、安定した計算が可能な分子モデル２１を作成することができる。
【００９２】
さらに、この実施形態の作成方法では、切断粒子モデル１１と隣接する粒子モデル６とを結合するボンドモデル７が、その長手方向に沿った軸線回りで回転し難い二重結合ボンドモデル７Ｂ（図４に示す）であっても、第１水素粒子モデル２２Ａと第２水素粒子モデル２２Ｂとを容易に離間させることができる。したがって、この実施形態の作成方法は、安定した計算が可能な分子モデル２１を、確実に作成することができる。これにより、高分子材料の物理現象を、コンピュータ１で確実に評価することができるため、高分子材料の開発及び製造に役立つ。
【００９３】
以上、本発明の特に好ましい実施形態について詳述したが、本発明は図示の実施形態に限定されることなく、種々の態様に変形して実施しうる。
【実施例】
【００９４】
分子動力学法（ＭＤ法）とＯＮＩＯＭ法とを組み合わせた伸長シミュレーションにおいて、高分子材料の全原子モデルから、ＱＭ計算対象として取り出した分子モデルが作成された（実施例１〜実施例３及び比較例）。実施例１〜３では、図６に示した処理手順に基づいて、第１水素粒子モデルと、第２水素粒子モデルとの間の距離が、予め定められた基準距離よりも大きくなるように、第１水素粒子モデルと、第２水素粒子モデルとを互いに離間させる離間工程が実施された。
【００９５】
実施例１の離間工程では、図１１に示した手順に基づいて、第１水素粒子モデル又は第２水素粒子モデルを回転させる工程が実施され、分子モデルが作成された。
【００９６】
実施例２の離間工程では、図１３の手順に基づいて、第１ポテンシャルのパラメータσ0（ファンデルワールス半径に相当）を徐々に大きくしながら、第１水素粒子モデル、及び、第２水素粒子モデルを対象とするＭＤ計算を行う工程が実施されて、分子モデルが作成された。
【００９７】
実施例３の離間工程では、図１４に示した手順に基づいて、第２ポテンシャルに基づいて、第１水素粒子モデル、及び、第２水素粒子モデルを対象とするＭＤ計算を行う工程が実施され、分子モデルが作成された。
【００９８】
比較のために、上記非特許文献１と同様の手順に基づいて、分子モデルの末端に配される粒子モデル（切断粒子モデル）に、原子の原子価や混成軌道に応じて、水素粒子モデルを結合させて、分子モデルが作成された。
【００９９】
実施例１〜実施例３及び比較例で作成された分子モデルを用いて、ＯＮＩＯＭ法に基づくＱＭ計算とＭＤ計算とを組み合わせた伸長シミュレーションが１０回実施され、伸長シミュレーションの成功回数が集計された。共通仕様は、次のとおりである。
コンピュータ：ＳＧＩ社製のワークステーション
ＣＰＵ：１２コア
メモリ：６４ＧＢ
高分子材料モデル：
セルの１辺の長さＬａ：２．３nm
全原子モデルの個数：１個
全原子モデル：
ポリイソプレン（イソプレン分子：１００重合）に硫黄（１原子）を架橋
粒子モデル（原子）の個数：１３０１個
変曲点での距離Ｒｐ：３．０９Å
第１距離Ｒａ：２．６９Å
第２距離Ｒｂ：３．３５Å
伸長シミュレーション：
伸長方向：一軸方向（Ｚ軸方向）
伸長速度Ｖａ：２５０ｍ／sec
セルの変形率の上限：１１００％
実施例２：
ファンデルワールス半径を大きくする時間間隔：１fs
ファンデルワールス半径の増分値：０．０１nm
実施例３：
単位時間ステップ数：５０００ＭＤステップ
テストの結果を、表１に示す。
【０１００】
【表１】
2020086773.tif 000007
【０１０１】
テストの結果、実施例１〜３の成功回数は、比較例の成功回数に比べて、大幅に増加した。したがって、実施例１〜３は、比較例に比べて、安定した計算が可能な分子モデルを作成することができた。
【０１０２】
実施例１では、１回の伸張シミュレーションにおいて、第１水素粒子モデル又は第２水素粒子モデルを回転させても、第１水素粒子モデルと、第２水素粒子モデルとを互いに離間させることができず、計算に失敗した。一方、実施例２及び実施例３では、全ての伸張シミュレーションにおいて失敗することはなかった。したがって、実施例２及び実施例３は、実施例１に比べて、安定した計算が可能な分子モデルを、より確実に作成することができた。
【０１０３】
実施例３は、実施例２に比べて、分子モデルを作成するのに必要な時間が２０％であった。したがって、実施例３は、分子モデルを効率よく作成することができた。
【符号の説明】
【０１０４】
Ｓ１  第１分子モデル及び第２分子モデルを取り出す工程
Ｓ２  第１分子モデルの切断粒子モデルに、第１水素粒子モデルを結合する工程
Ｓ３  第２分子モデルの切断粒子モデルに、第２水素粒子モデルを結合する工程
Ｓ４  第１水素粒子モデルと、第２水素粒子モデルとを互いに離間させる工程
【図１】
2020086773.tif 000008
【図２】
2020086773.tif 000009
【図３】
2020086773.tif 000010
【図４】
2020086773.tif 000011
【図５】
2020086773.tif 000012
【図６】
2020086773.tif 000013
【図７】
2020086773.tif 000014
【図８】
2020086773.tif 000015
【図９】
2020086773.tif 000016
【図１０】
2020086773.tif 000017
【図１１】
2020086773.tif 000018
【図１２】
2020086773.tif 000019
【図１３】
2020086773.tif 000020
【図１４】
2020086773.tif 000021
【図１５】
2020086773.tif 000022
【図１６】
2020086773.tif 000023