二つ石原石山に対する カルマンフィルタ有限要素法の適用

Presentation on theme: "二つ石原石山に対する カルマンフィルタ有限要素法の適用"— Presentation transcript:

1 二つ石原石山に対する カルマンフィルタ有限要素法の適用
加藤　有祐 中央大学大学院　応用力学研究室　Ｍ１

2 カルマンチームの歴史 カルマンチームの歴史はたぶん９５年くらいから始まります。 見てわかるようにshallwを用いた潮流の推定が 年度
研究者 数値解析例 研究テーマ 1995年度 稲本 --- カルマンフィルタ有限要素法の開発 1996年度 近藤 潮流推定問題 カルマンフィルタ有限要素法を用いた潮流現象の推定 1997年度 1998年度 須磨 地温推定問題 カルマンフィルタを用いた地中温度の最適推定 西脇 拡張カルマンフィルタを用いた入射波の同定 1999年度 藤本 カルマンフィルタ有限要素法に対する領域分割法の適用 2000年度 船越 潮流・濃度推定問題 カルマンフィルタ有限要素法を用いた御宿漁港における汚染物質濃度の推定 2001年度 2002年度 米川 AIC（赤池情報量基準）を考慮したカルマンフィルタ有限要素法 2003年度 中尾 土壌汚染の推定問題 カルマンフィルタ有限要素法を用いた土壌汚染の推定 2004年度 須賀 陸境界条件を考慮したカルマンフィルタ有限要素法 2005年度 脇田 開境界条件を考慮したカルマンフィルタ有限要素法 2005年度(B4) 加藤（有) 非圧縮粘性流れ場における流況推定問題 非圧縮粘性流れに対するカルマンフィルタ有限要素法の適用 2006年度(M1) 弾性波動推定問題 弾性波動伝播問題に対するカルマンフィルタ有限要素法の適用 カルマンチームの歴史はたぶん９５年くらいから始まります。 見てわかるようにshallwを用いた潮流の推定が 2006年度 加藤（拓） 弾性波動推定問題 弾性波動伝播問題に対するカルマンフィルタ有限要素法の適用 2007年度 尾島 － 未定

3 カルマンフィルタとは？ 予測 濾波 平滑 カルマンフィルター 宇宙工学 制御工学 通信工学 統計学 土木工学 経済学 Filtering
noise ノイズで乱された観測値から 高精度な解析解を推定する 時間方向の推定 空間方向の推定 宇宙工学　制御工学　 通信工学 統計学 　 土木工学　 経済学 予測 　　観測区間より将来の状態を推定 濾波 現在の観測区間の状態を推定 　 本研究では有限要素法とカルマンフィルターの２つの手法を用います。 カルマンフィルタというのは、１９６０年代にＲ．Ｅ．Ｋａｌｍａｎと Ｒ．Ｓ．Ｂｕｃｙが発表したもので、逐次法のひとつです。 当初は宇宙工学の分野で、人工衛星やミサイルなどの軌道を推定する 目的で使われていました。 それが現在では制御工学や通信工学、土木工学、経済学にまで発展して 使われています。 このカルマンフィルタは、ノイズで乱された観測値から 高精度な解析解を推定することができます。 つまり、ある観測されたデータを確率過程に基づいた カルマンフィルタというフィルタリングをかけることで、 ゴミを除去するフィルタのようにノイズを取り除くことが出来ます。 これには３種類ありまして、一つ目は予測です。 現在の観測区間よりも将来の状態を推定することができます。 二つ目は濾波で、現在の観測区間の状態を推定できます。 また平滑というのもあり、過去の状態を推定を指します。 本研究では濾波を用います。 しかしこのカルマンフィルタは大きな欠点があり、 時間方向に対する推定しかできません。 つまり空間方向の推定を行うことができないことになります。 このままでは 平滑 　　観測区間より過去の状態を推定

4 カルマンフィルタとは？ KALMAN FILTER 機械的誤差 人為的誤差

