特許 6579904 緩和スペクトルを算出する方法、物性値を算出する方法、装置、及びプログラム。

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

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

(11)【特許番号】6579904

(24)【登録日】2019年9月6日

(45)【発行日】2019年9月25日

(54)【発明の名称】緩和スペクトルを算出する方法、物性値を算出する方法、装置、及びプログラム。

(51)【国際特許分類】

   G06F 17/50 20060101AFI20190912BHJP

   G01N 33/44 20060101ALI20190912BHJP

   G06F 17/17 20060101ALI20190912BHJP

【ＦＩ】

   G06F17/50 612G

   G01N33/44

   G06F17/17

【請求項の数】9

【全頁数】12

(21)【出願番号】特願2015-203869(P2015-203869)

(22)【出願日】2015年10月15日

(65)【公開番号】特開2017-76280(P2017-76280A)

(43)【公開日】2017年4月20日

【審査請求日】2018年8月20日

(73)【特許権者】

【識別番号】000003148

【氏名又は名称】ＴＯＹＯ ＴＩＲＥ株式会社

(74)【代理人】

【識別番号】110000729

【氏名又は名称】特許業務法人 ユニアス国際特許事務所

(72)【発明者】

【氏名】狩野 康人

【審査官】 堀井 啓明

(56)【参考文献】

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

【文献】 特開２００９−２０４４０３（ＪＰ，Ａ）

【文献】 特開平０９−０６７４１３（ＪＰ，Ａ）

【文献】 特開平１０−３２６２７２（ＪＰ，Ａ）

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

Ｇ０６Ｆ１７／５０

Ｇ０１Ｎ３３／４４

Ｇ０６Ｆ１７／１７

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

【請求項1】

コンピュータが実行する方法であって、
予め設定された解析条件下において予め設定された分子モデルデータを用いた分子動力学計算を行い、各解析時点における応力データを算出するステップと、
前記応力データに基づき各解析時点毎の緩和弾性率を算出するステップと、
前記緩和弾性率を補間関数で近似した近似値を算出するステップと、
前記緩和弾性率の近似値を微分処理して緩和スペクトルを算出するステップと、
を含む緩和スペクトルを算出する方法。

【請求項2】

前記補間関数は、Ｂスプライン関数である、請求項１に記載の方法。

【請求項3】

前記緩和弾性率の近似値から緩和スペクトルを算出する演算は、式（２）に表されるＳｃｈｗａｒｚｌ−Ｓｔａｖｅｒｍａｎによる近似式を用いており、次数ｋは２次以上である、請求項１又は２に記載の方法。

【数1】

Ｇは緩和弾性率、τは緩和時間である。

【請求項4】

請求項１〜３のいずれかに記載の方法であって、
前記緩和スペクトルに基づき、貯蔵弾性率、損失弾性率及び損失正接のうち少なくとも１つの物性値を算出するステップを更に含む、物性値を算出する方法。

【請求項5】

予め設定された解析条件下において予め設定された分子モデルデータを用いた分子動力学計算の実行結果に基づき、各解析時点における応力データを算出する応力算出部と、
前記応力データに基づき各解析時点毎の緩和弾性率を算出する緩和弾性率算出部と、
前記緩和弾性率を補間関数で近似した近似値を算出する近似部と、
前記緩和弾性率の近似値を微分処理して緩和スペクトルを算出する緩和スペクトル算出部と、
を備える緩和スペクトルを算出する装置。

【請求項6】

前記補間関数は、Ｂスプラインである、請求項５に記載の装置。

【請求項7】

前記緩和弾性率の近似値から緩和スペクトルを算出する演算は、式（２）に表されるＳｃｈｗａｒｚｌ−Ｓｔａｖｅｒｍａｎによる近似式を用いており、次数ｋは２次以上である、請求項５又は６に記載の装置。

【数2】

Ｇは緩和弾性率、τは緩和時間である。

【請求項8】

請求項５〜７に記載の装置であって、
前記緩和スペクトルに基づき、貯蔵弾性率、損失弾性率及び損失正接のうち少なくとも１つの物性値を算出する物性値算出部を更に備える、物性値を算出する装置。

