(57)【要約】  流体、及び流体中に浮遊する複数の粒子を含むシミュレーション対象の流動層に含まれる粒子を仮想的に拡大し、個数を減らす粗視化を行う。粗視化の前後で、流動層の流れに関する無次元量及び熱輸送に関する無次元量が変化しない条件の下で、粒子及び流体の物性値、粒子及び流体に定義される物理量を変換する。変換後の物性値及び物理量を用いて流動層の挙動をシミュレーションする。これにより、粒子数が増えても計算負荷の増大を抑制することが可能になる。
2019181541.tif 000002
【特許請求の範囲】
【請求項１】
流体、及び前記流体中の複数の粒子を含むシミュレーション対象の流動層に含まれる前記粒子を仮想的に拡大し、個数を減らす粗視化を行い、
粗視化の前後で、流動層の流れに関する無次元量及び熱輸送に関する無次元量が変化しない条件の下で、前記粒子及び前記流体の物性値、前記粒子及び前記流体に定義される物理量を変換して、
変換後の前記物性値及び前記物理量を用いて流動層の挙動をシミュレーションするシミュレーション方法。
【請求項２】
粗視化の前後で変化しない流れに関する無次元量は、粒子レイノルズ数、アルキメデス数、及びフルード数である請求項１に記載のシミュレーション方法。
【請求項３】
粗視化の前後で変化しない熱輸送に関する無次元量は、プラントル数、粒子ヌセルト数、及びビオ数である請求項１または２に記載のシミュレーション方法。
【請求項４】
粗視化の前後において、前記流体の温度及び前記粒子の温度を不変とする条件の下で前記物性値及び前記物理量を変換する請求項１または２に記載のシミュレーション方法。
【請求項５】
粗視化の前後において、前記粒子と前記流体との間の熱伝達係数を不変とする条件の下で前記物性値及び前記物理量を変換する請求項１または２に記載のシミュレーション方法。
【請求項６】
粗視化の前後において、全体の前記粒子の顕熱を不変とする条件の下で前記物性値及び前記物理量を変換する請求項１または２に記載のシミュレーション方法。
【請求項７】
流体、及び前記流体中の複数の粒子を含むシミュレーション対象の流動層の前記流体及び前記粒子の物性値、前記流体及び前記粒子に定義される物理量の初期条件を取得するシミュレーション条件取得部と、
前記粒子を拡大する拡大率を取得する拡大率取得部と、
前記粒子を拡大しても流動層の流れに関する無次元量及び熱輸送に関する無次元量が変化しない条件の下で、前記シミュレーション条件取得部で取得された前記物理量の初期条件及び前記物性値を変換して、変換後の前記物性値及び前記物理量を用いて流動層の挙動をシミュレーションする演算部と
を有するシミュレーション装置。
【請求項８】
流体、及び前記流体中の複数の粒子を含むシミュレーション対象の流動層の前記流体及び前記粒子の物性値、前記流体及び前記粒子に定義される物理量の初期条件を取得する機能と、
前記粒子を拡大する拡大率を取得する機能と、
前記粒子を拡大しても流動層の流れに関する無次元量及び熱輸送に関する無次元量が変化しない条件の下で、取得された前記物理量の初期条件及び前記物性値を変換して、変換後の前記物性値及び前記物理量を用いて流動層の挙動をシミュレーションする機能と
をコンピュータに実現させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、シミュレーション方法、シミュレーション装置、及びプログラムに関する。
【背景技術】
【０００２】
粒子の挙動を解析する離散要素法（ＤＥＭ）と流体の流れ場を解析する数値流体力学（ＣＦＤ）とを連成させることで、固体粒子を流体中に浮遊させた状態の流動層の挙動を解析する手法が公知である（非特許文献１、２）。非特許文献１には、粒子数が増えたときの計算時間の増大を抑制するシミュレーション方法が提案されている。具体的には、粒子を拡大して粒子数を減らす処理（粗視化）を行い、粗視化の前後で支配方程式が同じになるように物性値や物理量を変換し、粗視化後の流動層についてシミュレーションを実行する。非特許文献２には、流動層内の熱輸送を評価する手法が提案されている。
【先行技術文献】
【非特許文献】
【０００３】
【非特許文献１】鷲野  公彰、許  志宏、川口  寿裕、辻  裕、「流動層のＤＥＭ計算における相似則モデル」、粉体工学会誌  第４４巻第３号、２００７年、第１９８頁〜第２０５頁
【非特許文献２】A.V. Patil, E.A.J.F. Peters, J.A.M. Kuipers, "Comparison of CFD-DEM heat transfer simulations with infrared/visual measurements", Chemical Engineering Journal 277 (2015) 388-401
【発明の概要】
【発明が解決しようとする課題】
【０００４】
非特許文献１には、熱輸送に関するパラメータについて何ら述べられていない。すなわち、非特許文献１に記載された方法は、冷間状態（温度変化せず、通常は室温状態）での流動層の挙動のシミュレーションに適用できるが、熱輸送が生じ得る熱間状態の流動層のシミュレーションに適用することはできない。熱間状態の流動層のシミュレーションに非特許文献２に記載の方法を適用すると、粒子数の増大に伴って計算負荷が大きくなってしまう。【０００５】
本発明の目的は、熱輸送が生じ得る流動層のシミュレーションにおいて、粒子数が増えても計算負荷の増大を抑制することが可能なシミュレーション方法、シミュレーション装置、及びプログラムを提供することである。
【課題を解決するための手段】
【０００６】
本発明の一観点によると、
流体、及び前記流体中の複数の粒子を含むシミュレーション対象の流動層に含まれる前記粒子を仮想的に拡大し、個数を減らす粗視化を行い、
粗視化の前後で、流動層の流れに関する無次元量及び熱輸送に関する無次元量が変化しない条件の下で、前記粒子及び前記流体の物性値、前記粒子及び前記流体に定義される物理量を変換して、
変換後の前記物性値及び前記物理量を用いて流動層の挙動をシミュレーションするシミュレーション方法が提供される。【０００７】
本発明の他の観点によると、
流体、及び前記流体中の複数の粒子を含むシミュレーション対象の流動層の前記流体及び前記粒子の物性値、前記流体及び前記粒子に定義される物理量の初期条件を取得するシミュレーション条件取得部と、
前記粒子を拡大する拡大率を取得する拡大率取得部と、
前記粒子を拡大しても流動層の流れに関する無次元量及び熱輸送に関する無次元量が変化しない条件の下で、前記シミュレーション条件取得部で取得された前記物理量の初期条件及び前記物性値を変換して、変換後の前記物性値及び前記物理量を用いて流動層の挙動をシミュレーションする演算部と
を有するシミュレーション装置が適用される。【０００８】
本発明のさらに他の観点によると、
流体、及び前記流体中の複数の粒子を含むシミュレーション対象の流動層の前記流体及び前記粒子の物性値、前記流体及び前記粒子に定義される物理量の初期条件を取得する機能と、
前記粒子を拡大する拡大率を取得する機能と、
前記粒子を拡大しても流動層の流れに関する無次元量及び熱輸送に関する無次元量が変化しない条件の下で、取得された前記物理量の初期条件及び前記物性値を変換して、変換後の前記物性値及び前記物理量を用いて流動層の挙動をシミュレーションする機能と
をコンピュータに実現させるためのプログラムが提供される。
【発明の効果】
【０００９】
粒子を粗視化して個数を減らすことにより、計算負荷を軽減させることができる。粗視化後の流動層の流れ及び熱輸送に関してシミュレーションを行った結果は、粗視化前の流動層の流れ及び熱輸送の状況を反映したものになる。このため、粗視化前の流動層の挙動を予測することができる。
【図面の簡単な説明】
【００１０】
【図１】図１Ａは、シミュレーション対象の流動層の一例を示す模式図であり、図１Ｂは、シミュレーション対象の粗視化後の流動層の一例を示す模式図である。
【図２】図２は、粒子及びガスの物性値、粒子及びガスに関して定義される種々の物理量について、本明細書で用いる記号及び粗視化の係数の一覧を示す図表である。
【図３】図３は、本実施例によるシミュレーション装置のブロック図である。
【図４】図４は、本実施例によるシミュレーション方法のフローチャートである。
【図５】図５は、実際に行ったシミュレーションのシミュレーション領域を示す斜視図である。
【図６】図６は、粗視化した流動層のシミュレーションにより求めた粗視化粒子の位置及び温度を、時系列で示す図である。
【図７】図７Ａ及び図７Ｂは、シミュレーション結果から求めた粒子の平均温度の時間変化を示すグラフである。
【図８】図８は、他の実施例によるシミュレーション方法で適用される変換則を示す図表である。
【発明を実施するための形態】
【００１１】
図１Ａ〜図７Ｂを参照して、実施例によるシミュレーション方法及び装置について説明する。
図１Ａは、シミュレーション対象の流動層の一例を示す模式図である。シミュレーション対象の領域１０内に複数の粒子１１を配置し、下方から上方に向かって領域１０内にガス１２を導入することにより形成される流動層の挙動をシミュレーションする。粒子１１の直径をＤｐ１で表す。本実施例では、粒子１１の各々を拡大して、その個数を減らす（以下、粗視化という。）ことにより、計算負荷を軽減させる。
【００１２】
図１Ｂは、シミュレーション対象の粗視化後の流動層の一例を示す模式図である。粒子１１を拡大して仮想的な粒子２１を得る。仮想的な粒子２１は、シミュレーション対象の領域２０内に配置される。粗視化後の領域２０の寸法は、粗視化前の領域１０の領域の寸法と同一である。仮想的な粒子２１の直径をＤｐ２で表す。拡大率Ｋを、粗視化前の粒子１１の直径に対する粗視化後の仮想的な粒子２１の直径の比と定義する。拡大率Ｋは以下の式で定義される。
【数１】
2019181541.tif 000003
【００１３】
粗視化後の粒子２１が配置された領域２０内に下方から上方に向かってガス２２を導入することにより形成される粗視化後の流動層について、数値流体力学（ＣＦＤ）と離散要素法（ＤＥＭ）とを連成させた解析を行う。粗視化後の仮想的な流動層と、粗視化前の実際の流動層とが相似則を満たすように、粗視化に際し、粒子１１及びガス１２の物性値及び種々の物理量を変換する。
【００１４】
次に、図２を参照して、粒子１１及びガス１２の物性値及び種々の物理量の変換則について説明する。
【００１５】
図２は、粒子及びガスの物性値、粒子及びガスに関して定義される種々の物理量について、本明細書で用いる記号及び粗視化の係数の一覧を示す図表である。粗視化前の実際の物性値及び物理量に粗視化の係数を乗ずることにより、粗視化後の流動層に関する物性値及び物理量が得られる。本明細書において、例えば式（１）に示したように、粗視化前の物性値及び物理量を表す記号には、下付きの添え字「１」を付し、粗視化後の物性値及び物理量を表す記号には、下付きの添え字「２」を付す。
【００１６】
流動層の流れに関する無次元量として、粒子レイノルズ数Ｒｅｐ、アルキメデス数Ａｒｐ、及びフルード数Ｆｒが挙げられる。これらの無次元量は以下の式で定義される。
【数２】
2019181541.tif 000004
ここで、ｇは重力加速度である。太字のＶ及びＵは、ベクトルであることを意味する。ボイド率εは、充填された粒子の全質量をＭ、粒子が充填された領域の見かけの体積をＶＡとして、以下の式で定義される。
【数３】
2019181541.tif 000005
【００１７】
粗視化の前後で、流動層の流れに関する無次元量である粒子レイノルズ数Ｒｅｐ、アルキメデス数Ａｒｐ、及びフルード数Ｆｒが変化しないという条件を設定する。さらに、ボイド率εが変化しないという条件、及びガス粘性係数μが変化しないという条件の下で、粗視化前後の物性値及び物理量の変換則を求めると、以下の変換則が得られる。
【数４】
2019181541.tif 000006
【００１８】
ガス密度ρｆ２の変換則から、ガス圧力ｐについて、以下の変換則が得られる。
【数５】
2019181541.tif 000007
【００１９】
粗視化の前後で粒子が充填された領域の見かけの体積ＶＡが変化せず、粒子の個数が粗視化によって１／Ｋ３に減少すると仮定すると、以下の変換則が得られる。
【数６】
2019181541.tif 000008
【００２０】
粒子質量流量ｍｐドットは、流路面積をＡとして以下の式で定義される。
【数７】
2019181541.tif 000009
この式から、以下の変換則が導出される。
【数８】
2019181541.tif 000010
【００２１】
さらに、熱輸送に関する無次元量についても、粗視化の前後で変化しないという条件を付す。熱輸送に関する無次元量として、プラントル数Ｐｒ、粒子ヌセルト数Ｎｕｐ、ビオ数Ｂｉが挙げられる。プラントル数Ｐｒ、粒子ヌセルト数Ｎｕｐ、ビオ数Ｂｉは以下の式で定義される。
【数９】
2019181541.tif 000011
ここで、Ｌｐは粒子の特徴長さであり、Ｌｐ＝Ｄｐ／６で定義することができる。
【００２２】
物性値の温度依存性を簡単化するために、粗視化の前後で粒子温度Ｔｐ及びガス温度Ｔが変化しないと仮定する。さらに、粒子熱伝達係数ｈも、粗視化の前後で変化しないと仮定する。この仮定の下で、以下の変換則が得られる。
【数１０】
2019181541.tif 000012
【００２３】
上述の仮定のみでは、粒子比熱ｃの変換則が定まらない。本実施例では、粒子比熱ｃの変換則を決定するために、粗視化の前後で粒子全体の顕熱Ｑｐ，ａｌｌが変化しないという仮定を導入する。粒子全体の顕熱Ｑｐ，ａｌｌは、粒子の個数をＮｐ、流動層に導入するガス温度Ｔに対する粒子の初期温度の差をΔＴｐとして、以下の式で定義される。
【数１１】
2019181541.tif 000013
【００２４】
粒子の個数Ｎｐは、粗視化によって約１／Ｋ３に減少するため、粒子全体の顕熱Ｑｐ，ａｌｌが粗視化の前後で不変と仮定すると、以下の変換則が得られる。
【数１２】
2019181541.tif 000014
【００２５】
粒子表面の伝熱量Ｑドットは、以下の式で定義される。
【数１３】
2019181541.tif 000015
この定義から、伝熱量Ｑドットについて以下の変換則が得られる。
【数１４】
2019181541.tif 000016
粒子表面の熱流束ｑドットについては、以下の変換則が得られる。
【数１５】
2019181541.tif 000017
【００２６】
図３は、本実施例によるシミュレーション装置のブロック図である。本実施例によるシミュレーション装置は、処理装置３０、入力装置３８、及び出力装置３９を含む。処理装置３０は、シミュレーション条件取得部３１、拡大率取得部３２、演算部３３、及び出力制御部３４を含む。
【００２７】
図３に示す各ブロックは、ハードウェア的には、コンピュータの中央処理ユニット（ＣＰＵ）をはじめとする素子や機械装置で実現することができ、ソフトウェア的にはコンピュータプログラム等によって実現することができる。図３では、ハードウェア及びソフトウェアの連携によって実現される機能ブロックが示されている。従って、これらの機能ブロックは、ハードウェア及びソフトウェアの組み合わせによって、種々の態様で実現することが可能である。
【００２８】
処理装置３０は入力装置３８及び出力装置３９と接続される。入力装置３８は、処理装置３０で実行される処理に関係するユーザからのコマンド及びデータの入力を受ける。入力装置３８として、例えばユーザが操作を行うことにより入力を行うキーボードやマウス、インターネット等のネットワークを介して入力を行う通信装置、ＣＤ、ＤＶＤ等の記録媒体から入力を行う読取装置等を用いることができる。
【００２９】
シミュレーション条件取得部３１は、入力装置３８を介してシミュレーション条件を取得する。シミュレーション条件には、シミュレーションに必要な種々の情報が含まれる。例えば、シミュレーション対象の粒子及びガスの物性値、粒子及びガスに関する物理量の初期条件、境界条件等が含まれる。拡大率取得部３２は、入力装置３８を介して拡大率Ｋ（図２）を取得する。
【００３０】
演算部３３は、シミュレーション条件及び拡大率Ｋに基づいて、粗視化前の物性値及び物理量に粗視化の係数（図２）を乗じることにより、粒子及びガスの粗視化後の物性値及び物理量の初期条件を算出する。粗視化後の物性値及び物理量の初期条件に基づき、ＣＦＤとＤＥＭとを連成させた流動層のシミュレーションを行う。
【００３１】
出力制御部３４は、シミュレーション結果を出力装置３９に出力する。例えば、粒子の位置及び温度の変動、ガスの温度分布の変動を、出力装置３９の表示画面に図形で表示する。
【００３２】
図４は、本実施例によるシミュレーション方法のフローチャートである。まず、シミュレーション条件取得部３１（図３）がシミュレーション条件を取得し（ステップＳ１）、拡大率取得部３２（図３）が拡大率Ｋ（図２）を取得する（ステップＳ２）。
【００３３】
その後、演算部３３（図３）は、シミュレーション条件として入力された物性値及び物理量の初期値を、粗視化後の値に変換する（ステップＳ３）。さらに、変換後の物性値及び物理量に基づいてシミュレーションを実行する（ステップＳ４）。シミュレーションが終了すると、出力制御部３４（図３）がシミュレーション結果を出力する（ステップＳ５）。
【００３４】
次に、図５〜図７Ｂを参照して、本実施例によるシミュレーション方法を用いて実際にシミュレーションを行った結果について説明する。このシミュレーションの対象は、非特許文献２に記載されているものと同一である。
【００３５】
図５は、シミュレーション領域４０を示す斜視図である。シミュレーション領域４０は、幅８ｃｍ、厚さ１．５ｃｍ、高さ２５ｃｍの直方体である。シミュレーション領域４０内に、直径１ｍｍの複数のガラス粒子を充填し、シミュレーション領域４０の底面からシミュレーション領域４０内にガスを導入する。粒子密度ρｐを２５００ｋｇ／ｍ３とした。粒子比熱ｃを８４０Ｊ／ｋｇ／Ｋとし、ガス定圧比熱ｃｐ，ｆを１０１０Ｊ／ｋｇ／Ｋとし、ガス粘性係数μを２．０×１０−５Ｐａ・ｓとした。シミュレーション領域４０内に充填する粒子の質量の合計を７５ｇとした。粒子の初期温度よりも低温のガスをシミュレーション領域４０内に導入した。ガスの流速を１．２０ｍ／ｓとした場合（流速が遅い場合）と、１．５４ｍ／ｓとした場合（流速が速い場合）とについてシミュレーションを行った。
【００３６】
拡大率Ｋを２にして粗視化した流動層と、元の流動層との２つについてシミュレーションを行った。
【００３７】
図６は、粗視化した流動層のシミュレーションにより求めた粗視化粒子の位置及び温度を、時系列で示す図である。図６の左から１番目、２番目、３番目、及び４番目の図は、それぞれ冷却開始時点、冷却開始からの経過時間がｔ、２ｔ、及び３ｔの流動層の状態を示す。各粒子の濃さは粒子の温度を表しており、温度が高いほど濃く表されている。ガスの流入によって粒子が流動し、時間の経過とともに粒子の温度が低下していることがわかる。
【００３８】
図７Ａ及び図７Ｂは、シミュレーション結果から求めた粒子の平均温度の時間変化を示すグラフである。横軸は冷却開始からの経過時間を任意単位で表し、縦軸は粒子の平均温度を、初期温度を基準とした相対値で表す。図７Ａは、ガス流速が遅い場合を示し、図７Ｂは、ガス流速が速い場合を示している。グラフ中の破線は粗視化前の流動層のシミュレーション結果を示し、実線は粗視化後の流動層のシミュレーション結果を示している。参考のために、非特許文献２に示された実験結果による粒子の温度変化を丸記号で示している。
【００３９】
図７Ａ及び図７Ｂに示したシミュレーション結果から、本実施例による方法で粗視化してシミュレーションを行っても、シミュレーション結果は実験結果とよく一致していることが確認できる。ガスの流速を速くすると、粒子の温度低下が速くなることも確認できる。このように、本実施例による粗視化の手法は、温度変化を伴う流動層の挙動のシミュレーションに適用することが可能である。
【００４０】
粗視化を行うことにより、シミュレーションに要する計算時間は、粗視化前の流動層のシミュレーションと比べて約１／３になった。このように、粗視化により計算負荷を低減させることができる。
【００４１】
次に、図８を参照して他の実施例による物性値及び物理量の変換則について説明する。
図８は、本実施例によるシミュレーション方法で適用される変換則を示す図表である。以下、図２に示した変換則と比較しながら説明する。
【００４２】
流動層の流れに関する無次元量、及び熱輸送に関する無次元量を、粗視化の前後で変化させない点は、図２に示した実施例の場合と同様である。粒子温度Ｔｐ、ガス温度Ｔ、粒子熱伝達係数ｈが、粗視化の前後で変化しない点についても、図２に示した実施例の場合と同様である。
【００４３】
図２に示した実施例では、粗視化の前後でガス粘性係数μが変化しないと仮定したが、本実施例では、粗視化の前後で粒子密度ρｐ及びガス密度ρｆが変化しないと仮定する。この仮定の下で、粒子全体の顕熱Ｑｐ，ａｌｌが粗視化の前後で変化しないと仮定すると、粒子比熱ｃも、粗視化の前後で変化しない。さらに、ガス圧力ｐも、粗視化の前後で変化しない。
【００４４】
図８に示した実施例では、粒子質量ｍｐ、ガス粘性係数μ、粒子比熱ｃ、ガス定圧比熱ｃｐ，ｆ、粒子質量流量ｍｐドットの変換則が、図２に示した変換則と異なっている。図８に示した変換則を用いて粒子及びガスの物性値と物理量を変換してシミュレーションを行ってもよい。
【００４５】
上述の各実施例は例示であり、異なる実施例で示した構成の部分的な置換または組み合わせが可能であることは言うまでもない。複数の実施例の同様の構成による同様の作用効果については実施例ごとには逐次言及しない。さらに、本発明は上述の実施例に制限されるものではない。例えば、種々の変更、改良、組み合わせ等が可能なことは当業者に自明であろう。
【符号の説明】
【００４６】
１０  シミュレーション対象領域
１１  粒子
１２  流体
２０  粗視化後のシミュレーション対象領域
２１  粗視化後の粒子
２２  粗視化後の流体
３０  処理装置
３１  シミュレーション条件取得部
３２  拡大率取得部
３３  演算部
３４  出力制御部
３８  入力装置
３９  出力装置
４０  シミュレーション領域
【図１】
2019181541.tif 000018
【図２】
2019181541.tif 000019
【図３】
2019181541.tif 000020
【図４】
2019181541.tif 000021
【図５】
2019181541.tif 000022
【図６】
2019181541.tif 000023
【図７】
2019181541.tif 000024
【図８】
2019181541.tif 000025
【国際調査報告】
2019181541.tif 000026
2019181541.tif 000027
2019181541.tif 000028
2019181541.tif 000029