5 カルマンフィルター有限要素法は時間方向だけでなく カルマンフィルター有限要素法は限られた点から
カルマンフィルター + 有限要素法 カルマンフィルター有限要素法は時間方向だけでなく 空間方向にも推定できる カルマンフィルター有限要素法は限られた点から 全体を推定できる

6 カルマンフィルターの基礎方程式 システムモデル 観測モデル : 状態ベクトル : 観測ノイズ : システムノイズ : 観測ベクトル
As a state equation of Kalman filter , it is applied state space model. State space model consist system model and observation model. : 状態遷移行列 : 観測行列 : 駆動行列

7 ＫＦ－ＦＥＭ 有限要素方程式 これらはカルマンフィルタのアルゴリズムです。 Ｋはカルマンゲイン、Ｐは推定誤差共分散、
Γは予測誤差協分散、ｘ＊は推定値で、ｘ＾は最適推定値です。 この赤い部分は状態遷移行列ですが、このマトリックスに 有限要素方程式を挿入することで、カルマンフィルタ有限要素法 となります。

8 弾性体にカルマンフィルタ有限要素法を適用する
目的 弾性体にカルマンフィルタ有限要素法を適用する 　　　実問題（二つ石）への適用 弾性体におけるＫＦ－ＦＥＭの 有効性の検証 本研究は、流体の運動を表すナビエ・ストークス方程式にカルマンフィルタ有限要素法を 適用することを目的とします。 現在まで、カルマンフィルタ有限要素法での推定には浅水長波方程式を主に用いられて いますが、ナビエストークス方程式は様々な問題点から適用例がありませんでした。 本研究では、そういった点から、ナビエ・ストークス方程式にカルマンフィルタ有限要素法を 適用することができるかどうかを検証します。

9 支配方程式 :全応力 :物体力 :密度 :加速度 :ひずみ :変位 応力の釣り合い方程式 ひずみー変位方程式 応力ーひずみ方程式
State equation :ひずみ 応力ーひずみ方程式 :変位 弾性応力ーひずみ行列 As a state equation, Balance of stress equation, Strain – displacement equation Stress-strain equation.

10 有限要素方程式 空間方向の離散化 時間方向の離散化
The spacial discretization is carried out applied the garlkn method Considering the effect of damping This finite elemen equation is obtained.

11 Finite element equation
有限要素方程式 Finite element equation This equation is finite element equation. S,L,R is expresses these.

12 状態遷移行列 F mat 状態遷移行列 F This matrix of finite element equation shows State transition matrix.

13 数値解析例１ 5m 10m 5m Nodes : 2577 5m Elements : 12008
Numerical Example 5m As a numerical example, This figure shows the finite element mesh. The total number of nodes and elements are 2577 and 12008, respectively. The image of this mesh is tunnel. Nodes : 2577 5m Elements : 12008

14 観測点 Observation point 2 Observation point 1 Estimation point
The upper surface assumed the ground. The other surfaces is slip conditions. As a observation point, this point and this point are chosen. And estimation point is set to this point. Estimation point

15 外力 弾性係数 密度 ポアソン比 Δt 1000KN/m2 Each parameter is expressed these.
External force Each parameter is expressed these. As a external force, Uniformly distributed load is added to tunnel face. This external force is image of blasting excavation. 1000KN/m2

16 観測データ X-direction Y-direction Z-direction Acceleration[m/sec2]
Time[s] 観測データ X-direction Acceleration[m/sec2] Time[s] Time[s] Observation data Y-direction Z-direction Acceleration[m/sec2] Acceleration[m/sec2] As a observation data, These figures show the analysis acceleration at this observation point using finite element method. This figure show the x-direction,y-direction,z-direction. In this research, add noise to observation data in order to ascertain the effect of Kalman filter finite element method. This research uses these observation data. Time[s] Time[s] Time[s]

17 結果＜加速度＞ X-direction Y-direction Z-direction Acceleration[m/sec2]
Time[s] X-direction Acceleration[m/sec2] Time[s] Result –acc- Acceleration[m/sec2] Time[s] Y-direction Z-direction As a result of acceleration, Red line is estimation value. Blue point is analysis value using FEM. Two line is in good agreement.

