JP2020086835A - ゴム材料のシミュレーション方法、及び、ゴム材料の製造方法 - Google Patents

Abstract

【課題】 ゴム材料の破壊特性をより精度良く評価する。【解決手段】 物性値に差異を有する第１高分子成分と第２高分子成分とを含むゴム材料の破壊現象のシミュレーション方法である。この方法は、第１高分子成分を第１粗視化粒子で表現した第１高分子鎖モデルを入力する工程Ｓ１１と、第２高分子成分を第２粗視化粒子で表現した第２高分子鎖モデルを入力する工程Ｓ２２と、第１高分子鎖モデル及び第２高分子鎖モデルを、セルに配置したゴム材料モデルを定義する工程Ｓ１６と、第１粗視化粒子間の第１相互作用ポテンシャルを定義する工程Ｓ２１と、第２粗視化粒子間の第２相互作用ポテンシャルを定義する工程Ｓ２２と、コンピュータが、ゴム材料モデルの変形を計算して、第１相互作用ポテンシャル、及び、第２相互作用ポテンシャルに基づく界面での破壊現象を計算する工程とを含む。【選択図】図３

Description

本発明は、物性値が異なる複数の高分子成分を含むゴム材料の破壊現象を、コンピュータを用いてシミュレーションするための方法などに関する。

下記非特許文献１は、ホモポリマーＡ、Ｂを含むポリマーブレンドについて、それらの界面での破壊特性を評価するためのシミュレーション方法を提案している。下記非特許文献１の方法では、先ず、ホモポリマーＡを、ビーズＡを用いてモデル化したビーズスプリングモデルＡが設定される。次に、ホモポリマーＢを、ビーズＢを用いてモデル化したビーズスプリングモデルＢが設定される。ビーズスプリングモデルＡのビーズＡ、Ａ間、及び、ビーズスプリングモデルＢのビーズＢ、Ｂ間には、相互作用ポテンシャルがそれぞれ定義されている。

次に、下記非特許文献１の方法では、予め定められた単位セルに、複数のビーズスプリングモデルＡ、Ｂが配置される。単位セルには、ビーズスプリングモデルＡとビーズスプリングモデルＢとがそれぞれ別の領域に配置され、それらの領域間で界面が形成されている。そして、下記非特許文献１の方法では、単位セルの変形が計算され、ポリマーブレンドの破壊エネルギー等が計算される。

Takeshi Aoyagi著、「Coarse-Grained Molecular Dynamics Study」、Nihon Reoroji Gakkaishi Vol.37、No.2、p75-79、2009年

上記非特許文献１の方法では、ポリマーモデルについて引張試験を行った場合、主として、ビーズスプリングモデルＡが配置された領域と、複数のビーズスプリングモデルＢが配置された領域との界面で破壊が生じる現象が計算される。しかしながら、実際のポリマーブレンドを変形させると、ホモポリマーＡ、Ｂのうち、耐変形性を示す物性値（例えば、弾性率）が相対的に小さいポリマー内で破壊が観察される。したがって、上記非特許文献１の方法では、ポリマーブレンドの破壊特性を精度よく評価するには、さらなる改善の余地があった。

発明者らは、鋭意研究を重ねた結果、ホモポリマーＡ、Ｂの物性値の差異を考慮して、ビーズＡ、Ａ間の相互作用ポテンシャル、及び、ビーズＢ、Ｂ間の相互作用ポテンシャルを互いに異なるものとして定義することが、ポリマーブレンドの破壊特性を評価する上で有効であることを知見した。

本発明は、以上のような実状に鑑み案出されたもので、ゴム材料の破壊特性をより精度良く評価することが可能なゴム材料のシミュレーション方法などを提供することを主たる目的としている。

本発明は、予め定められた物性値について差異を有する第１高分子成分と第２高分子成分とを少なくとも含むゴム材料の破壊現象を、コンピュータを用いてシミュレーションするための方法であって、前記第１高分子成分の分子鎖構造を複数の第１粗視化粒子を用いて表現した第１高分子鎖モデルを前記コンピュータに入力する工程と、前記第２高分子成分の分子鎖構造を複数の第２粗視化粒子を用いて表現した第２高分子鎖モデルを前記コンピュータに入力する工程と、前記第１高分子鎖モデルと前記第２高分子鎖モデルとの間に界面が形成されるように、前記第１高分子鎖モデル及び前記第２高分子鎖モデルを、前記ゴム材料の一部に対応する仮想空間であるセルに配置したゴム材料モデルを定義する工程と、前記第１粗視化粒子間の相互作用を定義するための第１相互作用ポテンシャルを定義する工程と、前記第２粗視化粒子間の相互作用を定義するための第２相互作用ポテンシャルを、前記物性値の差異に基づいて、前記第１相互作用ポテンシャルとは異なるものとして定義する工程と、前記コンピュータが、前記ゴム材料モデルの変形を計算して、前記第１相互作用ポテンシャル、及び、前記第２相互作用ポテンシャルに基づく前記界面での破壊現象を計算する工程とを含むことを特徴とする。

本発明に係る前記ゴム材料のシミュレーション方法において、前記第１相互作用ポテンシャル及び前記第２相互作用ポテンシャルは、相互作用の強さの定義するための強度パラメータを含み、前記第２相互作用ポテンシャルを定義する工程は、前記第２相互作用ポテンシャルに、前記第１相互作用ポテンシャルとは異なる前記強度パラメータを定義する工程を含んでもよい。

本発明に係る前記ゴム材料のシミュレーション方法において、前記第１相互作用ポテンシャル及び前記第２相互作用ポテンシャルは、以下の式で定義されてもよい。
ここで、
ｒ_ij：粒子間の距離
ｒ_ij,_cut：カットオフ距離
ε_ij：強度パラメータ
σ_ij：粒子の直径に相当

本発明に係る前記ゴム材料のシミュレーション方法において、前記強度パラメータを定義する工程は、予め定められた複数種類の前記強度パラメータε_ij毎に、前記物性値を計算する工程と、前記物性値の差異に基づいて、前記複数種類の前記強度パラメータε_ijから、前記第２相互作用ポテンシャルの前記強度パラメータε_ijを選択する工程とを含んでもよい。

本発明に係る前記ゴム材料のシミュレーション方法において、前記物性値は、弾性率を含んでもよい。

本発明の製造方法は、請求項１ないし５のいずれかに記載の前記シミュレーション方法において、前記性能が良好であると評価された前記ゴム材料モデルに基づいて、前記ゴム材料を製造する工程を含むことを特徴とする。

本発明のゴム材料のシミュレーション方法は、第１高分子成分の分子鎖構造を複数の第１粗視化粒子を用いて表現した第１高分子鎖モデルをコンピュータに入力する工程と、第２高分子成分の分子鎖構造を複数の第２粗視化粒子を用いて表現した第２高分子鎖モデルを前記コンピュータに入力する工程とを含んでいる。さらに、本発明の前記シミュレーションは、前記第１高分子鎖モデルと前記第２高分子鎖モデルとの間に界面が形成されるように、前記第１高分子鎖モデル及び前記第２高分子鎖モデルを、前記ゴム材料の一部に対応する仮想空間であるセルに配置したゴム材料モデルを定義する工程を含んでいる。

本発明の前記シミュレーション方法は、前記第１粗視化粒子間の相互作用を定義するための第１相互作用ポテンシャルを定義する工程と、前記第２粗視化粒子間の相互作用を定義するための第２相互作用ポテンシャルを、前記第１高分子成分と前記第２高分子成分との物性値の差異に基づいて、前記第１相互作用ポテンシャルとは異なるものとして定義する工程とを含んでいる。これにより、本発明の前記シミュレーション方法は、前記物性値の差異を考慮することが可能な前記ゴム材料モデルを作成することができる。

