説明

係数算出装置、係数算出方法、及び係数算出プログラム

【課題】簡便にひずみエネルギー密度関数の係数を算出することができる係数算出装置、係数算出方法、及び係数算出プログラムを提供すること。
【解決手段】本実施形態に係る係数算出装置は、大ひずみ粘弾性(Large Strain Viscoelasticity)モデルにおいて、超弾性モデルの構成式に使用されるひずみエネルギー密度関数の係数を材料試験の測定時間に対応する測定データから推定する係数算出装置であって、応力緩和分が除去された超弾性モデルの応力−ひずみ曲線の計算データを逐次的に算出する逐次算出手段と、応力−ひずみ曲線の計算データに基づいて超弾性モデルのひずみエネルギー密度関数の係数を算出する手段と、を備えるものである。

【発明の詳細な説明】
【技術分野】
【0001】
本発明は係数算出装置、係数算出方法、及び係数算出プログラムに関し、特に詳しくは、超弾性モデルのひずみエネルギー密度関数の係数を算出する係数算出装置、係数算出方法、及び係数算出プログラムに関する。
【背景技術】
【0002】
ゴムやエラストマーのような樹脂材料は大変形が可能であり、超弾性材料或いは超弾性体と呼ばれる。この超弾性材料の力学構成式(以下、超弾性モデル)は、下の式(1)に示すように、ひずみエネルギー密度関数Wを右Cauchy−Kirchhoff変形テンソルCで偏微分することにより、第2Piola−Kirchhoff応力テンソルSが得られるというものである。
【0003】

【0004】
しかし、右Cauchy−Kirchhoff変形テンソルCは、その固有値である主ストレッチを対角成分に持つテンソルであり、ひずみエネルギー密度関数Wは固有値の関数であることから、超弾性モデルの構成式は時間の変化を含まない。そのため、変形速度が変わると応力−ひずみ曲線も変化すると言う時間依存特性を表現することができない。
【0005】
これに対し、時間や温度依存特性を考慮可能な線形粘弾性モデルに超弾性モデルを組み込んだ大ひずみ粘弾性(Large Strain Viscoelasticity)モデルが提案されている(非特許文献1、2)。非特許文献2によれば、大ひずみ粘弾性モデルの構成式は下の式(2)である。
【0006】

