説明

多孔質体の焼結過程のシミュレーション方法およびプログラム

【課題】多孔質体の焼結過程を最初から最後まで正確にシミュレーションすることができるシミュレーション方法を提供する。
【解決手段】シミュレーション装置1は、各面の四方が円柱で囲まれた立方体の形状を有するユニットセルの集合体を多孔質体のモデルとして仮定したSchererモデル、またはM−Sモデルを用いて、多孔質体の焼結速度を算出し、算出した焼結速度に基づいて、塑性変形の項、焼結収縮の項、および自重の項を含む有限要素方程式の解を算出し、算出した有限要素方程式の解に基づいて、多孔質体の変化量を算出することを所定の計算時間刻みによって所定の終了条件を満たすまで繰り返す。

【発明の詳細な説明】
【技術分野】
【0001】
本発明は、多孔質体の焼結過程のシミュレーション方法等に関する。
【背景技術】
【0002】
金属やセラミックの粉末成形体を焼結する方法が、機械部品の製造や光ファイバの製造など、様々な製品の製造方法において用いられている。
【0003】
金属またはセラミックの粉末成形体の焼結をシミュレーションするための従来技術としては、例えば、非特許文献1〜3が存在する。
非特許文献1では、Levy−Misesの式に従う材料であるという条件下で、中央に立方体状の空隙がある立方体をユニットセルとし、焼結体をユニットセルの集合体としたモデルを利用して、焼結体の塑性変形の基礎式を導出している。
非特許文献2では、非特許文献1で導出された基礎式を用いて、有限要素法により、多孔質材料の変形形状を示している。
非特許文献3では、多孔質材料の構成式に焼結応力を導入し、焼結過程における粉末成形体の変形挙動を表す構成式を提案している。また、提案した構成式を有限要素法に組み入れ、自重の影響を考慮して焼結のシミュレーションを行っている。
【先行技術文献】
【非特許文献】
【0004】
【非特許文献1】大矢根守哉、島進、鴻野雄一郎著、「粉末焼結体の塑性基礎式」、日本金属学会論文集(第一部)第39巻第317号、p86〜94、1973年
【非特許文献2】島進、中西利介著、「剛塑性有限要素法による粉末の静水圧成形の解析」、塑性と加工(日本塑性加工学会誌)第29巻第325号、p139〜144、1988年
【非特許文献3】品川一成著、「焼結過程の有限要素シミュレーション」、日本機械学会論文集(A編)第62巻第593号、p240〜252、1996年
【発明の概要】
【発明が解決しようとする課題】
【0005】
しかしながら、非特許文献1〜3などの従来の技術では、完全に緻密な状態での密度に対する焼結前の多孔質体の密度の比である初期相対密度が高い多孔質体のみを対象としており、初期相対密度が低い多孔質体の焼結には適用することができない。
例えば、非特許文献3では、初期相対密度が0.43〜0.61の場合について検討している。
【0006】
本発明は、前述した問題点に鑑みてなされたもので、その目的とすることは、初期相対密度が低い場合においても、多孔質体の焼結過程を正確にシミュレーションすることができるシミュレーション方法等を提供することである。
【課題を解決するための手段】
【0007】
前述した目的を達成するために第1の発明は、多孔質体の焼結過程のシミュレーションを行うシミュレーション方法であって、各面の四方が円柱で囲まれた立方体の形状を有するユニットセルの集合体と仮定した第1の焼結モデルを用いて、多孔質体の焼結速度を算出する第1の算出ステップと、前記第1の算出ステップにおいて算出した焼結速度に基づいて、塑性変形の項、焼結収縮の項、および自重の項を含む有限要素方程式の解を算出し、算出した有限要素方程式の解に基づいて、多孔質体の変化量を算出する第2の算出ステップと、を所定の計算時間刻みによって所定の終了条件を満たすまで繰り返すことを特徴とするシミュレーション方法である。
第1の発明によって、初期相対密度が0.1程度から多孔質体の焼結過程をシミュレーションすることができる。
【0008】
第1の発明は、前記第2の算出ステップの後、所定の計算時間刻みごとに、前記第2の算出ステップにおいて算出した多孔質体の変化量に基づいて、前記第1の焼結モデルのパラメータを修正することが望ましい。
これによって、多孔質体の焼結過程を正確にシミュレーションすることができる。
【0009】
第1の発明における前記第1の算出ステップは、多孔質体の相対密度が所定の値を超えると、前記第1の焼結モデルに代えて、互いに独立で同じサイズの空孔が存在することを仮定した第2の焼結モデルを用いて、多孔質体の焼結速度を算出することが望ましい。
これによって、相対密度が1.0近くになるまで、多孔質体の焼結過程をシミュレーションすることができる。
【0010】
第1の発明は、前記第2の算出ステップの後、所定の計算時間刻みごとに、前記第2の算出ステップにおいて算出した多孔質体の変化量に基づいて、前記第2の焼結モデルのパラメータを修正することが望ましい。
これによって、多孔質体の焼結過程を正確にシミュレーションすることができる。
【0011】
第1の発明における前記第2の算出ステップにおいて算出する多孔質体の変化量は、例えば、少なくとも多孔質体の形状、密度分布、粒子径分布、応力分布、歪み分布、または温度分布のいずれか一つ以上である。
これらは、いずれも工程条件の最適化に必要なデータである。
【0012】
第2の発明は、コンピュータに第1の発明のシミュレーション方法を実行させるためのプログラムである。
【発明の効果】
【0013】
本発明により、初期相対密度が低い場合においても、多孔質体の焼結過程を正確にシミュレーションすることができるシミュレーション方法等を提供することができる。
【図面の簡単な説明】
【0014】
【図1】シミュレーション装置1のハードウエア構成図
【図2】シミュレーション装置1のソフトウエア構成図
【図3】2次元軸対称の物体の一例を示す図
【図4】Schererモデルのユニットセルを示す図
【図5】シミュレーション処理の流れを示すフローチャート
【図6】実プロセスにおける多孔質体の形状の変化を示す図
【図7】シミュレーション装置1による要素分割を示す図
【図8】シミュレーション装置1によって計算された加熱後の形状と温度分布を示す図
【図9】実験によって計測された実プロセスの形状とシミュレーション装置1によって計算された形状の比較を示す図
【図10】シミュレーション装置1によって計算された加熱過程での形状と密度分布を示す図
【0015】
以下図面に基づいて、本発明の実施形態を詳細に説明する。
最初に、図1、図2を参照しながら、シミュレーション装置1の概要を説明する。
図1は、シミュレーション装置1のハードウエア構成図である。尚、図1のハードウエア構成は一例であり、用途、目的に応じて様々な構成を採ることが可能である。
【0016】
シミュレーション装置1は、制御部11、記憶部12、メディア入出力部13、通信制御部14、入力部15、表示部16、周辺機器I/F部17等が、バス18を介して接続される。
【0017】
制御部11は、CPU(Central Processing Unit)、ROM(Read Only Memory)、RAM(Random Access Memory)等で構成される。
【0018】
CPUは、記憶部12、ROM、記録媒体等に格納されるプログラムをRAM上のワークメモリ領域に呼び出して実行し、バス18を介して接続された各装置を駆動制御し、シミュレーション装置1が行う後述する処理を実現する。
ROMは、不揮発性メモリであり、コンピュータのブートプログラムやBIOS等のプログラム、データ等を恒久的に保持している。
RAMは、揮発性メモリであり、記憶部12、ROM、記録媒体等からロードしたプログラム、データ等を一時的に保持するとともに、制御部11が各種処理を行う為に使用するワークエリアを備える。
【0019】
記憶部12は、HDD(ハードディスクドライブ)であり、制御部11が実行するプログラム、プログラム実行に必要なデータ、OS(オペレーティングシステム)等が格納される。プログラムに関しては、OS(オペレーティングシステム)に相当する制御プログラムや、後述する処理をコンピュータに実行させるためのアプリケーションプログラムが格納されている。
これらの各プログラムコードは、制御部11により必要に応じて読み出されてRAMに移され、CPUに読み出されて各種の手段として実行される。
【0020】
メディア入出力部13(ドライブ装置)は、データの入出力を行い、例えば、CDドライブ(−ROM、−R、−RW等)、DVDドライブ(−ROM、−R、−RW等)等のメディア入出力装置を有する。
通信制御部14は、通信制御装置、通信ポート等を有し、コンピュータとネットワーク間の通信を媒介する通信インタフェースであり、ネットワークを介して、他のコンピュータ間との通信制御を行う。ネットワークは、有線、無線を問わない。
【0021】
入力部15は、データの入力を行い、例えば、キーボード、マウス等のポインティングデバイス、テンキー等の入力装置を有する。
入力部15を介して、コンピュータに対して、操作指示、動作指示、データ入力等を行うことができる。
表示部16は、CRTモニタ、液晶パネル等のディスプレイ装置、ディスプレイ装置と連携してコンピュータのビデオ機能を実現するための論理回路等(ビデオアダプタ等)を有する。
【0022】
周辺機器I/F(インタフェース)部17は、コンピュータに周辺機器を接続させるためのポートであり、周辺機器I/F部17を介してコンピュータは周辺機器とのデータの送受信を行う。周辺機器I/F部17は、USBやIEEE1394やRS−232C等で構成されており、通常複数の周辺機器I/Fを有する。周辺機器との接続形態は有線、無線を問わない。
バス18は、各装置間の制御信号、データ信号等の授受を媒介する経路である。
【0023】
図2は、シミュレーション装置1のソフトウエア構成図である。シミュレーション装置1は、プリプロセッサ21、ソルバー22、ポストプロセッサ23を備える。
【0024】
プリプロセッサ21は、有限要素法によってシミュレーションを行うための前処理を行い、条件データ24を生成する。条件データ24は、解析対象の形状、メッシュ分割(要素分割)、材料定義、境界条件定義などシミュレーションに必要なデータである。
【0025】
ソルバー22は、条件データ24に基づいて、多孔質体の焼結速度を算出し、ついで有限要素法によってシミュレーションを行い、結果データ25を生成する。ソルバー22は、時間ステップごとに、解析対象に応じて定式化された有限要素方程式を解く。有限要素方程式の定式化とソルバー22の処理の詳細は後述する。
【0026】
ポストプロセッサ23は、結果データ25を表示部16に表示する。ポストプロセッサ23は、変形図、コンター図、ベクトル図、グラフ、アニメーションなど様々な表示を行う。また、ポストプロセッサ23は、結果データ25の時間ステップごとの解をリスト形式などで表示する。
【0027】
次に、有限要素方程式の定式化について説明する。
本発明の実施の形態においては、粘塑性有限要素法を用いる。以下では、焼結前の多孔質体から焼結後の焼結成形体までを総称して、多孔質体と呼ぶこととする。
【0028】
焼結工程では、多孔質体に焼結変形が生じて収縮する。多孔質体の収縮は、密度分布、粒度分布、温度分布、自重などの影響によって不均一になる。不均一な収縮によって、多孔質体には粘塑性変形も生じる。すなわち、全歪速度は、次式の通り、焼結収縮歪速度と塑性歪速度の和で表せる。
【数1】

