流動解析方法、流動解析装置、及び流動解析プログラム
【課題】計算コストを抑えるとともに固定壁の境界条件について任意形状に対応可能であり、さらに計算精度の良い流動解析方法、流動解析装置、及び流動解析プログラムを提供する。
【解決手段】粒子法を用いて流体の挙動を解析する流動解析方法であって、流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する距離判定ステップと、距離判定ステップにより壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の壁境界の曲率を算出する曲率算出ステップと、距離判定ステップにより壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、壁境界を挟んで対称に壁を表す複数の第２粒子を曲率算出ステップにより算出された曲率に応じた拡大率で複製して配置する壁粒子配置ステップと、壁粒子配置ステップにより配置された複数の第２粒子を境界条件として粒子間の相互作用に基づき流体の挙動を解析する解析ステップとを備える。
【解決手段】粒子法を用いて流体の挙動を解析する流動解析方法であって、流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する距離判定ステップと、距離判定ステップにより壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の壁境界の曲率を算出する曲率算出ステップと、距離判定ステップにより壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、壁境界を挟んで対称に壁を表す複数の第２粒子を曲率算出ステップにより算出された曲率に応じた拡大率で複製して配置する壁粒子配置ステップと、壁粒子配置ステップにより配置された複数の第２粒子を境界条件として粒子間の相互作用に基づき流体の挙動を解析する解析ステップとを備える。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、粒子法を用いて連続体である流体の挙動を解析する流動解析方法、流動解析装置、及び流動解析プログラムに関する。
【背景技術】
【０００２】
従来、流動解析方法の一つとして、均一な大きさの複数の粒子により流体を表現し、粒子間の相互作用の計算により複数の粒子の挙動をコンピュータを用いて解析する方法、いわゆる粒子法を利用した流動解析方法が知られている。この粒子法は、有限体積法や有限要素法の計算において、手間のかかる格子生成が不要であるという利点を有する。この粒子法を利用した流動解析方法は、空気や水等の自由表面を有する物体、特に大きく変形したり、分裂，合体したりする物体の挙動解析に適しており、原子力分野を始めとしてプラント，航空，船舶，環境等の様々な技術分野において広く利用されている。
【０００３】
粒子法を用いて流動解析を行うに際し、一般的に固定壁の境界条件は、位置を固定した複数の粒子により表現される。固定壁粒子には圧力変数が与えられており、流体粒子が壁に集まってきた場合に固定壁粒子の圧力は上昇し、流体粒子をはね返す。図８は、従来の流動解析方法における粒子法による壁粒子を用いた境界条件を示す概念図である。図８において、四角の枠線（境界）内に記載された複数の粒子２５は、流体を示すものとする。また、四角の枠線（境界）外に記載された複数の粒子２６は、壁粒子を示すものとする。図８に示すように、等しい大きさの複数の粒子２６からなる壁粒子は、格子状に規則正しく配置され、壁境界を表現する。このような手法を従来手法１と呼ぶこととする。従来手法１は、固定壁による流体粒子に与える力（すなわち壁境界の影響）を比較的正確に計算することができるため、流動解析の計算精度が良いという利点を有する。
【０００４】
例えば特許文献１には、通常のマクロな流れが解析できるような非統計的な粒子法を用いて流動を支配する物理過程を粒子間相互作用で表現し、これを用いた簡単な計算手法によって、これまで解析不可能であった砕波等の流動解析を短時間で高精度に行う流動解析装置及び流動解析方法が記載されている。特許文献１に記載の流動解析方法においても、流体を収容する容器（すなわち固定壁）は、座標位置を固定して流速０の複数粒子の集まりとして表現する旨が記載されており、従来手法１を利用したものである。
【０００５】
また、特許文献２には、汎用性が高く、濡れ性モデルを統一的に取り扱うことができ、空間解像度が粗くても精度よく解析を行うことができる流体解析装置及び流体解析方法が記載されている。この特許文献２における流体解析方法は、表面張力モデルにおいてエネルギ関係式を用いて壁の粒子からの作用力を算出する点に特徴があるが、壁を固定した複数粒子により表現する点においては特許文献１と同様に一般的な手法（従来手法１）を用いている。
【０００６】
壁境界を表現する別の手法として、壁境界付近の流体粒子を境界に対して対象となるように複製することにより、動的に壁粒子を生成する手法がある。図９は、従来の流動解析方法における粒子法による境界条件を境界付近の流体粒子を複製して生成した壁粒子により表現する場合の概念図である。図９において、四角の枠線（境界）内に記載された複数の粒子２７ａは、流体を示すものとする。また、四角の枠線（境界）外に記載された複数の粒子２７ｂは、壁粒子を示すものであり、境界に対して対称な位置の粒子２７ａを複製して生成されたものである。このような手法を従来手法２と呼ぶこととする。
【０００７】
従来手法２は、境界に流体粒子が隣接する場合にのみ壁粒子が生成されるため、流体の動きに応じて壁粒子の分布状況も動的に変化する。すなわち、従来手法２は、境界に流体粒子が隣接しない部分において壁粒子を生成する必要が無く、計算コストを抑えられる（計算時間が短い）という利点を有する。なお、従来手法２は、非特許文献１に記載されている。
【０００８】
さらに、壁境界を表現する別の手法として、壁粒子を一切用いずに壁からの反発力を計算する手法がある。図１０は、従来の流動解析方法における粒子法による境界条件を壁からの反発力の存在により表現する場合の概念図である。図１０において、四角の枠線（境界）内に記載された複数の粒子は、流体を示すものとする。図１０に示すように、境界内の複数の流体粒子は、壁に接触した際に壁に対してめり込んだ量に応じた壁からの反発力を受ける。このような手法を従来手法３と呼ぶこととする。
【０００９】
従来手法３は、壁に接触した流体粒子がペナルティ法で計算された壁境界からの反発力を受けるものとして行う解析手法であるが、壁粒子を一切用いないため、計算コストを抑えられる（計算時間が短い）という利点を有する。また、従来手法３は、壁境界の形状として図１０に示す矩形形状のような単純な形状のみならず、任意形状（ＣＡＤデータ）を取り扱うことが可能であるといった利点も有する。なお、従来手法３は、非特許文献２，３に記載されている。
【先行技術文献】
【特許文献】
【００１０】
【特許文献１】特開平７−３３４４８４号公報
【特許文献２】特開２００８−１１１６７５号公報
【非特許文献】
【００１１】
【非特許文献１】Ｓ．Ｊ．Ｃｕｍｍｉｎｓ ａｎｄ Ｍ．Ｒｕｄｍａｎ，Ａｎ ＳＰＨ ｐｒｏｊｅｃｔｉｏｎ ｍｅｔｈｏｄ，Ｊｏｕｒｎａｌ ｏｆ Ｃｏｍｐｕｔａｔｉｏｎａｌ Ｐｈｙｓｉｃｓ，１５２，５８４−６０７（１９９９）．
【非特許文献２】原田隆宏，越塚誠一，ＳＰＨにおける壁境界計算手法の改良，情報処理学会論文誌，４８，４，１８３８−１８４６（２００７）．
【非特許文献３】原田隆宏，越塚誠一，島崎克教，ＭＰＳ法における壁境界モデルの改良、日本計算工学会論文集，２００８，Ｐａｐｅｒ Ｎｏ．２００８０００６（２００８）．
【発明の概要】
【発明が解決しようとする課題】
【００１２】
しかしながら、上述した従来手法１は、等しい大きさの複数の粒子からなる壁粒子を格子状に規則正しく配置して壁境界を表現するため、複雑形状に不適切であるという欠点がある。この問題点について図１１を用いて説明する。図１１は、従来の流動解析方法において従来手法１により境界条件を表現する際の問題点を示す概念図である。図１１に示すように、境界が例えば滑らかな曲面を描いた円形形状等である場合に、格子状に規則正しく配置した壁粒子は、滑らかなカーブを表現することが困難であり、境界がいびつになるという問題点がある。壁粒子の大きさを小さなものにすれば滑らかなカーブに近い形状を表現することも可能であるが、当然のことながら多くの壁粒子を必要とし、計算コストは高く（計算時間は長く）なる。さらに、そもそも従来手法１は、従来手法２，３に比して、壁部分全体について壁粒子を配置する必要があるため、多くの壁粒子を必要とし、計算コストが高いという欠点を有する。
【００１３】
従来手法２は、壁境界付近の流体粒子を境界に対して対称となるように複製して壁粒子とするものであるため、コーナー部や複雑形状に対応が困難であり汎用性が低いという欠点を有する。例えば図９に示すように、従来手法２による壁粒子は、直線的な壁（境界）を模擬することは容易であるが、コーナー部には壁が存在しないこととなってしまい、所望する解析結果が得られないことも考えられる。
【００１４】
図１２は、従来手法２の方法を発展させたものであり、境界に対して対称となるように流体粒子を１度複製して壁粒子を生成した後に、さらに別の境界に対して対称となる位置に壁粒子を複製し、コーナー部の隙間を埋めるように壁粒子の配置を試みた場合の概念図である。図１２に示すように、流体粒子２８ａは、境界に対称な位置に壁粒子２８ｂを複製させる。さらに、壁粒子２８ｂは、別の境界に対称な位置に壁粒子２８ｃを複製させることにより、境界を表現する壁のコーナー部を形成する。同様に、流体粒子２９ａは、境界に対称な位置に壁粒子２９ｂを複製させる。さらに、壁粒子２９ｂは、別の境界に対称な位置に壁粒子２９ｃを複製させる。
【００１５】
しかしながら、図１２に示す手法は、左下に示すように１つの流体粒子３０が境界に対称な位置に壁粒子３１，３２を複製させる場合が考えられる。この場合に、壁粒子３１，３２の各々により境界に対称な位置に複製された壁粒子は、複製先において位置が重なってしまうという現象が生じうる。したがって、図１２に示す手法は、同じ位置に複数の壁粒子が存在する場合があり、流動解析の結果が不正確になる可能性がある。
【００１６】
従来手法３は、上述したように壁粒子を必要としないため、計算コストを抑えるとともに任意形状にも対応可能であるという利点を有するが、ペナルティ法により計算された反発力が壁に対して法線方向の力のみであるという欠点があり、計算精度が悪いという問題点がある。実際の固定壁は圧力勾配を有し、流体に対して法線方向のみならず接線方向の力も作用するからである。したがって、従来手法３を用いて行った流動解析は、計算精度が悪いため、場合によっては流体の挙動自体に悪影響を及ぼす可能性も考えられ、信頼性が低いという欠点がある。
【００１７】
本発明は上述した従来技術の問題点を解決するもので、計算コストを抑えるとともに固定壁の境界条件について任意形状に対応可能であり、さらに計算精度の良い流動解析方法、流動解析装置、及び流動解析プログラムを提供することを課題とする。
【課題を解決するための手段】
【００１８】
本発明に係る流動解析方法は、上記課題を解決するために、粒子法を用いて流体の挙動を解析する流動解析方法であって、流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する距離判定ステップと、前記距離判定ステップにより前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の前記壁境界の曲率を算出する曲率算出ステップと、前記距離判定ステップにより前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、前記壁境界を挟んで対称に壁を表す複数の第２粒子を前記曲率算出ステップにより算出された曲率に応じた拡大率で複製して配置する壁粒子配置ステップと、前記壁粒子配置ステップにより配置された前記複数の第２粒子を境界条件として粒子間の相互作用に基づき流体の挙動を解析する解析ステップとを備えることを特徴とする。
【発明の効果】
【００１９】
本発明によれば、計算コストを抑えるとともに固定壁の境界条件について任意形状に対応可能であり、さらに計算精度の良い流動解析を行うことができる。
【図面の簡単な説明】
【００２０】
【図１】本発明の実施例１の形態の流動解析装置の構成を示すブロック図である。
【図２】本発明の実施例１の形態の流動解析装置における壁境界計算部の詳細な構成を示すブロック図である。
【図３】本発明の実施例１の形態の流動解析装置の動作を示すフローチャート図である。
【図４】本発明の実施例１の形態の流動解析方法における格子幅が大きい場合の法線方向を示す模式図である。
【図５】本発明の実施例１の形態の流動解析方法における格子幅が小さい場合の法線方向を示す模式図である。
【図６】本発明の実施例１の形態の流動解析方法における粒子法による境界条件を境界付近の流体粒子を曲率に応じた拡大率で複製した壁粒子により表現する場合の概念図である。
【図７】本発明の実施例１の形態の流動解析方法における粒子法による境界条件を境界付近の流体粒子を曲率に応じた拡大率で複製した壁粒子により表現する場合の概念図である。
【図８】従来の流動解析方法における粒子法による壁粒子を用いた境界条件を示す概念図である。
【図９】従来の流動解析方法における粒子法による境界条件を境界付近の流体粒子を複製して生成した壁粒子により表現する場合の概念図である。
【図１０】従来の流動解析方法における粒子法による境界条件を壁からの反発力の存在により表現する場合の概念図である。
【図１１】従来の流動解析方法において従来手法１により境界条件を表現する際の問題点を示す概念図である。
【図１２】従来の流動解析方法において従来手法２を発展させた方法により境界条件を表現する際の問題点を示す概念図である。
【発明を実施するための形態】
【００２１】
以下、本発明の流動解析方法、流動解析装置、及び流動解析プログラムの実施の形態について図面を参照して詳細に説明する。
【実施例１】
【００２２】
図１は、本発明の実施例１の流動解析装置１の構成を示すブロック図である。図１を参照して、流動解析装置１の構成を説明する。本発明の実施形態となる流動解析装置１は、粒子法を用いて流体の挙動を解析する装置であり、パーソナルコンピュータ，ワークステーション，汎用コンピュータ等の公知の情報処理装置により構成され、図１に示すように、ＲＡＭ（Ｒａｎｄｏｍ Ａｃｃｅｓｓ Ｍｅｍｏｒｙ）２、ＲＯＭ（Ｒｅａｄ Ｏｎｌｙ Ｍｅｍｏｒｙ）３、初期設定データベース（ＤＢ）４、及びＣＰＵ（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）５を主な構成要素として備える。ＲＡＭ２は、ＣＰＵ５が実行する流動解析装置１の起動プログラムや流動解析プログラム等のコンピュータプログラム及び各種データを一時的に格納するワーキングエリアを提供する。ＲＯＭ３は、コンピュータプログラム及びこのコンピュータプログラムを実行するために必要な各種データを記憶する。なおＲＯＭ３は、磁気的又は光学的な記憶媒体若しくは半導体メモリ等のＣＰＵ５により読み取り可能な記憶媒体を含む構成を有し、ＲＯＭ３に記憶されているコンピュータプログラムや各種データの一部又は全部は電気通信回線を介してダウンロードされるように構成してもよい。
【００２３】
初期設定ＤＢ４は、粒子法により流動解析を実行するために必要な各種パラメータを格納する。具体的には、このパラメータには、粒子法における粒子ｉの粒子径ｒｉの初期値ｒiini、位置（重心の座標）、速度、温度、圧力、流体の物性（例えば粘度、密度、熱伝導率、比熱、潜熱等）、時間刻み幅Δｔ、計算終了時間Ｔ等の初期条件に関する情報が含まれる。ＣＰＵ５は、ＲＯＭ３に記憶されているコンピュータプログラムをＲＡＭ２内にロードし、読み出されたコンピュータプログラムを実行することにより、流動解析装置１全体の処理動作を制御する。具体的には、ＣＰＵ５は、ＲＯＭ３に記憶された流動解析プログラム及び流動解析プログラムを実行するために必要なデータをＲＡＭ２へロードし、流動解析プログラムに従って粒子法を利用した流動解析を実行する。本実施形態では、ＣＰＵ５は、ＲＯＭ３に記憶されている流動解析プログラムを実行することにより、壁境界計算部１０、外力項計算部１１、粘性項計算部１２、粒子移動追跡部１３、圧力項計算部１４、粒子データ更新部１５、及び終了判定部１６として機能する。壁境界計算部１０以外の各部における機能の詳細については後述する。
【００２４】
図２は、本実施例の流動解析装置１におけるＣＰＵ５内部の壁境界計算部１０の詳細な構成を示すブロック図であり、図２に示すように距離判定部１７、曲率算出部１８、及び壁粒子配置部１９から構成される。
【００２５】
距離判定部１７は、流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する。この「所定距離」は、予め距離判定部１７に設定されているものでもよいし、ユーザにより設定可能としてもよい。
【００２６】
曲率算出部１８は、距離判定部１７により壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の壁境界の曲率を算出する。ここで、曲率算出部１８は、予め壁境界からの距離を定めた距離関数に基づいて壁境界の曲率を算出する。具体的な算出方法は後述する。
【００２７】
壁粒子配置部１９は、距離判定部１７により壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、壁境界を挟んで対称に壁を表す複数の第２粒子を曲率算出部１８により算出された曲率に応じた拡大率で複製して配置する。
【００２８】
また、本実施例における壁粒子配置部１９は、配置した複数の第２粒子（壁粒子）の各々が有する圧力値を対応する第１粒子が有する圧力値と同じ値に設定する。
【００２９】
壁境界計算部１０を除くＣＰＵ５（すなわち外力項計算部１１、粘性項計算部１２、粒子移動追跡部１３、圧力項計算部１４、粒子データ更新部１５、及び終了判定部１６）は、本発明の解析部に対応し、壁粒子配置部１９により配置された複数の第２粒子（壁粒子）を境界条件として粒子間の相互作用に基づき流体の挙動を解析する。
【００３０】
次に、上述のように構成された本実施の形態の作用を説明する。図３は、本実施例の流動解析装置１の動作を示すフローチャート図である。流動解析装置１は、以下に示す流動解析処理を実行することにより流体の挙動を解析する。以下、図３に示すフローチャートを参照して、この流動解析処理を実行する際の流動解析装置１の動作及び本発明の流動解析方法について詳しく説明する。
【００３１】
図３に示すフローチャートは、流動解析装置１に対しオペレータが流動解析処理の開始を指示したタイミングで開始となり、流動解析処理はステップＳ１の処理に進む。
【００３２】
ステップＳ１の処理では、ＣＰＵ５が、初期設定ＤＢ４からＲＡＭ２内に解析対象となる流体の物性，時間刻み幅ΔＴ，計算終了時間Ｔ等の計算条件パラメータ、及び粒子の座標位置，速度，温度等の初期条件パラメータを読み出し、読み出されたパラメータに基づいて流体を複数の粒子により表現する。これにより、ステップＳ１の処理は完了し、流動解析処理はステップＳ２の処理に進む。
【００３３】
ステップＳ２の処理は、壁境界の計算を行うものである。まず、壁境界計算部１０内の距離判定部１７は、流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する（距離判定ステップ）。これは、各粒子について壁粒子として複製するか否かを判断するための処理であり、壁境界から離れた位置にある粒子は、壁粒子として複製されない。
【００３４】
を計算する。ここで、ｉは、粒子番号を示すものとする。本実施例において、この距離の計算は、空間に距離関数を定義することにより行われる。すなわち、距離判定部１７は、計算空間を格子状に区切り、各格子点において最も近い壁境界までの距離をあらかじめ計算した上で記憶しておく。粒子位置における壁までの最小距離は、その距離関数を補間することにより求められる。ここで、距離関数について、詳細に説明する。
【００３５】
本実施例で使用する距離関数は、ある任意の境界面からの距離を空間上に定義したものである。一般的には、計算領域を覆うようにして正方格子状に領域を分割し、その各格子点において壁境界からの距離を計算したものを用いる。壁の外側を正の値、内側を負の値とする。壁境界が動かない場合、距離関数は変化することがないため、前処理により計算しておくことで解析の計算コストを大幅に下げることが可能になる。
【００３６】
【数１】
【００３７】
【００３８】
より計算コストを抑えることができる。ここで、粒子ｉは、以下に示す８点で作られる立方体の中に存在しているとする。１：（Ｉ，Ｊ，Ｋ）、２：（Ｉ＋１，Ｊ，Ｋ）、３：（Ｉ＋１，Ｊ＋１，Ｋ）、４：（Ｉ，Ｊ＋１，Ｋ）、５：（Ｉ，Ｊ，Ｋ＋１）、６：（Ｉ＋１，Ｊ，Ｋ＋１）、７：（Ｉ＋１，Ｊ＋１，Ｋ＋１）、８：（Ｉ，Ｊ＋１，Ｋ＋１）。
【００３９】
格子点（Ｉ，Ｊ，Ｋ）は、（２）式により表される。
【数２】
【００４０】
【数３】
【００４１】
このようにして流体を表す複数の第１粒子の各々と壁境界との間の距離を求めた距離判定部１７は、当該距離が所定距離以下であるか否かを判定し、判定結果を出力する。
【００４２】
次に、曲率算出部１８は、距離判定部１７により壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の壁境界の曲率を算出する（曲率算出ステップ）。
【００４３】
【数４】
【００４４】
この法線方向は、粒子を複製する方向となる。次に、曲率算出部１８は、法線方向の発散を計算することにより、（５）式に示すように曲率ｋを求めることができる。曲率は、複製後の粒子の大きさを決定するのに用いられる。
【数５】
【００４５】
曲率算出部１８は、上述したように距離関数を利用して法線方向及び曲率を算出している。したがって、曲率算出部１８は、曲率算出ステップにおいて予め壁境界からの距離を定めた距離関数に基づいて壁境界の曲率を算出するといえる。
【００４６】
法線方向と曲率に関しては，粒子上の距離関数を用いて計算する方法と、格子点上において前処理により計算しておき、それを粒子上に線型補間する方法の２つが考えられる。１つ目の方法では従来の粒子法の定式化より、空間微分と発散を計算することにより求めることができる。２つ目の方法では以下のように法線方向と曲率を前計算しておく。格子点（Ｉ，Ｊ，Ｋ）における法線方向は、（６）式により表される。
【数６】
【００４７】
また、（６）式における〈∇ｄ〉Ｉ，Ｊ，Ｋは、（７）式のように表される。
【数７】
【００４８】
さらに、格子点（Ｉ，Ｊ，Ｋ）における曲率ｋ（Ｉ，Ｊ，Ｋ）は、（８）式のように表される。
【数８】
【００４９】
これらの格子点上の値を補間することにより、粒子上の値が計算される。
【００５０】
次に、法線方向と解像度の関係について述べる。上述した（６）式及び（７）式を使用して法線方向を算出する際に、角近辺における法線方向は、解像度すなわち格子の幅の影響を受ける。図４は、格子幅が大きい場合の法線方向を示す模式図である。また、図５は、格子幅が小さい場合の法線方向を示す模式図である。図４及び図５における太線は、壁境界面である。また、矢印は、法線方向を表す。ここで、点線において囲った点の法線方向を比較する。格子幅が大きい場合における法線方向は、図４に示すように斜め方向を指す。一方、格子幅が小さい場合における法線方向は、図５に示すように真横方向を指す。これは、法線方向を（６）式及び（７）式を用いて計算する際に、隣の距離関数を参照するため、角の近辺にある格子点では法線方向が斜め方向に算出されるためである。したがって、格子幅が大きいほど法線方向が斜め方向に計算される領域は広くなる。この格子幅を粒子径と同程度にすることにより、後述する複製された壁粒子が重ならないような法線方向が計算される。
【００５１】
次に、壁粒子配置部１９は、距離判定部１７により壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、壁境界を挟んで対称に壁を表す複数の第２粒子を複製する。ここで、壁粒子配置部１９は、複数の第２粒子を複製する際に、曲率算出部１８により算出された曲率に応じた拡大率で複製して配置する（壁粒子配置ステップ）。
【００５２】
粒子ｉの拡大率ｒｉは、粒子上の曲率ｋｉ及び粒子径Ｌｉを用いて、（９）式のように表される。
【数９】
【００５３】
【数１０】
【００５４】
壁粒子配置部１９は、（１１）式に示すように、拡大率が大きい場合には壁境界から離れた位置に複製した粒子を配置し、拡大率が小さい場合には壁境界に近い位置に複製した粒子を配置することにより、隙間の少ない壁を形成させる。
【００５５】
なお、壁粒子配置部１９は、複製した粒子を壁粒子（第２粒子）として配置する際に、配置した複数の第２粒子の各々が有する圧力値を対応する第１粒子（複製元の粒子）が有する圧力値と同じ値に設定する。設定した圧力値は、ステップＳ７で圧力項を計算する際に使用される。
【００５６】
図６及び図７は、本実施例の流動解析方法における粒子法による境界条件を境界付近の流体粒子を曲率に応じた拡大率で複製した壁粒子により表現する場合の概念図である。
【００５７】
壁粒子配置部１９は、例えば図６に示すように、壁境界との間の距離が所定距離以下であると判断された流体の粒子２０ａに対応して、壁境界を挟んで対称に壁粒子２０ｂを複製して配置する。この場合において、粒子２０ａ近傍の壁境界は直線的な壁であるため、曲率ｋ＝０である。したがって、拡大率ｒ＝１となり、壁粒子配置部１９は、粒子２０ａと同じ大きさの粒子２０ｂを複製し、壁粒子として配置する。
【００５８】
また、壁粒子配置部１９は、図６に示すように、壁境界との間の距離が所定距離以下であると判断された流体の粒子２１ａ，２２ａに対応して、壁境界を挟んで対称に壁粒子２１ｂ，２２ｂを複製して配置する。この場合において、粒子２１ａ，２２ａ近傍の壁境界は曲率が負（凹形状）であるため、拡大率ｒ＞１となり、壁粒子配置部１９は、粒子２１ａ，２２ａを大きく拡大して粒子２１ｂ，２２ｂを複製し、壁粒子として配置する。すなわち、曲率が負の場合には、流体の領域に比して壁の領域が広いため、壁粒子配置部１９は、壁の隙間を無くすように大きな壁粒子を配置する。
【００５９】
また、壁粒子配置部１９は、図７に示すように、壁境界との間の距離が所定距離以下であると判断された流体の粒子２３ａ，２４ａに対応して、壁境界を挟んで対称に壁粒子２３ｂ，２４ｂを複製して配置する。この場合において、粒子２３ａ，２４ａ近傍の壁境界は曲率が正（凸形状）であるため、拡大率ｒ＜１となり、壁粒子配置部１９は、粒子２３ａ，２４ａを小さく縮小して粒子２３ｂ，２４ｂを複製し、壁粒子として配置する。すなわち、曲率が正の場合には、流体の領域に比して壁の領域が狭いため、壁粒子配置部１９は、壁の隙間を無くすように小さな壁粒子を配置する。
【００６０】
なお、壁粒子配置部１９は、第１粒子と壁境界を挟んで対称に壁を表す第２粒子を複製するものであるが、ここで言う「対称」とは、算出した「法線方向」上に壁境界を挟んで「対称」の意である。したがって、図６に示す粒子２２ａ，２２ｂの組や、図７に示す粒子２４ａ，２４ｂの組は、一見、点線で示した壁境界に対して対称でないように見えるが、算出した法線方向上に対称に配置されており、法線方向が壁境界に対して斜め方向となる理由は図４，５を用いて説明したとおりである。また、壁粒子から壁境界までの距離は、拡大率に応じて変化する。
【００６１】
壁境界の計算が終わり、複数の第２粒子が壁粒子として配置されると、ステップＳ２の処理は完了し、流動解析処理はステップＳ３の処理に進む。
【００６２】
ステップＳ３以降のステップにおいて、本発明の流動解析装置は、壁粒子配置ステップにより配置された複数の第２粒子を境界条件として粒子間の相互作用に基づき流体の挙動を解析する（解析ステップ）。すなわち、ステップ３以降の動作は、全て解析ステップである。
【００６３】
ステップＳ３の処理では、外力項計算部１１が、重力等の各粒子に直接的に作用する外力項を計算する。粒子法における外力項の計算方法は本願発明の出願時点で既に公知であるので詳細な説明は省略する。これにより、ステップＳ３の処理は完了し、流動解析処理はステップＳ４の処理に進む。
【００６４】
ステップＳ４の処理では、粘性項計算部１２が、粘性項を計算し、粒子ｉが有する速度分布を近傍の粒子ｊに分配することにより、粘性による粒子の速度の拡散を計算する。具体的には、粘性項は以下の（１２）式により表される。
【数１１】
【００６５】
（１２）式中におけるパラメータｕ、ｔ、及びνは、それぞれ速度ベクトル、計算時刻、動粘性係数を示す。この（１２）式を粒子法で離散化すると、粘性項は陽的には（１３）式、陰的には（１４）式のように表される。
【数１２】
【００６６】
（１３）式及び（１４）式においてパラメータｄ、λ、ｎ０、ｗｉｊ、及びｋは、それぞれ次元数、拡散を解析解に一致させるための係数、基準粒子密度、粒子ｉ，ｊ間の重み、及びタイムステップを示す。
【００６７】
粒子ｉの速度変化Δｕは、以下の（１５）式により表されるので、これを（１４）式に代入することにより、以下の（１６）式に示す連立一次方程式が導かれる。粘性項計算部１２は、この（１６）式を解くことにより粘性項を陰的に解く。
【数１３】
【００６８】
このように粘性項を陰的に解くことにより、タイムステップｋを大きくすることができるので、発散することなく計算を安定的に行うことができる。これにより、このステップＳ４の処理は完了し、流動解析処理はステップＳ５の処理に進む。
【００６９】
ステップＳ５の処理では、粒子移動追跡部１３が、ステップＳ３及びステップＳ４の計算結果に基づいて粒子の速度を仮更新する。これにより、ステップＳ５の処理は完了し、流動解析処理はステップＳ６の処理に進む。
【００７０】
ステップＳ６の処理では、粒子移動追跡部１３が、ステップＳ５の処理により仮更新された速度による時間刻み幅Δｔの間での粒子の移動を追跡し、追跡結果に基づいて粒子の座標位置を仮更新する。これにより、ステップＳ６の処理は完了し、流動解析処理はステップＳ７の処理に進む。
【００７１】
ステップＳ７の処理では、圧力項計算部１４が、ステップＳ６の処理により算出された座標位置における各粒子の粒子数密度を算出し、算出された粒子数密度が所定値であるか否かを判別する。判別の結果、粒子数密度が所定値でない場合、圧力項計算部１４は、粒子数密度が所定値になるように粒子の圧力を修正することにより粒子の座標位置及び流速を修正する。そして圧力計算部１４は、修正された圧力によって生じる圧力項を算出する。なお、圧力計算部１４は、壁粒子の圧力については、ステップＳ２の壁粒子配置ステップにおいて設定された圧力値（すなわち複製元の流体粒子が有する圧力値と同じ圧力値）を用いて計算を行う。これにより、ステップＳ７の処理は完了し、流動解析処理はステップＳ８の処理に進む。
【００７２】
ステップＳ８の処理では、粒子データ更新部１５が、ステップＳ７の処理の算出結果に基づいてＲＡＭ２内に格納されている粒子の速度を更新する。これにより、ステップＳ８の処理は完了し、流動解析処理はステップＳ９の処理に進む。
【００７３】
ステップＳ９の処理では、粒子データ更新部１５が、ステップＳ７の処理の算出結果に基づいてＲＡＭ２内に格納されている粒子の座標位置を更新する。これにより、ステップＳ９の処理は完了し、流動解析処理はステップＳ１０の処理に進む。
【００７４】
ステップＳ１０の処理では、終了判定部１６が、解析時間ｔを時間刻み幅ΔＴだけインクリメントする。これにより、ステップＳ１０の処理は完了し、流動解析処理はステップＳ１１の処理に進む。
【００７５】
ステップＳ１１の処理では、終了判定部１６が、解析時間ｔが計算終了時刻Ｔであるか否かを判別する。判別の結果、解析時間ｔが計算終了時刻Ｔである場合、終了判定部１６は、ＲＡＭ２内に格納されている粒子データをＲＯＭ３の所定ファイルにコピーし、一連の流動解析処理を終了する。一方、解析時間ｔが計算終了時刻Ｔでない場合には、終了判定部１０は流動解析処理をステップＳ２の処理に戻す。
【００７６】
上述のとおり、本発明の実施例１の形態に係る流動解析装置１（流動解析方法）によれば、計算精度が良く、且つ計算コストを抑えるとともに固定壁の境界条件について任意形状（ＣＡＤデータ）に対応することができる。
【００７７】
すなわち、本発明の流動解析装置１は、等しい大きさの複数の粒子からなる壁粒子を格子上に規則正しく配置して壁境界を表現する従来手法１と異なり、壁境界の曲率に応じた拡大率で境界面付近の流体粒子を複製して壁粒子とするため、任意形状に対応可能という利点を有する。同様に、コーナー部や複雑形状に対応が困難な従来手法２と異なり、本発明の流動解析装置１は、コーナー部等にも隙間無く壁粒子を配置することができ、任意形状に対応可能であるのみならず、コーナー部等の壁粒子の影響も考慮して計算精度の良い解析結果を得ることができる。
【００７８】
さらに、本発明の流動解析装置１は、従来手法２と同様に境界に流体粒子が隣接しない部分において壁粒子を生成する必要が無く、計算コストを抑えられるという利点を有する。
【００７９】
また、ペナルティ法により計算された反発力が壁に対して法線方向の力のみであると従来手法３と異なり、本発明の流動解析装置１は、複製した壁粒子を配置することにより固定壁の圧力勾配も考慮に入れた計算を行うため、壁から流体に対する様々な方向の力を考慮に入れた計算精度の高い解析結果を得ることができる。
【００８０】
さらに、本発明の流動解析装置１は、複製した壁粒子を配置する際に、配置した複数の第２粒子（壁粒子）の各々が有する圧力値を対応する第１粒子（流体粒子）が有する圧力値と同じ値に設定するので、壁粒子の圧力について計算する必要が無く、計算コストを抑えることができるという利点を有する。
【００８１】
また、本発明の流動解析装置１は、壁の曲率を計算するために距離関数を使用しているので、予め各格子点における壁境界からの距離を計算したものを利用することができ、解析の計算コストを大幅に下げることが可能である。
【００８２】
また本実施形態における流動解析装置１は、モデム等の通信用部材を備えてインターネットを含む電気通信回線と接続可能なように構成してもよい。ここで、本発明の流動解析プログラムは、粒子法を用いて流体の挙動を解析する流動解析プログラムであって、流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する処理と、壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の壁境界の曲率を算出する処理と、壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、壁境界を挟んで対称に壁を表す複数の第２粒子を算出した壁境界の曲率に応じた拡大率で複製して配置する処理と、配置された複数の第２粒子を境界条件として粒子間の相互作用に基づき流体の挙動を解析する処理とをコンピュータに実行させるプログラムである。
【００８３】
当該流動解析プログラムは、電気通信回線からダウンロード等によって、ＲＡＭ２又はＲＯＭ３内に格納させてもよい。なお、この場合における通信ネットワークからダウンロードするためのダウンロードプログラムは、予めＲＯＭ３内に格納されているか、別の記憶媒体からＲＯＭ３内にインストールされるものとする。このように、本実施形態に基づいて当業者等によりなされる他の実施の形態、実施例及び運用技術等は全て本発明の範疇に含まれる。
【産業上の利用可能性】
【００８４】
本発明に係る流動解析方法、流動解析装置、及び流動解析プログラムは、粒子法を用いて連続体である流体の挙動を解析する流動解析装置等に利用可能である。
【符号の説明】
【００８５】
１…流動解析装置、２…ＲＡＭ、３…ＲＯＭ、４…初期設定ＤＢ、５…ＣＰＵ、１０…壁境界計算部、１１…外力項計算部、１２…粘性項計算部、１３…粒子移動追跡部、１４…圧力項計算部、１５…粒子データ更新部、１６…終了判定部、１７…距離判定部、１８…曲率算出部、１９…壁粒子配置部、２０ａ，２０ｂ，２１ａ，２１ｂ，２２ａ，２２ｂ，２３ａ，２３ｂ，２４ａ，２４ｂ，２５，２６，２７ａ，２７ｂ，２８ａ，２８ｂ，２８ｃ，２９ａ，２９ｂ，２９ｃ，３０，３１，３２…粒子。
【技術分野】
【０００１】
本発明は、粒子法を用いて連続体である流体の挙動を解析する流動解析方法、流動解析装置、及び流動解析プログラムに関する。
【背景技術】
【０００２】
従来、流動解析方法の一つとして、均一な大きさの複数の粒子により流体を表現し、粒子間の相互作用の計算により複数の粒子の挙動をコンピュータを用いて解析する方法、いわゆる粒子法を利用した流動解析方法が知られている。この粒子法は、有限体積法や有限要素法の計算において、手間のかかる格子生成が不要であるという利点を有する。この粒子法を利用した流動解析方法は、空気や水等の自由表面を有する物体、特に大きく変形したり、分裂，合体したりする物体の挙動解析に適しており、原子力分野を始めとしてプラント，航空，船舶，環境等の様々な技術分野において広く利用されている。
【０００３】
粒子法を用いて流動解析を行うに際し、一般的に固定壁の境界条件は、位置を固定した複数の粒子により表現される。固定壁粒子には圧力変数が与えられており、流体粒子が壁に集まってきた場合に固定壁粒子の圧力は上昇し、流体粒子をはね返す。図８は、従来の流動解析方法における粒子法による壁粒子を用いた境界条件を示す概念図である。図８において、四角の枠線（境界）内に記載された複数の粒子２５は、流体を示すものとする。また、四角の枠線（境界）外に記載された複数の粒子２６は、壁粒子を示すものとする。図８に示すように、等しい大きさの複数の粒子２６からなる壁粒子は、格子状に規則正しく配置され、壁境界を表現する。このような手法を従来手法１と呼ぶこととする。従来手法１は、固定壁による流体粒子に与える力（すなわち壁境界の影響）を比較的正確に計算することができるため、流動解析の計算精度が良いという利点を有する。
【０００４】
例えば特許文献１には、通常のマクロな流れが解析できるような非統計的な粒子法を用いて流動を支配する物理過程を粒子間相互作用で表現し、これを用いた簡単な計算手法によって、これまで解析不可能であった砕波等の流動解析を短時間で高精度に行う流動解析装置及び流動解析方法が記載されている。特許文献１に記載の流動解析方法においても、流体を収容する容器（すなわち固定壁）は、座標位置を固定して流速０の複数粒子の集まりとして表現する旨が記載されており、従来手法１を利用したものである。
【０００５】
また、特許文献２には、汎用性が高く、濡れ性モデルを統一的に取り扱うことができ、空間解像度が粗くても精度よく解析を行うことができる流体解析装置及び流体解析方法が記載されている。この特許文献２における流体解析方法は、表面張力モデルにおいてエネルギ関係式を用いて壁の粒子からの作用力を算出する点に特徴があるが、壁を固定した複数粒子により表現する点においては特許文献１と同様に一般的な手法（従来手法１）を用いている。
【０００６】
壁境界を表現する別の手法として、壁境界付近の流体粒子を境界に対して対象となるように複製することにより、動的に壁粒子を生成する手法がある。図９は、従来の流動解析方法における粒子法による境界条件を境界付近の流体粒子を複製して生成した壁粒子により表現する場合の概念図である。図９において、四角の枠線（境界）内に記載された複数の粒子２７ａは、流体を示すものとする。また、四角の枠線（境界）外に記載された複数の粒子２７ｂは、壁粒子を示すものであり、境界に対して対称な位置の粒子２７ａを複製して生成されたものである。このような手法を従来手法２と呼ぶこととする。
【０００７】
従来手法２は、境界に流体粒子が隣接する場合にのみ壁粒子が生成されるため、流体の動きに応じて壁粒子の分布状況も動的に変化する。すなわち、従来手法２は、境界に流体粒子が隣接しない部分において壁粒子を生成する必要が無く、計算コストを抑えられる（計算時間が短い）という利点を有する。なお、従来手法２は、非特許文献１に記載されている。
【０００８】
さらに、壁境界を表現する別の手法として、壁粒子を一切用いずに壁からの反発力を計算する手法がある。図１０は、従来の流動解析方法における粒子法による境界条件を壁からの反発力の存在により表現する場合の概念図である。図１０において、四角の枠線（境界）内に記載された複数の粒子は、流体を示すものとする。図１０に示すように、境界内の複数の流体粒子は、壁に接触した際に壁に対してめり込んだ量に応じた壁からの反発力を受ける。このような手法を従来手法３と呼ぶこととする。
【０００９】
従来手法３は、壁に接触した流体粒子がペナルティ法で計算された壁境界からの反発力を受けるものとして行う解析手法であるが、壁粒子を一切用いないため、計算コストを抑えられる（計算時間が短い）という利点を有する。また、従来手法３は、壁境界の形状として図１０に示す矩形形状のような単純な形状のみならず、任意形状（ＣＡＤデータ）を取り扱うことが可能であるといった利点も有する。なお、従来手法３は、非特許文献２，３に記載されている。
【先行技術文献】
【特許文献】
【００１０】
【特許文献１】特開平７−３３４４８４号公報
【特許文献２】特開２００８−１１１６７５号公報
【非特許文献】
【００１１】
【非特許文献１】Ｓ．Ｊ．Ｃｕｍｍｉｎｓ ａｎｄ Ｍ．Ｒｕｄｍａｎ，Ａｎ ＳＰＨ ｐｒｏｊｅｃｔｉｏｎ ｍｅｔｈｏｄ，Ｊｏｕｒｎａｌ ｏｆ Ｃｏｍｐｕｔａｔｉｏｎａｌ Ｐｈｙｓｉｃｓ，１５２，５８４−６０７（１９９９）．
【非特許文献２】原田隆宏，越塚誠一，ＳＰＨにおける壁境界計算手法の改良，情報処理学会論文誌，４８，４，１８３８−１８４６（２００７）．
【非特許文献３】原田隆宏，越塚誠一，島崎克教，ＭＰＳ法における壁境界モデルの改良、日本計算工学会論文集，２００８，Ｐａｐｅｒ Ｎｏ．２００８０００６（２００８）．
【発明の概要】
【発明が解決しようとする課題】
【００１２】
しかしながら、上述した従来手法１は、等しい大きさの複数の粒子からなる壁粒子を格子状に規則正しく配置して壁境界を表現するため、複雑形状に不適切であるという欠点がある。この問題点について図１１を用いて説明する。図１１は、従来の流動解析方法において従来手法１により境界条件を表現する際の問題点を示す概念図である。図１１に示すように、境界が例えば滑らかな曲面を描いた円形形状等である場合に、格子状に規則正しく配置した壁粒子は、滑らかなカーブを表現することが困難であり、境界がいびつになるという問題点がある。壁粒子の大きさを小さなものにすれば滑らかなカーブに近い形状を表現することも可能であるが、当然のことながら多くの壁粒子を必要とし、計算コストは高く（計算時間は長く）なる。さらに、そもそも従来手法１は、従来手法２，３に比して、壁部分全体について壁粒子を配置する必要があるため、多くの壁粒子を必要とし、計算コストが高いという欠点を有する。
【００１３】
従来手法２は、壁境界付近の流体粒子を境界に対して対称となるように複製して壁粒子とするものであるため、コーナー部や複雑形状に対応が困難であり汎用性が低いという欠点を有する。例えば図９に示すように、従来手法２による壁粒子は、直線的な壁（境界）を模擬することは容易であるが、コーナー部には壁が存在しないこととなってしまい、所望する解析結果が得られないことも考えられる。
【００１４】
図１２は、従来手法２の方法を発展させたものであり、境界に対して対称となるように流体粒子を１度複製して壁粒子を生成した後に、さらに別の境界に対して対称となる位置に壁粒子を複製し、コーナー部の隙間を埋めるように壁粒子の配置を試みた場合の概念図である。図１２に示すように、流体粒子２８ａは、境界に対称な位置に壁粒子２８ｂを複製させる。さらに、壁粒子２８ｂは、別の境界に対称な位置に壁粒子２８ｃを複製させることにより、境界を表現する壁のコーナー部を形成する。同様に、流体粒子２９ａは、境界に対称な位置に壁粒子２９ｂを複製させる。さらに、壁粒子２９ｂは、別の境界に対称な位置に壁粒子２９ｃを複製させる。
【００１５】
しかしながら、図１２に示す手法は、左下に示すように１つの流体粒子３０が境界に対称な位置に壁粒子３１，３２を複製させる場合が考えられる。この場合に、壁粒子３１，３２の各々により境界に対称な位置に複製された壁粒子は、複製先において位置が重なってしまうという現象が生じうる。したがって、図１２に示す手法は、同じ位置に複数の壁粒子が存在する場合があり、流動解析の結果が不正確になる可能性がある。
【００１６】
従来手法３は、上述したように壁粒子を必要としないため、計算コストを抑えるとともに任意形状にも対応可能であるという利点を有するが、ペナルティ法により計算された反発力が壁に対して法線方向の力のみであるという欠点があり、計算精度が悪いという問題点がある。実際の固定壁は圧力勾配を有し、流体に対して法線方向のみならず接線方向の力も作用するからである。したがって、従来手法３を用いて行った流動解析は、計算精度が悪いため、場合によっては流体の挙動自体に悪影響を及ぼす可能性も考えられ、信頼性が低いという欠点がある。
【００１７】
本発明は上述した従来技術の問題点を解決するもので、計算コストを抑えるとともに固定壁の境界条件について任意形状に対応可能であり、さらに計算精度の良い流動解析方法、流動解析装置、及び流動解析プログラムを提供することを課題とする。
【課題を解決するための手段】
【００１８】
本発明に係る流動解析方法は、上記課題を解決するために、粒子法を用いて流体の挙動を解析する流動解析方法であって、流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する距離判定ステップと、前記距離判定ステップにより前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の前記壁境界の曲率を算出する曲率算出ステップと、前記距離判定ステップにより前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、前記壁境界を挟んで対称に壁を表す複数の第２粒子を前記曲率算出ステップにより算出された曲率に応じた拡大率で複製して配置する壁粒子配置ステップと、前記壁粒子配置ステップにより配置された前記複数の第２粒子を境界条件として粒子間の相互作用に基づき流体の挙動を解析する解析ステップとを備えることを特徴とする。
【発明の効果】
【００１９】
本発明によれば、計算コストを抑えるとともに固定壁の境界条件について任意形状に対応可能であり、さらに計算精度の良い流動解析を行うことができる。
【図面の簡単な説明】
【００２０】
【図１】本発明の実施例１の形態の流動解析装置の構成を示すブロック図である。
【図２】本発明の実施例１の形態の流動解析装置における壁境界計算部の詳細な構成を示すブロック図である。
【図３】本発明の実施例１の形態の流動解析装置の動作を示すフローチャート図である。
【図４】本発明の実施例１の形態の流動解析方法における格子幅が大きい場合の法線方向を示す模式図である。
【図５】本発明の実施例１の形態の流動解析方法における格子幅が小さい場合の法線方向を示す模式図である。
【図６】本発明の実施例１の形態の流動解析方法における粒子法による境界条件を境界付近の流体粒子を曲率に応じた拡大率で複製した壁粒子により表現する場合の概念図である。
【図７】本発明の実施例１の形態の流動解析方法における粒子法による境界条件を境界付近の流体粒子を曲率に応じた拡大率で複製した壁粒子により表現する場合の概念図である。
【図８】従来の流動解析方法における粒子法による壁粒子を用いた境界条件を示す概念図である。
【図９】従来の流動解析方法における粒子法による境界条件を境界付近の流体粒子を複製して生成した壁粒子により表現する場合の概念図である。
【図１０】従来の流動解析方法における粒子法による境界条件を壁からの反発力の存在により表現する場合の概念図である。
【図１１】従来の流動解析方法において従来手法１により境界条件を表現する際の問題点を示す概念図である。
【図１２】従来の流動解析方法において従来手法２を発展させた方法により境界条件を表現する際の問題点を示す概念図である。
【発明を実施するための形態】
【００２１】
以下、本発明の流動解析方法、流動解析装置、及び流動解析プログラムの実施の形態について図面を参照して詳細に説明する。
【実施例１】
【００２２】
図１は、本発明の実施例１の流動解析装置１の構成を示すブロック図である。図１を参照して、流動解析装置１の構成を説明する。本発明の実施形態となる流動解析装置１は、粒子法を用いて流体の挙動を解析する装置であり、パーソナルコンピュータ，ワークステーション，汎用コンピュータ等の公知の情報処理装置により構成され、図１に示すように、ＲＡＭ（Ｒａｎｄｏｍ Ａｃｃｅｓｓ Ｍｅｍｏｒｙ）２、ＲＯＭ（Ｒｅａｄ Ｏｎｌｙ Ｍｅｍｏｒｙ）３、初期設定データベース（ＤＢ）４、及びＣＰＵ（Ｃｅｎｔｒａｌ Ｐｒｏｃｅｓｓｉｎｇ Ｕｎｉｔ）５を主な構成要素として備える。ＲＡＭ２は、ＣＰＵ５が実行する流動解析装置１の起動プログラムや流動解析プログラム等のコンピュータプログラム及び各種データを一時的に格納するワーキングエリアを提供する。ＲＯＭ３は、コンピュータプログラム及びこのコンピュータプログラムを実行するために必要な各種データを記憶する。なおＲＯＭ３は、磁気的又は光学的な記憶媒体若しくは半導体メモリ等のＣＰＵ５により読み取り可能な記憶媒体を含む構成を有し、ＲＯＭ３に記憶されているコンピュータプログラムや各種データの一部又は全部は電気通信回線を介してダウンロードされるように構成してもよい。
【００２３】
初期設定ＤＢ４は、粒子法により流動解析を実行するために必要な各種パラメータを格納する。具体的には、このパラメータには、粒子法における粒子ｉの粒子径ｒｉの初期値ｒiini、位置（重心の座標）、速度、温度、圧力、流体の物性（例えば粘度、密度、熱伝導率、比熱、潜熱等）、時間刻み幅Δｔ、計算終了時間Ｔ等の初期条件に関する情報が含まれる。ＣＰＵ５は、ＲＯＭ３に記憶されているコンピュータプログラムをＲＡＭ２内にロードし、読み出されたコンピュータプログラムを実行することにより、流動解析装置１全体の処理動作を制御する。具体的には、ＣＰＵ５は、ＲＯＭ３に記憶された流動解析プログラム及び流動解析プログラムを実行するために必要なデータをＲＡＭ２へロードし、流動解析プログラムに従って粒子法を利用した流動解析を実行する。本実施形態では、ＣＰＵ５は、ＲＯＭ３に記憶されている流動解析プログラムを実行することにより、壁境界計算部１０、外力項計算部１１、粘性項計算部１２、粒子移動追跡部１３、圧力項計算部１４、粒子データ更新部１５、及び終了判定部１６として機能する。壁境界計算部１０以外の各部における機能の詳細については後述する。
【００２４】
図２は、本実施例の流動解析装置１におけるＣＰＵ５内部の壁境界計算部１０の詳細な構成を示すブロック図であり、図２に示すように距離判定部１７、曲率算出部１８、及び壁粒子配置部１９から構成される。
【００２５】
距離判定部１７は、流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する。この「所定距離」は、予め距離判定部１７に設定されているものでもよいし、ユーザにより設定可能としてもよい。
【００２６】
曲率算出部１８は、距離判定部１７により壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の壁境界の曲率を算出する。ここで、曲率算出部１８は、予め壁境界からの距離を定めた距離関数に基づいて壁境界の曲率を算出する。具体的な算出方法は後述する。
【００２７】
壁粒子配置部１９は、距離判定部１７により壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、壁境界を挟んで対称に壁を表す複数の第２粒子を曲率算出部１８により算出された曲率に応じた拡大率で複製して配置する。
【００２８】
また、本実施例における壁粒子配置部１９は、配置した複数の第２粒子（壁粒子）の各々が有する圧力値を対応する第１粒子が有する圧力値と同じ値に設定する。
【００２９】
壁境界計算部１０を除くＣＰＵ５（すなわち外力項計算部１１、粘性項計算部１２、粒子移動追跡部１３、圧力項計算部１４、粒子データ更新部１５、及び終了判定部１６）は、本発明の解析部に対応し、壁粒子配置部１９により配置された複数の第２粒子（壁粒子）を境界条件として粒子間の相互作用に基づき流体の挙動を解析する。
【００３０】
次に、上述のように構成された本実施の形態の作用を説明する。図３は、本実施例の流動解析装置１の動作を示すフローチャート図である。流動解析装置１は、以下に示す流動解析処理を実行することにより流体の挙動を解析する。以下、図３に示すフローチャートを参照して、この流動解析処理を実行する際の流動解析装置１の動作及び本発明の流動解析方法について詳しく説明する。
【００３１】
図３に示すフローチャートは、流動解析装置１に対しオペレータが流動解析処理の開始を指示したタイミングで開始となり、流動解析処理はステップＳ１の処理に進む。
【００３２】
ステップＳ１の処理では、ＣＰＵ５が、初期設定ＤＢ４からＲＡＭ２内に解析対象となる流体の物性，時間刻み幅ΔＴ，計算終了時間Ｔ等の計算条件パラメータ、及び粒子の座標位置，速度，温度等の初期条件パラメータを読み出し、読み出されたパラメータに基づいて流体を複数の粒子により表現する。これにより、ステップＳ１の処理は完了し、流動解析処理はステップＳ２の処理に進む。
【００３３】
ステップＳ２の処理は、壁境界の計算を行うものである。まず、壁境界計算部１０内の距離判定部１７は、流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する（距離判定ステップ）。これは、各粒子について壁粒子として複製するか否かを判断するための処理であり、壁境界から離れた位置にある粒子は、壁粒子として複製されない。
【００３４】
を計算する。ここで、ｉは、粒子番号を示すものとする。本実施例において、この距離の計算は、空間に距離関数を定義することにより行われる。すなわち、距離判定部１７は、計算空間を格子状に区切り、各格子点において最も近い壁境界までの距離をあらかじめ計算した上で記憶しておく。粒子位置における壁までの最小距離は、その距離関数を補間することにより求められる。ここで、距離関数について、詳細に説明する。
【００３５】
本実施例で使用する距離関数は、ある任意の境界面からの距離を空間上に定義したものである。一般的には、計算領域を覆うようにして正方格子状に領域を分割し、その各格子点において壁境界からの距離を計算したものを用いる。壁の外側を正の値、内側を負の値とする。壁境界が動かない場合、距離関数は変化することがないため、前処理により計算しておくことで解析の計算コストを大幅に下げることが可能になる。
【００３６】
【数１】
【００３７】
【００３８】
より計算コストを抑えることができる。ここで、粒子ｉは、以下に示す８点で作られる立方体の中に存在しているとする。１：（Ｉ，Ｊ，Ｋ）、２：（Ｉ＋１，Ｊ，Ｋ）、３：（Ｉ＋１，Ｊ＋１，Ｋ）、４：（Ｉ，Ｊ＋１，Ｋ）、５：（Ｉ，Ｊ，Ｋ＋１）、６：（Ｉ＋１，Ｊ，Ｋ＋１）、７：（Ｉ＋１，Ｊ＋１，Ｋ＋１）、８：（Ｉ，Ｊ＋１，Ｋ＋１）。
【００３９】
格子点（Ｉ，Ｊ，Ｋ）は、（２）式により表される。
【数２】
【００４０】
【数３】
【００４１】
このようにして流体を表す複数の第１粒子の各々と壁境界との間の距離を求めた距離判定部１７は、当該距離が所定距離以下であるか否かを判定し、判定結果を出力する。
【００４２】
次に、曲率算出部１８は、距離判定部１７により壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の壁境界の曲率を算出する（曲率算出ステップ）。
【００４３】
【数４】
【００４４】
この法線方向は、粒子を複製する方向となる。次に、曲率算出部１８は、法線方向の発散を計算することにより、（５）式に示すように曲率ｋを求めることができる。曲率は、複製後の粒子の大きさを決定するのに用いられる。
【数５】
【００４５】
曲率算出部１８は、上述したように距離関数を利用して法線方向及び曲率を算出している。したがって、曲率算出部１８は、曲率算出ステップにおいて予め壁境界からの距離を定めた距離関数に基づいて壁境界の曲率を算出するといえる。
【００４６】
法線方向と曲率に関しては，粒子上の距離関数を用いて計算する方法と、格子点上において前処理により計算しておき、それを粒子上に線型補間する方法の２つが考えられる。１つ目の方法では従来の粒子法の定式化より、空間微分と発散を計算することにより求めることができる。２つ目の方法では以下のように法線方向と曲率を前計算しておく。格子点（Ｉ，Ｊ，Ｋ）における法線方向は、（６）式により表される。
【数６】
【００４７】
また、（６）式における〈∇ｄ〉Ｉ，Ｊ，Ｋは、（７）式のように表される。
【数７】
【００４８】
さらに、格子点（Ｉ，Ｊ，Ｋ）における曲率ｋ（Ｉ，Ｊ，Ｋ）は、（８）式のように表される。
【数８】
【００４９】
これらの格子点上の値を補間することにより、粒子上の値が計算される。
【００５０】
次に、法線方向と解像度の関係について述べる。上述した（６）式及び（７）式を使用して法線方向を算出する際に、角近辺における法線方向は、解像度すなわち格子の幅の影響を受ける。図４は、格子幅が大きい場合の法線方向を示す模式図である。また、図５は、格子幅が小さい場合の法線方向を示す模式図である。図４及び図５における太線は、壁境界面である。また、矢印は、法線方向を表す。ここで、点線において囲った点の法線方向を比較する。格子幅が大きい場合における法線方向は、図４に示すように斜め方向を指す。一方、格子幅が小さい場合における法線方向は、図５に示すように真横方向を指す。これは、法線方向を（６）式及び（７）式を用いて計算する際に、隣の距離関数を参照するため、角の近辺にある格子点では法線方向が斜め方向に算出されるためである。したがって、格子幅が大きいほど法線方向が斜め方向に計算される領域は広くなる。この格子幅を粒子径と同程度にすることにより、後述する複製された壁粒子が重ならないような法線方向が計算される。
【００５１】
次に、壁粒子配置部１９は、距離判定部１７により壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、壁境界を挟んで対称に壁を表す複数の第２粒子を複製する。ここで、壁粒子配置部１９は、複数の第２粒子を複製する際に、曲率算出部１８により算出された曲率に応じた拡大率で複製して配置する（壁粒子配置ステップ）。
【００５２】
粒子ｉの拡大率ｒｉは、粒子上の曲率ｋｉ及び粒子径Ｌｉを用いて、（９）式のように表される。
【数９】
【００５３】
【数１０】
【００５４】
壁粒子配置部１９は、（１１）式に示すように、拡大率が大きい場合には壁境界から離れた位置に複製した粒子を配置し、拡大率が小さい場合には壁境界に近い位置に複製した粒子を配置することにより、隙間の少ない壁を形成させる。
【００５５】
なお、壁粒子配置部１９は、複製した粒子を壁粒子（第２粒子）として配置する際に、配置した複数の第２粒子の各々が有する圧力値を対応する第１粒子（複製元の粒子）が有する圧力値と同じ値に設定する。設定した圧力値は、ステップＳ７で圧力項を計算する際に使用される。
【００５６】
図６及び図７は、本実施例の流動解析方法における粒子法による境界条件を境界付近の流体粒子を曲率に応じた拡大率で複製した壁粒子により表現する場合の概念図である。
【００５７】
壁粒子配置部１９は、例えば図６に示すように、壁境界との間の距離が所定距離以下であると判断された流体の粒子２０ａに対応して、壁境界を挟んで対称に壁粒子２０ｂを複製して配置する。この場合において、粒子２０ａ近傍の壁境界は直線的な壁であるため、曲率ｋ＝０である。したがって、拡大率ｒ＝１となり、壁粒子配置部１９は、粒子２０ａと同じ大きさの粒子２０ｂを複製し、壁粒子として配置する。
【００５８】
また、壁粒子配置部１９は、図６に示すように、壁境界との間の距離が所定距離以下であると判断された流体の粒子２１ａ，２２ａに対応して、壁境界を挟んで対称に壁粒子２１ｂ，２２ｂを複製して配置する。この場合において、粒子２１ａ，２２ａ近傍の壁境界は曲率が負（凹形状）であるため、拡大率ｒ＞１となり、壁粒子配置部１９は、粒子２１ａ，２２ａを大きく拡大して粒子２１ｂ，２２ｂを複製し、壁粒子として配置する。すなわち、曲率が負の場合には、流体の領域に比して壁の領域が広いため、壁粒子配置部１９は、壁の隙間を無くすように大きな壁粒子を配置する。
【００５９】
また、壁粒子配置部１９は、図７に示すように、壁境界との間の距離が所定距離以下であると判断された流体の粒子２３ａ，２４ａに対応して、壁境界を挟んで対称に壁粒子２３ｂ，２４ｂを複製して配置する。この場合において、粒子２３ａ，２４ａ近傍の壁境界は曲率が正（凸形状）であるため、拡大率ｒ＜１となり、壁粒子配置部１９は、粒子２３ａ，２４ａを小さく縮小して粒子２３ｂ，２４ｂを複製し、壁粒子として配置する。すなわち、曲率が正の場合には、流体の領域に比して壁の領域が狭いため、壁粒子配置部１９は、壁の隙間を無くすように小さな壁粒子を配置する。
【００６０】
なお、壁粒子配置部１９は、第１粒子と壁境界を挟んで対称に壁を表す第２粒子を複製するものであるが、ここで言う「対称」とは、算出した「法線方向」上に壁境界を挟んで「対称」の意である。したがって、図６に示す粒子２２ａ，２２ｂの組や、図７に示す粒子２４ａ，２４ｂの組は、一見、点線で示した壁境界に対して対称でないように見えるが、算出した法線方向上に対称に配置されており、法線方向が壁境界に対して斜め方向となる理由は図４，５を用いて説明したとおりである。また、壁粒子から壁境界までの距離は、拡大率に応じて変化する。
【００６１】
壁境界の計算が終わり、複数の第２粒子が壁粒子として配置されると、ステップＳ２の処理は完了し、流動解析処理はステップＳ３の処理に進む。
【００６２】
ステップＳ３以降のステップにおいて、本発明の流動解析装置は、壁粒子配置ステップにより配置された複数の第２粒子を境界条件として粒子間の相互作用に基づき流体の挙動を解析する（解析ステップ）。すなわち、ステップ３以降の動作は、全て解析ステップである。
【００６３】
ステップＳ３の処理では、外力項計算部１１が、重力等の各粒子に直接的に作用する外力項を計算する。粒子法における外力項の計算方法は本願発明の出願時点で既に公知であるので詳細な説明は省略する。これにより、ステップＳ３の処理は完了し、流動解析処理はステップＳ４の処理に進む。
【００６４】
ステップＳ４の処理では、粘性項計算部１２が、粘性項を計算し、粒子ｉが有する速度分布を近傍の粒子ｊに分配することにより、粘性による粒子の速度の拡散を計算する。具体的には、粘性項は以下の（１２）式により表される。
【数１１】
【００６５】
（１２）式中におけるパラメータｕ、ｔ、及びνは、それぞれ速度ベクトル、計算時刻、動粘性係数を示す。この（１２）式を粒子法で離散化すると、粘性項は陽的には（１３）式、陰的には（１４）式のように表される。
【数１２】
【００６６】
（１３）式及び（１４）式においてパラメータｄ、λ、ｎ０、ｗｉｊ、及びｋは、それぞれ次元数、拡散を解析解に一致させるための係数、基準粒子密度、粒子ｉ，ｊ間の重み、及びタイムステップを示す。
【００６７】
粒子ｉの速度変化Δｕは、以下の（１５）式により表されるので、これを（１４）式に代入することにより、以下の（１６）式に示す連立一次方程式が導かれる。粘性項計算部１２は、この（１６）式を解くことにより粘性項を陰的に解く。
【数１３】
【００６８】
このように粘性項を陰的に解くことにより、タイムステップｋを大きくすることができるので、発散することなく計算を安定的に行うことができる。これにより、このステップＳ４の処理は完了し、流動解析処理はステップＳ５の処理に進む。
【００６９】
ステップＳ５の処理では、粒子移動追跡部１３が、ステップＳ３及びステップＳ４の計算結果に基づいて粒子の速度を仮更新する。これにより、ステップＳ５の処理は完了し、流動解析処理はステップＳ６の処理に進む。
【００７０】
ステップＳ６の処理では、粒子移動追跡部１３が、ステップＳ５の処理により仮更新された速度による時間刻み幅Δｔの間での粒子の移動を追跡し、追跡結果に基づいて粒子の座標位置を仮更新する。これにより、ステップＳ６の処理は完了し、流動解析処理はステップＳ７の処理に進む。
【００７１】
ステップＳ７の処理では、圧力項計算部１４が、ステップＳ６の処理により算出された座標位置における各粒子の粒子数密度を算出し、算出された粒子数密度が所定値であるか否かを判別する。判別の結果、粒子数密度が所定値でない場合、圧力項計算部１４は、粒子数密度が所定値になるように粒子の圧力を修正することにより粒子の座標位置及び流速を修正する。そして圧力計算部１４は、修正された圧力によって生じる圧力項を算出する。なお、圧力計算部１４は、壁粒子の圧力については、ステップＳ２の壁粒子配置ステップにおいて設定された圧力値（すなわち複製元の流体粒子が有する圧力値と同じ圧力値）を用いて計算を行う。これにより、ステップＳ７の処理は完了し、流動解析処理はステップＳ８の処理に進む。
【００７２】
ステップＳ８の処理では、粒子データ更新部１５が、ステップＳ７の処理の算出結果に基づいてＲＡＭ２内に格納されている粒子の速度を更新する。これにより、ステップＳ８の処理は完了し、流動解析処理はステップＳ９の処理に進む。
【００７３】
ステップＳ９の処理では、粒子データ更新部１５が、ステップＳ７の処理の算出結果に基づいてＲＡＭ２内に格納されている粒子の座標位置を更新する。これにより、ステップＳ９の処理は完了し、流動解析処理はステップＳ１０の処理に進む。
【００７４】
ステップＳ１０の処理では、終了判定部１６が、解析時間ｔを時間刻み幅ΔＴだけインクリメントする。これにより、ステップＳ１０の処理は完了し、流動解析処理はステップＳ１１の処理に進む。
【００７５】
ステップＳ１１の処理では、終了判定部１６が、解析時間ｔが計算終了時刻Ｔであるか否かを判別する。判別の結果、解析時間ｔが計算終了時刻Ｔである場合、終了判定部１６は、ＲＡＭ２内に格納されている粒子データをＲＯＭ３の所定ファイルにコピーし、一連の流動解析処理を終了する。一方、解析時間ｔが計算終了時刻Ｔでない場合には、終了判定部１０は流動解析処理をステップＳ２の処理に戻す。
【００７６】
上述のとおり、本発明の実施例１の形態に係る流動解析装置１（流動解析方法）によれば、計算精度が良く、且つ計算コストを抑えるとともに固定壁の境界条件について任意形状（ＣＡＤデータ）に対応することができる。
【００７７】
すなわち、本発明の流動解析装置１は、等しい大きさの複数の粒子からなる壁粒子を格子上に規則正しく配置して壁境界を表現する従来手法１と異なり、壁境界の曲率に応じた拡大率で境界面付近の流体粒子を複製して壁粒子とするため、任意形状に対応可能という利点を有する。同様に、コーナー部や複雑形状に対応が困難な従来手法２と異なり、本発明の流動解析装置１は、コーナー部等にも隙間無く壁粒子を配置することができ、任意形状に対応可能であるのみならず、コーナー部等の壁粒子の影響も考慮して計算精度の良い解析結果を得ることができる。
【００７８】
さらに、本発明の流動解析装置１は、従来手法２と同様に境界に流体粒子が隣接しない部分において壁粒子を生成する必要が無く、計算コストを抑えられるという利点を有する。
【００７９】
また、ペナルティ法により計算された反発力が壁に対して法線方向の力のみであると従来手法３と異なり、本発明の流動解析装置１は、複製した壁粒子を配置することにより固定壁の圧力勾配も考慮に入れた計算を行うため、壁から流体に対する様々な方向の力を考慮に入れた計算精度の高い解析結果を得ることができる。
【００８０】
さらに、本発明の流動解析装置１は、複製した壁粒子を配置する際に、配置した複数の第２粒子（壁粒子）の各々が有する圧力値を対応する第１粒子（流体粒子）が有する圧力値と同じ値に設定するので、壁粒子の圧力について計算する必要が無く、計算コストを抑えることができるという利点を有する。
【００８１】
また、本発明の流動解析装置１は、壁の曲率を計算するために距離関数を使用しているので、予め各格子点における壁境界からの距離を計算したものを利用することができ、解析の計算コストを大幅に下げることが可能である。
【００８２】
また本実施形態における流動解析装置１は、モデム等の通信用部材を備えてインターネットを含む電気通信回線と接続可能なように構成してもよい。ここで、本発明の流動解析プログラムは、粒子法を用いて流体の挙動を解析する流動解析プログラムであって、流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する処理と、壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の壁境界の曲率を算出する処理と、壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、壁境界を挟んで対称に壁を表す複数の第２粒子を算出した壁境界の曲率に応じた拡大率で複製して配置する処理と、配置された複数の第２粒子を境界条件として粒子間の相互作用に基づき流体の挙動を解析する処理とをコンピュータに実行させるプログラムである。
【００８３】
当該流動解析プログラムは、電気通信回線からダウンロード等によって、ＲＡＭ２又はＲＯＭ３内に格納させてもよい。なお、この場合における通信ネットワークからダウンロードするためのダウンロードプログラムは、予めＲＯＭ３内に格納されているか、別の記憶媒体からＲＯＭ３内にインストールされるものとする。このように、本実施形態に基づいて当業者等によりなされる他の実施の形態、実施例及び運用技術等は全て本発明の範疇に含まれる。
【産業上の利用可能性】
【００８４】
本発明に係る流動解析方法、流動解析装置、及び流動解析プログラムは、粒子法を用いて連続体である流体の挙動を解析する流動解析装置等に利用可能である。
【符号の説明】
【００８５】
１…流動解析装置、２…ＲＡＭ、３…ＲＯＭ、４…初期設定ＤＢ、５…ＣＰＵ、１０…壁境界計算部、１１…外力項計算部、１２…粘性項計算部、１３…粒子移動追跡部、１４…圧力項計算部、１５…粒子データ更新部、１６…終了判定部、１７…距離判定部、１８…曲率算出部、１９…壁粒子配置部、２０ａ，２０ｂ，２１ａ，２１ｂ，２２ａ，２２ｂ，２３ａ，２３ｂ，２４ａ，２４ｂ，２５，２６，２７ａ，２７ｂ，２８ａ，２８ｂ，２８ｃ，２９ａ，２９ｂ，２９ｃ，３０，３１，３２…粒子。
【特許請求の範囲】
【請求項１】
粒子法を用いて流体の挙動を解析する流動解析方法であって、
流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する距離判定ステップと、
前記距離判定ステップにより前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の前記壁境界の曲率を算出する曲率算出ステップと、
前記距離判定ステップにより前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、前記壁境界を挟んで対称に壁を表す複数の第２粒子を前記曲率算出ステップにより算出された曲率に応じた拡大率で複製して配置する壁粒子配置ステップと、
前記壁粒子配置ステップにより配置された前記複数の第２粒子を境界条件として粒子間の相互作用に基づき流体の挙動を解析する解析ステップと、
を備えることを特徴とする流動解析方法。
【請求項２】
前記壁粒子配置ステップは、配置した前記複数の第２粒子の各々が有する圧力値を対応する前記第１粒子が有する圧力値と同じ値に設定することを特徴とする請求項１記載の流動解析方法。
【請求項３】
前記曲率算出ステップは、予め前記壁境界からの距離を定めた距離関数に基づいて前記壁境界の曲率を算出することを特徴とする請求項１又は請求項２記載の流動解析方法。
【請求項４】
粒子法を用いて流体の挙動を解析する流動解析装置であって、
流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する距離判定部と、
前記距離判定部により前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の前記壁境界の曲率を算出する曲率算出部と、
前記距離判定部により前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、前記壁境界を挟んで対称に壁を表す複数の第２粒子を前記曲率算出部により算出された曲率に応じた拡大率で複製して配置する壁粒子配置部と、
前記壁粒子配置部により配置された前記複数の第２粒子を境界条件として粒子間の相互作用に基づき流体の挙動を解析する解析部と、
を備えることを特徴とする流動解析装置。
【請求項５】
粒子法を用いて流体の挙動を解析する流動解析プログラムであって、
流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する処理と、
前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の前記壁境界の曲率を算出する処理と、
前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、前記壁境界を挟んで対称に壁を表す複数の第２粒子を算出した前記壁境界の曲率に応じた拡大率で複製して配置する処理と、
配置された前記複数の第２粒子を境界条件として粒子間の相互作用に基づき流体の挙動を解析する処理と、
をコンピュータに実行させることを特徴とする流動解析プログラム。
【請求項１】
粒子法を用いて流体の挙動を解析する流動解析方法であって、
流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する距離判定ステップと、
前記距離判定ステップにより前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の前記壁境界の曲率を算出する曲率算出ステップと、
前記距離判定ステップにより前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、前記壁境界を挟んで対称に壁を表す複数の第２粒子を前記曲率算出ステップにより算出された曲率に応じた拡大率で複製して配置する壁粒子配置ステップと、
前記壁粒子配置ステップにより配置された前記複数の第２粒子を境界条件として粒子間の相互作用に基づき流体の挙動を解析する解析ステップと、
を備えることを特徴とする流動解析方法。
【請求項２】
前記壁粒子配置ステップは、配置した前記複数の第２粒子の各々が有する圧力値を対応する前記第１粒子が有する圧力値と同じ値に設定することを特徴とする請求項１記載の流動解析方法。
【請求項３】
前記曲率算出ステップは、予め前記壁境界からの距離を定めた距離関数に基づいて前記壁境界の曲率を算出することを特徴とする請求項１又は請求項２記載の流動解析方法。
【請求項４】
粒子法を用いて流体の挙動を解析する流動解析装置であって、
流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する距離判定部と、
前記距離判定部により前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の前記壁境界の曲率を算出する曲率算出部と、
前記距離判定部により前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、前記壁境界を挟んで対称に壁を表す複数の第２粒子を前記曲率算出部により算出された曲率に応じた拡大率で複製して配置する壁粒子配置部と、
前記壁粒子配置部により配置された前記複数の第２粒子を境界条件として粒子間の相互作用に基づき流体の挙動を解析する解析部と、
を備えることを特徴とする流動解析装置。
【請求項５】
粒子法を用いて流体の挙動を解析する流動解析プログラムであって、
流体を表す複数の第１粒子の各々と壁境界との間の距離が所定距離以下であるか否かを判定する処理と、
前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対して最も近い位置の前記壁境界の曲率を算出する処理と、
前記壁境界との間の距離が所定距離以下であると判断された複数の第１粒子の各々に対応して、前記壁境界を挟んで対称に壁を表す複数の第２粒子を算出した前記壁境界の曲率に応じた拡大率で複製して配置する処理と、
配置された前記複数の第２粒子を境界条件として粒子間の相互作用に基づき流体の挙動を解析する処理と、
をコンピュータに実行させることを特徴とする流動解析プログラム。
【図１】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【図１２】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【図１０】
【図１１】
【図１２】
【公開番号】特開２０１０−２４３２９３（Ｐ２０１０−２４３２９３Ａ）
【公開日】平成２２年１０月２８日（２０１０．１０．２８）
【国際特許分類】
【出願番号】特願２００９−９１３３９（Ｐ２００９−９１３３９）
【出願日】平成２１年４月３日（２００９．４．３）
【出願人】（０００００３０７８）株式会社東芝 (54,554)
【公開日】平成２２年１０月２８日（２０１０．１０．２８）
【国際特許分類】
【出願日】平成２１年４月３日（２００９．４．３）
【出願人】（０００００３０７８）株式会社東芝 (54,554)
[ Back to top ]