JP2005352818A - 不均質材料のシミュレーションモデル作成方法 - Google Patents

Abstract

【課題】不均質材料のシミュレーションモデル作成方法であって、実行可能なシミュレーションモデルを表す行列を生成する前処理に時間をかけず、従来に比べて前処理を含めた１回のシミュレーション演算時間を短縮させる。
【解決手段】同一形状の有限要素によって不均質材料の代表領域の有限要素モデルを作成し、有限要素モデルの一方の端部に位置する節点の挙動が、他方の端部に位置する対応する節点の挙動に滑らかに繋がるように周期境界条件を定める。次に、周期境界条件を有限要素モデルに付与して拘束される節点の自由度を消去した縮小化行列を、有限要素の材料定数をパラメータとして入力することにより演算可能な行列として生成する処理ルーチンを作成する。代表領域における材料相の分散配置に基づいて有限要素モデルの各有限要素の材料定数を定め、この材料定数の値を前記処理ルーチンに与えて演算可能な縮小化行列を生成する。
【選択図】図２

Description

本発明は、材料特性の異なる複数の材料相が分散配置された不均質材料のシミュレーションモデル作成方法に関する。

タイヤのトレッドゴム部材やビードフィラーゴム部材等に用いられるフィラー配合ゴムや補強部材のバルク特性はタイヤ性能に影響を与える。このため、タイヤ性能のシミュレーション演算及びゴム材料特性のシミュレーション演算は、タイヤ開発の上で重要な要素となっている。バルク特性とは、材料特性の異なる材料相が分散配置された不均質材料の塊（バルク）の状態におけるマクロ特性であり、不均質材料を均質な材料と見なした場合の等価な材料特性である。

上記フィラー配合ゴムのような材料特性の異なる複数の材料相が分散配置された不均質材料のバルク特性をシミュレーション演算により評価するには、例えば不均質材料の代表領域を直方体形状に切り出してシミュレーションモデルを作成することが行われている（特許文献１、２、３）。
そして、このシミュレーションモデルがあたかも上下、左右、奥行き方向に連続して無限に連なってバルク状を成しているように、シミュレーションモデルにおける直方体の対向する面に周期境界条件を付与してシミュレーション演算を行う。

一方、上記シミュレーションモデルでは、有限要素等の単位セルを形成するためのメッシュ分割を容易に行う点から、同一の直方体形状の単位セル（以降、ボクセルともいう）を用いることが行われている（非特許文献１）。ここで、ボクセルは複数の離散点によって形状が特徴付けられたものであり、有限要素モデルでは有限要素をいう。
すなわち、ボクセルで構成したシミュレーションモデルを作成し、この後、ボクセル単位毎に材料相を１つずつ対応させてシミュレーション演算可能なモデル、すなわち、シミュレーション演算可能な全体行列（行列成分が数値化されて各ボクセルを表す要素行列が加算された行列）を作成する。この後、バルク特性を再現するために付与された周期境界条件を用いてシミュレーション演算可能な行列を修正した縮小化行列が前処理において作成され、この縮小化行列を用いてシミュレーション演算（本計算）を行う。

前処理では、例えば、直方体形状の代表領域における、対向する一方の端面上のボクセルの離散点を、他方の端面上のボクセルの離散点の挙動に拘束させる関係式を用いて、拘束される離散点の自由度をシミュレーション可能な全体行列から消去して縮小化行列を作成する。

しかし、この前処理では、代表領域の境界部に位置するボクセル毎に拘束すべき離散点が存在するため、例えば縦、横、奥行き方向にそれぞれ２５６個配列してシミュレーションモデルを構成した場合、拘束すべき離散点の数は極めて多数になる。このため、前処理に要する処理時間は極めて長くなり、シミュレーション演算に要する時間に比べて長くなる場合もある。
このように前処理を含めた１回のシミュレーション計算に時間かかるため、材料相の分散配置を変更し、或いは材料相の占有比率を変更して繰り返しシミュレーション計算を行う場合、極めて長時間の計算を要し、不均質材料の早期開発にとって実用上好ましくない。

特開平９−１８０００２号公報 特開平９−２８８６８９号公報 特開平１０−１１６１３号公報 高野，座古，日本機械学会論文集（Ａ編），61巻、583号，1995，pp.905

そこで、本発明は、実行可能なシミュレーションモデルを表す行列を生成する際、前処理に時間をかけず、前処理を含めた１回のシミュレーション演算時間を短縮させることのできる不均質材料のシミュレーションモデル作成方法を提供することを目的とする。

