論文

アンサンブルカルマンフィルタFEMによる波動水槽内の実水位観測値を用いた流れ場推定解析
倉橋 貴彦1・齋藤 浄2・衞藤 俊彦3
1 長岡技術科学大学大学院

機械創造工学専攻（Department of Mechanical Engineering, Graduate School of Nagaoka University of Technology）
2 長岡技術科学大学大学院

機械創造工学専攻(現： コトブキシーティング株式会社)2（Department of Mechanical Engineering, Graduate School of Nagaoka University of Technology）
3 環境都市工学科（Department of Civil Engineering, National Institute of Technology, Nagaoka College）

Estimation analysis of flow field in wave tankbased on the ensemble Kalman filter FEM using practical observed water elevation
Takahiko KURAHASHI1, Kiyora SAITO2 and Toshihiko ETO3
AbstractIn this paper, we present the estimation analysis of flow field in wave tank based on the ensemble Kalmanfilter FEM. The ensemble Kalman filter FEM is the fusion method of the ensemble Kalman filter theory
and the FEM, and the characteristic of this method is the estimation of flow field can be carried out
considering the observed physical variables. In this study, the observed water elevation is utilized for theflow estimation analysis. As the governing equation, the shallow water equation is introduced, and the“Stream Upwind Petrov Gelerkin method” (“SUPG method”) is employed to discretize the governingequation.
Key Words : ensemble Kalman filter FEM, SUPG method, shallow water flow

１．はじめにのが主であり，本論文では，非線形モデルも対象としたアンサンブルカルマンフィルタ 2), 3)を適用する．

本論文では，長岡高専内の造波水槽を用いて得アンサンブルカルマンフィルタによる解析において，られた水位測定値を用いて，流れ場の推定を行う解システム方程式および観測方程式が用いられるが，析について紹介する．近年，データ同化というネー

本検討では，システム方程式において有限要素法ミングで注目を集めており，測定値を利用したシミ（FEM）に基づく有限要素方程式を用いる点が一つュレーションによる解析結果の高精度化に関する検の特徴である．支配方程式の空間方向の離散化手法討が実施されている．理論としては，1960年に開発の一つである有限要素法は，一般に線形モデルを対されたカルマンフィルタの理論 をベースとするも象に用いられる境界要素法に比べ解析に用いる要素1)

17３．数値実験の数は多くなるという欠点もあるが，境界要素法のように基本解を求める必要もないため，非線形モデルの解析においても比較的容易に解を求めることがアンサンブルカルマンフィルタ FEM を用いた数可能である．このようなことから，本研究では，シ値実験として，長岡高専の波動水槽（図-1）を用いステム方程式において有限要素方程式を用いたアンて測定した水位の値を用いて流れ場推定解析を行う．

サンブルカルマンフィルタFEMを用いる．システム

水槽の上流部において 0.5m 間隔に 3 点，水位変動量の測定点を設定し，それ以降 4 点，流れ場の推定解析において解析結果を出力する推定点 A～D を設定する（図-2）．また，有限要素法に基づき，流れ場の推定解析を行うため，図-3 に示す有限要素メッシュを用いる．その他
計算条件としては，時間増分量Δt を 0.01s，時間ステップ数を 6000，基準水深 h を 0.46m，重力加速度 g を 9.81m/s2，システムノイズの平均値を 0，分散を 10-6，観測ノイズの平均値を 0，分散を 10-3 に設定する．また，サンプル数は1000 とする．

方程式の誘導に際しては，浅水長波方程式に有限要素法の一つであるSUPG法 4) を適用し，離散化を行う．時間方向の離散化には後退オイラー法を適用する．本論文では，長岡高専の波動水槽を用いた実験により得られた観測水位の値を用いて，アンサンブルカルマンフィルタFEMを用いた数値実験により流れ場の推定解析を行う．

２．支配方程式
本論文では，支配方程式として，式(1)，(2)に示す浅水長波方程式を導入する．式(1)，(2)は総和規約により表示しており，本論文では平面2次元のモデルを対象とする．

𝑢̇ 𝑖 + 𝑢𝑗 𝑢𝑖,𝑗 + 𝑔𝜂,𝑖 − 𝜈(𝑢𝑖,𝑗 + 𝑢𝑗,𝑖 ),𝑗 = 0(1)𝜂̇ + {(ℎ + 𝜂)𝑢𝑖 },𝑖 = 0(2)ここに，ui，h, η，g，νはx,y方向の流速成分，基準水深，基準水深からの水位変動量，重力加速度，動粘性係数を示す．動粘性係数νは式(3)により与えられる．

図-1 波動水槽の写真𝑘𝜈 = 𝑙 𝑢∗ ℎTop view

(3)

Observation point

0.256ABC

Estimation point

D

式(3)において，klはカルマン定数であり本論文では0.41と与える．またu*は摩擦速度を示しており，式Flow direction(4)により計算することができる．

図-2 水位の測定点と推定点の位置図