【0029】
また、有限要素法では、次式の通り、歪速度は節点速度に比例する。
【数2】

【0030】
塑性変形は多孔質材料の構成式に従うとすると、次式が得られる。
【数3】

【0031】
式(1)と式(2)を式(3)に代入すると、応力は次式で表される。
【数4】

【0032】
仮想仕事の原理から、節点力は次式で表される。
【数5】

【0033】
式(4)を式(5)に代入すると、次式が得られる。
【数6】

【0034】
式(6)が、焼結工程をシミュレーションするための有限要素方程式である。式(6)の右辺第一項は塑性変形、第二項は焼結収縮、第三項は自重の影響を表している。
ソルバー22は、計算時間刻みごとに、式(6)の有限要素方程式を解き、各要素の節点速度を算出する。そして、ソルバー22は、計算時間刻みごとに、多孔質体の変化量(変位、応力、全歪、塑性歪、焼結歪、相対密度など)を算出する。
【0035】
図3は、2次元軸対称の物体の一例を示す図である。図3の座標系は、r−θ―zである。2次元軸対称の物体は、z軸を対称軸とすると、円周方向、すなわちθ方向には物理量が変化しない。
本発明の実施の形態では、説明を簡略化する為、式(6)の有限要素方程式の解析対象は、図3に示すように、解析領域、境界条件、物性などの物理量が周方向に変化しない2次元軸対称とする。従って、以下では、多孔質材料におけるr−z断面31を解析対象とする。但し、本発明の実施の形態に係る技術的思想は、3次元を基にしており、容易に3次元に拡張できることは言うまでもない。
【0036】
2次元軸対称の場合、式(6)のDマトリックスは次式で表される。3次元の場合の式(6)のDマトリックスを式(7´)に示す。
【数7】

