粒子法による流体剛体連成解析法による橋梁流出シミュレーションと V&V
Verification and Validation of bridge wash-out simulation
by a particle based fluid rigid body interaction method
大屋朋子・九州大学大学院 工学府建設システム工学専攻
Tomoko OYA, Kyushu Univ.

浅井光輝・九州大学大学院准教授
工学研究院社会基盤部門

伊津野和行・立命館大学教授
理工学部都市システム工学科

磯部大吾郎・筑波大学教授
システム情報系構造エネルギー工学域

Mitsuteru ASAI, Kyushu Univ.

Kazuyuki IZUNO, Ritsumeikan Univ.

Daigoro ISOBE, Tsukuba Univ.

Key Words: SPH , particle method , Tsunami , fluid-rigid interaction
2011 年 3 月 11 日に発生した東北地方太平洋沖地震に伴う津波によって東北地方の太平洋沿岸部は甚大な被害を受けた．中でも橋梁に関しては，上部構造の流失被害が各地で発生し，震災後の支援・復旧活動に大きな支障をきたした．現在，今後危惧される巨大地震に伴う津波に備えて橋梁の津波対策が積極的に議論されており，実寸大での橋梁流失実験が現実的に困難であるため数値シミュレーションは橋梁の津波対策の検討に有効な手段となりえる．本研究では，数値解析手法の粒子法の一つである SPH 法を用い数値解析を行った．また，数値シミュレーションによる橋梁被害の事前予測と対策検討に向けて，事前段階として津波が橋梁に作用する流体力の精度検証を行った．また，粒子法による流体剛体連成解析を用い，橋梁上部構造の流失挙動解析の妥当性確認を実施した．本論文は，解析結果と水理実験を比較することで，流体剛体連成解析の精度と妥当性を確認した．
1. はじめに2011 年 3 月 11 日東北地方太平洋沖地震に伴い発生した津波により，2 万人近い死者・行方不明者が出ただけではなく，多くの社会インフラ構造物も津波により被害を受けた．イン
フラ構造物の中でも橋梁に関しては，宮城県仙台市から太平洋沿岸を上る国道 45 号に沿って，上部構造の流失が多数見られた．この流失被害により，被災地への救援活動の遅れや，被災地の復旧・復興にも遅れが生じた．この被害を受け，橋梁構造物の重要性が再確認され，早急な橋梁への津波対策が必要となってくる．そこで，本研究では数値シミュレーションを用い，津波が橋梁に作用する流体力や橋梁上部構造の流失の挙動を解析的に検討した．まず橋梁に作用する津波流体力評価に向けた小型模型実験の結果と解析結果の比較から，流体力に関する精度検証を行った．また，流体剛体連成解析において定式化はすでにされているものの，解析精度の検証は不明瞭であったため，本研究では模型実験の結果と解析結果から流体剛体連成解析の妥当性を確認した．また，解析手法は自由表面形状変化などの大変形解析に適している粒子法の一つである SPH (Smoothed Particle Hydrodynamics)法を採用した．
2. 解析手法
2-1 SPH 法の基礎式
SPH 法は，圧縮性流体の解析手法として提案されたのちに，非圧縮性流体解析法(ISPH 法)へと拡張された．連続体を有限個の粒子に離散化し，対象とする粒子の影響半径内に存在する近傍粒子に重み関数をかけて，一種の重み付き平均としてjj近似する手法である．粒子 j の質量 m ，密度 ρ とし，重み関数を W，影響半径を h とすると，物理量 f(x)は以下の式で表される．

N mjf ( x)  jj 1 ρf ( x)W ( x - x j , h)[1]

流体運動の支配方程式として，ナビエ・ストークス方程式を解く．本解析では，SPH 法の中でも圧力を陰的に，速度を陽的に解き，非圧縮性流体に適した解法とされる安定化ISPH 法 1）を採用した．
2-2 流体剛体連成解析
本研究では粒子法による流体剛体連成解析を用い数値解析を行った．流体剛体連成解析では，流体はナビエ・ストークス方程式，剛体は並進方向と回転方向の運動方程式を解くことで剛体粒子の速度や位置を更新する．以下に，並進方向と回転方向の運動方程式を示す．

dT F  Mg  F  F  M dt
fedω

 M  M  M  I dt  (ω  Iω)
