JP2022022999A - 表面アルゴリズムを使用して物理過程をシミュレートするためのコンピュータシステム - Google Patents

Abstract

【課題】流体流の要素をシミュレートするためのコンピュータ実施技法を提供すること。【解決手段】これらの技法は、複数のボクセルの状態ベクトルをメモリに記憶するステップであって、状態ベクトルが、ボクセルにおける複数の可能な運動量状態のうちの特定の運動量状態に対応する複数のエントリを含む、ステップと、ボクセルのサイズおよび向きとは無関係にサイズ指定され方向付けられた少なくとも１つの表面の表現をメモリに記憶するステップと、状態ベクトルに対して相互作用演算を実行するステップであって、相互作用演算が、異なる運動量状態の要素間の相互作用をモデル化する、ステップと、表面の表現に対して表面相互作用演算を実行するステップであって、表面相互作用演算が、表面とボクセルの実質的にすべての要素との間の相互作用をモデル化する、ステップと、要素の新しいボクセルへの移動を反映するために、状態ベクトルに対して移動演算を実行するステップとを含む。【選択図】図１

Description

本説明は、物理過程、例えば流体流（fluid flow）をシミュレートすることに関する。

格子ボルツマン法（ＬＢＭ：Ｌａｔｔｉｃｅ Ｂｏｌｔｚｍａｎｎ Ｍｅｔｈｏｄ）は、車および飛行機などの様々な幾何学形状の周囲の多様で複雑な流体流をシミュレートするために使用される。ソルバ全体の中でも重要な部分の１つは、境界条件を処理するそのアルゴリズム、すなわち、固体壁の表面近くなどの流体計算領域の縁部においてＬＢＭ粒子分布の動力学を処理するアルゴリズムである。このようなアルゴリズムをどれだけ適切に定式化するかが、シミュレーションの結果として得られる精度に直接的な影響を与える。

米国特許第５，８４８，２６０号では、厳密な質量保存、および境界面全体にわたる運動量流束（momentum fluxes）の正確な適用を保証する、体積（volume）ベースの定式化のための新規な技法（‘２６０スリップアルゴリズム）について説明した。

米国特許第５，８４８，２６０号

特許取得済みの手法は大きな成功を収めているが、いくつかの制限がある。例えば、壁の近くで結果として得られる流体量（fluid quantities）は、所望の人工ノイズよりも高くなる。このような人工物は、有効な表面の粗さとして現れ、結果として、数値的散逸が増加し、所望の境界層よりも厚くなる。これは、特に非常に平滑な流線形の筐体の場合、シミュレーションの次のレベルの精度に悪影響を与える可能性がある。第２に、特許取得済みのプロセスの安定域の広さは十分ではなく、したがって、特許取得済みのプロセスは、より高速な流体流のシミュレーションに対するその能力に限界がある。

‘２６０特許に記載のプロセス、すなわちスリップアルゴリズムは、表面表現とともに、格子ボルツマン流体領域の体積表現に基づいている。‘２６０特許では、表面は、それぞれが比表面積および表面法線を有する表面要素の集合として説明されている。‘２６０特許では、境界条件は、隣接する表面要素と流体領域内の粒子との相互作用によって実現される。

しかしながら、本明細書に述べる新規なスリップアルゴリズムでは、境界条件は、次の３つの基本的なプロセスで定義されるように、すべての表面要素と流体領域内の粒子との相互作用によって提供される。

いわゆる平行四辺形／平行六面体に従って、境界面近くの流体領域からの流入分布（incoming distribution）関数（境界面に向かって移動する粒子に対応する）を収集する。表面動力学は、流入分布関数の集合を、境界を通過する質量および運動量流束（momentum fluxes）に関する特定の境界条件を満たす流出分布（outgoing distribution）関数（境界面から反射された粒子に対応する）の集合に変換することを含む。平行四辺形／平行六面体に従って、流出分布関数を、境界近くの流体領域に分布させる。

この手法は、前述の制限に大いに対処する。新しい手法の全体的なフレームワークは‘２６０特許と同様であるが、流入分布関数の集合を流出分布の集合に変換する表面動力学において大きな変化がある。より具体的には、表面を通過するゼロ運動量流束を実現するために、表面動力学で使用されるプロセスが大幅に変更される。厳密なゼロ接線運動量流束を達成することは、正確な境界条件アルゴリズムの基礎であり、その結果、特定の量の運動量流束を後で正確に追加することができる。この変更は、前述の長期にわたる制限を克服する助けとなる。

一態様によれば、流体流（fluid flow）の要素をシミュレートするためのコンピュータ実施方法は、複数のボクセルの状態ベクトルをメモリに記憶するステップであって、状態ベクトルが、ボクセルにおける複数の可能な運動量状態（momentum states）のうちの特定の運動量状態に対応する複数のエントリを含む、ステップと、ボクセルのサイズおよび向きとは無関係にサイズ指定され方向付けられた少なくとも１つの表面の表現をメモリに記憶するステップと、状態ベクトル（state vectors）に対して相互作用演算（interaction operations）を実行するステップであって、相互作用演算が、異なる運動量状態の要素間の相互作用をモデル化する、ステップと、表面の表現に対して表面相互作用演算（surface interaction operations）を実行するステップであって、表面相互作用演算が、表面とボクセルの実質的にすべての要素との間の相互作用をモデル化する、ステップと、要素の新しいボクセルへの移動を反映するために、状態ベクトルに対して移動演算（move operations）を実行するステップとを含む。

さらなる一態様によれば、流体流の要素をシミュレートするためのデータ処理システムは、複数のボクセルの状態ベクトルをメモリに記憶することであって、状態ベクトルが、ボクセルにおける複数の可能な運動量状態のうちの特定の運動量状態に対応する複数のエントリを含む、記憶することと、ボクセルのサイズおよび向きとは無関係にサイズ指定され方向付けられた少なくとも１つの表面の表現をメモリに記憶することと、状態ベクトルに対して相互作用演算を実行することであって、相互作用演算が、異なる運動量状態の要素間の相互作用をモデル化する、相互作用演算を実行することと、表面の表現に対して表面相互作用演算を実行することであって、表面相互作用演算が、表面とボクセルの実質的にすべての要素との間の相互作用をモデル化する、表面相互作用演算を実行することと、要素の新しいボクセルへの移動を反映するために、状態ベクトルに対して移動演算を実行することとをデータ処理システムに行わせるための命令を含む。

さらなる一態様によれば、非一過性のコンピュータ可読媒体は、流体流の要素をシミュレートするためのコンピュータプログラム製品を記憶し、コンピュータプログラム製品は、複数のボクセルの状態ベクトルをメモリに記憶することであって、状態ベクトルが、ボクセルにおける複数の可能な運動量状態のうちの特定の運動量状態に対応する複数のエントリを含む、記憶することと、ボクセルのサイズおよび向きとは無関係にサイズ指定され方向付けられた少なくとも１つの表面の表現をメモリに記憶することと、状態ベクトルに対して相互作用演算を実行することであって、相互作用演算が、異なる運動量状態の要素間の相互作用をモデル化する、相互作用演算を実行することと、表面の表現に対して表面相互作用演算を実行することであって、表面相互作用演算が、表面とボクセルの実質的にすべての要素との間の相互作用をモデル化する、表面相互作用演算を実行することと、要素の新しいボクセルへの移動を反映するために、状態ベクトルに対して移動演算を実行することとをデータ処理システムに行わせるための命令を含む。

上記の態様のうちの１つまたは複数は、本明細書に記載の特徴の中でもとりわけ、以下の特徴のうちの１つまたは複数を含み得る。

表面相互作用演算を実行するための命令は、ファセットと相互作用する少なくとも１つのボクセルの第１の集合から要素を収集し、収集された要素とすべてのファセットとの間の相互作用をモデル化して表面相互作用要素の集合を生成し、表面相互作用要素を、ファセットと相互作用する少なくとも１つのボクセルの第２の集合に散乱させるための命令を含む。状態ベクトルのエントリは、ボクセルの特定の運動量状態における単位体積当たりの要素の密度を表す。状態ベクトルは、特定の運動量状態における単位体積当たりの要素の密度を表すとともにより多くの可能な値を有する複数の整数および／または浮動小数点値のうちの１つまたは複数を含む。

表面相互作用演算は、整数値および／または浮動小数点値のうちの１つまたは複数を使用して実行され、表面相互作用演算は、実数を表す値を使用して実行される。表面相互作用演算は、整数および／または浮動小数点数を使用して実行される。表面相互作用演算は、表面と交差するボクセルを部分ボクセルとして表すことを含む。少なくとも１つの表面の表現は、ボクセルのサイズおよび向きとは無関係にサイズ指定され方向付けられるとともに少なくとも１つの表面を表す複数のファセットを含み、表面相互作用演算は、ファセットとファセットの近くの少なくとも１つのボクセルの要素との間の相互作用をモデル化することをさらに含む。

