(57)【特許請求の範囲】
【請求項１】
材料の粘性率およびヤング率と内部エネルギとの関係を取得する第１取得部と、
前記材料を使用する計算対象物をモデル化した各粒子の位置、密度、速度および内部エネルギの初期値を取得する第２取得部と、
前記第２取得部が取得した内部エネルギおよび前記第１取得部が取得した関係により得られる粘性率およびヤング率により該粘性率を補正した補正粘性率に基づいて、所定の時間が経過した後の前記各粒子の位置、密度、速度および内部エネルギを計算する第１計算部と
を備えるシミュレーション装置。
【請求項２】
前記補正粘性率は、前記粘性率に、前記ヤング率と定数との積を加算して定める請求項１に記載のシミュレーション装置。
【請求項３】
前記補正粘性率は、式（１４）に基づいて定める請求項１に記載のシミュレーション装置。
【数１】
0006544099.tif 000025
【請求項４】
前記第１取得部は、
複数の温度における前記材料の粘性率、ヤング率、比熱および潜熱を取得する第３取得部と、
前記温度のうちの第１温度における材料の内部エネルギを取得する第４取得部と、
前記第３取得部が取得した複数の温度における温度、比熱および潜熱と前記第４取得部が取得した内部エネルギとから、材料の内部エネルギと温度と比熱と潜熱との変化の関係に基づいて前記複数の温度それぞれにおける内部エネルギを計算する第２計算部と
により得られた関係を取得する請求項１から請求項３のいずれか一つに記載のシミュレーション装置。
【請求項５】
材料の粘性率およびヤング率と内部エネルギとの関係を取得し、
前記材料を使用する計算対象物をモデル化した各粒子の位置、密度、速度および内部エネルギの初期値を取得し、
取得した内部エネルギの初期値および取得した関係により得られる粘性率およびヤング率により該粘性率を補正した補正粘性率に基づいて、所定の時間が経過した後の前記各粒子の位置、密度、速度および内部エネルギを計算する処理をコンピュータに実行させるシミュレーションプログラム。
【請求項６】
材料の粘性率およびヤング率と内部エネルギとの関係を取得し、
前記材料を使用する計算対象物をモデル化した各粒子の位置、密度、速度および内部エネルギの初期値を取得し、
取得した内部エネルギの初期値および取得した関係により得られる粘性率およびヤング率により該粘性率を補正した補正粘性率に基づいて、所定の時間が経過した後の前記各粒子の位置、密度、速度および内部エネルギを計算する処理をコンピュータに実行させるシミュレーション方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、シミュレーション装置、シミュレーションプログラムおよびシミュレーション方法に関する。
【背景技術】
【０００２】
溶融した金属を鋳型に流し込んだ状態で固化させて、鋳物を製作する鋳造が使用されている。溶融した金属の流動および冷えて固まる固化の過程をＳＰＨ（Smoothed Particle Method）法を用いて計算するシミュレーション方法が開示されている（非特許文献１）。
【０００３】
鋳型の中で固化した部分については、剛体の運動方程式を使用して解析することにより、シミュレーションに要する時間を短縮する方法が開示されている（特許文献１）。溶融金属の放射冷却の効果を考慮したシミュレーション方法が開示されている（特許文献２）。
【０００４】
弾性体の運動をヤング率を用いてシミュレートするシミュレーション方法が開示されている（非特許文献２）。
【先行技術文献】
【特許文献】
【０００５】
【特許文献１】特開２０１４−１４６３０２号公報
【特許文献２】特開２０１４−２１１７９８号公報
【非特許文献】
【０００６】
【非特許文献１】Paul W. Cleary、Extension of SPH to predict feeding, freezing, and defect creation in low pressure die casting、Applied Mathematical Modelling、２０１０年、３４巻、ｐ．３１８９−３２０１、
【非特許文献２】Yoichi Kawashima and Yuzuru Sakai、Large Deformation Analysis of Hyperelastic Materials Using SPH Method、e-Journal of Soft Materials、２００７年、３巻、ｐ．２１−２８
【発明の概要】
【発明が解決しようとする課題】
【０００７】
非特許文献１および特許文献２の方法では、溶融金属が固化する際にヤング率が変化することを考慮したシミュレーションを行うことができない。特許文献１では、液体と固体の２種類の状態のいずれかに分けて計算を行うため、両者の中間の状態においてヤング率が変化することを考慮したシミュレーションを行うことができない。そのため、これらの方法を使用する場合には、たとえば鋳型の内部で溶融金属の流動が停止する位置、鋳物の表面の窪みの発生等の現象のシミュレーションの精度が低い。
【０００８】
非特許文献２の方法では、ヤング率の影響を考慮する為に、細かい時間刻みで計算を行う必要がある。そのため、シミュレーションに必要な計算量が多くなる。
【０００９】
一つの側面では、精度の高いシミュレーションを少ない計算量で行うシミュレーション方法等を提供することを目的とする。
【課題を解決するための手段】
【００１０】
一態様では、本発明のシミュレーション装置は、材料の粘性率およびヤング率と内部エネルギとの関係を取得する第１取得部と、前記材料を使用する計算対象物をモデル化した各粒子の位置、密度、速度および内部エネルギの初期値を取得する第２取得部と、前記第２取得部が取得した内部エネルギおよび前記第１取得部が取得した関係により得られる粘性率およびヤング率により該粘性率を補正した補正粘性率に基づいて、所定の時間が経過した後の前記各粒子の位置、密度、速度および内部エネルギを計算する第１計算部とを備える。
【発明の効果】
【００１１】
一つの側面では、精度の高いシミュレーションを少ない計算量で行うシミュレーション方法等を提供することができる。
【図面の簡単な説明】
【００１２】
【図１】シミュレーション装置の構成を示す説明図である。
【図２】合金の粘性率μおよびヤング率σの温度依存性を示す説明図である。
【図３】粒子の挙動を示す説明図である。
【図４】物性値ＤＢのレコードレイアウトを示す説明図である。
【図５】粒子ＤＢのレコードレイアウトを示す説明図である。
【図６】シミュレーションプログラムの処理の流れを示すフローチャートである。
【図７】ｎ＋１のｘ、ｖ、ｕ、ρを算出するサブルーチンの処理の流れを示すフローチャートである。
【図８】重み関数の変化を示す説明図である。
【図９】実施の形態２の基本物性値ＤＢのレコードレイアウトを示す説明図である。
【図１０】実施の形態２のプログラムの処理の流れを示すフローチャートである。
【図１１】実施の形態３のシミュレーション装置の動作を示す機能ブロック図である。
【図１２】実施の形態４のシミュレーション装置の構成を示す説明図である。
【発明を実施するための形態】
【００１３】
［実施の形態１］
図１は、シミュレーション装置１０の構成を示す説明図である。シミュレーション装置１０は、ＣＰＵ（Central Processing Unit）１２、主記憶装置１３、補助記憶装置１４、通信部１５、入力部１６、表示部１７およびバスを備える。本実施の形態のシミュレーション装置１０は、汎用のパーソナルコンピューター、タブレット等の情報機器を使用する。
【００１４】
ＣＰＵ１２は、本実施の形態に係るプログラムを実行する演算制御装置である。ＣＰＵ１２には、一または複数のＣＰＵまたはマルチコアＣＰＵ等が使用される。ＣＰＵ１２は、バスを介してシミュレーション装置１０を構成するハードウェア各部と接続されている。
【００１５】
主記憶装置１３は、ＳＲＡＭ（Static Random Access Memory）、ＤＲＡＭ（Dynamic Random Access Memory）、フラッシュメモリ等の記憶装置である。主記憶装置１３には、ＣＰＵ１２が行う処理の途中で必要な情報およびＣＰＵ１２で実行中のプログラムが一時的に保存される。
【００１６】
補助記憶装置１４は、ＳＲＡＭ、フラッシュメモリ、ハードディスクまたは磁気テープ等の記憶装置である。補助記憶装置１４には、ＣＰＵ１２に実行させるプログラム、物性値ＤＢ３１および粒子ＤＢ３２等、プログラムの実行に必要な各種情報が保存される。
【００１７】
通信部１５は、図示しないインターネットまたはイントラネット等のネットワークとの通信を行うインターフェイスである。
【００１８】
入力部１６は、マウス、キーボード、タッチパネル、ペンタブレット、マイク等の機器であり、ユーザによる操作をシミュレーション装置１０が受け付ける際に使用する。表示部１７は、ディスプレイ、プリンタ、プロッタ等の機器であり、シミュレーション結果等を表示する。
【００１９】
図２は、合金の粘性率μおよびヤング率σの温度依存性を示す説明図である。本実施の形態のシミュレーション装置１０は、たとえば溶融した金属を鋳型に流し込んだ状態で固化させて、鋳物を製作する過程のシミュレーションに使用する。鋳物の素材には、アルミニウム合金、銅合金、亜鉛合金、マグネシウム合金、鉄合金等の合金が使用される。アルミニウム合金は、主成分であるアルミニウムに微量のケイ素、マグネシウム、銅などが添加されている。銅合金等他の合金も同様に、主成分に微量の元素が添加され、鋳造に適した特性を有する材料になっている。
【００２０】
図２Ａは、鋳造に使用する合金の粘性率μの温度依存性を示す図である。横軸は温度を、縦軸は粘性率μを示す。図２Ｂは、鋳造に使用する合金のヤング率σの温度依存性を示す図である。横軸は温度を、縦軸はヤング率σを示す。縦軸および横軸の数値は合金の成分によって異なる。
【００２１】
ここで、粘性率μおよびヤング率σについて説明する。粘性率μは、粘性の大小を表す値である。粘性とは、流体内部での流れの速度を一様にする力である。粘性率μの単位は、パスカル秒である。ヤング率σは、引張応力または圧縮応力と、その応力によって生じる歪との比を表す値である。ヤング率σの単位は、パスカルである。粘性率μおよびヤング率σは、材料によって異なる。粘性率μおよびヤング率σは、実験により測定することができる。また、素材の主要な特性として、材料メーカー、公的機関等から提供される場合もある。
【００２２】
合金は、低温では固体、高温では液体である。合金は、温度に応じて液体と固体との間で可逆的に変化する。純物質であれば、融点を境界にして固体と液体との間で明確な相転移が生じ、粘性率μおよび弾性率σも急激に変化する。しかし、複数の元素が混ざった合金では図に示すように固相線温度Ｔｓと液相線温度ＴＬとに挟まれた中間状態が存在する。固相線温度Ｔｓよりも低い温度では合金は固体であり、液相線温度ＴＬよりも高い温度では合金は液体である。
【００２３】
固体の合金では、粘性率μおよび弾性率σは大きな値である。一方、液体の合金では粘性率μは小さい値であり、弾性率σはほぼゼロである。固相線温度Ｔｓと液相線温度ＴＬとの間では、粘性率μおよび弾性率σは固体と液体の中間の値である。
【００２４】
図３は、粒子の挙動を示す説明図である。図３を使用して本実施の形態で使用するシミュレーションのモデルについて説明する。本実施の形態では、ＳＰＨ法を使用する。ＳＰＨ法は、連続体の挙動を解析するシミュレーション手法の一種であり、流体解析に適している。
【００２５】
粒子５１Ａ、粒子５１Ｂ、粒子５１Ｃは、シミュレーション対象の合金をモデル化した粒子を示す。粒子５１Ａ、粒子５１Ｂおよび粒子５１Ｃを区別せずに説明する場合には、粒子５１と記載する。１個の鋳物を鋳造する際の挙動をシミュレーションする際には、数十万個から数百万個の粒子５１を使用する。粒子５１Ａはｉ番目の粒子、粒子５１Ｂはｉ＋１番目の粒子、粒子５１Ｃはｉ＋２番目の粒子である。
【００２６】
それぞれの粒子５１の状態を、位置ベクトルｘ、速度ベクトルｖ、内部エネルギｕおよび密度ρで表現する。位置ベクトルｘは、粒子５１の位置を示す座標であり、３個の要素を含む１次元ベクトルである。速度ベクトルｖは、各座標軸の方向に対する粒子５１の速度を示す、３個の要素を含む１次元ベクトルである。内部エネルギｕは、粒子５１の内部に含まれるエネルギの量を示すスカラー量である。密度ρは、粒子５１の密度を示すスカラー量である。以下の説明では、下付添字は粒子５１の番号を、上付添字は計算回数を示す。たとえば、位置ベクトルｘinは、ｉ番目の粒子５１Ａのｎ回目の位置ベクトルを意味する。また、粒子５１の番号および計算回数が明確である場合および区別する必要が無い場合には、添字の記載を省略する場合がある。
【００２７】
粒子５１Ａのｎ回目の位置ベクトルｘin、速度ベクトルｖin、内部エネルギｕinおよび密度ρinが、ステップ時間ΔＴ後に変化する値を計算する。計算結果を位置ベクトルｘin+1、速度ベクトルｖin+1、内部エネルギｕin+1および密度ρin+1と示す。計算方法については後述する。この計算を、所定の条件を満足するまですべての粒子５１について繰り返すことにより、鋳物を鋳造する際の挙動をシミュレーションする。所定の条件とは、たとえばすべての粒子５１の内部エネルギｕが閾値以下になること、またはあらかじめ定めた回数の計算を終了することである。
【００２８】
適切なステップ時間ΔＴは、シミュレーションを行う対象のモデルに応じて定まる。ステップ時間ΔＴが大きすぎる場合には、シミュレーションの精度が低くなる。また、条件によっては結果が発散して正常にシミュレーションを終了することができなくなる。一方ステップ時間ΔＴが小さすぎる場合には、シミュレーションに要する時間が長くなる。本実施の形態においては、ステップ時間ΔＴに１マイクロ秒を使用する。
【００２９】
図４は、物性値ＤＢ３１のレコードレイアウトを示す説明図である。物性値ＤＢ３１は、材料の内部エネルギｕ、温度Ｔおよび材料の物性値を関連づけるＤＢである。物性値ＤＢ３１は、番号ｋフィールド、内部エネルギｕフィールド、温度Ｔフィールド、粘性率μフィールド、基準密度ρｓフィールド、比熱Ｃｖフィールド、潜熱ｑフィールド、ヤング率σフィールドおよび熱伝導率κフィールドを含む。物性値ＤＢ３１は、１つの内部エネルギｕの値に対して１つのレコードを有する。
【００３０】
番号ｋフィールドには、物性値レコードの通し番号が連番で記録されている。物性値レコードは内部エネルギｕフィールドに記録された値の昇順で配置されている。内部エネルギｕフィールドには、１ｋｇの材料が有する内部エネルギｕが記録されている。内部エネルギｕは、材料の有するエネルギの、所定の状態を基準とした相対値である。内部エネルギｕの単位はジュール毎キログラムである。温度Ｔフィールドには、内部エネルギｕの値に対応する温度Ｔが記録されている。温度Ｔの単位はケルビンである。粘性率μフィールドには、内部エネルギｕの値に対応する粘性率μが記録されている。前述の通り、粘性率μは、粘性、すなわち流体内部での流れの速度を一様にする力の大きさである。粘性率μの単位はパスカル秒である。基準密度ρｓフィールドには、内部エネルギｕの値に対応する材料の基準密度ρｓが記録されている。基準密度ρｓは、外圧を加えない場合の材料の密度である。基準密度ρｓの単位はキログラム毎立法メートルである。
【００３１】
比熱Ｃｖフィールドには、内部エネルギｕの値に対応する比熱Ｃｖが記録されている。比熱は、１キログラムの材料の温度Ｔを１ケルビン上昇させる際に要するエネルギの量である。比熱の単位はジュール毎キログラム毎ケルビンである。潜熱ｑフィールドには、内部エネルギｕの値に対応する潜熱ｑが記録されている。潜熱ｑは、１キログラムの材料を０ケルビンから内部エネルギｕに対応する状態に変化させる際に、材料の相変化に要するエネルギの量である。相変化とは、材料の状態がたとえば固体から液体に変化することを意味する。潜熱ｑの単位はジュール毎キログラムである。
【００３２】
ヤング率σフィールドには、内部エネルギｕの値に対応するヤング率σが記録されている。前述の通り、ヤング率σは、引張応力または圧縮応力と、その応力によって生じる歪との比を表す値である。ヤング率σの単位は、パスカルである。熱伝導率κフィールドには、内部エネルギｕの値に対応する熱伝導率κが記録されている。熱伝導率κは、材料中に温度の勾配がある場合にその勾配に沿って運ばれる熱流束の大きさを示す。熱伝導率κの単位はワット毎メートル毎ケルビンである。熱流束は単位時間に単位面積を横切る熱量である。熱流束の単位はワット毎平方メートルである。
【００３３】
図５は、粒子ＤＢ３２のレコードレイアウトを示す説明図である。粒子ＤＢ３２は、シミュレーション対象の合金をモデル化した粒子５１の番号、計算したステップ数、および粒子の状態を関連付けるＤＢである。粒子ＤＢ３２は、ステップ数ｎフィールド、粒子番号ｉフィールド、位置ベクトルｘフィールド、速度ベクトルｖフィールド、内部エネルギｕフィールドおよび密度ρフィールドを含む。
【００３４】
ステップ数ｎフィールドには、繰り返し計算を行ったステップ数ｎが記録されている。粒子番号ｉフィールドには、粒子５１の番号が記録されている。位置ベクトルｘフィールドには、ステップ数ｎおよび粒子番号ｉに対応する粒子の位置ベクトルｘinが記録されている。速度ベクトルｖフィールドには、ステップ数ｎおよび粒子番号ｉに対応する粒子の速度ベクトルｖinが記録されている。内部エネルギｕフィールドには、ステップ数ｎおよび粒子番号ｉに対応する粒子の内部エネルギｕinが記録されている。密度ρフィールドには、ステップ数ｎおよび粒子番号ｉに対応する粒子の密度ρinが記録されている。
【００３５】
ステップ数ｎが１である場合の各粒子レコードには、シミュレーションの初期条件が記録されている。ステップ数ｎが２である場合の各粒子レコードには、本実施の形態のシミュレーションを１サイクル行って得た、ステップ時間ΔＴ後の各粒子の状態が記録されている。以後、ステップ数ｎがｎである場合の各粒子レコードには、ステップ時間ΔＴとｎの積の時間後の各粒子の状態が記録されている。
【００３６】
たとえば、粒子ＤＢ３２に基づいて、ステップ数ｎごとの速度ベクトルｖinの時間変化を可視化することにより、鋳物を製作する際の液体金属の流動状態のシミュレーション結果を可視化することができる。また、シミュレーション終了時の位置ベクトルｘinの分布を可視化することにより、鋳物表面の窪みおよび内部の空洞のシミュレーション結果を可視化することができる。
【００３７】
図６は、シミュレーションプログラムの処理の流れを示すフローチャートである。図６を使用して、本実施の形態のシミュレーションプログラムの処理の流れを説明する。
【００３８】
ＣＰＵ１２は、粒子ＤＢ３２からシミュレーションを行うモデルの位置ベクトルｘ、速度ベクトルｖ、内部エネルギｕおよび密度ρの初期値を取得する（ステップＳ５０２）。前述の通り初期値は粒子ＤＢ３２のステップ数ｎフィールドが１のレコードに記録されている。
【００３９】
ＣＰＵ１２は、ステップ数ｎを１に設定する（ステップＳ５０３）。ＣＰＵ１２は、ｎ＋１のｘ、ｖ、ｕ、ρを算出するサブルーチンを起動する（ステップＳ５０４）。ｎ＋１のｘ、ｖ、ｕ、ρを算出するサブルーチンは、ステップ時間ΔＴ経過後の位置ベクトルｘ、速度ベクトルｖ、内部エネルギｕおよび密度ρを算出するサブルーチンである。ｎ＋１のｘ、ｖ、ｕ、ρを算出するサブルーチンの処理の流れは後述する。
【００４０】
ＣＰＵ１２はｎに１を加算する（ステップＳ５０５）。ＣＰＵ１２は、粒子ＤＢ３２にレコードを追加して、ステップ数ｎおよびｎ＋１のｘ、ｖ、ｕ、ρを算出するサブルーチンで得た計算結果を記録する（ステップＳ５０６）。
【００４１】
ＣＰＵ１２は、計算が終了したか否かを判定する（ステップＳ５０７）。たとえば、すべての粒子５１の内部エネルギｕが閾値以下になった場合、またはあらかじめ定めた回数の計算を終了した場合等、所定の条件が満たされた場合に、ＣＰＵ１２は計算が終了したと判定する。計算が終了したと判定した場合は（ステップＳ５０７でＹＥＳ）、ＣＰＵ１２は処理を終了する。計算が終了していないと判定した場合は（ステップＳ５０７でＮＯ）、ＣＰＵ１２はステップＳ５０４に戻る。
【００４２】
図７は、ｎ＋１のｘ、ｖ、ｕ、ρを算出するサブルーチンの処理の流れを示すフローチャートである。ｎ＋１のｘ、ｖ、ｕ、ρを算出するサブルーチンは、ステップ時間ΔＴ経過後の位置ベクトルｘ、速度ベクトルｖ、内部エネルギｕおよび密度ρを算出するサブルーチンである。図７を使用して、ｎ＋１のｘ、ｖ、ｕ、ρを算出するサブルーチンの処理の流れを説明する。
【００４３】
ＣＰＵ１２は、物性値ＤＢ３１を参照して、各粒子５１のｎ回目の内部エネルギｕに対応する熱伝達率κ（ｕ）、基準密度ρｓ（ｕ）、粘性率μ（ｕ）およびヤング率σ（ｕ）を算出する（ステップＳ５２１）。具体的には、ＣＰＵ１２は、ｉ番目の粒子５１の内部エネルギｕinをキーとして物性値ＤＢ３１の内部エネルギｕフィールドを検索する。内部エネルギｕinと一致するレコードがある場合には、そのレコードに記録されている熱伝達率κ、基準密度ρｓ、粘性率μおよびヤング率σを、内部エネルギｕinに対応する熱伝達率κ（ｕin）、基準密度ρｓ（ｕin）、粘性率μ（ｕin）およびヤング率σ（ｕin）に代入する。内部エネルギｕ（ｕin）が一致するレコードが無い場合には、近接する内部エネルギｕが記録されているレコードを使用して、補間により内部エネルギｕinに対応する熱伝達率κ（ｕin）、基準密度ρｓ（ｕin）、粘性率μ（ｕin）およびヤング率σ（ｕin）を算出する。補間には、線形補間その他任意の補間手法を使用することができる。
【００４４】
ＣＰＵ１２は、式（１）に基づいてステップ時間ΔＴの半分の時間が経過した後の各粒子５１の位置ベクトルｘを算出する（ステップＳ５２２）。
【数１】
0006544099.tif 000002
【００４５】
なお、以後に説明する数式において同一の記号は同一の意味で使用する。したがって、一度説明した記号の意味は、２回目以降は記載を省略する。
【００４６】
ＣＰＵ１２は、エネルギ保存則を離散化した式（２）に基づいてｎ＋１番目の時刻の各粒子５１の内部エネルギｕを算出する（ステップＳ５２３）。ここで、ｎ＋１番目の時刻は、ｎ番目の時刻からステップ時間ΔＴが経過した後の時刻である。
【数２】
0006544099.tif 000003
【００４７】
重み関数Ｗは、離れた位置に存在する粒子５１同士の影響の程度の重み付けを行う関数である。本実施の形態においては、重み関数Ｗに式（３）に示すスプライン関数を使用する。図８は、重み関数の変化を示す説明図である。
【数３】
0006544099.tif 000004
【００４８】
粒子の影響半径ｈは、初期状態の平均粒子間隔の２倍から３倍の値の定数である。定数βは、重み関数Ｗの全空間積分量を１に調整する値であり、３次元のシミュレーションを行う場合にはπｈ３、２次元のシミュレーションを行う場合は０．７πｈ２である。ここでπは円周率を示す。図８は、３次元のシミュレーションを行う場合のＷを示している。
【００４９】
なお、重み関数Ｗは式（３）に示すスプライン関数に限定しない。式（４）を満足する任意の関数Ｗ（ｒ）を重み関数Ｗに使用することができる。
【数４】
0006544099.tif 000005
【００５０】
ＣＰＵ１２は、ステップＳ５２３で算出したｎ＋１番目の時刻の各粒子５１の内部エネルギｕに対応する粘性率μ（ｕ）およびヤング率σ（ｕ）を算出する（ステップＳ５２４）。算出する方法は、ステップＳ５２１と同様である。
【００５１】
ＣＰＵ１２は、式（５）に基づいてｎ＋１番目の時刻における各粒子の補正粘性率μ´を算出する（ステップＳ５２５）。
【数５】
0006544099.tif 000006
【００５２】
補正粘性率μ´は、粘性率μにヤング率σの影響を加えた値である。本実施の形態においては、定数Ｃｏに１０００を使用する。
【００５３】
補正粘性率μ´の意味について説明する。本実施の形態では、粘性率μの代わりに補正粘性率μ´を使用することにより、擬似的にヤング率σの影響を考慮したシミュレーションを行う。すなわち、溶融した金属が冷えて上昇した材料の粘性率μに、ヤング率σと定数Ａ／Ｃｏとを積算した値を加算することにより、図２を用いて説明したヤング率σの変化の影響をシミュレーション結果に反映することができる。補正粘性率μ´を使用する方法については、後述する。
【００５４】
ＣＰＵ１２は、運動量保存則を離散化した式（６）に基づいてｎ＋１番目の時刻における各粒子の速度ベクトルｖを算出する（ステップＳ５２６）。
【数６】
0006544099.tif 000007
【００５５】
式（６）中の平均圧力ｐの定義を式（７）に示す。
【数７】
0006544099.tif 000008
【００５６】
式（７）中の圧力ｐは、式（８）により求めることができる。
【数８】
0006544099.tif 000009
【００５７】
式（６）中のΠについて説明する。運動量保存則を離散化した場合は、式（６）中のΠは、式（９）のように定義される。
【数９】
0006544099.tif 000010
【００５８】
ヤング率σの影響を考慮するために、式（９）中の粘性率μを前述の補正粘性率μ´に置き換えると、式（１０）を得ることができる。
【数１０】
0006544099.tif 000011
【００５９】
式（６）中のΠに式（１０）および式（５）を代入することにより、ヤング率σの影響を考慮したシミュレーションを行うことができる。
【００６０】
式（６）に戻り、ｎ＋１番目の時刻における各粒子の速度ベクトルｖを算出する方法について説明する。式（６）に、式（７）、式（８）および式（１０）を代入することにより、ｎ＋１番目の時刻における各粒子の速度ベクトルｖを未知数とする連立一次方程式を得ることができる。共役勾配法等のアルゴリズムを使用して連立一次方程式を解く事により、ｎ＋１番目の時刻における各粒子の速度ベクトルｖを求めることができる。連立一次方程式を解くアルゴリズムは従来から使用されているので、説明を省略する。
【００６１】
ＣＰＵ１２は、質量保存則を離散化した式（１１）に基づいてｎ＋１番目の時刻の各粒子５１の密度ρを算出する（ステップＳ５２７）。
【数１１】
0006544099.tif 000012
【００６２】
ＣＰＵ１２は、式（１２）に基づいてｎ＋１番目の時刻の各粒子５１の位置ベクトルｘを算出する（ステップＳ５２８）。
【数１２】
0006544099.tif 000013
【００６３】
ＣＰＵ１２は、ｎ＋１番目の時刻の各粒子５１の位置ベクトルｘ、速度ベクトルｖ、内部エネルギｕおよび密度ρを補助記憶装置１４に記憶する（ステップＳ５２９）。ＣＰＵ１２は、以上で処理を終了する。
【００６４】
本実施の形態によると、ヤング率σの影響を考慮した精度の高いシミュレーションを少ない計算量で行うことができる。
【００６５】
計算量について説明する。特許文献２の方法でヤング率σの影響を考慮したシミュレーションを行う場合のステップ時間ΔＴの目安は、０．１マイクロ秒以下である。これよりも大きなステップ時間ΔＴを使用する場合には、シミュレーションの精度が大幅に低下するか、発散して解が得られない状態になる。
【００６６】
一方、補正粘性率μ´を使用してヤング率σを疑似的に考慮した本実施の形態においては、１マイクロ秒程度のステップ時間ΔＴを使用することができる。したがって、特許文献２の方法と比較した場合には、本実施の形態では同等の精度のシミュレーションを十分の一以下の計算量で行うことができる。
【００６７】
シミュレーションプログラムは、図示しないネットワークを介して通信部１５に接続された大型計算機で実行されても良い。物性値ＤＢ３１および粒子ＤＢ３２は、図示しないネットワークを介して通信部１５に接続されたサーバに記憶されても良い。
【００６８】
シミュレーションを行う対象は、溶融金属の流動および固化の過程に限定しない。粘性と弾性とを兼ね備えた材料、すなわち粘弾性体の解析に、本実施の形態のシミュレーション装置を使用することができる。
【００６９】
［実施の形態２］
本実施の形態は、物性値ＤＢ３１を作成する方法に関する。なお、実施の形態１と共通する部分の説明は省略する。
【００７０】
図９は、実施の形態２の基本物性値ＤＢのレコードレイアウトを示す説明図である。基本物性値ＤＢは、材料の絶対温度と物性値とを関連づけるＤＢである。基本物性値ＤＢは、温度Ｔフィールド、粘性率μフィールド、基準密度ρｓフィールド、比熱Ｃｖフィールド、潜熱ｑフィールド、ヤング率σフィールド、熱伝導率κフィールドおよび備考フィールドを含む。基本物性値ＤＢは、１つの温度について１つのレコードを有する。基本物性値レコードは温度Ｔフィールドに記録された値の昇順で配置されている。基本物性値ＤＢは、補助記憶装置１４に記憶されている。
【００７１】
温度Ｔフィールドには、温度Ｔが記録されている。温度Ｔの単位はケルビンである。粘性率μフィールドには、温度Ｔに対応する粘性率μが記録されている。基準密度ρｓフィールドには、温度Ｔに対応する材料の基準密度ρｓが記録されている。比熱Ｃｖフィールドには、温度Ｔに対応する比熱Ｃｖが記録されている。潜熱ｑフィールドには、温度Ｔに対応する潜熱ｑが記録されている。ヤング率σフィールドには、温度Ｔに対応するヤング率σが記録されている。熱伝導率κフィールドには、温度Ｔに対応する熱伝導率κが記録されている。備考フィールドには、備考が記録されている。
【００７２】
基本物性値ＤＢに記録されている各物性値は材料の基本的な物性値であり、他のシミュレーションでも使用する場合が多いデータである。本実施の形態のプログラムは、基本物性値ＤＢを使用して実施の形態１で使用した物性値ＤＢ３１を作成するプログラムである。
【００７３】
図１０は、実施の形態２のプログラムの処理の流れを示すフローチャートである。図１０を使用して、本実施の形態のプログラムの処理の流れを説明する。
【００７４】
ＣＰＵ１２は、図９を使用して説明した基本物性値ＤＢを読み込む（ステップＳ５４１）。ＣＰＵ１２は、１番目のレコードに対応する内部エネルギｕの値ｕ（１）を補助記憶装置１４から取得する（ステップＳ５４２）。ｕ（１）には、任意の定数、たとえばゼロを使用することができる。ｕ（１）は、本実施の形態のプログラム中に定数として設定されていても良い。
【００７５】
ＣＰＵ１２は、カウンタｋを２に設定する（ステップＳ５４３）。ＣＰＵ１２は、式（１３）を使用してｋ番目の内部エネルギｕ（ｋ）を算出する（ステップＳ５４４）。式（１３）の右辺第２項から第４項に含まれる値は、基本物性値ＤＢに記録されている値である。
【数１３】
0006544099.tif 000014
【００７６】
ＣＰＵ１２は、ステップＳ５４１で読み込んだすべてのレコードの処理が完了したか否かを判定する（ステップＳ５４５）。完了していないと判定した場合は（ステップＳ５４５でＮＯ）、ＣＰＵ１２はカウンタｋに１を加算する（ステップＳ５４６）。その後、ＣＰＵ１２はステップＳ５４４に戻る。
【００７７】
すべてのレコードの処理が完了したと判定した場合は（ステップＳ５４５でＹＥＳ）、ＣＰＵ１２はステップＳ５４４で算出した内部エネルギｕと基本物性値ＤＢに記録されている各データとを組み合わせて、物性値ＤＢ３１を作成する。物性値ＤＢ３１は、図４を使用して説明したＤＢと同一のＤＢである。ＣＰＵ１２は、物性値ＤＢ３１を補助記憶装置１４に出力する（ステップＳ５４７）。ＣＰＵ１２は、その後処理を終了する。
【００７８】
本実施の形態によると、他のシミュレーションでも使用する基本的な物性値を記録した汎用の基本物性値ＤＢを使用して、ヤング率σの影響を考慮した精度の高いシミュレーションを行う事ができる。
【００７９】
本実施の形態のプログラムは、図７で説明した実施の形態１のサブルーチンのステップＳ５２１の前に挿入しても良い。この場合には、シミュレーションに必要な温度範囲の内部エネルギのみを計算しても良い。また、本実施の形態のプログラムをあらかじめ実行して、補助記憶装置１４に物性値ＤＢ３１を記憶しておいても良い。更に、基本物性値ＤＢに記録されていない温度範囲のデータを補間した物性値ＤＢ３１を補助記憶装置１４に記憶しておいても良い。
【００８０】
［実施の形態３］
図１１は、実施の形態３のシミュレーション装置の動作を示す機能ブロック図である。シミュレーション装置１０は、ＣＰＵ１２による制御に基づいて以下のように動作する。第１取得部４１は、物性値ＤＢ３１より材料の粘性率およびヤング率と内部エネルギとの関係を取得する。第２取得部４２は、粒子ＤＢ３２より材料を使用した計算対象物をモデル化した各粒子の位置、密度、速度および内部エネルギの初期値を取得する。計算部４３は、流体の方程式中の粘性率に補正粘性率を使用して、所定の時間が経過した後の前記各粒子の位置、密度、速度および内部エネルギを計算して粒子ＤＢ３２に記録する。
【００８１】
ここで、流体の方程式は、各粒子の位置、密度、速度および内部エネルギの間の関係を示す方程式である。また、補正粘性率は、第２取得部４２が取得した内部エネルギに対して第１取得部４１が取得した関係に基づいて関係付けられた粘性率を、第２取得部４２が取得した内部エネルギに対して第１取得部４１が取得した関係に基づいて関係付けられたヤング率を用いて補正した値である。
【００８２】
［実施の形態４］
実施の形態４は、汎用のコンピュータとプログラム３８とを組み合わせて動作させることにより、シミュレーション装置１０を実現する形態に関する。図１２は、実施の形態４のシミュレーション装置１０の構成を示す説明図である。図１２を使用して、本実施の形態の構成を説明する。なお、実施の形態１と共通する部分の説明は省略する。
【００８３】
本実施の形態のシミュレーション装置１０は、ＣＰＵ１２、主記憶装置１３、補助記憶装置１４、通信部１５、入力部１６、表示部１７、読取部２７およびバスを備える。シミュレーション装置１０は、汎用のパソコン等の情報処理装置である。
【００８４】
プログラム３８は、可搬型記録媒体２９に記録されている。ＣＰＵ１２は、読取部２７を介してプログラム３８を読み込み、補助記憶装置１４に保存する。またＣＰＵ１２は、シミュレーション装置１０内に実装されたフラッシュメモリ等の半導体メモリ２８に記憶されたプログラム３８を読出しても良い。さらに、ＣＰＵ１２は、通信部１５および図示しないネットワークを介して接続される図示しない他のサーバコンピュータからプログラム３８をダウンロードして補助記憶装置１４に保存しても良い。
【００８５】
プログラム３８は、入力装置１０の制御プログラムとしてインストールされ、主記憶装置１３にロードして実行される。これにより、情報処理装置は上述したシミュレーション装置１０として機能する。
【００８６】
各実施例で記載されている技術的特徴（構成要件）はお互いに組合せ可能であり、組み合わせすることにより、新しい技術的特徴を形成することができる。
今回開示された実施の形態はすべての点で例示であって、制限的なものでは無いと考えられるべきである。本発明の範囲は、上記した意味では無く、特許請求の範囲によって示され、特許請求の範囲と均等の意味および範囲内でのすべての変更が含まれることが意図される。
【００８７】
（付記１）
材料の粘性率およびヤング率と内部エネルギとの関係を取得する第１取得部と、
前記材料を使用する計算対象物をモデル化した各粒子の位置、密度、速度および内部エネルギの初期値を取得する第２取得部と、
前記第２取得部が取得した内部エネルギおよび前記第１取得部が取得した関係により得られる粘性率およびヤング率により該粘性率を補正した補正粘性率に基づいて、所定の時間が経過した後の前記各粒子の位置、密度、速度および内部エネルギを計算する第１計算部と
を備えるシミュレーション装置。
【００８８】
（付記２）
前記補正粘性率は、前記粘性率に、前記ヤング率と定数との積を加算して定める付記１に記載のシミュレーション装置。
【００８９】
（付記３）
前記補正粘性率は、式（１４）に基づいて定める付記１に記載のシミュレーション装置。
【数１４】
0006544099.tif 000015
【００９０】
（付記４）
前記第１取得部は、
複数の温度における前記材料の粘性率、ヤング率、比熱および潜熱を取得する第３取得部と、
前記温度のうちの第１温度における材料の内部エネルギを取得する第４取得部と、
前記第３取得部が取得した複数の温度における温度、比熱および潜熱と前記第４取得部が取得した内部エネルギとから、材料の内部エネルギと温度と比熱と潜熱との変化の関係に基づいて前記複数の温度それぞれにおける内部エネルギを計算する第２計算部と
により得られた関係を取得する付記１から付記３のいずれか一つに記載のシミュレーション装置。
【００９１】
（付記５）
前記第２計算部は、式（１５）に基づいて内部エネルギを計算する付記４に記載のシミュレーション装置。
【数１５】
0006544099.tif 000016
【００９２】
（付記６）
前記第１計算部は、式（１６）から式（２２）を使用して計算を行う付記１から付記３のいずれか一つに記載のシミュレーション装置。
【数１６】
0006544099.tif 000017
【数１７】
0006544099.tif 000018
【数１８】
0006544099.tif 000019
【数１９】
0006544099.tif 000020
【数２０】
0006544099.tif 000021
【数２１】
0006544099.tif 000022
【数２２】
0006544099.tif 000023
【００９３】
（付記７）
前記重み関数は、スプライン関数である付記６に記載のシミュレーション装置。
【００９４】
（付記８）
前記重み関数は、式（２３）で表される付記６に記載のシミュレーション装置。
【数２３】
0006544099.tif 000024
【００９５】
（付記９）
材料の粘性率およびヤング率と内部エネルギとの関係を取得し、
前記材料を使用する計算対象物をモデル化した各粒子の位置、密度、速度および内部エネルギの初期値を取得し、
取得した内部エネルギの初期値および取得した関係により得られる粘性率およびヤング率により該粘性率を補正した補正粘性率に基づいて、所定の時間が経過した後の前記各粒子の位置、密度、速度および内部エネルギを計算する処理をコンピュータに実行させるシミュレーションプログラム。
【００９６】
（付記１０）
材料の粘性率およびヤング率と内部エネルギとの関係を取得し、
前記材料を使用する計算対象物をモデル化した各粒子の位置、密度、速度および内部エネルギの初期値を取得し、
取得した内部エネルギの初期値および取得した関係により得られる粘性率およびヤング率により該粘性率を補正した補正粘性率に基づいて、所定の時間が経過した後の前記各粒子の位置、密度、速度および内部エネルギを計算する処理をコンピュータに実行させるシミュレーション方法。
【符号の説明】
【００９７】
１０  シミュレーション装置
１２  ＣＰＵ
１３  主記憶装置
１４  補助記憶装置
１５  通信部
１６  入力部
１７  表示部
２７  読取部
２８  半導体メモリ
２９  可搬型記録媒体
３１  物性値ＤＢ
３２  粒子ＤＢ
３８  プログラム
４１  第１取得部
４２  第２取得部
４３  計算部
【図１】
0006544099.tif 000026
【図２】
0006544099.tif 000027
【図３】
0006544099.tif 000028
【図４】
0006544099.tif 000029
【図５】
0006544099.tif 000030
【図６】
0006544099.tif 000031
【図７】
0006544099.tif 000032
【図８】
0006544099.tif 000033
【図９】
0006544099.tif 000034
【図１０】
0006544099.tif 000035
【図１１】
0006544099.tif 000036
【図１２】
0006544099.tif 000037