JP5012094B2 - 粘弾性特性決定方法、シミュレーションプログラムおよびシミュレーション装置、および、粘弾性特性を備えたシミュレーションモデルを用いるシミュレーション方法 - Google Patents
粘弾性特性決定方法、シミュレーションプログラムおよびシミュレーション装置、および、粘弾性特性を備えたシミュレーションモデルを用いるシミュレーション方法 Download PDFInfo
- Publication number
- JP5012094B2 JP5012094B2 JP2007055901A JP2007055901A JP5012094B2 JP 5012094 B2 JP5012094 B2 JP 5012094B2 JP 2007055901 A JP2007055901 A JP 2007055901A JP 2007055901 A JP2007055901 A JP 2007055901A JP 5012094 B2 JP5012094 B2 JP 5012094B2
- Authority
- JP
- Japan
- Prior art keywords
- value
- elastic modulus
- frequency
- viscoelastic
- simulation
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Expired - Fee Related
Links
Classifications
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
- Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation
Landscapes
- Tires In General (AREA)
Description
本発明は、材料の粘弾性の周波数分布から、時間の指数関数の和で表わされた時間領域の粘弾性特性を決定する方法に関する。また、本発明は、この時間領域の粘弾性特性を入力として、時間に依存する物理特性を算出するシミュレーション方法に関する。
タイヤに用いられるゴム材料は、一般に粘弾性特性を有し、粘弾性材料として扱われている。ここで、粘弾性特性とは、材料に応力又は歪みを与えたとき、そのとき生じる歪み又は応力が時間とともに変化する特性をいう。
近年、この粘弾性特性を有するゴム材料で構成されるタイヤの転がり抵抗、および、ゴム材料中に含まれるコンパウンドがタイヤの特性に与える影響の評価が、シミュレーションにより行われている。これらのシミュレーションは、一般的なシミュレーションツールであるＦＥＭで行なわれている。このＦＥＭでのシミュレーションにおいて、粘弾性特性は、時間領域で定められた粘弾性特性として入力される必要がある。
近年、この粘弾性特性を有するゴム材料で構成されるタイヤの転がり抵抗、および、ゴム材料中に含まれるコンパウンドがタイヤの特性に与える影響の評価が、シミュレーションにより行われている。これらのシミュレーションは、一般的なシミュレーションツールであるＦＥＭで行なわれている。このＦＥＭでのシミュレーションにおいて、粘弾性特性は、時間領域で定められた粘弾性特性として入力される必要がある。
しかし、一般に、材料の粘弾性特性は、周波数分布の形で測定される。
このため、上述したシミュレーションを行うために、材料の粘弾性の周波数分布（マスターカーブ）から、シミュレーション入力用の時間領域の粘弾性特性が決定される必要がある。非特許文献１および非特許文献２には、材料の粘弾性の周波数分布から、入力用の時間領域の粘弾性特性を決定する方法（ツール）が、種々提案されている。
このため、上述したシミュレーションを行うために、材料の粘弾性の周波数分布（マスターカーブ）から、シミュレーション入力用の時間領域の粘弾性特性が決定される必要がある。非特許文献１および非特許文献２には、材料の粘弾性の周波数分布から、入力用の時間領域の粘弾性特性を決定する方法（ツール）が、種々提案されている。
非特許文献１に記載された方法では、以下の手順で入力用の時間領域の粘弾性特性が決定されている。第１に、入力用の時間領域の粘弾性特性はProny級数の形で表わされ、このProny級数を具体的に定める緩和弾性定数と緩和時定数の初期値が設定される。次に、このProny級数を周波数領域で規定する入力用の貯蓄弾性率および入力用の損失弾性率が求められる。また別途、材料の粘弾性の周波数分布から材料の貯蓄弾性率および材料の損失弾性率が取得される。
次に、設定した貯蓄弾性率の設定値と取得した材料の貯蓄弾性率の取得値との差、および、設定した損失弾性率の設定値と取得した材料の損失弾性率の取得値との差の、周波数毎の、二乗和が算出される。最後に、算出された二乗和が判定条件の閾値より小でなければ、上記の緩和弾性定数と緩和時定数が変更され、上記のProny級数への緩和弾性定数と緩和時定数の設定以降のフローが繰り返される。また、算出された二乗和が判定条件の閾値より小ならば、上記の緩和弾性定数と緩和時定数を用いて、このProny級数が決定される。すなわち、入力用の時間領域の粘弾性特性が決定される。同時に、入力用の貯蓄弾性率の値、入力用の損失弾性率の値、入力用の損失弾性率の値を入力用の貯蓄弾性率の値で除算して定まる入力用の損失正接の値も決定される。
しかしながら、この入力用の損失正接の値は、以下の点で、材料の損失弾性率の値を材料の貯蓄弾性率の値で除算して定まる材料の損失正接の値と異なっている。取得した材料の損失正接の値は、一般に周波数に対してピークを持つ分布となり、そのピーク以下の周波数領域では、周波数の増加に対し、その値が、増加し続ける。これに対し、決定された損失正接の値は、周波数の増加に対し、低周波数ではその値は増加するが、その後急激に低下し、さらに増加する（変動する）変化がみられる。
タイヤのころがり抵抗のシミュレーションモデルに、決定された損失正接の値を入力として与え、シミュレーションを行うと、このシミュレーションモデルにより算出される転がり抵抗は、上記損失正接の周波数の増加に対する変化と同様の結果が得られる。実際のころがり抵抗は、タイヤの転動速度の増加に対し、その値が単調増加することから、このシミュレーションモデルにより算出される転がり抵抗は、実際のころがり抵抗を再現しない。
非特許文献２に記載された方法では、二乗和の算出以外の部分は、非特許文献１に記載された方法と同様の方法で、取得された材料の粘弾性の周波数分布から、時間領域の粘弾性特性が決定されている。ここで、二乗和の算出には、設定した貯蓄弾性率の設定値と取得した材料の貯蓄弾性率の取得値との差を材料の貯蓄弾性率の取得値で除算した値、および、設定した損失弾性率の設定値と取得した材料の損失弾性率の取得値との差を材料の損失弾性率の取得値で除算した値が用いられている。
非特許文献２に記載された方法により決定される損失正接の値も、周波数の増加に対し、低周波数ではその値は増加するが、その後急激に低下する（変動する）変化がみられる。
タイヤのころがり抵抗のシミュレーションモデルに、非特許文献２に記載された方法により決定された損失正接の値を入力として与え、シミュレーションを行うと、このシミュレーションモデルにより算出される転がり抵抗は、上記損失正接の周波数の増加に対する変化と同様の結果が得られる。このシミュレーションモデルにより算出される転がり抵抗も、実際のころがり抵抗を再現しない。
Ａｎａｌｙｓｉｓ Ｕｓｅｒ’ｓ Ｍａｎｕａｌ Ｖｏｌｕｍｅ ＩＩＩ：Ｍａｔｅｒｉａｌｓ（ＡＢＡＱＵＳ Ｖｅｒｓｉｏｎ ６．６）
第２４回Ｍｅｃｈ Ｄ＆Ａセミナー資料（株式会社メカニカルデザイン ２００５年１１月）
そこで、本発明は、材料の粘弾性の周波数分布から、時間の指数関数の和で表された時間領域の粘弾性特性の決定方法であって、上述したような損失正接の値の変動を抑えた、正確な時間領域の粘弾性特性の決定方法を提供することを目的とする。また、本発明は、数値解析可能な粘弾性特性を備えたシミュレーションモデルを用いて行なうシミュレーション方法であって、時間に依存する物理特性を正確に算出するシミュレーション方法を提供することを目的とする。
そこで、本発明は、上記問題点を解決するために、材料の粘弾性の周波数分布から、時間の指数関数の和で表された時間領域の粘弾性特性を定めるパラメータ値を算出し、前記時間領域の粘弾性特性を決定する粘弾性特性決定方法であって、前記時間領域の粘弾性特性を定めるパラメータ値の初期値として所定の値を設定する第１のステップと、前記材料の粘弾性の周波数分布から、貯蓄弾性率の値および損失弾性率の値を取得する第２のステップと、前記時間領域の粘弾性特性を定めるパラメータ値により設定された貯蓄弾性率の設定値、前記パラメータ値により設定された損失弾性率の設定値、前記第２のステップで取得した貯蓄弾性率の取得値および前記第２のステップで取得した損失弾性率の取得値から、前記貯蓄弾性率の設定値と取得値との誤差と、前記損失弾性率の設定値と取得値との誤差と、前記貯蓄弾性率の設定値の周波数に対する勾配の値と前記貯蓄弾性率の取得値の周波数に対する勾配の値との誤差と、前記損失弾性率の設定値の周波数に対する勾配の値と前記損失弾性率の取得値の周波数に対する勾配の値との誤差とを周波数毎に求め、この周波数毎の誤差を用いて、前記材料の粘弾性の周波数分布と前記パラメータ値により設定された設定粘弾性の周波数分布との分布誤差を算出する第３のステップと、前記分布誤差が収束判定条件の閾値より小ならば、前記パラメータ値を用いて前記時間領域の粘弾性特性を決定し、前記分布誤差が前記収束判定条件の閾値より小でなければ、前記パラメータ値を修正して前記第３のステップに戻る第４のステップとを有することを特徴とする粘弾性特性決定方法を提供する。
前記第３のステップでは、前記貯蓄弾性率の設定値と取得値、前記損失弾性率の設定値と取得値、前記貯蓄弾性率の設定値の周波数に対する勾配の値と前記貯蓄弾性率の取得値の周波数に対する勾配の値および前記損失弾性率の設定値の周波数に対する勾配の値と前記損失弾性率の取得値の周波数に対する勾配の値は、前記貯蓄弾性率および前記損失弾性率を対数変換した値を用いて求められる値であり、前記周波数分布は、周波数を対数変換した対数周波数分布であることが好ましい。
また、前記損失弾性率の設定値を前記貯蓄弾性率の取得値で除算して定まる損失正接の設定値、および、前記損失弾性率の取得値を前記貯蓄弾性率の取得値で除算して定まる損失正接の取得値から、前記損失正接の設定値の周波数に対する勾配の値と前記損失正接の取得値の周波数に対する勾配の値との誤差を周波数毎に求め、この誤差を前記分布誤差に含めて用いることが好ましい。
また、前記第３のステップでは、前記損失正接の設定値の周波数に対する勾配の値と前記損失正接の取得値の周波数に対する勾配の値は、前記損失正接を対数変換しない真数値を用いて求められる値であり、前記貯蓄弾性率の設定値と取得値、前記損失弾性率の設定値と取得値、前記貯蓄弾性率の設定値の周波数に対する勾配の値と前記貯蓄弾性率の取得値の周波数に対する勾配の値および前記損失弾性率の設定値の周波数に対する勾配の値と前記損失弾性率の取得値の周波数に対する勾配の値は、前記貯蓄弾性率および前記損失弾性率を対数変換した値を用いて求められる値であり、前記周波数分布は、周波数を対数変換した対数周波数分布であることが好ましい。
また、前記損失正接の設定値および取得値は、前記損失正接を対数変換しない真数値であり、前記損失正接の設定値と取得値との誤差を周波数毎に求め、この誤差を前記分布誤差に含めて用いることが好ましい。
なお、前記分布誤差は、前記誤差のそれぞれに重みづけ係数を乗算したものであり、前記損失正接の設定値と取得値との誤差に乗算する重みづけ係数は、前記誤差に乗算される重みづけ係数の中で最も大きいことが好ましい。また、重み付け係数は、設定損失正接の設定値と取得値との誤差に乗算する重み付け係数は１〜１００であり、他の誤差に乗算する重み付け係数は０〜１未満であることが好ましい。分布誤差は、最小二乗法により算出されることが好ましい。粘弾性率は、せん断粘弾性率、引張粘弾性率、体積粘弾性率のいずれかであることが好ましい。時間領域の粘弾性特性は、Prony級数で表わされることが好ましく、Prony級数の展開項数は、取得された材料の粘弾性を示すデータ数の１／２以上であることが好ましい。材料の粘弾性の周波数分布のデータ数は、周波数が１桁変化する範囲において、３つ以上含まれることが好ましい。
さらに、本発明は、数値解析可能な粘弾性特性を備えたシミュレーションモデルを用いて行なうシミュレーション方法であって、前記シミュレーションモデルの材料特性の１つとして、上述した粘弾性特性決定方法のいずれかにより決定された、時間の指数関数の和で表された時間領域の粘弾性特性を入力するステップと、前記シミュレーションモデルを作成するステップと、前記シミュレーションモデルに対して時間に依存する変形解析を行い、前記シミュレーションモデルから時間に依存する物理特性を算出するステップとを有することを特徴とするシミュレーション方法を提供する。
前記シミュレーションモデルは、ゴム材料と、このゴム材料を補強するために前記ゴム材料中に配される、無機化合物充填材、金属材及び繊維材の少なくともいずれか１つの補強材とを再現した有限要素モデルであることが好ましい。前記シミュレーションモデルは、タイヤモデルであることが好ましい。
なお、前記シミュレーションモデルの材料特性の１つとして、上述のシミュレーション方法により算出された粘弾性特性を入力することが好ましい。前記時間に依存する物理特性は、応力緩和特性、クリープ特性、定ひずみ負荷特性、動的粘弾性特性、ひずみ応答特性、荷重応答特性のいずれかであることが好ましい。前記無機化合物充填材は、カーボンブラックあるいはシリカであることが好ましい。
さらに、本発明では、数値解析可能な粘弾性特性を備えたシミュレーションモデルを用いて、コンピュータにシミュレーションを行なわせるシミュレーションプログラムであって、前記シミュレーションモデルの材料特性の１つとして、上述した粘弾性特性決定方法のいずれかにより決定された、時間の指数関数の和で表された時間領域の粘弾性特性をコンピュータに入力させる手順と、前記シミュレーションモデルを、コンピュータの演算部に作成させる手順と、前記シミュレーションモデルに対して時間に依存する変形解析をコンピュータの演算部に実行させ、前記シミュレーションモデルから時間に依存する物理特性をコンピュータの演算部に算出させる手順とを有することを特徴とするシミュレーションプログラムを提供する。
さらに、本発明では、数値解析可能な粘弾性特性を備えたシミュレーションモデルを用いてシミュレーションを行なうシミュレーション装置であって、前記シミュレーションモデルの材料特性の１つとして、上述した粘弾性特性決定方法のいずれかにより決定された、時間の指数関数の和で表された時間領域の粘弾性特性を入力する手段と、前記シミュレーションモデルを作成する手段と、前記シミュレーションモデルに対して時間に依存する変形解析を行い、前記シミュレーションモデルから時間に依存する物理特性を算出する手段とを有することを特徴とするシミュレーション装置を提供する。
本発明の粘弾性特性決定方法は、貯蓄弾性率および損失弾性率の値に加えて、貯蓄弾性率の周波数に対する勾配および損失弾性率の周波数に対する勾配の値を用いて、時間領域の粘弾性特性を決定する。これにより、本発明の貯蓄弾性率、損失弾性率および損失正接の値は、初期値の貯蓄弾性率、損失弾性率および損失正接の値が、設定された周波数帯内で変動する場合であっても、材料の貯蓄弾性率、損失弾性率および損失正接の値に近付き、この変動の振幅を減少させる。このため、本発明の粘弾性特性決定方法は、損失正接の値の変動を抑えた、正確な時間領域の粘弾性特性を決定することができる。
また、本発明の粘弾性特性を備えたシミュレーションモデルを用いるシミュレーション方法は、上述の方法により決定された時間領域で定められた粘弾性特性を入力とする。これにより、本発明のシミュレーションモデルの粘弾性特性は変動の振幅が減少したものとなる。このため、本発明のシミュレーション方法は、時間に依存する物理特性を正確に算出することができる。
本発明について、以下詳細に説明する。
図１は、本発明のシミュレーション装置の実施態様である、時間に依存する物理特性が算出されるシミュレーション装置１０の構成を機能的に示したブロック図である。
シミュレーション装置１０は、ＣＰＵ（演算部）１２、メモリ１４、ＲＯＭ１６、Ｉ／Ｏボード１８を備えたコンピュータで構成される。また、シミュレーション装置１０は、メモリ１４あるいはＲＯＭ１６に記憶されたアプリケーションソフトウェアを呼び出して、パラメータ値設定部２０、材料粘弾性取得部２２、分布誤差算出部２４、粘弾性特性決定部２６およびシミュレーション演算部２８をサブルーチン生成し、それぞれのサブルーチンを実行させる装置である。
図１は、本発明のシミュレーション装置の実施態様である、時間に依存する物理特性が算出されるシミュレーション装置１０の構成を機能的に示したブロック図である。
シミュレーション装置１０は、ＣＰＵ（演算部）１２、メモリ１４、ＲＯＭ１６、Ｉ／Ｏボード１８を備えたコンピュータで構成される。また、シミュレーション装置１０は、メモリ１４あるいはＲＯＭ１６に記憶されたアプリケーションソフトウェアを呼び出して、パラメータ値設定部２０、材料粘弾性取得部２２、分布誤差算出部２４、粘弾性特性決定部２６およびシミュレーション演算部２８をサブルーチン生成し、それぞれのサブルーチンを実行させる装置である。
シミュレーション装置１０は、取得された材料の粘弾性の周波数分布から、時間の指数関数の和で表された時間領域の粘弾性特性を決定する装置であり、この決定した粘弾性特性を用いて、時間に依存する物理特性を算出するシミュレーションを行う装置である。
また、シミュレーション装置１０は、ディスプレイ３０、プリンタ３２およびマウス・キーボード３４と接続されている。また、シミュレーション装置１０では、シミュレーションに必要な条件等が、ディスプレイ３０に表示された入力画面からマウス・キーボード３４にて入力指示され、粘弾性特性の決定結果が、ディスプレイ３０やプリンタ３２に数値、グラフあるいは図によって表示される。また、シミュレーション装置１０は、Ｉ／Ｏボード１８が図示されないＣＡＤ／ＣＡＭ等のデータ供給システムおよび計測システムと接続され、構成される。
パラメータ値設定部２０は、時間の指数関数の和で表された時間領域の粘弾性特性を定めるパラメータ値である緩和弾性定数および緩和時定数の初期値を設定する部分である。
材料粘弾性率取得部２２は、メモリ１４にマスターカーブとして記憶された材料の粘弾性の周波数分布を呼び出すことによって、取得された材料の粘弾性の周波数分布から、貯蓄弾性率の取得値および損失弾性率の取得値を取得する部分である。
分布誤差算出部２４は、設定されたパラメータ値により、貯蓄弾性率の値および損失弾性率の値が設定される粘弾性特性設定部分と、材料粘弾性率取得部２２で取得された貯蓄弾性率の取得値、粘弾性特性設定部分で設定された貯蓄弾性率の設定値、材料粘弾性率取得部２２で取得された損失弾性率の取得値、および、粘弾性特性設定部分で設定された損失弾性率の設定値を用いて分布誤差が算出される部分とで構成されている。
粘弾性特性決定部２６は、分布誤差の値と収束判定条件の閾値との大小が比較される部分と、分布誤差の値が収束判定条件の閾値より小ならば、パラメータ値を用いて時間領域の粘弾性特性が決定される部分と、分布誤差の値が収束判定条件の閾値より小でないならば、パラメータ値が別の値に修正されて再び設定される部分とで構成されている。
シミュレーション演算部２８は、シミュレーションモデルの材料特性の１つとして時間領域の粘弾性特性が入力される部分と、シミュレーションモデルが作成される部分と、シミュレーションモデルに対して時間に依存する変形解析が行なわれ、シミュレーションモデルから時間に依存する物理特性が算出される部分とで構成されている。
このような構成のシミュレーション装置１０において、以下の実施形態に示すように、材料の粘弾性の周波数分布から、時間の指数関数の和で表された時間領域の粘弾性特性が決定される。また、数値解析可能な粘弾性特性を備えたシミュレーションモデルを用いて、時間に依存する物理特性を算出するシミュレーションが行われる。
なお、本実施形態では、時間領域の粘弾性特性は、せん断粘弾性特性である。また、シミュレーションモデルは、図３に示すようなポリマー相Ａ内に粒状フィラー相Ｂが不均質に分散配置されているミクロコンパウンドモデルであり、時間に依存する物理特性は、応力緩和特性である。これらに基づいて詳細を説明するが、本発明はこれには限定されるものではない。
なお、本実施形態では、時間領域の粘弾性特性は、せん断粘弾性特性である。また、シミュレーションモデルは、図３に示すようなポリマー相Ａ内に粒状フィラー相Ｂが不均質に分散配置されているミクロコンパウンドモデルであり、時間に依存する物理特性は、応力緩和特性である。これらに基づいて詳細を説明するが、本発明はこれには限定されるものではない。
図２は、本発明の実施態様である応力緩和特性のシミュレーションの一例を示すフローチャートである。以下、その流れを説明する。
まず、パラメータ値設定部２０にて、時間の指数関数の和で表された時間領域の粘弾性特性を定めるパラメータ値の初期値として、緩和弾性定数および緩和時定数に所定の値が設定される（ステップＳ１００）。
初期値は、デフォルト設定される値でも、オペレータによりマウス・キーボード３４から入力される値でも良い。
まず、パラメータ値設定部２０にて、時間の指数関数の和で表された時間領域の粘弾性特性を定めるパラメータ値の初期値として、緩和弾性定数および緩和時定数に所定の値が設定される（ステップＳ１００）。
初期値は、デフォルト設定される値でも、オペレータによりマウス・キーボード３４から入力される値でも良い。
次に、材料粘弾性率取得部２２にて、メモリ１４から呼び出された材料の粘弾性を示す特性（実数部、虚数部で表わされた周波数分布）が取得される（ステップＳ１０２）。
ここで、材料の粘弾性は、いずれの粘弾性でも良く、例えば、引張粘弾性または体積粘弾性であっても良いが、本実施形態のせん断粘弾性が好ましい。
取得される特性として、材料のせん断粘弾性の周波数分布（複素せん断弾性率）を下記式（１−１）で表わすことにより、下記式（１−２）に示すせん断貯蓄弾性率の値、および、下記式（１−３）に示すせん断損失弾性率の値が取得される。本実施形態では、取得されたせん断貯蓄弾性率の値は、貯蓄弾性率の取得値とされ、取得されたせん断損失弾性率の値は、損失弾性率の取得値とされる。また、損失弾性率の取得値を貯蓄弾性率の取得値で除算して定まる値は、損失正接の取得値とされる。
ここで、材料の粘弾性は、いずれの粘弾性でも良く、例えば、引張粘弾性または体積粘弾性であっても良いが、本実施形態のせん断粘弾性が好ましい。
取得される特性として、材料のせん断粘弾性の周波数分布（複素せん断弾性率）を下記式（１−１）で表わすことにより、下記式（１−２）に示すせん断貯蓄弾性率の値、および、下記式（１−３）に示すせん断損失弾性率の値が取得される。本実施形態では、取得されたせん断貯蓄弾性率の値は、貯蓄弾性率の取得値とされ、取得されたせん断損失弾性率の値は、損失弾性率の取得値とされる。また、損失弾性率の取得値を貯蓄弾性率の取得値で除算して定まる値は、損失正接の取得値とされる。
なお、材料の粘弾性の周波数分布は、計測装置からデータとして得られる。
さらになお、材料の粘弾性の周波数分布のデータ数は、周波数が１桁変化する範囲において、例えば周波数が１〜１０〔Ｈｚ〕の範囲において、３つ以上含まれることが好ましい。これにより、後述する、時間の指数関数の和で表された時間領域の粘弾性特性が精度良く決定できる。
さらになお、材料の粘弾性の周波数分布のデータ数は、周波数が１桁変化する範囲において、例えば周波数が１〜１０〔Ｈｚ〕の範囲において、３つ以上含まれることが好ましい。これにより、後述する、時間の指数関数の和で表された時間領域の粘弾性特性が精度良く決定できる。
次に、分布誤差算出部２４にて、時間の指数関数の和で表された時間領域の粘弾性特性が、パラメータ値を用いて設定される（ステップＳ１０４）。
ここで、時間領域の粘弾性特性の粘弾性は、いずれの粘弾性でも良く、例えば、引張粘弾性または体積粘弾性であっても良いが、本実施形態のせん断粘弾性が好ましい。
ここで、時間領域の粘弾性特性の粘弾性は、いずれの粘弾性でも良く、例えば、引張粘弾性または体積粘弾性であっても良いが、本実施形態のせん断粘弾性が好ましい。
本実施形態において、時間の指数関数の和で表された時間領域の粘弾性特性は、下記式（２）にて設定されるProny級数である。ここで緩和弾性定数および緩和時定数の値は、ステップ１００で初期値として設定された値、または、後述するステップ１１２で設定された修正した値である。
なお、Prony級数の展開項数は、前述したステップＳ１０２で取得された材料の粘弾性を示すデータ数の１／２以上とすることが好ましい。これにより、後述する、時間の指数関数の和で表された時間領域の粘弾性特性が精度良く決定できる。
なお、Prony級数の展開項数は、前述したステップＳ１０２で取得された材料の粘弾性を示すデータ数の１／２以上とすることが好ましい。これにより、後述する、時間の指数関数の和で表された時間領域の粘弾性特性が精度良く決定できる。
本実施形態において、時間領域の粘弾性特性を定める緩和弾性定数および緩和時定数の値により、下記式（３−１）に示すせん断貯蓄弾性率、および、下記式（３−２）に示すせん断損失弾性率が設定される。ここで、これらのせん断貯蓄弾性率およびせん断損失弾性率は、上記式（２）にて設定されるProny級数が、フーリエ変換されることにより設定されるものである。
本実施形態では、設定されたせん断貯蓄弾性率の値は、貯蓄弾性率の設定値とされ、設定されたせん断損失弾性率の値は、損失弾性率の設定値とされる。また、損失弾性率の設定値を貯蓄弾性率の設定値で除算して定まる値は、損失正接の設定値とされる。
本実施形態では、設定されたせん断貯蓄弾性率の値は、貯蓄弾性率の設定値とされ、設定されたせん断損失弾性率の値は、損失弾性率の設定値とされる。また、損失弾性率の設定値を貯蓄弾性率の設定値で除算して定まる値は、損失正接の設定値とされる。
次に、分布誤差算出部２４にて、誤差が求められ、分布誤差が算出される（ステップＳ１０６）。
ここで貯蓄弾性率の設定値と取得値との誤差は、この２つの値の差を周波数毎に定量化するものである。例えば、この２つの値の差であっても良いが、本実施形態の下記式（４）に示すように、この２つの値の差を蓄弾弾性率の取得値で除算した値であることが好ましい。また、この誤差は、この２つの値の差を蓄弾性率の設定値で除算した値であっても良い。
損失弾性率の設定値と取得値との誤差、損失正接の設定値と取得値との誤差、貯蓄弾性率の設定値の周波数に対する勾配の値と貯蓄弾性率の取得値の周波数に対する勾配の値との誤差、損失弾性率の設定値の周波数に対する勾配の値と損失弾性率の取得値の周波数に対する勾配の値との誤差、および、損失正接の設定値の周波数に対する勾配の値と損失正接の取得値の周波数に対する勾配の値との誤差も、上記貯蓄弾性率の設定値と取得値との誤差と同様の誤差であっても良い。
ここで貯蓄弾性率の設定値と取得値との誤差は、この２つの値の差を周波数毎に定量化するものである。例えば、この２つの値の差であっても良いが、本実施形態の下記式（４）に示すように、この２つの値の差を蓄弾弾性率の取得値で除算した値であることが好ましい。また、この誤差は、この２つの値の差を蓄弾性率の設定値で除算した値であっても良い。
損失弾性率の設定値と取得値との誤差、損失正接の設定値と取得値との誤差、貯蓄弾性率の設定値の周波数に対する勾配の値と貯蓄弾性率の取得値の周波数に対する勾配の値との誤差、損失弾性率の設定値の周波数に対する勾配の値と損失弾性率の取得値の周波数に対する勾配の値との誤差、および、損失正接の設定値の周波数に対する勾配の値と損失正接の取得値の周波数に対する勾配の値との誤差も、上記貯蓄弾性率の設定値と取得値との誤差と同様の誤差であっても良い。
また、分布誤差は、上記誤差を定量化し累積加算したいずれの分布誤差であっても良い。分布誤差は、例えば、各誤差の絶対値を取り累積加算した分布誤差でも良く、また各誤差の絶対値をｋ乗（ｋ＞０）し累積加算した分布誤差でも良いが、上記式（４）に示すように、各誤差を２乗し累積加算した分布誤差であることが好ましい。
式（４）に示す分布誤差は、貯蓄弾性率の項および損失弾性率の項に加え、損失正接の項、貯蓄弾性率の周波数に対する勾配の項、損失弾性率の周波数に対する勾配の項および損失正接の周波数に対する勾配の項で構成されている。この構成により、詳細は後述するように、正確な時間領域の粘弾性特性を決定することができる。
ここで、式（４）に用いられる、貯蓄弾性率の設定値と取得値、損失弾性率の設定値と取得値、貯蓄弾性率の設定値の周波数に対する勾配の値と貯蓄弾性率の取得値の周波数に対する勾配の値および損失弾性率の設定値の周波数に対する勾配の値と損失弾性率の取得値の周波数に対する勾配の値は、貯蓄弾性率および損失弾性率を対数変換した値を用いて求められる値であっても良い。
損失正接の設定値と取得値および損失正接の設定値の周波数に対する勾配の値と損失正接の取得値の周波数に対する勾配の値は、損失正接を対数変換しない真数値を用いて求められる値であっても良い。
なお、上述した周波数に対する勾配の値は、対数変換した周波数の値を用いて求められる値であっても良く、また、対数変換しない周波数の真数値を用いて求められる値であっても良い。
損失正接の設定値と取得値および損失正接の設定値の周波数に対する勾配の値と損失正接の取得値の周波数に対する勾配の値は、損失正接を対数変換しない真数値を用いて求められる値であっても良い。
なお、上述した周波数に対する勾配の値は、対数変換した周波数の値を用いて求められる値であっても良く、また、対数変換しない周波数の真数値を用いて求められる値であっても良い。
式（４）の累積加算は、周波数を対数変換した対数周波数に対して行っても良い。
なお、本実施形態では、各誤差を２乗し累積加算するために、Newton-Raphson法を用いたが、他の公知の非線形の最小二乗法の解法を用いても良い。
なお、本実施形態では、各誤差を２乗し累積加算するために、Newton-Raphson法を用いたが、他の公知の非線形の最小二乗法の解法を用いても良い。
本実施形態では、分布誤差（χ２）を求めるために、重み付け係数の値を、ｗ０＝３．４８、ｗ１＝１．２７、ｗ２＝５．８、ｗ３＝１．３６、ｗ４＝１．２５、ｗ５＝０．１４としたが、これらの値は予め設定されたものである。ｗ０、ｗ１、ｗ３〜ｗ５の値は、予め設定された０〜１未満の値であり、ｗ２は予め設定された１〜１００の値とすることが、より正確な分布誤差（χ２）が求められる点で好ましい。
また、本実施形態では、分布誤差（χ２）を求めるために、式（４）を用いたが、下記式（５）および下記式（６）を用いても良い。
ここで、下記式（５）は、式（４）の右辺に対して、損失正接の項（第３項）を削除した式である。また、下記式（６）は、式（４）の右辺に対して、損失正接の項（第３項）および損失正接の周波数に対する勾配の項（第６項）を削除した式である。
ここで、下記式（５）は、式（４）の右辺に対して、損失正接の項（第３項）を削除した式である。また、下記式（６）は、式（４）の右辺に対して、損失正接の項（第３項）および損失正接の周波数に対する勾配の項（第６項）を削除した式である。
次に、粘弾性特性決定部２６にて、分布誤差の値と収束判定条件の閾値との大小が比較される（ステップＳ１０８）。分布誤差の値が収束判定条件の閾値より小ならば、フローは、粘弾性特性決定のステップＳ１１０に進められる。分布誤差の値が収束判定条件の閾値より小でなければ、フローは、パラメータ値修正のステップＳ１１２に進められる。
ここで、収束判定条件の閾値は、後述するステップＳ１１０にて、時間の指数関数の和で表された時間領域の粘弾性特性を必要な精度で決定するいずれの閾値でも良い。精度を重視する場合には、閾値は小さく設定され、算出コストの節約または算出時間の短縮を重視する場合には、閾値は大きく設定される。
ここで、収束判定条件の閾値は、後述するステップＳ１１０にて、時間の指数関数の和で表された時間領域の粘弾性特性を必要な精度で決定するいずれの閾値でも良い。精度を重視する場合には、閾値は小さく設定され、算出コストの節約または算出時間の短縮を重視する場合には、閾値は大きく設定される。
分布誤差の値が収束判定条件の閾値より小のときは、設定された応力緩和定数の値および応力時定数の値が最終的なパラメータ値として決定され、粘弾性特性決定部２６にて、時間の指数関数の和で表された時間領域の粘弾性特性が決定される（ステップＳ１１０）。
ここで、時間領域の粘弾性特性の粘弾性は、いずれの粘弾性でも良く、例えば、引張粘弾性または体積粘弾性であっても良いが、本実施形態のせん断粘弾性が好ましい。
本実施形態では、時間領域の粘弾性特性は、式（２）で設定されるProny級数に、ステップＳ１１０で決定された緩和弾性定数の値および緩和時定数の値が代入されて、決定される。
ここで、時間領域の粘弾性特性の粘弾性は、いずれの粘弾性でも良く、例えば、引張粘弾性または体積粘弾性であっても良いが、本実施形態のせん断粘弾性が好ましい。
本実施形態では、時間領域の粘弾性特性は、式（２）で設定されるProny級数に、ステップＳ１１０で決定された緩和弾性定数の値および緩和時定数の値が代入されて、決定される。
分布誤差の値が収束判定条件の閾値より小でないときは、パラメータ値が別の値に修正されて再び設定される（ステップＳ１１２）。その後、フローは、粘弾性特性設定のステップＳ１０４に進められる。
ここで、修正されるパラメータ値は、分布誤差の値と収束判定条件の閾値との差を用いて与えられる値であり、公知の方法を用いて値が修正される。
ここで、修正されるパラメータ値は、分布誤差の値と収束判定条件の閾値との差を用いて与えられる値であり、公知の方法を用いて値が修正される。
次に、シミュレーション演算部２８にて、シミュレーションモデルの材料特性の１つとして時間領域の粘弾性特性が入力される（ステップＳ１１４）。
ここで、時間領域の粘弾性特性の粘弾性は、本発明により決定されるいずれの粘弾性でも良く、例えば、引張粘弾性または体積粘弾性であっても良いが、本実施形態では、せん断粘弾性が好ましい。
ここで、時間領域の粘弾性特性の粘弾性は、本発明により決定されるいずれの粘弾性でも良く、例えば、引張粘弾性または体積粘弾性であっても良いが、本実施形態では、せん断粘弾性が好ましい。
次に、シミュレーション演算部２８にて、シミュレーションモデルが作成される（ステップＳ１１６）。
ここで、シミュレーションモデルは、数値解析可能な粘弾性特性を備えたいずれのシミュレーションモデルでも良く、例えば、均一な粘弾性特性を備えた直方体のモデルでもよく、金属材または繊維材が補強材として配されたゴム材料のモデルでも良く、タイヤの形状を模し、各部分毎に異なる粘弾性特性を備えたモデルでも良いが、本実施形態では、図３に示すミクロコンパウンドモデルであることが好ましい。
図３に示すミクロコンパウンドモデルは、ポリマー相Ａ（図３中灰色領域）と、カーボンブラックやシリカ等の無機化合物充填補強材である粒系フィラー相Ｂ（図３中の黒色領域）とからなり、ポリマー相Ａ中にフィラー相Ｂが不均質に分散配置されているモデルである。灰色領域のポリマー相Ａの材料特性は、時間に応じて変化する粘弾性を有し、この粘弾性は、ステップＳ１１０にて決定された時間領域の粘弾性特性により規定されている。一方、フィラー相Ｂの材料特性は、弾性特性で規定されている。なお、図３では、３０個×３０個×３０個にメッシュ分割した例が示されている。
ここで、シミュレーションモデルは、数値解析可能な粘弾性特性を備えたいずれのシミュレーションモデルでも良く、例えば、均一な粘弾性特性を備えた直方体のモデルでもよく、金属材または繊維材が補強材として配されたゴム材料のモデルでも良く、タイヤの形状を模し、各部分毎に異なる粘弾性特性を備えたモデルでも良いが、本実施形態では、図３に示すミクロコンパウンドモデルであることが好ましい。
図３に示すミクロコンパウンドモデルは、ポリマー相Ａ（図３中灰色領域）と、カーボンブラックやシリカ等の無機化合物充填補強材である粒系フィラー相Ｂ（図３中の黒色領域）とからなり、ポリマー相Ａ中にフィラー相Ｂが不均質に分散配置されているモデルである。灰色領域のポリマー相Ａの材料特性は、時間に応じて変化する粘弾性を有し、この粘弾性は、ステップＳ１１０にて決定された時間領域の粘弾性特性により規定されている。一方、フィラー相Ｂの材料特性は、弾性特性で規定されている。なお、図３では、３０個×３０個×３０個にメッシュ分割した例が示されている。
次に、シミュレーション演算部２８にて、上記シミュレーションモデルに対し、境界条件が与えられ、時間に依存する変形解析が行なわれ、上記シミュレーションモデルから、時間に依存する物理特性がシミュレーションにて算出される（ステップＳ１１８）。
ここで、シミュレーションモデルに与えられる境界条件は、時間に依存する物理特性の算出に用いるいずれの境界条件でも良く、例えば、力、応力、変位、熱のいずれでも良いが、本実施例では歪であることが好ましい。また、シミュレーションで算出される物理特性は、時間に依存するいずれの物理特性でも良く、例えば、クリープ特性、定ひずみ負荷特性、動的粘弾性特性のいずれであっても良いが、本実施形態では、応力を物理量として算出した応力緩和特性とすることが好ましい。
ここで、シミュレーションモデルに与えられる境界条件は、時間に依存する物理特性の算出に用いるいずれの境界条件でも良く、例えば、力、応力、変位、熱のいずれでも良いが、本実施例では歪であることが好ましい。また、シミュレーションで算出される物理特性は、時間に依存するいずれの物理特性でも良く、例えば、クリープ特性、定ひずみ負荷特性、動的粘弾性特性のいずれであっても良いが、本実施形態では、応力を物理量として算出した応力緩和特性とすることが好ましい。
本実施形態では、図３に示すミクロコンパウンドモデルに、時刻０において一定値の歪が与えられ、時刻０以降保持される。時刻０以降に発生する応力を算出することにより、応力緩和特性がシミュレーションにより算出され、このフローは終了する。
図４は、本実施形態のシミュレーションにより算出される応力緩和特性（以下「本発明の応力緩和特性」）の一例の結果を示す図である。また、比較例として、非特許文献１（以下「ＡＢＡＱＵＳ」）に記載の方法により算出されるシミュレーション結果である応力緩和特性（以下「ＡＢＡＱＵＳの応力緩和特性」）も示す。
図４に示すように、本発明の応力緩和特性を示す曲線では、時間の経過に対し、応力はなだらかに減少し、応力の緩和量（減少量）は減少している。これに対し、ＡＢＡＱＵＳの応力特性を示す曲線では、図４中のＡにて示す時間帯において、応力は急激に減少し、応力の緩和量（減少量）は増加している。
図４に示すように、本発明の応力緩和特性を示す曲線では、時間の経過に対し、応力はなだらかに減少し、応力の緩和量（減少量）は減少している。これに対し、ＡＢＡＱＵＳの応力特性を示す曲線では、図４中のＡにて示す時間帯において、応力は急激に減少し、応力の緩和量（減少量）は増加している。
ここで、応力緩和特性は、時刻０において一定値の歪を与えて時刻０以降保持した時に、時刻０以降に発生する応力である。一般に、粘弾性特性を備えた材料の応力緩和特性を示す曲線は、時間の経過に対し、応力はなだらかに減少し、また、応力の緩和量（減少量）は減少する、という特徴を備える。
本発明の応力緩和特性を示す曲線は、一般的な応力緩和特性を示す曲線が備える特徴を良く示している。これに対し、ＡＢＡＱＵＳの応力特性を示す曲線は、上述したように図４中のＡにて示す時間帯において、一般的な応力緩和特性が備える特徴と異なる特徴を示している。
本発明の応力緩和特性を示す曲線は、一般的な応力緩和特性を示す曲線が備える特徴を良く示している。これに対し、ＡＢＡＱＵＳの応力特性を示す曲線は、上述したように図４中のＡにて示す時間帯において、一般的な応力緩和特性が備える特徴と異なる特徴を示している。
本発明の応力緩和特性と、ＡＢＡＱＵＳの応力緩和特性との差異について図５（ａ）〜（ｃ）および図６（ａ）〜（ｃ）を用いて説明する。図５（ａ）〜（ｃ）は、上述した式（４）を用いて決定された時間領域の粘弾性特性から定まる貯蓄弾性率、損失弾性率および損失正接（以下「本発明の貯蓄弾性率、損失弾性率および損失正接」）の値を示し、図６（ａ）〜（ｃ）は、ＡＢＡＱＵＳに記載の下記式（７）を用いて決定された時間領域の粘弾性特性から定まる貯蓄弾性率、損失弾性率および損失正接（以下「ＡＢＡＱＵＳの貯蓄弾性率、損失弾性率および損失正接」）の値を示す。
また、図５（ａ）〜（ｃ）および図６（ａ）〜（ｃ）は、材料から取得された貯蓄弾性率、損失弾性率および損失正接（以下「材料の貯蓄弾性率、損失弾性率および損失正接」）の値も同時に示す。
また、図５（ａ）〜（ｃ）および図６（ａ）〜（ｃ）は、材料から取得された貯蓄弾性率、損失弾性率および損失正接（以下「材料の貯蓄弾性率、損失弾性率および損失正接」）の値も同時に示す。
図５（ａ）〜（ｃ）に示されるように、本発明の貯蓄弾性率、損失弾性率および損失正接の値は、材料の貯蓄弾性率、損失弾性率および損失正接の値と良く一致している。特に、粘弾性特性に大きな影響を与える損失正接の値は、周波数１０のマイナス１乗〜１０の４乗〔Ｈｚ〕の範囲で、材料の損失正接の値とほぼ一致している。
これに対して、図６（ａ）に示すＡＢＡＱＵＳの貯蓄弾性率の値は材料の貯蓄弾性率の値と比較的良く一致しているものの、図６（ｂ）、（ｃ）に示すＡＢＡＱＵＳの損失弾性率および損失正接の値は、材料の損失弾性率および損失正接の値と一致しない部分がある。特に、粘弾性特性に大きな影響を与える損失正接の値は、図６（ｃ）に示すＢおよびＣの周波数帯において、材料の損失正接の値より大きく、図６（ｃ）に示すＤの周波数帯において、材料の損失正接の値より小さい。
損失正接の値の変化は、応力緩和特性に大きな影響を与えるので、図４中のＡの時間帯におけるＡＢＡＱＵＳの応力緩和特性の値の急激な減少は、図６（ｃ）中のＣの周波数帯において、ＡＢＡＱＵＳの損失正接の値が、材料の損失正接の値より大きいことによって引き起こされている。
これに対して、図６（ａ）に示すＡＢＡＱＵＳの貯蓄弾性率の値は材料の貯蓄弾性率の値と比較的良く一致しているものの、図６（ｂ）、（ｃ）に示すＡＢＡＱＵＳの損失弾性率および損失正接の値は、材料の損失弾性率および損失正接の値と一致しない部分がある。特に、粘弾性特性に大きな影響を与える損失正接の値は、図６（ｃ）に示すＢおよびＣの周波数帯において、材料の損失正接の値より大きく、図６（ｃ）に示すＤの周波数帯において、材料の損失正接の値より小さい。
損失正接の値の変化は、応力緩和特性に大きな影響を与えるので、図４中のＡの時間帯におけるＡＢＡＱＵＳの応力緩和特性の値の急激な減少は、図６（ｃ）中のＣの周波数帯において、ＡＢＡＱＵＳの損失正接の値が、材料の損失正接の値より大きいことによって引き起こされている。
本発明の貯蓄弾性率、損失弾性率および損失正接の値と、ＡＢＡＱＵＳの貯蓄弾性率、損失弾性率および損失正接の値との差異について、貯蓄弾性率を例にとり説明する。
ここで、本発明の貯蓄弾性率の値は、上述した式（４）を用いて決定された時間領域の粘弾性特性（以下「本発明の粘弾性特性」）により定められている。また、ＡＢＡＱＵＳの貯蓄弾性率の値は、上述した式（７）を用いて決定された時間領域の粘弾性特性（以下「ＡＢＡＱＵＳの粘弾性特性」）により定められている。
なお、、ＡＢＡＱＵＳの粘弾性特性は式（７）に示す分布誤差（χ２）が判定条件の閾値より小となったときのパラメータ値により決定された粘弾性特性である。
ここで、本発明の貯蓄弾性率の値は、上述した式（４）を用いて決定された時間領域の粘弾性特性（以下「本発明の粘弾性特性」）により定められている。また、ＡＢＡＱＵＳの貯蓄弾性率の値は、上述した式（７）を用いて決定された時間領域の粘弾性特性（以下「ＡＢＡＱＵＳの粘弾性特性」）により定められている。
なお、、ＡＢＡＱＵＳの粘弾性特性は式（７）に示す分布誤差（χ２）が判定条件の閾値より小となったときのパラメータ値により決定された粘弾性特性である。
はじめに、ＡＢＡＱＵＳの貯蓄弾性率の値の決定方法を、図７（ａ）および（ｂ）を用いて説明する。図７（ａ）に、初期値のパラメータ値により設定された貯蓄弾性率（以下「初期値の貯蓄弾性率」）の値が示され、図７（ｂ）に、ＡＢＡＱＵＳの貯蓄弾性率の値が示される。
ＡＢＡＱＵＳの粘弾性特性の決定に用いられる分布誤差を求めるための式（７）は、貯蓄弾性率および損失弾性率の項で構成されている。このため、周波数毎の貯蓄弾性率の設定値と取得値の差を二乗した値と、周波数毎の損失弾性率の設定値と取得値との差を二乗した値の累積和を小さくするように、パラメータ値は修正され、ＡＢＡＱＵＳの粘弾性特性は決められる。
このため、初期値の貯蓄弾性率の値が、図７（ａ）に示すように、設定された周波数帯内で変動する場合、ＡＢＡＱＵＳの貯蓄弾性率の値は材料の貯蓄弾性率の値に比較的近付くが、図７（ｂ）に示すようにこの変動の振幅を減少させない。
ＡＢＡＱＵＳの粘弾性特性の決定に用いられる分布誤差を求めるための式（７）は、貯蓄弾性率および損失弾性率の項で構成されている。このため、周波数毎の貯蓄弾性率の設定値と取得値の差を二乗した値と、周波数毎の損失弾性率の設定値と取得値との差を二乗した値の累積和を小さくするように、パラメータ値は修正され、ＡＢＡＱＵＳの粘弾性特性は決められる。
このため、初期値の貯蓄弾性率の値が、図７（ａ）に示すように、設定された周波数帯内で変動する場合、ＡＢＡＱＵＳの貯蓄弾性率の値は材料の貯蓄弾性率の値に比較的近付くが、図７（ｂ）に示すようにこの変動の振幅を減少させない。
次に、本発明の貯蓄弾性率の値の決定方法を、図７（ａ）、（ｃ）および図８（ａ）、（ｂ）を用いて説明する。図７（ａ）に、初期値の貯蓄弾性率の値が示され、図７（ｃ）に、本発明の貯蓄弾性率の値が示される。また図８（ａ）に、初期値のパラメータ値により設定された貯蓄弾性率の周波数に対する勾配（以下「初期値の貯蓄弾性率の勾配」）の値が示され、図８（ｂ）に、本発明の貯蓄弾性率の周波数に対する勾配（以下「本発明の貯蓄弾性率の勾配」）の値が示される。
本発明の粘弾性特性の決定に用いられる分布誤差を求めるための上記式（４）は、貯蓄弾性率および損失弾性率の項に加え、損失正接、貯蓄弾性率の周波数に対する勾配、損失弾性率の周波数に対する勾配および損失正接の周波数に対する勾配の項で構成されている。
貯蓄弾性率の項は、周波数毎の貯蓄弾性率の設定値と取得値の差を貯蓄弾性率の取得値で除算した値の項であり、損失弾性率の項は、周波数毎の損失弾性率の設定値と取得値の差を損失弾性率の取得値で除算した値の項である。また、貯蓄弾性率の周波数に対する勾配の項は、周波数毎の貯蓄弾性率の設定値の周波数に対する勾配の値と貯蓄弾性率の取得値の周波数に対する勾配の値との差を、貯蓄弾性率の取得値の周波数に対する勾配の値で除算した値の項である。また、損失弾性率の周波数に対する勾配の項は、周波数毎の損失弾性率の設定値の周波数に対する勾配の値と損失弾性率の取得値の周波数に対する勾配の値との差を、周波数毎の損失弾性率の取得値の周波数に対する勾配の値で除算した値の項である。
貯蓄弾性率の項は、周波数毎の貯蓄弾性率の設定値と取得値の差を貯蓄弾性率の取得値で除算した値の項であり、損失弾性率の項は、周波数毎の損失弾性率の設定値と取得値の差を損失弾性率の取得値で除算した値の項である。また、貯蓄弾性率の周波数に対する勾配の項は、周波数毎の貯蓄弾性率の設定値の周波数に対する勾配の値と貯蓄弾性率の取得値の周波数に対する勾配の値との差を、貯蓄弾性率の取得値の周波数に対する勾配の値で除算した値の項である。また、損失弾性率の周波数に対する勾配の項は、周波数毎の損失弾性率の設定値の周波数に対する勾配の値と損失弾性率の取得値の周波数に対する勾配の値との差を、周波数毎の損失弾性率の取得値の周波数に対する勾配の値で除算した値の項である。
貯蓄弾性率の項の値を二乗した値と、損失弾性率の項の値を二乗した値と、貯蓄弾性率の周波数に対する勾配の項の値を二乗した値と、損失弾性率の周波数に対する勾配の項の値を二乗した値との累積和を小さくするように、パラメータ値は修正され、本発明の粘弾性特性は決められる。
ここで、貯蓄弾性率の項の値を二乗した値と、損失弾性率の項の値を二乗した値との累積和を小さくするように、パラメータ値が修正されることにより、初期値の貯蓄弾性率の値は材料の貯蓄弾性率の値に近付く。
また、貯蓄弾性率の周波数に対する勾配の項の値を二乗した値と、損失弾性率の周波数に対する勾配の項の値を二乗した値との累積和を小さくするように、パラメータ値が修正されることにより、図８（ａ）に示す初期値の貯蓄弾性率の勾配の値は、材料の貯蓄弾性率の勾配の値に近付き、図８（ｂ）に示す本発明の貯蓄弾性率の勾配となる。
したがって、貯蓄弾性率の項の値を二乗した値と、損失弾性率の項の値を二乗した値に加え、貯蓄弾性率の周波数に対する勾配の項の値を二乗した値と、損失弾性率の周波数に対する勾配の項の値を二乗した値との累積和を小さくするように、パラメータ値が修正されることにより、初期値の貯蓄弾性率の値が、図７（ａ）に示すように、設定された周波数帯内で変動する場合であっても、本発明の貯蓄弾性率の値は図７（ｃ）に示すように、材料の貯蓄弾性率の値に近付き、この変動の振幅を減少させることができる。
さらに、式（４）は、損失正接の周波数に対する勾配の項も有している。この損失正接の周波数に対する勾配の項の値を二乗した値を、さらに加えた累積和を小さくするように、パラメータ値が修正されることにより、初期値の貯蓄弾性率の値は、さらに材料の貯蓄弾性率の値に近付くことができる。
さらに、式（４）は、損失正接の項も有している。この損失正接の周波数に対する勾配の項の値を二乗した値を、さらに加えた累積和を小さくするように、パラメータ値が修正されることにより、初期値の貯蓄弾性率の値は、さらに材料の貯蓄弾性率の値に近付くことができる。
さらに、式（４）は、重み付け係数（ｗ０〜ｗ５）も有している。各項の値を二乗した値にこの重み付け係数が乗算される。この乗算された値の累積和を小さくするように、パラメータ値が修正されることにより、初期値の貯蓄弾性率の値は、さらに材料の貯蓄弾性率の値に近付くことができる。
なお、本発明の貯蓄弾性率の値を例にとり、詳細に説明したが、本発明の損失弾性率および損失正接の値も、上述した構成により、本発明の貯蓄弾性率の値と同様に、材料の損失弾性率および損失正接の値に近付くことができる。
なお、本発明の貯蓄弾性率の値を例にとり、詳細に説明したが、本発明の損失弾性率および損失正接の値も、上述した構成により、本発明の貯蓄弾性率の値と同様に、材料の損失弾性率および損失正接の値に近付くことができる。
本発明は、第１に貯蓄弾性率および損失弾性率の値に加え、損失正接、貯蓄弾性率の周波数に対する勾配、損失弾性率の周波数に対する勾配および損失正接の周波数に対する勾配の値を用いて、時間領域の粘弾性特性を決定する。次に、本発明は、この時間領域の粘弾性特性を入力として、ミクロコンパウンドモデルに対し応力緩和特性を算出するシミュレーションを行う。
この構成により、上述したように、本発明の貯蓄弾性率、損失弾性率および損失正接の値は、初期値の貯蓄弾性率、損失弾性率および損失正接の値が、設定された周波数帯内で変動する場合であっても、材料の貯蓄弾性率、損失弾性率および損失正接の値に近付き、この変動の振幅を減少させる。
この作用により、本発明は、一般的な特性を説明する応力緩和特性をシミュレーションにより算出でき、また、この応力緩和特性の算出を可能とする時間領域の粘弾性特性の決定ができる。
この構成により、上述したように、本発明の貯蓄弾性率、損失弾性率および損失正接の値は、初期値の貯蓄弾性率、損失弾性率および損失正接の値が、設定された周波数帯内で変動する場合であっても、材料の貯蓄弾性率、損失弾性率および損失正接の値に近付き、この変動の振幅を減少させる。
この作用により、本発明は、一般的な特性を説明する応力緩和特性をシミュレーションにより算出でき、また、この応力緩和特性の算出を可能とする時間領域の粘弾性特性の決定ができる。
また、本発明は、タイヤの転がり抵抗を算出するシミュレーションに用いられても良い。
ステップＳ１１６にて、図９に示すタイヤモデル４０を作成し、ステップＳ１１８にて転動処理を施し、タイヤの転がり抵抗を算出する。
ステップＳ１１６にて、図９に示すタイヤモデル４０を作成し、ステップＳ１１８にて転動処理を施し、タイヤの転がり抵抗を算出する。
タイヤモデル４０は、時間領域の粘弾性特性を備えた複数の有限要素によって構成された３次元モデルである。
ここで、タイヤモデル４０のいずれの部分が、粘弾性特性を備えた複数の有限要素によって構成されていても良く、例えば、ビードフィラー部分あるいはタイヤ全体がこの有限要素によって構成されていても良いが、本シミュレーションでは、トレッド部分がこの有限要素によって構成されていことが好ましい。また、粘弾性特性は、本発明で決定されるいずれの粘弾性特性であっても良く、例えば、上述の実施形態のミクロコンパウンドモデルから算出される応力緩和特性としても良いが、本シミュレーションでは、本発明の粘弾性特性決定のステップＳ１１０で決定された時間領域の粘弾性特性を用いることもできる。
なお、トレッド部分は、直径方向６個×幅方向３６個×円周方向５６個にメッシュ分割されている。
ここで、タイヤモデル４０のいずれの部分が、粘弾性特性を備えた複数の有限要素によって構成されていても良く、例えば、ビードフィラー部分あるいはタイヤ全体がこの有限要素によって構成されていても良いが、本シミュレーションでは、トレッド部分がこの有限要素によって構成されていことが好ましい。また、粘弾性特性は、本発明で決定されるいずれの粘弾性特性であっても良く、例えば、上述の実施形態のミクロコンパウンドモデルから算出される応力緩和特性としても良いが、本シミュレーションでは、本発明の粘弾性特性決定のステップＳ１１０で決定された時間領域の粘弾性特性を用いることもできる。
なお、トレッド部分は、直径方向６個×幅方向３６個×円周方向５６個にメッシュ分割されている。
ステップＳ１１８で行なう転動処理は、例えば、別途作成された図示されないリムモデルにタイヤモデル４０が結合されて、内圧充填処理が施され、さらに別途作成された剛体路面モデルに対して接地処理が施され、所定の加重が負荷される。さらにタイヤモデル４０に転動処理が施される。
本シミュレーションでは、タイヤサイズ：２０５／６５Ｒ１５、タイヤ内圧：２００ｋＰａ、接地荷重：４ｋＮ、転動速度：１０ｋｍ／ｈ、転動周期：０．１１秒が境界条件として与えられる。
本シミュレーションでは、タイヤサイズ：２０５／６５Ｒ１５、タイヤ内圧：２００ｋＰａ、接地荷重：４ｋＮ、転動速度：１０ｋｍ／ｈ、転動周期：０．１１秒が境界条件として与えられる。
本シミュレーションにより算出された転がり抵抗（タイヤ進行方向軸反力）は、３１．８Ｎである。この３１．８Ｎは、タイヤモデルと剛体路面モデルとの接触摩擦抵抗より算出される転がり抵抗の３倍である。本発明は、粘弾性が考慮された（時間に依存した）転がり抵抗を得ることができる。
さらに、本実施形態では、せん断弾性率を用いる時間領域の粘弾性特性の決定を例として説明したが、本発明は、引張弾性率または体積弾性率を用いる時間領域の粘弾性特性の決定に適用されても良いのは勿論である。
また、本実施形態では、応力緩和特性の算出を例として説明したが、本発明は、クリープ特性、定ひずみ負荷特性、あるいは動的粘弾性特性等の時間に依存する物理特性を算出するシミュレーションに用いられても良いのは勿論である。
また、本実施形態では、応力緩和特性の算出を例として説明したが、本発明は、クリープ特性、定ひずみ負荷特性、あるいは動的粘弾性特性等の時間に依存する物理特性を算出するシミュレーションに用いられても良いのは勿論である。
なお、上述した方法がコンピュータ上で行なわれるときは、数値解析可能な粘弾性特性を備えたシミュレーションモデルを用いてシミュレーションを行なわせる以下の手順を備えたシミュレーションプログラムが用いられる。
（１）第１の手順において、本発明により決定された時間の指数関数の和で表された時間領域の粘弾性特性をコンピュータに入力させる。
（２）第２の手順において、シミュレーションモデルを、コンピュータの演算部に作成させる。
（３）第３の手順において、シミュレーションモデルに対して時間に依存する変形解析をコンピュータの演算部に実行させ、前記シミュレーションモデルから時間に依存する物理特性をコンピュータの演算部に算出させる。
（１）第１の手順において、本発明により決定された時間の指数関数の和で表された時間領域の粘弾性特性をコンピュータに入力させる。
（２）第２の手順において、シミュレーションモデルを、コンピュータの演算部に作成させる。
（３）第３の手順において、シミュレーションモデルに対して時間に依存する変形解析をコンピュータの演算部に実行させ、前記シミュレーションモデルから時間に依存する物理特性をコンピュータの演算部に算出させる。
以上、本実施形態では、せん断粘弾性特性を粘弾性特性の例にとり、応力緩和特性を時間に依存する物理特性の例として詳細を説明したが、本発明は上述の本実施形態に限定はされず、本発明の要旨を逸脱しない範囲において、各種の改良または変更を行っても良いのは、もちろんである。
１０ シミュレーション装置
１２ ＣＰＵ
１４ メモリ
１６ ＲＯＭ
１８ Ｉ／Ｏボード
２０ パラメータ値設定部
２２ 材料粘弾性率取得部
２４ 分布誤差算出部
２６ 粘弾性特性決定部
２８ シミュレーション演算部
３０ ディスプレイ
３２ プリンタ
３４ マウス・キーボード
４０ タイヤモデル
１２ ＣＰＵ
１４ メモリ
１６ ＲＯＭ
１８ Ｉ／Ｏボード
２０ パラメータ値設定部
２２ 材料粘弾性率取得部
２４ 分布誤差算出部
２６ 粘弾性特性決定部
２８ シミュレーション演算部
３０ ディスプレイ
３２ プリンタ
３４ マウス・キーボード
４０ タイヤモデル
Claims (11)
- 材料の粘弾性の周波数分布から、時間の指数関数の和で表された時間領域の粘弾性特性を定めるパラメータ値を算出し、前記時間領域の粘弾性特性を決定する粘弾性特性決定方法であって、
前記時間領域の粘弾性特性を定めるパラメータ値の初期値として所定の値を設定する第１のステップと、
前記材料の粘弾性の周波数分布から、貯蓄弾性率の値および損失弾性率の値を取得する第２のステップと、
前記時間領域の粘弾性特性を定めるパラメータ値により設定された貯蓄弾性率の設定値、前記パラメータ値により設定された損失弾性率の設定値、前記第２のステップで取得した貯蓄弾性率の取得値および前記第２のステップで取得した損失弾性率の取得値から、前記貯蓄弾性率の設定値と取得値との誤差と、前記損失弾性率の設定値と取得値との誤差と、前記貯蓄弾性率の設定値の周波数に対する勾配の値と前記貯蓄弾性率の取得値の周波数に対する勾配の値との誤差と、前記損失弾性率の設定値の周波数に対する勾配の値と前記損失弾性率の取得値の周波数に対する勾配の値との誤差とを周波数毎に求め、この周波数毎の誤差を用いて、前記材料の粘弾性の周波数分布と前記パラメータ値により設定された設定粘弾性の周波数分布との分布誤差を算出する第３のステップと、
前記分布誤差が収束判定条件の閾値より小ならば、前記パラメータ値を用いて前記時間領域の粘弾性特性を決定し、前記分布誤差が前記収束判定条件の閾値より小でなければ、前記パラメータ値を修正して前記第３のステップに戻る第４のステップとを有することを特徴とする粘弾性特性決定方法。 - 前記第３のステップでは、前記貯蓄弾性率の設定値と取得値、前記損失弾性率の設定値と取得値、前記貯蓄弾性率の設定値の周波数に対する勾配の値と前記貯蓄弾性率の取得値の周波数に対する勾配の値および前記損失弾性率の設定値の周波数に対する勾配の値と前記損失弾性率の取得値の周波数に対する勾配の値は、前記貯蓄弾性率および前記損失弾性率を対数変換した値を用いて求められる値であり、
前記周波数分布は、周波数を対数変換した対数周波数分布である請求項１に記載の粘弾性特性決定方法。 - 前記損失弾性率の設定値を前記貯蓄弾性率の取得値で除算して定まる損失正接の設定値、および、前記損失弾性率の取得値を前記貯蓄弾性率の取得値で除算して定まる損失正接の取得値から、前記損失正接の設定値の周波数に対する勾配の値と前記損失正接の取得値の周波数に対する勾配の値との誤差を周波数毎に求め、この誤差を前記分布誤差に含めて用いる請求項１に記載の粘弾性特性決定方法。
- 前記第３のステップでは、前記損失正接の設定値の周波数に対する勾配の値と前記損失正接の取得値の周波数に対する勾配の値は、前記損失正接を対数変換しない真数値を用いて求められる値であり、
前記貯蓄弾性率の設定値と取得値、前記損失弾性率の設定値と取得値、前記貯蓄弾性率の設定値の周波数に対する勾配の値と前記貯蓄弾性率の取得値の周波数に対する勾配の値および前記損失弾性率の設定値の周波数に対する勾配の値と前記損失弾性率の取得値の周波数に対する勾配の値は、前記貯蓄弾性率および前記損失弾性率を対数変換した値を用いて求められる値であり、
前記周波数分布は、周波数を対数変換した対数周波数分布である請求項３に記載の粘弾性特性決定方法。 - 前記損失正接の設定値および取得値は、前記損失正接を対数変換しない真数値であり、
前記損失正接の設定値と取得値との誤差を周波数毎に求め、この誤差を前記分布誤差に含めて用いる請求項３または４に記載の粘弾性特性決定方法。 - 前記分布誤差は、前記誤差のそれぞれに重みづけ係数を乗算したものであり、前記損失正接の設定値と取得値との誤差に乗算する重みづけ係数は、前記誤差に乗算される重みづけ係数の中で最も大きい請求項５に記載の粘弾性特性決定方法。
- 数値解析可能な粘弾性特性を備えたシミュレーションモデルを用いて行なうシミュレーション方法であって、
前記シミュレーションモデルの材料特性の１つとして、請求項１〜６のいずれか１項に記載の粘弾性特性決定方法により決定された、時間の指数関数の和で表された時間領域の粘弾性特性を入力するステップと、
前記シミュレーションモデルを作成するステップと、
前記シミュレーションモデルに対して時間に依存する変形解析を行い、前記シミュレーションモデルから時間に依存する物理特性を算出するステップとを有することを特徴とするシミュレーション方法。 - 前記シミュレーションモデルは、ゴム材料と、このゴム材料を補強するために前記ゴム材料中に配される、無機化合物充填材、金属材及び繊維材の少なくともいずれか１つの補強材とを再現した有限要素モデルである請求項７に記載のシミュレーション方法。
- 前記シミュレーションモデルは、タイヤモデルである請求項７または８に記載のシミュレーション方法。
- 数値解析可能な粘弾性特性を備えたシミュレーションモデルを用いて、コンピュータにシミュレーションを行なわせるシミュレーションプログラムであって、
前記シミュレーションモデルの材料特性の１つとして、請求項１〜６のいずれか１項に記載の粘弾性特性決定方法により決定された、時間の指数関数の和で表された時間領域の粘弾性特性をコンピュータに入力させる手順と、
前記シミュレーションモデルを、コンピュータの演算部に作成させる手順と、
前記シミュレーションモデルに対して時間に依存する変形解析をコンピュータの演算部に実行させ、前記シミュレーションモデルから時間に依存する物理特性をコンピュータの演算部に算出させる手順とを有することを特徴とするシミュレーションプログラム。 - 数値解析可能な粘弾性特性を備えたシミュレーションモデルを用いてシミュレーションを行なうシミュレーション装置であって、
前記シミュレーションモデルの材料特性の１つとして、請求項１〜６のいずれか１項に記載の粘弾性特性決定方法により決定された、時間の指数関数の和で表された時間領域の粘弾性特性を入力する手段と、
前記シミュレーションモデルを作成する手段と、
前記シミュレーションモデルに対して時間に依存する変形解析を行い、前記シミュレーションモデルから時間に依存する物理特性を算出する手段とを有することを特徴とするシミュレーション装置。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2007055901A JP5012094B2 (ja) | 2007-03-06 | 2007-03-06 | 粘弾性特性決定方法、シミュレーションプログラムおよびシミュレーション装置、および、粘弾性特性を備えたシミュレーションモデルを用いるシミュレーション方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2007055901A JP5012094B2 (ja) | 2007-03-06 | 2007-03-06 | 粘弾性特性決定方法、シミュレーションプログラムおよびシミュレーション装置、および、粘弾性特性を備えたシミュレーションモデルを用いるシミュレーション方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| JP2008213724A JP2008213724A (ja) | 2008-09-18 |
| JP5012094B2 true JP5012094B2 (ja) | 2012-08-29 |
Family
ID=39834274
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| JP2007055901A Expired - Fee Related JP5012094B2 (ja) | 2007-03-06 | 2007-03-06 | 粘弾性特性決定方法、シミュレーションプログラムおよびシミュレーション装置、および、粘弾性特性を備えたシミュレーションモデルを用いるシミュレーション方法 |
Country Status (1)
| Country | Link |
|---|---|
| JP (1) | JP5012094B2 (ja) |
Families Citing this family (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP4603082B2 (ja) | 2009-02-03 | 2010-12-22 | 株式会社ブリヂストン | ゴム材料の変形挙動予測装置及びゴム材料の変形挙動予測方法 |
| JP2010181342A (ja) * | 2009-02-06 | 2010-08-19 | Bridgestone Corp | ゴム材料の変形挙動予測装置及びゴム材料の変形挙動予測方法 |
Family Cites Families (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP2003083874A (ja) * | 2001-09-10 | 2003-03-19 | Ricoh Co Ltd | 粘弾性材料特性解析方法、システムおよび記録媒体 |
| JP4299733B2 (ja) * | 2003-07-16 | 2009-07-22 | 株式会社ブリヂストン | 粘弾性応答性能予測方法 |
| JP4710261B2 (ja) * | 2004-06-30 | 2011-06-29 | 横浜ゴム株式会社 | タイヤ挙動シミュレーション用処理装置の動作方法 |
| JP4594043B2 (ja) * | 2004-11-15 | 2010-12-08 | 住友ゴム工業株式会社 | ゴム材料のシミュレーション方法 |
-
2007
- 2007-03-06 JP JP2007055901A patent/JP5012094B2/ja not_active Expired - Fee Related
Also Published As
| Publication number | Publication date |
|---|---|
| JP2008213724A (ja) | 2008-09-18 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| EP1674848B1 (en) | Method of simulating rolling tire | |
| US10436669B2 (en) | Interpolation engine for analysis of time-varying load data signals | |
| JP4931430B2 (ja) | タイヤの温度分布予測方法とタイヤの温度分布予測計算プログラム | |
| JP4152338B2 (ja) | タイヤ設計方法およびプログラム | |
| JP4469172B2 (ja) | タイヤのシミュレーション方法 | |
| JP4285991B2 (ja) | タイヤ経時変化予測方法、タイヤ特性予測方法、タイヤ設計方法、タイヤ製造方法およびプログラム | |
| JP3897477B2 (ja) | 応力−ひずみ関係シミュレート方法およびスプリングバック量予測方法 | |
| WO2013042600A1 (ja) | Chabocheモデルを用いた応力－ひずみ関係シミュレーション方法、応力－ひずみ関係シミュレーションシステム、応力－ひずみ関係シミュレーションプログラム | |
| Cho et al. | Numerical investigation of tire standing wave using 3-D patterned tire model | |
| JP2005306174A (ja) | タイヤ性能予測方法 | |
| Cho et al. | Finite element estimation of hysteretic loss and rolling resistance of 3-D patterned tire | |
| Hernandez et al. | Numerical prediction of three-dimensional tire-pavement contact stresses | |
| JP5012094B2 (ja) | 粘弾性特性決定方法、シミュレーションプログラムおよびシミュレーション装置、および、粘弾性特性を備えたシミュレーションモデルを用いるシミュレーション方法 | |
| JP5533181B2 (ja) | タイヤのシミュレーション方法 | |
| JP4961955B2 (ja) | 構造体のシミュレーション方法 | |
| Yang et al. | Optimization of reinforcement turn-up effect on tyre durability and operating characteristics for racing tyre design | |
| JP4486420B2 (ja) | タイヤ経時変化予測方法、装置、プログラム及び媒体 | |
| JP5128853B2 (ja) | 空気入りタイヤの設計方法 | |
| JP5564889B2 (ja) | 粘弾性試験方法、タイヤのシミュレーション方法およびタイヤのシミュレーション装置 | |
| Wei | A finite element based approach to characterising flexible ring tire (FTire) model for extended range of operating conditions | |
| JP6312975B2 (ja) | タイヤの耐久性評価方法及びこれを用いた設計方法 | |
| Gipser | The FTire tire model family | |
| Yintao et al. | Finite element modeling for steel cord analysis in radial tires | |
| JP5304093B2 (ja) | タイヤの転がり抵抗のシミュレーション方法および装置 | |
| JP4913964B2 (ja) | 実測値の推定方法及びそのプログラム |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20100225 |
|
| A977 | Report on retrieval |
Free format text: JAPANESE INTERMEDIATE CODE: A971007 Effective date: 20120214 |
|
| TRDD | Decision of grant or rejection written | ||
| A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 Effective date: 20120508 |
|
| A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 |
|
| A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20120521 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20150615 Year of fee payment: 3 |
|
| R150 | Certificate of patent or registration of utility model |
Free format text: JAPANESE INTERMEDIATE CODE: R150 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20150615 Year of fee payment: 3 |
|
| LAPS | Cancellation because of no payment of annual fees |