JP6443304B2

JP6443304B2 - 粘弾性材料特性解析装置

Info

Publication number: JP6443304B2
Application number: JP2015220140A
Authority: JP
Inventors: 祐介小西
Original assignee: Toyota Motor Corp
Current assignee: Toyota Motor Corp
Priority date: 2015-11-10
Filing date: 2015-11-10
Publication date: 2018-12-26
Anticipated expiration: 2035-11-10
Also published as: JP2017091193A

Description

本発明は、粘弾性材料の特性を解析し、その特性をモデリングする装置に関する。

ゴム材料などの粘弾性材料の特性を有限要素法により解析するためのモデルとして、弾性要素と粘弾性要素のつながり（以下、「ライン」とよぶ）を並列配置した粘弾性材料構成則が用いられる。特性モデルを定義するためには、ラインごとに、粘弾性要素の減衰特性を示す緩和時間を求める必要がある（特許文献１、非特許文献１参照）。一般に、ゴム材料に対する調和振動試験では加振する振幅の大きさによって応力とひずみの相関関係が変化するため、振幅依存性を再現できる特性モデルを用いることが望ましい。特許文献１では、緩和時間を歪速度のべき関数として定義することにより、振幅依存性に対応して予測精度を高めている。

特開２０１５−０７５３８３号公報

J.C. Simo and T.J.R Hughes著「Computational Inelasticity」 Interdisciplinary Applied Mathematics volume 7, Springer Verlag(1998)

緩和時間を定義したあと、各ラインの剛性割合を求める。この剛性割合を求める処理を「最適化処理」とよぶ。本発明者は、粘弾性材料の特性を表現するための１以上のライン（以下、「基本ライン」とよぶ）とは別に、モデルの剛性に影響しない計算上の便宜のためのライン（以下、「調整ライン」とよぶ）を付加することで最適化処理を高速化できることに想到した。特に、調整ラインに含まれる粘弾性要素は、特性モデルの挙動に影響しないものであることが理想である。しかし、歪速度が高い場合には、調整ラインに含まれる粘弾性要素の影響を無視できなくなる可能性があることがわかった。

本発明は、本発明者の上記課題認識に基づいてなされたものであり、その主たる目的は、粘弾性材料の特性モデルの再現精度を維持しつつ、計算処理を高速化するための技術、を提供することを目的とする。

上記課題を解決するために、本発明のある態様の粘弾性材料特性解析装置は、粘弾性材料に変位を与えたときの振幅および弾性率のデータセットを取得し、前記振幅の対数値に対する前記弾性率の対数値の変化率からべき数を算出する計算部と、前記粘弾性材料の応力−ひずみ特性をシミュレートする特性モデルを前記べき数に基づいて定義するモデル生成部と、を備える。
前記特性モデルにおいては、弾性要素と粘弾性要素が直列接続される１以上の基本ラインと、弾性要素と粘弾性要素が直列接続される調整ラインとが並列接続されており、前記モデル生成部は、前記基本ラインに含まれる粘弾性要素の緩和時間を歪速度の前記べき数によるべき乗値に比例する値として設定し、前記調整ラインの粘弾性要素の緩和時間を歪速度の逆数値に比例する値として設定する。

本発明によれば、粘弾性材料の特性モデルの再現精度を維持しつつ、計算処理を高速化させやすくなる。

粘弾性材料構成則を模式的に示す図である。調和振動試験における試験条件を示す図である。べき数の算出方法を模式的に示すグラフである。粘弾性材料特性解析装置の機能ブロック図である。粘弾性材料の特性解析処理過程を示すフローチャートである。本実施形態における特性モデル（粘弾性材料構成則）の模式図である。調整ラインのべき数を共通べき数としたときの歪速度と緩和時間の関係を示すグラフである。調整ラインのべき数を共通べき数としたときの歪速度と緩和時間の関係を示すグラフである。調整ラインのべき数を−１としたときの歪速度と緩和時間の関係を示すグラフである。

調整ラインに関して説明する前に、まず、調整ラインを持たない基本的な特性モデルを中心として説明する。

図１は、粘弾性材料構成則（特性モデル）の基本構成を模式的に示す図である。
本図に示す特性モデルでは、弾性率Ｇ０の弾性要素と、弾性率Ｇｉ（ｉ＝１〜Ｎ、Ｎは自然数）の弾性要素と粘性係数ηｉの粘性要素とが直列接続された複数の粘弾性要素とが並列接続される。本実施形態においては、弾性率Ｇ０の弾性要素のみを含むラインも基本ラインとよぶ。すなわち、図１では、Ｎ本の基本ラインが並列接続されている。複数の基本ラインを組み合わせた特性モデルは、タイヤやゴムブッシュなどのゴム部品の動特性を表すためのモデルとして用いられる。