【0007】
式(2)において、tは現在の時刻、τは時間0からtまでの時間、α、α、α、αは、時刻t=0の緩和弾性率をG、K、Prony級数の係数をG、K、G、Kとするとα=G/G、α=G/G、α=K/K、α=K/K、である。また、Jは体積変化率である。
【0008】
式(2)の第1項の積分は偏差成分を、第2項の積分は体積成分を表す。式(2)において、偏差成分、体積成分とも、線形粘弾性による応力緩和率(式(2)の[]内)に超弾性モデルの応力の時間増分(式(2)のd/dτ()の項)を乗じたものを現在の時刻tまで積分するというモデルである。
【0009】
大ひずみ粘弾性モデルにおいて、超弾性モデルのひずみエネルギー密度関数を求めるには、超弾性モデル単体のひずみエネルギー密度関数の係数を求めるのと同様、単軸引張試験、等方2軸引張試験、純せん断試験、体積試験(以下、まとめて材料試験)の応力−ひずみ曲線の測定値から算出する。超弾性モデル単体の場合、非特許文献3に示すように、応力−ひずみ曲線の測定値をそのまま使用することが可能であるが、大ひずみ粘弾性モデルの超弾性モデルの場合には、応力−ひずみ曲線の測定値に線形粘弾性の応力緩和の影響が含まれているため、応力緩和の影響を除去した応力−ひずみ曲線のデータを作成する必要がある。
【0010】
特許文献1に開示された粘弾性材料特性解析システムは、粘弾性材料を変形させる複数の変形モードに対応する実験データ、及び応力緩和特性値データをモデルに取り込んで、材料定数を算出する。また、粘弾性材料特性解析の解析結果は、応力−ひずみ線図としてグラフ化される。
【0011】
特許文献2に開示された応力ひずみ解析システムは、負荷によって生じる塑性ひずみ、粘弾性ひずみ、粘性ひずみを用いて、応力−ひずみ曲線を表現し、応力ひずみ解析から求めた応力−ひずみ曲線で囲まれる面積を、解析モデルの要素分割によって算出する。
【先行技術文献】
【特許文献】
【0012】
【特許文献1】特開2003−083874号公報
【特許文献2】特開2003−270060号公報
【非特許文献】
【0013】
【非特許文献1】G.A. Holzapfel, "ON LARGE STRAIN VISCOELASTICITY:CONTINUUM FORMULATION AND FINITE ELEMENT APPLICATIONS TO ELASTOMERIC STRUCTURES", Int. J. Numer. Meth. Eng.,Vol. 39, pp. 3903−3926 (1996)
【非特許文献2】ANSYS. Inc, "ANSYS Release 9.0 Documentation ANSYS, Inc. Theory Manual 4.7.7 Large Strain Viscoelasticity"
【非特許文献3】SIMULIA社,"Abaqus Version 6.7 Theory Manual 4.6"Large−strain elasticity"
【発明の概要】
【発明が解決しようとする課題】
【0014】
ひずみエネルギー密度関数の係数は、汎用シミュレータを利用して、算出することができる。汎用シミュレータなどを用いた場合、図12に示すように、3つの処理装置を用いた方法で係数を算出することができる。以下の3つの処理装置は、例えば、パーソナルコンピュータ等の演算処理装置である。
【0015】
第1の処理装置100は、入力部101と、計算部102と、出力部104とを有している。入力部101は、材料試験の応力−ひずみ曲線のデータを入力する。計算部102は、超弾性モデルのひずみエネルギー密度関数の係数を計算する算出手段121を有している。出力部104は、算出した係数を出力する。
【0016】
第2の処理装置200は、入力部201と、計算部202と、出力部204とを有している。入力部201は、第1の処理装置100の出力結果から、材料試験の応力シミュレーションを行うための設定を入力する。計算部202は、入力された設定にしたがって材料試験の応力シミュレーションを行う。計算部202は、応力シミュレーションによって応力―ひずみ曲線を算出する算出手段122を有している。出力部204は、シミュレーション結果から得られた応力−ひずみ曲線のデータを出力する。
【0017】
第3の処理装置300は、入力部301と、計算部302と、計算部402と、出力部304とを有している。入力部301は、第1の処理装置100の入力部101に入力した材料試験の応力−ひずみ曲線のデータと、第2の処理装置200の出力部204から出力された応力−ひずみ曲線のデータを入力する。計算部302は、2つの応力−ひずみ曲線のデータから、各ひずみでの応力緩和分を計算する算出手段123を有している。また、計算部402は、算出した応力緩和分を入力部101に入力した材料試験の応力から除去して応力緩和分を除いた応力−ひずみ曲線のデータを算出する算出手段124を有している。出力部304は、応力緩和分を除いた応力−ひずみ曲線のデータを出力する。
【0018】
出力部304から出力された応力−ひずみ曲線のデータは、第1の処理装置100の入力データとなり、計算部302で応力緩和率が許容範囲に入るまで計算が繰り返される。そのため、繰り返し計算の最初に第1の処理装置100の入力部101に入力されるのは材料試験の測定値である。また、計算部202にはANSYS(登録商標)などの汎用シミュレータが使用され、計算部102には汎用シミュレータに付属する係数算出プログラムなどが使用される。さらに、計算部302には、汎用表計算プログラムが使用される。
【0019】
このような構成を有する処理装置は次のように動作する。図13は、当該処理装置の動作フローを示す図である。
【0020】
図13に示すように、入力部101には試験データを入力する(ステップA101)。算出手段121は、入力部101から入力された材料試験の応力−ひずみ曲線のデータから線形最小2乗法や非線形最小2乗法により超弾性モデルのひずみエネルギー密度関数の係数を算出する(ステップA102)。
【0021】
入力部201は、応力シミュレーションに必要な設定を入力する(ステップA201)。算出手段122は、入力部201から入力された設定に従って有限要素法により材料試験の応力シミュレーションを行い、その応力−ひずみ曲線を計算する(ステップA202)。そして、出力部204は算出手段122で算出された応力−ひずみ曲線のデータを出力する。
【0022】
入力部301は、2つの応力−ひずみ曲線を入力する(ステップA301)。算出手段123は、入力部101に入力した応力−ひずみ曲線のデータと出力部204から出力された応力−ひずみ曲線のデータとを比較して、応力緩和率を計算する(ステップA302)。次に、応力緩和率が十分小さいかどうかを判断する(ステップA303)。十分小さい場合、算出手段121により算出した超弾性モデルのひずみエネルギー密度関数の係数を求める解とする。十分小さくない場合、算出手段124により、入力部101に入力した応力−ひずみ曲線のデータにおいて、応力緩和分を除去した応力を算出し、応力を新たな応力とした応力−ひずみ曲線のデータを算出する(ステップA402)。さらに、新たに作成した応力−ひずみ曲線のデータを入力部101に入力して、上述の作業を繰り返す。
【0023】
上記手法の第1の問題点は、処理算装置が別になっているため、処理装置間のデータの入出力に人手が掛かることである。その理由は、大ひずみ粘弾性の超弾性モデルのひずみエネルギー密度関数の係数の算出には、材料試験の応力−ひずみ曲線から応力緩和分を除去する必要があり、前記応力緩和分を除去した超弾性モデルの応力−ひずみ曲線のデータを算出する装置が存在しないためである。第2の問題点は、第1の問題点に加えて繰り返し計算を要することである。これにより、より多くの人手が掛かることになる。この理由も、応力緩和率を除去した超弾性モデルの応力−ひずみ曲線のデータを算出する装置が存在しないためである。従って、上記の算出方法では、係数を算出するための処理が煩雑になってしまうという問題点がある。
【0024】
また、特許文献1は、単に、実験データ、及び応力緩和特性値データを取り込むのみであり、ひずみが数%である粘弾性材料を解析対象としている。そのため、ひずみが数百%に及ぶ大ひずみ粘弾性材料、超弾性材料を扱うことはできない。
【0025】
さらに、特許文献2は、グラフの面積を求めてひずみエネルギーの値自体を算出しているが、大ひずみ粘弾性モデルの構成式の係数を算出することはできない。
【0026】
本発明は、このような問題点を鑑みてなされてものであって、簡便にひずみエネルギー密度関数の係数を算出することができる係数算出装置、係数算出方法、及び係数算出プログラムを提供することを目的とする。
【課題を解決するための手段】
【0027】
本発明にかかる係数算出装置は、大ひずみ粘弾性(Large Strain Viscoelasticity)モデルにおいて、超弾性モデルの構成式に使用されるひずみエネルギー密度関数の係数を材料試験の測定時間に対応する測定データから推定する係数算出装置であって、応力緩和分が除去された超弾性モデルの応力−ひずみ曲線の計算データを逐次的に算出する逐次算出手段と、前記応力−ひずみ曲線の計算データに基づいて超弾性モデルのひずみエネルギー密度関数の係数を算出する手段と、を備えるものである。
【0028】
本発明にかかる係数算出方法は、大ひずみ粘弾性(Large Strain Viscoelasticity)モデルにおいて、超弾性モデルの構成式に使用されるひずみエネルギー密度関数の係数を材料試験の測定時間に対応する測定データから推定する係数算出方法であって、応力緩和分が除去された超弾性モデルの応力−ひずみ曲線の計算データを逐次的に算出する逐次算出ステップと、前記応力−ひずみ曲線の計算データに基づいて超弾性モデルのひずみエネルギー密度関数の係数を算出するステップと、を備えるものである。
【0029】
本発明にかかる係数算出プログラムは、大ひずみ粘弾性(Large Strain Viscoelasticity)モデルにおいて、超弾性モデルの構成式に使用されるひずみエネルギー密度関数の係数を材料試験の測定時間に対応する測定データから推定する処理をコンピュータに実行させる係数算出プログラムであって、コンピュータに対して、応力緩和分が除去された超弾性モデルの応力−ひずみ曲線の計算データを逐次的に算出させる逐次算出ステップと、前記応力−ひずみ曲線の計算データに基づいて超弾性モデルのひずみエネルギー密度関数の係数を算出させるステップと、を備えるものである。
【発明の効果】
【0030】
本発明により、簡便にひずみエネルギー密度関数の係数を算出することができる係数算出装置、係数算出方法、及び係数算出プログラムを提供することができる。
【図面の簡単な説明】
【0031】
【図1】本発明にかかる係数算出装置の構成を示すブロック図である。
【図2】実施の形態1にかかる係数算出装置の構成を示すブロック図である。
【図3】実施の形態1にかかる係数算出方法を示すフローチャートである。
【図4】実施の形態2にかかる係数算出装置の構成を示す図である。
【図5】実施例で使用する単軸引張試験結果である。
【図6】実施例で使用するProny級数の緩和時間及び係数である。
【図7】実施例で使用するシフトファンクションの係数である。
【図8】実施例の計算の途中結果で、Prony級数の緩和時間と時刻0の緩和弾性率で図6の係数を除したものである。
【図9】実施例の計算の途中結果で、応力緩和分を除去した超弾性モデルの応力−ひずみ曲線のデータである。
【図10】実施例の計算によって得られた、Mooney−Rivlinの係数を示す表である。
【図11】本実施例の計算結果から商用シミュレータで得られた応力―ひずみ曲線と実測の材料試験結果の比較を表す図である。
【図12】商用シミュレータを用いた係数算出装置の構成を示すブロック図である。
【図13】商用シミュレータを用いた係数算出方法を示すフローチャートである。
【発明を実施するための形態】
【0032】
以下では、本発明を適用した具体的な実施の形態について、図面を参照しながら詳細に説明する。各図面において、同一要素には同一の符号が付されており、説明の明確化のため、必要に応じて重複説明は省略する。
【0033】
時間及び温度依存特性を考慮可能な線形粘弾性モデルに超弾性モデルを組み込んだ大ひずみ粘弾性(Large Strain Viscoelasticity)モデルにおいて、超弾性モデルの構成式に使用されるエネルギー密度関数の係数を材料試験の測定値から推定する。すなわち、大ひずみ粘弾性モデルで解析を行う際に、材料試験の測定データに基づいて、ひずみエネルギー密度関数の各係数を算出する。
【0034】
図1に本発明にかかる係数算出装置10の構成を示す。係数算出装置10は、応力緩和分が除去された超弾性モデルの応力−ひずみ曲線の計算データを逐次的に算出する手段11と、応力−ひずみ曲線の計算データに基づいて超弾性モデルのひずみエネルギー密度関数の係数を算出する手段12と、を備えている。
【0035】
係数を算出するために、まず、測定値のデータを係数算出装置に入力する。すなわち、試験片に対して引張試験などを行い、応力に対するひずみを測定する。これにより、応力―ひずみ曲線の測定データを取得することができる。そして、時間と、応力と、ひずみとが組となった測定データを取得する。この測定データにおける応力は応力緩和分を含んでいる。そのため、手段11は、逐次的な算出により応力緩和分が除去された超弾性モデルの応力−ひずみ曲線の計算データを算出する。
【0036】
すなわち、手段11は、入力された応力―ひずみ曲線の測定データを用いて、当該測定データに該当する時間における応力―ひずみ曲線の応力(計算データ)を算出する。そして、手段11は、算出された計算データと、次の時間に該当する測定データと、及び時間差分とを用いて次の時間に該当する計算データを逐次的に算出する。手段11は、式(2)の積分を行う代わりに、隣接する時間における時間差分と測定データとから逐次的に計算データを算出することにより、応力緩和がない状態での超弾性モデルの応力の計算データを求めることができる。
【0037】
手段12は、応力−ひずみ曲線の計算データに基づいて超弾性モデルのひずみエネルギー密度関数の係数を算出する。最小2乗法などによって、例えば、Mooney−Rivlinモデルの係数を算出する。
【0038】
このように大ひずみ粘弾性モデルにおいて解析を行う際に、ひずみエネルギー密度関数の係数を用いずに、測定データと測定時間の差分とから該当する時間の超弾性モデルの応力―ひずみ曲線の計算データを逐次的に算出している。すなわち、応力―ひずみ曲線から直接に応力緩和分を除去した計算データを求める。そして、応力緩和分が除去された応力―ひずみ曲線の計算データからエネルギー密度関数の係数を算出している。このようにすることで、材料固有の力学的特性を数学的に表現する構成式の一つであるエネルギー密度を簡便に算出することができる。
【0039】
実施形態1.
本実施形態にかかる係数算出方法、及び装置について、図2を用いて説明する。図2は、本実施形態にかかる係数算出装置の構成を示すブロック図である。図2を参照すると、係数算出装置は、キーボードやマウス等の入力部1と、係数の算出を行う計算部2と、計算時にデータを保持するための記憶部3と、ディスプレイや印刷機あるいはファイル装置等の出力部4とを有している。
【0040】
計算部2は、入力データ読み込み手段21と、シフトファクターを算出する手段22(以下、第一算出手段22)と、時刻0の緩和弾性率でProny級数の係数を除したα、αを算出する手段23(以下、第二算出手段23)と、応力緩和分を除去した超弾性モデルの応力−ひずみ曲線のデータを逐次的に計算する手段24(以下、第三算出手段24)と、ひずみエネルギー密度関数の係数算出手段25(以下、第四算出手段25)とから構成される。また、計算部2に接続された記憶部3には、3つの記憶部31〜33が設けられている。
【0041】
なお、これらの算出手段や記憶部は、物理的に同一の構成であってもよい。係数算出装置は、データの入出力が可能なパーソナルコンピュータ等の演算処理装置である。例えば、CPU(Central Processing Unit)、ROM(Read Only Memory)、RAM(Random Access Memory)、記憶ディスク、通信用のインターフェースなどを有し、係数算出プログラムにしたがって各種演算処理を実行する。RAMや記憶ディスクなどが各記憶部となる。また、CPUなどのプロセッサが各算出手段となる。そして、プロセッサなどは、予め格納されたコンピュータプログラムを実行することで、それぞれの処理を行う。従って、これらの算出手段は物理的に同一のプロセッサなどで構成されていても良い。例えば、記憶ディスクなどに記憶されている演算プログラムをCPUによって実行されると、入力されている測定値や設定値をRAMなどに読み込んで、各種演算を行う。そして、その計算過程のデータや計算結果のデータをRAMや記憶ディスクなどに書き込む。このように、演算処理装置に予め記録されている係数算出プログラムにしたがって、下記の演算処理が実行される。
【0042】
続いて、上述した構成のそれぞれの動作の概略を以下に説明する。
【0043】
入力データ読み込み手段21は、入力部1によって与えられた数値データを読み込み、記憶部3に書き込む。これにより、数値データが記憶部3に保存される。入力部1から与えられる数値データは、例えば、以下の4つである。
【0044】
(A)シフトファンクションの係数、(B)Prony級数の緩和時間とそれに対応する係数、(C)材料試験の測定値(試験データ)、(D)(C)の試験温度、である。(A)〜(D)の数値データは、入力部1によって入力される。また、これらは、作業者によって入力されてもよいし、他の装置から転送されてもよい。
【0045】
(A)、(D)のデータは第一記憶部31に、(B)、(C)のデータは第二記憶部32に記憶される。第一記憶部31は、入力されたシフトファンクションの係数と、試験温度とを設定値として格納する。また、第二記憶部32は、Prony級数の緩和時間とそれに対応する係数を設定値として格納し、材料試験の測定値(試験データ)を測定データとして格納する。
【0046】
ただし、(A)で、シフトファンクションの係数を入力する代わりに、シフトファクターaを入力しても良い。その場合、後述の第一算出手段22、及び(D)の試験温度の入力は不要となる。また、(B)で、Prony級数の代わりに時刻0の緩和弾性率でProny級数の係数を除したα、αを読み込んでも良い。この場合、後述のα、αiを算出する第二算出手段23は不要となる。
【0047】
また、(C)材料試験の測定値は、試験結果のひずみε(k=1・・・M)と、ひずみに対応する応力σと、時間tと、の組のデータ(以下、測定データ)である。したがって、測定データとして、ひずみと応力と時間とから構成されるデータセットがM組用意される。(C)について、材料試験の負荷を等ひずみ速度で与えた場合、各ひずみに対応する時間tの代わりに、ひずみ速度を与えてもよい。この場合、計算部2がtを計算しても良い。さらに、試験結果をひずみではなくストレッチとしても良い。公称ひずみεと主ストレッチλはλ=1+εの関係があるためである。主ストレッチλが公称ひずみεに換算される。
【0048】
第一算出手段22は、第一記憶部31からシフトファンクションの係数、試験温度を読み出す。そして、入力データ読み込み手段21は、これらのデータをシフトファンクションの式に代入することにより算出する。すなわち、第一算出手段22は、シフトファンクションの係数、及び試験温度からシフトファンクションを算出する。シフトファンクションの式は、記憶部3に予め登録されている。
【0049】
第二算出手段23は時刻0の緩和弾性率でProny級数の係数を除したα、αを算出する。そのため、第二算出手段23は第二記憶部32からProny級数の係数G、G(i=1・・・n)を読み取る。そして、第二算出手段23は、各Prony級数の係数を時刻0の緩和弾性率G=G+G+G+・・・+Gで除して、α=G/G、α=G/G(i=1・・・n)によりα、αを得る。
【0050】
第三算出手段24は、超弾性モデルの応力−ひずみ曲線の計算データを逐次的に算出する。そのため、第三算出手段24は、第一算出手段22で算出したシフトファクターaと、第二記憶部32から読み取ったProny級数の緩和時間と、第二算出手段23で算出したα、αと、測定データと、から測定データの各応力σに対応する超弾性モデルの応力を算出する。その際、第三算出手段24は、超弾性モデルの応力の初期値と、次の時間の測定データの応力σと、時間差分とから、逐次的に算出することで、応力緩和分を除去した応力−ひずみ曲線のデータを算出する。超弾性モデルの応力−ひずみ曲線のデータの逐次的算出は、次のようにして行われる。但し、以下では式(2)の偏差成分(式(2)の第1項)に対して述べるが、体積成分(式(2)の第2項)についても同様にして求めることができる。
【0051】
式(2)の第1項の積分項を展開すると、下の式(3)のようになる。但し、試験温度によりProny級数はシフトすることから、式(3)は、シフトファクターaも考慮した式となっている。