上記目的を達成するために、本発明は、材料特性の異なる複数の材料相が分散配置された不均質材料のシミュレーションモデルをコンピュータを用いて作成する不均質材料のシミュレーションモデル作成方法であって、複数の離散点によって形状が特徴付けられた同一形状の単位セルを、この単位セルを少なくとも二方向に沿って連続的に複数個配置することによって前記不均質材料の代表領域の離散化モデルを作成するステップと、前記単位セルのそれぞれの配置方向の一方の端部に位置する単位セルの離散点の挙動が、他方の端部に位置する対応する離散点の挙動に滑らかに繋がるように、前記一方の端部に位置する離散点の挙動を前記他方の端部に位置する対応する離散点の挙動に拘束する周期境界条件を定めるステップと、シミュレーション計算を行うための前記離散化モデルを表した全体行列から、前記周期境界条件を前記離散化モデルに付与して拘束される前記離散点の自由度を消去した縮小化行列を、前記単位セルの材料定数をパラメータとして入力することにより演算可能な行列として生成する処理ルーチンを作成するステップと、前記代表領域における材料相の分散配置に基づいて、前記離散化モデルの単位セル毎に材料定数を定めるとともに、定められた材料定数の値を前記処理ルーチンに与えて前記処理ルーチンを実行させることにより、縮小化行列の対応する行列成分に数値を与えてシミュレーション演算可能な縮小化行列を生成するステップと、を有することを特徴とする不均質材料のシミュレーションモデル作成方法を提供する。

前記不均質材料の代表領域及び前記単位セルの形状は、例えば、矩形形状又は直方体形状である。
前記シミュレーション演算可能な縮小化行列を生成するステップにおいて、例えば、前記代表領域の形状が前記離散化モデルの輪郭形状に一致しない場合、前記離散化モデルの少なくとも一方向における形状寸法を定数倍して前記離散化モデルを拡縮することにより、前記離散化モデルの輪郭形状を代表領域の形状に一致させて前記離散化モデルの単位セル毎に材料定数を定め、前記処理ルーチンを実行させる。
また、前記シミュレーション演算可能な縮小化行列を生成するステップにおいて、前記離散化モデルの各単位セルにそれぞれ１つの材料相を割り当て、割り当てられた材料相の材料定数の値を前記処理ルーチンに与えて実行させるのが好ましい。前記離散化モデルの各単位セル毎に材料定数を定める際、単位セルの対応領域に複数の材料相が占める場合、この単位セルの対応する領域に占める材料相の体積占有率又は面積占有率の最も大きい材料相の材料定数の値を前記処理ルーチンに与えて実行させるのが好ましい。

また、前記シミュレーション演算可能な縮小化行列を生成するステップにおいて、前記離散化モデルの単位セル毎に材料定数を定める際、単位セルの対応領域に複数の材料相が占める場合、この単位セルの対応領域に占める材料相の体積占有率又は面積占有率を重み係数として材料定数の加重平均値を求め、この加重平均値を前記処理ルーチンに与えて実行させるのも同様に好ましい。

前記離散化モデルは、例えば前記単位セルを有限要素とした有限要素モデルである。
あるいは、前記離散化モデルは、メッシュフリー法により離散化された部分モデルを有するものでああってもよい。この場合、前記部分モデルは、前記複数の材料相のうち剛性が最も低い材料相の配置部分のモデルを少なくとも含むのが好ましい。

なお、前記代表領域は、複数の輪郭線又は輪郭面が角度を成して接続された角部を有する多角形または多面体形状を有し、前記離散化モデルは、前記角部に対応する位置に前記単位セルの離散点を有し、この離散点の挙動は前記縮小化行列を生成する際拘束される。

前記縮小化行列は、例えば、前記材料相の弾性率を前記パラメータとして各行列成分が数値化される剛性行列であり、前記不均質材料の歪み・応力解析を前記シミュレーション計算によって行う。その際、前記複数の材料相における弾性率のうち、最大の弾性率は、最小の弾性率の１００倍以上である。
前記不均質材料には、カーポン又はシリカ等の粒状の補強材が材料相として分散配置されているのが好ましい。前記不均質材料には、例えば材料特性の異なる複数のエラストマーが材料相として分散配置される。さらに、不均質材料の空隙部に気体や液体が充たされた部分があり、この部分を前記材料相としてもよい。

本発明では、不均質材料の離散化モデルを表した全体行列から、周期境界条件を離散化モデルに付与して拘束される離散点の自由度を消去した縮小化行列を演算可能な行列として生成する処理ルーチンを予め作成しておく。この処理ルーチンでは、単位セルの材料定数をパラメータとして入力することにより演算可能な縮小化行列とする。このため、不均質材料の代表領域のシミュレーションモデルの演算可能な縮小化行列を前記処理ルーチンを用いて短時間に作成することができる。