【請求項9】

請求項１〜４のいずれかに記載の方法をコンピュータに実行させるプログラム。

【発明の詳細な説明】

【技術分野】

【0001】

本開示は、分子動力学計算により緩和スペクトルを算出する方法、物性値（貯蔵弾性率、損失弾性率及び損失正接など）を算出する方法、装置及びプログラムに関する。

【背景技術】

【0002】

貯蔵弾性率、損失弾性率及び損失正接等の物性値の計算は、次のように行われる。
（１）分子動力学計算を実行し、応力データを得る。
（２）グリーン・久保公式に基づき応力データから緩和弾性率を算出する。
（３）Ｓｃｈｗａｒｚｌ−Ｓｔａｖｅｒｍａｎによる近似式を用い、緩和弾性率に基づき緩和スペクトルを算出する。
（４）緩和スペクトルをＳｔｉｅｌｔｊｅｓ変換し、貯蔵弾性率、損失弾性率及び損失正接を算出する。

【0003】

例えば、特許文献１には、高分子モデルを用いて分子動力学計算を行い、その結果からグリーン・久保公式に基づき緩和弾性率を算出することが開示されている。緩和弾性率Ｇ（ｔ）は、式（１）を用いて算出できる。

【数1】

σ（ｔ）は時刻ｔのせん断応力テンソル、Ｖは体積、Ｔは温度、Ｋ_Ｂはボルツマン定数である

【先行技術文献】

【特許文献】

【0004】

【特許文献1】特開２０１４−２０３２６２号公報

【発明の概要】

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

【0005】

緩和スペクトルＨ（τ）を算出するためのＳｃｈｗａｒｚｌ−Ｓｔａｖｅｒｍａｎによる近似式は、式（２）に表される。式（２）の次数ｋを無限大にする計算は現実的に不可能であるので、式（３）に示す１次の式が利用しやすい。すなわち、緩和スペクトルの算出には、式（３）に示すように微分演算が含まれる。τは緩和時間である。

【数2】

【数3】

【0006】

緩和スペクトルＨ（τ）は、貯蔵弾性率、損失弾性率及び損失正接等の物性値を算出する基礎となるデータであるので、正しく計算されていることが望まれる。しかし、緩和スペクトルＨ（τ）がどの程度正しいかを見積もることは難しい。そこで、式（４）を用いて緩和スペクトルＨ（τ）をラプラス変換すれば、緩和弾性率Ｇ（ｔ）に戻すことができることを利用する。式（４）を用いて復元した緩和弾性率Ｇ（ｔ）_復元と元の緩和弾性率Ｇ（ｔ）が同じになれば、緩和スペクトルＨ（τ）が正しく算出されていることになる。図５に示すように、Ｇ（ｔ）とＧ（ｔ）_復元との残差の最小二乗和を算出して評価することが考えられる。

【数4】

【0007】

しかしながら、図５に示すように、短時間領域と長時間領域において両者が大きくことなっていることが分かる。すなわち、緩和スペクトルの算出精度が十分であるとはいえない。

【0008】

緩和スペクトルＨ（τ）を精度良く求めるためのアプローチの一つとして、Ｓｃｈｗａｒｚｌ−Ｓｔａｖｅｒｍａｎによる近似式（２）の次数を上げることが考えられる。しかしながら、図６Ａに示す１次から図６Ｂに示す２次になれば改良されるが、図６Ｃに示す３次になれば、逆に悪化してしまう場合があることが判明した。図６Ａに示す例では、式（２）の１次による算出結果の残差二乗和が0.0455である。図６Ｂに示す例では、式（２）の２次による算出結果の残差二乗和が0.0247である。図６Ｃに示す例では、式（２）の３次による算出結果の残差二乗和が0.829であった。３次で悪化している。

【0009】

上記した緩和スペクトルの算出精度が十分でないという認識、精度を向上させる課題、その解決手段に言及している文献は見つかっていない。

【0010】

本開示は、上記新たな知見に基づく課題に着目してなされたものであって、その目的は、算出精度を向上させた緩和スペクトル及びそれを用いた物性値を算出する方法、装置及びプログラムを提供することである。

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