ここで、ひずみは時間と対応するため、式(2)の右辺の第1項のd/dτ以降の()内(超弾性モデルの第2Piola−Kirchhoff応力)をS(t)と表記する。
【0052】
ここで、α、α(i=1・・・n)は、時刻0の緩和弾性率でProny級数の係数を除したもので、前記のα、α(i=1・・・n)である。また、tは現在の時刻、τは時間0からtまでの時刻、τiはαiに対応する緩和時間である。
【0053】
式(3)のΣの中の項を下の式(4)のようにs(t)と置く。

【0054】
式(4)において、s(t)をsすると、時刻tk+1=t+Δtでは次の式(5)のように表すことができる。

【0055】
式(5)の第2項の積分に中点則を適用すると、次の式(6)のようにsを逐次的に表すことができる。

【0056】
式(3)でt=tとし、式(6)を式(3)に代入して整理すると、時刻tの超弾性モデルの第2Piola−Kirchhoff応力の値S(t)は式(6)のsikを用いて次の式(7)のように逐次的に求めることができる。

【0057】
測定データは公称応力―公称ひずみのデータとして与えられることが多い。式(2)が応力緩和率(式(2)の[]内)に超弾性モデルの時間増分(式(2)のd/dτ()の項)を乗じたものを積分しているため、d/dτ()の第2Piola−Kirchhoff応力Sを公称応力である第1Piola−Kirchhoff応力Πに置換することにより、式(2)は、公称応力を求める式となる。当該公称応力を求める式に対して、上述した式(3)乃至(7)と同様に式変形することにより、式(6)、(7)に対応する公称応力の逐次計算式(6')、(7')を以下の通り得ることができる。