要素は流体の粒子を表し、ファセットは流体が流れる少なくとも１つの表面を表す。

表面相互作用演算の実行は、流入分布とボルツマン分布との間の差分の計算することと、すべての状態ベクトルの差分から合計運動量を決定することと、決定された差分に基づいて流出分布を生成することとを含む。

上記の態様のうちの１つまたは複数は、本明細書に開示される利点のうちの１つまたは複数を提供し得る。

固体壁面に関する境界条件を処理することに加えて、スリップアルゴリズムは、２つの異なる格子ボルツマンソルバの領域間および２つの異なる参照フレームの領域間の流体境界を処理するために適用することもできる。固体壁とは異なり、境界を横切る質量流束は一般に、ゼロではない。また、運動量流束は、ある領域から別の領域への流れの対流を表す追加の項を含む。これらはすべて、固体表面のための以前のスリップアルゴリズム（‘２６０）の非自明であるが単純な拡張である。

好ましい実施形態の以下の詳細な説明、および特許請求の範囲から、本発明の他の特徴および利点が明らかになる。

新規な表面動力学変換を含む、流体流をシミュレートするためのシステムを示す図である。 表面動力学変換を伴う格子ボルツマンモデルシミュレーションの定式化のための動作を示す流れ図である。 表面動力学変換を伴う格子ボルツマンモデルを使用したシミュレーション動作を示す流れ図である。 動力学変換を示す流れ図である。 ユークリッド空間で表される２つのＬＢＭモデルの速度成分を示す図である（先行技術）。 ユークリッド空間で表される２つのＬＢＭモデルの速度成分を示す図である（先行技術）。 補正されたＣＡＤ図面を使用する物理過程シミュレーションシステムが従う手順の流れ図である。 マイクロブロックの斜視図である（先行技術）。 図１のシステムによって使用される格子構造の図である（先行技術）。 図１のシステムによって使用される格子構造の図である（先行技術）。 可変分解能技法を示す図である（先行技術）。 可変分解能技法を示す図である（先行技術）。 粒子の動きを示す図である（先行技術）。 表面のファセットの影響を受ける領域を示す図である（先行技術）。 新規な表面動力学の流れ図である。 表面動力学を実行するための手順の流れ図である（先行技術）。

本発明の１つまたは複数の実施形態の詳細が、添付の図面および以下の説明に記載されている。本発明の他の特徴、目的、および利点が、説明および図面、ならびに特許請求の範囲から明らかになる。

流体流をシミュレートするための方法の１つは、いわゆる格子ボルツマンモデル（ＬＢＭ）である。ＬＢＭベースの物理過程シミュレーションシステムにおいて、流体流は、分布関数値によって表され、分布関数の時間発展を記述するよく知られた格子ボルツマン方程式（以下の方程式１を参照）を使用して、離散速度の集合で評価される。分布関数は、流動過程と衝突過程の２つの過程を含む。

表面動力学では、正確なゼロ接線運動量流束の境界条件が達成されるように、「流入」粒子分布の集合を「流出」粒子分布の集合にマッピング／変換する。上記のように、‘２６０特許と大きく異なるのは、この特有の部分である。特有のアルゴリズムの違いを以下の表面動力学のステップで説明する前に、最初に以下の表面動力学の概念の概要について説明する。

各表面要素について、

が収集ステップから取得された流入粒子分布関数を表すとする。この関数は、表面に向かう方向を指す速度値ｃ_iを有する（すなわち、

である）表面要素の近くの空間内の単位体積当たりの粒子数に対応し、式中、

は、流れ領域に向かう方向を指す表面要素の単位法線ベクトルを示す。流入分布の集合全体は、表面の方向を指しているすべての可能なｃ_i値を含む（すなわち、

）。

表面動力学では、「流入」分布関数の既知値の値に基づいて、表面要素の「流出」分布関数の値を求める。「流出」分布関数

は、表面から離れる方向を指す速度値ｃ_iを有する（すなわち、

である）表面近くの空間内の単位体積当たりの粒子数に対応する。

象徴的には、表面動力学は次の論理変換、

と定義され、上記の方程式１の左辺は、収集ステップを通して取得した流入分布関数の集合であり、方程式１の右辺は、表面動力学ステップを介して値を求めた流出分布関数の集合である。境界の表面を通る２つの基本的な流束がある。これらは、以下の数式で定義される質量流束および運動量流束である。
質量流束は、

(方程式2)
であり、運動量流束は、

(方程式3)
であり、式中、Ａは表面要素の面積である。

流体境界条件は、一般に、質量Ｍおよび運動量Ｐが、流体ソルバの境界面上で指定された値に従うように定義される。具体的には、粒子は固体壁を貫通できないので、ゼロ質量流束の条件、Ｍ＝０が適用される。表面を通るゼロ接線運動量流束は、

を有することになる。

固体壁では、これは摩擦のない境界条件を有することと同じである。ここでは、

は、表面要素に対して接線方向にある任意のベクトルであり、したがって、

である。

スリップアルゴリズムにおける表面動力学は、指定された値による２つの流束条件が満たされるように、流入分布関数の値を流出分布の値に変換する関係である。ゼロ質量流束条件とゼロ接線運動量流束条件を達成した後、指定された値による全体の質量流束条件および運動量流束条件を満たすための残りのプロセスは、‘２６０特許と新規な手法との間で同じであるので、そのプロセスについて、ここでは説明しない。

２つの流束条件の他にも多くの分布関数があるので、流出分布の値は一意ではない。実際、‘２６０特許および現在の手法は、流入分布関数を流出分布関数とリンクする関係が異なるため異なる値を与える。

ゼロ接線運動量流束のための表面動力学アルゴリズム
ゼロ質量流束およびゼロ接線運動量流束を満たす‘２６０特許における特有な関係は、次の定式化、

によって与えられ、式中、ｉ＊は、ｉに対するパリティ格子速度方向であり、すなわち、ｃ_i＝－ｃ_i*である。平衡分布

および

の計算は、密度および速度のサンプル値に基づく。計算の完全な記述は、‘２６０特許において指定されている。

例えば、等温ＬＢＭソルバの場合、平衡分布は、方程式５、

(方程式5)
によって与えられ、式中、ｗ_iおよびＴ₀は、格子ボルツマンモデルの知られている定数である。ここでは、密度ρおよび（

を伴う）速度ｕは、表面近くの流体領域からサンプリングされる。また、ρは、次の関係、

(方程式6)
を満たすことによってゼロ質量流束条件が満たされるように再スケーリングされる。

式（方程式４）を（方程式５）とともに運動量流束の定義（方程式３）に代入すると、

であると証明することができる。したがって、接線方向に沿った運動量流束、

は、ゼロである。

ここでは、量ｐ＝ρΤ₀は、表面の圧力値を示す。より一般的な境界条件は、ここで説明するゼロ質量流束条件およびゼロ接線運動量流束条件に適切な質量流束および運動量流束を追加することによって実現することができる。

‘２６０特許の定式化における１つの明確な特徴は、方程式（４）に示されているように、流出分布

の、流入分布

への直接の１対１の依存性である。現在はこれが、前述の欠陥または前述の欠陥のうちの少なくともいくつかの主な原因であると理解されている。

次に図１を参照すると、格子ボルツマン（ＬＢ）ベースのシミュレーションを実行するシステム１０は、新規な表面動力学変換３４ｂを含む。この実装におけるシステム１０は、クライアントサーバベースまたはクラウドベースのアーキテクチャに基づいており、大規模並列コンピューティングシステム１２（スタンドアロンまたはクラウドベース）として実装されたサーバシステム１２と、ネットワーク１５を介して結合されたクライアントシステム１４とを含む。サーバシステム１２は、メモリ１８、バスシステム１１、インターフェース２０（例えば、ユーザインターフェース／ネットワークインターフェース／ディスプレイまたはモニタインターフェースなど）、および処理デバイス２４を含む。メモリ１８には、メッシュ準備エンジン３２およびシミュレーションエンジン３４がある。

図１は、メモリ１８内のメッシュ準備エンジン３２を示しており、メッシュ準備エンジンは、サーバ１２とは異なるシステム上で実行されるサードパーティアプリケーションとすることができる。メッシュ準備エンジン３２がメモリ１８内で実行するのか、サーバ１２とは異なるシステム上で実行するのかにかかわらず、メッシュ準備エンジン３２は、ユーザ提供のメッシュ定義３０を受け取り、メッシュ準備エンジン３２は、メッシュを準備し、準備されたメッシュをシミュレーションエンジン３４に送信する（かつ／または記憶する）。

シミュレーションエンジン３４は、新規な表面動力学変換３４ｂを含む衝突相互作用モジュール３４ａと、境界モジュール３４ｃと、移流粒子衝突相互作用モジュール（advection particle collision interaction module）３４ｄと含む。システム１０は、２Ｄおよび／または３Ｄメッシュ（デカルトおよび／または曲線）、座標系、ならびにライブラリを格納するデータリポジトリ３８にアクセスする。