本発明の前記シミュレーション方法は、前記コンピュータが、前記ゴム材料モデルの変形を計算して、前記第１相互作用ポテンシャル、及び、前記第２相互作用ポテンシャルに基づく前記界面での破壊現象を計算する工程を含んでいる。これにより、本発明の前記シミュレーション方法は、前記第１高分子成分と前記第２高分子成分とを含むゴム材料の破壊特性を、精度良く評価することができる。

ゴム材料のシミュレーション方法を実行するためのコンピュータの一例を示す斜視図である。 ゴム材料のシミュレーション方法及び製造方法の処理手順の一例を示すフローチャートである。 ゴム材料モデル入力工程の処理手順の一例を示すフローチャートである。 第１高分子鎖モデル、第２高分子鎖モデル、及び、架橋剤粒子モデルが配置されたセルの一例を示す概念図である。 第１高分子鎖モデル及び架橋剤粒子モデルの一例を示す概念図である。 分布計算工程の処理手順の一例を示すフローチャートである。 第１高分子成分の分布、第２高分子成分の分布、及び、架橋剤の分布の一例を示すグラフである。 高分子成分の分布に基づいて第１高分子鎖モデル、及び、第２高分子鎖モデルが配置されたセルを部分的に拡大した概念図である。 （ａ）、（ｂ）は、一対の粗視化粒子を連結する工程の一例を説明する概念図である。 第２相互作用定義工程の処理手順の一例を示すフローチャートである。 強度パラメータ定義工程の処理手順の一例を示すフローチャートである。 物性値計算工程の処理手順の一例を示すフローチャートである。 高分子鎖モデルが配置されたセルの一例を示す概念図である。 強度パラメータε_ij毎に計算された応力−歪み曲線の一例を示すグラフである。 物性値の差異（ヤング率の比Ｅ／Ｅ_1.0）と、強度パラメータε_ijとの関係を示すグラフである。 （Ａ）は、ゴム材料モデル（歪：０．０）を示す図、（Ｂ）は、ゴム材料モデル（歪：０．１２）を示す図、（Ｃ）は、ゴム材料モデル（歪：０．３０）を示す図、（Ｄ）は、ゴム材料モデル（歪：０．６０）を示す図である。 （Ａ）は、比較例のゴム材料モデル（歪：０．０）を示す図、（Ｂ）は、比較例のゴム材料モデル（歪：０．１２）を示す図、（Ｃ）は、比較例のゴム材料モデル（歪：０．３０）を示す図、（Ｄ）は、比較例のゴム材料モデル（歪：０．６０）を示す図である。 実施例のゴム材料モデル、及び、比較例のゴム材料モデルの応力−歪み曲線を示すグラフである。

以下、本発明の実施の一形態が図面に基づき説明される。
本実施形態のゴム材料のシミュレーション方法（以下、単に「シミュレーション方法」ということがある。）では、ゴム材料の破壊現象が、コンピュータを用いてシミュレーション（計算）される。

図１は、ゴム材料のシミュレーション方法を実行するためのコンピュータの一例を示す斜視図である。コンピュータ１は、本体１ａ、キーボード１ｂ、マウス１ｃ、及び、ディスプレイ装置１ｄを含んでいる。本体１ａには、例えば、演算処理装置（ＣＰＵ）、ＲＯＭ、作業用メモリ、磁気ディスクなどの記憶装置、及び、ディスクドライブ装置１ａ１、１ａ２が設けられている。記憶装置には、本実施形態のシミュレーション方法を実行するためのソフトウェア等が予め記憶されている。

本実施形態のゴム材料は、第１高分子成分と第２高分子成分とを少なくとも含んでいる。第１高分子成分及び第２高分子成分としては、例えば、天然ゴム、ブタジエンゴム、ポリイソプレンゴム、又は、スチレンブタジエンゴム等が例示されるが、これらに限定されるわけではない。

第１高分子成分及び第２高分子成分は、例えば、上記の高分子成分の中から、互いに異なる高分子成分がそれぞれ選択される。本実施形態では、予め定められた物性値について差異を有するように、第１高分子成分及び第２高分子成分がそれぞれ選択される。物性値については、特に限定されない。本実施形態の物性値には、耐変形性を示す弾性率（ヤング率）が含まれる。本実施形態において、第１高分子成分がブタジエンゴムであり、第２高分子成分がスチレンブタジエンゴムである場合が例示される。

本実施形態のゴム材料には、架橋剤が含まれている。架橋剤としては、硫黄である場合が例示される。

図２は、ゴム材料のシミュレーション方法及び製造方法の処理手順の一例を示すフローチャートである。本実施形態のシミュレーション方法では、先ず、コンピュータ１に、ゴム材料モデルが入力される（ゴム材料モデル入力工程Ｓ１）。図３は、ゴム材料モデル入力工程Ｓ１の処理手順の一例を示すフローチャートである。図４は、第１高分子鎖モデル３、第２高分子鎖モデル４、及び、架橋剤粒子モデル５が配置されたセル１５の一例を示す概念図である。図５は、第１高分子鎖モデル３及び架橋剤粒子モデル５の一例を示す概念図である。図４は、一つの第１高分子鎖モデル３、一つの第２高分子鎖モデル４、及び、２つの架橋剤粒子モデル５を代表して示している。

本実施形態のゴム材料モデル入力工程Ｓ１では、先ず、図４及び図５に示されるように、第１高分子成分（本実施形態では、ブタジエンゴム）の分子鎖構造を、複数の第１粗視化粒子７を用いて表現した第１高分子鎖モデル３が、コンピュータ１に入力される（工程Ｓ１１）。本実施形態の第１高分子鎖モデル３は、粗視化モデル（本実施形態では、Kremer-Grestモデル）として定義されている。本実施形態の第１高分子鎖モデル３は、複数の第１粗視化粒子７と、隣接する第１粗視化粒子７、７間を結合する第１結合鎖モデル８とを含んで構成されている。

第１粗視化粒子７は、第１高分子成分（本実施形態では、ブタジエンゴム）の分子鎖のモノマー（図示省略）又はモノマーの一部分をなす構造単位を置換したものである。これにより、第１高分子鎖モデル３には、複数個（例えば、１０〜５０００個）の第１粗視化粒子７が設定される。第１粗視化粒子７は、後述の分子動力学計算において、運動方程式の質点として取り扱われる。このため、第１粗視化粒子７には、例えば、質量、体積、粒子径Ｄ１、又は、電荷などのパラメータが定義される。

本実施形態の第１結合鎖モデル８は、第１粗視化粒子７、７間に、伸びきり長が設定された第３相互作用ポテンシャルＰ３によって定義される。本実施形態の第３相互作用ポテンシャルＰ３は、下記式（２）の非調和ポテンシャルＵ^ch（ｒ）によって定義される。

ここで、各定数及び変数は、次のとおりである。
ｒ：粒子モデル間の距離
ｋ：粒子モデル間のばね定数
Ｒ₀：伸びきり長
なお、距離ｒ、及び、伸びきり長Ｒ₀は、第１粗視化粒子７の中心７ｃ（図５に示す）の座標に基づいて設定される。

非調和ポテンシャルＵ^ch（ｒ）の各定数及び各変数については、例えば、論文１（ Kurt Kremer & Gary S. Grest 著 「Dynamics of entangled linear polymer melts: A molecular-dynamics simulation」、J. Chem Phys. vol.92, No.8, 15 April 1990）に基づいて適宜設定される。これにより、工程Ｓ１１では、第１粗視化粒子７が伸縮自在に拘束された直鎖状の第１高分子鎖モデル３を定義することができる。第１高分子鎖モデル３は、コンピュータ１に記憶される。

次に、本実施形態のゴム材料モデル入力工程Ｓ１では、図４に示されるように、第２高分子成分（本実施形態では、スチレンブタジエンゴム）の分子鎖構造を、複数の第２粗視化粒子１１を用いて表現した第２高分子鎖モデル４が、コンピュータ１に入力される（工程Ｓ１２）。本実施形態の第２高分子鎖モデル４は、第１高分子鎖モデル３と同様に、粗視化モデルとして定義されている。本実施形態の第２高分子鎖モデル４は、複数の第２粗視化粒子１１と、隣接する第２粗視化粒子１１、１１間を結合する第２結合鎖モデル１２とを含んで構成されている。

