JP2021144626A - 流体シミュレーション方法および流体シミュレーションプログラム - Google Patents

Abstract

【課題】粒子法により、流体の一部が主領域から特定の流出量で流出するシミュレーションを行う際に、課された流出境界条件通りに粒子法シミュレーションを実施する方法を提案する。【解決手段】あらかじめ設定された主領域から流体が流出する系で、粒子法を用いて上記流体の挙動を計算するに際し、課された流出境界条件通りの粒子法シミュレーションを実現させる方法であって、上記主領域に助走領域を接続し、上記助走領域内の粒子の速度を、上記主領域への前記流体の流入境界条件と流出境界条件からなる境界条件および初期条件とから求められる上記主領域内の理論粒子数と、計算途上における上記主領域内の全粒子数との差に応じて変更する。【選択図】なし

Description

本発明は、流体シミュレーション方法および流体シミュレーションプログラムに関する。

流体の流動挙動のシミュレーションは、主に、空間を格子状にメッシュ分割して計算を行う格子法ＣＦＤ（Computational Fluid Dynamics）計算ソフトウェアを用いて行われている。現在、格子法ＣＦＤに基づく数多くの市販ソフトウェアが存在し、広く用いられているが、格子法ＣＦＤは、自由界面の取扱いが不得意とされている。

これに対して、粒子法は、流体を多数の粒子で表現し、各粒子の運動方程式を解いて流動挙動をシミュレーションする方法であり、空間を格子状にメッシュ分割する必要がなく、自由界面の取扱いが得意な手法とされている（例えば、非特許文献１参照）。

従来、粒子法は、シミュレーションを行う際の計算機負荷が大きいことが実用上の難点とされてきたが、昨今の計算機ハードウェアの発達により実用化の段階に入り、幾つかの市販ソフトウェアも利用可能となっている。

粒子法は、今後、自由界面を含んだ流体系のシミュレーション計算に広く用いられることが期待されているが、その計算手法の特徴から、流体流動を所与の流入条件および流出条件、とりわけ与えられた流出条件通りの計算を実現することが原理的に難しいという問題点がある。

すなわち、粒子法による流体の流動挙動のシミュレーションで、系全体についてシミュレーションすると計算機の負荷が高くなりすぎる場合には、系内の一部の領域を計算対象領域として設定し、計算対象領域内の流体についてのみシミュレーションするのが一般的である。

また、流体のシミュレーションでは、計算対象領域への流体の流入や計算対象領域からの流体の流出が発生しうる。粒子法では、計算対象領域内に特定の流速（流入量）で流体を流入させることは容易に行うことができる。しかし、例えば、計算対象領域に流体を流入させる流入境界の圧力値が一定となるような流入境界条件を課してシミュレーションを行うことは、困難である。

上述のような流入境界条件を課したシミュレーションの試みとして、特許文献１には、流体の流入口の近傍に粒子生成面（流入境界）を形成し、粒子生成面に応じて粒子を削除する粒子消失面を設定してシミュレーションを行い、流体の流入方向に逆流する粒子がこの粒子消失面を超えた場合に粒子を削除することにより、流入口での圧力上昇を回避する方法について記載されている。

特開２０１８−１８１６４号公報

田中良和他、「自由水面が大変形する局所急変流解析のための粒子法流体解析の統合環境の開発」、農業土木学会論文集、No.246、pp95(2006)

一方、流出境界条件についても、特定の流出量で流体を流出させるような流出境界条件を課したシミュレーションを行うことは、粒子法を用いたシミュレーションでは困難である。これは、粒子法では、計算対象領域内の粒子は、シミュレーションの結果として流出境界に到達するため、流出境界に条件として課された流出量に対応する数の粒子が到達するとは限らないためである。その結果、たとえば、貯水タンクに一定量の水が流入し、流入した水と同量の水が同時に排出される系で、貯水タンク内の水面の挙動をシミュレーションで解明しようとすると、ほとんど動かないはずの水面が大きく変動する等、妥当なシミュレーション結果が得られない現象が発生する。