次に図２を参照すると、物理オブジェクトの表現について流体流をシミュレートするためのプロセス４０が示されている。本明細書に述べる例では、物理オブジェクトは翼である。しかしながら、物理オブジェクトは任意の形状とすることができ、具体的には、平面および／または曲面を有することができるので、翼の使用は単に例示的なものである。プロセス４０は、シミュレートされている物理オブジェクトためのメッシュ（またはグリッド）を、例えば、クライアントシステム１４から受け取るか４２、またはデータリポジトリ３８から取得する。他の実施形態では、外部システムまたはサーバ１２のいずれかが、ユーザ入力に基づいて、シミュレートされている物理オブジェクトのためのメッシュを生成する。プロセスでは、取得されたメッシュから幾何学量を事前計算し４４、取得されたメッシュに対応する事前計算された幾何学量を使用して動力学格子ボルツマンモデルシミュレーションを実行する４６。格子ボルツマンモデルシミュレーションは、表面動力学変換、境界モデリング、およびＬＢＭメッシュ内の次のセルへの粒子の移流を含む、粒子分布の発展のシミュレーション４６を含む。

図３を参照すると、シミュレーションプロセス４６は、例えば、新規な表面動力学変換３４ｂ（図１）に適合するように修正された格子ボルツマンプロセス（ＬＢｐ）に従って、粒子分布の発展をシミュレートする。プロセス４６（図２を参照）は、衝突演算を実行し４６ａ、すべての粒子のすべての衝突演算の結果を動力学変換し４６ｂ、その後に境界モデリング４６ｃ、およびＬＢＭ空間内の次のセルへの粒子の移流４６ｄが続く。

次に図４を参照すると、粒子の時間発展をシミュレートするシミュレーションプロセス４６は、一般に従来の衝突演算４６ａ、境界モデリング４６ｃ、および移流４６ｄのプロセスを実行する。しかしながら、‘２６０特許に記載されている流入粒子分布の集合の流出粒子分布の集合への変換は、流出粒子分布の各集合を定式化するために、すべての流入粒子分布４６ｂの結果を取る動力学変換プロセスに置き換えられる。動力学変換４６ｂは、隣接するメッシュ位置から流入分布の集合を収集するのではなく、衝突演算から、メッシュ位置から分布のすべての流入集合を収集すること５２を含む。この新しい定式化は、上記の‘２６０特許における１対１の依存関係を回避し、その依存関係を、各「流出」分布関数がすべての「流入」分布に依存するように置き換える５４。

方程式４の代わりに、新規な手法で与えられた新しい関係が、以下の数学的定式化、

によって定義され、式中、δｆ_iは、

(方程式8)
によって定義される。量ｈは、

(方程式9)
によって定義される。

（方程式９）において、平衡分布は、（方程式５）と同じ方法で定義される。

の構築では、（方程式６）によって与えられるのと同じ質量流束制約を満たすために密度再スケーリングを適用する以外に、新しい定式化において平衡運動量流束の接線部分に追加の制約、

(方程式10)
が必要であり、式中、単位ベクトル

は、サンプリングされた速度方向ｕに平行である。後者は、表面要素に投影されるので

、したがって、

である。制約（方程式１０）は、ｕの大きさを再スケーリングすることによって実現することができる。

‘２６０特許の手法と同様に、（方程式７）～（方程式９）によって定義される新規な手法も、ゼロ質量流束およびゼロ接線運動量流束の境界条件を達成することが示され得る。

一見より複雑な手順であるが、新しい定式化における主な違いは、流出分布と流入分布との間の１対１の関係を回避することである。実際、（方程式７）～（方程式９）に示すように、流出分布

は、（方程式９）で表されるそれらのモーメント合計を通じたすべての流入分布に依存する。

新しい定式化は、前述の欠点を大幅に克服し、ノイズの少ない表面を提供し、境界近くでより高い流体速度をサポートする可能性がある。

詳細な例
以下の図７に述べる手順では、シミュレーション空間を構成するための識別された空き空間を有するＣＡＤ図面を使用する、流れシミュレーションプロセスについて説明する。前の図５および図６、ならびに図８、図９Ａ～図９Ｂ、図１０および図１１、ならびに図１２、図１３、および図１５では、これらの図が上記の参照特許に出現するので、これらの図はそれぞれ、先行技術としてラベル付けされている。

しかしながら、上記の参照特許には本明細書に記載の動力学変換プロセス４６ｂは記載されていないので、上記の特許に出現する図は、動力学変換プロセス４６ｂを使用する流れシミュレーションに対して加えられる変更を考慮していない。

シミュレーション空間のモデル化
ＬＢＭベースの物理過程シミュレーションシステムにおいて、流体流は、離散速度ｃ_iの集合で評価される分布関数値ｆ_iによって表される。分布関数の動力学は、方程式Ｉ．１、

方程式(I.1)
によって支配され、式中、

は、平衡分布関数として知られており、

と定義される。

方程式（Ｉ．１）は、分布関数ｆ_iの時間発展を表すよく知られた格子ボルツマン方程式である。左辺は、いわゆる「流動過程」による分布の変化を表す。流動過程は、流体のポケットが、あるメッシュ位置で出発し、次いで複数の速度ベクトルのうちの１つに沿って次のメッシュ位置に移動するときの過程である。その時点で、「衝突因子」、すなわち、出発している流体のポケットに対する近傍の流体のポケットの影響を計算する。流体は別のメッシュ位置に移動することしかできず、したがって、すべての速度のすべての成分が共通の速度の倍数になるように、速度ベクトルを適切に選択する必要がある。

第１の方程式の右辺は、流体のポケット間の衝突による分布関数の変化を表す前述の「衝突演算子」である。衝突演算子の特定の形態は、Ｂｈａｔｎａｇａｒ、Ｇｒｏｓｓ、およびＫｒｏｏｋ（ＢＧＫ）演算子によるものである。衝突演算子は、分布関数を、「平衡」形態である第２の方程式によって与えられる規定値になるようにする。

ＢＧＫ演算子は、衝突の詳細に関係なく、分布関数が、衝突

を通じて｛ｆ^eq（ｘ，ｖ，ｔ）｝によって与えられる明確に定義された局所平衡に近づくという物理的論拠に従って構築され、式中、パラメータτは、衝突を通じた平衡までの特性緩和時間を表す。粒子（例えば、原子または分子）を扱う場合、緩和時間は通常、定数と見なされる。

このシミュレーションから、質量ρおよび流体速度ｕなどの従来の流体変数は、以下の方程式（Ｉ．３）における単純和として得られる。

方程式(I.3)
式中、ρ、ｕ、およびＴはそれぞれ、流体の密度、速度、および温度であり、Ｄは離散化された速度空間の次元である（必ずしも物理空間の次元と等しいとは限らない）。

対称性を考慮することにより、速度値の集合は、構成空間にまたがるときに一定の格子構造を形成するように選択される。このような離散系の動力学は、形態
ｆ_i（ｘ＋ｃ_i，ｔ＋１）－ｆ_i（ｘ，ｔ）＝Ｃ_i（ｘ，ｔ）
を有するＬＢＥに従い、式中、衝突演算子は通常、上記のようなＢＧＫ形態を取る。平衡分布形態を適切に選択することにより、格子ボルツマン方程式が正しい水力学および熱水力学を生み出すことを理論的に示すことができる。すなわち、ｆ_i（ｘ，ｔ）から導出される水力学的モーメントは、巨視的極限ではナビエ－ストークス方程式に従う。これらのモーメントは、上記の方程式（Ｉ．３）によって定義される。

ｃ_iとｗ_iの集合値は、ＬＢＭモデルを定義する。ＬＢＭモデルは、スケーラブルなコンピュータプラットフォームに効率的に実装され、時間非定常流および複雑な境界条件に対して優れたロバスト性を伴って実行することができる。

ボルツマン方程式から流体系に対する巨視的な運動方程式を取得する標準的な技法は、完全なボルツマン方程式の逐次近似を取るチャップマン－エンスコッグ法である。流体系では、密度のわずかな乱れが音速で伝わる。気体系では一般に、音速は温度によって決定される。流れの圧縮率の影響の重要性は、特性速度と音速との比によって測定され、この比は、マッハ数として知られている。

流体流のシミュレーションを実施するための動力学変換４６ｂを含む、ＬＢＭベースのシミュレーションシステムの一般的な議論を以下に示す。ＬＢＭベースの物理過程シミュレーションシステムのさらなる説明については、参照として組み込まれる上記の米国特許を参照されたい。

