土木学会論文集 B2（海岸工学） Vol. 69，No. 2，2013，I_546-I_550
3次元LBM数値移動床の開発と掃流粒子群の挙動シミュレーション Modeling 3D Numerical Movable Bed based on LBM and its Application to the Analysis of Bed-Load Layer １
２
３
賀 露露 ・岡島雅史 ・喜岡 渉 ・北野利一
４
Lulu HE, Masashi OKAJIMA, Wataru KIOKA and Toshikazu KITANO The three-dimensional lattice Boltzmann method (LBM) is used to simulate the behavior of bed-load layer in open channel shear flow. The 3D LBM for particle-fluid suspensions originally developed by Ladd (1994) has been extended to free-surface flow and to higher Reynolds number regime. The SGS turbulence term based on Smagorinsky model is also included in the LBM flow simulation. The present simulation predicts well the behavior of each particle in the bed-load layer both in saltation and sheetflow regimes. Results of the present simulation for the height and length of saltation, the bed-load velocity and the existing probability density agree fairly well with the experimental data from Gotoh et al. (1994).
1. はじめに 流砂・漂砂現象を砂（粒子）の輸送問題としてとらえ， 流体・粒子間の相互干渉を正確に理解するために，実 験・数値解析の両面から研究が進められている．しかし，
一粒径の掃流砂の運動形態と底質に働く掃流力の関係を 定量的に調べることによった．
2. Ladd の LBM モデルとその拡張 （1）基本アルゴリズム
流体中の砂粒子間相互の干渉については，数値解析に複
Ladd（1994a, b）の LBM モデルに従って，格子上の仮
雑かつ膨大な計算を要する等の理由から，比較的研究例
想粒子分布の時間発展により流動を表現し，流体中を運
が少ない（例えば Feng ・ Michaelides，2003）．粒子・流
動する砂粒子（移動固体粒子）に link-bounce-back 法を適
体間の相互作用系としてのシートフロー漂砂の数値解析
用することによって，流体と移動固体粒子の相互作用力
は，従来より固液二相流系の観点から数多く行われてい
を運動量の交換に基づき計算する．巨視的な変数として
るものの，固液間や粒子間の 3 次元的相互作用を考慮し
密度 ρ ，流速 u に加えてモーメントフラックス ∏ を速度
た研究例は非常に少なく，3 次元個別要素法に基づく粒
分布関数から求める必要があることから，図-1 に示す
状体モデルを用いた後藤ら（2002），原田・後藤（2006）
3 次元 19 速度格子モデルを用いる．粒子の速度ベクトル
など限られたものになっている．2 次元モデルでは全て
e（i=0, 2, ..., 18）は，それぞれ0（i=0），e（i=1, 2, 3, 4, 5, 6）， i
の粒子が同一面内で接触するために運動の自由度が低下
√2e（i=7, 8, …, 18）で与えられる．ここで，e= ∆ x/ ∆ t
し，粒子運動による抵抗が相対的に増大するため，同一
（∆ x, ∆ tはそれぞれ空間，時間の格子間隔）である．
平面外での接触を表現できる 3 次元モデルによる個々粒
式（1）の格子ボルツマン方程式は，仮想粒子の並進
子の移動抵抗の評価が，特に混合砂における粒径階間の
と衝突の 2 つの過程を表しており，並進過程において粒
干渉を適切に評価する上で，不可欠とされている（後藤
子は速度に応じた方向の隣接した格子点へと移動する．
ら，2002） ．
衝突過程においては，粒子分布が割合 ∆ i で局所平衡状態
本研究は，Ladd（1994a, b）の固体粒子間空隙におけ る流体の圧縮性を考慮した格子ボルツマン法（Lattice Boltzmann Method，以下 LBM）をベースとして，自由表 面流下における掃流過程の 3 次元シミュレーションモデ ルに拡張し，3 次元数値移動床としての適用性を検討す るものである．モデルの検証は，開水路流れにおける単 １ 学生会員 ２ ３ フェロー ４ 正会員
名古屋工業大学大学院工学研究科社会 工学専攻 修（工） 中日本高速道路株式会社 Ph.D. 名古屋工業大学大学院教授 博（工） 名古屋工業大学大学院准教授 M.S.
図-1
格子形状
I_547
3 次元LBM数値移動床の開発と掃流粒子群の挙動シミュレーション
へ再配分される． …（1） f（x,t）は粒子分布関数であり，時間 t において，位置 x i の格子点上に存在する粒子の i 方向の速度ベクトルを持 つ粒子の数を示している．粒子分布関数を用いて，流体 場の密度 ρ ，流速 u，モーメントフラックス ∏ および局 所平衡分布関数 fieq はそれぞれ以下のように表される．
図-2
固体粒子表面境界
xb は境界点の位置 xb=x+1/2eb∆ t である．この bounce-back ………………………………………………………（2）
法では，固液間での流体の移流を許容していため，固体 粒子に作用する流体力が厳密に算定できない．そこで， Aidun ら（1998）が提案した固体粒子内部への流体の移 流を許さないアルゴリズムを採用して，固液界面のみで
…（3）
運動量交換が生じるようにした． 固体粒子に働く力 F とトルク T は，固体粒子の速度 U， 次ステップの粒子分布関数 fi*=fi+∆ i は，
，
p=ρ cs2 として．次式のように書くことができる．
角速度 Ω に依存する項と独立した項の和として半ステッ プ毎に求める．粒子速度や角速度に依存しない項につい てのみ示すと次のようである． ……… （10）
…………………………………………………（4） ここで， … （11） ………………（5）
求められた F，T を用いて，次式より固体粒子の速度 U， 角速度 Ω を算定し，粒子境界点を ∆ t後へ移動させる．
はそれぞれ， …（6）
…………………… （12）
……………………………（7）
…………………… （13）
添え字 i は速度ベクトル u の識別子で，1 は単位テンソル である．なお，式（6）の λ ， λ v はそれぞれせん断粘性
ここに，m は質量，Iは慣性モーメントである． （2）乱流項の追加 式（1）の BGK 近似の衝突項に ∆ i おける λ に，従来の
係数と体積粘性係数に関する固有値を示す． 固体粒子の表面は，格子点の間に新たな境界点を設置
粘性項に加えて Smagorinsky モデルによる SGS 乱流項を
することにより定義する（図-2）．粒子表面の境界条件は，
対応させた．すなわち，動粘性係数 ν に SGS 渦動粘性係
link-bounce-back 条件に従って与え，粒子表面外側の格子
数 vt を加えて vtotal = v + vt と表すと，λ は次式のように与え
点 b について，粒子表面から流入する粒子分布関数を以
られる．
下のように与える． ………………………… （14） …（8） SGS 渦動粘性係数は，格子スケールを ∆ x とすると歪速度 ただし，fb*（x,t）は位置 x，時間 t，における eb 方向の衝 ~ 突後の粒子分布関数であり，b は b と逆方向，すなわち eb~ =−eb である．ub は，質量の中心R，固体粒子速度U，角 速度 Ω を用いて次式で表される．
テンソルにSより次式で表される． ……………………………… （15） ここに，C は Smagorinsky 定数である．数値シミュレーシ ョンでは，乱流モデルにおけるフィルター長を 2∆ x 以下
………………………（9）
と仮定して，C=0.1 とした．歪速度︱S︱は LBM モデルで
I_548
土木学会論文集B2（海岸工学），Vol. 69，No. 2，2013
は次のように求めることができる． ……… （16） ただし，Qij は次式で与えられる．
………………… （17） （3）潤滑近似の導入 個々の粒子が運動により重なり合うのを防ぐために， 粒子同士が接近しようとするときに粒子間で働く斥力を Ngueyen ・ Ladd（2002）が提案した潤滑近似のアルゴリ ズムに従って算定した．なお，2 粒子間相互の潤滑近似 の計算におけるカットオフ距離 hN は，hN /∆ x = 0.2とした．
図-3
（4）自由表面流への拡張 この Ladd の LBM モデルを自由表面流に適用するため
自由表面・流速ベクトル・粒子のスナップショット （τ *=0.16）
件を課して掃流砂計算を行った．底面に接した砂粒子は，
に，巨視的変数に流体充填率 ε を加え，00.2 では τ * によらずほぼ固定されている．
τ *=0.23 について後藤ら（1993）の実験値と比較すると， 計算値のピーク値の方が若干大きく，ピークの位置もや や高くなっているものの，鉛直方向存在確率密度の重心 については実験結果とよく一致している．
5. おわりに 3 次元 LBM モデルを自由表面流下の掃流粒子群の数値 図-7
砂粒子の水平方向移動速度
シミュレーションモデルに拡張し，小規模移動床実験条 件での従来の実験や数値解析の結果と比較することによ り，その適用性を確認した．この LBMモデルにおいては， 実験値等を参考に半経験的に決める係数値を基本的に入 力値に必要としておらず，数値計算上の閾値（潤滑近似 のカットオフ距離 h N，セル充填率の閾値 κ ）を除けば， 乱流モデルにおける Smagorinsky 定数についてのみであ る．本数値モデルは並列計算にきわめて適したアルゴリ ズムであるため，大型模型実験規模の 3 次元数値移動床 への拡張も比較的容易である．
図-8
砂粒子の鉛直方向存在確率密度分布
次元混相流シミュレーションの結果に基づき提案した近 似式 Ls/d=8.32ln（τ *）+25.0，Hs/d=0.77ln（τ *）+2.74 とよく 一致している． 砂粒子の水平移動速度分布の τ * に対する変化につい て，τ *=0.23 における後藤ら（1994）の実験結果とともに 図-7 に示す．粒子間衝突がほとんど発生しない τ *=0.06の ケースを除いて，掃流層の上層部では単調増加分布とな っておらず，複数の変曲点を有する複雑な分布形状とな っている．この速度分布の蛇行域は，後藤ら（1993）が 示したように，中層部の砂粒子群との衝突を繰り返しな がら円弧状に大きく跳躍するサルテーション運動領域に 対応している．実験値と比較すると特に中層から底層の シートフロー領域で水平移動速度を小さく予測している が，底層部の砂粒子速度の欠損部分から表層付近のサル テーション移動速度に至るまで比較的よく再現されて いる． 次いで，掃流層 z の位置での移動砂粒子数をカウント し，全移動砂粒子数で除して鉛直方向の存在確率密度 P.D.F.を求めた．図-8 は，水路中央部の掃流層各層にお ける P.D.F.のシールズ数 τ * に対する変化の様子を示した ものである．τ * の増加に伴い，ピークの幅はやや広くな り，底層からピークに向かって急峻な立ち上がりと，上
参 考 文 献 後藤仁志・原田英治・酒井哲郎（2002）： 3 次元数値移動床に よるシートフロー層の鉛直分級過程の数値解析，海岸工 学論文集，第49 巻，pp.471-475. 後藤仁志・辻本哲郎・中川博次（1994）：流体・粒子相互作用 系としての掃流層の数値解析，土木学会論文集，第 485 号 /П-26，pp11-19. 後藤仁志・辻本哲郎・中川博次（1993）： Saltation からシート フローへの遷移過程における砂粒子群の運動特性，海岸 工学論文集，第40 巻，pp.326-330. 原田英治・後藤仁志（2006）：三次元数値移動床による混合粒 径シートフロー漂砂の分級過程の解析，土木学会論文集B， Vol.62, No.1, pp.128-138. Aidun, C. K., Y. Lu and E.-J. Ding (1998): Direct analysis of particulate suspensions with inertia using the discrete Boltzmann equation, J. Fluid Mech., 373, pp.287-311. Feng Z.-G. and E. E. Michaelides (2003): Fluid-particle interactions and resuspension in simple shear flow, J. Hydraul. Eng., 129, pp.985-994. Körner, C., T. Pohl, U. Rüde, N. Thürey and T. Zeiser (2005) : Parallel lattice Boltzmann methods for CFD applications, Numerical Solution of Partial Differential Equations on Parallel Computers, 51 of LNCSE, pp.439-465. Ladd, A. J. C. (1994a) : Numerical simulations of particulate suspensions via a discretized Boltzmann equation. Part 1. Theoretical foundation, J. Fluid Mech., 271, pp.285-309. Ladd, A. J. C. (1994b) : Numerical simulations of particulate suspensions via a discretized Boltzmann equation. Part 2. Numerical results, J. Fluid Mech., 271, pp. 311-339. Nguyen, N.-Q. and A. J. C. Ladd (2002) : Lubrication corrections for lattice-Boltzmann simulations of particle suspensions, Physical Review E, 66, 046708.