このように、粒子法では、計算対象領域から特定の流出量で流体が流出するようなシミュレーションを行うことは困難であり、こうしたシミュレーションを可能とするような課された流出境界条件通りに粒子法シミュレーションを実施する良い方法はこれまで提案されていない。

そこで、本発明の目的は、粒子法により、流体の一部が主領域から特定の流出量で流出するシミュレーションを行う際に、課された流出境界条件通りに粒子法シミュレーションを実施する方法を提案することにある。

上記課題を解決する本発明は、以下の通りである。
［１］あらかじめ設定された主領域から流体が流出する系で、粒子法を用いて前記流体の挙動を計算するに際し、課された流出境界条件通りの粒子法シミュレーションを実現させる方法であって、前記主領域に助走領域を接続し、前記助走領域内の粒子の速度を、前記主領域への前記流体の流入境界条件と流出境界条件からなる境界条件および初期条件とから求められる前記主領域内の理論粒子数と、計算途上における前記主領域内の全粒子数との差に応じて変更する、流体シミュレーション方法。
［２］あらかじめ設定された主領域から流体が流出する系で、粒子法を用いて前記流体の挙動を計算するプログラムであって、コンピュータに、前記主領域に接続する助走領域を設定する処理、および当該助走領域内の粒子の速度を、前記主領域内の理論粒子数と、前記計算途上における主領域内の全粒子数の差に応じて変更する処理を実行させるための流体シミュレーションプログラム。

本発明によれば、粒子法を用いて、流体の一部が計算対象領域から特定の流出量で流出するシミュレーションを行う際に、課された流出境界条件通りに粒子法シミュレーションを実施することができる。

実施例に使用した流体シミュレーションのモデルを示す図である。 助走領域の設定および助走領域に流入した粒子の速度制御を説明する図である。

（流体シミュレーション方法）
以下、図面を参照して本発明の実施形態について説明する。本発明による流体シミュレーション方法は、あらかじめ設定された主領域から流体が流出する系で、粒子法を用いて、水やマッハ数０．３以下の空気といった非圧縮性流体の一部が上記主領域から特定の流出量で流出する流動挙動をシミュレーションするに際し、主領域の流出境界からの非圧縮性流体の流出量を流出境界条件として課してシミュレーションを実施する方法である。ここで、上記主領域に助走領域を接続し、当該助走領域内の粒子の速度を、上記主領域への上記流体の流入境界条件と流出境界条件とからなる境界条件および初期条件から求められる理論粒子数と、計算途上における上記主領域の全粒子数との差に応じて変更することを特徴とする。

上述のように、粒子法では、計算対象領域内の粒子は、シミュレーションの結果として流出境界に到達するものであるため、流体の一部が主領域から特定の流出量で流出するような流出境界条件を課してシミュレーションを行うことは困難である。

本発明者らは、上述のような課された流出境界条件通りのシミュレーションを粒子法で実現する方法を鋭意検討した。その過程で、主領域の流出境界に到達した粒子を単純に消滅させる方法や、主領域内の流出境界に隣接して助走領域を設定し、助走領域に流入した粒子の速度を課された流出境界条件の流出量に対応する速度に固定する方法などを試みた。しかしながら、これらの方法では、後述する実施例に示すように、適切な計算結果を得ることが出来なかった。

しかし、主領域の流出境界に助走領域を接続し、助走領域内の粒子の速度を、適宜変えることにより、適切な計算結果を得ることができることを見出し、本発明を完成させたのである。以下、本発明における助走領域の設定方法および粒子の速度の変更方法について説明する。

まず、例として、図１に示す、供給管２から貯蔵部１に一定の流入量で非圧縮性流体（以下、特に断りがない限り、非圧縮性流体を「流体」と称する。）が供給され、排出管３から流入量と同じ流出量で流体が排出されるような流体の流動挙動をシミュレーションする場合について考える。この場合、現実の貯蔵では、流体の流入量と流出量とが同じであるため、貯蔵部１内の流体の上面の位置は変化しない。