第２粗視化粒子１１は、第２高分子成分（本実施形態では、スチレンブタジエンゴム）の分子鎖のモノマー（図示省略）又はモノマーの一部分をなす構造単位を置換したものである。これにより、第２高分子鎖モデル４には、複数個（例えば、１０〜５０００個）の第２粗視化粒子１１が設定される。第２粗視化粒子１１は、後述の分子動力学計算において、運動方程式の質点として取り扱われる。このため、第２粗視化粒子１１には、第１粗視化粒子７と同様のパラメータが定義される。

第２結合鎖モデル１２は、第２粗視化粒子１１、１１間に、伸びきり長が設定された第４相互作用ポテンシャルＰ４によって定義される。本実施形態の第４相互作用ポテンシャルＰ４は、上記式（２）の非調和ポテンシャルＵ^ch（ｒ）によって定義される。非調和ポテンシャルＵ^ch（ｒ）の各定数及び各変数については、第１高分子鎖モデル３の各定数及び各変数と同様に設定される。これにより、工程Ｓ１２では、第２粗視化粒子１１が伸縮自在に拘束された直鎖状の第２高分子鎖モデル４を定義することができる。第２高分子鎖モデル４は、コンピュータ１に記憶される。

次に、本実施形態のゴム材料モデル入力工程Ｓ１は、図４及び図５に示されるように、架橋剤をモデル化した架橋剤粒子モデル５が、コンピュータ１に入力される（工程Ｓ１３）。本実施形態の架橋剤粒子モデル５は、一つの独立した粒子モデルとして設定される。架橋剤粒子モデル５は、分子動力学計算において、運動方程式の質点として取り扱われる。このため、架橋剤粒子モデル５には、例えば、質量、体積、粒子径Ｄ２、又は、電荷などのパラメータが定義される。

架橋剤粒子モデル５の粒子径Ｄ２については、適宜設定することができる。本実施形態の粒子径Ｄ２は、第１粗視化粒子７及び第２粗視化粒子１１の粒子径Ｄ１と同一に設定されている。架橋剤粒子モデル５は、コンピュータ１に記憶される。

次に、本実施形態のゴム材料モデル入力工程Ｓ１では、図４に示されるように、セル１５が、コンピュータ１に入力される（工程Ｓ１４）。本実施形態のセル１５は、解析対象のゴム材料の一部に対応する仮想空間である。セル１５は、互いに向き合う三対の平面１５ａ、１５ｂを有する直方体として定義されている。各平面１５ａ、１５ｂには、周期境界条件が定義されている。このようなセル１５は、例えば、一方の平面１５ａから出て行った第１高分子鎖モデル３、第２高分子鎖モデル４、及び、架橋剤粒子モデル５の一部が、反対側の平面１５ｂから入ってくるように計算することができる。したがって、一方の平面１５ａと、反対側の平面１５ｂとが連続している（繋がっている）ものとして取り扱うことができる。

セル１５の一辺の各長さＬ１ａ、Ｌ１ｂ及びＬ１ｃは、適宜設定することができる。本実施形態の長さＬ１ａ、Ｌ１ｂ及びＬ１ｃは、第１高分子鎖モデル３、及び、第２高分子鎖モデル４の拡がりを示す量である慣性半径（図示省略）の２倍以上が望ましい。これにより、セル１５は、分子動力学計算において、周期境界条件による自己のイメージとの衝突の発生を防ぐことができるため、第１高分子鎖モデル３、及び、第２高分子鎖モデル４の空間的拡がりを適切に計算することができる。また、セル１５の大きさは、例えば１気圧で安定な体積に設定される。これにより、セル１５は、解析対象のゴム材料の少なくとも一部の体積を定義することができる。セル１５は、コンピュータ１に記憶される。

次に、本実施形態のゴム材料モデル入力工程Ｓ１では、コンピュータ１が、解析対象のゴム材料に含まれる成分（第１高分子成分、第２高分子成分、及び、架橋剤）の分布を求める（分布計算工程Ｓ１５）。図６は、分布計算工程Ｓ１５の処理手順の一例を示すフローチャートである。

本実施形態の分布計算工程Ｓ１５では、先ず、第１高分子成分、第２高分子成分、及び、架橋剤のχパラメータが求められる（工程Ｓ６１）。χパラメータは、Flory-Huggins理論に基づくものであり、第１高分子成分、第２高分子成分、及び、架橋剤の相溶性を表すものである。χパラメータの値が小さいほど、良溶媒であることを示している。第１高分子成分と第２高分子成分とによって形成される界面は、相溶性が高い（即ち、χパラメータが小さい）ほど、互いに混じり合った混合状態になりやすい。逆に、相溶性が低い（χパラメータが大きい）ほど、第１高分子成分及び第２高分子成分が互いに混じり合い難くなり、相分離した界面が形成される。

本実施形態の工程Ｓ６１では、ブタジエンゴム間の第１χパラメータ、スチレンブタジエンゴム間の第２χパラメータ、ブタジエンゴムとスチレンブタジエンゴムとの間の第３χパラメータが求められる。さらに、工程Ｓ６１では、架橋剤間の第４χパラメータ、ブタジエンゴムと架橋剤との間の第５χパラメータ、及び、スチレンブタジエンゴムと架橋剤との間の第６χパラメータが求められる。各χパラメータは、従来と同様に、（株）ＪＳＯＬ社製のソフトマテリアル総合シミュレーター（Ｊ−ＯＣＴＡ）に含まれるＳＵＳＨＩモデラを用いて、原子団寄与法、モンテカルロ法、分子動力学法、又は、ＱＳＰＲ法に基づいて計算することができる。各χパラメータは、コンピュータ１に記憶される。

次に、本実施形態の分布計算工程Ｓ１５では、解析対象のゴム材料に含まれる成分（第１高分子成分、第２高分子成分、及び、架橋剤）の分布が求められる（工程Ｓ６２）。本実施形態では、第１高分子成分、第２高分子成分、及び、架橋剤の配置のエントロピーを考慮した静的なＳＣＦ法（ Self-Consistent Field Method ）を用いて、ゴム材料に含まれる第１高分子成分（ブタジエンゴム）の分布、第２高分子成分（スチレンブタジエンゴム）の分布、及び、架橋剤の分布がそれぞれ求められる。ＳＣＦ法は、例えば、（株）ＪＳＯＬ社製のソフトマテリアル総合シミュレーター（Ｊ−ＯＣＴＡ）に含まれるＳＵＳＨＩを用いて計算することができる。

ＳＣＦ法の入力値としては、第１高分子成分（ブタジエンゴム）の分子量、第２高分子成分（スチレンブタジエンゴム）の分子量、架橋剤の分子量、及び、χパラメータ（本実施形態では、第１χパラメータ〜第６χパラメータ）が用いられる。これにより、工程Ｓ６２では、解析対象のゴム材料に含まれる第１高分子成分（ブタジエンゴム）の分布（濃度分布）、第２高分子成分（スチレンブタジエンゴム）の分布（濃度分布）、及び、架橋剤の分布（濃度分布）がそれぞれ求められる。図７は、第１高分子成分の分布、第２高分子成分の分布、及び、架橋剤の分布の一例を示すグラフである。このグラフでは、第１高分子成分、第２高分子成分及び架橋剤の体積分率と、それらのゴム材料での空間位置との関係を示している。ゴム材料に含まれる成分の分布は、コンピュータ１に入力される。