図５を参照すると、第１のモデル（２Ｄ－１）２００は、２１の速度を含む２次元モデルである。これらの２１の速度のうち、１つの速度（２０５）は、移動していない粒子を表し、３組の４つの速度は、格子のｘ軸またはｙ軸のいずれかに沿って正方向または負方向のいずれかに、正規化速度（ｒ）（２１０～２１３）、正規化速度の２倍（２ｒ）（２２０～２２３）、または正規化速度の３倍（３ｒ）（２３０～２３３）のいずれかで移動する粒子を表し、２組の４つの速度は、ｘ格子軸およびｙ格子軸の両方に対して正規化速度（ｒ）（２４０～２４３）または正規化速度の２倍（２ｒ）（２５０～２５３）で移動している粒子を表す。

図６を参照すると、第２のモデル（３Ｄ－１）２６０が示されており、このモデルは、各速度が図１０の矢印のうちの１つによって表される３９の速度を含む３次元モデルである。これらの３９の速度のうち、１つの速度は、移動していない粒子を表し、３組の６つの速度は、格子のｘ軸、ｙ軸、またはｚ軸に沿って正方向または負方向のいずれかに正規化速度（ｒ）、正規化速度の２倍（２ｒ）、または正規化速度の３倍（３ｒ）のいずれかで移動している粒子を表し、８つの速度は、ｘ、ｙ、ｚの３つの格子軸すべてに対して正規化速度（ｒ）で移動している粒子を表し、１２の速度は、ｘ、ｙ、ｚの２つの格子軸に対して正規化速度の２倍（２ｒ）で移動している粒子を表す。

１０１の速度を含む３Ｄ－２モデル、および３７の速度を含む２Ｄ－２モデルなど、より複雑なモデルを使用することもできる。１０１の速度のから成る３次元モデル３Ｄ－２の場合、１つの速度は、移動していない粒子（グループ１）を表し、３組の６つの速度は、格子のｘ軸、ｙ軸、またはｚ軸に沿って正方向または負方向のいずれかに正規化速度（ｒ）、正規化速度の２倍（２ｒ）、または正規化速度の３倍（３ｒ）のいずれかで移動している粒子（グループ２、４、および７）を表し、３組の８つの速度は、ｘ、ｙ、ｚの３つの格子軸すべてに対して正規化速度（ｒ）、正規化速度の２倍（２ｒ）、または正規化速度の３倍（３ｒ）で移動している粒子（グループ３、８、および１０）を表し、１２の速度は、ｘ、ｙ、ｚ格子軸のうちの２つに対して正規化速度の２倍（２ｒ）で移動している粒子（グループ６）を表し、２４の速度は、ｘ、ｙ、ｚ格子軸のうちの２つに対して正規化速度（ｒ）および正規化速度の２倍（２ｒ）で移動しており、残りの軸に対して移動していない粒子（グループ５）を表し、２４の速度は、ｘ、ｙ、ｚ格子軸のうちの２つに対して正規化速度（ｒ）で移動しており、残りの軸に対して正規化速度の３倍（３ｒ）で移動している粒子（グループ９）を表す。

３７の速度から成る２次元モデル２Ｄ－２の場合、１つの速度は、移動していない粒子（グループ１）を表し、３組の４つの速度は、格子のｘ軸またはｙ軸のいずれかに沿って正方向または負方向のいずれかに正規化速度（ｒ）、正規化速度の２倍（２ｒ）、または正規化速度の３倍（３ｒ）のいずれかで移動している粒子（グループ２、４、および７）を表し、２組の４つの速度は、ｘ格子軸およびｙ格子軸の両方に対して正規化速度（ｒ）または正規化速度の２倍（２ｒ）で移動している粒子を表し、８つの速度は、ｘ格子軸およびｙ格子軸のうちの一方に対して正規化速度（ｒ）で移動しており、他方の軸に対して正規化速度の２倍（２ｒ）で移動している粒子を表し、８つの速度は、ｘ格子軸およびｙ格子軸のうちの一方に対して正規化速度（ｒ）で移動しており、他方の軸に対して正規化速度の３倍（３ｒ）で移動している粒子を表す。

上記のＬＢＭモデルは、２次元と３次元の両方での流れの数値シミュレーションのための、特定のクラスの効率的でロバストな離散速度運動モデルを提供する。この種類のモデルは、特定の離散速度の集合、およびそれらの速度に関連する重みを含む。速度は、特に格子ボルツマンモデルとして知られる種類の離散速度モデルの正確で効率的な実装を容易にする、速度空間内のデカルト座標のグリッド点と一致する。このようなモデルを使用して、流れを高い充実度でシミュレートすることができる。

図７を参照すると、流体流などの物理過程をシミュレートするための手順２７０による物理過程シミュレーションシステムの動作が説明されている。流れのシミュレーションの前に、上記のＣＡＤ図面を使用して、シミュレーション空間がボクセルの集合体としてモデル化される（２７２）。シミュレーション空間は、コンピュータ支援設計（ＣＡＤ）プログラム、およびＣＡＤで生成された図面のギャップ補正処理を使用して生成される。例えば、ＣＡＤプログラムを使用して、風洞に配置された翼を描画することができる。

格子の分解能は、シミュレートされるシステムのレイノルズ数に基づいて選択することができる。レイノルズ数は、流れの粘度（ν）、流れの中の物体の特性長（Ｌ）、および流れの特性速度（ｕ）に関連する。
Ｒｅ＝ｕＬ／ν 方程式（Ｉ．４）

物体の特性長は、物体の大きいスケールの特徴を表す。例えば、マイクロデバイスの周囲の流れがシミュレートされる場合、マイクロデバイスの高さを特性長と見なすことができる。物体の小さい領域（例えば、自動車のサイドミラー）の周囲の流れが関心対象である場合、シミュレーションの分解能を上げるか、または分解能を上げた領域を関心領域の周囲で使用することができる。格子の分解能が上がるにつれて、ボクセルの寸法は減少する。

状態空間は、ｆ_i（ｘ，ｔ）として表され、式中、ｆ_iは、時間ｔにおいて３次元ベクトルｘで示される格子サイトにある、状態ｉの要素または粒子の単位体積当たりの数（すなわち、状態ｉの粒子の密度）を表す。既知の時間増分の場合、粒子の数は単にｆ_i（ｘ）と示される。格子サイトのすべての状態の組合せは、ｆ（ｘ）と示される。

状態の数は、各エネルギーレベル内の可能な速度ベクトルの数によって決定される。速度ベクトルは、ｘ、ｙ、およびｚの３次元を有する空間内の整数の線形速度である。多重種シミュレーションでは、状態の数が増加する。

各状態ｉは、特定のエネルギーレベル（すなわち、エネルギーレベル０、１、または２）における異なる速度ベクトルを表す。各状態の速度ｃ_iは、次のように、３次元の各次元における「速度」を用いて示される。
ｃ_i＝（ｃ_ix，ｃ_iv，ｃ_iz） 方程式（Ｉ．５）

エネルギーレベル０の状態は、どの次元においても移動していない停止した粒子を表し、すなわち、ｃ_stopped＝（０，０，０）である。エネルギーレベル１の状態は、３次元のうちの１次元において±１の速度を有し、他の２次元において０の速度を有する粒子を表す。エネルギーレベル２の状態は、３次元すべてにおいて±１の速度を有する粒子、または３次元のうちの１次元において±２の速度を有し、他の２次元において０の速度を有する粒子を表す。

３つのエネルギーレベルの可能な順列をすべて生成すると、合計３９の可能な状態（１つのエネルギー０状態、６つのエネルギー１状態、８つのエネルギー３状態、６つのエネルギー４状態、１２のエネルギー８状態、および６つのエネルギー９状態）が得られる。

各ボクセル（すなわち、各格子サイト）は、状態ベクトルｆ（ｘ）によって表される。状態ベクトルは、ボクセルのステータスを完全に定義し、３９個のエントリを含む。３９個のエントリは、１つのエネルギー０状態、６つのエネルギー１状態、８つのエネルギー３状態、６つのエネルギー４状態、１２のエネルギー８状態、および６つのエネルギー９状態に対応する。システムは、この速度集合を使用することによって、達成された平衡状態ベクトルのためのマクスウェル－ボルツマン統計を作り出すことができる。

処理効率のために、ボクセルは、マイクロブロックと呼ぶ２×２×２の体積にグループ化される。マイクロブロックは、ボクセルの並列処理を可能にして、データ構造に関連するオーバーヘッドを最小限に抑えるように編成される。マイクロブロック内のボクセルの省略表記は、Ｎ_i（ｎ）と定義され、式中、ｎは、マイクロブロック内の格子サイトの相対位置を表し、ｎ∈｛０，１，２，．．．，７｝である。

図８に、マイクロブロックを示す。

図９Ａおよび図９Ｂを参照すると、表面Ｓ（図９Ａ）を、シミュレーション空間（図９Ｂ）内にファセットＦ_αの集合体として表している。
Ｓ＝｛Ｆ_α｝ 方程式（Ｉ．６）
式中、αは、特定のファセットを挙げるインデックスである。ファセットは、ボクセル境界に制限されないが、ファセットが比較的少数のボクセルに影響を与えるように、ファセットのサイズは、典型的にはファセットに隣接するボクセルのサイズと同程度またはそれよりもわずかに小さい。ファセットには、表面動力学を実施する目的で特性が割り当てられる。具体的には、各ファセットＦ_αは、単位法線（ｎ_α）、表面積（Ａ_α）、中心位置（ｘ_α）、およびファセットの表面動力学特性を表すファセット分布関数（ｆ_i（α））を有する。全エネルギー分布関数ｑ_i（α）は、ファセットとボクセルの相互作用の流れ分布と同じ方法で扱われる。