【0037】
式(7)における定数A、m、nは、一般的にA=2.5、m=0.5、n=2.5が採用されている。しかしnの値については、本発明では初期相対密度がたとえば0.3以下、さらには0.2以下のものを対象としているため、n=2.5を採用すると良好な計算結果が得られない。そのため本発明では、0.5≦n<2.5が採用される。
【0038】
次に、式(6)における焼結収縮歪速度(焼結速度)を算出するための焼結モデルについて説明する。
【0039】
本発明の実施の形態では、非特許文献4「G.W.Scherer,“Sintering of low−density glasses:I,Theory”,J.American Ceramic Soc.,Vol.60,P.236−239,1977」において提案されている多孔質体のモデル(以下、「Schererモデル」と呼ぶ。)を焼結モデルとして採用する。
Schererモデルでは初期相対密度が、たとえば0.3以下という非常に低い相対密度を表現できるモデルである。
【0040】
図4は、Schererモデルのユニットセルを示す図である。
図4に示すように、Schererモデルのユニットセルは、各面の四方が円柱で囲まれ、円柱の中心軸を立方体の辺とした立方体の形状を有する。従って、各ユニットセルは、中心軸が異なる12個(=立方体の辺の数)の円柱の一部を含む。
Schererモデルの焼結速度は、次式で表される。
【数8】