【0011】

本開示は、上記目的を達成するために、次のような手段を講じている。

【0012】

すなわち、本開示の緩和スペクトルを算出する方法は、
予め設定された解析条件下において予め設定された分子モデルデータを用いた分子動力学計算を行い、各解析時点における応力データを算出するステップと、
前記応力データに基づき各解析時点毎の緩和弾性率を算出するステップと、
前記緩和弾性率を補間関数で近似した近似値を算出するステップと、
前記緩和弾性率の近似値を微分処理して緩和スペクトルを算出するステップと、
を含む。

【0013】

本開示の緩和スペクトルを算出する装置は、
予め設定された解析条件下において予め設定された分子モデルデータを用いた分子動力学計算の実行結果に基づき、各解析時点における応力データを算出する応力算出部と、
前記応力データに基づき各解析時点毎の緩和弾性率を算出する緩和弾性率算出部と、
前記緩和弾性率を補間関数で近似した近似値を算出する近似部と、
前記緩和弾性率の近似値を微分処理して緩和スペクトルを算出する緩和スペクトル算出部と、
を備える。

【0014】

このように、離散データである緩和弾性率を補間関数で近似するので、緩和弾性率のバラツキが平滑化され、微分処理による精度の悪化を抑制し、精度を向上させることが可能となる。

【図面の簡単な説明】

【0015】

【図1】本開示の緩和スペクトル及び物性値を算出する装置を示すブロック図。

【図2】本開示の装置で実行される算出処理ルーチンを示すフローチャート。

【図3】緩和弾性率及びそれを補間関数で近似した近似値を示す図。

【図4A】Ｂスプライン関数で近似した緩和弾性率Ｇ（ｔ）を１次の式（３）で緩和スペクトルＨ（τ）に変換し、式（４）のラプラス変換で緩和弾性率Ｇ（ｔ）に戻したものを示す図。

【図4B】Ｂスプライン関数で近似した緩和弾性率Ｇ（ｔ）を２次の式（５）で緩和スペクトルＨ（τ）に変換し、式（４）のラプラス変換で緩和弾性率Ｇ（ｔ）に戻したものを示す図。

【図4C】Ｂスプライン関数で近似した緩和弾性率Ｇ（ｔ）を３次の式（６）で緩和スペクトルＨ（τ）に変換し、式（４）のラプラス変換で緩和弾性率Ｇ（ｔ）に戻したものを示す図。

【図5】緩和弾性率Ｇ（ｔ）をそのまま１次の式（３）で緩和スペクトルＨ（τ）に変換し、式（４）のラプラス変換で緩和弾性率Ｇ（ｔ）に戻したものを示す図。

【図6A】緩和弾性率Ｇ（ｔ）をそのまま１次の式（３）で緩和スペクトルＨ（τ）に変換し、式（４）のラプラス変換で緩和弾性率Ｇ（ｔ）に戻したものを示す図。

【図6B】緩和弾性率Ｇ（ｔ）をそのまま２次の式（５）で緩和スペクトルＨ（τ）に変換し、式（４）のラプラス変換で緩和弾性率Ｇ（ｔ）に戻したものを示す図。

【図6C】緩和弾性率Ｇ（ｔ）をそのまま３次の式（６）で緩和スペクトルＨ（τ）に変換し、式（４）のラプラス変換で緩和弾性率Ｇ（ｔ）に戻したものを示す図。

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

【0016】

以下、本開示の一実施形態を、図面を参照して説明する。

【0017】

［緩和スペクトル及び物性値を算出する装置］
本実施形態の装置は、予め設定された分子モデルデータを用いた分子動力学計算を行い、緩和スペクトルを算出する装置である。また、装置は、算出した緩和スペクトルを用いて物性値（貯蔵弾性率、損失弾性率及び損失正接）を算出する。勿論、緩和スペクトルのみを算出するように構成してもよいし、貯蔵弾性率、損失弾性率及び損失正接のうち少なくとも１つを算出するように構成してもよい。

【0018】

