(57)【要約】【課題】格子ボルツマン法による流体解析を行なうにあたり、格子間隔を広げても実際の流路との差が小さい解析結果を得ることができる解析を行なう。【解決手段】多孔体１に隣接する格子点上の仮想流体粒子について並進過程における粒子分布を求める際に、仮想流体粒子が多孔体１に反射される反射率または多孔体１を透過する透過率を含んだ離散速度分布関数を用いる。【選択図】図４
2010224653.tif 000002
【特許請求の範囲】
【請求項１】
計算処理空間に複数の格子点を設け、前記格子点における仮想流体粒子の挙動を求める格子ボルツマン法を用いて、多孔体が配置された流路の流体解析を行う流体解析方法であって、
前記多孔体に隣接する前記格子点上において、前記格子ボルツマン法における並進過程時の離散速度分布関数として、前記仮想流体粒子が前記多孔体に反射される反射率または前記仮想流体粒子が前記多孔体を透過する透過率を含んだ関数を用いることにより流体解析を行う、
ことを特徴とする流体解析方法。
【請求項２】
請求項１記載の流体解析方法であって、前記反射率または前記透過率を、前記多孔体の浸透係数を用いて算出することを特徴とする流体解析方法。
【請求項３】
請求項２記載の流体解析方法であって、κを前記浸透係数、ｄを前記多孔体の厚さ、δxを前記複数の格子点の間隔、Ｒを前記反射率（０＜Ｒ＜１）、１−Ｒを前記透過率、νを動粘性係数、ｃsを前記仮想流体粒子の音速とした時に、以下の数式
2010224653.tif 000018
に基づいて前記反射率または前記透過率を算出することを特徴とする流体解析方法。
【請求項４】
請求項３記載の流体解析方法であって、前記多孔体はガスフィルタのフィルタ壁であり、前記ガスフィルタ中の排出ガス流路の流体解析を行うことを特徴とする流体解析方法。
【請求項５】
計算処理空間に複数の格子点を設け、前記格子点における仮想流体粒子の挙動を求める格子ボルツマン法を用いて、多孔体が配置された流路の流体解析を行う演算部を備えた流体解析装置であって、
前記演算部は、前記多孔体に隣接する前記格子点上において、前記格子ボルツマン法における並進過程時の離散速度分布関数として、前記仮想流体粒子が前記多孔体に反射される反射率または前記仮想流体粒子が前記多孔体を透過する透過率を含んだ関数を用いて流体解析を行う、
ことを特徴とする流体解析装置。
【請求項６】
計算処理空間に複数の格子点を設け、前記格子点における仮想流体粒子の挙動を求める格子ボルツマン法を用いて、多孔体が配置された流路の流体解析を行う流体解析装置に流体解析を実行させるプログラムであって、
前記流体解析装置に、前記多孔体に隣接する前記格子点上において、前記格子ボルツマン法における並進過程時の離散速度分布関数として、前記仮想流体粒子が前記多孔体に反射される反射率と、前記仮想流体粒子が前記多孔体を透過する透過率とを含んだ関数を用いて流体解析を行わせる、
ことを特徴とするプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、格子ボルツマン法を用いた流体解析方法、装置、プログラムに関する。
【背景技術】
【０００２】
流体解析の手法の一つとして格子ボルツマン法が知られている。格子ボルツマン法を用いた流体解析においては、まず、コンピュータ等の情報処理装置において格子によって区画された計算空間を作成する。次に、格子上を動き回る仮想的な粒子である仮想流体粒子を導入し、格子点における仮想流体粒子の分布（粒子の数）の変化を求めることにより流体解析を行う。
【０００３】
仮想流体粒子の分布の変化を求める上で、仮想流体粒子の挙動について並進過程（格子点から他の格子点に移動する過程）と衝突過程（仮想流体粒子同士が衝突する過程）の２つの過程を導入し、仮想流体粒子は並進過程と衝突過程を繰り返すと仮定する。更に、時間刻み（タイムステップ）δｔを予め定め、並進過程において、時間刻みδｔごとに仮想流体粒子は必ず格子点上に存在し、格子線の途中でとどまることはできないものと仮定する。また、衝突は瞬時に行われ、かつ格子点上でのみ起こるものと仮定する。
【０００４】
さらに、並進過程と衝突過程における仮想流体粒子の分布を示す離散速度分布関数をそれぞれ定める。この離散速度分布関数を用いて、ある格子点について、並進過程における仮想流体粒子の分布と衝突過程における仮想流体粒子の分布を時間刻みδｔごとに算出する。この演算を定められた範囲の格子点（例えば全格子点）について行い、仮想流体粒子の分布の変化を算出することにより流体解析を行っている。
【０００５】
この流体解析において、フィルタ等の多孔体を含む流路について流体解析を行う場合がある。この場合、図１０に示すように計算空間上に多孔体１が配置される。従来技術においては、多孔体１と重なる格子点と多孔体１とは重ならない格子点とで衝突過程の離散速度分布関数を異ならせることにより、多孔体１を含む流路の流体解析を行っている。
【０００６】
すなわち、多孔体１の内部にある流体の挙動はブリンクマン方程式に従うことが知られている。一方、多孔体１の外部の流体の挙動はナビエ＝ストークス方程式に従うことが知られている。そこで、非特許文献１においては、多孔体１と重なっている格子点における衝突過程ではブリンクマン方程式に従うように離散速度分布関数のパラメータを定め、他方、多孔体１と重なっていない格子点についてはナビエ＝ストークス方程式に従うように離散速度分布関数のパラメータを定めることにより流体解析を行っている。
【先行技術文献】
【非特許文献】
【０００７】
【非特許文献１】林、久保「ＤＰＦフィルタにおける煤煙粒子の濾過に関するコンピュータシュミレーションの研究（Ｃｏｍｐｕｔｅｒ  ｓｉｍｕｌａｔｉｏｎ  ｓｔｕｄｙ  ｏｎ  ｆｉｌｔｒａｔｉｏｎ  ｏｆ  ｓｏｏｔ  ｐａｒｔｉｃｌｅｓ  ｉｎ  ｄｉｅｓｅｌ  ｐａｒｔｉｃｕｌａｔｅ  ｆｉｌｔｅｒ）」、コンピュータ、計算とそのアプリケーション（Ｃｏｍｐｕｔｅｒｓ  ａｎｄ  Ｍａｔｈｅｍａｔｉｃｓ  ｗｉｔｈ  Ａｐｐｌｉｃａｔｉｏｎｓ）、平成２０年、第５５巻、ｐ．１４５０−１４６０
【発明の概要】
【発明が解決しようとする課題】
【０００８】
格子点間隔を広く取った場合、多孔体１と重なる格子点が少なくなり、多孔体１の影響が正確に反映されない状態で解析が行われる可能性が生じる。例えば、図１１に示す様に、多孔体１は流路を遮るように配置されているにも関わらず、流入側（図１１では左側）から流出側（右側）に移動する間に多孔体１と重なる格子点を経由しない経路２、つまり、計算上は多孔体１に遮られない経路２が生ずる可能性がある。このような状態で解析を行うと、実際の流路との差が大きい解析結果が出力されることになる。一方、格子点間隔を狭く設定し、多孔体１と重なり合う格子点の数を多くすると、計算の負荷が大きくなるという別の問題が生ずる。
【０００９】
そこで、本発明は、格子点間隔を広げても誤差の小さい解析結果を得ることのできる流体解析を行うことを目的とする。
【課題を解決するための手段】
【００１０】
本発明は、計算処理空間に複数の格子点を設け、前記格子点における仮想流体粒子の挙動を求める格子ボルツマン法を用いて、多孔体が配置された流路の流体解析を行う流体解析方法であって、前記多孔体に隣接する前記格子点上において、前記格子ボルツマン法における並進過程時の離散速度分布関数として、前記仮想流体粒子が前記多孔体に反射される反射率または前記仮想流体粒子が前記多孔体を透過する透過率を含んだ関数を用いることにより流体解析を行うことを特徴とする。
【００１１】
また、本発明は、前記反射率または前記透過率を、前記多孔体の浸透係数を用いて算出することを特徴とする。
【００１２】
また、本発明は、κを前記浸透係数、ｄを前記多孔体の厚さ、δxを前記複数の格子点の間隔、Ｒを前記反射率（０＜Ｒ＜１）、１−Ｒを前記透過率、νを動粘性係数、ｃsを前記仮想流体粒子の音速とした時に、以下の数式
2010224653.tif 000003
に基づいて前記反射率または前記透過率を算出することを特徴とする。
【００１３】
また、本発明は、請求項３記載の流体解析方法であって、前記多孔体はガスフィルタのフィルタ壁であり、前記ガスフィルタ中の排出ガス流路の流体解析を行うことを特徴とする。
【００１４】
さらに、本発明は、計算処理空間に複数の格子点を設け、前記格子点における仮想流体粒子の挙動を求める格子ボルツマン法を用いて、多孔体が配置された流路の流体解析を行う演算部を備えた流体解析装置であって、前記演算部は、前記多孔体に隣接する前記格子点上において、前記格子ボルツマン法における並進過程時の離散速度分布関数として、前記仮想流体粒子が前記多孔体に反射される反射率または前記仮想流体粒子が前記多孔体を透過する透過率を含んだ関数を用いて流体解析を行うことを特徴とする。
【００１５】
さらに、本発明は、計算処理空間に複数の格子点を設け、前記格子点における仮想流体粒子の挙動を求める格子ボルツマン法を用いて、多孔体が配置された流路の流体解析を行う流体解析装置に流体解析を実行させるプログラムであって、前記流体解析装置に、前記多孔体に隣接する前記格子点上において、前記格子ボルツマン法における並進過程時の離散速度分布関数として、前記仮想流体粒子が前記多孔体に反射される反射率と、前記仮想流体粒子が前記多孔体を透過する透過率とを含んだ関数を用いて流体解析を行わせることを特徴とする。
【発明の効果】
【００１６】
本発明によれば、多孔体に隣接する格子点で仮想流体粒子の挙動を算出する際、仮想流体粒子が多孔体に撥ね返される反射率または多孔体を透過する透過率を含んだ離散速度分布関数を用いる。これにより、格子点間隔を拡大した状態でも誤差の少ない流体解析結果を得ることが可能となる。
【図面の簡単な説明】
【００１７】
【図１】本発明にかかる流体解析装置のシステム構成図である。
【図２】計算空間に多孔体を配置した図である。
【図３】多孔体に隣接する格子点の並進過程を示す図である。
【図４】多孔体を跨ぐ仮想流体粒子および多孔体に撥ね返される仮想流体粒子の挙動を示す図である。
【図５ａ】従来技術と本発明の流体解析方法を示す図である。
【図５ｂ】従来技術と本発明の流体解析方法を示すずである。
【図６】反射率を求めるためのグラフを示す図である。
【図７】演算部の演算処理を示すフローチャートである。
【図８】本発明にかかる流体解析装置の機能ブロック図である。
【図９ａ】本発明による流体解析結果を示す図である。
【図９ｂ】本発明による流体解析結果を示す図である。
【図１０】従来技術における流体解析方法を示す図である。
【図１１】従来技術における流体解析方法を示す図である。
【発明を実施するための形態】
【００１８】
本発明の実施の形態における流体解析装置の構成を図１に示す。流体解析装置は、図１に示すように、演算部３、入力部４、表示部５及び記憶部６を含んで構成される。
【００１９】
演算部３は、いわゆるＣＰＵであり、記憶部６に予め記憶させてある流体解析プログラムを実行することによって後述する流体解析処理を行う。入力部４は、キーボード等の入力手段を含んで構成され、流体解析処理の対象となる多孔体１の大きさ、配置及び特性、計算空間の格子点の設定、全格子点の仮想流体粒子分布（各仮想流体粒子の離散速度分布関数ｆｉ（ｘ、０））、境界条件、演算時間の刻み、終了時刻等の入力をユーザから受け付ける。入力された情報は記憶部６に格納及び保持される。表示部５は、ディスプレイやプリンタ等の情報出力手段を含んで構成され、流体解析処理の対象となる計算空間、設定されたパラメータ、処理経過及び処理結果等の情報を表示出力する。記憶部６は、半導体メモリやハードディスク等の記憶手段を含んで構成され、流体解析プログラム、流体解析処理のために設定された情報、解析処理での演算結果等を格納及び保持する。記憶部６は、演算部３から適宜アクセス可能である。
【００２０】
次に、演算部３が行う流体解析処理について説明する。本発明では従来の格子ボルツマン法における多孔体１での流体粒子の挙動を求める演算過程の一部を新たな演算過程に代えて流体解析を行う。
【００２１】
格子ボルツマン法の説明にあたり、まず衝突過程から説明する。ここで、従来の衝突過程では、多孔体１と重なる格子点における流体の解析にはブリンクマン方程式を適用し、多孔体１と重ならない格子点における流体の解析にはナビエ＝ストークス方程式を適用したが、後述するように本発明においては多孔体１に格子点を重ねなくても流体解析を行うことができ、ナビエ＝ストークス方程式を適用するだけで解析が可能である。そのため、本実施形態においては図２に示すように格子点間隔を多孔体１の幅よりも大きく取って多孔体１と重なる格子点が無い様に多孔体１を配置している。以下に、多孔体１とは重ならない格子点に対してナビエ＝ストークス方程式が適用される衝突過程について説明する。
【００２２】
時刻ｔにおいて、位置ベクトルｘで示される格子点上に存在する仮想流体粒子のうち、格子点間の移動を表すベクトルである離散速度ベクトルがｅiである仮想流体粒子の分布を表す離散速度分布関数ｇｉ（ｘ，ｔ）は、下記数式１によって表される。
【００２３】
【数１】
2010224653.tif 000004
【００２４】
ここで、右辺のｆｉ（ｘ、ｔ）は衝突過程の前の並進過程における仮想流体粒子の分布を示している。また、離散速度ベクトルｅiの接尾辞（サフィックス）ｉは、三次元の流体解析における粒子が移動可能な方向を示す。例えば、各格子点に静止及び１４方向へ移動可能な場合には、下記数式２で示されるように、０から１４までの値を採る。
【００２５】
【数２】
2010224653.tif 000005
【００２６】
また、τは単一時間緩和係数と呼ばれる無次元の係数であり、τδｔ[ｓｅｃ]は仮想流体粒子が平衡状態に達するまでの時間を表す。単一時間緩和係数τは０．５＜τ＜∞の範囲で任意の値を設定することができる。また、ｆｉ（ｅｑ）（ｘ、ｔ）はある空間領域において平衡状態に達したときの仮想流体粒子の分布関数である。つまり、数式１は仮想流体粒子が衝突するたびに仮想流体粒子の分布が平衡状態に向かうことを表している。
【００２７】
ここで、ｆｉ（ｅｑ）（ｘ、ｔ）について説明する。前述したとおり、多孔体１と重なっていない格子点においては、仮想流体粒子の挙動はナビエ＝ストークス方程式に従う。ナビエ＝ストークス方程式に従うように種々のパラメータを定めると、ｆｉ（ｅｑ）（ｘ、ｔ）は下記数式３のようになる。
【００２８】
【数３】
2010224653.tif 000006
ここで、ｗｉは重み係数を表し、接尾辞ｉによって下記数式４のような値をとる。
【００２９】
【数４】
2010224653.tif 000007
【００３０】
また、ｐは圧力であり、ｕは流速ベクトルを表す。ここで、数式３について、重み係数ｗi、最小速さｃ、および離散速度ベクトルｅiは予め与えられているから、ｆi（ｅｑ）（ｘ、ｔ）は圧力ｐと流速ベクトルｕを求めることにより導出される。圧力ｐと流速ベクトルｕはそれぞれ数式５、６のように求めることができる。
【００３１】
【数５】
2010224653.tif 000008
【００３２】
【数６】
2010224653.tif 000009
【００３３】
つまり、圧力ｐおよび流速ベクトルｕは、衝突過程の前の並進過程における仮想流体粒子の分布から導くことができる。そこで、並進過程における離散速度分布関数について、以下に説明する。
【００３４】
本発明においては、並進過程において、多孔体１に隣接した格子点と多孔体１に隣接していない格子点とで離散速度分布関数の適用条件を変えている。これにより、仮想流体粒子が多孔体１に撥ね返され、または多孔体１を透過する過程を演算に反映させている。ここで、多孔体１に隣接する格子点とは、多孔体１と重なる格子線の両端にある格子点を指す。
【００３５】
まず、多孔体１に隣接していない格子点について、並進過程における離散速度分布関数ｆｉ（ｘ，ｔ）について以下の数式７が適用される。
【００３６】
【数７】
2010224653.tif 000010
【００３７】
ここで、右辺は並進過程の前の衝突過程における仮想流体粒子の分布を示している。
【００３８】
一方、多孔体１に隣接する格子点の並進過程における離散速度分布関数ｆｉ（ｘ，ｔ）は、離散速度ベクトルｅiに応じて場合分けされる。この場合分けを図３を使って説明する。格子点ｘは多孔体１に隣接している。この場合において、ｘとｘ−ｅ１との間、ｘとｘ−ｅ２との間、ｘとｘ−ｅ３との間の様に、格子点ｘ−ｅｉとｘとの間に多孔体１が存在しないときには前述の数式７が適用される。一方、格子点ｘとｘ−ｅ４との間のように、格子点ｘ−ｅｉとｘとの間に多孔体１が存在するとき、下記の数式８が適用される。
【００３９】
【数８】
2010224653.tif 000011
【００４０】
ここで、右辺のｇｉ（ｘ，ｔ）およびｇｉ（ｘ−ｅｉ，ｔ）はともに並進過程の前の衝突過程における仮想流体粒子の分布を示している。また、ｆ、ｇの接尾辞のうち、−ｉはｉとは反対方向に進む仮想流体粒子を表す。また、Ｒは反射率（０＜Ｒ＜１）であり、仮想流体粒子が多孔体１に撥ね返される比率を示している。一方、（１−Ｒ）は透過率であり、仮流流体粒子が多孔体１を通過する比率を示している。
【００４１】
ここで、数式８について、図４を用いて説明する。数式８の１番目の式において、右辺第１項は格子点ｘにいた仮想流体粒子のうち、隣の格子点に移動する際に多孔体１に撥ね返されて元の格子点ｘに戻った仮想流体粒子を表し、右辺第２項は多孔体１を挟んで格子点ｘの隣の格子点にいた仮想流体粒子のうち、多孔体１を通り抜けて格子点ｘに移動した仮想流体粒子を表している。
【００４２】
このように、本発明では多孔体１に隣接する格子点上の並進過程について、反射率または透過率を含んだ離散速度分布関数を用いることにより、多孔体１を含む流路の流体解析を正確に行うことができる。前述した図１１を例に取ると、反射率または透過率を含んだ離散速度分布関数が適用される格子点を経由せずに流入側から流出側に移動できる経路は採り得ない。このように格子点間隔を広く取っても、流れを遮るように多孔体１が配置された流路の構造を正確に反映した状態で解析を行うことができる。
【００４３】
さらに、本発明では多孔体１に隣接する格子点に対して反射率または透過率を含んだ離散速度分布関数が適用されることから、仮想流体粒子が多孔体１に反射されまたは多孔体１を透過する過程を演算に含めるために格子点を多孔体１に重ねなくても流体解析を行うことができる。したがって、前述した図２のように多孔体１に重なる格子点が無い様に多孔体１を配置し、ナビエ＝ストークス方程式のみを用いて流体解析を行うことができる。従来技術においてはナビエ＝ストークス方程式に加えて多孔体１に重なる格子点についてブリンクマン方程式を用いていたが、本発明ではブリンクマン方程式を適用しなくても流体解析が行えるから従来技術に比べて解析手法が簡潔になる。なお、多孔体１の形状や多孔体１の配置の都合上、多孔体１と重なる格子点が存在する場合には、その格子点についてはブリンクマン方程式が適用される。
【００４４】
次に、数式８における反射率Ｒの導出について説明する。本発明者らは下記の数式９にもとづいて反射率Ｒを導出できることを見出した。このことにつき、以下に説明する。
【００４５】
【数９】
2010224653.tif 000012
【００４６】
ここで、κは多孔体１の浸透係数、ｄは多孔体１の厚さ、νは動粘性係数（粘性係数μを密度ρで割ったもの）を表している。またｃsは仮想流体粒子の音速を示し、非圧縮性流体においては定数として扱われ、例えば、２次元の正方格子で解析を行う場合、
2010224653.tif 000013
である。
【００４７】
反射率Ｒの導出にあたり、図５ａのように、流体の流れを遮るように無限長の多孔体１を直線状に格子点上に配置した流路を考える。次に、同じ流路について図５ｂのように多孔体１の位置を微少量ずらして多孔体１と重なる格子点を無くした流路を考える。
【００４８】
図５ａの格子点における仮想流体粒子の離散分布速度関数ｆｉ（ｘ、ｔ）およびｇｉ（ｘ，ｔ）は従来技術であるブリンクマン方程式により導くことができる。ここで、図５ａと図５ｂは同じ流路を表しているから両図の同一地点における仮想流体粒子の分布は同一となることが理解される。そこで、ブリンクマン方程式を用いて図５ａの多孔体１内の格子点７における仮想流体粒子の分布を算出し、算出された仮想流体粒子の分布は図５ｂにおいて図５ａの格子点７に対応する格子点８における仮想流体粒子の分布とほぼ等しいと仮定し、算出された仮想流体粒子分布を数式８に代入することにより格子点８における反射率Ｒを算出する。
【００４９】
一方、ブリンクマン方程式により求めた仮想流体粒子分布に基づいて、数式５により圧力ｐを求め、また数式６により速度ｕを算出する。このｐと下記数式１０にもとづいて△ｐ〜を算出し、またｕと下記数式１１に基づいてｕ〜を算出して△ｐ〜／ｕ〜の値を求める。
【００５０】
【数１０】
2010224653.tif 000014
【００５１】
【数１１】
2010224653.tif 000015
【００５２】
この結果、仮想流体粒子分布がある値を取ったときの反射率Ｒと△ｐ〜／ｕ〜との組（Ｒ，△ｐ〜／ｕ〜）ができる。仮想流体粒子分布を種々変更して（Ｒ，△ｐ〜／ｕ〜）の組み合わせを複数得る。次に、図６に示すように△ｐ〜／ｕ〜とＲとの関係をプロットしたグラフを作成する。このプロットに沿う関数をフィッティング等により求めると、下記の数式１２を得る。
【００５３】
【数１２】
2010224653.tif 000016
【００５４】
ここで、ブリンクマン方程式による流体解析においては、流体の挙動は下記数式１３に示すダルシー則に従うことが知られている。
【００５５】
【数１３】
2010224653.tif 000017
【００５６】
数式１２と数式１３から、数式９が得られる。数式９においてδｘ、ｃs、ｄは予め定められており、κも実験により求めることができる。さらにνは密度ρによって変化するが、非圧縮性流体においてはρは圧力ｐのみを変数とする関数により導出可能であり、数式５でｐを算出することでνを導くことができる。これらの値を数式９に代入することにより反射率Ｒを導出することができる。Ｒが求められたことにより、数式８において並進過程における仮想流体粒子の分布を算出することができる。
【００５７】
以上のようにして、衝突過程における仮想流体粒子の分布と、並進過程における仮想流体粒子の分布を繰り返し算出し、各格子点における仮想流体粒子の分布の変化を求めることにより、流路の流体解析が行われる。
【００５８】
以上、本発明に係る流体解析の手法について説明した。次に、流体解析処理の具体的な手順につき、図７のフローチャートを参照して説明する。演算部３は、記憶部６に記憶されている流体解析処理プログラムを実行することによって以下に示す処理を実行する。これによって、後述するように流体解析装置９は図８に示す各手段として機能する。なお、以下の処理は、演算時間刻み毎に計算空間に設定された全格子点に位置する全仮想流体粒子に対して行う。計算空間に存在する総ての粒子についての処理が終了すると次の演算時間刻みの演算へ移行する。
【００５９】
まず演算部３は記憶部６に記憶された多孔体１の形状、初期条件、境界条件等を取得する（Ｓ１）。ここで、初期条件は、ｔ＝０における全格子点の仮想流体粒子分布（各仮想流体粒子の離散速度分布関数ｆｉ（ｘ、０））を含む。また、境界条件は、流路の壁において仮想流体粒子は壁を突き抜ける方向には移動しない等の条件を含む。
【００６０】
次に演算部３は前時刻ｔにおける仮想流体粒子分布ｆｉ（ｘ，ｔ）を数式５に代入して圧力ｐを得る。また、ｆｉ（ｘ，ｔ）を数式６に代入して速度ｕを得る（Ｓ２）。ここで得た圧力ｐと速度ｕを数式３に代入して平衡状態における仮想流体粒子分布ｆｉ（ｅｑ）（ｘ，ｔ）を得る（Ｓ３）。更にｆｉ（ｘ，ｔ）とｆｉ（ｅｑ）（ｘ，ｔ）とを数式１に代入することにより、現時刻ｔにおける衝突過程の仮想流体粒子分布ｇｉ（ｘ，ｔ）を得る（Ｓ４）。
【００６１】
演算部３は衝突過程の演算に続き並進過程の演算を行う。まず演算部３は演算の対象となる格子点ｘと、隣接する格子点ｘ−ｅｉとの間に多孔体１が存在するか否かを確認する（Ｓ５）。多孔体１が存在しなければｇｉ（ｘ−ｅｉ，ｔ）を数式７に代入することにより並進過程の仮想流体粒子分布ｆｉ（ｘ，ｔ＋δｔ）を得る（Ｓ６）。
【００６２】
一方、演算の対象となる格子点ｘと隣接する格子点ｘ−ｅｉとの間に多孔体１が存在する場合は演算部３は数式９に基づいて反射率を算出する（Ｓ７）。この反射率とＳ４で算出したｇ−ｉ（ｘ，ｔ）、ｇｉ（ｘ−ｅｉ，ｔ）を数式８に代入して並進過程の仮想流体粒子分布ｆｉ（ｘ，ｔ＋δｔ）を得る（Ｓ８）。
【００６３】
以上の衝突過程と並進過程の仮想流体粒子の分布の算出を、記憶部６に記憶された終了時刻まで行う。もしまだ終了時刻に達していない場合はＳ２に戻って演算を繰り返す（Ｓ９）。終了時刻に達したときは演算を終了し、表示部５に解析結果を出力する（Ｓ１０）
【００６４】
以上説明した流体解析処理を実行する流体解析装置９の機能ブロック図を図８に示す。流体解析装置９は、Ｓ１を実行する入力データ取得手段１０、Ｓ２からＳ４までを実行する衝突過程演算手段１１、Ｓ５からＳ８までを実行する並進過程演算手段１２、Ｓ９を実行する終了時刻確認手段１３、Ｓ１０を実行する解析結果出力手段１４を含んで構成される。
【００６５】
本発明による流体解析の結果を図９ｂに示す。図９ａは入力部４から入力した多孔体１を含む流路、図９ｂはその流路について本発明によるシミュレーションを行った結果である。このように、本発明に係る流体解析においては流路の構造が正確に反映された状態で解析が行われることが理解される。
【００６６】
以上、演算部本発明にかかる流体解析について説明した。この流体解析は種々の多孔体や流体について行うことが可能である。たとえば、多孔体である排出ガスフィルタを含む流路において排出ガスの挙動を解析したり、当該ガスフィルタ内の排出ガスの挙動を解析することができる。
【符号の説明】
【００６７】
１  多孔体、２  経路、３  演算部、４  入力部、５  表示部、６  記憶部、９  流体解析装置、１０  入力データ取得手段、１１  衝突過程演算手段、１２  並進過程演算手段、１３  終了時刻確認手段、１４  解析結果出力手段。
【図１】
2010224653.tif 000019
【図２】
2010224653.tif 000020
【図３】
2010224653.tif 000021
【図４】
2010224653.tif 000022
【図５ａ】
2010224653.tif 000023
【図５ｂ】
2010224653.tif 000024
【図６】
2010224653.tif 000025
【図７】
2010224653.tif 000026
【図８】
2010224653.tif 000027
【図９ａ】
2010224653.tif 000028
【図９ｂ】
2010224653.tif 000029
【図１０】
2010224653.tif 000030
【図１１】
2010224653.tif 000031