【0041】
Schererモデルでは、焼結の進行につれて、ユニットセルの辺の長さが短くなり、また円柱の半径が大きくなることで、多孔質体の相対密度の変化を表現する。しかしながら、相対密度が約0.942の状態になると、12個の円柱全てが接触してしまうので、それ以上の状態をこのモデルでは表現できない。そこで、相対密度が0.942以上の状態では、非特許文献5「J. K. Mackenzie and R.
Shuttleworth, “A phenomenonlogical theory of sintering”, Proc. Phys. Soc.(London),B62,P.833,1949」で提案されているモデル(以下、「M−Sモデル」と呼ぶ。)を焼結モデルとして採用する。
M−Sモデルは、焼結が進行した段階で、互いに独立で同じサイズの空孔が流動体内に存在することを仮定した焼結モデルである。
【0042】
M−Sモデルの焼結速度は、次式で表される。
【数9】

【0043】
焼結変形では、周りの変形の影響を受けるため、SchererモデルやM−Sモデルで算出された量を有限要素法の要素ごとに変化させることはできない。そこで、次の計算時間刻みにおけるシミュレーションのために、式(6)の有限要素方程式を解くことで得られる多孔質体の変化量を用いて、焼結モデルのパラメータを有限要素法の要素毎に修正する。
【0044】
Schererモデルにおける相対密度は、次式で表される。
【数10】

【0045】
ソルバー22は、式(6)の有限要素方程式を解くことで得られる相対密度ρ’*を式(11)に代入し、Schererモデルのパラメータであるユニットセルの辺の長さlを修正する。
【0046】
また、M−Sモデルのパラメータである空孔の数nは、Schererモデルのパラメータで表すと「1/l」である。そこで、次式により、空孔の数nを修正する。
【数11】