ここで、図１の系で、粒子法によって流体の挙動を計算してみる。その時、貯蔵部１内のみを主領域として計算対象とすれば、計算機の負荷の点で有利である。そこで、供給管２と貯蔵部１との接続境界４を流入境界として設定するとともに、排出管３と貯蔵部１との接続境界５を流出境界として設定した。そして、非圧縮性流体を水として、流入境界から特定の流入量で流体（粒子）を流入させ、流出境界に到達した粒子を単純に消滅させる条件でシミュレーションを行ったところ、後述する実施例に示すように、水位が上昇してしまうという結果を得た。
つまり、流体の流入量と流出量とが同じであるにもかかわらず、単に粒子法を適用して計算しただけでは、現実の状態は再現できない。

そこで、本発明では、例えば、図２に示すように、主領域である貯蔵部１に接続された排出管３のもう一方の端部を流出境界６として設定し、排出管３の接続境界５と流出境界６との間の領域を助走領域７として設定する。助走領域７にも粒子が充填されており、その粒子は速度Ｖで移動する。そして、主領域（すなわち、貯蔵部１）と助走領域７とを合わせて計算対象領域とし、計算対象領域全体を粒子法で計算する。

なお、助走領域７は、課された流出境界通りのシミュレーションを実現するための工夫として設ける仮想の領域であり、実際に存在する排出管３とは異なるものである。例えば、助走領域７の長さは、後述するように計算条件に応じて適宜調整されるものであり、実際の排出管３の長さとは異なるが、助走領域７の内径は排出管３の内径と同一に設定されるため、説明の便宜上、排出管３として説明している。この時、助走領域７内の粒子の速度Ｖを、適宜変更する。

上記助走領域７内の粒子の速度Ｖは、例えば、ある時刻の、粒子数差に応じて計算することにより実現できる。ここで、粒子数差とは、単位時間当たりの流体の主領域への流入条件と主領域からの流出条件および初期条件とから求められる理論粒子数と、計算によって得られた主領域内の全粒子数との差であり、ΔＮであらわす。なお、流出境界６から流出する粒子数と、助走領域７から流出する粒子数とは、同義である。よって、以降、両者を区別しないことがある。

上記の、貯蔵部１へ流体が流入し、同時に貯蔵部１から流体が流出する時の流体の挙動を調べる場合、貯蔵部１の容量が決まり、初期条件である貯蔵部１に最初に蓄えられていた流体の量、境界条件である流体の流入条件（単位時間あたりの流入量）と流出条件（単位時間あたりの流出量）が設定されれば、任意の時刻ｔにおける貯蔵部１内の流体の量が定まる。ここで流体の体積と粒子数との関係を決める（例えば７リットルの流体を３０万個の粒子で表現する）と、任意の時刻ｔにおける主領域内に存在する理論粒子数を求めることができる。流体の流入量と流出量とが同じ場合には、理論粒子数は不変であり、流入量が流出量よりも多ければ、その差に応じて、時間とともに理論粒子数は増加していき、逆に流入量が流出量よりも少なければ、その差に応じて、時間とともに理論粒子数は減少する。ここで、理論粒子数をＮ₀とし、時刻ｔでの理論粒子数をＮ₀（ｔ）とする。

一方、計算によって得られた時刻ｔにおける主領域内の全粒子数をＮ（ｔ）とする。すなわち、時刻ｔでの粒子数差は、以下の（１）式で与えられる。
ΔＮ（ｔ）＝Ｎ（ｔ）− Ｎ₀（ｔ） （１）

次に、助走領域７内の、時刻ｔでの粒子の速度Ｖは、すべて同じであり、すべての粒子が同じ方向に移動するとする。そして、時刻ｔからΔｔだけ進んだ時刻（ｔ＋Δｔ）での助走領域７内の粒子の速度Ｖ（ｔ＋Δｔ）を、上記式（１）を用いて、下記式（２）のように設定する。
Ｖ（ｔ＋Δｔ）＝Ｖ（ｔ）×｛１＋αΔＮ（ｔ）／Ｎ₀｝ （２）
ここで、αは定数であり、計算対象に応じて適切に設定される。

