(57)【要約】  物理プロセスのシミュレーションは、ボクセルのための状態ベクトルであって、モデルに対応し、ボクセルで生じ得る運動量状態である特定の運動量状態に対応するエントリーを具備する状態ベクトルをコンピュータアクセス可能メモリ内に記憶することを含む。作用操作は、状態ベクトルで実行される。作用操作は、モデルにしたがって、異なる運動量状態の要素間の作用をモデル化する。状態ベクトルで実行される移動操作は、モデルに従い、要素の動きを新たなボクセルに反映する。モデルは、０．１より大きいクヌーセン数である高クヌーセン数流れをシミュレーションするよう適用されている。
2010500654.tif 000002
【特許請求の範囲】
【請求項１】
コンピュータが読み取り可能な媒体で実現され、実行時に物理プロセスをシミュレーションするコンピュータプログラムプロダクトであって、コンピュータに
ボクセルのための状態ベクトルであって、あるモデルに対応し、ボクセルで生じ得る運動量状態である特定の運動量状態に対応する複数のエントリーを具備する状態ベクトルをコンピュータがアクセス可能なメモリ内に記憶することと、
前記モデルに従って異なる運動量状態の要素間の作用をモデリングする作用操作を前記状態ベクトルに実行することと、
移動操作を前記状態ベクトルに実行し、前記モデルに従って要素の動きを新たなボクセルに反映させることと、
を実行させるよう構成されているコンピュータプログラムプロダクトであり、
前記モデルは、クヌーセン数が０．１より大きい高クヌーセン数流れをシミュレーションするよう適用されていることを特徴とするコンピュータプログラムプロダクト。
【請求項２】
請求項１に記載のコンピュータプログラムプロダクトであって、コンピュータが少なくとも一つの表面の表現を前記メモリ内に記憶するようさらに構成されていることを特徴とするコンピュータプログラムプロダクト。
【請求項３】
請求項１に記載のコンピュータプログラムプロダクトであって、表面及びその表面近傍の少なくとも一つのボクセルの要素との間の作用をモデル化する表面の作用操作をコンピュータが実行するようさらに構成されており、
（ａ）前記要素は、前記表面に関連する接線運動量を有し、
（ｂ）前記表面の作用操作は、前記要素における前記接線運動量の少なくとも一部を保持し、
（ｃ）保持された前記接線運動量の一部は、摩擦パラメータに対応し、
（ｄ）前記摩擦パラメータは、前記表面近傍の圧力の変化に基づいて変化することを特徴とするコンピュータプログラムプロダクト。
【請求項４】
請求項１に記載のコンピュータプログラムプロダクトであって、前記モデルは、特定の速度ベクトルを有する多数の粒子を表す各状態を用いて、異なる状態に対応するように特定のボクセルでの粒子を示すことを特徴とするコンピュータプログラムプロダクト。
【請求項５】
請求項４に記載のコンピュータプログラムプロダクトであって、
前記モデルは、重みと各状態を結びつけ、
前記コンピュータが作用操作を実行することは、前記コンピュータが、状態ベクトルに対応する状態と関連する重みを用いてその状態ベクトルに作用操作を実行することを含むコンピュータプログラムプロダクト。
【請求項６】
請求項４に記載のコンピュータプログラムプロダクトであって、前記モデルは、速度が０である多数の粒子を表す第１状態グループを含むことを特徴とするコンピュータプログラムプロダクト。
【請求項７】
請求項６に記載のコンピュータプログラムプロダクトであって、前記モデルは、単位速度で移動する粒子と関連する状態を具備する第２状態グループを含むことを特徴とするコンピュータプログラムプロダクト。
【請求項８】
請求項７に記載のコンピュータプログラムプロダクトであって、前記第２状態グループは、格子軸に対して、異なる方向に移動する粒子に対応する状態を含むことを特徴とするコンピュータプログラムプロダクト。
【請求項９】
請求項７に記載のコンピュータプログラムプロダクトであって、前記モデルは、前記単位速度の２倍で移動する粒子と関連する状態を具備する第３状態グループを含むことを特徴とするコンピュータプログラムプロダクト。
【請求項１０】
請求項９に記載のコンピュータプログラムプロダクトであって、前記モデルは、前記単位速度の３倍で移動する粒子と関連する状態を具備する第４状態グループを含むことを特徴とするコンピュータプログラムプロダクト。
【請求項１１】
請求項１０に記載のコンピュータプログラムプロダクトであって、
前記第２状態グループは、第１サブグループを含み、その第１サブグループは、単一の格子軸に沿った前記単位速度での粒子の動きに相当する状態を含んでおり、
前記第３状態グループは、第２サブグループ及び第３サブグループを含み、その第２サブグループは、単一の格子軸に沿った前記単位速度の２倍での粒子の動きに相当する状態を含んでおり、その第３サブグループは、２つの格子軸に沿った前記単位速度での粒子の動きに相当する状態を含んでおり、
前記第４状態グループは、第４サブグループを含み、その第４サブグループは、単一の格子軸に沿った前記単位速度の３倍での粒子の動きに相当する状態を含んでいることを特徴とするコンピュータプログラムプロダクト。
【請求項１２】
請求項１１に記載のコンピュータプログラムプロダクトであって、
前記第４状態グループは、第５サブグループを含み、その第５サブグループは、前記格子軸の一つに沿った前記単位速度、かつ、前記格子軸の一つ以外の軸の一つに沿った前記単位速度の２倍での粒子の動きに相当する状態を含んでおり、
前記モデルは、前記単位速度の４倍で移動する粒子と関連する状態を具備する第５状態グループを含み、
前記第５状態グループは、第６サブグループ及び第７サブグループを含み、その第６サブグループは、２つの格子軸に沿った前記単位速度の２倍での粒子の動きに相当する状態を含んでおり、その第７サブグループは、前記格子軸の一つに沿った前記単位速度、かつ、前記格子軸の一つ以外の軸の一つに沿った前記単位速度の３倍での粒子の動きに相当する状態を含んでいることを特徴とするコンピュータプログラムプロダクト。
【請求項１３】
請求項１２に記載のコンピュータプログラムプロダクトであって、
前記第２状態グループは、第１サブグループを含み、その第１サブグループは、単一の格子軸に沿った前記単位速度での粒子の動きに相当する状態を含んでおり、
前記第３状態グループは、第２サブグループを含み、その第２サブグループは、単一の格子軸に沿った前記単位速度の２倍での粒子の動きに相当する状態を含んでおり、
前記第４状態グループは、第３サブグループ及び第４サブグループを含み、その第３サブグループは、単一の格子軸に沿った前記単位速度の３倍での粒子の動きに相当する状態を含んでおり、その第４サブグループは、前記格子軸の全てに沿った前記単位速度での粒子の動きに相当する状態を含んでおり、
前記モデルは、前記単位速度の４倍で移動する粒子と関連する状態を具備する第５状態グループを含み、
前記第５状態グループは、第５サブグループを含み、その第５サブグループは、３つの格子軸の内の２つの軸に沿った前記単位速度の２倍での粒子の動きに相当する状態を含んでいることを特徴とするコンピュータプログラムプロダクト。
【請求項１４】
請求項１３に記載のコンピュータプログラムプロダクトであって、
前記第４状態グループは、第６サブグループを含み、その第６サブグループは、前記３つの格子軸の内の２つの軸に沿って前記単位速度の３倍で動き、かつ、前記３つの格子軸の２つ以外の軸には動いていない粒子の動きに相当する状態を含んでおり、
前記モデルは、前記単位速度の６倍で移動する粒子と関連する状態を具備する第６状態グループを含み、
前記第６状態グループは、第７サブグループを含み、その第７サブグループは、前記３つの格子軸の全てに沿った前記単位速度の３倍での粒子の動きに相当する状態を含んでおり、
前記モデルは、前記単位速度の５倍で移動する粒子と関連する状態を具備する第７状態グループを含み、
前記第７状態グループは、第８サブグループ及び第９サブグループを含み、その第８サブグループは、前記３つの格子軸の２つに沿った前記単位速度での粒子の動きに相当する状態を含んでおり、その第９サブグループは、前記３つの格子軸の２つ以外の軸に対して前記単位速度の３倍での粒子の動きに相当する状態を含んでいることを特徴とするコンピュータプログラムプロダクト。
【請求項１５】
物理プロセスをシミュレーションするための方法であって、コンピュータで実行されるその方法は、
ボクセルのための状態ベクトルであって、あるモデルに対応し、ボクセルで生じ得る運動量状態である特定の運動量状態に対応する複数のエントリーを具備する状態ベクトルをコンピュータがアクセス可能なメモリ内に記憶することと、
前記モデルに従って異なる運動量状態の要素間の作用をモデリングする作用操作を前記状態ベクトルに実行することと、
移動操作を前記状態ベクトルに実行し、前記モデルに従って要素の動きを新たなボクセルに反映させることと、
を含み、
前記モデルは、クヌーセン数が０．１より大きい高クヌーセン数流れをシミュレーションするよう適用されていることを特徴とする方法。
【請求項１６】
請求項１５に記載の方法であって、少なくとも一つの表面の表現を前記メモリ内に記憶することをさらに含むことを特徴とする方法。
【請求項１７】
請求項１５に記載の方法であって、表面及びその表面近傍の少なくとも一つのボクセルの要素との間の作用をモデル化する表面の作用操作を実行することをさらに含んでおり、
（ａ）前記要素は、前記表面に関連する接線運動量を有し、
（ｂ）前記表面の作用操作は、前記要素における前記接線運動量の少なくとも一部を保持し、
（ｃ）保持された前記接線運動量の一部は、摩擦パラメータに対応し、
（ｄ）前記摩擦パラメータは、前記表面近傍の圧力の変化に基づいて変化することを特徴とする方法。
【請求項１８】
請求項１５に記載の方法であって、前記モデルは、特定の速度ベクトルを有する多数の粒子を表す各状態を用いて、異なる状態に対応するように特定のボクセルでの粒子を示すことを特徴とする方法。
【請求項１９】
請求項１８に記載の方法であって、
前記モデルは、重みと各状態を結びつけ、
作用操作を実行することは、状態ベクトルに対応する状態と関連する重みを用いてその状態ベクトルに作用操作を実行することを含む方法。
【請求項２０】
請求項１８に記載の方法であって、前記モデルは、速度が０である多数の粒子を表す第１状態グループを含むことを特徴とする方法。
【請求項２１】
請求項２０に記載の方法であって、前記モデルは、単位速度で移動する粒子と関連する状態を具備する第２状態グループを含むことを特徴とする方法。
【請求項２２】
請求項２１に記載の方法であって、前記第２状態グループは、格子軸に対して、異なる方向に移動する粒子に対応する状態を含むことを特徴とする方法。
【請求項２３】
請求項２１に記載の方法であって、前記モデルは、前記単位速度の２倍で移動する粒子と関連する状態を具備する第３状態グループを含むことを特徴とする方法。
【請求項２４】
請求項２３に記載の方法であって、前記モデルは、前記単位速度の３倍で移動する粒子と関連する状態を具備する第４状態グループを含むことを特徴とする方法。
【請求項２５】
請求項２４に記載の方法であって、
前記第２状態グループは、第１サブグループを含み、その第１サブグループは、単一の格子軸に沿った前記単位速度での粒子の動きに相当する状態を含んでおり、
前記第３状態グループは、第２サブグループ及び第３サブグループを含み、その第２サブグループは、単一の格子軸に沿った前記単位速度の２倍での粒子の動きに相当する状態を含んでおり、その第３サブグループは、２つの格子軸に沿った前記単位速度での粒子の動きに相当する状態を含んでおり、
前記第４状態グループは、第４サブグループを含み、その第４サブグループは、単一の格子軸に沿った前記単位速度の３倍での粒子の動きに相当する状態を含んでいることを特徴とする方法。
【請求項２６】
請求項２５に記載の方法であって、
前記第４状態グループは、第５サブグループを含み、その第５サブグループは、前記格子軸の一つに沿った前記単位速度、かつ、前記格子軸の一つ以外の軸の一つに沿った前記単位速度の２倍での粒子の動きに相当する状態を含んでおり、
前記モデルは、前記単位速度の４倍で移動する粒子と関連する状態を具備する第５状態グループを含み、
前記第５状態グループは、第６サブグループ及び第７サブグループを含み、その第６サブグループは、２つの格子軸に沿った前記単位速度の２倍での粒子の動きに相当する状態を含んでおり、その第７サブグループは、前記格子軸の一つに沿った前記単位速度、かつ、前記格子軸の一つ以外の軸の一つに沿った前記単位速度の３倍での粒子の動きに相当する状態を含んでいることを特徴とする方法。
【請求項２７】
請求項２６に記載の方法であって、
前記第２状態グループは、第１サブグループを含み、その第１サブグループは、単一の格子軸に沿った前記単位速度での粒子の動きに相当する状態を含んでおり、
前記第３状態グループは、第２サブグループを含み、その第２サブグループは、単一の格子軸に沿った前記単位速度の２倍での粒子の動きに相当する状態を含んでおり、
前記第４状態グループは、第３サブグループ及び第４サブグループを含み、その第３サブグループは、単一の格子軸に沿った前記単位速度の３倍での粒子の動きに相当する状態を含んでおり、その第４サブグループは、前記格子軸の全てに沿った前記単位速度での粒子の動きに相当する状態を含んでおり、
前記モデルは、前記単位速度の４倍で移動する粒子と関連する状態を具備する第５状態グループを含み、
前記第５状態グループは、第５サブグループを含み、その第５サブグループは、３つの格子軸の内の２つの軸に沿った前記単位速度の２倍での粒子の動きに相当する状態を含んでいることを特徴とする方法。
【請求項２８】
請求項２７に記載の方法であって、
前記第４状態グループは、第６サブグループを含み、その第６サブグループは、前記３つの格子軸の内の２つの軸に沿って前記単位速度の３倍で動き、かつ、前記３つの格子軸の２つ以外の軸には動いていない粒子の動きに相当する状態を含んでおり、
前記モデルは、前記単位速度の６倍で移動する粒子と関連する状態を具備する第６状態グループを含み、
前記第６状態グループは、第７サブグループを含み、その第７サブグループは、前記３つの格子軸の全てに沿った前記単位速度の３倍での粒子の動きに相当する状態を含んでおり、
前記モデルは、前記単位速度の５倍で移動する粒子と関連する状態を具備する第７状態グループを含み、
前記第７状態グループは、第８サブグループ及び第９サブグループを含み、その第８サブグループは、前記３つの格子軸の２つに沿った前記単位速度での粒子の動きに相当する状態を含んでおり、その第９サブグループは、前記３つの格子軸の２つ以外の軸に対して前記単位速度の３倍での粒子の動きに相当する状態を含んでいることを特徴とする方法。
【請求項２９】
物理プロセスをシミュレーションするコンピュータシステムであって、そのシステムは、
ボクセルのための状態ベクトルであって、あるモデルに対応し、ボクセルで生じ得る運動量状態である特定の運動量状態に対応する複数のエントリーを具備する状態ベクトルをコンピュータがアクセス可能なメモリ内に記憶することと、
前記モデルに従って異なる運動量状態の要素間の作用をモデリングする作用操作を前記状態ベクトルに実行することと、
移動操作を前記状態ベクトルに実行し、前記モデルに従って要素の動きを新たなボクセルに反映させることと、
を実行するよう構成されており、
前記モデルは、クヌーセン数が０．１より大きい高クヌーセン数流れをシミュレーションするよう適用されていることを特徴とするシステム。
【請求項３０】
請求項２９に記載のシステムであって、少なくとも一つの表面の表現を前記メモリ内に記憶するようさらに構成されていることを特徴とするシステム。
【請求項３１】
請求項２９に記載のシステムであって、表面及びその表面近傍の少なくとも一つのボクセルの要素との間の作用をモデル化する表面の作用操作を実行するようさらに構成されており、
（ａ）前記要素は、前記表面に関連する接線運動量を有し、
（ｂ）前記表面の作用操作は、前記要素における前記接線運動量の少なくとも一部を保持し、
（ｃ）保持された前記接線運動量の一部は、摩擦パラメータに対応し、
（ｄ）前記摩擦パラメータは、前記表面近傍の圧力の変化に基づいて変化することを特徴とするシステム。
【請求項３２】
請求項２９に記載のシステムであって、前記モデルは、特定の速度ベクトルを有する多数の粒子を表す各状態を用いて、異なる状態に対応するように特定のボクセルでの粒子を示すことを特徴とするシステム。
【請求項３３】
請求項３２に記載のシステムであって、
前記モデルは、重みと各状態を結びつけ、
作用操作を実行することは、状態ベクトルに対応する状態と関連する重みを用いてその状態ベクトルに作用操作を実行することを含むシステム。
【請求項３４】
請求項３２に記載のシステムであって、前記モデルは、速度が０である多数の粒子を表す第１状態グループを含むことを特徴とするシステム。
【請求項３５】
請求項３４に記載のシステムであって、前記モデルは、単位速度で移動する粒子と関連する状態を具備する第２状態グループを含むことを特徴とするシステム。
【請求項３６】
請求項３５に記載のシステムであって、前記第２状態グループは、格子軸に対して、異なる方向に移動する粒子に対応する状態を含むことを特徴とするシステム。
【請求項３７】
請求項３６に記載のシステムであって、前記モデルは、前記単位速度の２倍で移動する粒子と関連する状態を具備する第３状態グループを含むことを特徴とするシステム。
【請求項３８】
請求項３７に記載のシステムであって、前記モデルは、前記単位速度の３倍で移動する粒子と関連する状態を具備する第４状態グループを含むことを特徴とするシステム。
【請求項３９】
請求項３８に記載のシステムであって、
前記第２状態グループは、第１サブグループを含み、その第１サブグループは、単一の格子軸に沿った前記単位速度での粒子の動きに相当する状態を含んでおり、
前記第３状態グループは、第２サブグループ及び第３サブグループを含み、その第２サブグループは、単一の格子軸に沿った前記単位速度の２倍での粒子の動きに相当する状態を含んでおり、その第３サブグループは、２つの格子軸に沿った前記単位速度での粒子の動きに相当する状態を含んでおり、
前記第４状態グループは、第４サブグループを含み、その第４サブグループは、単一の格子軸に沿った前記単位速度の３倍での粒子の動きに相当する状態を含んでいることを特徴とするシステム。
【請求項４０】
請求項３９に記載のシステムであって、
前記第４状態グループは、第５サブグループを含み、その第５サブグループは、前記格子軸の一つに沿った前記単位速度、かつ、前記格子軸の一つ以外の軸の一つに沿った前記単位速度の２倍での粒子の動きに相当する状態を含んでおり、
前記モデルは、前記単位速度の４倍で移動する粒子と関連する状態を具備する第５状態グループを含み、
前記第５状態グループは、第６サブグループ及び第７サブグループを含み、その第６サブグループは、２つの格子軸に沿った前記単位速度の２倍での粒子の動きに相当する状態を含んでおり、その第７サブグループは、前記格子軸の一つに沿った前記単位速度、かつ、前記格子軸の一つ以外の軸の一つに沿った前記単位速度の３倍での粒子の動きに相当する状態を含んでいることを特徴とするシステム。
【発明の詳細な説明】
【技術分野】
【０００１】
本開示は、流体流量のような物理プロセスのコンピュータシミュレーションに関する。
【背景技術】
【０００２】
高レイノルズ数流れは、ナビエ−ストーク微分方程式の離散化解を生成することによってシミュレーションが行われる。その離散化は、巨視的物理量（例えば、密度、温度、流速）を表す変化量に基づく多くの空間的な離散位置のそれぞれで高精度浮動小数点演算操作を実行することによって行われる。他のアプローチは、微分方程式を、格子ガス（又はセル）オートマトンとして一般的に知られているものに置き換える。その場合、ナビエ−ストーク方程式を解くことによって得られる巨視的レベルのシミュレーションは、巨視的レベルのモデルに置き換えられる。その巨視的レベルのモデルは、格子上の位置を移動する粒子に操作を実行する。
【０００３】
従来の格子ガスシミュレーションは、ショートベクトルのビットによって表される粒子を用いて各格子位置での限られた数の粒子を仮定する。各ビットは、特定の方向に移動する粒子を表す。例えば、ベクトルにおける１ビットは、特定の方向に沿って移動する粒子の存在（１が設定される）又は非存在（０が設定される）を示してもよい。そのようなベクトルは、６ビットとしてもよく、例としては、値１１００００は、Ｘ軸の反対方向に沿って移動する２つの粒子、そして、粒子はＹ軸及びＸ軸に沿って移動していないことを示す。衝突ルールのセットは、各位置での粒子間の衝突の挙動を制御する（例えば、１１００００ベクトルは、００１１００ベクトルとなってもよく、これは、Ｘ軸に沿って移動する２つの粒子間の衝突がＹ軸に沿って離れる２つの粒子をもたらすことを示す）。ルールは、状態ベクトルを参照テーブルに適用することによって実装される。これは、ビット上で置換を行うことによって行われる（例えば、１１００００から００１１００へ変換する）。すると、粒子は、隣接している位置に移動する（例えば、Ｙ軸に沿って移動する２つの粒子は、Ｙ軸に沿って隣接している左及び右の位置に移動しようとする）。
【０００４】
改善された格子ガス技術において、各格子位置での状態ベクトルは、粒子のエネルギー及び移動方向の変動を提供するよう、より多くのビットを含み（例えば、亜音速流れについては５４ビット）、そして、完全な状態ベクトルのサブセットを含む衝突ルールが用いられる。さらに改善されたシステムにおいては、単一粒子は、各格子位置又は各ボクセル（本開示ではこれら２つの用語が交互に用いられる）での各運動量の状態で存在することが可能であることにとどまらない。例えば、８ビットの実装において、０乃至２５５の粒子が、特定のボクセルで特定の方向に移動することができる。ビットのセットの代わりとなる状態ベクトルは、整数のセットである（例えば、８ビットバイトのセットは、０から２５５の範囲における整数を有する）。それぞれの整数は、所定の状態における粒子の数を表す。
【０００５】
最近になって、格子ボルツマン法（ＬＢＭ）は、流体のメゾスコピックな表現を使用して、複雑な形態における不安定な３Ｄ圧縮性乱流を従来の計算流体力学（ＣＦＤ）アプローチで可能なものより深いレベルでシミュレーションを行う。ＬＢＭ法の概要を以下に示す。
【０００６】
ボルツマンレベルのメゾスコピック表現
流体系が、いわゆる「メゾスコピック」レベルで運動論的方程式によって表すことができることは、統計物理学において周知である。このレベルでは、個々の粒子の詳細な運動を決定する必要はない。代わりに、流体の特性が、単一粒子の位相空間を用いて定義される粒子の分布関数ｆ＝ｆ（ｘ，ｖ，ｔ）によって表される。ここで、ｘは、空間座標であり、ｖは、粒子の速度についての座標系である。質量密度や、流速、温度などのような典型的な流体力学的物理量は、粒子の分布関数についての単純なモーメントである。粒子の分布関数についての力学は、以下のボルツマン方程式を遵守する。
2010500654.tif 000003
方程式（１）
【０００７】
ここで、Ｆ（ｘ，ｔ）は、（ｘ，ｔ）における外部の体積力又は自己無撞着的に生成される体積力を表す。衝突項Ｃは、様々な速度及び位置における粒子の作用を表す。衝突項Ｃのために特定の形態を明確にすることなく、上記ボルツマン方程式が、周知である（本来、ボルツマン方程式で導かれるような）希薄気体の状況だけでなく全ての流体系に適用可能であることは重要である。
【０００８】
一般的には、Ｃは、二点相関関数の複雑な複数次元の積分を含む。分布関数ｆを用いた閉鎖系を形成する目的だけではなく、効率的な計算の目的に関しても、最も便利かつ物理的に安定した形態の一つは、ＢＧＫオペレータである。ＢＧＫオペレータは、物理的根拠にしたがって構成され、どんなに衝突が小さいものであっても、分布関数は、特定の局所平衡に近づける。その局所平衡は、｛ｆeq（ｘ，ｖ，ｔ）｝によって以下の衝突を介して与えられる。
2010500654.tif 000004
方程式（２）
【０００９】
ここで、パラメータτは、衝突を介した平衡に対する特有の緩和時間を表す。粒子（例えば、原子や分子など）を処理する際、緩和時間は、大体一定に保たれる。「合成物」の（水による運動）表現において、この緩和時間は、ひずみ速度や乱流運動エネルギーなどのような流体力学的変数の関数である。したがって、乱流は、乱流粒子（「渦巻き」）のガスとして局所的に決まっている固有性を用いて表すことが可能である。
【００１０】
ボルツマン−ＢＧＫ方程式の数値的解法は、計算上、ナビエ−ストークス方程式の解法を越えるいくつかの有利な点を有する。第１に、その方程式において、複雑な非線形項又は高次空間導関数が無いことが即座に分かるであろう。したがって、移流の不安定性に関する問題はほとんど存在しない。記述のこのレベルで、圧力を扱う必要が無いことから、方程式は、局所的であり、アルゴリズムの並列化のための重要な有利性を提供する。線形移流オペレータの他の望ましい特徴は、二次空間導関数を用いた拡散性のオペレータが無いという事実とともに、どのように粒子が固体表面と正確に作用するかを近似させるかという方法において、流体の偏微分方程式（ｐａｒｔｉａｌ  ｄｉｆｆｅｒｅｎｔｉａｌ  ｅｑｕａｔｉｏｎｓ、「ＰＤＥ」）のための数学的条件よりむしろ、非滑り面又は滑り面などのような、物理的境界条件の理解のし易さである。直接的な利益の一つは、固体表面上の接触についての動きの問題を扱わないことである。その利益は、シミュレーションソフトウエアに基づく格子ボルツマンが、複雑な乱流空気力学のシミュレーションを実行可能とすることに役立つ。加えて、有限表面粗さなどのような境界による所定の物理的特徴は、力に組み込むことが可能である。さらに、ＢＧＫ衝突オペレータは、非常に局所的である。さらに、矛盾のない体積力の計算は、近傍の情報のみを通じて達成することが可能である。したがって、ボルツマン−ＢＧＫ方程式の計算は、並列処理のために効果的に適用され得る。
【００１１】
格子ボルツマンの公式化
ボルツマン連続体方程式を解決することは、位置及び速度の位相空間における微分積分方程式の数値的計算を必要とする点で重要な課題を示す。位置だけでなく速度の位相空間が離散化され得ることが順守されるとき顕著な単純化が行われる。その単純化は、ボルツマン方程式の解法のために効率的な数値アルゴリズムをもたらす。水力学的物理量は、単純合計の項において説明され得る。その項は、近傍の情報に依存する。格子ボルツマン方程式（ｌａｔｔｉｃｅ  Ｂｏｌｔｚｍａｎｎ  ｅｑｕａｔｉｏｎ）の公式化が、速度ｖ（∈｛ｃi，ｉ. . .，ｂ｝）の個別のセットで粒子の発生を規定する格子ガスモデルに基づいていたとしても、この方程式は、ボルツマン連続体方程式の離散化として第１原理から体系的に導かれ得る。結果として、ＬＢＥは、格子ガスアプローチに関連する周知の問題に悩まされない。したがって、位相空間における連続体の分布関数ｆ（ｘ，ｖ，ｔ）を扱う代わりに、離散速度インデックスをラベルされた添え字（下付文字）を有する離散型分布ｆi（ｘ，ｔ）の有限集合を探知すれば足りる。巨視的記述の代わりに、この運動論的方程式を扱うことの重要な利点は、問題の場所によって系の増加した位相空間が補正される。対象性の考慮に起因して、速度のセットは、それが配置空間にわたっているとき所定の格子構造を形成するこのような方法で選択される。そのような離散系の力学は、形式ｆi（ｘ＋ｃi，ｔ＋１）−ｆi（ｘ，ｔ）＝Ｃi（ｘ，ｔ）を有するＬＢＥに従う。ここで、衝突オペレータは、上述のようなＢＧＫ形式をとる。平衡分布形式の適切な選択によって、格子ボルツマン方程式が正確な流体力学及び熱水力学をもたらすことは、論理的に示され得る。すなわち、ｆi（ｘ，ｔ）から導かれる水力学モーメントは、巨視的極限においてナビエ−ストークス方程式に従う。これらモーメントは以下のように示される。
2010500654.tif 000005
方程式（３）
ここで、ρ、ｕ、Ｔは、それぞれ、流体密度、速度、温度であり、そして、Ｄは、離散速度空間における次元である（物理空間の次元とは全く等しくない）。
【図面の簡単な説明】
【００１２】
【図１】ＬＢＭモデルの速度成分を示す。
【図２】ＬＢＭモデルの速度成分を示す。
【図３】物理プロセスシミュレーションシステムによる手順のフローチャートである。
【図４】マイクロブロックの斜視図である。
【図５Ａ】図３のシステムが用いる格子構造の説明図である。
【図５Ｂ】図３のシステムが用いる格子構造の説明図である。
【図６】可変分解能技術の説明。
【図７】可変分解能技術の説明。
【図８】表面の小面によって影響を受ける領域の説明。
【図９】ボクセルから表面への粒子の動きの説明。
【図１０】表面から表面への粒子の動きの説明。
【図１１】表面動力学を実行する手順のフローチャートである。
【図１２】異なるサイズのボクセル間における作用の説明
【図１３】可変分解能条件下で小面との作用をシミュレーションする手順のフローチャート。
【発明を実施するための形態】
【００１３】
Ａ．シミュレーション空間のモデル化
物理プロセスシミュレーションシステムは、改良された格子ボルツマン法（ＬＢＭ）モデルを用いる。そのモデルは、システムが正確かつ効果的に高クヌーセン数流れをシミュレーションすることを可能にする。高クヌーセン数流れは、分子平均自由行程（すなわち、粒子が他の粒子と衝突することなく移動するよう期待される距離）と特性長さ基準（characteristic length scale）（すなわち、ボクセルサイズ）の比率が高いところの流れである。高クヌーセン数流れは、非常に低い密度のガスで発生し得る。その非常に低い密度のガスは、例として、高地や、マイクロデバイスなどで見られるような非常に小さい空間規模を有する表面の周りの流体中で見られる。
【００１４】
ＬＢＭベースのシミュレーションにおいて、流体流量は、離散速度ｃiのセットで評価される分布関数値ｆiによって表すことが可能である。分布関数の力学は、方程式１によって決定する。ここで、ｆi（０）は、方程式（４）で定義されている平衡分布関数として既知である。
2010500654.tif 000006
方程式（４）
【００１５】
上記の一次方程式は、周知の格子ボルツマン方程式である。その方程式は、分布方程式ｆiの時間発展を表す。左辺は、いわゆる「流動プロセス」に起因する分布の変化を表す。流動プロセスは、流体のポケットが、グリッド位置で始まって、速度ベクトルの一つに沿って次のグリッド位置へ移動するときのことを示す。その時、我々は、「衝突要因」すなわち、流体のポケットが始まることによる近くの流体のポケットへの影響を算出する。流体は、他のグリッド位置への移動のみが可能である。これにより、全体の速度の全要素が共通のスピードの倍数となるよう速度ベクトルの適切な選択が必要となる。
【００１６】
一次方程式の右辺は、前述の「衝突オペレータ」である。その衝突オペレータは、流体のポケット間の衝突に起因する分布関数の変化を表す。ここで用いられる衝突オペレータの特定の形態は、Ｂｈａｔｎａｇａｒ−Ｇｒｏｓｓ−Ｋｒｏｏｋ（ＢＧＫ）に起因する。それは、分布関数を二次方程式によって付与される規定値に近づける。その場合、平衡形態となる。
【００１７】
二次方程式は、流体系が平衡状態に到達しているとき、分布関数が何になる傾向があるかを示す。格子ボルツマン方程式における平衡状態の形態は、巨視的な流体力学的挙動を判断するのに重要な役割を果す。
【００１８】
このシミュレーションから、質量ρ及び流体速度ｕなどのような従来の流体変数は、方程式（３）において単純な合計として得られる。ここで、ｃi及びｗiについての集約値はＬＢＭモデルを定義する。
【００１９】
ＬＢＭモデルは、拡張可能なコンピュータプラットフォーム上で効果的に実装することができ、一定ではない時間の流れ及び複雑な境界条件のために強度な安定性を伴って実行することが可能である。
【００２０】
ボルツマン方程式から流体系のための動作についての巨視的方程式を得る標準的な技術は、チャップマン・エンスコッグ（Chapman-Enskog）法である。その場合、クヌーセン数が小さいと仮定すると、完全なボルツマン方程式の逐次近似法が必要とされる。得られる巨視的方程式は、正確さの順にそれぞれ、オイラー方程式、ナビエ−ストークス方程式、バーネット方程式、スーパーバーネット方程式である。ナビエ−ストークスの近似のレベルは、多くの流体工学の問題にとって十分なものであると考えられる一方、それは、クヌーセン数が小さくない（例えば、Ｋｎ＞０．１）場合には十分とはならない。このような状況では、バーネト及びスーパーバーネットの近似のレベルが、しばしば用いられるが、これら近似に関する理論上及び実施上の困難性の両方が存在する。
【００２１】
従前のＬＢＭ系では、その系が採用したモデルによって導入されている離散的数値結果が、高クヌーセン数において重要となることから、高クヌーセン数流れを正確にシミュレーションすることは困難であった。以下で説明するようなモデルを採用することによって、これらの問題は、回避が可能となり、高クヌーセン数での物理流量について正確な数値で再生されたものを導くＬＢＭシミュレーションが、達成され得る。
【００２２】
以下で説明するような表１及び表２は、４つの代表的なモデル（３Ｄ−１、２Ｄ−１、３Ｄ−２、２Ｄ−２、）を説明する。それらモデルは、高クヌーセン数流れをシミュレーションするのに用いることが可能である。
【００２３】
表１は、バーネット順（Burnett order）での高クヌーセン数動力学的モデルの離散速度及び離散重さである。対称なグループにおけるポイントの数、速度ｃi、重さｗiを一覧する。添え字ＦＳは、完全な対称性を有する順列を意味する。全ての速度は、音速に対して正規化されている。単位速度あたりの大きさは、
2010500654.tif 000007
である。
2010500654.tif 000008
【００２４】
表２は、バーネット順での高クヌーセン数動力学的モデルの個別の離散速度及び離散重さである。対称なグループにおけるポイントの数、速度ｃi、重さｗiを一覧する。添え字ＦＳは、完全な対称性を有する順列を意味する。全ての速度は、音速に対して正規化されている。単位速度あたりの大きさは、ｒ＝１．１９６９７９７７０３９３０７４３５９０である。
2010500654.tif 000009
【００２５】
流体系において、密度の微小擾乱は、音速で移動する。ガス系において、音速は、一般に、温度によって決定する。流体内の圧縮率の効果についての重要性は、特性速度と音速との比率によって測定される。これは、マッハ数として既知である。上述の表における所定のモデルでは、全ての速度は、音速の単位で与えられる。これら速度の方向及び大きさの選択は、巨視的な挙動の正確な複製を高クヌーセン数での格子ボルツマン方程式によって保証する。添え字ＦＳは、完全は対称性を伴った置換が行われることを意味する。例えば、２次元のベクトルにおいて、（ｒ，０）＿ＦＳは、（ｒ，０）、（−ｒ，０）、（０，ｒ）、（０，−ｒ）によって示される４つの速度を意味する。
【００２６】
また、図１で示すように、第１モデル（２Ｄ−１）１００は、２１の速度を有する２次元モデルである。これら２１の速度の一つ（１０５）は、移動しない粒子を表す。４つの速度についての３つのセットは、それぞれ、正規化されている速度（ｒ）（１１０乃至１１３）で、正規化されている速度の２倍（２ｒ）（１２０乃至１２３）で、正規化されている速度の３倍（３ｒ）（１３０乃至１３３）で、格子におけるｘ軸又はｙ軸のどちらかに沿って正方向又は負方向のいずれかの方向に移動する粒子を表す。４つの速度についての２つのセットは、正規化されている速度（ｒ）（１４０乃至１４３）又は正規化されている速度の２倍（２ｒ）で、格子のｘ軸又はｙ軸の両方に対して移動する粒子を表す。
【００２７】
また、図２に示すように、第２モデル（３Ｄ−１）２００は、３９の速度を含んでいる３次元モデルである。この場合、各速度は、図２の矢印の一つによって表される。これら３９の速度の一つは、移動しない粒子を表す。６つの速度についての３つのセットは、それぞれ、正規化されている速度（ｒ）、正規化されている速度の２倍（２ｒ）、正規化されている速度の３倍（３ｒ）のいずれかの速度で、格子におけるｘ軸、ｙ軸、又はｚ軸に沿って正方向又は負方向のどちらかの方向に移動する粒子を表す。８つの速度は、正規化されている速度（ｒ）で、格子のｘ軸、ｙ軸、又はｚ軸の３つ全てに対して移動する粒子を表す。１２の速度は、正規化されている速度の２倍（２ｒ）で、格子のｘ軸、ｙ軸、又はｚ軸の２つに対して移動する粒子を表す。
【００２８】
第３モデル（３Ｄ−２）及び第４モデル（２Ｄ−２）は、より複雑なモデルであり、３Ｄ−２は、１０１の速度を含み、２Ｄ−２は、３７の速度を含む。それら速度は、それぞれ表１及び表２内に記載されているように各軸に沿っている要素によって、より明確に表現される。
【００２９】
３次元モデル３Ｄ−２に関して、１０１の速度の一つは、移動しない粒子を表す（グループ１）。６つの速度についての３つのセットは、それぞれ、正規化されている速度（ｒ）、正規化されている速度の２倍（２ｒ）、正規化されている速度の３倍（３ｒ）のいずれかで、格子におけるｘ軸、ｙ軸、又はｚ軸に沿って正方向又は負方向のどちらかの方向に移動する粒子を表す（グループ２、グループ４、グループ７）。８つの速度についての３つのセットは、それぞれ、正規化されている速度（ｒ）で、正規化されている速度の２倍（２ｒ）で、正規化されている速度の３倍（３ｒ）で、格子のｘ軸、ｙ軸、ｚ軸の３つ全てに対して移動する粒子を表す（グループ３、グループ８、グループ１０）。１２の速度は、正規化されている速度の２倍（２ｒ）で、格子のｘ軸、ｙ軸、又はｚ軸の２つに対して移動する粒子を表す（グループ６）。２４の速度は、格子のｘ軸、ｙ軸、又はｚ軸の２つに対して正規化されている速度（ｒ）及び正規化されている速度の２倍（２ｒ）で移動し、そして、残りの軸に対しては移動しない粒子を表す（グループ５）。２４の速度は、格子のｘ軸、ｙ軸、又はｚ軸の２つに対して正規化されている速度（ｒ）で移動し、そして、残りの軸に対して正規化されている速度の３倍（３ｒ）で移動する粒子を表す（グループ９）。
【００３０】
２次元モデル２Ｄ−２に関して、３７の速度の一つは、移動しない粒子を表す（グループ１）。４つの速度についての３つのセットは、それぞれ、正規化されている速度（ｒ）で、正規化されている速度の２倍（２ｒ）で、正規化されている速度の３倍（３ｒ）で、格子におけるｘ軸又はｙ軸のどちらかに沿って正方向又は負方向のいずれかの方向に移動する粒子を表す（グループ２、グループ４、グループ７）。４つの速度についての２つのセットは、正規化されている速度（ｒ）又は正規化されている速度の２倍（２ｒ）で、格子のｘ軸又はｙ軸の両方に対して移動する粒子を表す。８つの速度は、格子のｘ軸及びｙ軸の一つに対して正規化されている速度（ｒ）で移動し、そして、他の軸に対して正規化されている速度の２倍（２ｒ）で移動する粒子を表す。８つの速度は、格子のｘ軸及びｙ軸の一つに対して正規化されている速度（ｒ）で移動し、そして、他の軸に対して正規化されている速度の３倍（３ｒ）で移動する粒子を表す。
【００３１】
上述のＬＢＭモデルは、２次元及び３次元両方における高クヌーセン数流れの数値シミュレーションのために、効果的かつ強固な離散速度動力学的モデルの特定のクラスを提供する。この種類のモデルは、それら速度に関連する離散した速度及び重さの特定のセットを含む。速度は、速度空間におけるデカルト座標の格子ポイントと一致する。その速度空間は、正確かつ効率的な離散速度モデルの実装を促進する。とりわけ、格子ボルツマンモデルとして既知である。そのようなモデルを用いることで、高クヌーセン数を伴う流体を、非常に正確にシミュレーションすることができる。
【００３２】
図３を参照すると、物理プロセスシミュレーションシステムが、手順３００によって動作し、流体流量などのような物理プロセスをシミュレーションする。シミュレーションに先立ち、一群のボクセルのようにシミュレーション空間がモデル化される（ステップ３０２）。通常、シミュレーション空間は、コンピュータを使ったデザイン（ＣＡＤ）プログラムを用いて生成される。例として、ＣＡＤプログラムは、風洞に配置されているマイクロデバイスを描画するよう用いられ得る。したがって、ＣＡＤプログラムによって生成されるデータは、適切な分解能を有する格子構造を加えるよう、さらに、シミュレーション空間内でオブジェクト及び表面を明らかにするよう処理が実行される。
【００３３】
格子の分解能は、シミュレーションを行うシステムのレイノルズ数に基づいて選択され得る。レイノルズ数は、流量の粘性（ｖ）、流量内のオブジェクトの特性長さ（Ｌ）、そして、流量の特性速度（ｕ）に関連し、以下で表される。
Ｒｅ＝ｕＬ／ｖ
方程式（５）
【００３４】
オブジェクトの特性長さは、オブジェクトについて長さの基準である特徴を表す。例えば、マイクロデバイスの周りの流量がシミュレーションされる場合、マイクロデバイスの高さは、特性長さとなるよう考慮され得る。オブジェクトの小さい領域の周り（例えば、自動車のサイドミラーなど）の流量が対象であるとき、シミュレーションの分解能を上昇させることができ、また、上昇させた分解能の場所を対象の領域の周辺に用いることができる。格子の分解能が上昇すると、ボクセルの大きさは減少する。
【００３５】
状態空間は、ｆi（ｘ，ｔ）として表される。ここで、ｆiは、時間ｔにおける３次元ベクトルｘによって示される格子位置での状態ｉの単位体積あたりの要素又は粒子の数である（すなわち、状態ｉにおける粒子の密度）。既知の時間増分に関しては、粒子の数は、単にｆi（ｘ）として示される。格子位置の全ての状態の組み合わせは、ｆ（ｘ）として示される。
【００３６】
状態の数は、各エネルギーレベルの範囲で起こり得る速度ベクトルの数によって決定する。速度ベクトルは、３次元（ｘ、ｙ、ｚ）空間内の整数で表現される線速度からなる。多くの種類のシミュレーションのために、状態の数は増加する。
【００３７】
各状態ｉは、特定のエネルギーレベル（すなわち、エネルギーレベルが０や１又は２など）で異なる速度ベクトルを表す。各状態の速度ｃiは、以下のように３次元それぞれの速度で表される。
ｃi＝（ｃi,x，i,y，ｃi,z）
方程式（６）
【００３８】
エネルギーレベルが０である状態は、どの次元にも移動していない停止している粒子を表す。すなわち、ｃstopped＝（０，０，０）である。エネルギーレベルが１の状態は、３つの次元のうちの一つの速度が±１、他の２つの次元の速度が０である粒子を表す。エネルギーレベルが２である状態は、３つの次元全てにおける速度が±１であるか、又は３つの次元のうちの一つ速度が±２かつ他の２つの次元の速度が０であるかのどちらかである粒子を表す。
【００３９】
３つのエネルギーレベルの起こり得る順列の全てを生成すると、考え得る状態の合計は３９となる（エネルギーが０の状態が一つ、エネルギーが１の状態が六つ、エネルギーが３の状態が八つ、エネルギーが４の状態が六つ、エネルギーが８の状態が１２、エネルギーが９の状態が六つ）。
【００４０】
各ボクセル（すなわち、各格子位置）は、状態ベクトルｆ（ｘ）によって表される。状態ベクトルは、ボクセルのステータスを十分に定義し、３９のエントリーを含む。３９のエントリーは、エネルギーが０の状態である一つ、エネルギーが１の状態である六つ、エネルギーが３の状態である八つ、エネルギーが４の状態である六つ、エネルギーが８の状態である１２、エネルギーが９の状態が六つからなる。この速度のセットを用いることによって、システムは、達成される平衡状態ベクトルのためにマクスウェルボルツマン統計を生成することが可能となる。
【００４１】
処理効率に関して、ボクセルは、マイクロブロックと呼ばれる２×２×２の立体に分類される。マイクロブロックは、ボクセルの並列処理を可能にするよう、さらに、データ構造に関連するオーバーヘッドを最小にするよう構成されている。マイクロブロックにおけるボクセルのための省略表記は、Ｎi（ｎ）として定義される。ここで、ｎは、マイクロブロック内の格子位置についての関連する位置を表し、ｎ∈｛０、１、２・・・７｝となる。マイクロブロックを図４に示す。
【００４２】
図５Ａ及び図５Ｂを参照すると、表面Ｓ（図５Ａ）は、小面Ｆαの一群としてシミュレーション空間内に示される（図５Ｂ）。
Ｓ＝｛Ｆα｝
方程式（７）
ここで、αは、特定の小面を示すインデックスである。小面は、ボクセル境界に制限されないが、通常、小面に隣接しているボクセルのサイズとほぼ同じ又はそれよりわすかに小さい大きさである。これにより、小面は、比較的少数のボクセルに影響を及ぼす。小面は、表面動力学を実行する目的のための特性を有している。具体的には、各小面Ｆαは、単位法線（ｎα）、表面積（Ａα）、中心位置（ｘα）、小面分布関数（ｆi（α））を有する。その小面分布関数は、小面の表面動力学的特性を表す。
【００４３】
図６を参照すると、分解能の異なるレベルをシミュレーション空間の異なる領域で用いることが可能であり、処理効率を向上させる。通常、オブジェクト６５５の周りの領域６５０は、最も重要である。従って、その領域６５０は、最も高い分解能を用いてシミュレーションが行われる。粘性の効果は、オブジェクトからの距離とともに減少することから、減少した分解能のレベルは、領域６６０、６６５（すなわち、ボクセルの体積は大きくなる）をシミュレーションするために用いられる。その領域６６０、６６５は、オブジェクト６５５から距離をあけて位置している。同様に、図７に示す通り、分解能の低いレベルは、オブジェクト７７５ほど顕著な特徴を有さない領域７７０をシミュレートするために用いることができる。一方、最高レベルの分解能は、オブジェクト７７５で最も顕著な特徴を有する領域７８０をシミュレーションするために用いる。中心からはなれた領域７８５は、最も低いレベルの分解能及び最も大きいボクセルを用いてシミュレーションされる。
【００４４】
Ｂ．小面によって影響を受けるボクセルの識別
図３を参照すると、一度、シミュレーション空間がモデル化されると（ステップ３０２）、一又は二以上の小面によって影響を受けるボクセルが識別される（ステップ３０４）。ボクセルは、さまざまな形で小面から影響を受ける可能性がある。まず、一又は二以上の小面が交差しているボクセルは、交差していないボクセルと比較して体積が小さいことから影響を受けやすい。これは、小面及びその小面によって表現される表面の下にある物質が、ボクセルの部分を占有することから生じる。部分要因Ｐf（ｘ）は、小面からは影響を受けないボクセルの部分を示す（すなわち、流体又は他の物質によって占有される部分に関して、流量がシミュレーションされる）。交差していないボクセル関しては、Ｐf（ｘ）は１である。
【００４５】
粒子を小面に移動させ、さらに、粒子を小面から受けとることによって一又は二以上の小面と作用するボクセルは、小面から影響を受けるボクセルとして識別される。小面が交差している全てのボクセルは、小面から粒子を受けとる少なくとも一つの状態、そして、粒子を小面に移動させる少なくとも一つの状態を含む。ほとんどの場合、追加のボクセルは、そのような状態を含む。
【００４６】
図８を参照すると、速度ベクトルｃiが０ではない各状態iに関して、小面Ｆαは、平行六面体Ｇiαによって定義されている領域から粒子を受けとり、その領域へ粒子を移動させる。その平行六面体Ｇiαは、速度ベクトルｃiと小面の単位法線ｎαとのベクトルドット積の大きさによって定義される高さ（｜ｃiｎi｜）、そして、小面の表面積Ａαによって定義される底部を有する。これにより、平行六面体Ｇiαの体積Ｖiαは、以下の方程式で表される。
Ｖiα＝｜ｃiｎα｜Ａα
方程式（８）
【００４７】
状態の速度ベクトルが小面に向けられているとき（｜ｃiｎi｜＜０）、小面Ｆαは、体積Ｖiαから粒子を受けとり、状態の速度ベクトルが小面から遠ざかる方に向けられているとき（｜ｃiｎi｜＞０）、小面Ｆαは、粒子を領域へ移動させる。内側の角のような出っ張っていない特徴の周辺で見られる状況である、他の小面が平行六面体Ｇiαの一部を占有するとき、以下で説明するように、この表現は、変更しなければならない。
【００４８】
小面Ｆαの平行六面体Ｇiaは、複数のボクセルの一部または全てと重なる可能性がある。ボクセル又はそれの一部の数は、ボクセルのサイズに対する小面のサイズ、状態のエネルギー、及び、格子構造に対する小面の方向に依存する。影響を受けているボクセルの数は、小面のサイズとともに増加する。したがって、上述のように、小面のサイズは通常、小面近くに配置されているボクセルのサイズとほぼ同じ又は小さくなるよう選択される。
【００４９】
平衡六面体Ｇiαによって重なるボクセルの一部Ｎ（ｘ）は、Ｖiα（ｘ）として定義される。この項を用いることで、ボクセルＮ（ｘ）と小面Ｆαとの間を移動する状態ｉの粒子における流束Γiα（ｘ）は、ボクセルと重なる領域の体積（Ｖiα（ｘ））をボクセルにおける状態ｉの粒子の密度（Ｎi（ｘ））に乗じることで等しくなる。
Γiα（ｘ）＝Ｎi（ｘ）Ｖiα（ｘ）
方程式（９）
【００５０】
平行六面体Ｇが一又は二以上の小面によって交差しているときは、以下の条件が正しくなる。
Ｖiα＝ΣＶα（ｘ）（ｘ）  ＋  ΣＶiα（β）（ｆi）
方程式（１０）
ここで、第１の和は、Ｇiαによって重なっている全てのボクセルを定義し、第２項は、Ｇiαと交差する全ての小面を定義する。他の小面と平行六面体が交差しないとき、この式は、以下のように変換される。
Ｖiα＝ΣＶiα（ｘ）
方程式（１１）
【００５１】
Ｃ．シミュレーションの実行
一又は二以上の小面によって影響を受けているボクセルが、一度識別されると（ステップ３０４）、タイマーが初期化され、シミュレーションを開始する（ステップ３０６）。シミュレーションの各時間増分の間、ボクセルからボクセルへの粒子の動きが、移流段階でシミュレーションされる（ステップ３０８乃至３１６）。その移流段階は、粒子と表面の小面との作用を明らかにする。次に、衝突段階は、各ボクセル内の粒子の作用をシミュレーションする（ステップ３１８）。その後、タイマーには、一定時間量が追加される。一定時間量が追加されたタイマーがシミュレーション完了を示していない場合（ステップ３２２）、移流段階及び衝突段階（ステップ３０８乃至３２０）は繰り返される。一定時間量が追加されたタイマーがシミュレーション完了を示している場合（ステップ３２２）、シミュレーション結果が、記憶及び／又は表示される（ステップ３２４）。
【００５２】
１．表面についての境界条件
表面との作用を正確にシミュレーションするために、各小面は、４つの境界条件を満たさなければならない。第１に、小面が受けとる粒子の質量すべてを合わせた質量は、その小面が移動させる粒子の質量すべてを合わせたものと等しくしなければならない（すなわち、小面への正味質量の流束は、０にしなければならない）。第２に、小面が受けとる粒子のエネルギーすべてを合わせたエネルギーは、その小面が移動させる粒子のエネルギーすべてを合わせたものと等しくしなければならない（すなわち、小面への正味エネルギーの流束は、０にしなければならない）。これら２つの条件は、各エネルギーレベル（すなわち、エネルギーレベルが１及び２）での正味質量の流速が０になるよう要求することによって満たすことが可能となる。
【００５３】
他の２つの境界条件は、小面と作用する粒子の正味運動量に関連する。ここでは滑り面と呼ばれるような表面摩擦が無い表面に関して、接線での正味運動量流束は、０としなければならず、法線での正味運動量流束は、小面での局部圧力に等しくしなければならない。したがって、小面の法線ｎαと垂直である、受けとられる全てを合わせた運動量と移動した全てを合わせた運動量の要素（すなわち、接線要素）は、等しくなければならず、小面の法線ｎαと平行である、受けとられる全てを合わせた運動量と移動した全てを合わせた運動量の要素間の相違（すなわち、法線要素）は、小面での居部圧力に等しくなければならない。非滑り面の表面については、表面の摩擦は、小面が受けとる粒子における接線運動量全てを合わせた運動量に応じて、摩擦量に関連する要因により、小面が移動させた粒子における接線運動量全てを合わせた運動量を減少させる。
【００５４】
２．ボクセルから小面への収集
粒子と表面との間の作用をシミュレーションする第１ステップとして、粒子は、ボクセルから集められて小面にもたらされる（ステップ３０８）。上述の通り、ボクセルＮ（ｘ）と小面Ｆαとの間の状態ｉの粒子における流束は、以下のように示される。
Γiα（ｘ）＝Ｎi（ｘ）Ｖiα（ｘ）
方程式（１２）
【００５５】
これより、小面Ｆαに向けられている各状態ｉ（ｃiｎα＜０）に関して、ボクセルが小面Ｆαにもたらす粒子の数は、以下のように示される。
2010500654.tif 000010
方程式（１３）
【００５６】
Ｖiα（ｘ）が０値ではない場合、ボクセルのみを集約しなければならない。上述の通り、小面のサイズは、Ｖiα（ｘ）がわずかなボクセルのみのために０ではない値となるよう選択される。Ｖiα（ｘ）とＰf（ｘ）は、非整数値とすることができるため、Γα（ｘ）は、実数として記憶及び処理が行われる。
【００５７】
３．小面から小面への移動
次に、粒子が小面間で移動する（ステップ３１０）。小面Ｆαの入流状態（ｃiｎα＜０）の間の平行六面体Ｇiαが、他の小面Ｆβと交差する場合、すると、小面Ｆαが受けとる状態ｉの粒子における一部は、Ｆβから得られることとなる。特に、小面Ｆαは、小面Ｆβが生成する状態ｉの粒子の一部を、前回の時間増分の間、受けとろうとする。この関係は図１０に示してある。ここでは、小面Ｆβと交差する平行六面体Ｇiαの一部１０００は、小面Ｆαと交差する平行六面体Ｇiβの一部１００５に相当する。上述の通り、交差している部分をＶiα（β）として示す。この項を用いることで、小面Ｆβと小面Ｆαとの間の状態ｉの粒子における流束を以下のように表すことができる。
Γiα（β，ｔ−１）＝Γi（β）Ｖiα（β）／Ｖiα
方程式（１４）
ここで、Γi（β，ｔ−１）は、前回の時間増分の間に小面Ｆβによってもたらされる状態ｉの粒子の量である。これによると、小面Ｆαに向けられている各状態ｉ（ｃiｎα＜０）に関して、小面Ｆαに他の小面がもたらす粒子の数は、以下のように表される。
2010500654.tif 000011
方程式（１５）
そして、小面への状態ｉの粒子における全流束は、以下のように表される。
2010500654.tif 000012
方程式（１６）
【００５８】
小面分布関数とも呼ばれる、小面に対しての状態ベクトルＮ（α）は、ボクセル状態ベクトルの５４のエントリーに対応する５４のエントリーを有する。小面分布関数Ｎ（α）のインプット状態は、粒子の流束を体積Ｖiαで除算して得られるこれらの状態に等しくなる。
ｃiｎα＜０である場合、
Ｎi（α）＝ΓiIN（α）／Ｖiα
方程式（１７）
【００５９】
小面分布関数は、小面からアウトプット流束を生成するためのシミュレーションツールであり、必ずしも実際の粒子を表すとは限らない。正確なアウトプット流束を生成するために、値は、分布関数の他の状態に割当てられる。外部の状態は、内部の状態を取り込むために上述した技術を用いて取り込まれる。
Ｎi（α）＝ΓiOTHER（α）／Ｖ
方程式（１８）
ｃiｎα≧０である場合、ΓiOTHER（α）は、ΓiIN（α）を生成するための上記技術を用いて決定されるが、この入流状態（ｃiｎα＜０）以外の状態に（ｃiｎα≧０）その技術を適用する。別のアプローチでは、ΓiOTHER（α）は、前回のステップからΓiOUT（α）の値を用いて生成され得る。これにより、以下のように表される。
ΓiOTHER（α，ｔ）＝ΓiOUT（α，ｔ−１）
方程式（１９）
【００６０】
平衡状態に関して（ｃiｎα＝０）、Ｖiα及びＶiα（ｘ）ともに０となる。Ｎi（α）についての式において、Ｖiα（ｘ）は分数の分子部（ΓiOTHER（α）ついての式により）で表され、Ｖiαは、分数の分母部（Ｎi（α）ついての式により）で表される。したがって、平衡状態についてのＮｉ（α）は、Ｖiα及びＶiα（ｘ）が０に近づくにつれてＮi（α）の極限として決定される。
【００６１】
速度が０である状態の値（すなわち、静止状態及び状態が（０，０，０，２）、（０，０，０，−２）である場合）は、温度及び圧力の初期条件に基づいてシミュレーション開始時に初期化される。そして、これら値は時間とともに調整される。
【００６２】
４．小面での表面動力学の実行
次に、表面動力学は、上述の４つの境界条件を満たすよう各小面のために実行される（ステップ３１２）。小面のために表面動力学を実行するための手順が、図１１に示されている。まず、小面Ｆαの法線での運動量全てを合わせた運動量は、その小面で粒子全ての運動量Ｐ（α）を決定することによって以下のように決められる。
2010500654.tif 000013
方程式（２０）
全てのｉに関して、これにより、法線の運動量Ｐｎ（α）は、以下のように決定する。
Ｐn（α）＝ｎα・Ｐ（α）
方程式（２１）
【００６３】
法線のこの運動量は、プッシング／プリング技術を用いて排除され（ステップ１１１０）、Ｎn-（α）を生成する。この技術によって、粒子は、法線の運動量のみに影響を及ぼす方向における状態の間で移動する。米国特許番号５，５９４，６７１に開示されているプッシング／プリング技術は、参照によって組み込まれる。
【００６４】
したがって、Ｎn-（α）の粒子は、衝突して、ボルツマン分布Ｎn-β（α）をもたらす（ステップ１１１５）。流体力学を実行することに関して以下に示すように、ボルツマン分布は、衝突ルールのセットをＮn-（α）に適用することによって達成され得る。
【００６５】
小面Ｆαに対する出流束分布は、入流束分布及びボルツマン分布に基づいて決定される（ステップ１１２０）。第１に、入流束分布Γi（α）とボルツマン分布との間の相違は、以下のように決定される。
ΔΓi（α）＝ΓiIN（α）−Ｎn-βi（α）Ｖiα
方程式（２２）
【００６６】
この相違を用いることで、出流束分布は、
ΓiOUT（α）＝Ｎn-βi（α）Ｖiα−．Δ．Γi*（α）
方程式（２３）
ｎαｃi＞０である場合、ｉ＊は、状態ｉとは逆方向である状態である。例えば、状態ｉが（１，１，０，０）である場合、すると、状態ｉ＊は、（−１，−１，０，０）となる。表面摩擦及び他の要因を明らかにすることで、出流束分布は、以下の通りさらに改善することが可能である。
2010500654.tif 000014
方程式（２４）
ｎαｃi＞０である場合、Ｃfは表面摩擦の関数であり、ｔiαはｎαと垂直である第１接線ベクトルであり、ｔ2αはｎα及びＴ1αの両方に垂直である第２接線ベクトルである。そして、ΔＮj,1及びΔＮj,2は、状態ｉのエネルギー（ｊ）及び示されている接線ベクトルに対応する分布関数である。この分布関数は、以下のように決定される。
2010500654.tif 000015
方程式（２５）
この場合、ｊは、エネルギーレベルが１の状態である１、そして、エネルギーレベルが２の状態である２と等しい。
【００６７】
Γiout（α）のための方程式の各項の関数は、以下の通りである。第１項及び第２項は、法線での運動量流束の流束境界条件であって衝突がボルツマン分布を生成することに有効である範囲の条件を表すが、接線運動量の流束異常も含む。第３項及び第４項は、この異常を補正する。その異常は、不適切な衝突を原因として離散性効果又は非ボルツマン構造を起すことがある。最後に、第３項は、表面上の接線運動量流束における所望の変化を表すよう表面摩擦の特定量を加える。摩擦係数Ｃfの生成は、以下に示す。ベクトル操作を含む全ての項は、シミュレーションを開始する前に算出可能な幾何学的な要素であることに留意すべきである。
【００６８】
これにより、接線速度は以下のように決定する。
ｕi（α）＝（Ｐ（α）−Ｐn（α）ｎα）／ρ
方程式（２６）
ここで、ρは、小面分布の密度である。
2010500654.tif 000016
方程式（２７）
【００６９】
前述同様、入流束分布とボルツマン分布との間の相違は以下のように決定される。
ΔΓi（α）＝ΓiIN（α）−Ｎn-βi（α）Ｖiα
方程式（２８）
【００７０】
出流束分布は以下のようになる。
ΓiOUT(α)=Ｎn-βi(α)Ｖiα−ΔΓi*(α)+Ｃf(ｎαｃi)［Ｎn-βi*(α)−Ｎn-βi(α)]Ｖiα
方程式（２９）
これは、従来技術によって決定される出流束分布の２つの第１線に相当するが、接線の異常な流束の補正を必要としない。
【００７１】
どちらのアプローチを用いても、結果として生じる流束分布は、全ての運動量流束の条件を満たす。すなわち、以下のようになる。
2010500654.tif 000017
方程式（３０）
ここで、ｐαは、小面Ｆαでの平衡圧であり、粒子を小面にもたらすボクセルの平均化した密度値及び温度に基づく。また、ｕαは、その小面での平均速度である。
【００７２】
質量境界条件及びエネルギー境界条件を満たすことを保証するために、インプットエネルギーとアウトプットエネルギーとの間の相違は、以下の通り、各エネルギーレベルｊ毎に測定される。
2010500654.tif 000018
方程式（３１）
ここでインデックスｊは、状態ｉのエネルギーを意味する。このエネルギーの相違は、差分項を生成するよう用いられる。
2010500654.tif 000019
方程式（３２）
ｃjiｎα＞０である場合、この差分項は、出流束を変更するよう用いられる。これにより、流束は以下のようになる。
ΓαjiOUTf＝ΓαjiOUT＋δΓαji
方程式（３３）
ｃjiｎα＞０である場合、この動作は、接線運動量が変化しないでいる間、質量流束及びエネルギー流束を訂正する。流量が、小面の近傍でほぼ均一であって、さらに、平衡化状態に近い場合、この調整はわずかなものとなる。調整後に結果として生じる法線の運動量流束は、平衡圧である値に若干変更される。その変更は、近傍の平均特性さらに、近傍における非均一特性又は非平衡特性による訂正に基づく。
【００７３】
５．ボクセルからボクセルへの移動
再度、図３を参照すると、粒子は、３次元の直線格子に沿ってボクセル間を移動する（ステップ３１４）。このボクセルからボクセルへの移動は、小面と作用しないボクセル上（すなわち、表面近くに位置していないボクセル）で実行される動きのみの動作である。典型的なシミュレーションにおいて、大多数のボクセルは、表面と作用するほど表面近くに位置していない。
【００７４】
個々の状態それぞれは、整数倍の速度で移動する粒子を表す。さらに、その移動は、ｘ、ｙ、及びｚからなる３つの次元それぞれにおける格子に沿って移動する。整数倍の速度の整数は、０、±１、そして、±２を含む。速度の符号は、方向を示しており、その方向は、粒子が対応する軸に沿って移動する向きである。
【００７５】
表面と作用しないボクセルに関して、動きの動作は、計算的には非常に単純である。状態の集合全ては、その状態の現在のボクセルからその状態の目的のボクセルへ時間増分ごとの間に移動する。同時に、目的のボクセルにおける粒子は、あるボクセルから、その粒子の目的のボクセルへ移動する。例えば、＋１ｘ及び＋１ｙの方向（１，０，０）に移動するエネルギーレベル１の粒子は、それの現在のボクセルからｘ方向には＋１を超えて移動し、他の方向のボクセルへの移動は０である。粒子は、最終的にはそれの目的のボクセルに移動する前と同じ同じ状態（１，０，０）で到達する。ボクセル内の作用は、他の粒子及び表面との局所的な作用に基づいて粒子数をそのような状態のために変化させようとする。もしそうでない場合、粒子は、格子に沿って同一の速度及び方向で移動し続けようとする。
【００７６】
動きの動作は、一又は二以上の表面と作用するボクセルのためにわずかに複雑となる。これは、小面に移動させられる一又は二以上のわずかな粒子が生じる可能性がある。そのようなわずかな粒子の小面への移動は、ボクセルに残り続けるわずかな粒子を生じさせる。これらわずかな粒子は、小面によって占有されるボクセルへ移動する。例えば、図９を参照すると、ボクセル９０５についての状態ｉの粒子の部分９００が、小面９１０へ移動したとき（ステップ３０８）、残りの部分９１５がボクセル９２０へ移動する。小面９１０は、そのボクセル９２０に位置し、そして、状態ｉの粒子は、そのボクセル９２０から小面９１０へ向けられる。したがって、状態集合が２５に等しく、そして、Ｖiα（ｘ）が０．２５に等しい（すなわち、ボクセルの４分の１が平行六面体Ｇｉαと交差する）場合、すると、６．２５の粒子は小面Ｆαへ移動しようとし、１８．７５の粒子は、小面Ｆαによって占有されるボクセルへ移動しようとする。複数の小面が単一のボクセルと交差することから、一又は二以上の小面によって占有されるボクセルＮ（ｆ）へ移動する状態ｉの粒子の数は、以下のように示される。
2010500654.tif 000020
方程式（３４）
ここで、Ｎ（ｘ）はソースボクセルである。
【００７７】
６．小面からボクセルへの散乱
次に、各小面から外に向かう粒子は、ボクセルへ散乱する（ステップ３１６）。基本的に、このステップは、前述の集めるステップを逆転させたものである。粒子は、ボクセルから小面へこのステップによって移動する。状態ｉにおける小面ＦαからＮ（ｘ）へ移動する粒子の数は、以下のように表される。
2010500654.tif 000021
方程式（３５）
ここで、Ｐf（ｘ）は、一部のボクセルの体積の減少を表す。これにより、各状態ｉに関して、小面からボクセルＮ(x)に向けられている粒子の総数は、以下のように表される。
2010500654.tif 000022
方程式（３６）
【００７８】
小面からボクセルへ粒子を散乱させた後、そのボクセルと周囲のボクセルから移動する粒子とを結びつけ、結果を整数化する。特定のボクセルにおける特定の方向は、アンダーフロー（負となる）かオーバーフロー（８ビットで実装した場合は２５５を超える）のどちらとなる可能性がある。これは、質量、運動量、そして、エネルギーの増加か損失のどちらかを、これら量が許容範囲の値に適合するよう切り捨てられた後にもたらす。そのような現象から保護するために、境界外の質量、運動量、エネルギーは、問題となる状態の切捨てより前に集められる。ある状態が属しているエネルギーに関して、増加した値（アンダーフローに起因）又は損失した値（オーバーフローに起因）に等しい質量は、同じエネルギーを有する無作為に（又は順に）選択される状態に加えられる。そして、その状態は、オーバーフロー又はアンダーフローに影響を受けやすい。この質量及びエネルギーの増加に起因する増加した運動量は、集められて、切捨てに基づいて運動量に追加される。質量のカウンターが０に到達したとき、質量及びエネルギーの両方は、同じエネルギー状態に質量の追加のみを行うことによって訂正される。最後に、運動量は、運動量の蓄積機構（accumulator）が０に戻るまでプッシング／プリング技術を用いて訂正される。
【００７９】
７．流体力学の実行
最後に、流体力学が実行される（ステップ３１８）。このステップは、ミクロ動力学作用又はイントラボクセル作用と呼ばれる。同様に、移流の手順は、インターボクセル作用と呼ばれる。以下に開示するミクロ動学作用は、小面で粒子を衝突させるよう用いることができ、ボルツマン分布を発生させる。
【００８０】
格子ボルツマン方程式モデルにおいて、流体力学は、ＢＧＫ衝突モデルとして既知である特定の衝突オペレータによって保証されている。この衝突モデルは、実際の流体系における分布の動力学に近似させる。この衝突プロセスは、方程式１及び方程式２の右辺によって表すことが可能である。移流のステップの後、流体系の保存量、特に密度、運動量、そして、エネルギーは、方程式３を用いて分布関数から得られる。これらの量により、方程式（２）のｆeqによって示される平衡分布関数は、方程式（４）で十分に特定される。速度ベクトルのセットｃiと重量の選択肢の両方は、方程式２とともに表１にリストされている。方程式２は、巨視的な挙動が正確な流体力学方程式に従うことを保証する。
【００８１】
Ｄ．可変分解能
図１２を参照すると、可変分解能（図６、図７、及び上述してあるような）は、異なるサイズのボクセル（以下、粗いボクセル１２００及び細かいボクセル１２０５とする）を用いる。（２つの異なるサイズを有するボクセルについての以下の説明において、当然のことながら、開示されている技術は、追加の分解能のレベルを提供するよう３又はそれ以上の異なるサイズのボクセルに適用することが可能である。）粗いボクセルと細かいボクセルとの接触面は、可変分解能（ＶＲ）接触面１２１０と呼ばれる。
【００８２】
可変分解能が、表面又はその近くで用いられているとき、小面は、ＶＲ接触面の両側でボクセルと相互に作用し得る。これら小面は、ＶＲ接触小面１２１５（ＦαIC）又はＶＲ細小面１２２０（ＦαIF）に分類される。ＶＲ接触面１２１５は、ＶＲ接触面の粗くなっている側に位置する小面であって、細かいボクセルにまで伸びる粗い平行六面体１２２５を有する小面である。（粗い平行六面体のために、ｃiは、粗いボクセルの寸法にしたがって合わせてあり、さらに、細かい平行六面体のために、ｃiは、細かいボクセルの寸法にしたがって合わせてある。）ＶＲ細小面１２２０は、ＶＲ接触面の細かくなっている側に位置している小面であって、粗いボクセルにまで伸びる細かい平行六面体１２３０を有する小面である。また、接触面に関連する処理には、粗い小面１２３５（ＦαC）及び細かい小面１２４０（ＦαF）との作用を含めてもよい。
【００８３】
ＶＲ小面の両方のタイプに関して、表面動力学は、細かい基準で実行され、上述の通りに動作する。しかしながら、粒子がＶＲ小面より移動、そして、ＶＲ小面へ移動する方法に関して、ＶＲ小面は、他の小面とは異なる。
【００８４】
ＶＲ小面との作用は、図１３に示してある可変分解能手順１３００を用いて操作される。この手順の多くのステップは、非ＶＲ小面との作用に関して上述した同様なステップを用いて実行される。手順１３００は、粗いものへの時間帯のステップ（すなわち、粗いボクセルに対応する時間間隔）の間、実行される。そして、その粗いものへの時間帯（可変分解能手順）は、２つの段階を含む。また、その２つの段階それぞれは、細かいものへの時間帯のステップに対応する。小面の表面動力学は、細かいものへの時間帯それぞれの間に実行される。したがって、ＶＲ接触小面ＦαICは、完全に同じサイズかつ方向である２つの細かい小面としてみなされ、それぞれ、黒い小面ＦαICbと赤い小面ＦαICrと呼ばれる。黒い小面ＦαICbは、粗いものへの時間帯のステップ内の第１の細かいものへの時間帯のステップと関連しており、一方、赤い小面ＦαICrは、粗いものへの時間帯のステップ内の第２の細かいものへの時間帯のステップと関連している。
【００８５】
最初に、粒子は、表面から表面への第１移流段階で小面間を移動する（ステップ１３０２）。粒子は、Ｖ-αβの重み要素を用いて黒い小面ＦαICbから粗い小面ＦβCへ移動する。そのＶ-αβは、細かい平行六面体のブロック化されていない部分（図１２の１２４５参照）を差し引いた粗い平行六面体のブロック化されていない部分（図１２の１２２５参照）の体積に相当する。その粗い平行六面体は、小面Ｆαから伸びており、小面Ｆβの後方にあり、その細かい平行六面体は、小面Ｆαから伸びており、小面Ｆβの後方にある。細かいボクセルについてのｃiの大きさは、粗いボクセルについてのｃiの大きさの２分の１である。上述の通り、小面Ｆαについての平行六面体の体積は、以下のように定義される。
Ｖiα＝｜ｃiｎα｜Ａα
方程式（３７）
【００８６】
したがって、小面の表面積Ａαは、粗い平行六面体と細かい平行六面体との間で変化しないことから、さらに、単位法線ｎαは、常に一つの大きさであることから、小面に対応する細かい平行六面体の体積は、その小面についての対応する粗い平行六面体の体積の２分の１である。
【００８７】
粒子は、Ｖαβの重み要素を用いて粗い小面ＦαCから黒い小面ＦβICbへ移動する。そのＶαβは、細かい平行六面体のブロック化されていない部分の体積に相当する。その細かい平行六面体は、小面Ｆαから伸びており、小面Ｆβの後方にある。
【００８８】
粒子は、Ｖαβの重み要素を用いて赤い小面ＦαICrから粗い小面ＦβCへ移動し、さらに、Ｖ-αβの重み要素を用いて粗い小面ＦαCから赤い小面ＦβICrへ移動する。
【００８９】
粒子は、Ｖαβの重み要素を用いて赤い小面ＦαICrから黒い小面ＦβICbへ移動する。この段階において、黒から赤への移流は、発生しない。さらに、黒い小面及び赤い小面は、連続した時間帯のステップに対応することから、黒から黒への移流（また、赤から赤への移流）は絶対に発生しない。同じ理由から、この段階における粒子は、Ｖαβの重み要素を用いて赤い小面ＦαICrから細かい小面ＦβIF又はＦβFへ移動し、さらに、同じ重み要素を用いて細かい小面ＦαIF又はＦαFから黒い小面ＦαICbへ移動する。
【００９０】
最後に、粒子は、同じ重み要素を用いて細かい小面ＦαIF又はＦαFから他の細かい小面ＦβIF又はＦβFへ移動し、さらに、ＶCαβの重み要素を用いて粗い小面ＦαCから他の粗い小面ＦCへ移動する。そのＶCαβは、粗い平行六面体のブロック化されていない部分の体積に相当する。その粗い平行六面体は、小面Ｆαから伸びており、小面Ｆβの後方にある。
【００９１】
粒子が表面間を移動した後、粒子は、第１の集める段階（ステップ１３０４乃至１３１０）においてボクセルから集められる。粒子は、細かい小面ＦαFのために細かいボクセルから細かい平行六面体を用いて集められ（ステップ１３０４）、また、粗い小面ＦαCのために粗いボクセルから粗い平行六面体を用いて集められる（ステップ１３０６）。次に、粒子は、黒い小面ＦαIRbのために、そして、細かいＶＲ小面ＦαIFのために、粗いボクセル及び細かいボクセルの両方から細かい平行六面体を用いて集められる（ステップ１３０８）。最後に、粒子は、赤い小面ＦαIRrのために粗いボクセルから、粗い平行六面体と細かい平行六面体との間の相違を用いて集められる（ステップ１３１０）。
【００９２】
次に、細かいボクセル又はＶＲ小面と作用する粗いボクセルは、細かいボクセルのコレクション内で急激に増加する（ステップ１３１２）。粗いものへの時間帯の単独のステップ内で、粒子を細かいボクセルに送ろうとする粗いボクセルの状態が、増加する。例えば、小面と交差しない粗いボクセルの適切な状態は、図４のマイクロブロックのように方向付けされている８つの細かいボクセルに増加する。一又は二以上の小面と交差する粗いボクセルの適切な状態は、全て、及び／又は一部の細かいボクセルのコレクションに発展する。その細かいボクセルは、小面と交差していない粗いボクセルの一部に相当する。それの増加に起因する粗いボクセルと細かいボクセルについての粒子密度Ｎi（ｘ）は、等しくなるが、細かいボクセルは、部分要因Ｐfを有することができる。そのＰfは、粗いボクセルの部分要因及び他の細かいボクセルの部分要因とは異なる。
【００９３】
したがって、表面動力学は、細かい小面ＦαIF及びＦαFのために実行され（ステップ１３１４）、さらに、黒い小面ＦαICbのために実行される（ステップ１３１６）。動力学は、図１１及び上記示した手順を用いて実行される。
【００９４】
次に、粒子は、細かいボクセル間で移動する（ステップ１３１８）。その細かいボクセルは、実際の細かいボクセル、そして、粗いボクセルの増加に起因する細かいボクセルを含む。一度、粒子が移動すると、粒子は、細かい小面ＦαIF及びＦαFから細かいボクセルへ散乱する（ステップ１３２０）。
【００９５】
また、粒子は、黒い小面ＦαICbから細かいボクセル（粗いボクセルを増加させることに起因する細かいボクセルも含む）へ散乱する（ステップ１３２２）。表面が存在しないときにボクセルが粒子を受けとろうとする場合、粒子は細かいボクセルへ散乱する。特に、ボクセルが実際の細かいボクセル（粗いボクセルの増加に起因する細かいボクセルは含まない）であるときや、ボクセルＮ（ｘ）を超える一単位速度であるボクセルＮ（ｘ＋ｃi）が、実際の細かいボクセルであるとき、又は、ボクセルＮ（ｘ）を超える一単位速度であるボクセルＮ（ｘ＋ｃi）が、粗いボクセルの増加に起因する細かいボクセルであるとき、粒子は、ボクセルＮ（ｘ）へ散乱する。
【００９６】
最後に、細かいものへの時間帯の第１のステップは、細かいボクセル上で流体力学を実行することによって完了する（ステップ１３２４）。流体力学がボクセルのために実行されると、そのボクセルには、粗いボクセルを増加させたことに起因する細かいボクセルは含まれない（ステップ１３１２）。
【００９７】
手順１３００は、細かいものへの時間帯の第２のステップの間、同じようなステップを実行する。最初に、粒子は、表面から表面への第２の移流段階において表面の間を移動する（ステップ１３２６）。粒子は、黒い小面から赤い小面へ、そして、黒い小面から細かい小面へ、細かい小面から赤い小面へ、さらに、細かい小面から細かい小面へ移動する。
【００９８】
粒子が表面間を移動した後、粒子は、第２の集める段階（ステップ１３２８乃至１３３０）において、ボクセルから集められる。粒子は、赤い小面ＦαIRrのために細かいボクセルから細かい平行六面体を用いて集められる（ステップ１３２８）。また、粒子は、細かい小面ＦαF及びＦαIFのために細かいボクセルから細かい平行六面体を用いて集められる（ステップ１３３０）。
【００９９】
したがって、上述の通り、表面動力学は、細かい小面ＦαIF及びＦαFのために実行され（ステップ１３３２）、そして、粗い小面ＦαCのために実行され（ステップ１３３４）、さらに、赤い小面ＦαICrのために実行される（ステップ１３３６）。
【０１００】
次に、粒子は、細かい分解能を用いてボクセル間を移動する（ステップ１３３８）。これにより、粒子は、細かいボクセル及び粗いボクセルを表す細かいボクセルへ、または、それら細かいボクセルセルより移動する。そして、粒子は、粗い分解能を用いてボクセル間を移動する（ステップ１３４０）。これにより、粒子は、粗いボクセルへ、そして、粗いボクセルより移動する。
【０１０１】
次に、全てを合わせるステップにおいて、粗いボクセルを表す細かいボクセル（すなわち、粗いボクセルの増加に起因する細かいボクセル）が、粗いボクセルと合体する間、粒子は、小面からボクセルへ散乱する（ステップ１３４２）。この全てを合わせるステップにおいて、粒子は、粗い平行六面体を用いて粗い小面から粗いボクセルへや、細かい平行六面体を用いて細かい小面から細かいボクセルへ、細かい平行六面体を用いて赤い小面から細かいボクセル又は粗いボクセルへ、さらに、粗い平行六面体と細かい平行六面体との間の相違を用いて黒い小面からから粗いボクセルへ散乱する。最後に、流体力学は、細かいボクセル及び粗いボクセルのために実行される（ステップ１３４４）。
【０１０２】
多くの実装が開示されているが、それでもやはり、多くの変更が特許請求の範囲の思想及び範囲から逸脱することなく構成させることが理解されるであろう。それに応じ、他の実装は、特許請求の範囲の技術的範囲内に含まれる。
【図１】
2010500654.tif 000023
【図２】
2010500654.tif 000024
【図３】
2010500654.tif 000025
【図４】
2010500654.tif 000026
【図５Ａ−５Ｂ】
2010500654.tif 000027
【図６】
2010500654.tif 000028
【図７】
2010500654.tif 000029
【図８】
2010500654.tif 000030
【図９】
2010500654.tif 000031
【図１０】
2010500654.tif 000032
【図１１】
2010500654.tif 000033
【図１２】
2010500654.tif 000034
【図１３】
2010500654.tif 000035
【国際調査報告】
2010500654.tif 000036