図１に示すように、装置１は、設定部１０と、分子動力学計算実行部１１と、応力算出部１２と、緩和弾性率算出部１３と、近似部１４と、緩和スペクトル算出部１５と、物性値算出部１６と、を有する。これら各部１０〜１６は、ＣＰＵ、メモリ、各種インターフェイス等を備えたパソコン等の情報処理装置において予め記憶されている図３に示す物性値算出処理ルーチンをＣＰＵが実行することによりソフトウェア及びハードウェアが協働して実現される。

【0019】

図１に示す設定部１０は、キーボードやマウス等の既知の操作部を介してユーザからの操作を受け付け、解析対象となる分子モデル（ゴムモデルなど）に関する情報の設定、温度・体積などの分子動力学計算に用いる各種解析条件の設定を実行し、これらをメモリに記憶する。

【0020】

図１に示す分子動力学計算実行部１１は、予め設定された分子モデルデータを用いて予め設定された解析条件下において分子動力学計算を行い、各粒子の座標及び速度の時系列データを算出する。一例として、ランジュバン動力学のアルゴリズムによる定温における分子動力学計算などが挙げられる。

【0021】

図１に示す応力算出部１２は、分子動力学計算実行部１１の算出結果（時系列データ）に基づき、各解析時間における応力データを算出する。具体的には、各解析時点における応力テンソルを算出する。

【0022】

図１に示す緩和弾性率算出部１３は、応力データに基づき各解析時点毎の緩和弾性率を算出する。具体的には、時系列データから、せん断応力テンソルの自己相関関数を計算し、式（１）に示すグリーン・久保公式に代入して緩和弾性率Ｇ（ｔ）を算出する。図３は、鎖長１００の高分子１００本で構成される高分子メルトモデルを用いて算出した緩和弾性率Ｇ（ｔ）をバツ印で示す。同図は縦軸及び横軸共にログスケールに設定している。同図に示すように、解析時点毎に緩和弾性率の値が存在する離散的な時系列データである。

【数5】

σ（ｔ）は時刻ｔのせん断応力テンソル、Ｖは体積、Ｔは温度、Ｋ_Ｂはボルツマン定数である。

【0023】

図１に示す近似部１４は、緩和弾性率Ｇ（ｔ）を補間関数で近似した近似値を算出する。補間関数としては、Ｂスプライン、ベジェ、非一様有理Ｂスプライン（Non-Uniform Rational Basis Spline）、スプライン、ラグランジュ補間が挙げられる。本実施形態では、３次のＢスプラインを用いているが、これに限定されず、２次のＢスプラインも使用できる。この中では最もＢスプラインが取り扱いしやすい。ベジェは多数の制御点を持つ場合には接続条件が煩雑になるからである。スプラインは、滑らかな曲線が得られるものの、弾性緩和率の値をそのまま制御点にすれば計算結果の振動（バラツキ）をダイレクトに拾ってしまうからである。ラグランジュ補間は制御点が多くなると振動性が出てしまうからである。

【0024】

図３は、緩和弾性率Ｇ（ｔ）をバツ印で示すと共に、緩和弾性率Ｇ（ｔ）を３次のＢスプライン関数で近似し、得られる近似値（Ｂスプライン関数の値）を実線で描いた図である。バツ印で示す離散的データである緩和弾性率Ｇ（ｔ）が平滑化され、バラツキが低減されていることが分かる。

【0025】

図１に示す緩和スペクトル算出部１５は、近似値を微分処理して緩和スペクトルＨ（τ）を算出する。具体的には、Ｓｃｈｗａｒｚｌ−Ｓｔａｖｅｒｍａｎによる近似式（２）を用い、緩和弾性率Ｇ（ｔ）の近似値（Ｂスプライン関数の値）を微分処理して緩和スペクトルＨ（τ）を算出する。式（３）は１次、式（５）は２次、式（６）は３次である。

【数6】

【数7】

【数8】

【0026】