次に、本実施形態のゴム材料モデル入力工程Ｓ１では、第１高分子成分、第２高分子成分及び架橋剤の分布に基づいて、コンピュータ１が、第１高分子鎖モデル３、第２高分子鎖モデル４、及び、架橋剤粒子モデル５をセル１５内に配置する（工程Ｓ１６）。本実施形態の工程Ｓ１６では、ＤＢＭＣ法（Density Biased Monte Carlo）を用いて、複数の第１高分子鎖モデル３、複数の第２高分子鎖モデル４、及び、複数の架橋剤粒子モデル５が、セル１５内に配置される。ＤＢＭＣ法は、例えば、上記したソフトマテリアル総合シミュレーター（Ｊ−ＯＣＴＡ）に含まれるＣＯＧＮＡＣを用いて計算することができる。

ＤＢＭＣ法の入力値としては、工程Ｓ６２で求められた第１高分子成分（ブタジエンゴム）の分布（濃度分布）、第２高分子成分（スチレンブタジエンゴム）の分布（濃度分布）、及び、架橋剤の分布（濃度分布）が用いられる。これにより、本実施形態の工程Ｓ１６では、解析対象のゴム材料に含まれる第１高分子成分の分布、第２高分子成分の分布、及び、架橋剤の分布（図７に示す）に基づいて、複数の第１高分子鎖モデル３、複数の第２高分子鎖モデル４、及び、複数の架橋剤粒子モデル５を、セル１５内に配置することができる。図８は、高分子成分の分布に基づいて第１高分子鎖モデル３、及び、第２高分子鎖モデル４が配置されたセル１５を部分的に拡大した概念図である。図８において、図４に示した架橋剤粒子モデル５、第１結合鎖モデル８及び第２結合鎖モデル１２を省略するとともに、第２高分子鎖モデル４を着色して表示している。

図８に示されるように、工程Ｓ１６では、図７に示した各高分子成分の分布に基づいて、複数の第１高分子鎖モデル３（ブタジエンゴム）と、複数の第２高分子鎖モデル４（スチレンブタジエンゴム）との界面構造が作成されたゴム材料モデル１６が定義される。

次に、本実施形態のゴム材料モデル入力工程Ｓ１では、図４に示されるように、第１高分子鎖モデル３の第１粗視化粒子７、及び、第２高分子鎖モデル４の第２粗視化粒子１１に、引力及び斥力が作用する相互作用ポテンシャルが定義される（工程Ｓ１７）。本実施形態の相互作用ポテンシャルは、一対の第１高分子鎖モデル３、３の第１粗視化粒子７、７間に定義される第５相互作用ポテンシャルＰ５（図示省略）、一対の第２高分子鎖モデル４、４の第２粗視化粒子１１、１１間に定義される第６相互作用ポテンシャルＰ６（図示省略）、及び、第１高分子鎖モデル３の第１粗視化粒子７と第２高分子鎖モデル４の第２粗視化粒子１１との間に定義される第７相互作用ポテンシャルＰ７を含んでいる。第５相互作用ポテンシャルＰ５〜第７相互作用ポテンシャルＰ７は、例えば、下記式（１）のＬＪポテンシャルＵ^LJ _ij（ｒ）で定義することができる。

ここで、
ｒ_ij：粒子間の距離
ｒ_ij,_cut：カットオフ距離
ε_ij：強度パラメータ
σ_ij：粒子の直径に相当
なお、距離ｒ_ij、カットオフ距離ｒ_ij,_cut、及び、粒子の直径に相当するσ_ijは、第１粗視化粒子７の中心７ｃ（図５に示す）の座標に基づいて定義される。

第５相互作用ポテンシャルＰ５〜第７相互作用ポテンシャルＰ７の各定数及び各変数については、例えば、論文２（Ali Makke, Michel Perez, Olivier Lame and Jean-Louis Barrat著、 「Mechanical testing of glassy and rubbery polymers in numerical simulations: Role of boundary conditions in tensile stress experiments 」J. Chem. Phys. vol.131, 014904 (2009)）に基づいて適宜設定される。なお、第５相互作用ポテンシャルＰ５〜第７相互作用ポテンシャルＰ７の各強度パラメータε_ijについては、適宜設定することができる。強度パラメータε_ijは、ブタジエンゴム及びスチレンブタジエンゴムのχパラメータと相関がある。したがって、第５相互作用ポテンシャルＰ５〜第７相互作用ポテンシャルＰ７の各強度パラメータε_ijは、ブタジエンゴム及びスチレンブタジエンゴムのχパラメータに基づいて同定されるのが望ましい。本実施形態の第５相互作用ポテンシャルＰ５〜第７相互作用ポテンシャルＰ７の各強度パラメータε_ijには、例えば、１．０に設定される。第５相互作用ポテンシャルＰ５〜第７相互作用ポテンシャルＰ７は、コンピュータ１に記憶される。

次に、本実施形態のゴム材料モデル入力工程Ｓ１では、架橋剤粒子モデル５に、引力及び斥力が作用する相互作用ポテンシャルが定義される（工程Ｓ１８）。本実施形態の相互作用ポテンシャルは、一対の架橋剤粒子モデル５、５間に定義される第８相互作用ポテンシャルＰ８、架橋剤粒子モデル５と第１高分子鎖モデル３の第１粗視化粒子７との間に定義される第９相互作用ポテンシャルＰ９、及び、架橋剤粒子モデル５と第２高分子鎖モデル４の第２粗視化粒子１１との間に定義される第１０相互作用ポテンシャルＰ１０を含んでいる。第８相互作用ポテンシャルＰ８〜第１０相互作用ポテンシャルＰ１０は、例えば、上記式（１）のＬＪポテンシャルＵ^LJ _ij（ｒ）を用いて定義することができる。

第８相互作用ポテンシャルＰ８〜第１０相互作用ポテンシャルＰ１０の強度パラメータε_ijについては、適宜設定することができる。強度パラメータε_ijは、架橋剤のχパラメータと相関がある。従って、強度パラメータε_ijは、架橋剤のχパラメータに基づいて同定されるのが望ましい。第８相互作用ポテンシャルＰ８〜第１０相互作用ポテンシャルＰ１０は、コンピュータ１に記憶される。

次に、本実施形態のゴム材料モデル入力工程Ｓ１では、コンピュータ１が、分子動力学計算に基づいて、第１高分子鎖モデル３及び第２高分子鎖モデル４を対象に構造緩和を計算する（工程Ｓ１９）。本実施形態の分子動力学計算では、例えば、セル１５について所定の時間、第１高分子鎖モデル３の第１粗視化粒子７、及び、第２高分子鎖モデル４の第２粗視化粒子１１が古典力学に従うものとして、ニュートンの運動方程式が適用される。そして、各時刻での第１粗視化粒子７及び第２粗視化粒子１１の動きが、単位時間ステップ毎に追跡される。

分子動力学計算では、セル１５において、圧力及び温度が一定、又は体積及び温度が一定に保たれる。これにより、工程Ｓ１９では、実際のゴム材料の分子運動に近似させて、第１高分子鎖モデル３、及び、第２高分子鎖モデル４の初期配置を精度よく緩和することができる。本実施形態では、第１高分子鎖モデル３、及び、第２高分子鎖モデル４の初期配置が十分に緩和されるまで行われる。構造緩和の計算は、例えば、上記Ｊ−ＯＣＴＡに含まれるＣＯＧＮＡＣを用いて処理することができる。

本実施形態の工程Ｓ１９では、第１粗視化粒子７、７間に定義される第３相互作用ポテンシャルＰ３（図５に示す）、及び、第２粗視化粒子１１、１１間に定義される第４相互作用ポテンシャルＰ４（図５に示す）が、非調和ポテンシャルＵ^ch（ｒ）として定義されている。これにより、工程Ｓ１９では、第１粗視化粒子７、７間、及び、第２粗視化粒子１１、１１間で、柔軟に伸縮することができるため、構造緩和をスムーズに行うことができる。さらに、本実施形態の第３相互作用ポテンシャルＰ３及び第４相互作用ポテンシャルＰ４の各定数及び各変数には、同一のものが採用されている。このため、セル１５において、第１高分子鎖モデル３の配置、及び、第２高分子鎖モデル４の配置に偏りが生じるのを防ぐことができる。