以下、本発明の不均質材料のシミュレーションモデル作成方法について、添付の図面に示される好適実施例を基に詳細に説明する。

図１は本発明の不均質材料のシミュレーションモデル作成方法を実施する処理装置１０の構成を機能的に示したブロック図である。
処理装置１０は、入力操作系１２、コンピュータ１４及びディスプレイ１６を有する。
入力操作系１２は、マウスやキーボードであり、各種情報をオペレータの指示により入力するデバイスである。

コンピュータ１４は、ＣＰＵ１８、メモリ２０を有し、この他に図示されないＲＯＭ等を有する。コンピュータ１４は、ＲＯＭ等に記憶されたコンピュータソフトウェアを実行することにより、処理ルーチン作成部２２、縮小化行列作成部２４及びシミュレーション演算部２６を機能的に形成し、本発明の不均質材料のシミュレーションモデル作成方法を実施する部分である。
ディスプレイ１６は、入力操作系１２を用いてオペレータが指示できるように入力画面を表示し、又後述する有限要素モデル及びシミュレーション演算結果を表示する部分である。

処理ルーチン作成部２２は、バルク特性を解析するためのシミュレーション演算可能な縮小化行列を作成するための処理ルーチンを作成する部分である。処理ルーチンは、シミュレーション計算を行うための、有限要素モデルを表した全体行列から、周期境界条件によって拘束される節点の自由度を消去した縮小化行列を、各有限要素の材料定数をパラメータとして入力することにより演算可能な行列として生成するサブプログラムである。
処理ルーチン作成部２２は、不均質材料の代表領域を再現した有限要素モデルを有限要素（ボクセル）によって作成し、作成した有限要素モデルに付与する周期境界条件を定めるとともに、さらに周期境界条件によって節点の自由度が縮小化された縮小化行列を生成する上記処理ルーチンを作成する。作成された処理ルーチンは、メモリ２０に記憶保持される。

縮小化行列作成部２４は、各有限要素に対して割り当てられた材料相に従って材料定数を定めるとともに、メモリ２０に記憶保持された処理ルーチンを呼び出し、定められた材料定数の値を処理ルーチンに入力することによりに、縮小化行列の行列成分に数値を割り当て、演算可能な縮小化行列を作成する部分である。縮小化行列は、各行列成分の値がメモリ２０内にアドレス管理されて割り当てられて記憶され、各有限要素毎に材料定数を定めることにより、対応する行列成分のメモリに数値化が振り分けられる。これにより縮小化行列はシミュレーション演算可能な行列となる。

シミュレーション演算部２６は、作成された演算可能な縮小化行列を用いてシミュレーション演算を行う。例えば、縮小化行列は、剛性行列であり、所定の位置に外力や強制変位を付与して応力・歪み解析を行って変形状態を再現したり、最大歪み分布を算出する。あるいは、縮小化行列は熱伝導行列であり、所定の位置に温度を付与して熱伝導を解析して温度分布や線膨脹係数を算出する。あるいは、外力として所定の物理量(変位、温度、濃度など)の勾配に相当する物体力分布を与え、さらに所定の応力分布を考慮した自己釣合い式を演算することで、弾性特性、熱伝導特性、熱弾性特性、粘弾性特性等を算出してもよい。

このような処理装置１０で実施される不均質材料のシミュレーションモデルの作成方法を具体的に説明する。
図２は不均質材料のシミュレーションモデルの作成方法の一例の流れを示すフローチャートである。図３は不均質材料の代表領域の一例の断面図である。この不均質材料は、第１のポリマー相Ａ（図３中白色領域）と、第２のポリマー相Ｂ（図３中黒領域）と、粒系フィラー相Ｃ（図３中の円形状の灰色領域）と、粒系フィラー相の周りを取り巻くフィラー・ポリマー境界相Ｄ（図３中の灰色領域）とが不均質に分散配置されている。すなわち、材料特性の異なる複数のエラストマーが材料相として分散配置され、さらに、粒状の補強材が材料相として分散配置されている。この不均質材料は、複数の材料相における弾性率のうち、例えば最大の弾性率は最小の弾性率の１００倍以上である。不均質材料を構成する材料相は、エラストマーやフィラー等の固体相であるばかりでなく、気体や液体が満たされた空隙相であってもよい。
このような不均質材料についてシミュレーションモデルが作成される。