図１０を参照すると、処理効率を向上させるために、シミュレーション空間の異なる領域内で異なるレベルの分解能を使用することができる。典型的には、物体３２２の周囲の領域３２０が最大関心対象であり、したがって、この領域が、最高分解能でシミュレートされる。粘度の影響は物体からの距離とともに減少するので、レベルが減少していく分解能（すなわち、拡大されたボクセル体積）を使用して、物体３２２からの距離が増加する間隔で位置する領域３２４、３２６をシミュレートする。

同様に、図１１に示すように、より低レベルの分解能を使用して、物体３４２のそれほど重要でない特徴の周囲の領域３４０をシミュレートすることができ、一方、最高レベルの分解能を使用して、物体３４２の最も重要な特徴（例えば、先端面および後端面）の周囲の領域３４４をシミュレートする。外側の領域３４６は、最低レベルの分解能および最大のボクセルを使用してシミュレートされる。

Ｃ．ファセットの影響を受けるボクセルの識別
再び図１１を参照すると、シミュレーション空間がモデル化されると（２７２）、１つまたは複数のファセットの影響を受けるボクセルが識別される（２７４）。ボクセルは、ファセットによって様々な形で影響を受ける場合がある。第１に、１つまたは複数のファセットが交差するボクセルは、交差していないボクセルと比較してボクセルの体積が減少しているという点で影響を受ける。これは、ファセットおよびファセットによって表される表面下にある材料がボクセルの一部を占有するために発生する。小数係数（fractional factor）Ｐ_f（ｘ）は、ファセットの影響を受けないボクセルの部分（すなわち、流れをシミュレートする対象となっている流体または他の材料によって占有される可能性のある部分）を示す。交差していないボクセルの場合、Ｐ_f（ｘ）は１に等しい。

粒子をファセットに移送するかまたはファセットから粒子を受け取ることによって１つまたは複数のファセットと相互作用するボクセルも、ファセットの影響を受けるボクセルとして識別される。ファセットが交差するボクセルはすべて、ファセットから粒子を受け取る少なくとも１つの状態と、ファセットに粒子を移送する少なくとも１つの状態とを含むことになる。ほとんどの場合、追加のボクセルもこのような状態を含む。

図１２を参照すると、各状態ｉが０以外の速度ベクトルｃ_iを有する場合、ファセットＦ_αは、平行六面体Ｇ_iαによって画定される領域から粒子を受け取るか、その領域に粒子を移送し、平行六面体Ｇ_iαは、速度ベクトルｃ_iとファセットの単位法線ｎ_αとのベクトルドット積（｜ｃ_iｎ_α｜）の大きさによって定義された高さと、ファセットの表面積Ａ_αによって定義される底面とを有し、その結果、平行六面体のＧ_iαの体積Ｖ_iαは、
Ｖ_iα＝｜ｃ_iｎ_α｜Ａ_α 方程式（Ｉ．７）
に等しい。

ファセットＦ_αは、状態の速度ベクトルがファセットの方に向いている（｜ｃ_iｎ_α｜＜０）場合は体積Ｖ_iαから粒子を受け取り、状態の速度ベクトルがファセットから離れる方に向いている（｜ｃ_iｎ_α｜＞０）場合はその領域に粒子を移送する。以下に説明するように、内部の隅部などの非凸面特徴の近傍で起こり得る状態である、別のファセットが平行六面体のＧ_iαの一部を占めるとき、この式は修正される。

ファセットＦ_αの平行六面体のＧ_iαは、複数のボクセルの一部または全部と重なることがある。ボクセルまたはその一部の数は、ボクセルのサイズに対するファセットのサイズ、状態のエネルギー、および格子構造に対するファセットの向きに依存する。影響を受けるボクセルの数は、ファセットのサイズとともに増加する。したがって、上述のように、ファセットのサイズは、通常、ファセットの近くに位置するボクセルのサイズと同程度またはそれよりも小さくなるように選択される。

平行六面体Ｇ_iαと重なるボクセルＮ（ｘ）の一部は、Ｖ_iα（ｘ）と定義される。この項を使用すると、ボクセルＮ（ｘ）とファセットＦ_αとの間を移動する状態ｉの粒子の流束Γ_iα（ｘ）は、ボクセル内の状態ｉの粒子の密度（Ｎ_i（ｘ））と、ボクセルと重なる領域の体積（Ｖ_iα（ｘ））との積に等しくなる。
Γ_iα（ｘ）＝Ｎ_i（ｘ）Ｖ_iα（ｘ） 方程式（Ｉ．８）

平行六面体のＧ_iαが１つまたは複数のファセットと交差する場合、次の条件、
Ｖ_iα＝ΣＶ_iα（ｘ）＋ΣＶ_iα（β） 方程式（Ｉ．９）
は真であり、式中、第１の和は、Ｇ_iαと重なるすべてのボクセルに相当し、第２の項はＧ_iαと交差するすべてのファセットに相当する。平行六面体Ｇ_iαが別のファセットと交差していない場合、この式は、
Ｖ_iα＝ΣＶ_iα（ｘ） 方程式（Ｉ．１０）
に縮小する。

Ｄ．シミュレーションの実行
１つまたは複数のファセットの影響を受けるボクセルが識別されると（２７４）、シミュレーションを開始するためにタイマが初期化される（２７６）。シミュレーションの各時間増分中に、ボクセルからボクセルへの粒子の移動が、粒子と表面ファセットとの相互作用に相当する移流段階（２７８～２８６）によってシミュレートされる。次に、衝突段階（２８８）は、各ボクセル内の粒子の相互作用をシミュレートする。その後、タイマが増分される（２９０）。増分されたタイマが、シミュレーションが完了したことを示さない場合（２９２）、移流段階および衝突段階（２７８～２９０）が繰り返される。増分されたタイマが、シミュレーションが完了したことを示す場合（２９２）、シミュレーションの結果が、記憶および／または表示される（２９４）。

１．表面の境界条件
表面との相互作用を正しくシミュレートするために、各ファセットは４つの境界条件を満たす。第１に、ファセットが受け取る粒子の合計質量は、ファセットが移送した粒子の合計質量に等しい（すなわち、ファセットへの正味質量流束は０に等しい）。第２に、ファセットが受け取る粒子の合計エネルギーは、ファセットが移送した粒子の合計エネルギーに等しい（すなわち、ファセットへの正味エネルギー流束は０に等しい）。これらの２つの条件は、各エネルギーレベル（すなわち、エネルギーレベル１および２）での正味質量流束が０に等しいことを要求することによって満たすことができる。

他の２つの境界条件は、ファセットと相互作用する粒子の正味運動量に関連する。本明細書では滑り面と呼ぶ表面摩擦のない表面の場合、正味接線運動量流束は０に等しく、正味法線運動量流束はファセットでの局所圧力に等しい。したがって、ファセットの法線ｎ_αに垂直な受け取った合計運動量の成分と、ファセットの法線ｎ_αに垂直な移送した合計運動量の成分と（すなわち、接線成分）は等しく、一方、ファセットの法線ｎ_αに平行な受け取った合計運動量の成分と、ファセットの法線ｎ_αに平行な移送した合計運動量の成分と（すなわち、法線成分）の差分がファセットでの局所圧力に等しい。非滑り面の場合、表面の摩擦によって、ファセットが移送する粒子の合計接線運動量が、ファセットが受け取る粒子の合計接線運動量と比較して、摩擦量に関連する係数だけ減少する。

２．ボクセルからファセットへの収集
粒子と表面との間の相互作用をシミュレートする際、粒子が、ボクセルから収集され、ファセットに提供される（２７８）。上述のように、ボクセルＮ（ｘ）とファセットＦ_αとの間の状態ｉの粒子の流束は、
Γ_iα（ｘ）＝Ｎ_i（ｘ）Ｖ_iα（ｘ） 方程式（Ｉ．１１）
である。

この式から、各状態ｉがファセットＦ_αの方に向いている（ｃ_iｎ_α＜０）場合、ボクセルによってファセットＦ_αに提供される粒子の数は、
Γ_iαV→F＝Σ_XΓ_iα（ｘ）＝Σ_XＮ_i（ｘ）Ｖ_iα（ｘ） 方程式（Ｉ．１２）
である。

Ｖ_iα（ｘ）が０以外の値を有するボクセルのみが合計される。上述のように、ファセットのサイズは、Ｖ_iα（ｘ）が少数のボクセルに対してのみ０以外の値を有するように選択される。Ｖ_iα（ｘ）およびＰ_f（ｘ）は非整数値を有することがあるので、Γ_α（ｘ）は実数として記憶および処理される。