本実施形態の工程Ｓ１９では、架橋剤粒子モデル５を対象とした構造緩和が計算されていない。これにより、架橋剤粒子モデル５は、セル１５において、実際の架橋剤の分布を維持することができる。

次に、本実施形態のゴム材料モデル入力工程Ｓ１は、コンピュータ１が、隣接する分子鎖モデル（第１高分子鎖モデル３又は第２高分子鎖モデル４）について、一対の粗視化粒子（第１粗視化粒子７又は第２粗視化粒子１１）を連結する（工程Ｓ２０）。本実施形態の工程Ｓ２０では、一対の粗視化粒子（第１粗視化粒子７又は第２粗視化粒子１１）を、架橋剤粒子モデル５を介して連結している。図９（ａ）、（ｂ）は、一対の粗視化粒子を連結する工程の一例を説明する概念図である。

図９（ａ）、（ｂ）に示されるように、工程Ｓ２０では、架橋剤粒子モデル５の位置から予め定められた範囲Ｒ（図９（ａ）に示す）内に存在する一対の粗視化粒子（本例では、第１粗視化粒子７、７）を連結している。範囲Ｒについては、適宜設定することができる。本実施形態の範囲Ｒは、例えば、架橋剤粒子モデル５の中心５ｃからの半径ｒｄが、０．８σ以内の領域に設定されている。

工程Ｓ２０では、図９（ｂ）に示されるように、一対の粗視化粒子（本例では、第１粗視化粒子７、７）が、架橋剤粒子モデル５を介して連結される。本実施形態では、一方の粗視化粒子（第１粗視化粒子７又は第２粗視化粒子１１）と架橋剤粒子モデル５との間、及び、他方の粗視化粒子（第１粗視化粒子７又は第２粗視化粒子１１）と架橋剤粒子モデル５との間が、第３結合鎖モデル１８を介して接続される。

第３結合鎖モデル１８は、第１１相互作用ポテンシャルＰ１１によって設定される。本実施形態の第１１相互作用ポテンシャルＰ１１は、上記式（２）の非調和ポテンシャルＵ^ch（ｒ）によって定義される。また、上記式（２）において、第１１相互作用ポテンシャルＰ１１の定数及び変数は、例えば、図４に示した第３相互作用ポテンシャルＰ３及び第４相互作用ポテンシャルＰ４の定数又は変数と同様に定義されうる。

このように、工程Ｓ２０では、隣接する高分子鎖モデル（本実施形態では、第１高分子鎖モデル３又は第２高分子鎖モデル４）の一対の粗視化粒子（本実施形態では、第１粗視化粒子７又は第２粗視化粒子１１）が、架橋剤粒子モデル５を介して連結される。これにより、解析対象の架橋されたゴム材料を再現したゴム材料モデル１６（図４に示す）を作成することができる。ゴム材料モデル１６は、コンピュータ１に記憶される。

図４に示されるように、第１高分子鎖モデル３及び第２高分子鎖モデル４について、第１粗視化粒子７、７間に定義される第３相互作用ポテンシャルＰ３、及び、第２粗視化粒子１１、１１間に定義される第４相互作用ポテンシャルＰ４は、非調和ポテンシャルＵ^ch（ｒ）によって定義されている。さらに、第３相互作用ポテンシャルＰ３及び第４相互作用ポテンシャルＰ４の各定数及び各変数には、同一の値が設定されている。このため、第１高分子鎖モデル３及び第２高分子鎖モデル４は、第１高分子成分及び第２高分子成分の物性値の差異が考慮されていない。したがって、工程Ｓ１１〜工程Ｓ２０を経て作成された相分離界面（即ち、第１高分子成分及び第２高分子成分が互いに分離した界面）のゴム材料モデル１６は、例えば、後述の変形（引張試験）の計算に用いられると、上記非特許文献１と同様に、第１高分子鎖モデル３が配置された領域と、第２高分子鎖モデル４が配置された領域との界面で破壊が生じる現象が計算される。

しかしながら、第１高分子成分及び第２高分子成分を含む相分離した実際のゴム材料を変形させると、第１高分子成分及び第２高分子成分のうち、耐変形性を示す物性値（例えば、弾性率）が相対的に小さい第１高分子成分又は第２高分子成分の一方（本実施形態では、第１高分子成分）内で破壊が観察される。

発明者らは、鋭意研究を重ねた結果、第１高分子成分及び第２高分子成分の物性値の差異を考慮して、第１粗視化粒子７、７間の相互作用、及び、第２粗視化粒子１１、１１間の相互作用を互いに異なるものとして定義することが、第１高分子成分及び第２高分子成分を含む実際のゴム材料の破壊特性を評価する上で有効であることを知見した。

このような知見に基づいて、本実施形態のシミュレーション方法（ゴム材料モデル入力工程Ｓ１）では、第１高分子成分と第２高分子成分との物性値の差異に基づいて、第２粗視化粒子１１、１１間の相互作用を定義するための第２相互作用ポテンシャルＰ２を、第１粗視化粒子７、７間の相互作用を定義するための第１相互作用ポテンシャルＰ１とは異なるものとして定義している。

本実施形態のゴム材料モデル入力工程Ｓ１では、第１粗視化粒子７、７間の相互作用を定義するための第１相互作用ポテンシャルＰ１（図４に示す）が定義される（工程Ｓ２１）。第１相互作用ポテンシャルＰ１は、相互作用の強さの定義するための強度パラメータε_ijを含み、上記式（１）のＬＪポテンシャルによって定義される。本実施形態の工程Ｓ２１では、上記式（２）で定義された第１粗視化粒子７、７間の第３相互作用ポテンシャルＰ３（図４に示す）を無効にし、上記式（１）で定義される第１相互作用ポテンシャルＰ１（図４に示す）が定義される。

第１相互作用ポテンシャルＰ１の各定数及び各変数については、例えば、上記論文２に基づいて適宜設定される。本実施形態の第１相互作用ポテンシャルＰ１の強度パラメータε_ijは、１．０に設定される。第１相互作用ポテンシャルＰ１は、コンピュータ１に記憶される。

次に、本実施形態のゴム材料モデル入力工程Ｓ１では、第２粗視化粒子１１、１１間の相互作用を定義するための第２相互作用ポテンシャルＰ２（図４に示す）が定義される（第２相互作用定義工程Ｓ２２）。第２相互作用ポテンシャルＰ２は、相互作用の強さの定義するための強度パラメータε_ijを含み、上記式（１）のＬＪポテンシャルによって定義される。第２相互作用定義工程Ｓ２２では、第１高分子成分と第２高分子成分との物性値の差異に基づいて、第１相互作用ポテンシャルＰ１（図４に示す）とは異なるものとして、第２相互作用ポテンシャルＰ２を定義している。図１０は、第２相互作用定義工程Ｓ２２の処理手順の一例を示すフローチャートである。

本実施形態の第２相互作用定義工程Ｓ２２では、先ず、上記式（２）で定義された第２粗視化粒子１１、１１間の第４相互作用ポテンシャルＰ４（図４に示す）を無効にする工程Ｓ３１、及び、第２粗視化粒子１１、１１間に、上記式（１）で定義される第２相互作用ポテンシャルＰ２（図４に示す）を定義する工程Ｓ３２が実施される。

本実施形態の工程Ｓ３２において、第２相互作用ポテンシャルＰ２の各定数及び各変数については、例えば、上記論文２に基づいて適宜設定される。なお、強度パラメータε_ijは、後述の強度パラメータ定義工程Ｓ３３において、第１相互作用ポテンシャルＰ１とは異なる強度パラメータε_ijが定義される。

次に、本実施形態の第２相互作用定義工程Ｓ２２では、第２相互作用ポテンシャルＰ２に、第１相互作用ポテンシャルＰ１とは異なる強度パラメータε_ijが定義される（強度パラメータ定義工程Ｓ３３）。図１１は、強度パラメータ定義工程Ｓ３３の処理手順の一例を示すフローチャートである。