まず、複数の離散点によって形状が特徴付けられた同一の直方体（６面体）形状のボクセル（単位セル）を、直交する三方向に沿って隣接しかつ連続的に複数個配置することによって不均質材料の代表領域の有限要素モデル（離散化モデル）を作成する（ステップＳ１０）。
不均質材料の代表領域は、例えば直方体形状とし互いに直交する３方向にメッシュ分割することにより、直方体形状の有限要素を連続的に配置し、例えば２５６個×２５６個×２５６個の有限要素により構成した有限要素モデルを作成する。
図４は、直方体形状の有限要素によって構成された有限要素モデルの一例を示している。

次に、有限要素のそれぞれの配置方向の一方の端部に位置する有限要素の節点（離散点）の挙動が、他方の端部に位置する節点（離散点）の挙動に滑らかに繋がるように周期境界条件が定められる（ステップＳ２０）。
次に、定められた周期境界条件を有限要素モデルに付与することにより、有限要素モデルの一方の端部に位置する節点を、他方の端部に位置する節点の挙動に拘束させて、拘束された節点の自由度を、有限要素モデルを表した全体行列から消去した縮小化行列を生成する処理ルーチンが作成される（ステップＳ３０）。

具体的には、各行列成分の値がメモリ２０の所定のアドレスに記憶保持されることによって演算可能な縮小化行列が作成されるが、処理ルーチンはこのアドレス毎に記憶保持する各行列成分の値を算出し、かつ記憶すべきメモリ２０のアドレスを指定する部分である。
上述したように縮小化行列は、有限要素モデルを表した全体行列に周期境界条件を付与されて拘束される節点の自由度が消去されたものである。縮小化行列は、有限要素モデルの各有限要素の行列成分を加算処理してできる全体行列から行列成分が修正されたものである。

この処理ルーチンでは、縮小化行列の行列成分に与える値を各有限要素の材料定数がパラメータとして入力可能に設定されており、このパラメータの値（材料定数の値）が入力されることで、縮小化行列の各行列成分に割り当てる値を自動的に算出し、この値をメモリ２０の各行列成分に応じたアドレスに割り当てて記憶保持させる。
図５（ａ）及び（ｂ）は有限要素モデルおよび各有限要素の行列成分を表した一例である。ここでは、簡易に説明するために有限要素Ａ〜Ｄ、節点１〜９からなる有限要素モデルで説明し、行列は剛性行列とする。また、各節点の自由度は１自由度系であるとする。

図５（ａ）中の有限要素Ａ〜Ｄのそれぞれは、図５（ｂ）に示すような剛性行列とする。
例えば有限要素Ａの１行２列の行列成分ｋ１_1２は、有限要素Ａにおける節点１と節点２との間の剛性を表し、これらの行列成分は有限要素の形状を規定することにより一意的に定められる。例えば、各有限要素の材料相を定めることで、ヤング率、剪断剛性及びポアソン比等の材料定数が定まり、各有限要素の剛性行列の行列成分の値が定まる。
このような剛性行列は、図５（ａ）に示す有限要素の配列に従って、剛性行列が加算処理されて下記[数１]に示す下記全体剛性マトリクスＭとなる。
なお、図５（ｂ）中のSymm.は、剛性行列の対角成分を中心として対称な行列成分を有することを意味する。例えば２行３列の行列成分ｋ２₂₂は３行２列の行列成分と同様のものとして定義される。

また、全体剛性行列Ｍの非ゼロの行列成分をＫ_ijで(i,j＝１〜８の自然数)表して下記[数２]のように表される。

そして、各節点の自由度をｕ_i（i＝１〜９）、外力ベクトルをＦ_i（i＝１〜９）として行列方程式が下記式（１）のように生成される。

一方、図５（ａ）に示す有限要素モデルにおける周期境界条件は、下記[数４]のように表される。節点２と節点８、及び節点４と節点６は、互いに有限要素モデルの端部に位置し、周期境界条件が付与されるべき対応節点である。したがって、節点２と節点８、及び節点４と節点６の自由度は同じ値に拘束される。
また、節点１，３，７，９は有限要素モデルの２方向の境界面の角部に位置するため、節点１と３、及び節点７と９、及び節点１と７、及び節点３及び９は同じ値に拘束される。このような角部に位置する節点の挙動が滑らかに繋がる周期境界条件は、各節点の値が０に拘束されることを意味する。

この周期境界条件のうち、節点１，３，７，９の周期境界条件（ｕ₁＝ｕ₃＝ｕ₇＝ｕ₉＝０）を式（１）に付与することで、下記式（２）となる。さらに、式（２）に節点２と節点８、及び節点４と節点６の周期境界条件を付与し、拘束する節点８及び節点６を消去することにより、下記式（３）のように表される。
この式（３）を有限要素Ａ〜Ｄの各行列成分で表すと下記式（４）のようになる。