【0047】
次に、図5を参照しながら、シミュレーション処理を説明する。
図5は、シミュレーション処理の流れを示すフローチャートである。図5に示す処理は、ソルバー22によるものである。
【0048】
シミュレーション装置1の制御部11は、条件データ24を入力する(S101)。条件データ24は、例えば、多孔質体の初期形状、要素形状、材料特性(初期密度、初期粒子径、比熱、粘度など)、境界条件、熱処理条件(加熱温度、加熱部位の移動速度など)、計算時間刻みの間隔、終了条件などである。
【0049】
制御部11は、多孔質体の温度分布を算出する(S102)。温度分布は、例えば、熱処理条件の関数として定義しても良いし、有限要素法による熱伝導解析を行うことで算出しても良い。
【0050】
制御部11は、S102において算出した温度分布に基づいて、粘度を算出する(S103)。粘度は、温度の関数として定義しておく。
【0051】
制御部11は、S103において算出した粘度に基づいて、焼結速度を算出する(S104)。制御部11は、相対密度が94.2%未満の状態であれば、式(8)、(9)に示すSchererモデルによって焼結速度を算出する。相対密度が0.942以上の状態であれば、式(10)に示すM−Sモデルによって焼結速度を算出する。
【0052】
制御部11は、S104において算出した焼結速度に基づいて、式(6)に示す有限要素方程式の解を算出する(S105)。更に、制御部11は、多孔質材料の変化量(変位、応力(全応力、塑性応力、焼結応力)、歪み(全歪、塑性歪、焼結歪)、密度など)を算出する。そして所望の頻度で結果データ25を記憶部12に保存することもできる(S106)。
【0053】
制御部11は、焼結モデルのパラメータを修正する(S107)。制御部11は、相対密度が0.942未満の状態であれば、式(11)によってSchererモデルのパラメータであるユニットセルの辺の長さを修正する。相対密度が0.942%以上の状態であれば、式(12)によってM−Sモデルのパラメータである空孔の数を修正する。
【0054】
制御部11は、計算時間刻みをインクリメントする(S108)。尚、S102において有限要素法による熱伝導解析を行う場合、熱伝導解析と塑性変形解析の計算時間刻みの間隔を異なるものとしても良い。
なお、上記計算時間刻みは、焼結の進行による変形の大きさや熱伝導による温度変化の大きさに依存するため、計算結果を吟味して決める必要がある。通常は0.1秒から100秒の間の値を取れば良い。
【0055】
制御部11は、終了条件を満たすかどうか確認する(S109)。終了条件は、例えば、累計時間、計算ステップ回数、相対密度、多孔質材料の温度など様々な値によって定義することができる。
終了条件を満たさない場合、制御部11は、S102から繰り返す。
終了条件を満たす場合、制御部11は、結果データ25を記憶部12に保存し、処理を終了する。
【0056】
以上の通り、シミュレーション装置1は、各面の四方が円柱で囲まれた立方体の形状を有するユニットセルの集合体を仮定したSchererモデル、またはM−Sモデルを用いて、多孔質体の焼結速度を算出し、算出した焼結速度に基づいて、塑性変形の項、焼結収縮の項、および自重の項を含む有限要素方程式の解を算出し、算出した有限要素方程式の解に基づいて、多孔質体の変化量を算出することを所定の計算時間刻みによって所定の終了条件を満たすまで繰り返す。これによって、多孔質体の焼結過程を最初から最後まで正確にシミュレーションすることができる。
【実施例1】
【0057】
次に、本発明の実施の形態に係る実施例1について説明する。実施例1では、初期相対密度が約0.2以下である多孔質体の焼結工程について、実プロセスと、シミュレーション装置1によるシミュレーション結果の比較を行った。
なお、以下で参照する図7、8、10は、カラー画像をグレースケール変換したものである。
【0058】
図6は、多孔質体を下端から順次加熱して焼結させた場合の実プロセスにおける多孔質体の形状の変化を示す図である。なお、多孔質体は、相対密度が0.157である中心部と、中心部の外周に位置し、相対密度が0.214である外周部からなり、r方向(図3参照)における中心部と外周部の外径比は約1/4である。図6に示すグラフは、多孔質体のr−z断面(図3参照)の形状の変化を示すものである。横軸がr軸、縦軸がz軸である。51は焼結前の初期形状である。55は、焼結後の形状である。
図6に示す焼結後の形状55は、焼結を途中で止めて、多孔質材料を炉から引き上げた状態を示している。従って、焼結後の形状55には、緻密化が終了して径が細くなった領域と、緻密化が途中の領域と、殆ど緻密化が進んでいない領域がある。このような実プロセスとシミュレーション装置1によるシミュレーション結果を比較することで、焼結を終えた最終的な形状のみならず、途中経過の形状の整合性についても評価することができる。
【0059】
表1は、シミュレーションの主な計算条件を示している。
【表1】