ここで、並列接続される基本ラインの剛性割合をγｉ（ｉ＝０〜Ｎ）、粘弾性要素における緩和時間をτｉ＝ηｉ／Ｇｉとすると、この粘弾性材料モデルにおける時刻ｔでの応力Ｓは、下記の式（１）、（２）により表される。

式（１）において、Ｓは、第２Ｐｉｏｌａ−Ｋｉｒｃｈｈｏｆｆ応力を示し、上付き記号を有するＳ○は、粘性力成分を除去した弾性成分のみの応力であることを示す。Ｊは、粘弾性材料の体積変化率を示す。体積変化率Ｊは、ある物質点における変形前と変形後の位置の線形変換の関係を示す変形勾配テンソルＦのデターミナント（ｄｅｔ）を用いて、Ｊ＝ｄｅｔ［Ｆ］により表される。演算子ＤＥＶは、右Ｃａｕｃｈｙ−ＧｒｅｅｎテンソルＣ＝ＦＴ・Ｆを用いて、下記式（３）で表される。なお、式（３）における［・］は、演算子ＤＥＶの演算対象となる変数を表す。

また、Ｑｉは、それぞれの粘弾性要素における粘性力を示し、式（１）におけるＱｉは式（２）に示す発展方程式により表される。式（２）において、

は、超弾性体におけるひずみポテンシャルエネルギーの偏差成分を表す。また、

は、体積成分を除去した修正右Ｃａｕｃｈｙ−Ｇｒｅｅｎテンソルであり、式（４）で表される。

上記の式（１）、（２）により表される第２Ｐｉｏｌａ−Ｋｉｒｃｈｈｏｆｆ応力Ｓは、積分因子ｅｘｐ（ｔ／τｉ）を用いた畳み込み積分形式を積分することで、下記の式（５）に示す積分形で表すことができる。

なお、式（５）における

は、ひずみポテンシャルエネルギーの体積成分である。また、ｇ（ｔ）は、緩和関数であり、上記の式（６）により表される。

ここで、時間ｔの関数として表される式（５）を下記の式（７）に変形することで、時刻ｔｎ＋１における第２Ｐｉｏｌａ−Ｋｉｒｃｈｈｏｆｆ応力Ｓｎ＋１を得ることができる。なお、計算ステップｎにおける関数を（・）ｎ、計算ステップｎ＋１における関数を（・）ｎ＋１と表記している。

式（７）における

は、中点定理を用いて［ｔｎ，ｔｎ＋１］の時間間隔で積分した近似解として得られる中間関数であり、下記式（８）で表される。

ここで、

は、下記式（９）、（１０）で定義される。

なお、式（９）は、Ｋｉｒｃｈｈｏｆｆ弾性応力

を、下記式（１１）

により定義することにより、下記式（１２）により表すこともできる。

また、Ｋｉｒｃｈｈｏｆｆ応力テンソル

は、第２Ｐｉｏｌａ−Ｋｉｒｃｈｈｏｆｆ応力Ｓｎ＋１を用いて下記の式（１３）により表すことができる。

したがって、式（７）、（１３）より、Ｋｉｒｃｈｈｏｆｆ応力テンソルを下記の式（１４）により表すことができる。

なお、式（１４）における演算子ｄｅｖは、下記の式（１５）により定義される。なお、式（１５）における（・）は、演算子ＤＥＶの演算対象となる変数を表す。

また、式（１４）における緩和関数ｇ^＊は、下記の式（１６）により定義される。

以上より、中間関数として

を計算ステップ毎に保持することにより、Ｋｉｒｃｈｈｏｆｆ応力テンソルの値を数値解析により得ることができる。なお、上述の式変形の詳細については、例えば、非特許文献１に記載される。

本実施形態では、図１に示す特性モデルにおいて、粘弾性要素の減衰特性を表す緩和時間τｉを歪速度に依存した変数とすることにより、振幅依存性の再現性を高めている。振幅依存性の精度を高める手法として、緩和時間τｉを歪速度のべき関数とすることが有効である（特許文献１）。本実施形態では、下記式（１７）に示す関係式を満たす緩和時間τｉを用いる。

式（１７）において、

は、偏差成分を除去したＧｒｅｅｎ−Ｌａｇｒａｎｇｅひずみテンソルであり、修正右Ｃａｕｃｈｙ−Ｇｒｅｅｎテンソル