本実施形態の強度パラメータ定義工程Ｓ３３では、先ず、予め定められた複数種類の強度パラメータε_ij毎に、物性値を計算する（物性値計算工程Ｓ４１）。図１２は、物性値計算工程Ｓ４１の処理手順の一例を示すフローチャートである。図１３は、高分子鎖モデル１７が配置されたセル２１の一例を示す概念図である。

本実施形態の物性値計算工程Ｓ４１では、先ず、図１３に示されるように、物性値を計算するための高分子鎖モデル１７が、コンピュータ１に入力される（工程Ｓ５１）。高分子鎖モデル１７は、図４に示した第１高分子鎖モデル３及び第２高分子鎖モデル４と同様に、複数の粗視化粒子１９と、隣接する粗視化粒子１９、１９間を結合する結合鎖モデル２０とを含んで構成されている。高分子鎖モデル１７には、図４に示した第１高分子鎖モデル３又は第２高分子鎖モデル４がそのまま用いられてもよいし、第１高分子鎖モデル３及び第２高分子鎖モデル４とは独立して定義されてもよい。

粗視化粒子１９は、後述の分子動力学計算において、運動方程式の質点として取り扱われる。このため、粗視化粒子１９には、図４に示した第１粗視化粒子７及び第２粗視化粒子１１と同様のパラメータが定義される。

結合鎖モデル２０は、粗視化粒子１９、１９間に設定される第１５相互作用ポテンシャルＰ１５によって定義される。第１５相互作用ポテンシャルＰ１５は、図４に示した第１相互作用ポテンシャルＰ１及び第２相互作用ポテンシャルＰ２と同様に、上記式（１）のＬＪポテンシャルによって定義される。第１５相互作用ポテンシャルＰ１５の各定数及び各変数については、例えば、強度パラメータε_ijを除いて、例えば、第１相互作用ポテンシャルＰ１及び第２相互作用ポテンシャルＰ２と同一のものを設定することができる。高分子鎖モデル１７は、コンピュータ１に記憶される。

次に、本実施形態の物性値計算工程Ｓ４１では、セル２１が、コンピュータ１に入力される（工程Ｓ５２）。本実施形態のセル２１は、図４に示したゴム材料モデル１６のセル１５とは独立して定義される。セル２１は、ゴム材料モデル１６のセル１５と同一の大きさに設定されても良いし、異なる大きさに設定されてもよい。セル２１は、コンピュータ１に記憶される。

次に、本実施形態の物性値計算工程Ｓ４１では、複数の高分子鎖モデル１７がセル２１内に配置され（工程Ｓ５３）、高分子鎖モデル１７の粗視化粒子１９に、引力及び斥力が作用する第１６相互作用ポテンシャルＰ１６が定義される（工程Ｓ５４）。本実施形態の第１６相互作用ポテンシャルＰ１６は、一対の高分子鎖モデル１７、１７の粗視化粒子１９、１９間にそれぞれ定義される。第１６相互作用ポテンシャルＰ１６は、第５相互作用ポテンシャルＰ５〜第７相互作用ポテンシャルＰ７（図４に示す）と同様に、上記式（１）のＬＪポテンシャルＵ^LJ _ij（ｒ）で定義することができる。

第１６相互作用ポテンシャルＰ１６の各定数及び各変数については、例えば、上記論文２に基づいて適宜設定される。本実施形態の第１６相互作用ポテンシャルＰ１６の強度パラメータε_ijは、１．０に設定される。第１６相互作用ポテンシャルＰ１６は、コンピュータ１に記憶される。

次に、本実施形態の物性値計算工程Ｓ４１では、コンピュータ１が、分子動力学計算に基づいて、高分子鎖モデル１７を対象に構造緩和を計算する（工程Ｓ５５）。分子動力学計算については、上述のとおりである。本実施形態では、高分子鎖モデル１７の初期配置が十分に緩和されるまで行われる。これにより、単一の高分子成分を含むゴム材料をモデル化したゴム材料モデル（以下、単に「単一ゴム材料モデル」ということがある。）２２が定義される

次に、本実施形態の物性値計算工程Ｓ４１は、コンピュータ１が、予め定められた複数種類の強度パラメータε_ij毎に、物性値を計算する（工程Ｓ５６）。工程Ｓ５６では、高分子鎖モデル１７の第１５相互作用ポテンシャルＰ１５に設定される複数種類の強度パラメータε_ij毎に、単一ゴム材料モデル２２の伸長計算をして、物性値が計算される。

複数種類の強度パラメータε_ijについては、例えば、第１高分子成分と第２高分子成分との物性値の差異に基づいて、適宜設定することができる。強度パラメータε_ijの一例としては、第１高分子成分と第２高分子成分との物性値に差異が無い状態の強度パラメータε_ij（例えば、１．０）から、第１高分子成分と第２高分子成分との物性値の差異が大きくなる強度パラメータε_ij（例えば、２．４）の間で複数設定（例えば、ε_ij＝１．０、１．２、１．４、１．６、１．８、２．０、２．２、及び、２．４）されるのが望ましい。なお、強度パラメータε_ijは、このような態様に限定されない。

工程Ｓ５６では、先ず、複数種類の強度パラメータε_ijから選択された一つの強度パラメータε_ijが、高分子鎖モデル１７の第１５相互作用ポテンシャルＰ１５に定義される。次に、工程Ｓ５６では、単軸引張り試験に基づいて、単一ゴム材料モデル２２を一方向に（例えば、Ｙ軸方向に０％〜５０％）伸長させる変形シミュレーションが実施され、その物性値（本実施形態では、弾性率）が計算される。変形シミュレーションについては、後述のゴム材料モデル１６を用いたゴム材料のシミュレーション方法と同様の手順で行うことができる。

工程Ｓ５６において、物性値の計算は、設定された複数の強度パラメータε_ij毎に実施される。図１４は、強度パラメータε_ij毎に計算された応力−歪み曲線の一例を示すグラフである。図１４では、強度パラメータε_ijが１．０、１．６、及び、２．４のみを代表して示している。図１５は、物性値の差異（ヤング率の比Ｅ／Ｅ_1.0）と、強度パラメータε_ijとの関係を示すグラフである。ヤング率の比Ｅ／Ｅ_1.0は、強度パラメータε_ij＝１．０のヤング率Ｅ_1.0と他の強度パラメータε_ijのヤング率Ｅとの割合を示している。複数種類の強度パラメータε_ij毎に計算された物性値は、コンピュータ１に記憶される。

次に、図１１に示されるように、本実施形態の強度パラメータ定義工程Ｓ３３は、物性値の差異に基づいて、複数種類の強度パラメータε_ij（例えば、ε_ij＝１．０、１．２、１．４、１．６、１．８、２．０、２．２、及び、２．４）から、第２相互作用ポテンシャルＰ２の強度パラメータε_ijを選択する（工程Ｓ４２）。工程Ｓ４２では、先ず、実際の第１高分子成分と第２高分子成分との物性値の差異が求められる。物性値（弾性率）の差異としては、第１高分子成分のヤング率Ｅ₁と第２高分子成分のヤング率Ｅ₂との比Ｅ₂／Ｅ₁が求められる。第１高分子成分のヤング率Ｅ₁と第２高分子成分のヤング率Ｅ₂は、例えば、ＪＩＳＫ６２５１「加硫ゴム及び熱可塑性ゴム−引張特性の求め方」に記載の試験方法に準拠して引張特性を求め、歪みに対して応力が直線的に変化する範囲（即ち、線形応答が得られた範囲）内の引張応力と、これに対応する歪みとの比から求めることができる。

次に、本実施形態の工程Ｓ４２では、図１５のグラフにおいて、実際の第１高分子成分と第２高分子成分との物性値の差異（ヤング率の比Ｅ₂／Ｅ₁）に対応する強度パラメータε_ijが選択される。図１５のグラフにおいて、例えば、物性値の差異（ヤング率の比Ｅ₂／Ｅ₁）が２．５の場合、強度パラメータε_ij＝２．０が選択される。