式（４）の左辺の剛性行列が、全体剛性行列Ｍに周期境界条件を付与して縮小化された縮小化剛性行列である。この縮小化剛性行列の行列成分は、各有限要素Ａ〜Ｄにおける所定の剛性行列の行列成分の加算和となっているので、各有限要素Ａ〜Ｄにおける剛性行列の行列成分の値、すなわち、有限要素Ａ〜Ｄの材料定数が定められて行列成分の値が既知となると縮小化行列の行列成分を求めることができる。
上記式（４）で表される縮小化行列の行列成分の値を記憶すべきメモリ２０のアドレスが縮小化剛性行列の行列成分毎に設定されており、処理ルーチンに、各有限要素の材料定数がパラメータとして入力されると、各有限要素における剛性行列の行列成分の値を算出し、この算出した値を縮小化行列の行列成分における対応するアドレスに割り当てて記憶保持させる処理を行うサブプログラムである。このようにして、各有限要素の行列成分の値が所定のアドレスに割り当てられ、既に記憶されている値に加算されて記憶される。このように、処理ルーチンは式（４）中の左辺に示す縮小化行列の各行列成分の値を算出して、演算可能とする縮小化行列を作成する。

例えば、下記［数８］のように縮小化行列の行列成分をa1〜a6とすると、このa1〜a6としてメモリ２０に記憶保持するアドレスが設定されており、このアドレスと各有限要素の剛性行列の行列成分との対応関係が下記［数９］のように定められる。

処理ルーチンの実行により、上記対応関係に従って各有限要素の剛性行列の各行列成分の値がメモリ２０に割り当てられる。
処理ルーチンはメモリに記憶保存される。すなわち、ステップＳ１０〜ステップＳ３０によって処理ルーチンが一度作成されると、この処理ルーチンが以降のステップＳ４０〜６０において呼び出されて用いられる。

次に、有限要素（ボクセル）毎に材料定数が設定される（ステップＳ４０）。
具体的には、不均質材料の代表領域が切り出され、この代表領域における材料相の分散配置に基づいて、有限要素モデルの有限要素毎に材料定数が定められる。
図６は、図３に示す代表領域に対して作成された有限要素モデルの一例を示す図である。
図６中の右側拡大図に示すように各有限要素単位に１つの材料相が占めるように設定する。
有限要素モデルの各有限要素毎に材料定数を定める際、図７（ａ）に示すように有限要素の対応領域に複数の材料相（図７（ａ）中の灰色の相及び白色の相）が占める場合、この有限要素の対応する領域に占める材料相の体積占有率又は面積占有率の最も大きい材料相の材料定数の値を定める。図７（ａ）中の領域Ｒは、灰色の相の占有率が白色の相の占有率に比べて大きいため、この領域Ｒの有限要素には、灰色の相の材料定数が割り当てられて材料定数が定められる。図７（ａ）のように複数の材料相が１つの有限要素の領域に占める場合、図７（ｂ）のように１つの材料相が有限要素単位で定められる。
なお、材料定数の設定は上記方法に限定されない。例えば、有限要素の対応領域に占める材料相の体積占有率又は面積占有率を重み係数として各材料相の材料定数の加重平均値を求め、この荷重平均値をこの有限要素における材料定数として定めてもよい。

次に、設定された材料定数の値を、各有限要素毎にステップ３０で作成された処理ルーチンのパラメータとして入力することで、処理ルーチンは上述したように、各有限要素における行列の行列成分の値を算出し、この算出した値を縮小化行列成分における対応するメモリ２０のアドレスに割り当てる（ステップＳ５０）。
こうして作成された縮小化行列は各行列成分に数値が付与された演算可能な行列であり、しかも周期境界条件が付与されている。この縮小化行列に所定の外力が付与されてシミュレーション演算が行われる（ステップＳ６０）。例えば縮小化行列が剛性行列の場合、外力は不均質材料に加わる力又は応力等であり、熱伝導行列の場合、外力は不均質材料に供給される熱量である。

シミュレーション演算は、付与された外力と縮小化行列とから各節点の挙動ｕ_iを算出することによって行われる（ステップＳ６０）。
図８（ａ），（ｂ）は、図６に示す有限要素モデルを用いて縮小化剛性行列を作成し、シミュレーション演算を行った結果を示す図である。
図８（ａ）は、所定の外力が付与されたときの有限要素モデルの変形状態を示した変形図であり、図８（ｂ）はその時の最大主歪み分布である。
図８（ａ）で白色、灰色、黒色領域は、それぞれ図３における第２のポリマー相、第１のポリマー相及びフイラー・ポリマー境界相、粒系フィラー相を示している。そして、有限要素モデルの境界面は外力により変形している。

