マトリックスゴムとフィラーとを含むゴム材料の損失正接tanδをコンピュータを用いて計算するゴム材料のシミュレーション方法であって、前記コンピュータに、前記マトリックスゴム及び前記フィラーをそれぞれ有限個の要素で分割した得たゴム材料モデルを入力するステップと、前記コンピュータが、入力された前記ゴム材料モデルの変形条件として、大小2種類の歪を付与するステップと、前記ゴム材料モデルの大小2種類の各歪状態の変形シミュレーションを行い、各要素について、それぞれの状態での垂直歪及びせん断歪を計算し、前記計算された歪を用いて各要素のエネルギーロスを計算するステップと、前記各要素のエネルギーロスの総和、前記大小2種類の歪の差から得られる歪振幅、及びゴム材料モデルの弾性率に基づいてゴム材料モデルの損失正接tanδを計算するステップとを含むことを特徴とするゴム材料のシミュレーション方法。前記コンピュータは、前記各要素の単位体積当たりのエネルギーロスを下記式1で計算する請求項1記載のゴム材料のシミュレーション方法。W=ΣπEεp/22tanδ...1ここで、Σは、垂直歪及びせん断歪の各成分の全ての和を示し、Eは各要素の貯蔵弾性率、εpは大小2種類の歪が作用したときの要素の歪変化量、tanδは各要素の損失正接である。前記コンピュータは、前記ゴム材料モデルの損失正接tanδを下記式2で計算する請求項1又は2記載のゴム材料のシミュレーション方法。tanδ=Wa/{πEaεa2}...2ここで、Waは各要素の単位体積当たりのエネルギーロスの総和、Eaはゴム材料モデルの弾性率、εaは大小2種類の歪の差から得られる歪振幅である。ところで、上述のシミュレーションでは、ゴム材料モデルの変形条件として、図12に示したような正弦波状に変化する歪aを連続的に与えることができない。このため、シミュレーションでは、図11aに示されるように、歪が、一定の時間増分Δtの間隔で離散的にかつ正弦波状に付与される。そして、この離散的に与えられる各歪に応じた応力が収束計算等によりそれぞれ求められる。応力の応答は、例えば図11bに示される。また、応力の変化を調べるために、離散的に得られた応力の各計算点に正弦波でカーブフィットを施して近似波形を求め、この近似波形に基づいて、例えば応力の振幅σ0等が求められていた。しかしながら、上述のようなシミュレーション方法では、応力の波形を精度良く得るためには、図11bに示されるように1周期中に応力の計算数を比較的密に行う必要があり、多くの計算時間が必要であった。また、応力と歪との位相差が安定するためには、複数の周期の計算が必要であるため、やはり多くの計算時間が必要であった。本発明は、以上のような問題点に鑑み案出なされたもので、各要素の歪変化量に基づくエネルギーロスに着目し、これを利用して、損失正接tanδを短時間でかつ精度良く計算しうるゴム材料のシミュレーション方法を提供することを主たる目的としている。本発明のうち請求項1記載の発明は、マトリックスゴムとフィラーとを含むゴム材料の損失正接tanδをコンピュータを用いて計算するゴム材料のシミュレーション方法であって、前記コンピュータに、前記マトリックスゴム及び前記フィラーをそれぞれ有限個の要素で分割した得たゴム材料モデルを入力するステップと、前記コンピュータが、入力された前記ゴム材料モデルの変形条件として、大小2種類の歪を付与するステップと、前記ゴム材料モデルの大小2種類の各歪状態の変形シミュレーションを行い、各要素について、それぞれの状態での垂直歪及びせん断歪を計算し、前記計算された歪を用いて各要素のエネルギーロスを計算するステップと、前記各要素のエネルギーロスの総和、前記大小2種類の歪の差から得られる歪振幅、及びゴム材料モデルの弾性率に基づいてゴム材料モデルの損失正接tanδを計算するステップとを含むことを特徴とする。また請求項2記載の発明は、前記コンピュータは、前記各要素の単位体積当たりのエネルギーロスを下記式1で計算する請求項1記載のゴム材料のシミュレーション方法である。W=ΣπEεp/22tanδ...1ここで、Σは、垂直歪及びせん断歪の各成分の全ての和を示し、Eは各要素の貯蔵弾性率、εpは大小2種類の歪が作用したときの要素の歪変化量、tanδは各要素の損失正接である。また請求項3記載の発明は、前記コンピュータは、前記ゴム材料モデルの損失正接tanδを下記式2で計算する請求項1又は2記載のゴム材料のシミュレーション方法である。tanδ=Wa/{πEaεa2}...2ここで、Waは各要素の単位体積当たりのエネルギーロスの総和、Eaはゴム材料モデルの弾性率、εaは大小2種類の歪の差から得られる歪振幅である。本発明では、マトリックスゴムとフィラーとを含むゴム材料の損失正接tanδをコンピュータを用いて計算するゴム材料のシミュレーション方法において、コンピュータに、マトリックスゴム及びフィラーをそれぞれ有限個の要素で分割した得たゴム材料モデルを入力するステップと、前記コンピュータが、入力された前記ゴム材料モデルの変形条件として、大小2種類の歪を付与するステップと、前記ゴム材料モデルの大小2種類の各歪状態において、それぞれ垂直歪とせん断歪の各成分を計算し、前記計算された歪を用いて各要素の単位体積当たりのエネルギーロスを計算するステップと、前記各要素のエネルギーロスの総和、前記大小2種類の歪の差から得られる歪振幅、及びゴム材料モデルの弾性率とに基づいてゴム材料モデルの損失正接tanδを計算するステップとを含むことを特徴とする。従って、本発明によれば、大小2つの歪の状態での各要素のエネルギーロスから、ゴム材料モデル全体のエネルギーロスを計算することができる。従って、より少ない計算量でゴム材料モデルの損失正接tanδを計算することができる。以下、本発明の実施の一形態が図面に基づき説明する。本発明では、マトリックスゴムとフィラーとを含むゴム材料の損失正接tanδをコンピュータを用いて計算するゴム材料のシミュレーション方法である。図1には、このようなシミュレーション方法を実施するためのコンピュータ装置1が示されている。該コンピュータ装置1は、本体1a、キーボード1b、マウス1c及びディスプレイ装置1dを含んで構成されている。本体1aの内部には、CPU、ROM、作業用メモリー及び磁気ディスク等の大容量記憶装置が設けられる。また、本体1aには、CD-ROMやフレキシブルディスクのドライブ装置1a1、1a2が設けられる。そして、前記大容量記憶装置には後述する本発明のシミュレーション方法を実行するための処理手順プログラムが記憶されている。図2には、本実施形態のシミュレーション方法の処理手順の一例が示される。本実施形態では、先ず、ゴム材料モデルが設定されるステップS1。本実施形態のゴム材料モデルは、例えばタイヤ用のゴム材料を模したものであって、このようなゴム材料は、例えば初期歪10、周波数10Hz、歪片振幅1、温度30°Cの測定条件で、損失正接tanδが0.1~0.4程度の値を示すものが一般的である。図3a及びそのX部拡大図であるbに視覚化して示されるように、ゴム材料モデル2は、例えばカーボンやシリカといったフィラーがマトリックスゴム中に分散して配されたフィラー入りのゴム材料がモデル化される。また、図3aのモデルは、繰り返し最小単位の微視構造としてのゴム材料モデル2の一例が視覚化して示されており、この実施形態において、微視構造としてのゴム材料モデル2は、例えば300nm×300nmの正方形である。該ゴム材料モデル2は、解析しようとするフィラー配合のゴム材料の微小領域が、有限個の小さな要素メッシュ2a、2b、2c...に置き換えられたものである。各要素2a、2b、2c...は、数値解析が可能に定義される。前記数値解析が可能とは、例えば有限要素法、有限体積法、差分法又は境界要素法といった数値解析法により、各要素ないし系全体についての変形計算が可能なことを意味する。具体的には、各要素2a、2b、2c...について、座標系における節点座標値、要素形状、材料特性などが定義される。各要素2a、2b、2c...には、例えば2次元平面としての三角形ないし四辺形の要素、3次元要素としては、例えば4ないし6面体の要素が好ましく用いられる。このゴム材料モデル2は、前記コンピュータ装置1にて取り扱い可能な数値データで構成される。本実施形態のゴム材料モデル2は、後述する変形シミュレーションにおいて平面ひずみ状態の解析、より詳しくはy軸方向の引張変形シミュレーションが行われる。つまり、本実施形態では、z方向にはひずみを持たない2次元のシミュレーションが行われる。また、本実施形態のゴム材料モデル2は、ゴムマトリックス部分がモデル化されたマトリックスモデル3と、このマトリックスモデル3の中に配置されかつ前記フィラーがモデル化されたフィラーモデル4と、各フィラーモデル4の周りを環状に取り囲みかつ前記マトリックスモデル3よりも硬い物性が定義された界面モデル5図3bに示すとを含む。前記マトリックスモデル3は、ゴム材料モデル2の主要部を構成し、かつ、例えば三角形ないし四辺形の複数個の要素を用いて表現されている。変形計算を行うために、マトリックスモデル3を構成する各要素には、例えば、その物性として応力と伸びとの関係を表す関数が定義される。本実施形態のゴム材料のシミュレーション方法では、ゴム弾性応答を表現するために、前記マトリックスモデル3及び界面モデル5のゴム部分は、いずれも分子鎖網目理論に基づいて計算が行われる。なお、分子鎖網目理論については、本件出願人らが提案している特開2010-205165号公報や特開2009-216612号公報に記載されており、これらの文献に記載の方法で本シミュレーションの計算過程に導入されている。従って、分子鎖網目理論については、ここでの詳細な説明は省略する。前記フィラーモデル4は、シリカからなるフィラーが四辺形の複数個の要素を用いてモデル化されたもので、全体として円形に形成されている。三次元モデルの場合、フィラーモデル4は、球形にモデル化されるのが望ましい。また、シリカは、直径約10~300nm程度であり、ゴムに比べて非常に硬い粒子からなる。このため、フィラーモデル4には、このような解析対象となるシリカの物性とほぼ等しい物性が設定される。即ち、本実施形態において、フィラーモデル4は、粘弾性体ではなく弾性体として取り扱われる。また、フィラーモデル4の粒子の個数は、例えば、解析対象のゴム材料のシリカ配合量などに基づいて適宜決定される。前記界面モデル5は、シリカからなるフィラーとマトリックスゴムとを化学的に結合させるシランカップリング剤の働きをシミュレーションに取り込むためにモデル化されたものである。本実施形態の界面モデル5は、フィラーモデル4の周りを小さい厚さtで環状に連続して取り囲むように設定されている。従って、界面モデル5の内周面はフィラーモデル4の外周面に接触しかつ固着されている。本実施形態では、界面モデル5の内周面とフィラーモデル4の外周面とは、互いに剥離しない条件が設定されるが、必要に応じて、予め定めた値以上の応力が生じたときに、フィラーモデル4と界面モデル5との境界を分離させるような条件が設定されても良い。なお、界面モデル5の外周面は、マトリックスゴムモデル3に接触しかつ固着されている。前記界面モデル5の厚さtは、特に限定されるものではないが、種々の実験結果等に鑑み、フィラーモデル4の直径の10~30程度、より好ましくは15~25程度に設定されるのが実際のゴム材料と整合する点で望ましい。また、前記界面モデル5にも、応力と伸びとの関係が定義される。実際のシリカの界面結合材の物性に鑑み、この界面モデル5はマトリックスゴムよりも硬い物性、即ち、マトリックスモデル3よりも伸び難く定義される。ただし、界面モデル5は、フィラーモデル4よりは軟らかいのは言うまでもない。また、マトリックスモデル3、フィラーモデル4及び界面モデル5には、それぞれ、これまでの実験結果等に基づいて、物性についてのパラメータ例えば、密度、弾性率、損失正接tanδ等が入力される。次に、上記マトリックスモデル3、フィラーモデル4及び界面モデル5の複合体であるゴム材料モデル2を変形させるための変形条件が設定されるステップS2。従来、ゴム材料からtanδを実測する場合、ゴム材料モデル2に図4に示されるような周期的な歪が与えられる。また、従来、シミュレーションでゴム材料モデル2のtanδを計算する場合、図11aに示したように、擬似的に、周期的な歪がゴム材料モデルに離散的に与えられる。本発明では、前記ゴム材料モデル2からtanδを計算するためのゴム材料モデルの変形条件として、大小2種類の歪が与えられる。設定される上記大小2種類の歪の値は、特に限定されるものではないが、物性値として本来得たいtanδの測定条件に基づいて決定することが望ましい。例えば、図5に仮想線で示されるように、本来得たいtanδの測定条件が、初期歪10、振幅歪±1の振幅歪であった場合、シミュレーションにおいては、前記振幅歪の最大ピーク歪εmax=11と、最小ピーク歪εmin=9との値が、上記大小2種類の歪の値として設定されるのが好ましい。これにより、上記測定条件により近い変形条件でシミュレーションを行うことが可能になる。次に本実施形態のシミュレーション方法では、上述のように設定されたゴム材料モデル2を用いて変形シミュレーションが行われるステップS3。コンピュータ装置1が行う変形シミュレーションの具体的な処理手順は、図6に示される。変形シミュレーションでは、先ずゴム材料モデル2の各種のデータがコンピュータ装置1に入力されるステップS31。入力されるデータには、各要素に定義された節点の位置や材料特性といった情報が含まれる。コンピュータ装置1では、入力されたデータに基づいて各要素の剛性マトリックスを作成しステップS32、しかる後、全体構造の剛性マトリックスを組み立てるステップS33。全体構造の剛性マトリックスには、既知節点の変位、節点力が導入されステップS34、剛性方程式の解析が行われる。そして、未知節点変位が決定されステップS35、先に定められた変形条件各ひずみに基づいた応力を計算するステップS36。この際、コンピュータ装置1は、先に述べたように、ゴム材料モデル2に、変形条件として2種類の歪の一つを与え、計算を行う。そして、コンピュータ装置1は、それらの結果を記憶装置等に出力するステップ37。ステップS38では、計算を終了させるか否かの判定がなされる。設定された2種類の歪について全て計算が終えていると判断された場合、図1のステップS4に戻る。他方、設定された2種類の歪について全て計算が終えていないと判断された場合には、新たなあ歪の条件でステップS32以降が繰り返される。このようなシミュレーション変形計算は、例えば有限要素法を用いたエンジニアリング系の解析アプリケーションソフトウエア例えば米国リバモアソフトウェアテクノロジー社で開発改良されたLS-DYNA等を用いて行うことができる。また、本シミュレーションは、均質化法漸近展開均質化法に基づいて行われる。均質化法は、図6に示されるように、図3に示した微視構造均質化法ではユニットセルとも呼ばれるを周期的に持っているゴム材料全体Mを表現するxIと、前記微視構造を表現するyIとの独立した2変数が用いられる。微視的スケールと巨視的スケールという異なる尺度の場におけるそれぞれ独立した変数を漸近展開することにより、図3に示した微視構造のモデル構造を反映させたゴム材料全体の平均的な力学応答を近似的に求めることができる。前記変形計算が行われると、コンピュータ装置1は、その結果から必要な物理量を取得するメモリー等に読む込みことができるステップS4。前記物理量としては、前記大小2種類の各歪条件における各要素の垂直歪及びせん断歪が含まれる。図8には、二次元のソリッド要素からなる要素2aを例に挙げ、aにはその垂直歪が示される。該垂直歪は、x軸方向の垂直歪ε1とy軸方向の垂直歪ε2とを含む。三次元の場合、垂直歪には、さらにz軸方向の歪が含まれる。また、図8b及びcには、せん断歪が示され、これには、x軸方向のせん断歪ε3とy軸方向のせん断歪ε4とが含まれる。三次元の場合、せん断歪には、さらにz軸方向の歪が含まれる。次に、ゴム材料モデル2のtanδが計算されるステップS5。tanδの計算の一例は、図9に示されている。本実施形態において、コンピュータ装置1は、各要素2a、2b、2c...の単位体積当たりのエネルギーロスWii=1~要素数を下記式1で計算するステップS51。Wi=ΣπEεp/22tanδi...1式1は、歪の変化に着目して実験によって得られた公知の理論式であり、式中、Σは要素の垂直歪及びせん断歪の各成分の全ての和、Eは要素の貯蔵弾性率、εpは大小2種類の歪εmax、εminが作用したときの要素の歪の変化量、tanδiは要素の損失正接をそれぞれ表している。次に、コンピュータ装置1は、各要素のエネルギーロスWiを計算し終えると、これらの総和Waを計算するステップS52。次に、コンピュータ装置1は、ゴム材料モデル2の損失正接tanδを下記式2で計算するステップS53。tanδ=Wa/{πEaεa2}...2式2は、上記式1の理論式を、ゴム材料モデルの系全体に適用したもので、Waは前記各要素の単位体積当たりのエネルギーロスWiの総和、Eaはゴム材料モデルの弾性率、εaは大小2種類の歪の差から得られる歪振幅{εmax-εmin/2}である。なお、ゴム材料モデルの弾性率Eaは、その応力σと歪εとの関係から、下式によって求める。Ea=σmax-σmin/εmax-εmin上記σは、全要素の応力と要素の体積から求めた平均の応力で、σmax-σminは、その最大と最小の差を意味している。このように、本発明によれば、大小2つの歪の状態での各要素のエネルギーロスから、ゴム材料モデル全体のエネルギーロスを計算し、さらにそこからtanδを計算することができる。従って、より少ない計算量でゴム材料モデルのtanδを計算することができる。図10には、図3に示したゴム材料モデル2を使用して、本発明のシミュレーション方法で求めたtanδと従来手法で求めたtanδとを比較したグラフが示される。これらのシミュレーションの各条件は次の通りである。ゴム材料モデルの初期伸長:10周波数:10Hzひずみの片振幅:1従来手法では、ゴム材料モデルの変形条件として、図11aで示したように、短い時間増分Δt1で徐々に変化する歪が複数回離散的に与えられた。なお、従来手法では、計算される応力が、歪に対する位相差が安定するまで繰り返し計算された。本実施例では、この計算に7周期を要した。他方、本発明では、ゴム材料モデルの変形条件として、図5に示したように、前記振幅歪の最大ピーク歪εmax=11と、最小ピーク歪εmin=9との2種類だけが与えられた。図10から明らかなように、本発明と従来手法とは、tanδに関して、正の相関を有していることが確認できる。また、tanδの計算に要した時間は、従来手法の役1/15であり、大幅な短縮が確認できた。従って、本発明によれば、より少ない計算時間で精度良くゴム材料モデルのtanδを計算予測することができる。