を用いて下記の式（１８）で表される。

また、

は、歪速度であり、歪テンソルの時間微分値を表す。また、

は、歪速度の大きさを表し、３次元の歪テンソルを用いる場合には、下記式（１９）により表される。

なお、べき数ｍｉおよび比例定数１／Ａｉは、本実施の形態で導入する値であり、調和振動試験の試験結果を用いて後述する関係式により導出される。

つづいて、式（１７）に示す緩和時間τｉの導入に際し、数値解析上必要となる式の変形について述べる。上述の式（１４）を用いてＫｉｒｃｈｈｏｆｆ応力を算出するためには、緩和時間τｉを含む式を変形する必要がある。具体的には、式（８）を下記の式（２０）に、式（９）を下記の式（２１）に変形すればよい。

ここで、緩和時間τｉの関数は、下記の式（２２）、（２３）により定義される。

なお、式（２０）の右辺第２項および式（２１）における緩和時間τｉとして式（２３）を用いるのは、中点定理を用いて［ｔｎ，ｔｎ＋１］の時間間隔で積分した近似解を用いたためである。

次に、本実施の形態に示す特性モデルにおける材料定数の同定方法について述べる。

図２は、調和振動試験における試験条件を示す図である。
本図は、横軸を時間ｔとして、予歪Ｅｐｒｅ、振幅ε、周波数ωのとした場合における粘弾性材料の試験体に加振される歪みＥ＝Ｅｐｒｅ＋εｓｉｎ（ωｔ）のグラフを示す。このような歪みＥを加振した場合における試験体の動的弾性率Ｇを測定することにより、下記式（２４）の関係を用いてべき数ｍｉを導出できる。

なお、式（２４）は、図２に示す調和振動試験の入力に対して、式（２）に式（１７）を代入した発展方程式である下記式（２５）を解くことにより得られる。

図３は、べき数ｍｉの算出方法を模式的に示すグラフである。
本図は、調和振動試験において加振する振幅εの対数値ｌｏｇ（ε）を横軸として、予歪Ｅｐｒｅおよび周波数ωを変化させたときの試験結果値である動的弾性率Ｇの対数値ｌｏｇ（Ｇ）をプロットしたものである。グラフ上の近似線Ｌ１は、予歪Ｅｐｒｅ１と周波数ω１に対して、振幅の値ε１〜ε３を変化させたときのプロット値に対応する傾きｌ１の近似線である。同様に、近似線Ｌ２、Ｌ３は、予歪Ｅｐｒｅ２、Ｅｐｒｅ３と周波数ω２、ω３に対して、振幅の値ε１〜ε３を変化させたときのプロット値に対応する傾きｌ２、ｌ３の近似線である。式（２４）の両辺に対数をとることにより、下記式（２６）を得ることができることから、この傾きｌｉからべき数ｍｉを得ることができる。

式（２６）より、べき数はｍｉ＝−１−ｌｉの関係式により得ることができる。なお、式（２６）におけるαは、図３に示す近似線の切片を示す。これにより、調和振動試験における加振周波数ωｉに対応するべき数ｍｉを得ることができる。なお、図３に示す振幅εの値は例示であり、べき数ｍｉを得るために４以上の振幅εの値を用いてもよいし、振幅εの値を２としてよい。

なお、べき数ｍｉの値が加振周波数ωｉによってそれほど変化しない場合には、それぞれの周波数に対応するｍｉを得る代わりに、それぞれのｍｉの平均値として得られる共通のべき数ｍを用いることとしてもよい。共通のべき数ｍを用いることにより、周波数ごとに異なるべき数ｍｉを用いる場合と比べて、計算負荷を減らすことができる。

また、比例定数１／Ａｉは、上記方法で得られたべき数ｍｉを用いて、下記式（２７）の関係から得ることができる。上述のべき数ｍｉにおける導出方法と同様、式（２７）における周波数ωｉは、調和振動試験における加振周波数であり、この周波数を変化させることにより周波数成分ωｉが異なるそれぞれの粘弾性要素に対応する比例定数１／Ａｉを得ることができる。

べき数ｍｉおよび比例定数１／Ａｉを同定したあと、剛性割合γｉを決めるための最適化処理を実行する。

図４は、粘弾性材料特性解析装置１００の機能ブロック図である。
粘弾性材料特性解析装置１００の各構成要素は、任意のコンピュータのＣＰＵ、メモリ、メモリにロードされた本図の構成要素を実現するプログラム、そのプログラムを格納するハードディスクなどの記憶ユニット、ネットワーク接続用インタフェースを中心にハードウェアとソフトウェアの任意の組み合わせによって実現される。そして、その実現方法、装置にはいろいろな変形例があることは、当業者には理解されるところである。以下説明する各図は、ハードウェア単位の構成ではなく、機能単位のブロックを示している。

