(57)【要約】【課題】解適合格子法及び多重格子法を、分木構造を介して接続し、高速かつ効率的にシミュレーション解析を行なうためのシミュレーション装置、シミュレーション方法及びシミュレーションプログラムを提供する。【解決手段】流体解析装置２０の制御部２１は、解適合格子法を用いて計算格子の細分化を行なう。そして、制御部２１は、ＶサイクルＡＭＧ法による計算処理を実行するための粗格子の生成処理を実行する。ここでは、制御部２１は、処理対象格子を特定し、系譜を用いて１つ上位のレベルの親格子を特定する。そして、この親格子から生成された他の子格子の状態を確認し、各子格子が同じレベルである場合には、分木構造に従って１つ上位のレベルの粗格子に統合し、計算格子を統合して計算手順の記録処理を実行する。そして、このように生成された粗格子を用いて、支配方程式の収束解を算出する。【選択図】図１
2008197921.tif 000002
【特許請求の範囲】
【請求項１】
解適合格子法及び多重格子法の演算処理手段と分木構造情報記憶手段とを備えたシミュレーション装置であって、
前記演算処理手段が、
モデルパラメータを取得し、初期の計算格子を生成し、
物理量の分布に基づいて、前記初期の計算格子を親格子として、分木法により細分化した子格子を関連付けた分木構造を生成して分木構造情報記憶手段に記録し、
前記分木構造情報記憶手段に記録された分木構造において、同じ親格子から生成された子格子であって、リーフ格子となるすべての計算格子を結合して粗格子を生成し、
解適合格子法により前記粗格子に基づいて物理量を計算するための計算式を生成し、
前記計算式に基づいて物理量を計算した計算結果を取得し、
前記分木構造情報記憶手段に記録された分木構造に基づいて細分化された細分格子を取得し、
前記細分格子に基づいて物理量を計算するための連立方程式に対して、前記計算結果を用いて物理量の計算処理を実行することを特徴とするシミュレーション装置。
【請求項２】
前記粗格子の生成は、同じ親格子から生成されたすべての子格子が揃ったときに結合して粗格子を生成することを特徴とする請求項１に記載のシミュレーション装置。
【請求項３】
前記粗格子の生成は、同じ親格子から生成された子格子であって、隣接する子格子が存在する場合には、これらの子格子を結合して粗格子を生成することを特徴とする請求項１に記載のシミュレーション装置。
【請求項４】
前記粗格子に基づいて物理量を計算するための計算式は、前記細分格子の計算式に基づいて代数的多重格子法により生成することを特徴とする請求項１〜３のいずれか一つに記載のシミュレーション装置。
【請求項５】
前記初期の計算格子の細分化は、異相界面を特定し、この異相界面を含む計算格子を細分化することを特徴とする請求項１〜４のいずれか一つに記載のシミュレーション装置。
【請求項６】
解適合格子法及び多重格子法の演算処理手段と分木構造情報記憶手段とを備えたシミュレーション装置を用いてシミュレーション解析を行なう方法であって、
前記演算処理手段が、
モデルパラメータを取得し、初期の計算格子を生成し、
物理量の分布に基づいて、前記初期の計算格子を親格子として、分木法により細分化した子格子を関連付けた分木構造を生成して分木構造情報記憶手段に記録し、
前記分木構造情報記憶手段に記録された分木構造において、同じ親格子から生成された子格子であって、リーフ格子となるすべての計算格子を結合して粗格子を生成し、
解適合格子法により前記粗格子に基づいて物理量を計算するための計算式を生成し、
前記計算式に基づいて物理量を計算した計算結果を取得し、
前記分木構造情報記憶手段に記録された分木構造に基づいて細分化された細分格子を取得し、
前記細分格子に基づいて物理量を計算するための計算式に対して、前記計算結果を用いて物理量の計算処理を実行することを特徴とするシミュレーション方法。
【請求項７】
解適合格子法及び多重格子法の演算処理手段と分木構造情報記憶手段とを備えたシミュレーション装置を用いてシミュレーション解析を行なうためのプログラムであって、
前記演算処理手段を、
モデルパラメータを取得し、初期の計算格子を生成し、
物理量の分布に基づいて、前記初期の計算格子を親格子として、分木法により細分化した子格子を関連付けた分木構造を生成して分木構造情報記憶手段に記録し、
前記分木構造情報記憶手段に記録された分木構造において、同じ親格子から生成された子格子であって、リーフ格子となるすべての計算格子を結合して粗格子を生成し、
解適合格子法により前記粗格子に基づいて物理量を計算するための計算式を生成し、
前記計算式に基づいて物理量を計算した計算結果を取得し、
前記分木構造情報記憶手段に記録された分木構造に基づいて細分化された細分格子を取得し、
前記細分格子に基づいて物理量を計算するための計算式に対して、前記計算結果を用いて物理量の計算処理を実行する手段として機能させることを特徴とするシミュレーションプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、解適合格子を用いたシミュレーション解析において、高速に計算結果を算出するためのシミュレーション装置、シミュレーション方法及びシミュレーションプログラムに関する。
【背景技術】
【０００２】
今日、多様な分野でシミュレーション解析が行われている。このシミュレーション解析では、解析空間に生成した計算格子を用いて、空間の物理法則を表現する方程式群を数値解析により解く。このようなシミュレーション計算では、解析する空間に計算格子を用意し、その計算格子上で、公知の支配方程式（質量保存式、運動量保存式等）を連立させて計算する。
【０００３】
この場合、解析精度を高めるためには、格子間隔を密にして計算する。しかし、必要以上に格子点を増やすと、計算量が多くなるため計算時間がかかり、また多くのメモリ容量が必要になり、計算負荷が大きくなる。
【０００４】
そこで、物体の周囲等のように物理量の変化の激しい場所については、精度良く解析するために局所的に格子間隔を密にすることがある。解析空間に生成された格子間隔の粗な主格子上に、格子間隔の密な副格子を配置する重畳格子法を利用する試みがなされている（例えば、特許文献１を参照。）。この文献記載の流体解析装置では、流れ場が複雑になると予想される場所の主格子上に、密な格子間隔の副格子を自動的に配置する。具体的には、副格子を主格子上に配置する際、副格子を所定の数用意し、主格子の各格子点における物理量の変化率を計算する。そして、その変化率の大きい順番に格子点に優先順位を付与する。次に、優先順位により格子点に副格子が配置されていない場合には、その格子点から変化率の所定の値までの範囲の格子点を覆うように副格子の領域を決定し、副格子の各格子上の変化率の値に応じて副格子の格子間隔を制御することにより、副格子を自動的に主格子上に配置する。
【０００５】
また、効率的な計算のために、解析空間を局所的に細分化するＡＭＲ法（Adaptive Mesh Refinement method）を利用する場合もある（例えば、特許文献２を参照。）。この文
献に記載された数値解析装置においては、初期格子データと分割しきい値、削除しきい値と物理量などの各種データやパラメータを記憶する。この数値解析装置の指針値計算手段は、密度などの物理量をもとに指針値を計算し、分割判定手段は分割しきい値と指針値とを比較する。削除判定手段は、削除しきい値と指針値を比較し、指針値が分割しきい値以上であれば、格子分割手段が要素を分割する。指針値が削除しきい値以下であれば、格子削除手段が詳細要素を削除し元のサイズの要素に戻す。更に、遷移要素作成手段は、要素の生成削除で生じた中間節点を解消するための遷移要素を生成する。中間節点が解消された格子について有限要素法により数値解析を行なう。
【０００６】
支配方程式は計算格子上の代数方程式（大規模連立方程式）に帰着され、これを解くことにより格子点の物理量が求められる。大規模連立方程式の解法としては一般に反復法が用いられる。反復法には、ＳＯＲ（Successive Over-Relaxation）法などの古典的反復法やＣＧ（Conjugate-Gradient）法のようなクリロフ部分空間法など様々な解法があるが、中でも多重格子法（Multi Grid：ＭＧ法）はその収束効率の高さから注目されている。多重格子法では、格子間隔（粗さ）の異なる複数の格子を用意し、各格子上で緩和法を用いることによって各格子のサイズに応じた残差成分（細かい格子における残差の高周波成分
や粗い格子における残差の低周波成分）を除去する。特に、代数的多重格子法（Algebraic Multi Grid：ＡＭＧ法）では、幾何的多重格子法に対して、与えられた連立一次方程式の係数行列から粗格子の係数行列を代数的に生成し、予め粗格子を作成する必要がない。
【特許文献１】特開平６−２２１９５６号公報（第１頁）
【特許文献２】特開２００５−５６２８８号公報（第１頁）
【発明の開示】
【発明が解決しようとする課題】
【０００７】
シミュレーション解析においては、上述のように解適合格子法や多重格子法が利用されている。しかし、解適合格子法や多重格子法は別個に検討されており、相互関係を考慮した検討は行なわれていない。
【０００８】
本発明は、上記課題を解決するためになされたものであり、その目的は、解適合格子法を用いることにより計算負荷を軽減しながら高精度の解析を実現するシミュレーション装置、シミュレーション方法及びシミュレーションプログラムを提供することにある。更に、解適合格子法の特徴を活かしながら多重格子法を用いることにより、高速かつ効率的にシミュレーション解析を実現することを目的とする。
【課題を解決するための手段】
【０００９】
上記問題点を解決するために、請求項１に記載の発明は、解適合格子法及び多重格子法の演算処理手段と分木構造情報記憶手段とを備えたシミュレーション装置であって、前記演算処理手段が、モデルパラメータを取得し、初期の計算格子を生成し、物理量の分布に基づいて、前記初期の計算格子を親格子として、分木法により細分化した子格子を関連付けた分木構造を生成して分木構造情報記憶手段に記録し、前記分木構造情報記憶手段に記録された分木構造において、同じ親格子から生成された子格子であって、リーフ格子となるすべての計算格子を結合して粗格子を生成し、解適合格子法により前記粗格子に基づいて物理量を計算するための計算式を生成し、前記計算式に基づいて物理量を計算した計算結果を取得し、前記分木構造情報記憶手段に記録された分木構造に基づいて細分化された細分格子を取得し、前記細分格子に基づいて物理量を計算するための計算式に対して、前記計算結果を用いて物理量の計算処理を実行することを要旨とする。
【００１０】
請求項２に記載の発明は、請求項１に記載のシミュレーション装置において、前記粗格子の生成は、同じ親格子から生成されたすべての子格子がリーフ格子である場合に、すべての子格子を結合して粗格子を生成することを要旨とする。
【００１１】
請求項３に記載の発明は、請求項１に記載のシミュレーション装置において、前記粗格子の生成は、同じ親格子から生成された子格子であって、隣接する子格子が存在する場合には、これらの子格子を結合して粗格子を生成することを要旨とする。
【００１２】
請求項４に記載の発明は、請求項１〜３のいずれか一つに記載のシミュレーション装置において、前記粗格子に基づいて物理量を計算するための計算式は、前記細分格子の計算式に基づいて代数的多重格子法により生成することを要旨とする。
【００１３】
請求項５に記載の発明は、請求項１〜４のいずれか一つに記載のシミュレーション装置において、前記初期の計算格子の細分化は、混相流における異相界面を特定し、この異相界面を含む計算格子を細分化することを要旨とする。
【００１４】
請求項６に記載の発明は、解適合格子法及び多重格子法の演算処理手段と分木構造情報記憶手段とを備えたシミュレーション装置を用いてシミュレーション解析を行なう方法であって、前記演算処理手段が、モデルパラメータを取得し、初期の計算格子を生成し、物
理量の分布に基づいて、前記初期の計算格子を親格子として、分木法により細分化した子格子を関連付けた分木構造を生成して分木構造情報記憶手段に記録し、前記分木構造情報記憶手段に記録された分木構造において、同じ親格子から生成された子格子であって、リーフ格子となるすべての計算格子を結合して粗格子を生成し、解適合格子法により前記粗格子に基づいて物理量を計算するための計算式を生成し、前記計算式に基づいて物理量を計算した計算結果を取得し、前記分木構造情報記憶手段に記録された分木構造に基づいて細分化された細分格子を取得し、前記細分格子に基づいて物理量を計算するための計算式に対して、前記計算結果を用いて物理量の計算処理を実行することを要旨とする。
【００１５】
請求項７に記載の発明は、解適合格子法及び多重格子法の演算処理手段と分木構造情報記憶手段とを備えたシミュレーション装置を用いてシミュレーション解析を行なうためのプログラムであって、前記演算処理手段を、モデルパラメータを取得し、初期の計算格子を生成し、物理量の分布に基づいて、前記初期の計算格子を親格子として、分木法により細分化した子格子を関連付けた分木構造を生成して分木構造情報記憶手段に記録し、前記分木構造情報記憶手段に記録された分木構造において、同じ親格子から生成された子格子であって、リーフ格子となるすべての計算格子を結合して粗格子を生成し、解適合格子法により前記粗格子に基づいて物理量を計算するための計算式を生成し、前記計算式に基づいて物理量を計算した計算結果を取得し、前記分木構造情報記憶手段に記録された分木構造に基づいて細分化された細分格子を取得し、前記細分格子に基づいて物理量を計算するための計算式に対して、前記計算結果を用いて物理量の計算処理を実行する手段として機能させることを要旨とする。
【００１６】
（作用）
請求項１、６又は７に記載の発明によれば、演算処理手段が、モデルパラメータを取得し、初期の計算格子を生成し、物理量の分布に基づいて、初期の計算格子を親格子として、分木法により細分化した子格子を関連付けた分木構造を生成して分木構造情報記憶手段に記録する。そして、分木構造情報記憶手段に記録された分木構造において、同じ親格子から生成された子格子であって、リーフ格子となるすべての計算格子を結合して粗格子を生成する。そして、粗格子に基づいて物理量を計算するための計算式を生成し、計算式に基づいて物理量を計算した計算結果を取得する。次に、分木構造情報記憶手段に記録された分木構造に基づいて細分化された細分格子を取得し、細分格子に基づいて物理量を計算するための計算式に対して、計算結果を用いて物理量の計算処理を実行する。これにより、解適合格子法及び多重格子法を、分木構造を介してシームレスに接続しているため、詳細かつ効率的なシミュレーション計算を行なうことができる。
【００１７】
請求項２に記載の発明によれば、粗格子の生成は、同じ親格子から生成されたすべての子格子が揃ったときに結合して粗格子を生成する。これにより、分木構造を利用して粗化することができる。
【００１８】
請求項３に記載の発明によれば、粗格子の生成は、同じ親格子から生成された子格子であって、隣接する子格子が存在する場合には、これらの子格子を結合して粗格子を生成する。これにより、粗化のための処理は複雑になるが、より早く収束解を取得することができる。
【００１９】
請求項４に記載の発明によれば、粗格子に基づいて物理量を計算するための計算式は、細分格子の計算式に基づいて代数的多重格子法により生成する。これにより、細分格子の計算式を計算するための行列式を用いて、粗格子の方程式を計算するための行列式を効率的に算出することができる。
【００２０】
請求項５に記載の発明によれば、初期の計算格子の細分化は、異相界面を特定し、この
異相界面を含む計算格子を細分化する。これにより、異相界面の解析を詳細に行なうことができる。
【発明の効果】
【００２１】
本発明によれば、計算負荷を軽減しながら、高精度の解析を実現できる。
【発明を実施するための最良の形態】
【００２２】
以下、本発明を具体化した一実施形態を、図１〜図１１を用いて説明する。本実施形態では、３次元空間において、混相流における異相界面挙動の解析を目的としたシミュレーション装置、シミュレーション方法及びシミュレーションプログラムとして説明する。
【００２３】
本実施形態では、３次元の解析空間を、ＡＭＲ法を利用して計算格子により細分化する。そして、異相界面が存在する場合には、計算格子を細分化する。例えば、図９に示すように、第１相５００（例えば水）に、第２相５０１（例えば気泡）が混在するような系の場合、異相界面（第１相５００と第２相５０１との界面）では計算格子を密に分割する。
【００２４】
本実施形態では、図１に示すように、シミュレーション装置としての流体解析装置２０を用いる。流体解析装置２０は、シミュレーション解析を実行するコンピュータシステムである。この流体解析装置２０は、ＣＰＵ、ＲＡＭ、ＲＯＭ等から構成された制御手段としての制御部２１、計算結果データ記憶部２２を備える。更に、流体解析装置２０は、キーボードやポインティングデバイス等から構成された入力部１１や、ディスプレイ等から構成された出力部１２を備える。
【００２５】
流体解析装置２０の制御部２１は、演算処理手段２１１及びメモリ２１２を含んで構成され、計算格子を用いて後述する処理（初期の計算格子生成段階、分木構造情報の記録段階、粗格子生成段階、連立方程式生成段階、粗格子に基づいて計算処理段階、細分格子取得段階、細分格子に基づいて計算段階の各処理）を行なう。
【００２６】
本実施形態では、計算格子の管理には、分木型データ構造を用いる。分木型データ構造では各計算格子が分岐点となり、各計算格子を細分化することによって分木を伸ばしていく。本実施形態では、３次元では計算格子を等方的に８等分する８分木を用いるが、２次元では等方的に４等分する４分木を用いる。基本になる計算格子を「ルート格子（root cell ）」、分割される計算格子を「親格子（parent cell ）」、分割された計算格子を「子格子（child cell）」、計算格子の世代を「レベル」、末端の計算格子を「リーフ格子（leaf cell ）」と呼ぶ。各計算格子は子格子と親格子へのポインタを持ち、分木を辿ることによって各計算格子間を移動する。
【００２７】
例えば、図７に示すように、このルート格子ＲＣが細分化されたレベルＬ１には、細分格子（子格子）が形成される。本実施形態では、計算格子は８分木構造を用いるため、ルート格子ＲＣを細分化した場合には、８個の子格子が生成される。レベルＬ１の計算格子の細分化を行なった場合には、レベルＬ２において８個の細分格子（孫格子）が生成され、更に、孫格子の細分化を行なった場合にはレベルＬ３において８個の細分格子（曾孫格子）が生成される。
【００２８】
計算を行なう際には、隣接格子など、任意の位置の計算格子を特定する必要がある。この場合、親子関係を辿って計算格子の位置関係を特定するために「系譜」を用いる。ここで、「系譜」とはルート格子から任意の計算格子まで、どの計算格子を通ってきたかという情報を記憶したものである。図８に示すように、３次元の計算格子においては、親格子内での子格子の位置（ｘｙｚ）を各方向に「０」又は「１」で表わす。レベル「ｎ」における計算格子Ｃ１は、レベル「ｎ＋１」において、８個の計算格子Ｃ２に細分化される。
この場合、各計算格子Ｃ２は、計算格子Ｃ１に対するレベル「ｎ＋１」において（０００）〜（１１１）により特定される。また、（１００）の計算格子Ｃ２は、レベル「ｎ＋２」において、８個の計算格子Ｃ３に細分化される。そして、各レベルでの得られた（ｘｙｚ）の値を並べることにより、ルート格子に対して各細分格子の位置を特定することができる。
【００２９】
メモリ２１２は分木構造情報記憶手段として機能し、モデルパラメータ、初期分木構造、計算格子及び物理量に関するデータが記憶される。
モデルパラメータデータ領域には、シミュレーション計算を実行するために必要な条件に関するデータが記録される。モデルパラメータとしては、タイムステップの繰り返し回数、反復法の繰り返し回数、収束判定値、流体の物性値（粘度、密度等）、解析領域の形状構造（例えば、ＣＡＤデータ）等がある。
【００３０】
初期分木構造データ領域には、３次元空間に配置された３次元オブジェクトの形状に基づいて生成された格子に関するデータが記録される。
計算格子データ領域には、計算の際に用いる格子に関するデータが記録される。ここでは、各ルート格子を特定するためのデータ、ルート格子から細分化されて生成された計算格子（子格子、孫格子等）を特定するための系譜データが記録される。更に、各計算格子を特定する系譜データに対して、粗格子生成のための順番に関するデータ（計算手順データ）が記録される。
【００３１】
物理量データ領域には、計算処理途上での３次元空間内の計算格子における物理量に関するデータが記録される。
また、計算結果データ記憶部２２には、各タイムステップにおける計算格子毎に、計算処理により算出された物理量（計算結果）に関するデータが格納される。
【００３２】
このように構成されたシステムを用いて行なうシミュレーション計算処理を、図２〜図６に従って説明する。ここでは、まず、計算処理の全体像（図２）を説明する。そして、この計算処理に用いられる初期格子の設定処理（図３）、連立方程式求解処理（図４、図５、図６）を説明する。
【００３３】
（計算処理の全体像）
まず、図２を用いて計算処理の全体像を説明する。
流体解析装置２０の制御部２１は、初期格子の設定処理を実行する（ステップＳ１−１）。ここでは、後述するように、演算処理手段２１１が、メモリ２１２に記憶された初期値（モデルパラメータや物理量初期値）を用いて、計算処理のための計算格子を生成する。
【００３４】
そして、流体解析装置２０の制御部２１は、非定常解析を開始する（ステップＳ１−２）。これにより、演算処理手段２１１は、タイムステップ識別子を初期値に設定し、後述する時間変化に伴う流体の動きのシミュレーション計算を実行する。
【００３５】
各タイムステップにおいては、流体解析装置２０の制御部２１は、各タイムステップにおける物理量を算出するための連立方程式（計算式）の作成処理を実行する（ステップＳ１−３）。具体的には、演算処理手段２１１は、この系における公知の支配方程式（質量保存式、運動量保存式等）に対して、計算格子の物理量（速度、圧力等）を未知数とする連立代数方程式を生成する。
【００３６】
次に、流体解析装置２０の制御部２１は、連立方程式の求解処理を実行する（ステップＳ１−４）。本実施形態では、演算処理手段２１１は、公知の前処理付きクリロフ部分空
間法により解の算出を行なう。ここでは、後述するように粗格子を生成するための計算手順を決定し、この粗格子を用いてＶサイクルＡＭＧ法をクリロフ部分空間法の前処理に適用する。
【００３７】
流体解析装置２０の制御部２１は、連立方程式の収束解として算出された物理量についてファイル出力処理を実行する（ステップＳ１−５）。具体的には、演算処理手段２１１は、計算を行なったタイムステップ識別子に関連付けて、計算格子毎に算出した物理量に関するデータを計算結果データ記憶部２２に記録する。
【００３８】
次に、流体解析装置２０の制御部２１は、分木構造の更新処理を実行する（ステップＳ１−６）。具体的には、演算処理手段２１１は、先のタイムステップにおいて算出した物理量の勾配を用いて各計算格子の細分化を行なう。本実施例では、演算処理手段２１１は、算出した物理量に応じて異相界面領域を検出する。異相界面追跡法にはＶＯＦ（Volume
of Fluid ）法を利用する。ＶＯＦ法は、計算格子毎に流体の体積占有率であるＶＯＦ関数を定義し、ＶＯＦ関数の輸送方程式を解くことにより界面を追跡する手法である。例えば、図１０に示すように、計算格子５１０において、所定の領域の物理量の勾配が大きい場合には、計算格子５２０のように細分化を行なう。更に、計算格子５２０において、所定の領域の物理量の勾配が更に大きい場合には、計算格子５３０のように細分化を行なう。
【００３９】
演算処理手段２１１は、この界面領域を含む計算格子を８分木により細分化する。細分化された計算格子に界面領域が含まれる場合には、この計算格子を再度、８分木により細分化する。これを、予め定められた階層（レベル）まで繰り返すことにより、分木構造を更新する。
【００４０】
流体解析装置２０の制御部２１は、タイムステップをインクリメントする（ステップＳ１−７）。具体的には、演算処理手段２１１は、タイムステップ識別子に対して所定値を加算して、後続のタイムステップのシミュレーション計算処理を実行する。ここでは、予め定められたタイムステップ（タイムステップ識別子の上限値）に達するまで、ステップＳ１−３〜Ｓ１−７の処理を繰り返す。
【００４１】
（初期格子の設定処理）
次に、図２を用いて、ステップＳ１−１における初期格子の設定処理について説明する。
【００４２】
まず、流体解析装置２０の制御部２１は、モデルパラメータの取得処理を実行する（ステップＳ２−１）。具体的には、演算処理手段２１１は、解析対象領域について３次元オブジェクトのＣＡＤデータを、入力部１１を介して取得する。また、演算処理手段２１１は、タイムステップの繰り返し回数、反復法の繰り返し回数、収束判定値、流体の物性値（粘度、密度等）等に関するデータを、入力部１１を介して取得する。そして、演算処理手段２１１は、取得したデータをモデルパラメータとしてメモリ２１２に格納する。
【００４３】
次に、流体解析装置２０の制御部２１は、初期格子生成処理を実行する（ステップＳ２−２）。具体的には、演算処理手段２１１は、３次元オブジェクトが配置された３次元空間を所定の大きさの初期格子を用いて細分化する。本実施形態では、計算格子の形状として直方体の直交格子を用いる。そして、演算処理手段２１１は、生成した計算格子の特定情報（計算格子の形状情報、３次元空間における各計算格子の位置情報）をメモリ２１２に格納する。
【００４４】
次に、流体解析装置２０の制御部２１は、境界格子の細分化処理を実行する（ステップ
Ｓ２−３）。ここでは、３次元オブジェクトの形状に対応させて計算格子を細分化する。まず、演算処理手段２１１は、初期格子では３次元オブジェクトの形状を表現できない領域を特定する。例えば、直交格子を用いる場合には、曲線部分や傾斜部分、微細構造等は初期格子を用いて表現することができない。このため、演算処理手段２１１は、初期格子を細分化して、３次元オブジェクトへのフィッティングを行なう。
【００４５】
具体的には、演算処理手段２１１は、初期格子を用いて表現できない境界領域を検出して、この領域を含む計算格子を８分木に従って８個の計算格子に細分化する。これにより、第１レベルの細分格子（子格子）が生成される。
【００４６】
更に、フィッティングができていない子格子（境界格子）が存在する場合には、この子格子を系譜に基づいて特定し、この子格子を８分木に従って８個の計算格子に細分化して、第２レベルの細分格子（孫格子）を生成する。このように、演算処理手段２１１は、順次、境界領域の計算格子を細分化する。この細分化は、モデルパラメータにおいて予め設定されたレベルまで繰り返す。
【００４７】
そして、流体解析装置２０の制御部２１は、初期分木構造の記憶処理を実行する（ステップＳ２−４）。具体的には、演算処理手段２１１は、初期格子に対して細分化した計算格子について、系譜（レベル情報と配置情報）に関するデータを、初期分木構造データとしてメモリ２１２に格納する。
【００４８】
次に、流体解析装置２０の制御部２１は、初期値の設定処理を実行する（ステップＳ２−５）。具体的には、演算処理手段２１１は、入力部１１を介して、流体の物理量の初期値を取得する。ここでは、３次元空間における座標に対して、物理量（流体解析の場合には流体の速度や圧力）に関するデータを取得する。そして、演算処理手段２１１は、この物理量を初期状態データとして計算結果データ記憶部２２に記録する。演算処理手段２１１は、取得した初期値を、その座標に応じて各計算格子に付与する。
【００４９】
次に、流体解析装置２０の制御部２１は、界面格子の細分化処理を実行する（ステップＳ２−６）。具体的には、演算処理手段２１１は、各計算格子に付与した初期値に応じて、異相界面の計算格子の細分化を行なう。ここでも、ステップＳ２−３の処理と同様に、演算処理手段２１１は、界面領域を検出して、この領域を含む計算格子を８分木に従って８個の計算格子に細分化する。演算処理手段２１１は、順次、界面領域の計算格子を細分化する。この細分化は、モデルパラメータにおいて予め設定されたレベルまで繰り返す。具体的には、演算処理手段２１１は、初期分木構造に対して細分化した計算格子について、系譜（レベル情報と配置情報）に関するデータを、計算格子データとしてメモリ２１２に格納する。
【００５０】
（連立方程式求解処理）
次に、図４を用いて、連立方程式求解処理（ステップＳ１−４）を説明する。この連立方程式求解処理においては、求解のための粗格子の生成処理（図５）、ＶサイクルＡＭＧ法による計算処理（図６）を実行する。
【００５１】
まず、流体解析装置２０の制御部２１は、離散式に基づいて連立方程式の係数行列と右辺の計算処理を実行する（ステップＳ３−１）。ここでは、系の支配方程式を計算格子で離散化する。具体的には、演算処理手段２１１は、計算格子毎に、隣接する計算格子との関係に基づいて、系の支配方程式（質量保存式、運動量保存式）について離散式を生成する。この場合、メモリ２１２に格納された物理量を用いる。係数行列を「Ａ」、定数項を「ｂ」として解くべき連立方程式「Ａｘ＝ｂ」が生成される。
【００５２】
次に、流体解析装置２０の制御部２１は、ＡＭＧ法による求解のための粗格子の生成処理を実行する（ステップＳ３−２）。この粗格子の生成処理は、すべての計算格子について、分木構造レベルを粗化した粗格子を生成するための処理である。この粗格子の生成処理を、図５を用いて説明する。
【００５３】
流体解析装置２０の制御部２１は、すべての計算格子毎に粗格子の生成処理を実行し（ステップＳ４−１）、演算処理手段２１１は、粗格子の生成処理の処理対象格子を順次、特定する。
【００５４】
そして、流体解析装置２０の制御部２１は、分木構造に従って１つ上位のレベルの粗格子に統合可能かどうかを判断する（ステップＳ４−２）。具体的には、演算処理手段２１１は、系譜を用いて、この計算格子の１つ上位のレベルの格子（親格子）を特定する。そして、演算処理手段２１１は、この親格子から生成された子格子（細分格子）の状態を確認し、これらのすべての細分格子が同じレベルのリーフ格子となる場合には、分木構造に従って１つ上位のレベルの粗格子に統合できると判断する。
【００５５】
上位のレベルの粗格子に統合可能な場合（ステップＳ４−２において「ＹＥＳ」の場合）、流体解析装置２０の制御部２１は、計算格子を統合して計算手順の記録処理を実行する（ステップＳ４−３）。具体的には、演算処理手段２１１は、求解計算時において、分木構造に従って、この親格子から生成された計算格子が同じレベルである場合には、分木構造に従って１つ上位のレベルの粗格子に統合するように、その計算手順を、メモリ２１２の計算格子データに含めて記録する。
【００５６】
一方、子格子が同じレベルになっていないため、上位のレベルの粗格子に統合できない場合（ステップＳ４−２において「ＮＯ」の場合）、流体解析装置２０の制御部２１は、この計算格子に関する統合処理は行なわない。
【００５７】
そして、流体解析装置２０の制御部２１は、以上の判断処理をすべての格子について、順次、繰り返し実行する（ステップＳ４−４）。
次に、流体解析装置２０の制御部２１は、計算格子の粗化のレベルを判定する（ステップＳ４−５）。ここでは、演算処理手段２１１は、計算格子が十分に粗いかどうかを確認する。具体的には、メモリ２１２に、粗格子の生成終了の条件として、計算格子の個数（終了判定個数）に関するデータを記憶させておき、演算処理手段２１１は、計算格子が終了判定個数以下になったかどうかを判定する。
【００５８】
まだ十分に粗くなっていないと判定された場合（ステップＳ４−５において「ＮＯ」の場合）には、流体解析装置２０の制御部２１は、再度、最初からステップＳ４−１からの処理を繰り返す。一方、十分に粗くなっていると判定された場合（ステップＳ４−５において「ＹＥＳ」の場合）には、流体解析装置２０の制御部２１は、求解のための粗格子の生成処理を終了する。ここで、ＡＭＧ法によりＶサイクルを実行する場合の段数が決定される。この段数はメモリ２１２に記憶させておく。
【００５９】
そして、図４に示すように、流体解析装置２０の制御部２１は、前処理付きクリロフ部分空間法処理を実行する（ステップＳ３−３）。本実施形態では、連立方程式の解法としては、ＢｉＣＧＳＴＡＢ法を利用する。前処理としては、ＡＭＧ法によりＶサイクルを実行し、前処理付きクリロフ部分空間法により誤差が、メモリ２１２に記憶された収束判定値以下になるまで計算を繰り返す。
【００６０】
次に、ＶサイクルＡＭＧ法による計算処理を、図６を用いて説明する。
流体解析装置２０の制御部２１は、緩和法を適用して近似解の計算処理を実行する（ス
テップＳ５−１）。緩和法では、異なる粗さの計算格子を複数、用意しておき、緩和プロセスの途中で粗い計算格子に移り、誤差を早く収束させる。計算格子のサイズと同程度の波長の誤差成分が最も減衰しやすいので、波長の長い誤差を減衰させるには、大きなサイズの粗格子が有効である。
【００６１】
次に、流体解析装置２０の制御部２１は、計算格子の粗さの判断処理を実行する（ステップＳ５−２）。具体的には、演算処理手段２１１は、予め設定された段数（例えば、５段）に達したかどうかを、粗格子のレベルに基づいて判定する。
【００６２】
最も粗い格子に達していないと判定された場合（ステップＳ５−２において「ＮＯ」の場合）には、流体解析装置２０の制御部２１は、次の粗格子の取得処理を実行する（ステップＳ５−３）。具体的には、演算処理手段２１１は、図５に示した粗格子の生成処理において算出した計算手順に従って、次回の計算に用いる粗格子を取得する。ここでは、例えば、図１１に示すように、構造格子６００、構造格子６１０、構造格子６２０、構造格子６３０、構造格子６４０、構造格子６５０のように、順次、計算格子の粗化を行なう。
【００６３】
次に、流体解析装置２０の制御部２１は、近似解に対する残差を計算し、粗格子への制限補間処理を実行する（ステップＳ５−４）。係数行列を「Ａ」、定数項を「ｂ」として解くべき連立方程式を「Ａｘ＝ｂ」とおく。近似解を「ｘ’」とし、そのときの残差を「ｒ’」とした場合、「ｒ’＝ｂ−Ａｘ’」となる。演算処理手段２１１は、この残差を粗格子に制限補間する。ここで、真の解「ｘ」と近似解「ｘ’」との差を「δｘ」とした場合、修正方程式「Ａδｘ＝ｒ’」が成り立つ。そして、ステップＳ５−１に戻り、修正方程式を用いて計算処理を繰り返す。
【００６４】
一方、最も粗い格子に達したと判定された場合（ステップＳ５−２において「ＹＥＳ」の場合）には、最も粗い計算格子（ｎ）上の修正方程式の解ｘ（ｎ）が算出される。そこで、流体解析装置２０の制御部２１は、細分格子の取得処理を実行する（ステップＳ５−５）。具体的には、図５に示した粗格子の生成処理において算出した計算手順を遡って、次回の計算に用いる細分格子を取得する。
【００６５】
流体解析装置２０の制御部２１は、解を細分格子に延長補間する（ステップＳ５−６）。具体的には、演算処理手段２１１は、解ｘ（ｎ）を一つ細かい計算格子（ｎ−１）に延長補間する。
【００６６】
そして、流体解析装置２０の制御部２１は、緩和法を適用して近似解の計算処理を実行する（ステップＳ５−７）。具体的には、演算処理手段２１１は、「ｘ（ｎ−１）＋ｘ（ｎ）」を初期値として、計算格子（ｎ−１）における修正方程式を緩和する。
【００６７】
そして、流体解析装置２０の制御部２１は、最も細かい格子かどうかの判定処理を実行する（ステップＳ５−８）。最も細かい格子に達していないと判定された場合（ステップＳ５−８において「ＮＯ」の場合）には、流体解析装置２０の制御部２１は、次の細分格子の取得処理を実行する（ステップＳ５−５）。ここでは、例えば、図１１に示すように、構造格子６５０、構造格子６４０、構造格子６３０、構造格子６２０、構造格子６１０、構造格子６００のように、順次、計算格子の細分化を行なう。
【００６８】
一方、最も細かい格子に達したと判定された場合（ステップＳ５−８において「ＹＥＳ」の場合）には、流体解析装置２０の制御部２１は、ＶサイクルＡＭＧ法による計算処理を終了する。
【００６９】
以上、本実施形態によれば、以下に示す効果を得ることができる。
・  上記実施形態では、ＡＭＲ法及びＡＭＧ法を、分木構造を介してシームレスに接続しているため、詳細かつ効率的なシミュレーション計算を行なうことができる。すなわち、解析空間全体を細分化したのでは、膨大なメモリ容量、計算時間が必要になるが、ＡＭＲ法を適用することにより、界面近傍のみを細分化して、計算負荷を軽減することができる。そして、ＡＭＧ法により収束解を計算する場合、分木構造を用いて効率的に粗化を行なうことができる。
【００７０】
・  上記実施形態では、流体解析装置２０の制御部２１は、分木構造に従って１つ上位のレベルの粗格子に統合可能かどうかを判断する（ステップＳ４−２）。そして、上位のレベルの粗格子に統合可能な場合（ステップＳ４−２において「ＹＥＳ」の場合）、流体解析装置２０の制御部２１は、計算格子を統合して計算手順の記録処理を実行する（ステップＳ４−３）。流体解析装置２０の制御部２１は、生成した計算格子毎に粗格子の生成処理を実行する（ステップＳ４−１〜Ｓ４−４）。このようにして生成された粗格子を用いてＶサイクルＡＭＧ法による計算処理を実行する（ステップＳ３−３）。これにより、ＡＭＲ法の適用時の分木構造（ここでは、８分木）をそのまま用いて、ＡＭＧ法における粗化を行なうことができる。
【００７１】
特に、演算処理手段２１１は、系譜を用いて、この計算格子の１つ上位のレベルの格子（親格子）を特定し、この親格子から生成された計算格子が同じレベルである場合には、分木構造に従って１つ上位のレベルの粗格子に統合する（ステップＳ４−３）。これにより、どのレベルであっても、同じ親格子から生成された各細分格子が同じレベルである場合には、分木構造に従って１つ上位のレベルの粗格子に統合して、解析空間全体を粗化することができる。従って、全体を通じて均一に波数の長い誤差を減衰させることができる。
【００７２】
・  上記実施形態では、ＶＯＦ関数の輸送方程式を解くことにより、効率的に自由表面の形状を決定し、流体体積の保存に関しても精度がよく計算することができる。
・  上記実施形態では、ＡＭＧ法を用いるため、細分格子の係数行列から粗格子の係数行列を自動生成することができる。
【００７３】
なお、上記実施形態は、以下の態様に変更してもよい。
○  上記実施形態では、流体解析に用いたが、これに限定されるものではなく、構造解析や電磁場解析等、各種解析のシミュレーションにも適用することが可能である。
【００７４】
○  上記実施形態では、Ｖサイクルを用いたが、反復方法はこれに限定されるものではなく、Ｗサイクル、Ｆサイクルを用いることができる。
○  上記実施形態では、直方体の格子を利用して離散化したが、離散化の形状はこれに限定されるものではなく、四面体等を用いることができる。特に、ＡＭＧ法の場合には、子格子を結合して粗化を行なう場合、親格子と同じ形状である必要はなく、多様な形状を利用することができる。
【００７５】
○  上記実施形態では、分木構造として８分木を用いたが、これに限定されるものではなく、他の細分化法を利用することも可能である。例えば、２次元の場合には４分木を用いることも可能である。
【００７６】
○  上記実施形態では、分木構造に従って１つ上位のレベルの粗格子に統合可能かどうかを判断する（ステップＳ４−２）。具体的には、演算処理手段２１１は、系譜を用いて、この計算格子の１つ上位のレベルの格子（親格子）を特定する。そして、演算処理手段２１１は、この親格子から生成された子格子（細分格子）の状態を確認する。各細分格子が同じレベルのリーフ格子である場合には、分木構造に従って１つ上位のレベルの粗格子
に統合できると判断する。これに代えて、系譜に基づいて、親格子から生成された子格子となる計算格子の相互関係を記録したテーブル（隣接格子テーブル）を生成して、メモリに保持させておくことも可能である。この隣接格子テーブルを用いることにより、親格子に戻ることなく細分格子間で隣接関係を判断できるので、計算格子の粗化段階において、同じ親格子から生成されて隣接する子格子を効率的に特定することができる。
【００７７】
○  上記実施形態では、流体解析装置２０の制御部２１は、分木構造に従って１つ上位のレベルの粗格子に統合可能かどうかを判断する（ステップＳ４−２）。ここでは、粗格子の生成は、同じ親格子から生成されたすべての子格子がリーフ格子となったときに結合して粗格子を生成する。これに代えて、同じ親格子から生成された子格子であって、隣接する子格子が存在する場合には、これらの子格子を結合して粗格子を生成することも可能である。具体的には、演算処理手段２１１は、系譜を用いて、この計算格子の１つ上位のレベルの格子（親格子）を特定する。次に、演算処理手段２１１は、この親格子から生成された他の子格子の状態を確認する。そして、隣接する計算格子が同じレベルである場合には、これらを粗格子として統合（粗化）する。これにより、計算処理は複雑になるが、波数の長い誤差伝播を早くして収束を早めることができる。
【図面の簡単な説明】
【００７８】
【図１】本発明の一実施形態のシステムの概略図。
【図２】本発明の一実施形態の処理手順の説明図。
【図３】本発明の一実施形態の処理手順の説明図。
【図４】本発明の一実施形態の処理手順の説明図。
【図５】本発明の一実施形態の処理手順の説明図。
【図６】本発明の一実施形態の処理手順の説明図。
【図７】八分木構造の説明図。
【図８】計算格子の細分化の説明図。
【図９】異相界面の説明図。
【図１０】計算格子の細分化の説明図。
【図１１】計算格子の粗化の説明図。
【符号の説明】
【００７９】
１１…入力部、１２…出力部、２０…流体解析装置、２１…制御部、２１１…演算処理手段、２１２…メモリ、２２…計算結果データ記憶部、Ｃ１，Ｃ２，Ｃ３…計算格子。
【図１】
2008197921.tif 000003
【図２】
2008197921.tif 000004
【図３】
2008197921.tif 000005
【図４】
2008197921.tif 000006
【図５】
2008197921.tif 000007
【図６】
2008197921.tif 000008
【図７】
2008197921.tif 000009
【図８】
2008197921.tif 000010
【図９】
2008197921.tif 000011
【図１０】
2008197921.tif 000012
【図１１】
2008197921.tif 000013