なお、不均質材料の切り出した代表領域の形状が有限要素モデルの輪郭形状に一致しないときは、有限要素モデルの少なくとも一方向における形状寸法を定数倍して有限要素モデルを拡縮することにより、有限要素モデルの輪郭形状を代表領域の形状に一致させて、有限要素モデルにおける材料定数の設定を行ってもよい。

このような不均質材料の代表領域を再現した有限要素モデルの作成は、予め処理ルーチンを作成し、この処理ルーチンを用いてシミュレーション演算可能な縮小化行列を効率よく作成することができる。
例えば、２５６個×２５６個×１個の３次元の有限要素の配置構成の有限要素モデルにおいて、材料定数の値を有限要素毎に定めることにより、周期境界条件が付与された縮小化行列を作成する処理ルーチンを作成しておき、この処理ルーチンを用いて演算可能な縮小化行列を作成する、ステップＳ４０及び５０に要する処理時間は５分であり、ステップＳ６０におけるシミュレーション演算に要する本計算の処理時間は１４分であり、合計処理時間は１９分である。一方、従来のように演算可能な全体行列を作成した後、周期境界条件を付与して縮小化行列を作成する場合、縮小化行列の作成に要する前処理時間は３０分となる。本計算に要する処理時間は同じ縮小化行列を用いるため上記処理時間と同様の１４分であり、合計処理時間は４４分となる。
特に、材料相の分散配置を種々替えて、あるいは材料相の代表領域における体積占有比率を替えてシミュレーション計算を繰り返し行う場合、演算可能な縮小化行列を繰り返し作成する必要がある。このため、処理ルーチンを何度も用いて縮小行列を短時間に作成することができる。