３．ファセットからファセットへの移動
次に、粒子は、ファセット間を移動する（２８０）。ファセットＦ_αが流入状態（ｃ_iｎ_α＜０）にある平行六面体Ｇ_iαが別のファセットＦ_βと交差する場合、ファセットＦ_αが受け取る状態ｉの粒子の一部はファセットＦβから流入することになる。具体的には、ファセットＦ_αは、前の時間増分中にファセットＦ_βによって作り出された状態ｉの粒子の一部を受け取ることになる。この関係を図１５に示しており、平行六面体Ｇ_iαのうちのファセットＦ_βと交差する部分３８０は、平行六面体Ｇ_iβのうちのファセットＦ_αと交差する部分３８２に等しい。上述のように、交差部分はＶ_iα（β）で表される。この項を使用すると、ファセットＦ_βとファセットＦ_αとの間の状態ｉの粒子の流束は、
Γ_iα（β，ｔ－１）＝Γ_i（β）Ｖ_iα（β）／Ｖ_iα 方程式（Ｉ．１３）
と記述することができ、式中、Γ_i（β，ｔ－１）は、前の時間増分中にファセットＦ_βによって作り出された状態ｉの粒子の測定量である。この式から、各状態ｉがファセットＦ_αの方に向いている（ｃ_iｎ_α＜０）場合、他のファセットによってファセットＦ_αに提供される粒子の数は、
Γ_iαF→F＝Σ_βΓ_iα（β）＝Σ_βΓ_i（β，ｔ－１）Ｖ_iα（β）／Ｖ_iα 方程式（Ｉ．１４）
であり、ファセット内への状態ｉの粒子の全流束は、
Γ_iIN（α）＝Γ_iαF→F＋Γ_iαF→F＝Σ_xＮ_i（ｘ）Ｖ_iα＋Σ_βΓ_i（β，ｔ－１）Ｖ_iα（β）／Ｖ_iα 方程式（Ｉ．１５）
である。

ファセット分布関数とも呼ぶファセットの状態ベクトルＮ（α）は、ボクセル状態ベクトルのＭ個のエントリに対応するＭエントリを有する。Ｍは離散格子速度の数である。ファセット分布関数Ｎ（α）の入力状態は、これらの状態になる粒子の流束を体積Ｖ_iαで除算したものに等しく設定され、ｃ_iｎ_α＜０の場合、
Ｎ_i（α）＝Γ_iIN（α）／Ｖ_iα 方程式（Ｉ．１６）
である。

ファセット分布関数は、ファセットからの出力流束を生成するためのシミュレーションツールであり、必ずしも実際の粒子を表すものではない。正確な出力流束を生成するために、分布関数の他の状態に値を割り当てる。上述した内向き状態を代入するための技法を使用して、外向き状態を代入すると、ｃ_iｎ_α≧０の場合、
Ｎ_i（α）＝Γ_iOTHER（α）／Ｖ_iα 方程式（Ｉ．１７）
であり、式中、Γ_iOTHER（α）は、上述したΓ_iIN（α）を生成するための技法を使用して求められるが、流入状態（ｃ_iｎ_α＜０）以外の状態（ｃ_iｎ_α≧０）にこの技法を適用する。代替の手法では、前のタイムステップからのΓ_iOUT（α）の値を使用して、Γ_iOTHER（α）を、
Γ_iOTHER（α，ｔ）＝Γ_iOUT（α，ｔ－１） 方程式（Ｉ．１８）
のように生成することができる。

平行状態（ｃ_iｎ_α＝０）の場合、Ｖ_iαとＶ_iα（ｘ）はいずれも０である。Ｎ_i（α）を求める式では、Ｖ_iα（ｘ）は、（Γ_iOTHER（α）を求める式からの）分子に現れ、Ｖ_iαは、（Ｎ_i（α）を求める式からの）分母に現れる。したがって、平行状態を求めるＮ_i（α）は、Ｖ_iαおよびＶ_iα（ｘ）が０に近づくときのＮ_i（α）の極限として決定される。速度が０の状態（すなわち、静止状態、ならびに状態（０，０，０，２）および状態（０，０，０，－２））の値は、シミュレーションの開始時に、温度および圧力の初期条件に基づいて初期化される。その後、これらの値は、時間の経過とともに調整される。

４．ファセット表面動力学の実行
次に、各ファセットが上記の境界条件を満たすように表面動力学を実行する（２８２）。図１４に、ファセットの表面動力学を実行するための手順３９０を示す。

２７８の間に、速度および密度を、ボクセルからファセットＦ_αまでサンプリングする。次いで、速度を、表面に沿って投影する（すなわち）

である。サンプリングした密度および投影した速度に基づいてボルツマン平衡分布を計算し３９２、密度および速度を、それぞれ方程式６および方程式１０で指定された制約を満たすようにスケーリングする。次いで、結果として得られる密度および速度を使用して、新しいボルツマン平衡分布を計算する（３９２）。流入分布と新しいボルツマン分布との間の差分（３９４）、およびすべての流入状態とそれに対応するボルツマン分布（３９６）との間の差分の合計運動量は、次のように計算される。

方程式(I.19)

運動量差分を、ファセットの接線方向に沿って投影する（３９８）。

上記の運動量差分およびボルツマン分布から、ゼロ接線流束

を満たすことによって、完全なスリップ境界条件を満たすように流出流束を計算する（３９９）。

方程式(I.21)

表面摩擦および他の因子を説明するために、流出流束分布をさらに、

と改良することができ、式中、Ｃ_fは表面摩擦係数である。完全な質量保存およびエネルギー保存に必要とされる様々なエネルギーレベルの格子に表面摩擦および補正を適用することに関する詳細な説明は、‘２６０特許に示されている。

５．ボクセルからボクセルへの移動
再び図７を参照すると、３次元直線格子に沿って、粒子をボクセル間で移動させる（２８４）。このボクセルからボクセルへの移動は、ファセットと相互作用しないボクセル（すなわち、表面の近くに位置していないボクセル）に対して実行される唯一の移動演算である。典型的なシミュレーションでは、表面と相互作用するのに十分なほど表面の近くに位置していないボクセルが、ボクセルの大多数を構成する。

各分離状態は、３次元のｘ、ｙ、およびｚのそれぞれにおいて、格子に沿って整数速度で移動する粒子を表す。整数速度は、０、±１、および±２を含む。速度の符号は、粒子が対応する軸に沿って移動している方向を示す。

表面と相互作用しないボクセルの場合、移動演算は、計算上極めて単純である。毎時間増分中に、状態の母集団全体は、その現在のボクセルからその移動先ボクセルに移動される。同時に、移動先ボクセルの粒子は、そのボクセルからそれ自体の移動先ボクセルに移動される。例えば、＋１ｘおよび＋１ｙ方向（１，０，０）に移動しているエネルギーレベル１の粒子は、その現在のボクセルからｘ方向に＋１および他の方向に０にあるボクセルに移動される。粒子は、最終的に、移動前と同じ状態（１，０，０）でその移動先ボクセルに到達する。ボクセル内の相互作用は、他の粒子および表面との局所的な相互作用に基づいて、この状態の粒子数を変更する可能性がある。そうでない場合、粒子は、同じ速度および方向で格子に沿って移動し続けることになる。

１つまたは複数の表面と相互作用するボクセルの場合、移動演算は、わずかに複雑になる。これにより、１つまたは複数の分画粒子がファセットに移送される可能性がある。このような分画粒子がファセットに移送されると、分画粒子がボクセルに残ることになる。これらの分画粒子は、ファセットが占有するボクセルに移送される。

図１３を参照すると、ボクセル３６２の状態ｉの粒子の一部３６０がファセット３６４に移動すると（２７８）、残りの部分３６６は、ファセット３６４が位置するボクセル３６８に移動し、状態ｉの粒子は、そのボクセル３６８から、ファセット３６４に方向付けられる。したがって、状態の母集団が２５に等しく、Ｖ_iα（ｘ）が０．２５に等しい（すなわち、ボクセルの４分の１が平行六面体Ｇ_iαと交差する）場合、６．２５個の粒子が、ファセットＦ_αに移動され、１８．７５個の粒子が、ファセットＦ_αによって占有されたボクセルに移動される。複数のファセットが単一のボクセルと交差する可能性があるので、１つまたは複数のファセットによって占有されたボクセルＮ（ｆ）に移送される状態ｉの粒子の数は、

であり、式中、Ｎ（ｘ）はソースボクセルである。

６．ファセットからボクセルへの散乱
次に、各ファセットからの流出粒子がボクセルに散乱される（２８６）。本質的に、この散乱は、粒子がボクセルからファセットに移動される収集の逆である。ファセットＦ_αからボクセルＮ（ｘ）に移動する状態ｉの粒子の数は、