このように、本発明では、接続境界５で主領域に助走領域７を接続し、助走領域７内の粒子の速度Ｖを、粒子数差ΔＮに応じて変更する。

時刻ｔからΔｔの間に、助走領域７内の粒子は、距離Ｌ＝Ｖ（ｔ＋Δｔ）×Δｔを一斉に移動する。そして、流出境界６に到達した粒子は、消滅する。その反面、接続境界５付近では、助走領域７内の粒子が移動しており、その空いた空間に、粒子法による計算の結果に従って、主領域から粒子が補充される。

上記助走領域７のサイズは、流体の性質、流体の排出量、流出境界の大きさ等により、適宜設定してよい。ただし、助走領域７が長すぎると、計算の負荷が大きくなりすぎ、計算時間が必要以上に長くなる。よって、計算の負荷と、得られる計算結果のバランスを見ながら設定する。これにより、特定の流出量で流出する流体の流動挙動をシミュレーションすることができる。

また、上記助走領域７の形状は、特に限定されるものではないが、単純な円筒または中空の角材であることが好ましい。というのも、粒子法では、境界も球で設定されるため、複雑な形状では、計算の負荷が大きくなってしまうからである。

以下、本発明による流体シミュレーション方法の具体的な手順を説明するが、これに限定されない。なお、下記の各ステップは、シミュレーション装置として構成するコンピュータにより行うことができる。シミュレーション装置は、シミュレーション装置として構成されたコンピュータにインストールされたプログラムが、コンピュータのＣＰＵに実行させる処理により実現される。

なお、演算速度を向上させるためにＧＰＵを搭載したコンピュータも存在するが、上記コンピュータのＣＰＵに実行させる処理の意味する中にＧＰＵを搭載したコンピュータによる処理が含まれていることは言うまでもない。

まず、ステップＳ１において、粒子法による流体シミュレーションを行うための入力データを取得する。このステップＳ１においては、例えば図１に示した貯蔵部１のような、流体が流動する領域の形状やサイズ、流体の物性値、初期条件として貯蔵部１にあらかじめ蓄えられている流体の量とその流体の流速、境界条件として流体の流入条件と流出条件を取得し（ステップＳ１１）、計算モデルに対して主領域を決定して流体としての粒子が主領域に流入する流入境界および流出境界を設定し（ステップＳ１２）、主領域の流出境界に対して流出方向に接して存在する助走領域７を設定し（ステップＳ１３）、主領域と助走領域をあわせて計算対象領域として、計算対象領域内の粒子数を設定し、各粒子の座標を設定する（ステップＳ１４）。

ここで、主領域は、流体の大部分が存在する領域であり、たとえば、貯水タンク等が該当する。また、流体の流入条件は、たとえば主領域へ流入境界から流入する流体の単位時間当たりの流入量であり、流体の流出条件は、例えば主領域から流出境界へ流出する単位時間当たりの流体の流出量である。

次に、ステップＳ２において、ステップＳ１にて取得した入力データに基づいて、粒子法によるシミュレーションを開始する。このステップＳ２においては、主領域の流入境界から粒子を流入させるか否かの判定を行い（ステップＳ２１）、粒子を流入させると判定した場合に所定の速度で粒子を流入させる（ステップＳ２２）。流入した粒子の速度は、シミュレーションにおける流体の流入量に対応する速度とすることができる。

続いて、本発明の助走領域７を含む計算モデル設定方法を採用した流体解析手法に基づいて、粒子の運動方程式（ナビエストークス方程式）を解いて粒子の速度を求める（ステップＳ２３）。ここで用いる流体解析手法としては、Smoothed Particle hydrodynamics（ＳＰＨ）法やMoving Particle Semi-implicit（ＭＰＳ）法などが挙げられ、これらの適切な手法を使用することができる。

次いで、ステップＳ２３で更新した速度に基づいて、各粒子の位置を１ステップ分の時間（Δｔ）だけ更新する（ステップＳ２４）。