図４Ａ〜Ｃは、Ｂスプライン関数で近似した緩和弾性率Ｇ（ｔ）を、１次の式（３）、２次の式（５）、３次の式（６）で緩和スペクトルＨ（τ）に変換し、それらを式（４）のラプラス変換で緩和弾性率Ｇ（ｔ）_復元に戻したものを示す。同４Ａに示す１次の式（２）による算出結果の残差二乗和が0.0324である。図４Ｂに示す２次の式（２）による算出結果の残差二乗和が0.0173である。図４Ｃに示す３次の式（２）による算出結果の残差二乗和が0.0133であった。次数が上がるほど残差が小さくなり、精度が向上している。図６Ａ〜Ｃに比べても残差が小さくなり、精度が向上している。精度が向上した原因について、発明者は、次のように考えている。従来では、緩和弾性率データにはバラツキが含まれており、緩和弾性率Ｇ（ｔ）から緩和スペクトルＨ（τ）を算出する際に行う微分処理によって、算出精度が悪化していた。本開示の方法では、緩和弾性率Ｇ（ｔ）を補間関数で近似しているので、緩和弾性率Ｇ（ｔ）の振動成分（バラツキ）が平滑化され、悪化の原因が低減されたためと考えている。

【0027】

したがって、式（２）に表されるＳｃｈｗａｒｚｌ−Ｓｔａｖｅｒｍａｎによる近似式を用いており、次数は１以上であればよい。より一層精度を向上させるのであれば２次以上がよい。１次から２次への向上代が最も大きく、２次から３次への向上代は小さくなる。よって、演算コストに対する精度の向上代を考慮すれば、２次が好ましい。

【0028】

図１に示す物性値算出部１６は、緩和スペクトルＨ（τ）に基づき、貯蔵弾性率Ｇ’（ω）、損失弾性率Ｇ’’（ω）及び損失正接ｔａｎδのうち少なくとも１つの物性値を算出する。具体的には、貯蔵弾性率Ｇ’（ω）を式（７）で、損失弾性率Ｇ’’（ω）を式（８）で算出する。損失正接ｔａｎδは式（９）で算出する。ωは周波数であり、δは歪と応力の位相差を示す。なお、物性値算出部１６は実装に応じて省略可能である。

【数9】

【0029】

［緩和スペクトル及び物性値を算出する装置］
上記装置１を用い、緩和スペクトルを算出する方法を説明する。また、緩和スペクトルを用いて物性値を算出する方法を、図２を参照して説明する。

【0030】

まず、ステップＳＴ１において、初期設定部１０は、解析対象となる分子モデルデータ、温度・体積などの分子動力学計算に用いる各種解析条件の設定を実行し、これらをメモリに記憶する。

【0031】

ステップＳＴ２において、図１に示す分子動力学計算実行部１１は、メモリに記憶されている分子モデルデータを用いて予め設定された解析条件下において分子動力学計算を行い、各粒子の座標及び速度の時系列データを算出する。

【0032】

次のステップＳＴ３において、応力算出部１２は、時系列データに基づき各解析時間における応力データを算出する。次のステップＳＴ４において、緩和弾性率算出部１３は、応力データに基づき緩和弾性率Ｇ（ｔ）を算出する。

【0033】

次のステップＳＴ５において、近似部１４は、緩和弾性率Ｇ（ｔ）を補間関数で近似した近似値を算出する。

【0034】

次のステップＳＴ６において、緩和スペクトル算出部１５は、緩和弾性率Ｇ（ｔ）の近似値を微分処理して緩和スペクトルＨ（τ）を算出する。

【0035】

次のステップＳＴ７において、物性値算出部１６は、緩和スペクトルＨ（τ）に基づき、貯蔵弾性率Ｇ’（ω）、損失弾性率Ｇ’’（ω）及び損失正接ｔａｎδのうち少なくとも１つの物性値を算出する。

【0036】

以上のように、本実施形態の緩和スペクトルを算出する方法は、
予め設定された解析条件下において予め設定された分子モデルデータを用いた分子動力学計算を行い、各解析時点における応力データを算出するステップ（ＳＴ２〜ＳＴ３）と、
応力データに基づき各解析時点毎の緩和弾性率Ｇ（ｔ）を算出するステップ（ＳＴ４）と、
緩和弾性率Ｇ（ｔ）を補間関数で近似した近似値を算出するステップ（ＳＴ５）と、
緩和弾性率Ｇ（ｔ）の近似値を微分処理して緩和スペクトルＨ（τ）を算出するステップ（ＳＴ６）と、
を含む。