𝑢∗ =𝑔𝑛∗2 √𝑢𝑘 𝑢𝑘
1ℎ3(4)

ここにn*はマニングの粗度係数を示している．式(1)，(2)を空間方向および時間方向に離散化するために，
図-3 有限要素メッシュ図（0～3m まで）

SUPG法および後退オイラー法を適用する．アンサンブルカルマンフィルタFEMによる計算の流れは文献5), 6)を参照して頂くこととする．
18本数値実験に用いた水位変動量の測定値をマニン

0.080.06

スライダーボリュームの値は 0.25 である．流れ場の0.04

Water elevation, [m]グの粗度係数を図-4 に示す．測定の際に設定した推定解析においてはマニングの粗度係数をn*=0.01m-1/3s，0.1m-1/3s，0.25m-1/3s と変え数値実験を行った．推定点 A～D における水位の推定値を図-5～図-8 に示す．どのマニングの粗度係数の設定値

0.020.01

00123456-0.02

0.1
0.25-0.04

においても水位の計算を適切に行えていることを確-0.06

認できる．実際の推定解析においては，流路底面の-0.08

Time, [s]状況等からマニングの粗度係数を設定することになるが，このパラメータにより流れ場の推定解析の結図-6 粗度係数を変えた場合の推定点 B における果も異なるため，その他の水位変動量等の計測値を水位変動量の推定値用いて適切に同定する必要もあると考えられる．
0.08

Water elevation, [m]

0.06
0.04
0.02
0.01

00123456-0.02

0.1
0.25

-0.04
-0.06
-0.08

Time, [s]図-7 粗度係数を変えた場合の推定点 C における水位変動量の推定値
図-4 水路入り口部における水位変動量の測定値

0.08

0.06

Water elevation, [m]0.08Water elevation, [m]

0.06
0.04
0.02
0.01

001

2

3-0.024560.04
0.02
0.01

001

2

3-0.02

4

560.1
0.25

-0.04

0.1
-0.06

0.25

-0.08

-0.04

Time, [s]

-0.06
-0.08

図-8 粗度係数を変えた場合の推定点 D におけるTime, [s]水位変動量の推定値
図-5 粗度係数を変えた場合の推定点 A における

また，マニングの粗度係数を 0.01 と設定した場合に水位変動量の推定値

おいて，時間 t=0.2, 0.4, 0.6s における水位変動量の分布図を図-9～図-11 に示す．結果より，波動が適切に伝播している様子を確認することができる．

19

推定解析を行った．支配方程式としては，非線形項を有する浅水長波方程式を導入し，空間方向に対してはSUPG法，時間方向に対しては後退オイラー法による離散化を行った．実際の測定水位を用いた場合においても流れ場の推定は良好に行えることを確認できた．また，マニングの粗度係数を変えた場合においても，流れ場の推定解析を行えることは確認できたが，実際には租度係数も推定しながら流れ場の解析を行える方が好ましいと考えられる．物性値の同定を入れた流れ場の推定解析は今後の課題とする．
図-9 粗度係数を 0.01 m-1/3s とした場合のt=0.2s における水面形の分布図

謝辞：本論文の解析結果は，九州大学情報基盤研
究開発センターの高性能アプリケーションサーバを使わせていただいた上で計算を行ったものである．
また，水位の測定においては，著者と同大学の小林雅輝氏にも協力頂いた．ご協力頂いた皆様に対して，ここに謝意を表す．

参考文献
1)

R. E. Kalman, “A New Approach to Linear Filtering andPrediction Problems”, Transactions of the ASME. J. BasicEng., Vol.82, No.1, pp.35-45, 1960.

2)

図-10 粗度係数を 0.01 m-1/3s とした場合のG. Evensen, “Sequential data assimilation with a nonlinear
quasi-geostrophic model using Monte Carlo methods tot=0.4s における水面形の分布図

forecast error statistics”, J. Geophys. Research, Vol.99,pp.10143-10162, 1994.
3)

G. Evensen, “The Ensemble Kalman filter: theoreticalformulation and practical implementation”, Ocean Dyn.,
Vol.53, pp.343-367, 2003.

4)

中山 司, “流れ解析のための有限要素法入門”, 東京大学出版会, 2008.

5)

T. Kurahashi and K. Saito, “Investigation of the numericalaccuracy of shallow water flow estimation based on theensemble Kalman filter using the SUPG FEM”, Int. J. ofCompt. Fld. Dyn., Vol.32, pp.146-157, 2018.

6)

倉橋 貴彦, 小林 雅輝, 衛藤 俊彦, “アンサンブルカルマ

図-11 粗度係数を 0.01 m-1/3s とした場合のンフィルタ FEM に基づく浅水流推定解析に対する時

t=0.6s における水面形の分布図

変サンプル数の適用”, 長岡工業高等専門学校紀要,
Vol.54, pp.1-5, 2018.
(2019. 8. 1 受付)

４．おわりに
本論文では，波動水槽を用いた測定水位を用いてアンサンブルカルマンフィルタFEMによる流れ場の20