fe[2][3]

ここで，M は剛体の質量，g は重力加速度，Ff は剛体に作用する流体力，Fe は剛体間の接触力，Mf は流体力によるモーメント，Me は接触力によるモーメント，I は慣性モーメントとする．
並進方向と回転方向の運動方程式から算出された並進速度
T，回転速度 ω を用い，剛体の位置と速度を更新する．以下に剛体粒子の位置，速度更新の式を示す．

xin1  rg n1  R n1ri0un1

Tn1ωn1n1 0iR r[4][5]

ここで， rg は剛体中心の移動を表すベクトル，R は重心周りの剛体の回転を表す回転行列とする．

日本機械学会〔№16-4〕第29回計算力学講演会CD-ROM論文集〔2016.9.22-24,名古屋〕

3. 解析手法の精度検証と妥当性確認
3-1 流体力の精度検証
橋梁に作用する流体力評価に向けた小型模型実験の結果と解析結果を比較することにより，本研究の解析手法における流体力の精度検証を行った．水路に矩形模型を設置し，ゲートを急開することにより，貯水部に溜めた水を放流し，模型に作用する流体力(模型に対する水平抗力，揚力)を計測する．
図-1 に抗力の実験結果と解析結果を示す．解析は実験と同様に，波が模型に衝突した際に抗力は最大となる．実験結果と解析結果を比較すると，解析は実験よりも 0.3 秒程遅く波が到達していた．また，抗力について実験と計算値を比較す
ると，実験値は解析値の 6 割程であった．実験計測では，ロードセルでの誤差，計測開始時点での時間計測の遅れが生じているのではないかとの憶測から，解析値を 0.3 秒早めに，抗力を 6 割に補正することで水平抗力は実験値と類似した傾向を得た．
図-2 には揚力の実験結果と解析結果を示す．実験と同様，初めに波が模型を下から上へと押し上げ，その後模型の上に乗った水により，下向きの力が作用する結果となった．なお
実験において，2.5～3.0 秒にかけて揚力が振動している要因
としては，片持ち状態で固定している模型が振動した結果であると推測している．(橋梁模型については，実験では片持ち梁的に模型を支持しているが，本研究では，小さな隙間をモデル化することが困難であることから，両端固定で模型をモデル化した)また，抗力と同様に解析は実験よりも 0.3 秒程遅く波が到達し，実験は解析の 6 割程の流体力しか受けていないことから，水平抗力時と同様な補正を加えた結果を同図に示した．抗力，揚力どちらにおいても，補正した解析結果と実験結果を比較すると，概ね近い結果となった．
3-2 剛体移動の妥当性確認
流体剛体連成解析における剛体移動の妥当性を確認するた
め，簡単な矩形モデルを用い，実験との比較検証を行った．
図-3 には立命館大学で実施した水理模型実験の概略図を示す．この実験では，ゲートを急開することで段波を発生させ，橋梁上部構造を模擬した模型に段波を衝突させ，この際の上部構造の運動を計測した．模型運動の計測には 3D モーションキャプチャーシステムを利用し，耐水機能のある自発光マーカーを模型上部に設置し，模型が流失する際の回転角度(奥行き方向)と流失の挙動を記録した(図-4)．なお，初期の貯水高さは 350mm とし，実験結果と比較を行った．
解析条件は粒子径を 0.5cm に設定したことで，総粒子数は約 200 万となった．また，時間増分は 0.0005sec で実時間は2.5 秒程度である．図-5 に数値解析と水理実験の回転角度の結果を示す．各実験結果にばらつきがあるため，4 回の実験結果と数値解析による計算結果を比較した．実験結果と同様に，数値解析の結果も津波が衝突した後，回転角度が正から負へと移行する傾向が見られた．

図-1 抗力図-2 揚力図-3 実験，解析の概略図

図-4 流失の挙動

4. おわりに本研究では，粒子法による流体剛体連成解析の精度検証と妥当性確認を行った．水理実験の結果と解析結果を比較することにより，流体力の精度や，剛体移動の妥当性を確認することができた．
参考文献
1) M. Asai, AM. Aly, Y. Sonoda and Y. Sakai, A stabilized
incompressible SPH method by relaxing the density invariancecondition, Int. J. for Applied Mathematics, Vol.2012(2012), Article
ID 139583

図-5 回転角度