粘弾性材料特性解析装置１００は、計算部１０２、モデル格納部１０４およびモデル生成部１０６を含む。計算部１０２は、粘弾性材料の振動実験結果を取得し、べき数ｍｉ等の各種材料定数を同定する。モデル生成部１０６は、各材料定数に基づいて、粘弾性材料の特性モデルを生成し、モデル格納部１０４に格納する。

図５は、粘弾性材料の特性解析処理過程を示すフローチャートである。
まず、粘弾性材料の振動実験を行い、計算部１０２は予歪Ｅｐｒｅ、周波数ω、動的弾性率Ｇ等の各種データを取得する（Ｓ１２）。計算部１０２は、これらのデータに基づいてべき数ｍｉや材料定数Ａｉを上述の計算式に基づいて算出する（Ｓ１４）。モデル生成部１０６は、剛性割合γｉを求めるために最適化処理を実行する（Ｓ１６）。最適化処理は既知のアルゴリズムの応用でよい。本実施形態における最適化処理では、ヒューリスティック法の一種である解適応焼きなまし法により剛性割合を求める。通常、この最適化処理は解が収束するまでに時間がかかる。

図６は、本実施形態における特性モデル１１０（粘弾性材料構成則）の模式図である。
特性モデル１１０においては、弾性率Ｇ０の弾性要素のみを含む基本ライン１１２−０と、弾性率Ｇｉの弾性要素および粘性係数ηｉ（緩和時間τｉ）の粘性要素の双方を含む基本ライン１１２ｉ（１１２−１〜１１２−Ｎ）に加えて、弾性率Ｇｘの弾性要素および粘性係数ηｘ（緩和時間τｘ）の粘性要素を含む調整ライン１１４も並列接続される。調整ライン１１４を追加することにより、全基本ライン１１２の剛性割合と調整ライン１１４の剛性割合の合計値を１とすることでソルバに制約をかけることができるため、最適化処理を高速化できる。本発明者の実験では調整ライン１１４を導入することにより最適化処理に要する時間を約６割カットできる。

調整ライン１１４は最適化処理の高速化のために便宜的に導入されるものであるから、特性モデル１１０への影響を極力抑える必要がある。本実施形態においては、材料定数（ｍｉやＡｉなど）が未知の状態で、シミュレートしたい粘弾性材料の振動実験を行って実験データを取得し、各基本ライン１１２および調整ライン１１４のべき数ｍｉ、比例定数１／Ａｉを同定する。最適化処理を簡潔にするため、すべての基本ライン１１２のべき数ｍｉは同一であるとしてもよい。たとえば、上述のように、各基本ライン１１２のｍｉの平均値として得られる共通のべき数ｍ（以下、「共通べき数ｍ」とよぶ）をすべての基本ライン１１２に適用してもよい。

調整ライン１１４のべき数ｍｘは、本実施形態においては−１に固定される。理由については後述する。調整ライン１１４の比例定数１／Ａｘに含まれる材料定数Ａｘは、１００×Ａｐに設定される。Ａｐは、基本ライン１１２−０〜１１２−Ｎのうち、もっとも大きな材料定数Ａである。定数Ａｘを大きくすると、式（１７）からも明らかなように調整ライン１１４の緩和時間τｘは非常に小さくなる。緩和時間τｘが小さいことは、粘弾性要素がほとんど効かないことを意味する。定数Ａｘを大きくすることで、調整ライン１１４は特性モデル１１０の全体としての特性にほとんど影響しなくなる。歪速度が小さいときには、緩和時間τｘは特に小さくなる。

図７および図８は、調整ラインのべき数を共通べき数としたときの歪速度と緩和時間の関係を示すグラフである。
横軸は歪速度Ｅ’の対数値を示し、縦軸は緩和時間τの対数値を示す。図７には４種類の基本ライン１１２に対応するグラフ１１２ａ〜１１２ｄと、調整ライン１１４に対応するグラフ１１４が示されている。基本ライン１１２ａ〜１１２ｄの材料定数Ａはそれぞれ０．０１，０．１，１，１０であるとする。また、調整ライン１１４の材料定数Ａは１０００である。式（１７）から明らかなように、各グラフの傾きはべき数ｍを示す。図７，図８では基本ライン１１２と同じく、調整ライン１１４のべき数ｍｘも共通べき数ｍとしているため、すべてのグラフの傾きは同一である。調整ライン１１４は材料定数Ａが大きいため、基本ライン１１２に比べて緩和時間τが小さい。