18 結果＜速度＞ X-direction Y-direction Z-direction velocity[m/sec]
Time[s] velocity[m/sec] Time[s] Result –vel- velocity[m/sec] Time[s] Y-direction Z-direction

19 結果＜変位＞ X-direction Z-direction Y-direction displacement[m]
Time[s] displacement[m] Result –dis- displacement[m] Z-direction Time[s] Y-direction Time[s]

20 数値解析例２＜箕ノ輪山＞

21 数値解析例２＜二つ石原石山＞ Observed points Blasting area

22 数値解析例２＜二つ石原石山＞ Nodes : 5056 Elements : 24961 As a numerical example,
This model is kinowa mountain. Kinowa mountain located in miyagi prefecture. This figure shows the finite element mesh. The total number of nodes and elements are 5056 and 24961, respectively. Magnify the image for display. Nodes : 5056 Elements : 24961

23 外力 密度 ポアソン比 Δt Each parameter is expressed these. As a external force,
Uniformly distributed load is added to this face. In fact, blasting powder that is called unfo is used. This external force is image of blasting excavation.

24 Layer CM CL Young’s modulus CM C L
This figure shows the image of layer. This green part is CM layer. This red part is CL layer.

25 観測点 観測点 推定点 This point and this point are observation points.
Accelerometer set to these points.

26 観測データ X-direction Y-direction Z-direction Time[s] Time[s] Time[s]
Acceleration[m/sec2] Time[s] Result –vel- Y-direction Z-direction Acceleration[m/sec2] Acceleration[m/sec2] As a result of acceleration at this point, Red line is analysis value using FEM. Blue line is observed data at this point. Compared with two lines, accuracy of calculation falls in comparison with this point. This reason depend on the Damping coefficient and young’s modules. So, Require ample examination of these Parameter. この点での解析解より精度が落ちている。 その理由として、この点より、距離や位置が荷重面より離れているため、 減衰係数や弾性係数の影響を受けている。 さらに適切な係数を考慮する必要がある。 Time[s] Time[s]

27 推定結果＜加速度＞ X-direction 100 100 Y-direction Z-direction Time[s] Time[s]
Acceleration[m/sec2] Time[s] Result –acc- 100 100 Y-direction Z-direction Acceleration[m/sec2] Acceleration[m/sec2] As a result of acceleration at this point, Red line is analysis value using FEM. Blue line is observed data at this point. Compared with two lines, Amplitude is almost same. Focus attention on observed value, occurred small amplitude vibration. But purpose of this research, peak amplitude, structure behavior are important. Therefore, good results are obtained in this case. Time[s] Time[s]

28 弾性体にＫＦ－ＦＥＭを適用することができた トンネルモデルにおいて観測誤差を取り除くことができた 二つ石原石山において精度の良い推定ができた
おわりに 弾性体にＫＦ－ＦＥＭを適用することができた トンネルモデルにおいて観測誤差を取り除くことができた 二つ石原石山において精度の良い推定ができた In my research, many problems is solved.

29 仮定 システム方程式の線形性（linearity) システム及び観測雑音の白色性（whiteness) 雑音のガウス性（gaussian)
最小２乗規範（quadratic criteria)

30 Future Work Application to the Actual Model
Identification of Young’s modulus Using the KF-FEM 実問題に適用したい。 そのために たくさんの節点数や要素数で計算できるように工夫する。 そしてカルマンフィルターを用いたヤング係数の同定を行う。 I want to grapple with a actual problem in the future. However, the kalman filter uses a lot of calculation capacity. Therefore, it is necessary to devise it to take capacity a lot. I want to identify young’s modulus using KF-FEM.

31 Trace is diagonal matrix sum.
It’s mean error variance. If this trace converge, I obtain optimal kalman gain.

32 Displacement F matrix

33 Newmark βmethod

34 Velocity F matrix

35 Newmark βmethod

36