【0037】

本実施形態の緩和スペクトルを算出する装置１は、
予め設定された解析条件下において予め設定された分子モデルデータを用いた分子動力学計算の実行結果に基づき、各解析時点における応力データを算出する応力算出部１２と、
応力データに基づき各解析時点毎の緩和弾性率Ｇ（ｔ）を算出する緩和弾性率算出部１３と、
緩和弾性率Ｇ（ｔ）を補間関数で近似した近似値を算出する近似部１４と、
緩和弾性率Ｇ（ｔ）の近似値を微分処理して緩和スペクトルＨ（τ）を算出する緩和スペクトル算出部１５と、
を備える。

【0038】

このように、離散データである緩和弾性率Ｇ（ｔ）を補間関数で近似するので、緩和弾性率Ｇ（ｔ）のバラツキが平滑化され、微分処理による精度の悪化を抑制し、精度を向上させることが可能となる。

【0039】

本実施形態において、補間関数はＢスプライン関数である。
この構成によれば、緩和弾性率Ｇ（ｔ）の各値をそのまま制御点に使用できるので、実装が容易である。また、滑らかな曲線が得られ、且つ、他の補間関数に比べてバラツキを抑制できるので、他の補間関数に比べて、精度を向上させることが可能となる。

【0040】

本実施形態において、緩和弾性率Ｇ（ｔ）の近似値から緩和スペクトルＨ（τ）を算出する演算は、式（２）に表されるＳｃｈｗａｒｚｌ−Ｓｔａｖｅｒｍａｎによる近似式を用いており、次数ｋは２次以上が好ましい。

【数10】

Ｇは緩和弾性率である。
この構成によれば、１次よりも精度を向上させることが可能となる。

【0041】

本実施形態の方法において、緩和スペクトルＨ（τ）に基づき、貯蔵弾性率Ｇ’（ω）、損失弾性率Ｇ’’（ω）及び損失正接ｔａｎδのうち少なくとも１つの物性値を算出するステップ（ＳＴ７）を更に含む。
本実施形態の装置において、緩和スペクトルＨ（τ）に基づき、貯蔵弾性率Ｇ’（ω）、損失弾性率Ｇ’’（ω）及び損失正接ｔａｎδのうち少なくとも１つの物性値を算出する物性値算出部１６を更に備える。
この構成によれば、従来に比べて精度のよい緩和スペクトルＨ（τ）を用いるので、物性値の算出精度を向上させることが可能となる。

【0042】

本実施形態に係るコンピュータプログラムは、上記方法を構成する各ステップをコンピュータに実行させるプログラムである。このプログラムを実行することによっても、上記方法の奏する作用効果を得ることが可能となる。言い換えると、上記方法を使用しているとも言える。

【0043】

以上、本開示の実施形態について図面に基づいて説明したが、具体的な構成は、これらの実施形態に限定されるものでないと考えられるべきである。本開示の範囲は、上記した実施形態の説明だけではなく特許請求の範囲によって示され、さらに特許請求の範囲と均等の意味および範囲内でのすべての変更が含まれる。

【0044】

例えば、図１に示す各部１０〜１６は、所定プログラムをコンピュータのＣＰＵで実行することで実現しているが、各部を専用メモリや専用回路で構成してもよい。

【0045】

上記の各実施形態で採用している構造を他の任意の実施形態に採用することは可能である。各部の具体的な構成は、上述した実施形態のみに限定されるものではなく、本開示の趣旨を逸脱しない範囲で種々変形が可能である。

【符号の説明】

【0046】

１２…応力算出部
１３…緩和弾性率算出部
１４…近似部
１５…緩和スペクトル算出部
１６…物性値算出部

【図1】

【図2】

【図3】

【図4A】

【図4B】

【図4C】

【図5】

【図6A】

【図6B】

【図6C】