図８に示す特性影響領域１１６は、粘弾性材料のゴム特性への影響が大きな領域を概念的に示したものである。たとえば、歪速度が１のとき、基本ライン１１２ｄは特性影響領域１１６から外れているため、基本ライン１１２ｄはゴム特性にほとんど影響しない。一方、基本ライン１１２ｃは特性影響領域１１６に含まれているため、基本ライン１１２ｃの粘弾性要素はゴム特性に影響しやすい。調整ライン１１４は、粘弾性材料のゴム特性に一切影響を与えないのが理想である。調整ライン１１４の材料定数Ａを大きく設定することで、調整ライン１１４を特性影響領域１１６から外れるように調整できる。しかし、歪速度が大きくなると、粘弾性材料に含まれる粘弾性要素の影響が無視しがたくなる。図８においては、歪速度が１０００付近になると、調整ライン１１４のグラフが特性影響領域１１６に包含されている。すなわち、歪速度が大きいときには、わずかな緩和時間を有する粘弾性要素であっても、粘弾性材料のゴム特性に影響を及ぼすことがある。

図９は、調整ラインのべき数を−１としたときの歪速度と緩和時間の関係を示すグラフである。
調整ライン１１４Ａはべき数ｍｘを共通べき数ｍとしたときのグラフであり、調整ライン１１４Ｂはべき数ｍｘを−１に固定したときのグラフである。上述のように、べき数ｍｘを基本ライン１１２と同じく共通べき数ｍとした場合、歪速度が大きくなるときに調整ライン１１４Ａと特性影響領域１１６が重なってしまう。そこで本実施形態においては、調整ライン１１４Ｂのようにべき数ｍｘを−１に固定している。粘弾性要素のべき数は−１〜０の範囲にあるため、べき数ｍｘを最小値にしているともいえる。

調整ライン１１４Ｂのべき数ｍｘを最小値にすることで、調整ライン１１４Ｂの傾きを最大化し、調整ライン１１４Ｂを特性影響領域１１６から外す。べき数ｍｘ＝−１と設定するときに傾きが最大になるので最も影響抑制しやすいが、少なくとも、共通べき数ｍよりも小さなべき数ｍｘを設定することでも効果がある。たとえば、共通べき数ｍが−０．５のときには、少なくとも、調整ライン１１４のべき数ｍｘを−０．５未満にするだけでも調整ライン１１４のゴム特性への影響を抑制できる。

すぐり（切欠部分）を有するサスペンションブッシュなど複雑なゴム製品は、局所的に応力が集中することもある。また、ゴム材料の静的な応力−ひずみ特性は非線形特性を示すこともある。実験前にゴムの変形を定量的に予測するのは難しいため、調整ライン１１４がゴム特性に影響を及ぼすかどうかを事前に判断するのは難しい。調整ライン１１４のべき数ｍｘを共通べき数ｍとする場合には、振動実験の結果によって調整ライン１１４の緩和時間τｘが左右されてしまう。本実施形態においては、調整ライン１１４のべき数ｍｘを実験結果に影響されない定数、特に、最小値の−１に設定することで、調整ライン１１４のゴム特性に対する影響を抑制している。

以上、本発明を実施の形態をもとに説明した。実施の形態は例示であり、それらの各構成要素や各処理プロセスの組合せにいろいろな変形例が可能なこと、またそうした変形例も本発明の範囲にあることは当業者に理解されるところである。

１００粘弾性材料特性解析装置、１０２計算部、１０４モデル格納部、１０６モデル生成部、１１０特性モデル、１１２基本ライン、１１４調整ライン、１１６特性影響領域。

Claims

粘弾性材料に変位を与えたときの振幅および弾性率のデータセットを取得し、前記振幅の対数値に対する前記弾性率の対数値の変化率からべき数を算出する計算部と、
前記粘弾性材料の応力−ひずみ特性をシミュレートする特性モデルを前記べき数に基づいて定義するモデル生成部と、を備え、
前記特性モデルにおいては、弾性要素と粘弾性要素が直列接続される１以上の基本ラインと、弾性要素と粘弾性要素が直列接続される調整ラインとが並列接続されており、
前記モデル生成部は、前記基本ラインに含まれる粘弾性要素の緩和時間を歪速度の前記べき数によるべき乗値に比例する値として設定し、前記調整ラインの粘弾性要素の緩和時間を歪速度の逆数値に比例する値として設定することを特徴とする粘弾性材料特性解析装置。