【0060】
図7は、シミュレーション装置1による要素分割を示す図である。図7に示す要素は、四角形4節点の要素タイプである。図7に示す多孔質体の形状は、図6に示す実プロセスの焼結前の初期形状51と一致している。
【0061】
図8は、シミュレーション装置1によって計算された加熱後の形状と温度分布を示す図である。図8の原図(カラー画像)では、温度分布が色分けで示されている。
図8に示す多孔質体の形状と、図6に示す実プロセスの焼結後の形状55とを比較すると、両者の形状がほぼ一致していることが分かる。
【0062】
図9は、実プロセスにおける焼結後の形状55とシミュレーション装置1によって計算された焼結後形状の比較を示す図である。図9の横軸は、長手位置、すなわち図6のz軸方向の位置を示している。図9の縦軸は、多孔質材料の半径、すなわち図6のr軸方向の位置を示している。
両者の形状を比較すると、シミュレーション装置1のシミュレーション結果が、実プロセスの焼結後形状を良く再現できていることが分かる。
【実施例2】
【0063】
次に、本発明の実施の形態に係る実施例2について説明する。実施例2では、初期相対密度が0.3である多孔質体の焼結工程のシミュレーションを行った。
なお、多孔質体は、相対密度が1である中心材の外周に形成されており、r方向における中心材と外周部の外径比は約1/4である。
【0064】
図10は、焼結工程のシミュレーション結果を示す図である。図10の原図(カラー画像)では、密度分布が色分けで示されている。
図10は、加熱温度が1460℃の条件でシミュレーションを行い、形状と密度の時間変化を示しており、左に行くほど時間が経過した形状を示している。なお、計算時間刻みは、5秒とした。図10を見ると、加熱時間と共に加熱部位が上昇することによって、緻密化される部分が上方に移動する様子が分かる。また、最終的な焼結形成体の長さは0.1m弱縮んだことが分かる。この縮み量は実プロセスと同程度の量であり、このように中心材の外周に形成された多孔質材料においても妥当な結果が得られていると言える。
【0065】
本発明は、初期相対密度が低い場合(たとえば0.3以下)に特に有効であり、光ファイバの製造工程における多孔質母材の焼結工程などに適用可能である。
【0066】
以上、添付図面を参照しながら、本発明に係るシミュレーション装置等の好適な実施形態について説明したが、本発明はかかる例に限定されない。当業者であれば、本願で開示した技術的思想の範疇内において、各種の変更例又は修正例に想到し得ることは明らかであり、それらについても当然に本発明の技術的範囲に属するものと了解される。
【符号の説明】
【0067】
1………シミュレーション装置
11………制御部
12………記憶部
13………メディア入出力部
14………通信制御部
15………入力部
16………表示部
17………周辺機器I/F部
18………バス
21………プリプロセッサ
22………ソルバー
23………ポストプロセッサ
24………条件データ
25………結果データ

【特許請求の範囲】
【請求項1】
多孔質体の焼結過程のシミュレーションを行うシミュレーション方法であって、
各面の四方が円柱で囲まれた立方体の形状を有するユニットセルの集合体と仮定した第1の焼結モデルを用いて、多孔質体の焼結速度を算出する第1の算出ステップと、
前記第1の算出ステップにおいて算出した焼結速度に基づいて、塑性変形の項、焼結収縮の項、および自重の項を含む有限要素方程式の解を算出し、算出した有限要素方程式の解に基づいて、多孔質体の変化量を算出する第2の算出ステップと、
を所定の計算時間刻みによって所定の終了条件を満たすまで繰り返すことを特徴とするシミュレーション方法。
【請求項2】
前記第2の算出ステップの後、所定の計算時間刻みごとに、前記第2の算出ステップにおいて算出した多孔質材料の変化量に基づいて、前記第1の焼結モデルのパラメータを修正することを特徴とする請求項1に記載のシミュレーション方法。
【請求項3】
前記第1の算出ステップは、多孔質体の相対密度が所定の値を超えると、前記第1の焼結モデルに代えて、互いに独立で同じサイズの空孔が存在することを仮定した第2の焼結モデルを用いて、多孔質体の焼結速度を算出することを特徴とする請求項1または2に記載のシミュレーション方法。
【請求項4】
前記第2の算出ステップの後、所定の計算時間刻みごとに、前記第2の算出ステップにおいて算出した多孔質体の変化量に基づいて、前記第2の焼結モデルのパラメータを修正することを特徴とする請求項3に記載のシミュレーション方法。
【請求項5】
前記第2の算出ステップにおいて算出する多孔質体の変化量は、少なくとも多孔質体の形状、密度分布、または温度分布のいずれか一つであることを特徴とする請求項1から4のいずれか1に記載のシミュレーション方法。
【請求項6】
コンピュータに請求項1から5のいずれか1に記載のシミュレーション方法を実行させるためのプログラム。
【請求項7】
請求項6に記載のプログラムを搭載することを特徴とする装置。

【図1】
image rotate

【図2】
image rotate

【図3】
image rotate

【図5】
image rotate

【図4】
image rotate

【図6】
image rotate

【図7】
image rotate

【図8】
image rotate

【図9】
image rotate

【図10】
image rotate