であり、式中、Ｐ_f（ｘ）は、部分ボクセルの体積減少に相当する。この式から、状態ｉごとにファセットからボクセルＮ（ｘ）に方向付けられる粒子の総数は、

である。

ファセットからボクセルに粒子を散乱させ、それらを周囲のボクセルから移流した粒子と組み合わせ、結果を整数化した後、一定のボクセルにおける一定の方向が、アンダーフローする（負になる）か、またはオーバーフローする（８ビットの実装形態では２５５を超える）可能性がある。その結果、質量、運動量、およびエネルギーの量が値の許容範囲に収まるように切り捨てられて、質量、運動量、およびエネルギーが増加または減少することになる。このような事象から保護するために、問題のある状態を切り捨てる前に、範囲外の質量、運動量、およびエネルギーが蓄積される。状態が属するエネルギーの場合、（アンダーフローにより）増加した値、または（オーバーフローにより）減少した値に等しい量の質量が、同じエネルギーを有するとともにそれ自体がオーバーフローまたはアンダーフローの影響を受けないランダムに（または逐次に）選択された状態に追加される。質量およびエネルギーのこの追加から生じる追加の運動量は、蓄積され、切り捨てによる運動量に追加される。同じエネルギー状態にのみ質量を追加することによって、質量カウンタが０に達したときに質量とエネルギーの両方が補正される。最終的に、運動量は、運動量アキュムレータが０に戻るまで、プッシュ／プル技法を使用して補正される。

７．流体動力学の実行
流体動力学を実行する（２８８）、図７。このステップは、マイクロ動力学またはボクセル内演算と呼ばれることがある。同様に、移流手順は、ボクセル間演算と呼ばれることがある。以下に説明するマイクロ動力学演算を使用して、ファセットで粒子を衝突させ、ボルツマン分布を作り出すこともできる。

流体動力学は、ＢＧＫ衝突モデルとして知られる特定の衝突演算子によって格子ボルツマン方程式モデルにおいて保証される。この衝突モデルは、実際の流体系における分布の動力学を模倣する。衝突過程については、方程式１および方程式２の右辺によって適切に説明することができる。移流ステップの後、流体系、具体的には密度、運動量、およびエネルギーの保存量を、方程式３を使用して分布関数から取得する。これらの量から、方程式（２）のｆ^eqで示された平衡分布関数は、方程式（４）で完全に指定される。速度ベクトルの集合ｃ_i、重み、両方の選択肢は、表１にリスト化されており、方程式２とともに、巨視的挙動が正しい水力学方程式に従うことを保証する。

可変分解能
（米国特許出願公開第２０１３／０１５１２２１号に記載の）可変分解能を使用することもでき、様々なサイズのボクセル、例えば、粗いボクセルおよび細かいボクセルなど使用することになる。

本明細書に記載の主題および機能動作の実施形態は、デジタル電子回路、有形に具体化されたコンピュータソフトウェアもしくはファームウェア、コンピュータハードウェア（本明細書に開示された構造およびそれらの構造的均等物を含む）において、またはこれらのうちの１つまたは複数の組合せにおいて実装することができる。本明細書に記載の主題の実施形態は、１つまたは複数のコンピュータプログラム（すなわち、データ処理装置による実行のために、またはデータ処理装置の動作を制御するために、有形の非一過性プログラムキャリア上に符号化されたコンピュータプログラム命令の１つまたは複数のモジュール）として実装することができる。コンピュータ記憶媒体は、機械可読記憶デバイス、機械可読記憶基板、ランダムもしくはシリアルアクセスメモリデバイス、またはそれらのうちの１つまたは複数の組合せとすることができる。

「データ処理装置」という用語は、データ処理ハードウェアを指し、例として、プログラマブルプロセッサ、コンピュータ、または複数のプロセッサもしくはコンピュータを含む、データを処理するためのあらゆる種類の装置、デバイス、および機械を包含する。装置はまた、専用論理回路（例えば、ＦＰＧＡ（フィールドプログラマブルゲートアレイ）またはＡＳＩＣ（特定用途向け集積回路））とするか、またはそれをさらに含むことができる。ハードウェアに加えて、装置は、コンピュータプログラム用の実行環境を作り出すコード（例えば、プロセッサファームウェア、プロトコルスタック、データベース管理システム、オペレーティングシステム、またはそれらの１つまたは複数の組合せを構成するコード）を任意選択で含むことができる。

プログラム、ソフトウェア、ソフトウェアアプリケーション、モジュール、ソフトウェアモジュール、スクリプト、またはコードとも呼ばれることも記載されることもあるコンピュータプログラムは、コンパイル型もしくはインタープリタ型言語、または宣言型もしくは手続き型言語を含む、任意の形態のプログラミング言語で記述することができ、スタンドアロンプログラムとして、またはモジュール、コンポーネント、サブルーチン、もしくはコンピューティング環境での使用に適した別のユニットとしてなど、任意の形態で展開することができる。コンピュータプログラムは、ファイルシステム内のファイルに対応してもよいが、対応する必要はない。プログラムは、他のプログラムまたはデータ（例えば、マークアップ言語ドキュメント内、当該のプログラム専用の単一ファイル内、または複数の協調ファイル（例えば、１つまたは複数のモジュール、サブプログラム、もしくはコードの一部を格納するファイル）内に格納された１つまたは複数のスクリプト）を保持するファイルの一部に格納することができる。コンピュータプログラムは、プログラムが１つのコンピュータ上で実行されるように、または１つのサイトに配置されているか、もしくは複数のサイトにわたって分散されて、データ通信ネットワークによって相互接続されている複数のコンピュータ上で実行されるように展開することができる。

コンピュータプログラムの実行に適したコンピュータは、汎用もしくは専用マイクロプロセッサもしくはその両方、または任意の他の種類の中央処理ユニットに基づくことができる。コンピュータプログラム命令およびデータを記憶するのに適したコンピュータ可読媒体には、例として、半導体メモリデバイス（例えば、ＥＰＲＯＭ、ＥＥＰＲＯＭ、およびフラッシュメモリデバイス）、磁気ディスク（例えば、内部ハードディスクまたはリムーバブルディスク）、光磁気ディスク、ならびにＣＤ－ＲＯＭおよびＤＶＤ－ＲＯＭディスクを含む、媒体およびメモリデバイス上のすべての形態の不揮発性メモリを含む。プロセッサおよびメモリは、専用論理回路によって補完されるか、または専用論理回路を組み込まれてもよい。

本明細書に記載の主題の実施形態は、バックエンドコンポーネントを（例えば、データサーバとして）含むか、ミドルウェアコンポーネント（例えば、アプリケーションサーバ）を含むか、もしくはフロントエンドコンポーネント（例えば、ユーザが本明細書に記載の主題の実装と対話できるグラフィカルユーザインターフェースもしくはウェブブラウザを有するクライアントコンピュータ）を含むコンピューティングシステムで、または１つまたは複数のそのようなバックエンドコンポーネント、ミドルウェアコンポーネント、もしくはフロントエンドコンポーネントの任意の組合せで実装することができる。システムのコンポーネントは、任意の形態または媒体のデジタルデータ通信（例えば、通信ネットワーク）によって相互接続することができる。通信ネットワークの例には、ローカルエリアネットワーク（ＬＡＮ）およびワイドエリアネットワーク（ＷＡＮ）（例えば、インターネット）が含まれる。

コンピューティングシステムは、クライアントおよびサーバを含むことができる。クライアントおよびサーバは一般に、互いに離れており、典型的には通信ネットワークを通して対話する。クライアントとサーバの関係は、それぞれのコンピュータ上で実行するとともにクライアント－サーバ関係を有するコンピュータプログラムによって生じる。いくつかの実施形態において、サーバは、（例えば、データを表示し、ユーザデバイスと対話するユーザからのユーザ入力を受信する目的で）データ（例えば、ＨＴＭＬページ）を、クライアントとして機能するユーザデバイスに送信する。ユーザデバイスで生成されたデータ（例えば、ユーザ対話の結果）は、ユーザデバイスからサーバで受信することができる。

主題の特定の実施形態について説明してきた。他の実施形態も、以下の特許請求の範囲内にある。例えば、特許請求の範囲に記載の動作は、異なる順序で実行することができ、それでもなお望ましい結果を達成することができる。一例として、添付の図に示すプロセスは、望ましい結果を達成するために、図示の特定の順序または逐次的順序を必ずしも必要としない。場合によっては、マルチタスクおよび並列処理が有利である場合がある。

１０ システム
１１ バスシステム
１２ 大規模並列コンピューティングシステム／サーバシステム
１５ ネットワーク
１８ メモリ
２０ インターフェース
２４ 処理デバイス
３０ ユーザ提供のメッシュ定義
３２ メッシュ準備エンジン
３４ シミュレーションエンジン
３４ａ 衝突相互作用モジュール
３４ｂ 表面動力学変換
３４ｃ 境界モジュール
３４ｄ 移流粒子衝突相互作用モジュール
３８ データリポジトリ

