(57)【要約】【課題】多相流及び格子ボルツマン方式に関するＳｈａｎ−Ｃｈｅｎ法などの様々なタイプの数値シミュレーション手法を用いて実施する。【解決手段】本方法は、格子速度セットにおいて粒子間の衝突を引き起こす移送を用いて流体容積中の粒子移送をシミュレートする段階と、粒子移送に関する分布関数を生成する段階とを含み、分布関数は熱力学的ステップ及び粒子衝突ステップを含み、熱力学的ステップは粒子衝突ステップからは実質的に独立して切り離されている。【選択図】図１３
2019053752.tif 000002
【特許請求の範囲】
【請求項１】
格子速度セットにおいて流体容積中での粒子の移送をシミュレートする段階であって、前記移送が前記粒子間の衝突を引き起こす段階と、
前記粒子の移送に関する分布関数を生成する段階と、
を含む方法であって、
前記分布関数は、熱力学的ステップ及び粒子衝突ステップを含み、前記熱力学的ステップは、前記粒子衝突ステップからは実質的に独立して切り離されている、方法。
【請求項２】
前記分布関数は、移流ステップをさらに含み、前記熱力学的ステップは、前記熱力学的ステップで前記粒子衝突ステップを強化するのではなく、前記熱力学的ステップで前記移流ステップを強化することによって、前記分布関数の一部に含まれる、請求項１に記載の方法。
【請求項３】
前記熱力学的ステップは、前記移送中の前記流体容積の温度を含む、請求項１に記載の方法。
【請求項４】
前記分布関数を生成する段階は、
特定時間ｔにおける前記流体容積中の特定位置ｘでの前記衝突に対する衝突後分布関数ｆi′（ｘ，ｔ）を決定する段階であって、ここで、ｆi′（ｘ，ｔ）＝ｆi（ｘ，ｔ）＋  Ｃi（ｘ，ｔ）であり、Ｃiは衝突演算子、ｆiは前記衝突前の前記粒子に関する分布関数である、段階と、
分画部分ｇi（ｘ，ｔ）を前記衝突後分布関数ｆi′（ｘ，ｔ）から差し引いて移流前粒子密度分布関数ｆi′′（ｘ，ｔ）＝  ｆi′（ｘ，  ｔ）  −ｇi（ｘ，  ｔ）を得る段階であって、前記粒子の一部ｆi′（ｘ，ｔ）−ｇi（ｘ，ｔ）は流体容積中の別の位置に移流し、ｇi（ｘ，ｔ）は移流しない粒子の分布関数を表す、段階と、
時間ｔ＋Δｔにおいて前記粒子の一部が前記流体容積中の他の位置への移流をシミュレートする段階であって、前記他の位置は（ｘ＋ｃiΔｔ）で表され、ｃiは衝突前の前記粒子の速度ベクトルであり、Δｔは特定の時間ｔと別の時点との間の間隔である、段階と、  前記移流シミュレーションに基づいて、前記移流粒子の密度分布関数
2019053752.tif 000081
を得る段階であって、
2019053752.tif 000082
であり、
2019053752.tif 000083
は位置ｘ−ｃiから移流した位置ｘにおける粒子の分布関数である、段階と、
前記先に差し引いた部分ｇi（ｘ，ｔ）を前記密度分布関数に戻して加えて移流後密度分布関数
2019053752.tif 000084
を生成する段階と、
時間ｔ＋Δｔにおける位置ｘでの前記粒子の質量、運動量、及び温度を計算する段階と、
前記計算された温度、質量、及び運動量を用いてｇi（ｘ，ｔ＋Δｔ）を決定する段階と、
差分ｇi（ｘ，ｔ＋Δｔ）−ｇi（ｘ，ｔ）を運動状態ｆi（ｘ，ｔ）＋Ｃi（ｘ，ｔ）に加える段階と、
を含む、請求項１に記載の方法。
【請求項５】
ｇiは、
2019053752.tif 000085
に従って定義され、
ρは流体密度であり、
Ｔ0は一定格子温度であり、
Ｐは前記流体容積中の圧力であり、
Ｔは前記計算された温度であり、
ｗiは一定の重み因子である、
請求項４に記載の方法。
【請求項６】
前記生成される分布関数は、
2019053752.tif 000086
であり、
ｘは前記容積中の特定位置であり、
ｔは特定の第１の時点であり、
ｉは前記集合での格子速度のインデックス数であり、
ｃiは衝突前の前記粒子の速度ベクトルであり、
Ｃiは衝突演算子であり、
Δｔは前記第１の時点と第２の時点との間の間隔であり、
ｇiは前記熱力学的ステップであり、
ｆiは時間ｔにおける位置ｘでの前記粒子に対する分布関数である、
請求項４に記載の方法。
【請求項７】
静止状態を
2019053752.tif 000087
に一致するように修正することによって、前記質量を保存する段階をさらに含む、請求項４に記載の方法。
【請求項８】
前記粒子の衝突ステップは、等温平衡分布関数を含む、請求項１に記載の方法。
【請求項９】
前記分布関数は、時間ｔ＋Δｔにおける位置ｘ＋ｃiΔｔでの前記粒子に対する分布関数であり、
2019053752.tif 000088
で表され、
ｘは前記容積中の特定位置であり、
ｔは特定の第１の時点であり、
ｉは前記セットでの格子速度のインデックス数であり、
Ｃiは衝突演算子であり、
ｃiは衝突前の前記粒子の速度ベクトルであり、
Δｔは前記第１の時点と第２の時点との間の間隔であり、
ｇiは前記熱力学的ステップであり、
ｆiは時間ｔにおける位置ｘでの前記粒子に対する分布関数である、
請求項１に記載の方法。
【請求項１０】
ｇiは、
2019053752.tif 000089
に従って定義され、
ρは流体密度であり、
Ｔ0は一定格子温度であり、
Ｐは前記流体容積中の圧力であり、
Ｔは前記流体の実際の温度であり、
ｗiは一定の重み因子である、
請求項９に記載の方法。
【請求項１１】
前記格子速度セットは、格子ボルツマン法に基づく、請求項１に記載の方法。
【請求項１２】
１又は２以上の処理装置で実行可能な命令を記憶する１又は２以上の機械可読ハードウェア記憶装置であって、
格子速度セットにおいて流体容積中での粒子の移送をシミュレートする段階であって、前記移送が前記粒子間の衝突を引き起こす段階と、
前記粒子の移送に関する分布関数を生成する段階と、
を含む方法であって、
前記分布関数は、熱力学的ステップ及び粒子衝突ステップを含み、前記熱力学的ステップは、前記粒子衝突ステップからは実質的に独立して切り離されている、方法を含む、１又は２以上の機械可読ハードウェア記憶装置。
【請求項１３】
前記分布関数は、移流ステップをさらに含み、前記熱力学的ステップは、前記熱力学的ステップで前記粒子衝突ステップを強化するのではなく、前記熱力学的ステップで前記移流ステップを強化することによって、前記分布関数の一部に含まれる、請求項１２に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項１４】
前記熱力学的ステップは、前記移送中の前記流体容積の温度を含む、請求項１２に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項１５】
前記分布関数を生成する段階は、
特定時間ｔにおける前記流体容積中の特定位置ｘでの前記衝突に対する衝突後分布関数ｆi′（ｘ，ｔ）を決定する段階であって、ここで、ｆi′（ｘ，ｔ）＝ｆi（ｘ，ｔ）＋  Ｃi（ｘ，ｔ）であり、Ｃiは衝突演算子、ｆiは前記衝突前の前記粒子に関する分布関数である、段階と、
分画部分ｇi（ｘ，ｔ）を前記衝突後分布関数ｆi′（ｘ，ｔ）から差し引いて移流前粒子密度分布関数ｆi′′（ｘ，ｔ）＝  ｆi′（ｘ，  ｔ）  −ｇi（ｘ，  ｔ）を得る段階であって、前記粒子の一部ｆi′（ｘ，ｔ）−ｇi（ｘ，ｔ）は流体容積中の別の位置に移流し、ｇi（ｘ，ｔ）は移流しない粒子の分布関数を表す、段階と、
時間ｔ＋Δｔにおいて前記粒子の一部が前記流体容積中の他の位置への移流をシミュレートする段階であって、前記他の位置は（ｘ＋ｃiΔｔ）で表され、ｃiは衝突前の前記粒子の速度ベクトルであり、Δｔは特定の時間ｔと別の時点との間の間隔である、段階と、  前記移流シミュレーションに基づいて、前記移流粒子の密度分布関数
2019053752.tif 000090
を得る段階であって、
2019053752.tif 000091
であり、
2019053752.tif 000092
は位置ｘ−ｃiから移流した位置ｘにおける粒子の分布関数である、段階と、
前記先に差し引いた部分ｇi（ｘ，ｔ）を前記密度分布関数に戻して加えて移流後密度分布関数
2019053752.tif 000093
を生成する段階と、
時間ｔ＋Δｔにおける位置ｘでの前記粒子の質量、運動量、及び温度を計算する段階と、
前記計算された温度、質量、及び運動量を用いてｇi（ｘ，ｔ＋Δｔ）を決定する段階と、
差分ｇi（ｘ，ｔ＋Δｔ）−ｇi（ｘ，ｔ）を運動状態ｆi（ｘ，ｔ）＋Ｃi（ｘ，ｔ）に加える段階と、
を含む、請求項１２に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項１６】
ｇiは、
2019053752.tif 000094
に従って定義され、
ρは流体密度であり、
Ｔ0は一定格子温度であり、
Ｐは前記流体容積中の圧力であり、
Ｔは前記計算された温度であり、
ｗiは一定の重み因子である、
請求項１５に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項１７】
前記生成される分布関数は、
2019053752.tif 000095
であり、
ｘは前記容積中の特定位置であり、
ｔは特定の第１の時点であり、
ｉは前記集合での格子速度のインデックス数であり、
ｃiは衝突前の前記粒子の速度ベクトルであり、
Ｃiは衝突演算子であり、
Δｔは前記第１の時点と第２の時点との間の間隔であり、
ｇiは前記熱力学的ステップであり、
ｆiは時間ｔにおける位置ｘでの前記粒子に対する分布関数である、
請求項１５に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項１８】
静止状態を
2019053752.tif 000096
に一致するように修正することによって、前記質量を保存する段階をさらに含む、請求項１５に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項１９】
前記粒子の衝突ステップは、等温平衡分布関数を含む、請求項１２に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項２０】
前記分布関数は、時間ｔ＋Δｔにおける位置ｘ＋ｃiΔｔでの前記粒子に対する分布関数であり、
2019053752.tif 000097
で表され、
ｘは前記容積中の特定位置であり、
ｔは特定の第１の時点であり、
ｉは前記セットでの格子速度のインデックス数であり、
Ｃiは衝突演算子であり、
ｃiは衝突前の前記粒子の速度ベクトルであり、
Δｔは前記第１の時点と第２の時点との間の間隔であり、
ｇiは前記熱力学的ステップであり、
ｆiは時間ｔにおける位置ｘでの前記粒子に対する分布関数である、
請求項１２に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項２１】
ｇiは、
2019053752.tif 000098
に従って定義され、
ρは流体密度であり、
Ｔ0は一定格子温度であり、
Ｐは前記流体容積中の圧力であり、
Ｔは前記流体の実際の温度であり、
ｗiは一定の重み因子である、
請求項２０に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項２２】
前記格子速度セットは、格子ボルツマン法に基づく、請求項１２に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項２３】
１又は２以上の処理装置と、
前記１又は２以上の処理装置で実行可能な命令を記憶する１又は２以上の機械可読ハードウェア記憶装置と、
を含むシステムであって、前記１又は２以上の機械可読ハードウェア記憶装置は、
格子速度セットにおいて流体容積中での粒子の移送をシミュレートする段階であって、前記移送が前記粒子間の衝突を引き起こす段階と、
前記粒子の移送に関する分布関数を生成する段階と、
を含む方法であって、
前記分布関数は、熱力学的ステップ及び粒子衝突ステップを含み、前記熱力学的ステップは、前記粒子衝突ステップからは実質的に独立して切り離されている方法を含む、システム。
【請求項２４】
前記分布関数は、移流ステップをさらに含み、前記熱力学的ステップは、前記熱力学的ステップで前記粒子衝突ステップを強化するのではなく、前記熱力学的ステップで前記移流ステップを強化することによって、前記分布関数の一部に含まれる、請求項２３に記載のシステム。
【請求項２５】
前記熱力学的ステップは、前記移送中の前記流体容積の温度を含む、請求項２３に記載のシステム。
【請求項２６】
前記分布関数を生成する段階は、
特定時間ｔにおける前記流体容積中の特定位置ｘでの前記衝突に対する衝突後分布関数ｆi′（ｘ，ｔ）を決定する段階であって、ここで、ｆi′（ｘ，ｔ）＝ｆi（ｘ，ｔ）＋  Ｃi（ｘ，ｔ）であり、Ｃiは衝突演算子、ｆiは前記衝突前の前記粒子に関する分布関数である、段階と、
分画部分ｇi（ｘ，ｔ）を前記衝突後分布関数ｆi′（ｘ，ｔ）から差し引いて移流前粒子密度分布関数ｆi′′（ｘ，ｔ）＝  ｆi′（ｘ，  ｔ）  −ｇi（ｘ，  ｔ）を得る段階であって、前記粒子の一部ｆi′（ｘ，ｔ）−ｇi（ｘ，ｔ）は流体容積中の別の位置に移流し、ｇi（ｘ，ｔ）は移流しない粒子の分布関数を表す、段階と、
時間ｔ＋Δｔにおいて前記粒子の一部が前記流体容積中の他の位置への移流をシミュレートする段階であって、前記他の位置は（ｘ＋ｃiΔｔ）で表され、ｃiは衝突前の前記粒子の速度ベクトルであり、Δｔは特定の時間ｔと別の時点との間の間隔である、段階と、  前記移流シミュレーションに基づいて、前記移流粒子の密度分布関数
2019053752.tif 000099
を得る段階であって、
2019053752.tif 000100
であり、
2019053752.tif 000101
は位置ｘ−ｃiから移流した位置ｘにおける粒子の分布関数である、段階と、
前記先に差し引いた部分ｇi（ｘ，ｔ）を前記密度分布関数に戻して加えて移流後密度分布関数
2019053752.tif 000102
を生成する段階と、
時間ｔ＋Δｔにおける位置ｘでの前記粒子の質量、運動量、及び温度を計算する段階と、
前記計算された温度、質量、及び運動量を用いてｇi（ｘ，ｔ＋Δｔ）を決定する段階と、
差分ｇi（ｘ，ｔ＋Δｔ）−ｇi（ｘ，ｔ）を運動状態ｆi（ｘ，ｔ）＋Ｃi（ｘ，ｔ）に加える段階と、
を含む、請求項２３に記載のシステム。
【請求項２７】
ｇiは、
2019053752.tif 000103
に従って定義され、
ρは流体密度であり、
Ｔ0は一定格子温度であり、
Ｐは前記流体容積中の圧力であり、
Ｔは前記計算された温度であり、
ｗiは一定の重み因子である、
請求項２６に記載のシステム。
【請求項２８】
前記生成される分布関数は、
2019053752.tif 000104
であり、
ｘは前記容積中の特定位置であり、
ｔは特定の第１の時点であり、
ｉは前記集合での格子速度のインデックス数であり、
ｃiは衝突前の前記粒子の速度ベクトルであり、
Ｃiは衝突演算子であり、
Δｔは前記第１の時点と第２の時点との間の間隔であり、
ｇiは前記熱力学的ステップであり、
ｆiは時間ｔにおける位置ｘでの前記粒子に対する分布関数である、
請求項２６に記載のシステム。
【請求項２９】
前記操作は、静止状態を
2019053752.tif 000105
に一致するように修正することによって、前記質量を保存する段階をさらに含む、請求項２６に記載のシステム。
【請求項３０】
前記粒子の衝突ステップは、等温平衡分布関数を含む、請求項２３に記載のシステム。
【請求項３１】
前記分布関数は、時間ｔ＋Δｔにおける位置ｘ＋ｃiΔｔでの前記粒子に対する分布関数であり、
2019053752.tif 000106
で表され、
ｘは前記容積中の特定位置であり、
ｔは特定の第１の時点であり、
ｉは前記セットでの格子速度のインデックス数であり、
Ｃiは衝突演算子であり、
ｃiは衝突前の前記粒子の速度ベクトルであり、
Δｔは前記第１の時点と第２の時点との間の間隔であり、
ｇiは前記熱力学的ステップであり、
ｆiは時間ｔにおける位置ｘでの前記粒子に対する分布関数である、
請求項２３に記載のシステム。
【請求項３２】
ｇiは、
2019053752.tif 000107
に従って定義され、
ρは流体密度であり、
Ｔ0は一定格子温度であり、
Ｐは前記流体容積中の圧力であり、
Ｔは前記流体の実際の温度であり、
ｗiは一定の重み因子である、
請求項２３に記載のシステム。
【請求項３３】
前記格子速度セットは、格子ボルツマン法に基づく、請求項３２に記載のシステム。
【発明の詳細な説明】
【技術分野】
【０００１】
（優先権の主張）
本出願は、２０１３年７月３１日出願の米国仮特許出願第６１／８６０，３９２号に対して米国特許法第１１９条（ｅ）の下で優先権を主張するものであり、その開示内容全体は引用により本明細書に組み込まれている。
【背景技術】
【０００２】
格子ボルツマン法（ＬＢＭ）（又は熱格子ボルツマン法（ＴＬＢＭ））は、流体シミュレーションに関する計算流体力学（ＣＦＤ）の１つの態様である。ナビエ−ストークス方程式を解く代わりに、離散ボルツマン方程式を解いてＢｈａｔｎａｇａｒ−Ｇｒｏｓｓ−Ｋｒｏｏｋ（ＢＧＫ）などの衝突モデルを用いてニュートン流体の流れをシミュレートする。限定数の粒子にわたって流れ及び衝突過程をシミュレートすることにより、本質的な粒子相互作用は、より大きな集団にわたって適用可能な粘性流挙動の小世界を明示する。
【先行技術文献】
【特許文献】
【０００３】
【特許文献１】米国特許第５５９４６７１号明細書
【発明の概要】
【課題を解決するための手段】
【０００４】
一般に、本明細書は、格子速度セットにおいて、粒子間の衝突を引き起こす移送を用いて流体容積中の粒子移送をシミュレートして、粒子移送に関する分布関数を生成するための技術を記載し、分布関数は、熱力学的ステップ及び粒子衝突ステップを含み、熱力学的ステップは、粒子衝突ステップから実質的に独立して切り離されている。
【０００５】
いくつかの実施例では、分布関数は、移流ステップをさらに含み、熱力学的ステップは、熱力学的ステップを用いて粒子衝突ステップを強化するのではなく、熱力学的ステップを用いて移流ステップを強化することによって、分布の一部に含まれる。別の実施例では、熱力学的ステップは、移送中の流体容積の温度を含む。さらに別の実施例では、生成することは、特定時間ｔにおける流体容積中の特定位置ｘでの衝突に関する衝突後分布関数ｆi′（ｘ，ｔ）を決定するステップであって、
2019053752.tif 000003
であり、ここで、Ｃiは衝突演算子、ｆiは衝突前の粒子に関する分布関数である段階と、分画部分ｇi（ｘ，ｔ）を衝突後分布関数ｆi′（ｘ，ｔ）から差し引いて移流前の粒子密度分布関数
2019053752.tif 000004
を得る段階であって、粒子の一部ｆi′（ｘ，ｔ）−ｇi（ｘ，ｔ）が流体容積中の別の位置に移流し、ｇi（ｘ，ｔ）は移流しない粒子の分布関数を表す、段階と、時間ｔ＋Δｔにおいて粒子の一部の流体容積中の他の位置への移流をシミュレートする段階であって、他の位置が（ｘ＋ｃiΔｔ）として表され、ｃiが衝突前の粒子の速度ベクトルであり、Δｔが特定の時間ｔと別の時点との間の間隔である、段階と、移流シミュレーションに基づいて、移流する粒子の密度分布関数
2019053752.tif 000005
を得る段階であって、
2019053752.tif 000006
であり、
2019053752.tif 000007
は、位置ｘ−ｃiΔｔから移流した位置ｘにおける粒子の分布関数である、段階と、先に差し引いた部分ｇi（ｘ，ｔ）を密度分布関数へ戻して加えて移流後密度分布関数
2019053752.tif 000008
を生成する段階と、時間ｔ＋Δｔにおける位置ｘでの粒子の質量、運動量、及び温度を計算する段階と、計算された温度、質量、及び運動量を用いてｇi（ｘ，ｔ＋Δｔ）を決定する段階と、差分ｇi（ｘ，ｔ＋Δｔ）−ｇi（ｘ，ｔ）を運動状態ｆi（ｘ，ｔ）＋Ｃi（ｘ，ｔ）に加える段階と、を含む。
【０００６】
いくつかの実施例では、ｇiは以下のように定義され、
2019053752.tif 000009
ここで、ρは流体密度、Ｔ0は一定格子温度、Ｐは流体容積中の圧力、Ｔは計算された温度、及びｗiは一定重み因子である。生成された分布関数は以下の通りであり、
2019053752.tif 000010
ここで、ｘは容積中の特定位置、ｔは特定の第１の時点、ｉは格子速度セットのインデックス数、ｃiは衝突前の粒子の速度ベクトル、Ｃiは衝突演算子、Δｔは第１の時点と第２の時点との間の間隔、ｇiは熱力学的ステップ、及びｆiは時間ｔにおける位置ｘでの粒子に対する分布関数である。
【０００７】
いくつかの実施例では、本方法は、静止状態を以下の式に従って変更することによって質量を保存することを含む。
2019053752.tif 000011
別の実施例では、粒子の衝突ステップは、等温平衡分布関数を含む。他の実施例では、分布関数は、時間ｔ＋Δｔにおける位置ｘ＋ｃiΔｔでの分布関数であり、以下の通りであり、
2019053752.tif 000012
ここで、ｘは容積中の特定位置、ｔは特定の第１の時点、ｉは格子速度のセットのインデックス数、Ｃiは衝突演算子、ｃiは衝突前の粒子の速度ベクトルであり、Δｔは、第１の時点と第２の時点との間の間隔、ｇiは熱力学的ステップ、ｆiは時間ｔにおける位置ｘでの粒子に対する分布関数である。さらに別の実施例では、ｇiは以下のように定義され、
2019053752.tif 000013
ここで、ρは流体密度、Ｔ0は一定格子温度、Ｐは流体容積中の圧力、Ｔは計算された温度、及びｗiは一定の重み因子である。いくつかの実施例では、格子速度セットは、格子ボルツマン法に基づく。
【０００８】
前述の全て又は一部は、１又は２以上の不揮発性機械可読記憶媒体に記憶されかつ１又は２以上の処理装置で実行可能な命令を含むコンピュータプログラムとして実装することができる。前述の全て又は一部は、装置、方法、又は電子システムとして実装して前述の機能を実行することができ、電子システムは、１又は２以上の処理装置及び実行可能な命令を記憶するためのメモリを含むことができる。
【０００９】
１又は２以上の実施構成の詳細は、添付図面及び以下の説明で示される。他の特徴、目的、及び利点は、図面及び詳細な説明、並びに特許請求の範囲から明らになるはずである。
【００１０】
本システム、方式、及び技術は、多相流及び格子ボルツマン方式に関するＳｈａｎ−Ｃｈｅｎ法などの様々なタイプの数値シミュレーション手法を用いて実施することができる。格子ボルツマン方式に関する追加の情報を本明細書で以下に説明する。しかしながら、本明細書に説明するシステム及び技術は、格子ボルツマン方式を用いるシミュレーションに限定されず、他の数値シミュレーション手法に適用することができる。
【００１１】
システム及び技術は、格子ボルツマン方式を使用する格子ガスシミュレーションを使用して実施することができる。従来の格子ガスシミュレーションは、各格子部位での限られた数の粒子を仮定し、粒子は、ビットのショートベクトルによって表される。各ビットは、特定の方向に移動する粒子を表している。例えば、ベクトル中の１ビットは、特定の方向に沿って移動する粒子の存在（１に設定時）又は不在（０に設定時）を表すことができる。そのようなベクトルは、６ビットを有することができ、例えば、値１１００００は、２つの粒子がＸ軸に沿って反対方向に移動し、かつＹ軸及びＺ軸に沿って移動する粒子がないことを示している。１組の衝突規則が、各部位での粒子間の衝突の挙動を支配する（例えば、１１００００ベクトルは、００１１００ベクトルになることができ、Ｘ軸に沿って移動する２つの粒子間の衝突は、Ｙ軸に沿って離れる２つの粒子を生成したことを示す）。規則は、ビットに対して置換を実行する（例えば、１１００００を００１１００に変形する）ルックアップテーブルに状態ベクトルを供給することによって実施される。粒子は、その後に、隣接部位に移動される（例えば、Ｙ軸に沿って移動する２つの粒子は、Ｙ軸に沿って左右に隣接部位に移動されると考えられる）。
【００１２】
強化されたシステムにおいて、各格子部位での状態ベクトルは、粒子エネルギ及び移動方向の変動を与えるためにより多くのビットを含み（例えば、亜音速流れに対して５４ビット）、完全な状態ベクトルの部分集合を伴う衝突規則が使用される。更に強化されたシステムにおいて、１つよりも多い粒子が、各格子部位又はボクセル（これらの２つの用語は、本明細書を通して交換可能に使用される）で各運動量状態で存在することが許される。例えば、８ビット実装において、０〜２５５個の粒子は、特定のボクセルで特定の方向に移動している可能性がある。状態ベクトルは、１組のビットである代わりに、１組の整数であり（例えば、１組の８ビットバイトが、０〜２５５の範囲の整数を与える）、その各々は、与えられた状態にある粒子の数を表する。
【００１３】
更なる強化において、格子ボルツマン方法（ＬＢＭ）は、従来の計算流体力学（「ＣＦＤ」）手法で可能であるよりも深いレベルで複合的な幾何学形状における３Ｄ不安定圧縮性乱流過程をシミュレートするために流体のメゾスコピック表現を使用する。ＬＢＭ方法の概要を以下に与える。
【００１４】
ボルツマンレベルのメゾスコピック表現
流体システムは、いわゆる「メゾスコピック」レベルに対する運動方程式によって表すことができることが統計物理学で公知である。このレベルに対しては、個々の粒子の詳細な運動は決定する必要はない。これに代えて、流体の特性は、単一粒子位相空間、
2019053752.tif 000014
を使用して定義される粒子分布関数によって表され、ｘは、空間座標であり、νは、粒子速度座標である。質量、密度、流量速度、及び温度のような典型的な流体力学的な量は、粒子分布関数の単純なモーメントである。粒子分布関数の力学は、ボルツマン方程式に従う。
2019053752.tif 000015
方程式（１）
ここで、Ｆ（ｘ，ｔ）は、（ｘ，ｔ）で外部又は自己矛盾なく発生された物体力を表している。衝突項Ｃは、様々な速度及び位置の粒子の相互作用を表している。衝突項Ｃに対して特定の形態を指定することなく、上述のボルツマン方程式は、希薄ガスの公知の状況（ボルツマンによって本来構成されたような）だけでなく全ての流体システムに適用可能であることを強調することが重要である。
【００１５】
一般的に、Ｃは、二点相関関数の複雑な多次元積分を含む。分布関数ｆのみで閉じたシステムを形成するために、並びに効率的な計算のために、最も便利かつ物理的に一貫した形態の１つは、公知のＢＧＫ演算子である。ＢＧＫ演算子は、衝突の詳細がたとえ何であろうとも、分布関数は、衝突：
2019053752.tif 000016
方程式（２）
を通じて
2019053752.tif 000017
によって与えられる明確な局所平衡に近づくという物理的議論に従って構成され、ここで、パラメータτは、衝突を通じた平衡までの固有緩和時間を表している。粒子（例えば、原子又は分子）を扱うと、緩和時間は、典型的には定数として取られる。「混成」（水力−動力学）表現において、この緩和時間は、歪み率及び乱流運動エネルギなどのような流体力学的変数の関数である。すなわち、乱流は、局所的に決定された固有特性を有する乱流粒子（「渦」）のガスとして表すことができる。
【００１６】
ボルツマン−ＢＧＫ方程式の数値解は、ナビエ−ストークス方程式の解に優るいくつかの計算上の利点を有する。第１に、複雑な非線形項又は高次空間導関数が方程式にないことを直ちに認識することができ、すなわち、移流不安定性に関する問題がほとんどない。このレベルの記述では、方程式は、圧力を扱う必要がないので局所的であり、これは、アルゴリズム並列化に対してかなりの利点を提供する。２次空間導関数を有する拡散演算子がないという事実と共に線形移流演算子の別の望ましい特徴は、流体偏微分方程式（「ＰＤＥ」）の数学的な条件ではなく、実際に粒子が固体面と真に相互作用する方法を模倣するように非スリップ面又はスリップ面のような物理的境界条件を実現する際の容易さである。直接の利益の１つは、固体面上でインタフェースの移動を処理する問題がないということであり、これは、格子ボルツマンベースのシミュレーションソフトウエアが失敗なく複合乱流空気力学をシミュレートすることを可能にするのを補助する。これに加えて、有限粗度面のような境界からのある一定の物理特性も、力に組み込むことができる。更に、ＢＧＫ衝突演算子は、純粋に局所的であり、一方、自己矛盾のない物体力の計算は、近−隣接情報のみを通して達成することができる。その結果、ボルツマン−ＢＧＫ方程式の計算は、並列処理に実質的に適応させることができる。
【００１７】
格子ボルツマン方式
連続体ボルツマン方程式を解くことは、それが位置及び速度位相空間における積分微分方程式の数値評価を伴うということで有意な課題を表している。位置だけでなく速度位相空間も離散化することができることが観測された時に大規模な簡素化が行われ、それは、ボルツマン方程式の解のための効率的な数値アルゴリズムをもたらした。流体力学的な量は、高々最近隣接情報に依存する単純和の項で書くことができる。従来的には、格子ボルツマン方程式の方式は、速度の離散集合：
2019053752.tif 000018
に対する粒子の推移を規定する格子ガスモデルに基づいていたが、この方程式は、連続体ボルツマン方程式の離散化としての第１の原理から系統的に導出することができる。その結果、ＬＢＥは、格子ガス手法に関連付けられた公知の問題を被らない。従って、位相空間における連続体分布関数ｆ（ｘ，ν，ｔ）を取り扱う代わりに、離散速度指数をラベル付けする下付き文字を有する離散分布の有限集合
2019053752.tif 000019
を追跡することが必要なだけである。巨視的記述の代わりにこの動態方程式を取り扱う重要な利点は、システムの位相空間の増大が、問題の局所性によってオフセットされるということである。
【００１８】
対称性の考慮に起因して、速度値の集合は、それらが構成空間内で張られた時にある一定格子構造を形成するように選択される。そのような離散システムの力学は、形態：
2019053752.tif 000020
を有するＬＢＥに従い、ここで、衝突演算子は、上述のように通常はＢＧＫ形態を取る。平衡分布形態の適正な選択により、格子ボルツマン方程式は正しい流体力学及び熱−流体力学を生じさせることを理論的に示すことができる。すなわち、
2019053752.tif 000021
から導出された流体力学モーメントは、巨視的限界においてナビエ−ストークス方程式に従う。これらのモーメントは、
2019053752.tif 000022
方程式（３）
して定義され、ここで、ρ、ｕ、及びＴは、それぞれ、流体密度、速度、及び温度であり、Ｄは、離散化された速度空間の次元（物理空間次元に全く等しくない）である。
【００１９】
他の特徴及び利点は、図面及び特許請求の範囲を含めて、以下の説明から明らかとなる。
【図面の簡単な説明】
【００２０】
【図１】ＬＢＭモデルの速度成分を例示する。
【図２】ＬＢＭモデルの速度成分を例示する。
【図３】物理過程シミュレーションシステムが従う手順の流れ図である。
【図４】マイクロブロックの斜視図である
【図５Ａ−Ｂ】図３のシステムで使用する格子構造を例示する。
【図６】可変分解能技術を例示する。
【図７】可変分解能技術を例示する。
【図８】表面のファセットにより影響を受ける領域を例示する。
【図９】ボクセルから表面への粒子の移動を例示する。
【図１０】表面から表面への粒子の移動を例示する。
【図１１】表面動力学を実行するための手順の流れ図である。
【図１２】粒子衝突ステップから独立した熱力学ステップを用いて分布関数を生成する過程の流れ図である。
【図１３】粒子衝突ステップから独立した熱力学ステップを用いて分布関数を生成するシステムの構成要素のブロック図である。
【発明を実施するための形態】
【００２１】
Ａ．概要
本開示に適合するシステムは、温度結合を衝突ステップから取り除いて、温度を組み込んだ別個の熱力学的ステップとして温度結合をＬＢＭの中に導入することによって、温度をＬＢＭに結合する。熱力学的ステップ（ＬＢＭを用いて生成される）を分布関数に導入することによって、この分布関数は、高マッハ数及び高い温度範囲用途のシミュレーションを行うために用いることができる（温度範囲の制限を課する１次項を修正することによって温度をＬＢＭに結合しないという理由で）。
【００２２】
Ｂ．シミュレーション空間のモデル化
ＬＢＭベースの物理過程シミュレーションシステムでは、離散速度ｃiの集合で評価される分布関数値ｆiによって、流体の流れを表現することができる。分布関数の力学は、方程式４により支配され、ここにｆi（０）は平衡分布関数として知られ、
2019053752.tif 000023
方程式（４）
として定義され、ここで、
2019053752.tif 000024
である。
2019053752.tif 000025
方程式（５）
この方程式は、分布関数ｆiの時間発展を説明する、よく知られた格子ボルツマン方程式である。左辺はいわゆる「流動過程」による分布の変化を表す。流動過程は、流体ポケットが或るグリッド位置から出発して、それから速度ベクトルの１つに沿って次のグリッド位置へ移動する過程である。その地点で、「衝突演算子」つまり近接した流体ポケットの開始流体ポケットへの影響を計算する。流体は別のグリッド位置に移動することだけが可能なので、全速度の成分全体は共通した速度の倍数であるように速度ベクトルの適切な選択が必要である。
【００２３】
第１の方程式の右辺は、流体ポケット間の衝突による分布関数の変化を表す前述の「衝突演算子」である。本明細書で使用する衝突演算子の特別な形式は、Ｂｈａｔｎａｇａｒ、Ｇｒｏｓｓ及びＫｒｏｏｋ（ＢＧＫ）によるものである。衝突演算子は、分布関数を、「平衡」形である第２の式により与えられる規定値へ向かうように強いる。
【００２４】
このシミュレーションから、密度ρ及び流速ｕなどの従来の流体変数は、式（３）の単純和として得ることができる。ここで、ｃi及びｗiの集団的な値はＬＢＭモデルを規定する。ＬＢＭモデルは、大規模に実現可能なコンピュータ・プラットフォームで効率的に実行することができ、時間的に不安定な流れ及び複雑な境界条件に対して高い堅牢性で走らせることができる。
【００２５】
ボルツマン方程式から流体システムに対する巨視的な運動方程式を得るための標準的な技術は、Ｃｈａｐｍａｎ−Ｅｎｓｋｏｇ法であり、そこでは完全なボルツマン方程式の逐次近似がなされる。
【００２６】
流体システムでは、密度の小さな擾乱が音速で伝播する。気体系では、音速は一般に温度で決まる。流れにおける圧縮性の影響の重要性は、特徴的な速度と音速との比によって評価され、それはマッハ数として知られる。
【００２７】
図１を参照すると、第１モデル（２Ｄ−１）１００は２１個の速度を含む２次元モデルである。これら２１個の速度の内、１個（１０５）は動いていない粒子を表し、３組の４速度は、正規化速度（ｒ）（１１０−１１３）、正規化速度の２倍（２ｒ）（１２０−１２３）又は正規化速度の３倍（３ｒ）（１３０−１３３）で、格子のｘ軸又はｙ軸に沿って正方向又は負方向に動いている粒子を表し、２組の４速度は、正規化速度（ｒ）（１４０−１４３）又は正規化速度の２倍（２ｒ）（１５０−１５３）で、２組の４速度は、格子のｘ軸とｙ軸の両方に対して動いている粒子を表す。
【００２８】
さらに図２に図示するように、第２モデル（３Ｄ−１）２００は３９個の速度を含む３次元モデルであり、各速度は図２の矢印の内の１つにより表されている。これら３９個の速度の内、１個は動いていない粒子を表し、３組の６速度は、正規化速度（ｒ）、正規化速度の２倍（２ｒ）又は正規化速度の３倍（３ｒ）で、格子のｘ軸、ｙ軸又はｚ軸に沿って正方向又は負方向に動いている粒子を表し、８個は、正規化速度（ｒ）でｘ、ｙ、ｚの全３格子軸に対して動いている粒子を表し、１２個は、正規化速度の２倍（２ｒ）でｘ、ｙ、ｚの格子軸の内の２つに対して動いている粒子を表す。
【００２９】
例えば、３Ｄ−２は１０１個の速度を含み、２Ｄ−２モデルは３７個の速度を含むなど、さらに複雑なモデルも使用することができる。
【００３０】
３次元モデル３Ｄ−２に関しては、１０１個の速度の内、１個は動いていない粒子を表し（グループ１）、３組の６速度は、正規化速度（ｒ）、正規化速度の２倍（２ｒ）又は正規化速度の３倍（３ｒ）で、格子のｘ軸、ｙ軸又はｚ軸に沿って正方向又は負方向に動いている粒子を表し（グループ２、４及び７）、３組の８速度は、正規化速度（ｒ）、正規化速度の２倍（２ｒ）又は正規化速度の３倍（３ｒ）で、ｘ、ｙ、ｚの全３格子軸に対して動いている粒子を表し（グループ３、８及び１０）、１２個は、正規化速度の２倍（２ｒ）でｘ、ｙ、ｚの格子軸の内の２つに対して動いている粒子を表し（グループ６）、２４個は、正規化速度（ｒ）及び正規化速度の２倍（２ｒ）でｘ、ｙ、ｚの格子軸の内の２つに対して動いており、残る１軸に対しては動いていない粒子を表し（グループ５）、２４個は、正規化速度（ｒ）でｘ、ｙ、ｚの格子軸の内の２つに対して動いており、残る１軸に対しては正規化速度の３倍（３ｒ）で動いている粒子を表す（グループ９）。
【００３１】
２次元モデル２Ｄ−２に関しては、３７個の速度の内、１個は動いていない粒子を表し（グループ１）、３組の４速度は、正規化速度（ｒ）、正規化速度の２倍（２ｒ）又は正規化速度の３倍（３ｒ）で、格子のｘ軸又はｚ軸に沿って正方向又は負方向に動いている粒子を表し（グループ２、４及び７）、２組の４速度は、正規化速度（ｒ）又は正規化速度の２倍（２ｒ）で、格子のｘ軸とｙ軸の両方に対して動いている粒子を表し、８速度は、正規化速度（ｒ）で格子のｘ軸とｙ軸の内の１つに対して、他方の軸に対しては正規化速度の２倍（２ｒ）で動いている粒子を表し、８速度は、正規化速度（ｒ）で格子のｘ軸とｙ軸の内の１つに対して、他方の軸に対しては正規化速度の３倍（３ｒ）で動いている粒子を表す。
【００３２】
前記のＬＢＭモデルは、２次元及び３次元の流れの数値シミュレーション用の効率的で堅牢な、離散速度の動力学的モデルの特定のクラスを提供する。この種のモデルには、離散速度とこれらの速度に関係する重みとの特定の組が含まれる。速度は速度空間においてデカルト座標のグリッド点と一致し、これにより、離散速度モデル、特に格子ボルツマンモデルとして知られる種類の正確で効率的な実行が容易になる。このようなモデルを用いて、流れを高い忠実度でシミュレートすることができる。
【００３３】
図３を参照すると、物理過程シミュレーションシステムは、流体の流れなどの物理過程をシミュレートするために手順３００に従って作動する。シミュレーションに先立って、シミュレーション空間をボクセルの集合としてモデル化する（ステップ３０２）。通常、シミュレーション空間はコンピュータ支援設計（ＣＡＤ）を用いて生成される。例えば、ＣＡＤプログラムは、風洞に置かれた微小デバイスを描くために使用することができる。その後に、適切な分解能を備える格子構造を加えるために、及びシミュレーション空間内の物体及び表面に対処するために、ＣＡＤプログラムにより作成されたデータを処理する。
【００３４】
格子の分解能は、シミュレートするシステムのレイノルズ数に基づいて選択することができる。レイノルズ数は、流れの粘度（ν）、流れ中の物体の特徴的な長さ（Ｌ）、及び流れの特徴的な速度に関係し、
2019053752.tif 000026
方程式（６）
である。
【００３５】
物体の特徴的な長さは、物体の大規模な特徴を表す。例えば、微小デバイスの周りの流れがシミュレートされる場合には、微小デバイスの高さは、特徴的な長さであると見なすことができる。物体の小領域（例えば、自動車のサイドミラー）周りの流れが関心領域である場合には、シミュレーションの分解能を増大させることができ、又は関心領域の周りには分解能が増大した領域を使用することができる。ボクセルの寸法は、格子の分解能が増加するにつれて減少する。
【００３６】
状態空間は、ｆi（ｘ，ｔ）として表され、ｆiは、時間ｔにおいて３次元ベクトルｘにより示される格子部位で状態ｉにある単位容積当たりの要素又は粒子の数（つまり、状態ｉの粒子密度）を表す。既知の時間増分に対しては、粒子数を単にｆi（ｘ）と呼ぶ。１つの格子部位の全ての状態の組み合わせをｆ（ｘ）として示す。
【００３７】
状態の数は、各エネルギ準位内の可能な速度ベクトルの数により決定される。速度ベクトルは、３次元ｘ、ｙ、及びｚを備える空間において整数値の線速度から構成される。状態の数は、複数の種シミュレーションに対して増大される。
【００３８】
各状態ｉは、特定エネルギ準位での異なる速度ベクトルを表す（つまり、エネルギ準位０、１又は２）。各状態の速度ｃiは、以下のように３つの次元の各々におけるその「速度」で示す。
2019053752.tif 000027
方程式（７）
【００３９】
エネルギ準位０状態は、いずれの次元でも動いていない、つまりｃstopped＝（０、０、０）である静止粒子を表す。エネルギ準位１状態は、３つの次元の内の１つで±１の速度、他の２つの次元では速度０を有する粒子を表す。エネルギ準位２状態は、全ての３つの次元で±１の速度、又は３つの次元のうちの１つで±２の速度及び他の２つの次元で速度０を有する粒子を表す。
【００４０】
３つのエネルギ準位の可能な順列の全てを生成することにより、合計３９個の可能な状態（１個のエネルギ０状態、６個のエネルギ１状態、８個のエネルギ３状態、６個のエネルギ４状態、１２個のエネルギ８状態、及び６個のエネルギ９状態）が得られる。
【００４１】
各ボクセル（つまり、各格子部位）は、状態ベクトルｆ（ｘ）により表される。状態ベクトルはボクセルのステータスを完全に規定し、３９個のエントリを含む。３９個のエントリは、１個のエネルギ０状態、６個のエネルギ１状態、８個のエネルギ３状態、６個のエネルギ４状態、１２個のエネルギ８状態、及び６個のエネルギ９状態に対応する。この速度セットを使用することによって、システムは、実現された平衡状態ベクトルに関してマックスウェル−ボルツマン統計を生成することができる。
【００４２】
処理効率のために、ボクセルをマイクロブロックと呼ばれる２ｘ２ｘ２ボリュームにグループ化する。マイクロブロックは、ボクセルの並列処理を可能にし、かつデータ構造に関連づけされたオーバーヘッドを最小化するように構成される。マイクロブロック内のボクセルの略記法は、Ｎi（ｎ）として規定され、ｎは、マイクロブロック内の格子部位の相対位置を表し、ｎ∈｛０，１，２，．．．，７｝である。マイクロブロックを図４に示す。
【００４３】
図５Ａ及び５Ｂを参照すると、表面Ｓはシミュレーション空間においてファセットＦαの集合として表され、
2019053752.tif 000028
方程式（８）
ここで、αは、特定のファセットを列挙する指標である。ファセットは、ボクセルの境界に限定されず、通常は、ファセットが比較的小さな数のボクセルに影響を与えるように、ファセットを隣接するボクセルの大きさの程度に又はそれより僅かに小さくサイズ設定される。表面動力学を実行するために、ファセットに特性を割り当てる。特に、各ファセットＦαは、単位法線（ｎα）、表面積（Ａα）、中心位置（ｘα）、及びファセットの表面動力学特性を説明するファセット分布関数（ｆi（α））を有する。
【００４４】
図６を参照すると、処理効率を向上させるために、異なるレベルの分解能をシミュレーション空間の異なる領域で使用することができる。通常は、物体６５５の周りの領域６５０が最大関心領域なので、最高分解能でシミュレートされる。粘性の影響は物体からの距離と共に減少するので、低レベルの分解能（つまり、拡張されたボクセルボリューム）を使用して、物体６５５から増大する距離で離間する領域６６０、６６５をシミュレートする。同様に、図７に図示するように、より低レベルの分解能を使用して、物体７７５のより有意でない特徴部の周りの領域７７０をシミュレートすることができるが、最高レベルの分解能を使用して、最も有意な特徴部（例えば、先端面及び後端面）の周りの領域７８０をシミュレートする。周辺領域７８５は、最低レベル分解能と最大のボクセルとを用いてシミュレートする。
【００４５】
Ｃ．ファセットに影響を受けるボクセルの特定
図３を再び参照すると、シミュレーション空間をモデル化した状態で（ステップ３０２）、１又は２以上のファセットにより影響を受けるボクセルを特定する（ステップ３０４）。ボクセルは、多くの点でファセットの影響を受ける場合がある。最初に、１又は２以上のファセットと交差するボクセルは、交差しないボクセルと比べて低減された容積を有するという点で影響を受ける。これは、ファセット及びファセットにより表される表面下にある物質がボクセルの一部を占めるので発生する。分画因子Ｐf（ｘ）は、ファセットの影響を受けないボクセルの部分（つまり、流れがシミュレートされる流体又は他の物質が占める可能性がある部分）を示す。交差しないボクセルに対しては、Ｐf（ｘ）は１に等しい。
【００４６】
粒子をファセットに移動するか又は粒子をファセットから受け取ることによって１又は２以上のファセットと相互作用するボクセルも、ファセットにより影響を受けるボクセルとして特定される。ファセットと交差する全てのボクセルは、ファセットから粒子を受け取る少なくとも１つの状態と、粒子をファセットへ移動する少なくとも１つの状態とを含むことになる。ほとんどの場合に、付加的なボクセルは、このような状態を含むことになる。
【００４７】
図８を参照すると、０以外の速度ベクトルｃiを有する各状態ｉに対して、ファセットＦαは、平行六面体Ｇiαにより画定される領域から粒子を受け取るか又はこの領域へ粒子を移動するが、この平行六面体は速度ベクトルｃiとファセットの単位法線とのベクト
ルドット積（｜ｃiｎi｜）の大きさにより規定される高さとファセットの表面積Ａαにより規定される底面とを有するので、平行六面体Ｇiαの容積Ｖiαは次のようになる。
2019053752.tif 000029
方程式（９）
【００４８】
ファセットＦαは、状態の速度ベクトルがファセットに向けられた時に（｜ｃiｎi｜＜０）、容積Ｖiαから粒子を受け取り、状態の速度ベクトルがファセットから離れるように向けられた時に（｜ｃiｎi｜＞０））、その領域へ粒子を移動する。以下に説明するように、この式は、別のファセットが平行六面体Ｇｉαの一部、すなわち、内部のコーナのような非凸面特徴部の近くに発生する可能性がある状態を占有した時には修正する必要がある。
【００４９】
ファセットＦαの平行六面体Ｇiαは、複数のボクセルの内の一部又は全てと重なる場合がある。ボクセル又はその部分の数は、ボクセルのサイズに対するファセットのサイズ、状態のエネルギ、及び格子構造に対するファセットの向きに依存する。影響を受けるボクセルの数は、ファセットのサイズと共に増加する。従って、ファセットのサイズは、上述のように、通常はファセット近くに位置するボクセルのサイズと同程度に又はそれより小さいように選択される。
【００５０】
平行六面体Ｇiαと重なるボクセルＮ（ｘ）の一部分は、Ｖiα（ｘ）として定義する。この項を使用して、ボクセルＮ（ｘ）とファセットＦαとの間を移動する状態ｉ粒子の流束Γiα（ｘ）は、ボクセル（Ｎi（ｘ））における状態ｉ粒子の密度にボクセルと重なる領域の容積（Ｖiα（ｘ））を乗じたものに等しい。
2019053752.tif 000030
方程式（１０）
【００５１】
平行六面体Ｇiαが１又は２以上のファセットと交差する場合、以下の条件が正しく、
2019053752.tif 000031
方程式（１１）
ここで、第１の和は、Ｇiαと重なる全ボクセルに相当し、第２の項は、Ｇiαと交差する全てのファセットに相当する。平行六面体Ｇiαが別のファセットと交差しない場合、この式は以下に帰する。
2019053752.tif 000032
方程式（１２）
【００５２】
Ｄ．シミュレーションの実行
１又は２以上のファセットにより影響を受けるボクセルが特定された状態で（ステップ３０４）、シミュレーションを開始するためにタイマを初期化する（ステップ３０６）。シミュレーションの各時間増分中に、粒子のボクセルからボクセルへの移動は、粒子の表面ファセットとの相互作用に対応する移流段階（ステップ３０８−３１６）によりシミュレートされる。次いで、衝突段階（ステップ３１８）は、各ボクセル内の粒子の相互作用をシミュレートする。次に、タイマを増分する。増分されたタイマがシミュレーションの完了（ステップ３２２）を示さない場合、移流段階及び衝突段階（ステップ３０８−３２０）が繰り返される。増分されたタイマがシミュレーションの完了（ステップ３２２）を示す場合、シミュレーションの結果を記憶及び／又は表示する（ステップ３２４）。
【００５３】
１．表面に対する境界条件
表面との相互作用を正確にシミュレートするために、各ファセットは、４つの境界条件を満たす必要がある。第１に、ファセットが受け取る粒子の合計質量は、ファセットが移送する粒子の合計質量に等しくなければならない（つまり、ファセットに対する正味の質量流束は０でなければならない）。第２に、ファセットが受け取る粒子の合計エネルギは、ファセットが移送する粒子の合計エネルギに等しくなければならない（つまり、ファセットに対する正味のエネルギ流束は０でなければならない）。これら２つの条件は、各エネルギ準位（つまり、エネルギ準位１及び２）での正味質量流束が０に等しいことを必要とすることにより満たすことができる。
【００５４】
他の２つの境界条件は、ファセットと相互作用する粒子の正味運動量に関連する。本明細書では滑り面として言及する、表面摩擦のない表面に対しては、正味接線方向運動量流束は０でなければならず、正味法線方向運動量流束はファセットでの局所的な圧力に等しくなければならない。従って、受け取った合計運動量と移送した合計運動量の、ファセットの法線ｎαに対して垂直な成分（つまり、接線方向成分）は等しくなければならないが、受け取った合計運動量と移送した合計運動量の、ファセットの法線ｎαに対して平行な成分（つまり、法線方向成分）の差はファセットでの局所的な圧力に等しくなければならない。非滑り面に対しては、表面の摩擦は、ファセットが受け取る粒子の合計接線方向運動量に対するファセットが移送する粒子の合計接線方向運動量を摩擦量に関連する係数だけ減少させる。
【００５５】
２．ボクセルからファセットへの集まり
粒子と表面間の相互作用をシミュレートする場合の第１のステップとして、粒子は、ボクセルから集めてファセットへ供給する（ステップ３０８）。前述のように、ボクセルＮ（ｘ）とファセットＦαとの間の状態ｉ粒子の流束は、以下のようになる。
2019053752.tif 000033
方程式（１３）
【００５６】
この式から、ファセットＦα（ｃiｎα＜０）に向けられた各状態ｉに対して、ボクセルからファセットＦαに供給される粒子数は、以下のようになる。
2019053752.tif 000034
方程式（１４）
【００５７】
Ｖiα（ｘ）が０以外の値を有するボクセルのみを合計する必要がある。前述のように、ファセットのサイズは、Ｖiα（ｘ）が少数のボクセルに対してのみ０以外の値を有するように選択される。Ｖiα（ｘ）及びＰf（ｘ）が非整数値を取る場合があるので、Γα（ｘ）は、実数として記憶及び処理される。
【００５８】
３．ファセットからファセットへ移動
次に、粒子をファセット間で移動させる（ステップ３１０）。ファセットＦαの流入状態（ｃiｎα＜０）の平行六面体Ｇiαが別のファセットＦβと交差する場合、ファセットＦαが受け取る状態ｉの粒子の一部は、ファセットＦβから到来することになる。特に、ファセットＦαは、以前の時間増分中にファセットＦβにより生成される状態ｉ粒子の一部を受け取ることになる。この関係を図１０に示すが、ファセットＦβと交差する平行六面体Ｇiαの部分１０００は、ファセットＦαと交差する平行六面体Ｇiβの部分１００５に等しい。前述のように、交差する部分は、Ｖiα（β）として示す。この項を使用して
、ファセットＦβとファセットＦαとの間の状態ｉ粒子の流束は、以下のように説明することができ、
2019053752.tif 000035
方程式（１５）
ここで、Γi（β，ｔ−１）は、前回の時間増分中にファセットＦβによって生成された状態ｉ粒子の尺度である。この式から、ファセットＦα（ｃiｎα＜０）に向けられた各状態ｉに対して、他のファセットによりファセットＦαに供給される粒子数は、
2019053752.tif 000036
方程式（１６）
であり、ファセットへの状態ｉ粒子の全流束は、以下のようになる。
2019053752.tif 000037
方程式（１７）
【００５９】
ファセット分布関数とも呼ばれるファセットに対する状態ベクトルＮ（α）は、ボクセル状態ベクトルのＭエントリに対応するＭエントリを有する。Ｍは、離散格子速度の数である。ファセット分布関数Ｎ（α）の入力状態は、それらの状態への粒子の流束を容積Ｖiαで割ったものに等しく設定され、ｃiｎα＜０に対して以下のようになる。
2019053752.tif 000038
方程式（１８）
【００６０】
ファセット分布関数は、ファセットからの出力流束を生成するためのシミュレーションツールであり、必ずしも実際の粒子を表すものではない。正確な出力流束を生成するために、他状態の分布関数に値を割り当てる。外向きの状態は、内向きの状態をポピュレートするために前述の技術を使用してポピュレートされ、ｃiｎα≧０に対して
2019053752.tif 000039
方程式（１９）
となり、ここで、Γiother（α）は、ΓiIN（α）を生成するための前述の技術を用いるが、流入状態（ｃiｎα＜０）以外の状態（ｃiｎα≧０）にこの技術を適用して決定する。代替技術では、前回の時間ステップからΓiOUT（α）の値を用いて、Γiother（α）を生成することができるので、
2019053752.tif 000040
方程式（２０）
となる。
【００６１】
平行状態（ｃiｎα＝０）に対しては、Ｖiα及びＶiα（ｘ）は共に０である。Ｎi（α）の式において、Ｖｉα（ｘ）は分子に現れ（Γiother（α）の式から）、Ｖiαは、分
母に現れる（Ｎi（α）の式から）。従って、平行状態に対するＮi（α）は、Ｖiα及びＶiα（ｘ）が０に近づく時のＮi（α）の極限値として決定される。
【００６２】
０速度を有する状態（つまり、静止状態及び状態（０、０、０、２）及び（０、０、０、−２））の値は、温度と圧力に関する初期条件に基づいてシミュレーションの最初に初期化される。これらの値は、その後、時間と共に調整される。
【００６３】
４．ファセットの表面力学の実行
次に、各ファセットが前述の４つの境界条件を満たすように、表面力学を実行する。ファセットに対して表面力学を実行するための手順を図１１に示す。最初に、ファセットＦαに垂直な合計運動量は、ファセットにおける粒子の合計運動量Ｐ（α）を確定することによって、以下のように決定され（ステップ１１０５）、全てのｉに対して
2019053752.tif 000041
方程式（２１）
これから、法線方向運動量Ｐn（α）は、以下のように決定される。
2019053752.tif 000042
方程式（２２）
【００６４】
その後、この法線方向運動量は、Ｎn-（α）を生成するために押込／引抜技術（ステップ１１１０）を用いて除去する。この技術に従って、法線方向運動量にのみ影響を及ぼすように状態間で粒子を移動させる。この押込／引抜技術は米国特許第５，５９４，６７１号に記載されており、この特許は、引用により組み込まれている。
【００６５】
その後、Ｎn-（α）の粒子を衝突させてボルツマン分布Ｎn-β（α）を生成する（ステップ１１１５）。流体力学の実行に関して以下で説明するように、ボルツマン分布は、１組の衝突規則をＮn-（α）に適用することによって実現することができる。
【００６６】
その後、ファセットＦαに関する流出流束分布を流入流束分布及びボルツマン分布に基づいて決定する（ステップ１１２０）。最初に、流入流束分布Γi（α）とボルツマン分布との差を以下のように決定する。
2019053752.tif 000043
方程式（２３）
【００６７】
この差分を用いると、流出流束分布は、ｎαｃi＞０に対して
2019053752.tif 000044
方程式（２４）
であり、ここで、ｉ*は、状態ｉと反対の方向を有する状態である。例えば、状態ｉが（１、１、０、０）の場合、状態ｉ*は、（−１、−１、０、０）である。表面摩擦及び他の要因を考慮するために、流出流束分布を以下のようにさらに精緻化することができる。
ｎαｃi＞０に対して、
2019053752.tif 000045
方程式（２５）
ここで、Ｃfは、表面摩擦関数であり、ｔ1αは、ｎαに垂直な第１の接線ベクトルであり、ｔ2αは、ｎα及びｔ1αに垂直な第２の接線ベクトルであり、ΔＮj,1及びΔＮj,2は、状態ｉ及び表示される接線ベクトルのエネルギ（ｊ）に対応する分布関数である。分布関数は以下の方程式に従って決定する。
2019053752.tif 000046
方程式（２６）
ここで、ｊは、エネルギ準位１に対して１、エネルギ準位２に対して２に等しい。
【００６８】
ΓiOUT（α）の方程式の各項の関数は以下の通りである。第１及び第の２項は、衝突がボルツマン分布を生み出すことに有効な範囲で法線方向運動量流束の境界条件を実行するが、接線方向運動量流束偏差を含む。第４及び第５の項は、不十分な衝突による離散効果又は非ボルツマン構造に起因して生ずる場合があるこの偏差を補正する。最後に、第３の項は、表面上の接線方向運動量の所望の変化を強いるために特定の量の表面分画を追加する。摩擦係数Ｃfの生成を以下に説明する。ベクトル演算を伴う全ての項は、シミュレーションの開始に先立って計算することのできる幾何形状因子であることに留意されたい。
【００６９】
これにより、接線方向速度を以下のように決定する。
2019053752.tif 000047
方程式（２７）
ここで、ρは、ファセット分布の密度である。
2019053752.tif 000048
方程式（２８）
【００７０】
前と同様に、流入流束分布とボルツマン分布との差を以下のように決定する。
2019053752.tif 000049
方程式（２９）
【００７１】
従って、流出流束分布は、以下のようになる。
2019053752.tif 000050
方程式（３０）
この方程式は、前回の技術によって求めた流出流束分布の最初の２行に対応するが、特異な接線方向流束に対する補正を必要としない。
【００７２】
あらゆる手法を用いても、結果として生じる流束分布は、運動量流束の条件の全てを満たし、つまり、
2019053752.tif 000051
方程式（３１）
であり、ここで、ｐαは、ファセットＦαにおける平衡圧力であり、ファセットに粒子を与えるボクセルの平均密度及び温度値に基づいており、ｕαは、ファセットの平均速度である。
【００７３】
質量及びエネルギの境界条件が満たされることを保証するために、入力エネルギと出力エネルギとの差を、以下のように各エネルギ準位ｊに対して測定する。
2019053752.tif 000052
方程式（３２）
ここで、指標ｊは、状態ｉのエネルギを示す。その後、このエネルギ差を使用して差の項を生成する。
ｃjiｎα＞０に対して、
2019053752.tif 000053
方程式（３３）
この差の項は、流束が以下になるように流出流束を修正するために使用する。
ｃjiｎα＜０に対して、
2019053752.tif 000054
方程式（３４）
この演算は、質量及びエネルギの流束を補正するが、接線方向の運動量流束を不変のままにする。この調整は、流れがファセットの近傍でほぼ均一で平衡状態に近い場合は小さい。この調整の後に、結果として得られた法線方向運動量流束は、近傍の不均一性又は非平衡特性に起因する補正を加えて近傍の平均特性に基づいて平衡圧力である値に僅かに変更される。
【００７４】
５．ボクセルからボクセルへの移動
再び図３を参照すると、粒子を３次元直線格子に沿ってボクセル間で移動させる（ステップ３１４）。このボクセルからボクセルへの移動は、ファセットとは相互作用しないボクセル（つまり、表面近くに位置しないボクセル）で実行される唯一の移動演算である。典型的なシミュレーションでは、表面と相互作用するほどには近くに位置していないボクセルは、ボクセルの大多数を構成する。
【００７５】
異なる状態の各々は、３次元ｘ、ｙ及びｚの各々において整数速度で格子に沿って移動する粒子を表す。整数速度は、０、±１及び±２を含む。速度の符号は、粒子が対応する軸に沿って移動している方向を示す。
【００７６】
表面と相互作用しないボクセルに関して、移動演算は、計算的にはきわめて簡単である。各時間増分中に、状態の全粒子数を現在のボクセルから目的地ボクセルへ移動させる。同時に、目的地ボクセルの粒子をボクセルから固有の目的地ボクセルへ移動させる。例えば、＋１ｘ方向（１、０、０）に移動しているエネルギ準位１の粒子は、現在のボクセルからｘ方向に＋１、他の方向に対して０のボクセルに移動させる。粒子は、結局、移動（１、０、０）の前に有していたのと同じ状態で目的地ボクセルに行く。ボクセル内の相互作用は、他の粒子及び表面との局所的な相互作用に基づいて、その状態に対して粒子総数を変化させることになる。そうでない場合は、粒子は同じ速度及び同じ方向で格子に沿って移動し続けることになる。
【００７７】
移動演算は、１又は２以上の表面と相互作用するボクセルに対しては少し複雑となる。これにより、１又は２以上の分画粒子は、ファセットに移送される。このような分画粒子のファセットへの移送により分画粒子はボクセル内に残る。これらの分画粒子は、ファセットで占有されるボクセルへ移送される。例えば、図９を参照すると、ボクセル９０５の状態ｉ粒子の部分９００がファセット９１０に移動した場合（ステップ３０８）、残りの部分９１５は、ボクセル９２０へ移動させるが、そこにはファセット９１０が位置しており、状態ｉの粒子はそこからファセット９１０へ向けられる。従って、状態粒子数が２５に等しく、Ｖiα（ｘ）が０．２５に等しい場合（つまり、ボクセルの４分の１が平行六面体と交差する）、６．２５個の粒子は、ファセットＦαに移動させることになり、１８．７５個の粒子は、ファセットＦαにより占有されたボクセルへ移動させることになる。複数のファセットが単一のボクセルと交差することができるので、１又は２以上のファセットにより占有されたボクセルＮ（ｆ）に移送される状態ｉの粒子数は、以下のようになる。
2019053752.tif 000055
方程式（３５）
ここで、Ｎ（ｘ）は、ソースボクセルである。
【００７８】
６．ファセットからボクセルへの散乱
次に、各ファセットからの流出粒子をボクセルに散乱させる（ステップ３１６）。本質的に、このステップは、粒子をボクセルからファセットへ移動させた収集ステップの逆である。ファセットＦαからボクセルＮ（ｘ）へ移動する状態ｉの粒子数は、以下の通りである。
2019053752.tif 000056
方程式（３６）
ここで、Ｐf（ｘ）は、部分的なボクセルの容積減少に対処する。これから、各状態ｉに対して、ファセットからボクセルＮ（ｘ）に向けられる全粒子数は以下のようになる。
2019053752.tif 000057
方程式（３７）
【００７９】
粒子をファセットからボクセルまで散乱させ、周囲のボクセルから流入した粒子と結合して結果を整数化した後に、特定のボクセルの特定の方向がアンダーフローするか（負になる）、又はオーバーフロー（８ビット実装で２５５を超える）する可能性がある。これは、これらの量が許容範囲の値に適合するように切り捨てられた後に、質量、運動量、及びエネルギの取得又は損失をもたらす。このような事象の発生を防ぐために、範囲を越えた質量、運動量、及びエネルギが、問題のある状態の切り捨てに先立って蓄積される。状態が属するエネルギに関して、得られる（アンダーフローによる）又は失われる（オーバーフローによる）値に等しい量の質量は、同じエネルギを有しかつそれ自体オーバーフロー又はアンダーフローを受けないランダムに（又は連続的に）選択された状態に追加される。この質量及びエネルギの追加から生じる追加の運動量は蓄積され、切り捨てられた運動量に追加される。同じエネルギ状態にだけ質量を追加することにより、質量カウンタが０に到達した場合に質量及びエネルギの両方が補正される。最後に、運動量アキュムレータが０に戻るまで、運動量は、押込／引抜技術を利用して補正される。
【００８０】
７．流体力学の実行
最後に、流体力学を実行する（ステップ３１８）。このステップは、微視的動力学又はボクセル内演算と呼ぶことができる。同様に、移流手順は、ボクセル内演算と呼ぶことができる。また、以下に説明する微視的動力学演算を使用して、ファセットで粒子を衝突させてボルツマン分布を生成することができる。
【００８１】
流体力学は、ＢＧＫ衝突モデルとして知られる特定の衝突演算子によって格子ボルツマン方程式において保証される。この衝突モデルは、実際の流体システムの分布の動力学を模倣している。衝突過程は、方程式１の右辺及び式２によりうまく説明することができる。移流ステップの後で、流体システムの保存量、具体的には、密度、運動量、及びエネルギは、方程式３を用いて分布関数から取得する。これらの量から、方程式（２）のｆeqで説明される平衡分布関数は、方程式（４）により完全に特定される。速度ベクトルセットｃi及び重みの選択は表１に列挙され、方程式２と一緒に、巨視的挙動が正確な流体力学の方程式に従うことを保証する。
【００８２】
Ｅ．ハイブリッド熱格子ボルツマン法
格子ボルツマン法（ＬＢＭ）は、従来の計算流体力学（ＣＦＤ）の代わりに、非圧縮性極限近傍での広範な産業用途のために用いられる。熱的用途に関して、Ｈ定理の欠如に起因して、エネルギ保存則を用いたＬＢソルバーは非常に不安定であり、正確な結果を与えるにはより高次のモデルを必要とする。質量及び運動量に加えて、エネルギに対する境界条件をさらに適用するには、複雑なアルゴリズムが必要となる。これらの困難性を回避するために、大部分のＬＢＭ法は、エネルギに関する別の有限差分ソルバー又は格子ボルツマンスカラーのいずれかを使用する。便宜上、以下の方程式の番号付与は（１．０）で再開することになる。
【００８３】
標準等温ＬＢＭは以下の方程式により与えられ、
2019053752.tif 000058
方程式（１．０）
ここで、ｆi（ｘ，ｔ）は、（ｘ，ｔ）における速度値ｃiに対する粒子分布関数であり、Ｃi（ｘ，ｔ）は、粒子衝突演算子である。最も一般的に使用される衝突演算子は、例えば、以下のＢＧＫ形式を有する。
2019053752.tif 000059
方程式（１．１）
【００８４】
方程式（１．１）において、ｆi（ｘ，ｔ）は粒子分布関数であり、τは緩和時間である。標準平衡分布、ｆieqは、以下の通りである。
2019053752.tif 000060
ここで、定数の重み因子ｗi及び格子温度Ｔ0は、格子セットに基づく。重み因子の目的は、予測次数までのモーメント等方性を満たすことである。例えば、Ｄ３Ｑ１９に対して、Ｔ0値は（１／３）であり、重み因子は２次までのモーメント等方性を満たす。前記の等温ＬＢＭは状態方程式Ｐ＝ρＴ0をもたらす。本実施例では、Ｐは圧力（例えば、Ｎｍ-2で）、ρは密度（例えば、ｋｇｍ-3で）、及びＴ0は温度定数である。
【００８５】
ハイブリッド熱ＬＢＭソルバーの場合、温度をＬＢＭに戻す結合は、２次モーメントが実際の圧力をもたらすように平衡分布の１次項を修正することによって実行される。修正されたＬＢＭ方程式は以下の通りである。
2019053752.tif 000061
ここで、
2019053752.tif 000062
である。
【００８６】
方程式（１．４）では、Ｐは実際の状態方程式を満たす。理想気体ではＰ＝ρＲＴであり、Ｒは気体定数である。静止状態の平衡状態における因子βは、
2019053752.tif 000063
により与えられ、ｑは運動格子の総数である
【００８７】
前記の平衡状態の修正は、０次及び１次のモーメントを変えず、２次モーメントに関して正確な熱力学的圧力をもたらす。
2019053752.tif 000064
【００８８】
この手法は、従来のＬＢＭに対して最小限の修正しか必要とせず、２次モーメントとして得られる圧力は、全格子方向を含む。このアルゴリズムの大きな欠点は、最大次数の項（ｌｅａｄｉｎｇ  ｏｒｄｅｒ  ｔｅｒｍ）の変化に起因して該アルゴリズムが課する温度範囲の制限である。例えば、理想気体に関して、Ｄ３Ｑ１９格子セットに対する静止状態の分布は以下の通りである。
2019053752.tif 000065
【００８９】
ＲＴ＞１．５Ｔ0となる高温に対しては、静止状態の平衡状態は負となるので不安定となる。また、ＲＴ＜１．５Ｔ0となる低温に対しては、運動状態の正値範囲は減少し、同様に高い速度での不安定性につながる。この欠点に起因して、大部分のハイブリッド熱ＬＢＭモデルは、低温範囲及び低マッハ数の用途に制限される。
【００９０】
Ｆ．力スキーム（ｆｏｒｃｅ  ｓｃｈｅｍｅ）を用いた手法
別の手法において、等温ＬＢＭから生じる圧力と実際の状態方程式との差異は、以下のように計算される。
2019053752.tif 000066
【００９１】
前記の項Ｆは、外部物体力として加える。物体力の離散的格子作用を防ぐために、速度定義に対する修正を有する高次の力の項を用いる必要がある。付加される力の項は圧力勾配なので、平衡分布の正値範囲は悪影響を受けない。しかしながら、この付加的な力の項の計算は、ＬＢＭの等方性に影響を及ぼす有限差分計算を必要とする。壁面近傍での勾配の計算はやはり難しく、これは高熱流束シミュレーションに対して非ゼロの法線速度をもたらす可能性がある。
Ｇ．温度結合アルゴリズム
【００９２】
前記セクションで検討したように、従来の結合方法自体、長所及び短所を有する。平衡状態の最大次数を変更する場合（セクションＥ）、得られる圧力勾配は全格子方向を含み、実行するのが容易である。しかしながら、平衡状態が変わると正値範囲が減少するので、ＬＢソルバーの使用可能範囲が縮小する。物体力の項の導入は（セクションＦ）、平衡状態に影響を及ぼさないので、広範な速度及び温度に対して安定性を提供する。しかしながら、有限ステンシルを用いた勾配の計算は、ＬＢスキームの等方性に影響を及ぼし、グリッド依存性、格子配向依存性、及び正確な境界条件を適用する際の複雑性のようなエラーを発生させる。
【００９３】
本セクションでは、現行の方法の欠点を解消するための温度結合アルゴリズムを説明する。線形安定解析法は、平衡状態が負の場合には必ず衝突が数値的不安定性をもたらすことを示す。さらに、Ｃｈａｐｍａｎｎマルチスケール展開は、圧力勾配が移流過程の結果であることを示す。従って、温度結合を衝突ステップから取り除き、以下のように移流中に導入することが可能である。
2019053752.tif 000067
ここで、平衡分布は式（１．２）と同じであり、圧力結合項は以下の通りである。
2019053752.tif 000068
【００９４】
本実施例では、ρは流体密度、Ｔ0は一定格子温度、Ｐは流体容積中の圧力、Ｔは移送中の流体温度、及びｗiは一定重み因子である。前記方程式（３．１）において、第１の部分は衝突過程の数学的表現（つまり、ｆi（ｘ，ｔ）＋Ｃi（ｘ，ｔ））であり、Ｃiは衝突演算子である。この方程式の第２の部分は、熱力学的ステップ（例えば、ｇi（ｘ，ｔ＋Δｔ）−ｇi（ｘ，ｔ））である。マルチスケール解析を行うことによって、結合項の２次モーメントは、質量及び運動量保存の式における他の量に影響を及ぼすことなく正確な圧力勾配をもたらすことを示すことができる。
【００９５】
温度結合の別の方法は、以下のように時間変動を結合項に加えることである。
2019053752.tif 000069
【００９６】
方程式（３．１）と方程式（３．２）との相違は、新たな力の項ｇがｔ＋Δｔで計算される点にある。このことは、追加の時間微分項を質量保存則に導入する。これは静止状態を以下のように変更することによって簡単に取り除くことができる。
2019053752.tif 000070
【００９７】
移流中に温度を導入するためのこれらの技術は、従来のＬＢＭ法に対して最小限の変更を必要とするだけである。経時変化する温度結合のためのアルゴリズムである方程式（３．３）を図１２に関して以下に説明する。
【００９８】
図１２を参照すると、システムは粒子移送の分布関数を生成する処理１２００を実行し、熱力学的ステップは衝突ステップから切り離される。本実施例では、分布関数は移流ステップを含み、熱力学的ステップは、熱力学的ステップで粒子衝突ステップを強化することではなく、熱力学的ステップで移流ステップを強化することによって、分布関数の一部に含まれる。一般に、移流ステップは、粒子の移送を含む（例えば、一方の領域から他方の領域へ水平方向に）。本実施例では、熱力学的ステップは、移送中の流体容積の温度を含む。
【００９９】
作動時、システムは、格子速度セットにおいて、流体容積中の粒子の移送をシミュレートする（１２０２）。この移送は粒子間の衝突を引き起こす。また、システムは、粒子移送に関する分布関数を生成し（１２０３）、分布関数は熱力学的ステップ及び粒子衝突ステップを含み、熱力学的ステップは、粒子衝突ステップから実質的に独立して切り離されている。この生成過程は、動作１２０４、１２０６、１２０８、１２１０、１２１２、１２１４、１２１６、１２１８を含む。
【０１００】
図１２の例では、システムは、特定時間ｔにおける流体容積中の特定位置ｘでの衝突に対する衝突後分布関数ｆi′（ｘ，ｔ）を決定する（１２０４）。ここで、
2019053752.tif 000071
であり、ここで、Ｃiは衝突演算子、ｆiは衝突前の粒子に関する分布関数である。
【０１０１】
本実施例では、システムは、衝突後分布関数ｆi′（ｘ，ｔ）から分画部分ｇi（ｘ，ｔ）を差し引いて（１２０６）、移流前粒子密度関数
2019053752.tif 000072
を取得するようになっており、ｇi（ｘ，ｔ）は、移流中に同じ位置に維持する必要がある粒子の分画部分を表す（例えば、移流ステップ中に移流しない粒子部分に対する分布関数）。本実施例では、システムは、粒子を、以下に説明するように移流する第１の部分と、移流ステップ中に移流しない第２の部分とに分割する。
【０１０２】
システムは、粒子の第２の部分の、時間ｔ＋Δｔにおける流体の容量中の他の位置（例えば、隣接セル）への移流をシミュレートし（１２０８）、他の位置は（ｘ＋ｃiΔｔ）として表され、ｃiは衝突前の粒子の速度ベクトルであり、Δｔは特定時間ｔと別の時点との間の間隔である。１つの実施例では、粒子は、離散的な粒子速度方向に沿って隣接セルへ移流する。移流する粒子の一部は、粒子の第２の部分（例えば、ｆi′（ｘ，ｔ）−ｇi（ｘ，ｔ））である。システムは、移流のシミュレーションに基づいて移流粒子の密度分布関数
2019053752.tif 000073
を取得し、ここで
2019053752.tif 000074
である。
2019053752.tif 000075
は、以下に説明するように位置ｘ−ｃiから移流した位置ｘにおける分布関数である。本実施例では、ｆi〃（ｘ−ｃiΔｔ，ｔ）は、時間ｔにおける位置ｘ−ｃiΔｔでの移流前粒子密度分布関数を表す。項ｘ−ｃiΔｔは、粒子が位置ｘへ移流する別の位置（例えば、隣接セル内）を表す。この別の位置ｘ−ｃiΔｔは、元の位置ｘと、粒子が時間間隔Δｔで（速度値ｃiに基づいて）移動可能な距離とに基づく。位置ｘ−ｃiΔｔの粒子が位置ｘへ移流している場合、方程式
2019053752.tif 000076
は、時間ｔ＋Δｔにおける粒子位置ｘでの粒子の分布関数が、時間ｔにおける位置ｘ−ｃiΔｔでの移流前粒子密度分布関数と同じであることを表す。
【０１０３】
システムは、先に差し引いた部分ｇi（ｘ，ｔ）を密度分布関数
2019053752.tif 000077
（隣接セルからの移流）に戻して加えて（１２１２）、以下の移流後密度分布関数を生成するようになっており、
2019053752.tif 000078
これは、以下のように表すこともできる。
2019053752.tif 000079
【０１０４】
本実施例では、先に差し引いた部分ｇi（ｘ，ｔ）を密度分布関数に戻して加えることによって、質量が保存されることになる。さらに、前記の方程式の全項の位置ｘに対してｃiを付加することで理解できるように、この付加により方程式（３．１）が回復する。
【０１０５】
システムは、時間ｔ＋Δｔにおける位置ｘでの粒子の質量、運動量、及び温度を計算する（１２１４）。システムは、計算された温度、質量、及び運動量を用いてｇi（ｘ，ｔ＋Δｔ）を決定する（１２１６）。システムは、温度結合の時間変動を考慮するために、ｇi（ｘ，ｔ＋Δｔ）−ｇi（ｘ，ｔ）の差分を各運動状態ｆi（ｘ，ｔ＋Δｔ）に加える（１２１８）。この操作は、以下のように式（３．３）を回復させる。
2019053752.tif 000080
【０１０６】
前述の運動状態に加えられた粒子の数は静止状態から差し引かれて、方程式（３．４）のように質量を保存する。方程式（３．４）におけるｇ項の総和は、非静止状態からの時間ｔ及び時間ｔ＋Δｔのｇグループの粒子の間の粒子数の差分を表す。システムは、全体の質量保存を保証するように、差分を補うために静止状態に関する方程式に総和の項を導入する。システムは、操作１２０４、１２０６、１２０８、１２１０、１２１２、１２１４、１２１６、１２１８を繰り返す。１つの時間ステップｔから次の時間ステップｔ＋Δｔへの全サイクルは、以下同様に繰り返される。
【０１０７】
１つの実施例では、方程式（３．１）及び（３．３）は、位置ｘで速度ｃiを有する粒子が位置ｘ−ｃiから移流するのではなく、位置ｘ＋ｃiで速度ｃiを有する粒子が位置ｘから移流する場合に得られる。方程式（３．１）を得るために、ｇが位置ｘ＋ｃiΔｔにおいて時間ｔで計算される（操作１２１６で）以外、システムは処理１２００に含まれる操作を繰り返す。
【０１０８】
本明細書に説明される技術を用いると、平衡分布（ｆieq）は、分布関数（例えば、ｆi（ｘ＋ｃiΔｔ，ｔ＋Δｔ））では不変なので、安定範囲は、等温ＬＢＭのものと同じままである。これにより、高マッハ数及び高温度範囲用途のシミュレーションが可能となる。加えて、移流の結果として計算された圧力勾配は、全格子方向に関与し、これによりＬＢＭソルバーの等方性が維持される。有限差分近似法を用いて計算された力は付加されないので、局所的な運動量保存が保証される。速度定義は、力の方法とは違って正確であり、離散的な格子作用を除去するために代わりの速度定義を必要とする。本明細書に記載の技術は、大規模な圧力勾配計算及び複雑な境界条件処理を伴わないので、計算コストが低い。また、これらの技術は、壁近傍での正確な圧力勾配を保証する、より簡単な境界条件を有するので、正確な速度プロファイルをもたらす。これは、高温用途の熱流束の予測に関して重要である。これらの技術は、磁力、重力、多相用途での要素間力などの、保存力を伴う任意のＬＢＭ用途に容易に拡張することができる。
【０１０９】
図１３は、ネットワーク環境１３００の構成要素のブロック図である。ネットワーク環境１３００にはさらにシステム１３０２を含み（本明細書で説明する技術を実施するための）、そのシステムにはメモリ１３０４、バスシステム１３０６、及びプロセッサ１３０８が含まれる。メモリ１３０４は、ハードドライブと、ダイナミック・ランダム・アクセス・メモリ、機械可読ハードウェア記憶装置、機械可読媒体又は他タイプの固定機械可読記憶装置などのランダムアクセスメモリとを備えることができる。例えば、データバス及びマザーボードを含むバスシステム１３０６は、システム１３０２の構成要素間のデータ通信を確立し制御するために使用することができる。プロセッサ１３０８は、１又は２以上のマイクロプロセッサ及び／又は処理装置を含むことができる。一般に、プロセッサ１３０８は、データの受信と記憶、及びネットワーク上の通信（図示せず）が可能な任意の適切なプロセッサ及び／又は論理回路を含むことができる。
【０１１０】
システム１３０２は、サーバ、分散型計算システム、デスクトップ型コンピュータ、ラップトップ、携帯電話機、ラックマウント・サーバなどのデータ受信が可能な種々の計算装置のいずれでもよい。システム１３０２は、単一のサーバとすること、或いは同じ場所に又は別の場所にある一群のサーバとすることができる。図示するシステム１３０２は、入力／出力（「Ｉ／Ｏ」）インタフェース１３１０を経由してデータを受信することができる。Ｉ／Ｏインタフェース１３１０は、イーサネット（登録商標）、無線ネットワーク・インタフェース、光ファイバ・ネットワーク・インタフェース、モデムなどの、ネットワーク上のデータを受信可能な任意タイプのインタフェースとすることができる。システム１３０２は、速度モデル、シミュレーションデータなどを記憶するように構成することができるデータ格納部１３１２との通信のために設定される。
【０１１１】
実施形態は、デジタル電子回路、又はコンピュータハードウエア、ファームウェア、ソフトウェア、又はその組合せに実施することができる。本明細書に説明する技術の装置は、プログラマブルプロセッサによる実行のために機械可読媒体（例えば、ストレージデバイス）に有形に具現化又は記憶されたコンピュータプログラム製品に実施することができ、方法アクションは、入力データに対して演算して出力を生成することにより、本明細書に説明する技術の演算を実行する命令のプログラムを実行するプログラマブルプロセッサによって実行することができる。本明細書に説明する技術は、データ及び命令をデータストレージシステム、少なくとも１つの入力デバイス、及び少なくとも１つの出力デバイスから受信して、データ及び命令をそれらに送信するように結合された少なくとも１つのプログラマブルプロセッサを含むプログラマブルシステム上で実行可能である１つ又はそれよりも多くのコンピュータプログラムで実行することができる。各コンピュータプログラムは、高レベル手順又はオブジェクト指向プログラミング言語、又は必要に応じてアセンブリ又はマシン語で実行することができ、いずれにしても、言語は、コンパイル又は解釈された言語とすることができる。
【０１１２】
適切なプロセッサには、一例として、汎用及び専用マイクロプロセッサの両方が挙げられる。一般的に、プロセッサは、命令及びデータを読取専用メモリ及び／又はランダムアクセスメモリから受信することになる。一般的に、コンピュータは、データファイルを記憶する１つ又はそれよりも多くの大容量ストレージデバイスを含むことになり、そのようなデバイスには、内蔵ハードディスク及び取外し可能ディスクのような磁気ディスク、光磁気ディスク、光ディスクが挙げられる。コンピュータプログラム命令及びデータを有形に具現化するのに適するストレージデバイスには、一例として、ＥＰＲＯＭ、ＥＥＰＲＯＭ、及びフラッシュメモリデバイスのような半導体メモリデバイス、内部ハードディスク及び着脱式ディスクのような磁気ディスク、光磁気ディスク、及びＣＤ−ＲＯＭディスクを含む全ての形態の不揮発性メモリが挙げられる。上述のいずれも、ＡＳＩＣ（特定用途向け集積回路）によって補足されるか、又はＡＳＩＣに組み込むことができる。
【０１１３】
いくつかの実施を説明した。それにも関わらず、本発明の精神及び範囲から逸脱することなく様々な修正を行うことができることは理解されるであろう。従って、他の実施も以下の特許請求の範囲内である
【符号の説明】
【０１１４】
１３００  ネットワーク環境
１３０２  システム
１３０４  メモリ
１３０６  バスシステム
１３０８  処理装置
１３１０  インタフェース
１３１２  データ格納部
【図１】
2019053752.tif 000108
【図２】
2019053752.tif 000109
【図３】
2019053752.tif 000110
【図４】
2019053752.tif 000111
【図５Ａ−Ｂ】
2019053752.tif 000112
【図６】
2019053752.tif 000113
【図７】
2019053752.tif 000114
【図８】
2019053752.tif 000115
【図９】
2019053752.tif 000116
【図１０】
2019053752.tif 000117
【図１１】
2019053752.tif 000118
【図１２】
2019053752.tif 000119
【図１３】
2019053752.tif 000120
【手続補正書】
【提出日】平成30年12月3日(2018.12.3)
【手続補正１】
【補正対象書類名】特許請求の範囲
【補正対象項目名】全文
【補正方法】変更
【補正の内容】
【特許請求の範囲】
【請求項１】
所定の一定格子温度に対する温度の範囲にわたる流体の流れをシミュレートするための粒子に関する安定した分布関数を生成するコンピュータ実施方法であって、コンピュータによって
格子速度セットにおいて流体容積中での粒子の移送をシミュレートする段階であって、前記移送が前記粒子間の衝突を引き起こす、前記シミュレートする段階と、
前記粒子の移送に関する分布関数を生成する生成段階であって、
前記流体容積中の粒子の熱力学的特性から実質的に独立した、前記容積中の位置に関する衝突後分布関数を決定する段階と、
分画部分を前記衝突後分布関数から差し引く段階であって、前記粒子の一部は流体容積中の別の位置に移流する、前記差し引く段階と、
前記容積の熱力学的特性に従って、前記粒子の一部の移送に関する移送分布関数を計算する段階と、
前記衝突後分布関数と前記移送分布関数から前記分布関数を生成する段階と、を含む前記生成段階と、
を含む方法。
【請求項２】
前記分布関数は、熱力学的ステップ，粒子衝突ステップ及び移流ステップを含み、前記熱力学的ステップは、前記熱力学的ステップで前記粒子衝突ステップを強化するのではなく、前記熱力学的ステップで前記移流ステップを強化することによって、前記分布関数に含まれる、請求項１に記載の方法。
【請求項３】
前記熱力学的ステップは、前記移送中の前記流体容積の温度を含む、請求項２に記載の方法。
【請求項４】
特定時間ｔにおける前記流体容積中の特定位置ｘでの前記衝突に対する衝突後分布関数は、ｆi′（ｘ，ｔ）であり、ここで、ｆi′（ｘ，ｔ）＝ｆi（ｘ，ｔ）＋Ｃi（ｘ，ｔ）であり、Ｃiは衝突演算子、ｆiは前記衝突前の前記粒子に関する分布関数であり、
前記分画部分は、ｇi（ｘ，ｔ）であり、移流前粒子密度分布関数ｆi′′（ｘ，ｔ）＝ｆi′（ｘ，ｔ）−ｇi（ｘ，ｔ）を得るために、前記衝突後分布関数ｆi′（ｘ，ｔ）から差し引かれ、ｆi′（ｘ，ｔ）−ｇi（ｘ，ｔ）により示される前記粒子の一部は流体容積中の別の位置に移流し、ｇi（ｘ，ｔ）は移流しない粒子の分布関数を表し、
前記生成段階は、さらに、
時間ｔ＋Δｔにおいて前記粒子の一部が前記流体容積中の別の位置へ移流することをシミュレートする段階であって、前記別の位置は（ｘ＋ｃiΔｔ）で表され、ｃiは衝突前の前記粒子の速度ベクトルであり、Δｔは特定の時間ｔと別の時点との間の間隔である、段階と、
前記移流のシミュレーションに基づいて、移流粒子の密度分布関数
2019053752.tif 000121
を得る段階であって、
2019053752.tif 000122
であり、
2019053752.tif 000123
は位置ｘ−ｃiΔｔから移流した位置ｘにおける粒子の分布関数である、段階と、
先に差し引いた部分ｇi（ｘ，ｔ）を前記密度分布関数に戻して加えて移流後密度分布関数
2019053752.tif 000124
を生成する段階と、
時間ｔ＋Δｔにおける位置ｘでの前記粒子の質量、運動量、及び温度を計算する段階と、
計算された温度、質量、及び運動量を用いてｇi（ｘ，ｔ＋Δｔ）を決定する段階と、
差分ｇi（ｘ，ｔ＋Δｔ）−ｇi（ｘ，ｔ）を運動状態ｆi（ｘ，ｔ）＋Ｃi（ｘ，ｔ）に加える段階と、
を含む、請求項１に記載の方法。
【請求項５】
ｇiは、
2019053752.tif 000125
に従って定義され、
ρは流体密度であり、
Ｔ0は一定格子温度であり、
Ｐは前記流体容積中の圧力であり、
Ｔは前記計算された温度であり、
ｗiは一定の重み因子である、
請求項４に記載の方法。
【請求項６】
前記加える段階が為された後の分布関数は、
2019053752.tif 000126
であり、
ｘは前記容積中の特定位置であり、
ｔは特定の第１の時点であり、
ｉは前記格子速度セットでのインデックスであり、
ｃiは衝突前の前記粒子の速度ベクトルであり、
Ｃiは衝突演算子であり、
Δｔは前記第１の時点と第２の時点との間の間隔であり、
ｇiは熱力学的ステップであり、
ｆi（ｘ，ｔ）は時間ｔにおける位置ｘでの前記粒子に対する分布関数である、
請求項４に記載の方法。
【請求項７】
静止状態を
2019053752.tif 000127
に一致するように修正することによって、前記質量を保存する段階をさらに含む、請求項４に記載の方法。
【請求項８】
前記粒子衝突ステップは、等温平衡分布関数を含む、請求項２に記載の方法。
【請求項９】
前記分布関数は、時間ｔ＋Δｔにおける位置ｘ＋ｃiΔｔでの前記粒子に対する分布関数であり、
2019053752.tif 000128
で表され、
ｘは前記容積中の特定位置であり、
ｔは特定の第１の時点であり、
ｉは前記格子速度セットでのインデックスであり、
Ｃiは衝突演算子であり、
ｃiは衝突前の前記粒子の速度ベクトルであり、
Δｔは前記第１の時点と第２の時点との間の間隔であり、
ｇiは熱力学的ステップであり、
ｆi（ｘ，ｔ）は時間ｔにおける位置ｘでの前記粒子に対する分布関数である、
請求項１に記載の方法。
【請求項１０】
ｇiは、
2019053752.tif 000129
に従って定義され、
ρは流体密度であり、
Ｔ0は一定格子温度であり、
Ｐは前記流体容積中の圧力であり、
Ｔは流体の実際の温度であり、
ｗiは一定の重み因子である、
請求項９に記載の方法。
【請求項１１】
前記格子速度セットは、格子ボルツマン法に基づく、請求項１に記載の方法。
【請求項１２】
１又は２以上の処理装置で実行可能な命令を記憶する１又は２以上の機械可読ハードウェア記憶装置であって、
前記命令が実行されることにより所定の一定格子温度に対する温度の範囲にわたる流体の流れをシミュレートするための粒子に関する安定した分布関数を生成する方法が実行されるものであり、
前記方法は、
格子速度セットにおいて流体容積中での粒子の移送をシミュレートする段階であって、前記移送が前記粒子間の衝突を引き起こす、前記シミュレートする段階と、
前記粒子の移送に関する分布関数を生成する生成段階であって、
前記流体容積中の粒子の熱力学的特性から実質的に独立した、前記容積中の位置に関する衝突後分布関数を決定する段階と、
分画部分を前記衝突後分布関数から差し引く段階であって、前記粒子の一部は流体容積中の別の位置に移流する、前記差し引く段階と、
前記容積の熱力学的特性に従って、前記粒子の一部の移送に関する移送分布関数を計算する段階と、
前記衝突後分布関数と前記移送分布関数から前記分布関数を生成する段階と、を含む前記生成段階と、を含む、１又は２以上の機械可読ハードウェア記憶装置。
【請求項１３】
前記分布関数は、熱力学的ステップ，粒子衝突ステップ及び移流ステップを含み、前記熱力学的ステップは、前記熱力学的ステップで前記粒子衝突ステップを強化するのではなく、前記熱力学的ステップで前記移流ステップを強化することによって、前記分布関数に含まれる、請求項１２に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項１４】
前記熱力学的ステップは、前記移送中の前記流体容積の温度を含む、請求項１３に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項１５】
特定時間ｔにおける前記流体容積中の特定位置ｘでの前記衝突に対する衝突後分布関数は、ｆi′（ｘ，ｔ）であり、ここで、ｆi′（ｘ，ｔ）＝ｆi（ｘ，ｔ）＋Ｃi（ｘ，ｔ）であり、Ｃiは衝突演算子、ｆiは前記衝突前の前記粒子に関する分布関数であり、
前記分画部分は、ｇi（ｘ，ｔ）であり、移流前粒子密度分布関数ｆi′′（ｘ，ｔ）＝ｆi′（ｘ，ｔ）−ｇi（ｘ，ｔ）を得るために、前記衝突後分布関数ｆi′（ｘ，ｔ）から差し引かれ、ｆi′（ｘ，ｔ）−ｇi（ｘ，ｔ）により示される前記粒子の一部は流体容積中の別の位置に移流し、ｇi（ｘ，ｔ）は移流しない粒子の分布関数を表し、
前記生成段階は、さらに、
時間ｔ＋Δｔにおいて前記粒子の一部が前記流体容積中の別の位置へ移流することをシミュレートする段階であって、前記別の位置は（ｘ＋ｃiΔｔ）で表され、ｃiは衝突前の前記粒子の速度ベクトルであり、Δｔは特定の時間ｔと別の時点との間の間隔である、段階と、
前記移流のシミュレーションに基づいて、移流粒子の密度分布関数
2019053752.tif 000130
を得る段階であって、
2019053752.tif 000131
であり、
2019053752.tif 000132
は位置ｘ−ｃiΔｔから移流した位置ｘにおける粒子の分布関数である、段階と、
先に差し引いた部分ｇi（ｘ，ｔ）を前記密度分布関数に戻して加えて移流後密度分布関数
2019053752.tif 000133
を生成する段階と、
時間ｔ＋Δｔにおける位置ｘでの前記粒子の質量、運動量、及び温度を計算する段階と、
計算された温度、質量、及び運動量を用いてｇi（ｘ，ｔ＋Δｔ）を決定する段階と、
差分ｇi（ｘ，ｔ＋Δｔ）−ｇi（ｘ，ｔ）を運動状態ｆi（ｘ，ｔ）＋Ｃi（ｘ，ｔ）に加える段階と、
を含む、請求項１２に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項１６】
ｇiは、
2019053752.tif 000134
に従って定義され、
ρは流体密度であり、
Ｔ0は一定格子温度であり、
Ｐは前記流体容積中の圧力であり、
Ｔは前記計算された温度であり、
ｗiは一定の重み因子である、
請求項１５に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項１７】
前記加える段階が為された後の分布関数は、
2019053752.tif 000135
であり、
ｘは前記容積中の特定位置であり、
ｔは特定の第１の時点であり、
ｉは前記格子速度セットでのインデックスであり、
ｃiは衝突前の前記粒子の速度ベクトルであり、
Ｃiは衝突演算子であり、
Δｔは前記第１の時点と第２の時点との間の間隔であり、
ｇiは熱力学的ステップであり、
ｆi（ｘ，ｔ）は時間ｔにおける位置ｘでの前記粒子に対する分布関数である、
請求項１５に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項１８】
前記方法は、静止状態を
2019053752.tif 000136
に一致するように修正することによって、前記質量を保存する段階をさらに含む、請求項１５に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項１９】
前記粒子衝突ステップは、等温平衡分布関数を含む、請求項１３に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項２０】
前記分布関数は、時間ｔ＋Δｔにおける位置ｘ＋ｃiΔｔでの前記粒子に対する分布関数であり、
2019053752.tif 000137
で表され、
ｘは前記容積中の特定位置であり、
ｔは特定の第１の時点であり、
ｉは前記格子速度セットでのインデックスであり、
Ｃiは衝突演算子であり、
ｃiは衝突前の前記粒子の速度ベクトルであり、
Δｔは前記第１の時点と第２の時点との間の間隔であり、
ｇiは熱力学的ステップであり、
ｆi（ｘ，ｔ）は時間ｔにおける位置ｘでの前記粒子に対する分布関数である、
請求項１２に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項２１】
ｇiは、
2019053752.tif 000138
に従って定義され、
ρは流体密度であり、
Ｔ0は一定格子温度であり、
Ｐは前記流体容積中の圧力であり、
Ｔは流体の実際の温度であり、
ｗiは一定の重み因子である、
請求項２０に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項２２】
前記格子速度セットは、格子ボルツマン法に基づく、請求項１２に記載の１又は２以上の機械可読ハードウェア記憶装置。
【請求項２３】
１又は２以上の処理装置と、
前記１又は２以上の処理装置で実行可能な命令を記憶する１又は２以上の機械可読ハードウェア記憶装置と、
を含む、所定の一定格子温度に対する温度の範囲にわたる流体の流れをシミュレートするための粒子に関する安定した分布関数を生成するシステムであって、
前記命令が実行されることにより方法が実行されるものであり、
前記方法は、
格子速度セットにおいて流体容積中での粒子の移送をシミュレートする段階であって、前記移送が前記粒子間の衝突を引き起こす、前記シミュレートする段階と、
前記粒子の移送に関する分布関数を生成する生成段階であって、
前記流体容積中の粒子の熱力学的特性から実質的に独立した、前記容積中の位置に関する衝突後分布関数を決定する段階と、
分画部分を前記衝突後分布関数から差し引く段階であって、前記粒子の一部は流体容積中の別の位置に移流する、前記差し引く段階と、
前記容積の熱力学的特性に従って、前記粒子の一部の移送に関する移送分布関数を計算する段階と、
前記衝突後分布関数と前記移送分布関数から前記分布関数を生成する段階と、を含む前記生成段階と、を含む、システム。
【請求項２４】
前記分布関数は、熱力学的ステップ，粒子衝突ステップ及び移流ステップを含み、前記熱力学的ステップは、前記熱力学的ステップで前記粒子衝突ステップを強化するのではなく、前記熱力学的ステップで前記移流ステップを強化することによって、前記分布関数に含まれる、請求項２３に記載のシステム。
【請求項２５】
前記熱力学的ステップは、前記移送中の前記流体容積の温度を含む、請求項２４に記載のシステム。
【請求項２６】
特定時間ｔにおける前記流体容積中の特定位置ｘでの前記衝突に対する衝突後分布関数は、ｆi′（ｘ，ｔ）であり、ここで、ｆi′（ｘ，ｔ）＝ｆi（ｘ，ｔ）＋Ｃi（ｘ，ｔ）であり、Ｃiは衝突演算子、ｆiは前記衝突前の前記粒子に関する分布関数であり、
前記分画部分は、ｇi（ｘ，ｔ）であり、移流前粒子密度分布関数ｆi′′（ｘ，ｔ）＝ｆi′（ｘ，ｔ）−ｇi（ｘ，ｔ）を得るために、前記衝突後分布関数ｆi′（ｘ，ｔ）から差し引かれ、ｆi′（ｘ，ｔ）−ｇi（ｘ，ｔ）により示される前記粒子の一部は流体容積中の別の位置に移流し、ｇi（ｘ，ｔ）は移流しない粒子の分布関数を表し、
前記生成段階は、さらに、
時間ｔ＋Δｔにおいて前記粒子の一部が前記流体容積中の別の位置へ移流することをシミュレートする段階であって、前記別の位置は（ｘ＋ｃiΔｔ）で表され、ｃiは衝突前の前記粒子の速度ベクトルであり、Δｔは特定の時間ｔと別の時点との間の間隔である、段階と、
前記移流のシミュレーションに基づいて、移流粒子の密度分布関数
2019053752.tif 000139
を得る段階であって、
2019053752.tif 000140
であり、
2019053752.tif 000141
は位置ｘ−ｃiΔｔから移流した位置ｘにおける粒子の分布関数である、段階と、
先に差し引いた部分ｇi（ｘ，ｔ）を前記密度分布関数に戻して加えて移流後密度分布関数
2019053752.tif 000142
を生成する段階と、
時間ｔ＋Δｔにおける位置ｘでの前記粒子の質量、運動量、及び温度を計算する段階と、
計算された温度、質量、及び運動量を用いてｇi（ｘ，ｔ＋Δｔ）を決定する段階と、
差分ｇi（ｘ，ｔ＋Δｔ）−ｇi（ｘ，ｔ）を運動状態ｆi（ｘ，ｔ）＋Ｃi（ｘ，ｔ）に加える段階と、
を含む、請求項２３に記載のシステム。
【請求項２７】
ｇiは、
2019053752.tif 000143
に従って定義され、
ρは流体密度であり、
Ｔ0は一定格子温度であり、
Ｐは前記流体容積中の圧力であり、
Ｔは前記計算された温度であり、
ｗiは一定の重み因子である、
請求項２６に記載のシステム。
【請求項２８】
前記加える段階が為された後の分布関数は、
2019053752.tif 000144
であり、
ｘは前記容積中の特定位置であり、
ｔは特定の第１の時点であり、
ｉは前記格子速度セットでのインデックスであり、
ｃiは衝突前の前記粒子の速度ベクトルであり、
Ｃiは衝突演算子であり、
Δｔは前記第１の時点と第２の時点との間の間隔であり、
ｇiは熱力学的ステップであり、
ｆi（ｘ，ｔ）は時間ｔにおける位置ｘでの前記粒子に対する分布関数である、
請求項２６に記載のシステム。
【請求項２９】
前記方法は、静止状態を
2019053752.tif 000145
に一致するように修正することによって、前記質量を保存する段階をさらに含む、請求項２６に記載のシステム。
【請求項３０】
前記粒子衝突ステップは、等温平衡分布関数を含む、請求項２４に記載のシステム。
【請求項３１】
前記分布関数は、時間ｔ＋Δｔにおける位置ｘ＋ｃiΔｔでの前記粒子に対する分布関数であり、
2019053752.tif 000146
で表され、
ｘは前記容積中の特定位置であり、
ｔは特定の第１の時点であり、
ｉは前記格子速度セットでのインデックスであり、
Ｃiは衝突演算子であり、
ｃiは衝突前の前記粒子の速度ベクトルであり、
Δｔは前記第１の時点と第２の時点との間の間隔であり、
ｇiは熱力学的ステップであり、
ｆi（ｘ，ｔ）は時間ｔにおける位置ｘでの前記粒子に対する分布関数である、
請求項２３に記載のシステム。
【請求項３２】
ｇiは、
2019053752.tif 000147
に従って定義され、
ρは流体密度であり、
Ｔ0は一定格子温度であり、
Ｐは前記流体容積中の圧力であり、
Ｔは流体の実際の温度であり、
ｗiは一定の重み因子である、
請求項３１に記載のシステム。
【請求項３３】
前記格子速度セットは、格子ボルツマン法に基づく、請求項３２に記載のシステム。
【外国語明細書】
2019053752000001.pdf