続いて、助走領域７の終端の流出境界（図２の流出境界６）を粒子が通過したか否かの判定を行い、通過したと判定した場合には、粒子を削除する（ステップＳ２５）。

そして、１ステップ分の計算の終了し、計算結果を出力する（ステップＳ２６）。

続いて、予め定めたステップ分の処理を終了したか否かを判定し（ステップＳ２７）、終了していない場合には、上記ステップＳ２１〜ステップＳ２６までの処理を繰り返し、終了している場合には、シミュレーションを終了する。

（流体シミュレーションプログラム）
続いて、本発明に係る粒子法によるシミュレーションプログラムについて説明する。本発明に係るプログラムは、あらかじめ設定された主領域から流体が流出する系で、粒子法を用いて上記流体の挙動を計算するプログラムである。ここで、コンピュータに、上記主領域に接続する助走領域を設定する処理、および当該助走領域内の粒子の速度を、上記主領域内の理論粒子数と、上記計算途上における主領域内の全粒子数の差に応じて変更する処理を実行させることを特徴とする。

上述した本発明による方法で流体シミュレーションを行う装置を、コンピュータで構成することができる。コンピュータに、上記装置として機能させるための本発明による流体シミュレーションプログラムは、各コンピュータの内部または外部に備えられる記憶部に記憶することができる。そのような記憶部は、外付けハードディスクなどの外部記憶装置、或いはＲＯＭ又はＲＡＭなどの内部記憶装置で実現することができる。

各コンピュータに備えられる制御部は、中央演算処理装置（ＣＰＵ）などの制御で実現することができる。すなわち、ＣＰＵが、上記シミュレーション装置としての機能を実現するための処理内容が記述されたプログラムを、適宜、記憶部から読み込むことによって、流体シミュレーション装置としての機能をコンピュータ上で実現させることができる。

以下、本発明の実施例について説明するが、本発明は実施例に限定されない。

（比較例１）
図１に示した、貯蔵部１と、それに接続された供給管２、排出管３からなる系で、非圧縮性流体を水として、粒子法によりシミュレーションを行った。その際、貯蔵部１は貯水タンクとした。
貯水タンク１の大きさを、内部の直径が２０ｃｍ、高さが５０ｃｍの円筒とした。供給管２および排出管３の直径はともに１８ｍｍとし、それぞれ底面から５ｃｍの高さ位置に接続した。また、貯水タンク１に供給管２から一定の流入量で水を供給し、排出管から流入量と同じ流出量で水を排出する、とした。そして、供給管２と貯水タンク１との接続境界４を流入境界とし、貯水タンク１と排出管３との接続境界５を流出境界として設定した。また、貯水タンク１内を主領域とし、この主領域を計算対象領域として、貯水タンク１内の水の流れおよび水面の変動を粒子法で計算した。

また、初期条件として貯水タンク１における初期貯水容量は７リットルとし、境界条件として水の流入量および流出量は０．１５リットル／秒とした。そして、貯水タンク１内の７リットルの水を３０万個の粒子で表現した。これにより、理論的には、定常状態では毎秒６４００個の粒子が６０ｃｍ／秒の粒子速度で流入境界から流入し、流出境界から流出することになる。

粒子法の具体的な流体解析手法としてはＳＰＨ法を用い、各時間ステップにおいて各粒子について運動方程式（ナビエストークス方程式）を解いて粒子の速度を更新した。また、１ステップ分の時間Δｔに応じて粒子の位置を更新し、排出管３と貯水タンク１との接続境界５である流出境界に到達した粒子を消去した。このような計算を繰り返し行い、上記水の流動挙動のシミュレーションを行った。その結果、貯水タンク１の水量が時間とともに増大するというシミュレーション結果になった。