Claims (25)

1. 流体流の要素をシミュレートするためのコンピュータ実施方法であって、
   複数のボクセルの状態ベクトルをメモリに記憶するステップであって、前記状態ベクトルが、ボクセルにおける複数の可能な運動量状態のうちの特定の運動量状態に対応する複数のエントリを含む、ステップと、
   前記ボクセルのサイズおよび向きとは無関係にサイズ指定され方向付けられた少なくとも１つの表面の表現をメモリに記憶するステップと、
   前記状態ベクトルに対して相互作用演算を実行するステップであって、前記相互作用演算が、異なる運動量状態の要素間の相互作用をモデル化する、ステップと、
   前記表面の前記表現に対して表面相互作用演算を実行するステップであって、前記表面相互作用演算が、前記表面とボクセルの実質的にすべての要素との間の相互作用をモデル化する、ステップと、
   要素の新しいボクセルへの移動を反映するために、前記状態ベクトルに対して移動演算を実行するステップと
   を含む、コンピュータ実施方法。
2. 前記表面相互作用演算が、
   ファセットと相互作用する少なくとも１つのボクセルの第１の集合から要素を収集することと、
   前記収集された要素とすべてのファセットとの間の相互作用をモデル化して表面相互作用要素の集合を作り出すことと、
   前記表面相互作用要素を、前記ファセットと相互作用する少なくとも１つのボクセルの第２の集合に散乱させることと
   を含む、請求項１に記載の方法。
3. 状態ベクトルのエントリが、ボクセルの特定の運動量状態における単位体積当たりの要素の密度を表す、請求項１に記載の方法。
4. 状態ベクトルが、特定の運動量状態における単位体積当たりの要素の密度を表すとともにより多くの可能な値を有する複数の整数および／または浮動小数点値のうちの１つまたは複数を含む、請求項１に記載の方法。
5. 前記表面相互作用演算が、整数値および／または浮動小数点値のうちの１つまたは複数を使用して実行され、前記表面相互作用演算が、実数を表す値を使用して実行される、請求項１に記載の方法。
6. 前記表面相互作用演算が、整数および／または浮動小数点数を使用して実行される、請求項１に記載の方法。
7. 表面と交差するボクセルを部分ボクセルとして表すステップをさらに含む、請求項１に記載の方法。
8. 少なくとも１つの表面の前記表現が、前記ボクセルのサイズおよび向きとは無関係にサイズ指定され方向付けられるとともに少なくとも１つの表面を表す複数のファセットを含み、
   前記表面相互作用演算が、ファセットと前記ファセットの近くの少なくとも１つのボクセルの前記要素との間の相互作用をモデル化することを含む、請求項１に記載の方法。
9. 前記要素が、流体の粒子を表し、前記ファセットが、前記流体が流れる少なくとも１つの表面を表す、請求項１に記載の方法。
10. 表面相互作用演算を実行するステップが、
    流入分布とボルツマン分布との間の差分を計算することと、
    すべての状態ベクトルの差分から合計運動量を決定することと、
    前記決定された差分に基づいて流出分布を生成することと
    をさらに含む、請求項１に記載の方法。
11. 流体流の要素をシミュレートするためのデータ処理システムであって、
    複数のボクセルの状態ベクトルをメモリに記憶することであって、前記状態ベクトルが、ボクセルにおける複数の可能な運動量状態のうちの特定の運動量状態に対応する複数のエントリを含む、記憶することと、
    前記ボクセルのサイズおよび向きとは無関係にサイズ指定され方向付けられた少なくとも１つの表面の表現をメモリに記憶することと、
    前記状態ベクトルに対して相互作用演算を実行することであって、前記相互作用演算が、異なる運動量状態の要素間の相互作用をモデル化する、相互作用演算を実行することと、
    前記表面の前記表現に対して表面相互作用演算を実行することであって、前記表面相互作用演算が、前記表面とボクセルの実質的にすべての要素との間の相互作用をモデル化する、表面相互作用演算を実行することと、
    要素の新しいボクセルへの移動を反映するために、前記状態ベクトルに対して移動演算を実行することと
    を前記データ処理システムに行わせるための命令を含む、データ処理システム。
12. 前記表面相互作用演算を実行するための前記命令が、
    ファセットと相互作用する少なくとも１つのボクセルの第１の集合から要素を収集することと、
    前記収集された要素とすべてのファセットとの間の相互作用をモデル化して表面相互作用要素の集合を作り出すことと、
    前記表面相互作用要素を、前記ファセットと相互作用する少なくとも１つのボクセルの第２の集合に散乱させることと
    を行うための命令を含む、請求項１１に記載のデータ処理システム。
13. 状態ベクトルのエントリが、ボクセルの特定の運動量状態における単位体積当たりの要素の密度を表す、請求項１１に記載のデータ処理システム。
14. 状態ベクトルが、特定の運動量状態における単位体積当たりの要素の密度を表すとともにより多くの可能な値を有する複数の整数および／または浮動小数点値のうちの１つまたは複数を含む、請求項１１に記載のデータ処理システム。
15. 前記表面相互作用演算が、整数値および／または浮動小数点値のうちの１つまたは複数を使用して実行され、前記表面相互作用演算が、実数を表す値を使用して実行される、請求項１１に記載のデータ処理システム。
16. 前記表面相互作用演算が、整数および／または浮動小数点数を使用して実行される、請求項１１に記載のデータ処理システム。
17. 前記表面相互作用演算を実行するための前記命令が、
    表面と交差するボクセルを部分ボクセルとして表すための命令を含む、請求項１１に記載のデータ処理システム。
18. 少なくとも１つの表面の前記表現が、ボクセルのサイズおよび向きとは無関係にサイズ指定され方向付けられるとともに少なくとも１つの表面を表す複数のファセットを含み、前記命令が、
    ファセットと前記ファセットの近くの少なくとも１つのボクセルの前記要素との間の相互作用をモデル化するための命令をさらに含む、請求項１１に記載のデータ処理システム。
19. 前記要素が、流体の粒子を表し、前記ファセットが、前記流体が流れる少なくとも１つの表面を表す、請求項１１に記載のデータ処理システム。
20. 前記表面相互作用演算を実行するための前記命令が、
    流入分布とボルツマン分布との間の差分を計算することと、
    すべての状態ベクトルの差分から合計運動量を決定することと、
    前記決定された差分に基づいて流出分布を生成することと
    を行うための命令を含む、請求項１１に記載のデータ処理システム。
21. 流体流の要素をシミュレートするためのコンピュータプログラム製品を記憶する有形のコンピュータ可読媒体であって、前記コンピュータプログラム製品が、
    複数のボクセルの状態ベクトルをメモリに記憶することであって、前記状態ベクトルが、ボクセルにおける複数の可能な運動量状態のうちの特定の運動量状態に対応する複数のエントリを含む、記憶することと、
    前記ボクセルのサイズおよび向きとは無関係にサイズ指定され方向付けられた少なくとも１つの表面の表現をメモリに記憶することと、
    前記状態ベクトルに対して相互作用演算を実行することであって、前記相互作用演算が、異なる運動量状態の要素間の相互作用をモデル化する、相互作用演算を実行することと、
    前記表面の前記表現に対して表面相互作用演算を実行することであって、前記表面相互作用演算が、前記表面とボクセルの実質的にすべての要素との間の相互作用をモデル化する、表面相互作用演算を実行することと、
    要素の新しいボクセルへの移動を反映するために、前記状態ベクトルに対して移動演算を実行することと
    をデータ処理システムに行わせるための命令を含む、コンピュータ可読媒体。
22. 前記表面相互作用命令が、
    ファセットと相互作用する少なくとも１つのボクセルの第１の集合から要素を収集することと、
    前記収集された要素とすべてのファセットとの間の相互作用をモデル化して表面相互作用要素の集合を作り出すことと、
    前記表面相互作用要素を、前記ファセットと相互作用する少なくとも１つのボクセルの第２の集合に散乱させることと
    を行うための命令をさらに含む、請求項２１に記載のコンピュータ可読媒体。
23. 状態ベクトルが、特定の運動量状態における単位体積当たりの要素の密度を表すとともにより多くの可能な値を有する複数の整数および／または浮動小数点値のうちの１つまたは複数を含む、請求項２１に記載のコンピュータ可読媒体。
24. 前記表面相互作用演算が、整数値および／または浮動小数点値のうちの１つまたは複数を使用して実行され、前記表面相互作用演算が、実数を表す値を使用して実行される、請求項２１に記載のコンピュータ可読媒体。
25. 表面相互作用演算を実行することが、
    流入分布とボルツマン分布との間の差分を計算することと、
    すべての状態ベクトルの差分から合計運動量を決定することと、
    前記決定された差分に基づいて流出分布を生成することと
    を行うための命令をさらに含む、請求項２１に記載のコンピュータ可読媒体。

Images