次に、本実施形態の強度パラメータ定義工程Ｓ３３では、第２相互作用ポテンシャルＰ２（図４に示す）に、選択された強度パラメータ（本例では、ε_ij＝２．０）が定義される（工程Ｓ４３）。強度パラメータε_ij＝２．０は、第１相互作用ポテンシャルＰ１（図４に示す）の強度パラメータε_ij＝１．０とは異なるものである。これにより、第２相互作用定義工程Ｓ２２では、第２相互作用ポテンシャルＰ２を、第１高分子成分と第２高分子成分との物性値の差異（本実施形態では、ヤング率の比Ｅ₂／Ｅ₁）に基づいて、第１相互作用ポテンシャルＰ１とは異なるものとして定義することができる。したがって、本実施形態のシミュレーション方法（ゴム材料モデル入力工程Ｓ１）は、第１高分子成分と第２高分子成分との物性値の差異を考慮することが可能なゴム材料モデル１６（図４に示す）を作成することができる。ゴム材料モデル１６は、コンピュータ１に記憶される。

本実施形態のゴム材料モデル入力工程Ｓ１では、第２相互作用定義工程Ｓ２２の後に、コンピュータ１が、分子動力学計算に基づいて、図４に示した第１高分子鎖モデル３及び第２高分子鎖モデル４を対象に構造緩和を計算するのが望ましい。これにより、第２相互作用定義工程Ｓ２２で設定される第２相互作用ポテンシャルＰ２の変化によって、分子同士の凝集度合い、セル２１のサイズ、セル２１内の密度（即ち、セル２１内の粒子数／セル２１のサイズ）が変化しても、新たな第２相互作用ポテンシャルＰ２に基づいて平衡化された構造を作成することができるため、セル２１を適切な密度にすることができる。本実施形態では、NPT（N：粒子数、P：圧力、T：温度）一定の条件下で構造緩和が行われるのが望ましく、また、セル１５の圧力がゼロ（又は、ゼロに近い任意の値）になるように、構造緩和が計算されるのが望ましい。これにより、ポテンシャルの異なる構造間で初期応力が異なることを避けることができ、界面状態の異なる構造間での応力比較が可能になる。

本実施形態のゴム材料モデル入力工程Ｓ１では、架橋剤をモデル化した架橋剤粒子モデル５が、セル１５に配置されたが、このような態様に限定されない。本実施形態のゴム材料モデル入力工程Ｓ１では、架橋剤粒子モデル５（図４及び図５に示す）を省略してもよい。この場合、工程Ｓ２０では、架橋剤粒子モデル５を介さずに、一対の粗視化粒子（第１粗視化粒子７又は第２粗視化粒子１１）を連結するのが望ましい。また、ゴム材料に配合される充填剤（フィラー）をモデル化したフィラーモデル（図示省略）が、セル１５に配置されてもよい。これにより、充填剤が配合されたゴム材料を再現したゴム材料モデル１６を作成することができる。

次に、本実施形態のシミュレーション方法では、コンピュータ１が、図４に示したゴム材料モデル１６の変形を計算（以下、単に、「変形計算」ということがある。）する（工程Ｓ２）。ゴム材料モデル１６の変形計算としては、例えば、引張変形や圧縮変形など適宜選択することができる。本実施形態の変形計算は、引張変形が行われる。

工程Ｓ２では、Ｙ軸方向において、ゴム材料モデル１６の一端１６ａ及び他端１６ｂを離間させて、ゴム材料モデル１６の引張変形が計算される。本実施形態では、例えば、特開２０１６−８１２９７号公報に記載の方法と同様に、分子動力学計算の単位時間（ＭＤステップ）あたりの歪が１０^−５以上の速度Ｖ１で、ゴム材料モデル１６をアフィン変形に基づいて引張変形させている。また、工程Ｓ２では、第１高分子鎖モデル３、第２高分子鎖モデル４、及び、架橋剤粒子モデル５の熱運動（分子運動）が計算される。

次に、本実施形態のシミュレーション方法では、ゴム材料モデル１６の応力が計算される（工程Ｓ３）。工程Ｓ３では、工程Ｓ２で与えられた歪に基づいて、ゴム材料モデル１６の応力が単位時間毎に計算される。さらに、工程Ｓ３では、ゴム材料モデル１６の引張変形によって計算される第１粗視化粒子７、第２粗視化粒子１１、及び、架橋剤粒子モデル５の動きが、単位時間毎に追跡される。図１６（Ａ）は、ゴム材料モデル（歪：０．０）を示す図である。図１６（Ｂ）は、ゴム材料モデル（歪：０．１２）を示す図である。図１６（Ｃ）は、ゴム材料モデル（歪：０．３０）を示す図である。図１６（Ｄ）は、ゴム材料モデル（歪：０．６０）を示す図である。

工程Ｓ３では、ゴム材料モデル１６に与えられた歪、各相互作用ポテンシャル及び運動方程式に基づいて、第１高分子鎖モデル３、第２高分子鎖モデル４、及び、架橋剤粒子モデル５の熱運動が計算される。これらの熱運動が計算されることにより、セル１５には、第１高分子鎖モデル３、第２高分子鎖モデル４、及び、架橋剤粒子モデル５（図示省略）が存在しない小領域が形成される。このような小領域は、ゴム材料モデル１６に形成された空孔（ボイド）２４として定義される。

図１６（Ａ）〜（Ｄ）に示されるように、本実施形態の変形計算では、第２高分子鎖モデル４よりも小さい強度パラメータε_ijが定義された第１高分子鎖モデル３において、破壊が生じる（空孔２４が生じる）現象が計算された。これは、実際のゴム材料の破壊現象に合致している。したがって、本実施形態のシミュレーション方法では、ゴム材料モデル入力工程Ｓ１で作成されたゴム材料モデル１６が、変形計算を行うシミュレーションに用いられることにより、第１高分子成分と第２高分子成分とを含むゴム材料の破壊特性を、精度良く評価することができる。ゴム材料モデル１６の応力、第１高分子鎖モデル３、第２高分子鎖モデル４、及び、架橋剤粒子モデル５の各座標値は、コンピュータ１に記憶される。

次に、本実施形態のシミュレーション方法では、ゴム材料モデル１６の変形の計算結果に基づいて、ゴム材料の性能が評価される（工程Ｓ４）。本実施形態の工程Ｓ４では、例えば、ゴム材料モデル１６の変形結果に基づいて計算される応力−歪み曲線や、空孔（ボイド）２４の成長曲線等に基づいて、ゴム材料の強度が評価される。

本実施形態では、ゴム材料の性能（変形の計算結果）が良好であると評価された場合（工程Ｓ４において、「Ｙ」）、ゴム材料モデル１６に基づいて、ゴム材料が製造される（工程Ｓ５）。一方、ゴム材料の性能が良好でないと判断された場合（工程Ｓ４において、「Ｎ」）、ゴム材料の諸条件を変更し（工程Ｓ６）、かつ、図３に示した手順に基づいて、ゴム材料モデル１６を再作成して（工程Ｓ７）、及び、工程Ｓ１〜工程Ｓ４が再度実施される。このように、本実施形態のシミュレーション方法、及び、製造方法では、ゴム材料の性能が良好になるまで、ゴム材料の諸条件が変更されるため、良好な強度（破壊特性）を有するゴム材料を、効率よく製造することができる。

以上、本発明の特に好ましい実施形態について詳述したが、本発明は図示の実施形態に限定されることなく、種々の態様に変形して実施しうる。

予め定められた物性値について差異を有する第１高分子成分と第２高分子成分とを少なくとも含むゴム材料の数値解析用のゴム材料モデルが作成された（実施例、比較例）。