上述のような結果になった理由としては、実際の水の流動挙動では、貯水タンク１と排出管３との接続境界５では流体は貯水タンク１の外に排出されているため、上記接続境界５での圧力は貯水タンク１の壁面よりも低いと考えられる。しかし、接続境界５に到達した粒子を消去する手法では、上述のような境界面での圧力低下の効果を計算に反映できないために、境界面に到達する粒子数が理論的な数値である毎秒６４００個よりも少なくなった結果、貯水タンク１の水量が時間とともに増大したと考えられる。

（比較例２）
比較例１では、貯水タンク１内のみを計算対象領域としたが、比較例２では、排出管３の内部の領域のうち、貯水タンク１との接続境界５から５ｃｍまでの領域を助走領域として設定し、貯水タンク１内の主領域と助走領域とを計算対象領域とした。また、貯水タンク１との接続境界５から５ｃｍのところに到達した粒子を消去した。このように構成した理由は、実際の流れのように流出口に向けた流れを計算でも再現することによって、接続境界５での圧力が貯水タンク１の壁面よりも低い状態を再現するためである。

そして、助走領域に流入した粒子の速度を、理論的に求めた粒子速度と等しい速度である６０ｃｍ／秒に固定してシミュレーションを行った。その他の条件は比較例１と全て同じである。その結果、今度は貯水タンク１の水量が時間とともに減少するというシミュレーション結果になった。

このように正しく計算できない理由は下記と思われる。すなわち、粒子法は密度誤差が大きい計算手法であり、非圧縮条件を課して計算しているとはいっても、数パーセントの誤差ないしゆらぎは存在している。今回の計算条件では、設定した流出口位置の流出境界面を通過して来る粒子群の密度が平均密度よりも高かったために、理論的に求めた粒子速度と等しい速度６０ｃｍ／秒で吸引すると、流出境界面に到達する粒子数が理論的な数値である毎秒６４００個よりも多くなってしまい、貯水タンク１の水量が時間とともに減少してしまう結果となったと考えられる。

（比較例３）
上述の比較例１および比較例２のような試行の結果を踏まえ、接続境界５に隣接した助走領域７を設けて助走領域７の終端で粒子を消去する計算方法とすることによって、実際の流れのように流出口に向けた流れを計算でも再現することにより接続境界面での圧力が壁面よりも低い状態を再現するようにしつつ、接続境界面に到達する粒子群に密度誤差があっても定常的な流出流動計算が正しくできるように、助走領域７に入った粒子の粒子速度を制御する方法を考えた。

助走領域７の粒子速度を制御する方法として、計算中に何らかの物理量をセンシングして、定常状態であるとすると実現できる当該物理量の値との偏差をフィードバックさせることによって、助走領域７に入った粒子の粒子速度を制御することで、定常状態の流出流動計算を実現することを考えた。上記助走領域７に流入した粒子の速度制御は、例えば、流出境界６に単位時間当たりに到達する粒子数の偏差量に応じて行うことができる。具体的には、ある時刻ｔにおいて、単位時間当たりに流入境界から計算対象領域内に流入する粒子数をＮ_inとし、ある時刻ｔにおいて単位時間当たりに流出境界６から流出する粒子数をＮ_out（ｔ）とした時、それらの差ΔＮ₁（ｔ）を偏差量として求める。ΔＮ₁（ｔ）は下記の式で与えられる。
ΔＮ₁（ｔ）＝Ｎ_in（ｔ）−Ｎ_out（ｔ） （３）
上記式（１）を用いて、Δｔ後の速度Ｖ（ｔ＋Δｔ）は、下記式（４）のように制御することができる。
Ｖ（ｔ＋Δｔ）＝Ｖ（ｔ）×｛１＋βΔＮ₁（ｔ）／Ｎ_in｝ （４）
ここで、βは定数であり、計算対象に応じて適切に設定される。

比較例２と同様に、貯水タンク１である主領域と、貯水タンク１との接続境界５から５ｃｍまでの助走領域７とを計算対象領域として、流体のシミュレーションを行った。ただし、助走領域７内の粒子の速度を、式（３）および（４）に基づいて変更した。ここで、Ｎ_in＝６４００個／秒とした。その他の条件は比較例２と全て同じとした。