一方、有限要素の要素行列を有限要素の配列に従って加算処理して全体行列を作成し、この全体行列に周期境界条件を付与するために、ペナルティー数Ｋ（行列成分に比べて非常に大きな数）を所定の行列成分に加えてシミュレーション演算の本計算を行うこともできる。例えば[数２]に示す全体行列の節点４と節点６を拘束する場合、下記[数１０］に示すように４行４列、４行６列、６行４列及び６行６列の行列成分にそれぞれＫ、−Ｋ、Ｋ、−Ｋの値を加算する。

しかし、この場合、上述の縮小化行列のように拘束すべき節点の自由度を消去しないので、全体行列の行列サイズのまま本計算を行う。このため本計算に要する処理時間は、縮小化行列を用いて行う本計算の処理時間は２２分となり、上述したよう本発明における前処理及び本計算の合計処理時間１９分に比べて増大し、合計処理時間１９分に比べて長くなる。
このように、本発明では、縮小化行列を作成する処理時間及びシミュレーション演算処理の処理時間の合計処理時間は従来に比べて短く、効率よくシミュレーション演算を行うことができる。

なお、上記実施形態では、不均質材料の代表領域を３次元領域として、直方体形状（立方体形状を含む）の有限要素を用いて単位セルを構成した３次元有限要素モデルを作成するが、この単位セルは直方体形状に限定されず、４面体形状、６面体形状等、複数の同一形状の単位セルを、単位セルを互いに接合させて少なくとも二方向に沿って連続的に複数個配置することによって不均質材料の代表領域の空間を隙間無く占める形状であればどのようなものであってもよい。また、不均質材料の代表領域を２次元平面領域として有限要素を矩形形状とするものであってもよく、正３角形、正６角形等の２次元形状で構成された２次元の有限要素モデルを作成してもよい。また、ジグソーパズルの１ピースの形状のように、単位セルの輪郭部を互いに接合させて少なくとも二方向に沿って連続的に複数個配置することによって不均質材料の代表領域の領域を隙間無く占める形状であればどのようなものであってもよい。

さらに、単位セルは、有限要素に限らない。例えばメッシュフリー法で規定されるサポートの範囲を単位セルとして、メッシュフリー法によりモデル化してもよい。
メッシュフリー法とは、図９に示すように、複数の要素点（離散点）が連続して配置される状態において、注目する要素点５０ａ（図９中、黒丸）を中心として距離ρ₀ 内の範囲をサポートとし、この範囲において図９に示すようなスプライン関数、例えば４次のスプライン関数で重み関数ｗ（ｘ，ｙ）が設定されている。すなわち、各要素点において距離ρ₀内の範囲をサポートとして単位セルが設定されている。

この重み関数ｗ（ｘ，ｙ）を用いて、注目する要素点を中心としたサポート内の各要素点の位置座標（ｘ_i ，ｙ_i ）（ｉは１〜ｎの整数；ｎはサポート内の要素点の総数）から重み係数を求め、各要素点の位置座標と重み係数とから、物理量を表す内挿関数Ｎ（ｘ，ｙ）を下記式（５）によって設定する。さらに、下記式（６）に従ってサポート内の要素点における物理量の値φ_i（ｉ＝１〜ｎの整数）から物理量の近似関数φ^h （ｘ，ｙ）を定める。この設定方法は、サポート内において、物理量の値φ_iに対する誤差が最小化するように近似関数φ^h （ｘ，ｙ）を定めるものであり、つまり、移動最小二乗法による近似手法によって求められるものである。詳細は、「計算力学ハンドブック 第１巻 有限要素法（構造編集）」（日本機械学会、１９９８年、第３７７〜３７９頁）に記載されている。

この内挿関数を用いて変形する不均質材料の１つの材料相の応力場における物理量の近似関数を表し、この物理量の近似関数を、ガラーキン法により定式化されている、有限要素法で用いられる応力場の支配方程式に代入することで、有限要素法における各有限要素における剛性行列に対応した剛性行列を作成することができる。このように、上記式（５）で作成される内挿関数は、サポート内の要素点の位置情報と重み関数ｗ（ｘ，ｙ）で定まる重み係数によって定まるので、変形を受けて変位した要素点に応じて内挿関数が変化することを特徴とする。この点、モデルを設定した時点で、物理量を表す内挿関数が有限要素の形状に応じて一意的に定まる有限要素法とは異なる。

このように、メッシュフリー法により作成されるモデルでは、内挿関数はサポート内の要素点の位置情報と重み関数ｗ（ｘ，ｙ）で定まる重み係数によって定まるので、変形を受けて変位した要素点に応じて内挿関数が変化するため、大変形の計算において発散しにくいといった特徴を有する。このことから、メッシュフリー法により離散化される部分は、複数の材料相のうち剛性が最も低い材料相の配置部分を少なくとも含むようにし、残りの部分は、有限要素で構成するモデルであってもよい。これにより、剛性が１００倍以上異なるエラストマーとカーボンブラック等の粒状の補強材からなる不均質材料であっても、剛性の弱い材料相の局部的な大変形に対して発散することなくシミュレーション演算を行うことができる。

以上、本発明の不均質材料のシミュレーションモデル作成方法について詳細に説明したが、本発明は上記実施形態に限定されず、本発明の主旨を逸脱しない範囲において、種々の改良や変更をしてもよいのはもちろんである。

本発明の不均質材料のシミュレーションモデル作成方法を実施する処理装置１０の構成を機能的に示したブロック図である。 本発明の不均質材料のシミュレーションモデル作成方法の一例の流れを示すフローチャートである。 不均質材料の代表領域の一例の断面図である。 本発明の不均質材料のシミュレーションモデル作成方法において生成される有限要素モデルの一例を示す図である。 （ａ）及び（ｂ）は有限要素モデルおよび各有限要素の行列成分を説明する図である。 図３に示す代表領域に対して作成された有限要素モデルの一例を示す図である。 （ａ）及び（ｂ）は、本発明の不均質材料のシミュレーションモデル作成方法において材料定数の割り当てを説明する図である。 （ａ），（ｂ）は、シミュレーション演算を行った結果の例を示す図である。 本発明の不均質材料のシミュレーションモデル作成方法において用いられるメッシュフリー法を説明する図である。

符号の説明

１０ 処理装置
１２ 入力操作系
１４ コンピュータ
１６ ディスプレイ
１８ ＣＰＵ
２０ メモリ
２２ 処理ルーチン作成部
２４ 縮小化行列作成部
２６ シミュレーション演算部

Claims (14)

1. 材料特性の異なる複数の材料相が分散配置された不均質材料のシミュレーションモデルをコンピュータを用いて作成する不均質材料のシミュレーションモデル作成方法であって、
   複数の離散点によって形状が特徴付けられた同一形状の単位セルを、この単位セルを少なくとも二方向に沿って連続的に複数個配置することによって前記不均質材料の代表領域の離散化モデルを作成するステップと、
   前記離散化モデルにおける単位セルのそれぞれの配置方向の一方の端部に位置する単位セルの離散点の挙動が、他方の端部に位置する対応する離散点の挙動に滑らかに繋がるように、前記一方の端部に位置する離散点の挙動を前記他方の端部に位置する対応する離散点の挙動に拘束する周期境界条件を定めるステップと、
   シミュレーション計算を行うための前記離散化モデルを表した全体行列から、前記周期境界条件を前記離散化モデルに付与して拘束される前記離散点の自由度を消去した縮小化行列を、前記単位セルの材料定数をパラメータとして入力することにより演算可能な行列として生成する処理ルーチンを作成するステップと、
   前記代表領域における材料相の分散配置に基づいて、前記離散化モデルの単位セル毎に材料定数を定めるとともに、定められた材料定数の値を前記処理ルーチンに与えて前記処理ルーチンを実行させることにより、縮小化行列の対応する行列成分に数値を与えてシミュレーション演算可能な縮小化行列を生成するステップと、を有することを特徴とする不均質材料のシミュレーションモデル作成方法。
2. 前記不均質材料の代表領域及び前記単位セルの形状は、矩形形状又は直方体形状である請求項１に記載の不均質材料のシミュレーションモデル作成方法。
3. 前記シミュレーション演算可能な縮小化行列を生成するステップにおいて、前記代表領域の形状が前記離散化モデルの輪郭形状に一致しない場合、前記離散化モデルの少なくとも一方向における形状寸法を定数倍して前記離散化モデルを拡縮することにより、前記離散化モデルの輪郭形状を代表領域の形状に一致させて前記離散化モデルの単位セル毎に材料定数を定め、前記処理ルーチンを実行させる請求項１又は２に記載の不均質材料のシミュレーションモデル作成方法。
4. 前記シミュレーション演算可能な縮小化行列を生成するステップにおいて、前記離散化モデルの各単位セルにそれぞれ１つの材料相を割り当て、割り当てられた材料相の材料定数の値を前記処理ルーチンに与えて実行させる請求項１〜３のいずれか１項に記載の不均質材料のシミュレーションモデル作成方法。
5. 前記シミュレーション演算可能な縮小化行列を生成するステップにおいて、前記離散化モデルの各単位セル毎に材料定数を定める際、単位セルの対応領域に複数の材料相が占める場合、この単位セルの対応する領域に占める材料相の体積占有率又は面積占有率の最も大きい材料相の材料定数の値を前記処理ルーチンに与えて実行させる請求項４に記載の不均質材料のシミュレーションモデル作成方法。
6. 前記シミュレーション演算可能な縮小化行列を生成するステップにおいて、前記離散化モデルの単位セル毎に材料定数を定める際、単位セルの対応領域に複数の材料相が占める場合、この単位セルの対応領域に占める材料相の体積占有率又は面積占有率を重み係数として材料定数の加重平均値を求め、この加重平均値を前記処理ルーチンに与えて実行させる請求項１〜３のいずれか１項に記載の不均質材料のシミュレーションモデル作成方法。
7. 前記離散化モデルは、前記単位セルを有限要素とした有限要素モデルである請求項１〜６のいずれか１項に記載の不均質材料のシミュレーションモデル作成方法。
8. 前記離散化モデルは、メッシュフリー法により離散化された部分モデルを有する請求項１〜７のいずれか１項に記載の不均質材料のシミュレーションモデル作成方法。
9. 前記部分モデルは、前記複数の材料相のうち剛性が最も低い材料相の配置部分のモデルを少なくとも含む請求項８に記載の不均質材料のシミュレーションモデル作成方法。
10. 前記代表領域は、複数の輪郭線又は輪郭面が角度を成して接続された角部を有する多角形または多面体形状を有し、
    前記離散化モデルは、前記角部に対応する位置に前記単位セルの離散点を有し、この離散点の挙動は前記縮小化行列を生成する際拘束される請求項１〜９のいずれか１項に記載の不均質材料のシミュレーションモデル作成方法。
11. 前記縮小化行列は、前記材料相の弾性率を前記パラメータとして各行列成分が数値化される剛性行列であり、前記不均質材料の歪み・応力解析を前記シミュレーション計算によって行う請求項１〜１０のいずれか１項に記載の不均質材料のシミュレーションモデル作成方法。
12. 前記複数の材料相における弾性率のうち、最大の弾性率は、最小の弾性率の１００倍以上である請求項１１に記載の不均質材料のシミュレーションモデル作成方法。
13. 前記不均質材料には、粒状の補強材が材料相として分散配置されている請求項１〜１２のいずれか１項に記載の不均質材料のシミュレーションモデル作成方法。
14. 前記不均質材料には、材料特性の異なる複数のエラストマーが材料相として分散配置されている請求項１〜１３のいずれか１項に記載の不均質材料のシミュレーションモデル作成方法。

Images