実施例及び比較例では、第１高分子成分の分子鎖構造を複数の第１粗視化粒子を用いて表現した第１高分子鎖モデル、及び、第２高分子成分の分子鎖構造を複数の第２粗視化粒子を用いて表現した第２高分子鎖モデルが、コンピュータに入力された。第１高分子鎖モデルには、第１粗視化粒子間の相互作用を定義するための第１相互作用ポテンシャルが定義された。第２高分子鎖モデルには、第２粗視化粒子間の相互作用を定義するための第２相互作用ポテンシャルが定義された。さらに、実施例及び比較例では、架橋剤をモデル化した架橋剤粒子モデルが、コンピュータに入力された。

そして、実施例及び比較例では、図７に示した第１高分子成分の分布、第２高分子成分の分布、及び、架橋剤の分布に基づいて、第１高分子鎖モデル、第２高分子鎖モデル、及び、架橋剤粒子モデルがセル内に配置され、ゴム材料モデルが作成された。

実施例の第２相互作用ポテンシャルは、図１１に示した強度パラメータ定義工程の処理手順にしたがい、第１高分子成分と第２高分子成分との物性値の差異に基づいて、第１相互作用ポテンシャルとは異なるものとして定義された。一方、比較例の第２相互作用ポテンシャルは、物性値の差異を考慮せずに、第１相互作用ポテンシャルと同一のものが採用された。そして、実施例及び比較例のゴム材料モデルの引張変形が計算された。共通仕様は、次のとおりである。

第１、第２高分子鎖モデル：
第１、第２粗視化粒子の個数（鎖長）：２００個
各本数：４０本
第３、第４相互作用ポテンシャル：
粒子モデル間のばね定数ｋ：３０［ε／σ²］
伸びきり長Ｒ₀：１［σ］
架橋剤粒子モデル：２２９個
セル：
体積：１６２２９［σ³］
実施例の密度：０．９４［／σ³］
比較例の密度：０．８７［／σ³］
架橋密度：０．０１［σ³］
第５〜第１０相互作用ポテンシャル：
カットオフ距離ｒ_ij,_cut：２．５［σ］
強度パラメータε_ij：１．０［ε］
第１１相互作用ポテンシャル：
粒子モデル間のばね定数ｋ：３０［ε／σ²］
伸びきり長Ｒ₀：１［σ］
強度パラメータ定義工程：
高分子鎖モデル：
粗視化粒子の個数（鎖長）：２００個
分子数：４０本
第１５相互作用ポテンシャル：
粒子モデル間のばね定数ｋ：３０［ε／σ²］
伸びきり長Ｒ₀：１［σ］
強度パラメータε：１．０、１．２、１．４、１．６、
１．８、２．０、２．２、２．４［ε］
第１６相互作用ポテンシャル：
カットオフ距離ｒ_ij,_cut：２．５［σ］
強度パラメータε_ij：１．０［ε］
粒子の直径に相当σ_ij：１［σ］
第１相互作用ポテンシャル：
強度パラメータε_ij：１［ε］
カットオフ距離ｒ_ij,_cut：２．５［σ］
粒子の直径に相当σ_ij：１［σ］
第２相互作用ポテンシャル：
強度パラメータε_ij（実施例）：２［ε］
強度パラメータε_ij（比較例）：１［ε］
カットオフ距離ｒ_ij,_cut：２．５［σ］
粒子の直径に相当σ_ij：１［σ］
計算条件：
時間刻み（単位時間）：０．００６［τ］
合計時間：６０００［τ］
ひずみ速度：０．０００１［１／τ］
最大歪み：０．６
ポアソン比：０

図１７（Ａ）は、比較例のゴム材料モデル（歪：０．０）を示す図である。図１７（Ｂ）は、比較例のゴム材料モデル（歪：０．１２）を示す図である。図１７（Ｃ）は、比較例のゴム材料モデル（歪：０．３０）を示す図である。図１７（Ｄ）は、比較例のゴム材料モデル（歪：０．６０）を示す図である。比較例では、主として、第１高分子鎖モデル３が配置された領域と、第２高分子鎖モデル４が配置された領域との界面で破壊が生じる現象が計算され、実際のゴム材料の破壊現象に合致しなかった。

一方、実施例では、図１６（Ａ）〜（Ｄ）に示されるように、第２高分子鎖モデル４よりも小さい強度パラメータε_ijが定義された第１高分子鎖モデル３が配置された領域において、破壊が生じる（空孔２４が生じる）現象が計算されており、実際のゴム材料の破壊現象に合致した。

図１８は、実施例のゴム材料モデル、及び、比較例のゴム材料モデルの応力−歪み曲線を示すグラフである。実施例は、比較例に比べて、ヤング率が大きく計算された。これは、実際のゴム材料のヤング率の傾向に近似している。

このように、実施例のゴム材料モデルは、変形計算を行うシミュレーションに用いられることにより、第１高分子成分と第２高分子成分とを含むゴム材料の破壊特性を、精度良く評価できることを確認できた。

Ｓ１１ 第１高分子鎖モデルを入力する工程
Ｓ１２ 第２高分子鎖モデルを入力する工程
Ｓ２１ 第１相互作用ポテンシャルを定義する工程
Ｓ２２ 第２相互作用ポテンシャルを定義する工程

Claims (6)

1. 予め定められた物性値について差異を有する第１高分子成分と第２高分子成分とを少なくとも含むゴム材料の破壊現象を、コンピュータを用いてシミュレーションするための方法であって、
   前記第１高分子成分の分子鎖構造を複数の第１粗視化粒子を用いて表現した第１高分子鎖モデルを前記コンピュータに入力する工程と、
   前記第２高分子成分の分子鎖構造を複数の第２粗視化粒子を用いて表現した第２高分子鎖モデルを前記コンピュータに入力する工程と、
   前記第１高分子鎖モデルと前記第２高分子鎖モデルとの間に界面が形成されるように、前記第１高分子鎖モデル及び前記第２高分子鎖モデルを、前記ゴム材料の一部に対応する仮想空間であるセルに配置したゴム材料モデルを定義する工程と、
   前記第１粗視化粒子間の相互作用を定義するための第１相互作用ポテンシャルを定義する工程と、
   前記第２粗視化粒子間の相互作用を定義するための第２相互作用ポテンシャルを、前記物性値の差異に基づいて、前記第１相互作用ポテンシャルとは異なるものとして定義する工程と、
   前記コンピュータが、前記ゴム材料モデルの変形を計算して、前記第１相互作用ポテンシャル、及び、前記第２相互作用ポテンシャルに基づく前記界面での破壊現象を計算する工程とを含む、
   ゴム材料のシミュレーション方法。
2. 前記第１相互作用ポテンシャル及び前記第２相互作用ポテンシャルは、相互作用の強さの定義するための強度パラメータを含み、
   前記第２相互作用ポテンシャルを定義する工程は、前記第２相互作用ポテンシャルに、前記第１相互作用ポテンシャルとは異なる前記強度パラメータを定義する工程を含む、請求項１記載のゴム材料のシミュレーション方法。
3. 前記第１相互作用ポテンシャル及び前記第２相互作用ポテンシャルは、以下の式で定義される、請求項２記載のゴム材料のシミュレーション方法。
   ここで、
   ｒ_ij：粒子間の距離
   ｒ_ij,_cut：カットオフ距離
   ε_ij：強度パラメータ
   σ_ij：粒子の直径に相当
4. 前記強度パラメータを定義する工程は、予め定められた複数種類の前記強度パラメータε_ij毎に、前記物性値を計算する工程と、
   前記物性値の差異に基づいて、前記複数種類の前記強度パラメータε_ijから、前記第２相互作用ポテンシャルの前記強度パラメータε_ijを選択する工程とを含む、請求項３記載のゴム材料のシミュレーション方法。
5. 前記物性値は、弾性率を含む、請求項１ないし４のいずれかに記載のゴム材料のシミュレーション方法。
6. 請求項１ないし５のいずれかに記載の前記ゴム材料のシミュレーション方法において、前記性能が良好であると評価された前記ゴム材料モデルに基づいて、前記ゴム材料を製造する工程を含む、
   ゴム材料の製造方法。