【0058】
但し、Π(t)は時刻tでの超弾性モデルの応力の値、Π(t)は時刻tでの公称応力の値(=測定データの値)である。つまり、所定の測定時間である時刻tに対応する超弾性モデルの応力−ひずみ曲線の計算データΠ(t)は、式(7')により、時刻tに対応する測定データΠ(t)と、時刻tより短い測定時間である時刻tk−1に対応する超弾性モデルの応力−ひずみ曲線の計算データであるΠ(tk−1)と、時刻tと時刻tk−1との時間差分であるΔtとに基づいて、k=0から測定時間の最大値までについて逐次的に算出することができる。また、式(7')に含まれるπ(tk−1)は、式(6')により、超弾性モデルの応力−ひずみ曲線の計算データΠ(tk−1)及びΠ(tk−2)の差分に基づいて求めることができる。つまり、計算データΠ(t)は、時刻tより短い測定時間である時刻tk−1及び時刻tk−2に対応する計算データΠ(tk−1)及びΠ(tk−2)の差分から求めることができる。同様に、πは、式(7')により求められた計算データΠ(t)と及びΠ(tk−1)の差分に基づいて求めることで、逐次的に算出することができる。尚、πは、次の式(4')である。

【0059】
但し、式(6')はπの漸化式として解くことができる。式(6')を解くと、次の式(6'')となる。

さらに、式(6'')を式(7')に代入することで次の式(7'')を得る。

【0060】
式(7'')を用いると、所定の測定時間である時刻tの超弾性モデルの応力−ひずみ曲線の計算データΠHG(tk)は、所定の測定時間である時刻tより短い時間の超弾性モデルの応力−ひずみ曲線の計算データであるΠ(t)(j=1・・・k−1)に基づいて求めることができる。式(7'')を用いる場合、式(7')と式(6' ')を用いた場合のようにπk−1を求めるためにπk−2の値を保持しておく記憶領域が不要となる。
【0061】
測定データの(公称)応力σから超弾性モデルの応力Π(t)を求めるには、まず、測定データを時間に関して昇順で整列する。k=0を時刻t=0とし、材料試験では時刻0の応力が0であるためΠ(t)=0、π=0とする。次に、式(7')でΠ(t)にσkを代入してΠ(t)を求めた後、式(6')によりπを更新することをk=1から繰り返し、逐次的に全て応力のデータσ (k=1・・・M)を超弾性モデルの応力Π(t)(k=1・・・M)に変換する。但し、式(7'')を用いる場合、πを0とする必要もπk−1を更新する必要も無く、時刻0の超弾性モデルの応力Π(t)=0として、逐次超弾性モデルの応力Π(t)(k=1・・・M)へ変換することができる。
【0062】
また、第四算出手段25は、第三算出手段24で求めた超弾性モデルの応力−ひずみのデータ(ε)(k=1・・・M)からひずみエネルギー密度関数の係数を算出する。ここでは、線形最小2乗法や非線形最小2乗法を用いて、係数を算出する。これにより、超弾性モデルの構成式となるひずみエネルギー密度関数の係数を算出することができる。そして、第四算出手段25は係数を出力部4に出力する。あるいは、第四算出手段25は係数を記憶部3に書き込む。このように、簡便に係数を算出することができる。
【0063】
次に、図3を参照して本実施形態の全体動作について詳細に説明する。図3は、本実施の形態にかかる係数算出方法を示すフローチャートである。また、図3の各ステップには、図2で示した各算出部や各記憶部などが合わせて示されている。
【0064】
まず、入力部1から与えられた材料試験結果のひずみε、それに対応する応力σと時間t(k=1・・・M)、Prony級数の緩和時間τ(i=1・・・n)とそれに対応する係数G、G、シフトファンクションの係数C、C、・・・(個数はシフトファンクションに依存)、試験温度Teは、入力データ読み込み手段21を通じて読み込まれる。そして、読み込まれたシフトファンクションの係数と試験温度を第一記憶部31に記憶し、Prony級数と測定データを第二記憶部32に記憶する (図3のステップA1)。
【0065】
次に、第一記憶部31からシフトファンクションの係数と試験温度を読み取り、シフトファンクションの式からシフトファクターaを決定する(ステップA2)。シフトファンクションはT−N式やWLF式などが使用可能であり、例えばT−N式では2つの係数C、Cを必要とし、次の式からシフトファクターが求められる。

【0066】
シフトファクターaの算出後、Prony級数の緩和時間τ(i=1・・・n)とそれに対応する係数G、G、測定データε、σ、t(k=1・・・M)を第二記憶部32から読み込む(ステップA3)。次に、G=G+G+G+・・・+Gとして、α=G/G、α=G/G(i=1・・・n)により、α、α(i=1・・・n)を求める(ステップA4)。なお、ステップA3とステップA4は順序が逆でも良い。
【0067】
その後、測定データをtk−1<tとなるように整列する(ステップA5)。その際、k=0を時刻0とする。そのため、t=0である。また、時刻tの測定データの応力をσ = 0、超弾性モデルの応力の初期値をΠ(t)=0、式(4')のπの初期値π(i=1・・・n)を0とする。但し、式(7'')を用いる場合、πは使用しない。
【0068】
次に、k=1として、測定データの応力σ(k=1・・・M)を超弾性モデルの応力Π(t)に変換することを開始する(ステップA6)。
【0069】
k番目のデータに対し、前記シフトファクターa、Prony級数の係数を時刻0の緩和弾性率で除したα、α(i=1・・・n)とそれに対応する緩和時間τ、τ、k−1番目とk番目のデータの時間の差Δt=t−tk−1、k−1番目の超弾性モデルの応力Π(tk−1)、及びπk−1から、式(7')を用いて超弾性モデルの応力Π(t)を求め、求めた超弾性モデルのデータ(εk, Π(t))を記録部33に記録する(ステップA7)。
【0070】
次に、式(6')及び(7')を用いる場合、式(6')によりπを更新する(ステップA8)。式(7'')を用いる場合、ステップA8は不要である。さらにkの値が測定データ数に達したかを確認する(ステップA9)。測定データ数に達しない場合、kを1つ増やして(ステップ10)ステップA7へ戻る。
【0071】
一方、kが測定データ数に達して全てのデータの処理を終えた場合、記憶部33から超弾性モデルのひずみεと応力Π(t)(k=1・・・M)を読み取り、前記読み取った応力−ひずみ曲線のデータから超弾性モデルのひずみエネルギー密度関数の係数を算出する(ステップA11)。
【0072】
次に、本実施の形態の効果について説明する。
【0073】
本実施の形態では、測定データから逐次計算により応力緩和分を除去した応力−ひずみ曲線を算出し、前記応力−ひずみ曲線に対して大ひずみ粘弾性モデルの超弾性モデルのひずみエネルギー密度関数の係数を算出している。このため、係数算出工数を低減できる。また、演算処理装置間のデータのやり取り回数を低減することができるため、簡便に係数を算出することができる。
【0074】
また、本実施の形態では、ひずみエネルギー密度関数の係数を求めることなく、測定データと測定時間の差分とから該当する時間の超弾性モデルの応力―ひずみ曲線の計算データを逐次的に算出している。すなわち、応力緩和分を除去した構造式から逐次的に算出することで、測定データとその時間差分から順次、計算することが可能となる。よって、応力緩和分の誤差が小さくなるまで、繰り返し、ひずみエネルギー密度関数の係数、及び応力―ひずみ曲線を算出する必要がなくなる。応力シミュレーションによってひずみエネルギー密度関数の係数から応力―ひずみ曲線を算出するための商用のシミュレータが不用となるため、より簡便に係数を求めることができる。応力緩和分の誤差が小さくなるまでひずみエネルギー密度関数の係数を繰り返し算出することがなくなるため、計算時間を短縮することができる。よって、簡便に係数を算出することができ、利便性を向上することができる。
【0075】
実施の形態2.
図4を参照すると、本発明の実施形態2は、インターネット等のネットワーク5を通じてサーバコンピュータ6とクライアントコンピュータ7とが接続されたネットワークシステムにおいて、ネットワーク5を通じてクライアントコンピュータ7から係数算出要求をサーバコンピュータ6へ送信し、この要求を受けたサーバコンピュータ6がネットワーク5を通じて係数算出処理を行って求めた係数の値を要求元のクライアントコンピュータ7へ返却するサービスを有料または無料で提供する。
【0076】
サーバコンピュータ6には、図2の計算部2と記憶部3とが設けられている。クライアントコンピュータ7には、図2の入力部1と出力部4とが設けられている。クライアントコンピュータ7とサーバコンピュータ6がネットワーク5を通じて通信可能に接続される。クライアントコンピュータ7に設けられた入力部1及び出力部4と、サーバコンピュータ6に設けられた計算部2とがネットワーク5を通じて通信可能に接続される。従って、クライアントコンピュータ7の入力部1で入力された各データがサーバコンピュータ6の記憶部3に記憶される。また、サーバコンピュータ6の計算部2で計算されたデータが、クライアントコンピュータ7の出力部4に出力される。
【0077】
次に本実施の形態の動作を説明する。
【0078】
クライアントコンピュータ7のユーザは、ネットワーク5を通じてサーバコンピュータ6にクライアントコンピュータ7を接続し、材料試験結果のひずみε、それに対応する応力σと時間t(k=1・・・M)、Prony級数の緩和時間τ(i=1・・・n)とそれに対応する係数G、G、シフトファンクションの係数C、C、・・・、試験温度Tをサーバコンピュータ6へ送信する。
【0079】
サーバコンピュータ6は、係数算出要求を受信すると、それに含まれる材料試験結果のひずみε、それに対応する応力σと時間t(k=1・・・M)、Prony級数の緩和時間τ(i=1・・・n)とそれに対応する係数G、G、シフトファンクションの係数C、C、・・・、試験温度Tを記憶装置に記憶し、第1の実施の形態における計算部2と同様の処理を実行して係数を算出する。そして、サーバコンピュータ6は、ネットワーク5を通じて、それらの値を含む応答メッセージをクライアントコンピュータ7へ送信する。
【0080】
クライアントコンピュータ7は、サーバコンピュータ6から応答メッセージを受信すると、それに含まれる係数の値を記憶部に記憶する。記憶された係数の値は、各種のシミュレーションに利用される。このように、物理的に異なる処理装置を用いた場合でも、実施の形態1と同様の効果を得ることができる。
【0081】
実施例
次に、具体的な実施例を説明する。
【0082】
図5にエラストマー樹脂の単軸引張試験データ(試験温度T=125℃、ひずみ速度0.1%/secで一定)を示す。図6にProny級数の緩和時間τ(i=1・・・39)と対応する係数G、G、を示す。図5及び図6では、各データを2段で示している。また、シフトファクター算出に使用するシフトファンクションはT−N式とし、図7にシフトファンクションの係数を示す。さらに、算出する超弾性モデルのひずみエネルギー密度関数Wは下の式(9)のようにMooney−Rivlinモデルとし、その次数は3とする(9パラメータのMooney−Rivlinモデル)。算出する係数はAi,jである。但し、入力に使用する材料試験が単軸引張試験であるため、偏差成分しか求めることができないことから、式(9)は偏差成分のみを記す。なお、式(9)において、Iは、偏差右Cauchy−Green変形テンソルの不変量を示す。
【0083】

【0084】
入力データの読み込み手段21は、測定データのひずみεkと対応する応力σと時間t(k =1・・・65)、Prony級数の緩和時間τ(i=1・・・39)と対応する係数G、G、シフトファンクションの係数C、C、及び試験温度T=125℃を読み込み、シフトファンクションの係数と試験温度を記憶部31に、測定データとProny級数を第二記憶部32に記憶する(図3のステップA1)。
【0085】
次に、第一記憶部31からシフトファンクションの係数C、Cと試験温度Tを読み取り、T−N式からシフトファクターaを算出する。(ステップA2)。但し、T−N式ではC、Tを絶対温度で入力する必要があることから、絶対温度へ変換するために273℃を加えている。式(10)は、シフトファクターaの算出過程及び算出値である。
【0086】

【0087】
シフトファクターaの算出後、Prony級数の緩和時間τ(i=1・・・39)と対応する係数G、G、及び測定データのひずみεk、応力σ、時間t(k=1・・・65)を第二記憶部32から読み込む(ステップA3)。
【0088】
その後、上記のように、図7のデータを用いてα、αを算出する(ステップA4)。すなわち、α=G/G、α=G/G(i=1・・・n)によりα、αiを得る。算出したα、αを図8に示す。なお、図8では、データを2段で示している。
【0089】
次に、測定データを時間が早い順に整列する(ステップA5)。そして、測定データのカウンタであるkを1として、各測定データの応力から超弾性モデルの応力を逐次的に算出することを開始する(ステップA6)。
【0090】
まず、1番目のデータについて、式(7')又は式(7'')から超弾性モデルの応力Π(t)を求め、前記Π(t)を対応するひずみε1との組として記憶部33に記録する(ステップA7)。
【0091】
前述のステップA7で式(7')を用いた場合、次に式(6')によりπを計算する(ステップA8)。
【0092】
その後、ステップA9によりkの値が測定データ数65に達したかを確認する。k=1ゆえまだ測定データ数に達していないことから、ステップA10に進んでk=2としてからステップA7に戻る。
【0093】
測定データのカウンタkが測定データ数65に達するまで上記の操作を繰り返す。これによって算出した超弾性モデルの応力−ひずみ曲線のデータを図9に示す。なお、図9では、データを2段で示している。
【0094】
測定データのカウンタkが測定データ数65に達した場合、第三記憶部33から応力緩和分を除去した超弾性モデルの応力−ひずみ曲線のデータε、Π(k =1・・・65)を読み出し、線形最小2乗法によりMooney−Rivlinの係数Ai,jを算出する(ステップA11)。算出したMooney−Rivlinの係数を図10に示す。すなわち、本発明の手法で得た解を図10に示す。
【0095】
また、本発明の手法で得た係数を用いた際の商用シミュレータの応力−ひずみ曲線と材料試験結果との比較を図11に示す。図11から、本手法で得た係数は実用上十分な精度を有することが分かる。このように、簡便且つ、高精度に係数を算出することができる。
【0096】
以上、本発明の実施の形態および実施例について説明したが、本発明は以上の例に限定されず、その他各種の付加変更が可能である。本発明の実施例では、偏差成分の係数算出に適用したが、体積成分の算出にも同様に適用可能である。また、単軸引張試験の測定データから係数を算出したが、等方2軸引張試験、純せん断試験、体積試験等の他の材料試験の測定データを用いても同様に係数を算出できる。加えて、ひずみエネルギー密度関数をMooney−Rivlinとし、次数を3としたが、3以外の次数でも良く、さらに他のひずみエネルギー密度関数を用いても良い。例えば、Neo−HookeモデルやOgdenモデルやArruda−Boyceモデルのひずみエネルギー密度関数の係数を求めてもよい。
【0097】
また、本発明の実施例では、ひずみ速度が一定であることとしたが、これに限定されず、ひずみと時間の対応が測定できれば良い。
【0098】
尚、ここで、本発明と上述した特許文献1及び特許文献2との差異を述べる。特許文献1では、データを取り込むだけであり、自然法則(構成式)に基づいた数値処理を行っていない。そのため、本発明における応力緩和分の除去を行うことはできない。また、特許文献2では、グラフの面積を求めるだけであり、応力−ひずみ曲線にフィットするよう自然法則(構成式)の関数の係数を算出することはできない。
【0099】
本発明によれば、力学構成式の積分項に含まれる力学モデルの係数を、材料試験の測定値から推定する分野に利用できる。特に、大ひずみ粘弾性モデルの超弾性モデルのひずみエネルギー密度関数の係数を算出する用途に適用できる。
【0100】
さらに、本発明は上述した実施の形態のみに限定されるものではなく、既に述べた本発明の要旨を逸脱しない範囲において種々の変更が可能であることは勿論である。
【符号の説明】
【0101】
1 入力部
2 計算部
21 入力データ読み込み手段
22 第一算出手段(シフトファクタ算出手段)
23 第二算出手段(Prony級数からα、αを算出する手段)
24 第三算出手段(応力緩和分を除去した超弾性モデルの応力−ひずみ曲線のデータを逐次的に計算する手段)
25 第四算出手段(ひずみエネルギー密度関数の係数計算手段)
3 記憶部
31 第一記憶部(シフトファンクションの係数、試験温度記憶部)
32 第二記憶部(Prony級数と測定データ記憶部)
33 第三記憶部(超弾性モデルの応力と対応するひずみの記憶部)
4 出力部
5 ネットワーク
6 サーバコンピュータ
10 係数算出装置
11 手段
12 手段
100 第1の処理装置
101 入力部
102 計算部
121 係数算出部
104 出力部
200 第2の処理装置
201 入力部
202 計算部
204 出力部
300 第3の処理装置
301 入力部
302 計算部
304 出力部
402 計算部

【特許請求の範囲】
【請求項1】
大ひずみ粘弾性(Large Strain Viscoelasticity)モデルにおいて、超弾性モデルの構成式に使用されるひずみエネルギー密度関数の係数を材料試験の測定時間に対応する測定データから推定する係数算出装置であって、
応力緩和分が除去された超弾性モデルの応力−ひずみ曲線の計算データを逐次的に算出する逐次算出手段と、
前記応力−ひずみ曲線の計算データに基づいて超弾性モデルのひずみエネルギー密度関数の係数を算出する手段と、を備える係数算出装置。
【請求項2】
前記逐次算出手段は、所定の測定時間に対応する測定データ(A)及び前記所定の測定時間より短い測定時間に対応する超弾性モデルの応力−ひずみ曲線の計算データ(B)に基づいて、超弾性モデルの応力−ひずみ曲線の計算データを算出することを特徴とする請求項1に記載の係数算出装置。
【請求項3】
前記逐次算出手段は、前記測定データ(A)及び前記計算データ(B)に加え、前記所定の測定時間と前記所定の測定時間より短い測定時間との差分(C)に基づいて、超弾性モデルの応力−ひずみ曲線の計算データを算出することを特徴とする請求項2に記載の係数算出装置。
【請求項4】
前記逐次算出手段は、前記所定の測定時間より短い測定時間に対応する2つの計算データの差分から逐次的に算出される計算値を用いることを特徴とする請求項1乃至3のいずれか1項に記載の係数算出装置。
【請求項5】
シフトファクターを算出する手段、及びProny級数の係数を時刻0の緩和弾性率で除した係数を算出する手段の少なくとも一方の手段を有している請求項1乃至4のいずれか1項に記載の係数算出装置。
【請求項6】
大ひずみ粘弾性(Large Strain Viscoelasticity)モデルにおいて、超弾性モデルの構成式に使用されるひずみエネルギー密度関数の係数を材料試験の測定時間に対応する測定データから推定する係数算出方法であって、
応力緩和分が除去された超弾性モデルの応力−ひずみ曲線の計算データを逐次的に算出する逐次算出ステップと、
前記応力−ひずみ曲線の計算データに基づいて超弾性モデルのひずみエネルギー密度関数の係数を算出するステップと、を備える係数算出方法。
【請求項7】
前記逐次算出ステップは、所定の測定時間に対応する測定データ(A)及び前記所定の測定時間より短い測定時間に対応する超弾性モデルの応力−ひずみ曲線の計算データ(B)に基づいて、超弾性モデルの応力−ひずみ曲線の計算データを算出することを特徴とする請求項6に記載の係数算出方法。
【請求項8】
前記逐次算出ステップは、前記測定データ(A)及び前記計算データ(B)に加え、前記所定の測定時間と前記所定の測定時間より短い測定時間との差分(C)に基づいて、超弾性モデルの応力−ひずみ曲線の計算データを算出することを特徴とする請求項7に記載の係数算出方法。
【請求項9】
前記逐次算出ステップは、前記所定の測定時間より短い測定時間に対応する2つの計算データの差分から逐次的に算出される計算値を用いることを特徴とする請求項6乃至8のいずれか1項に記載の係数算出方法。
【請求項10】
シフトファクターを算出するステップ、及びProny級数の係数を時刻0の緩和弾性率で除した係数を算出するステップの少なくとも一方の手段を有している請求項6乃至9のいずれか1項に記載の係数算出方法。
【請求項11】
大ひずみ粘弾性(Large Strain Viscoelasticity)モデルにおいて、超弾性モデルの構成式に使用されるひずみエネルギー密度関数の係数を材料試験の測定時間に対応する測定データから推定する処理をコンピュータに実行させる係数算出プログラムであって、
コンピュータに対して、
応力緩和分が除去された超弾性モデルの応力−ひずみ曲線の計算データを逐次的に算出させる逐次算出ステップと、
前記応力−ひずみ曲線の計算データに基づいて超弾性モデルのひずみエネルギー密度関数の係数を算出させるステップと、を備える係数算出プログラム。
【請求項12】
前記逐次算出ステップは、所定の測定時間に対応する測定データ(A)及び前記所定の測定時間より短い測定時間に対応する超弾性モデルの応力−ひずみ曲線の計算データ(B)に基づいて、超弾性モデルの応力−ひずみ曲線の計算データを算出させることを特徴とする請求項6に記載の係数算出プログラム。
【請求項13】
前記逐次算出ステップは、前記測定データ(A)及び前記計算データ(B)に加え、前記所定の測定時間と前記所定の測定時間より短い測定時間との差分(C)に基づいて、超弾性モデルの応力−ひずみ曲線の計算データを算出させることを特徴とする請求項7に記載の係数算出プログラム。
【請求項14】
前記逐次算出ステップは、前記所定の測定時間より短い測定時間に対応する2つの計算データの差分から逐次的に算出される計算値を用いることを特徴とする請求項11乃至13のいずれか1項に記載の係数算出プログラム。
【請求項15】
コンピュータに対して、
シフトファクターを算出させるステップ、及びProny級数の係数を時刻0の緩和弾性率で除した係数を算出させるステップの少なくとも一方の手段を有している請求項11乃至14のいずれか1項に記載の係数算出プログラム。

【図1】
image rotate

【図2】
image rotate

【図3】
image rotate

【図4】
image rotate

【図5】
image rotate

【図6】
image rotate

【図7】
image rotate

【図8】
image rotate

【図9】
image rotate

【図10】
image rotate

【図11】
image rotate

【図12】
image rotate

【図13】
image rotate


【公開番号】特開2010−8396(P2010−8396A)
【公開日】平成22年1月14日(2010.1.14)
【国際特許分類】
【出願番号】特願2009−741(P2009−741)
【出願日】平成21年1月6日(2009.1.6)
【出願人】(000004237)日本電気株式会社 (19,353)
【Fターム(参考)】