式（４）のβの値を適宜変更してシミュレーションを行ったところ、流出量の時間平均はほぼＮ_out＝６４００個／秒に制御できたものの、流出量の時間変動が非常に大きく計算が不安定となることが分かった。βの値を種々変更しても計算の不安定さは解消されなかった。

（実施例１）
比較例３と同様に、貯水タンク１である主領域と助走領域７とを計算対象領域として、流体のシミュレーションを行った。ただし、助走領域７内の粒子の速度を、式（１）および（２）に基づいて変更した。

その結果α＝９．８であれば、流出量の時間変動はわずかにあるものの、平均流出量は所望の値である６４００個／秒であり、水面の変動もほとんどないという安定的なシミュレーション結果を得た。

比較例３で計算が不安定となる理由は、比較例２の項でも記載したように粒子法は密度誤差が大きい計算手法であるためと考えられる。すなわち、計算主領域である貯水タンク１内粒子の密度にもわずかに誤差はあり、助走領域７内粒子の密度誤差はさらに大きいものと考えられる。比較例３の方法はこの助走領域７内の粒子の密度誤差はないものとして当該領域内の粒子速度を制御するものではあるが、そのフィードバック制御が必ずしも助走領域７内の粒子の密度誤差を減少させるようには働かないためと考えられる。

一方、実施例１の方法は、そのフィードバック制御量を密度誤差変動の大きい助走領域７内の粒子に直接依存させずに、計算対象の大部分を占める主計算領域の粒子数誤差を減少させる方向にフィードバック制御量を決定するようにしているので、計算が安定すると同時に、主領域の粒子密度誤差も減少させる方向に働いていると考えられる。

上記は、貯水タンク１への水の流入量と貯水タンク１からの水の流出量は同じであるという定常流れ前提でシミュレーションを行ったが、貯水タンク１への水の流入量と貯水タンク１からの水の流出量が異なる場合であっても、適用できることは言うまでもない。実施例１の系では、α＝９．８であれば実際の定常流れ流体の挙動を再現できることが判明した。よって、貯水タンク１への水の流入量と水の貯水タンク１からの流出量が異なる場合であっても、αを適宜調整すれば、実際の流体の挙動を再現できることになる。

また、貯水タンク１への水の流入がない、すなわち理論粒子数が単純に減少していく場合でも、本発明が適用可能であるし、貯水タンク１への水の流入量と貯水タンク１からの水の流出量が、ある時点で変わる場合であっても、本発明が適用可能であり、変更された条件にあわせて理論粒子数と助走領域内の粒子の速度が正しく計算される。また、非圧縮性流体が水以外の場合であっても、本発明が適用可能である。非圧縮性流体として、各種の水溶液やコロイド溶液、マッハ数が０．３以下の空気などの気体が例示できる。

本発明によれば、本発明によれば、粒子法により、流体の一部が計算対象領域から特定の流出量で流出するシミュレーションを行う際の流出境界条件を設定することができる。

１ 貯蔵部（貯水タンク）
２ 供給管
３ 排出管
４，５ 接続境界
６ 流出境界
７ 助走領域

Claims (2)

1. あらかじめ設定された主領域から流体が流出する系で、粒子法を用いて前記流体の挙動を計算するに際し、課された流出境界条件通りの粒子法シミュレーションを実現させる方法であって、
   前記主領域に助走領域を接続し、
   前記助走領域内の粒子の速度を、前記主領域への前記流体の流入境界条件と流出境界条件からなる境界条件および初期条件とから求められる前記主領域内の理論粒子数と、計算途上における前記主領域内の全粒子数との差に応じて変更する、流体シミュレーション方法。
2. あらかじめ設定された主領域から流体が流出する系で、粒子法を用いて前記流体の挙動を計算するプログラムであって、
   コンピュータに、
   前記主領域に接続する助走領域を設定する処理、および当該助走領域内の粒子の速度を、前記主領域内の理論粒子数と、前記計算途上における主領域内の全粒子数の差に応じて変更する処理を実行させるための流体シミュレーションプログラム。

Images