この記事は、「数値計算 Advent Calendar 2018」アドベントカレンダー22日目の記事です。
目次
- 目的
- 有限要素法とGalerkin法
- Lax-Milgramの定理
- 今後の展開
- 参考文献
目的
数値計算手法の一つである有限要素法(FEM)の数学的な解釈について考えます.
FEMの中でも、メジャーなGalerkin法の概要と数学的な意味づけと簡単な例を紹介します.
有限要素法とGalerkin法
有限要素法（FEM）とは?
端的に言うと、「微分方程式を、近似的に解くための数値解析の方法」
もう少し詳しく言うと「メッシュを用いて空間的に離散化された場の中から弱形式化された偏微分方程式の解を見つける方法」
解析対象を細かく分割し、その分割した各々の要素について近似的に応力と変位その関係 (弾性,弾塑性) を求め、要素の集合体である連続体に対して成立する線形の方程式をマトリックス法により解く
Galerkin法とは?
端的に言うと「重み付き残差法（微分方程式の境界値問題の近似解法の一つ）において、基底関数と重み関数に同一の関数系を用いる」方法
説明のため、以下のような、1次元のPoisson方程式の問題を考えます.
問題
$b:(0,1) \rightarrow \mathbb{R}, \ u_D \in \mathbb{R}, \ p \in \mathbb{R}$ が与えられたとき,
\begin{eqnarray}
\left\{
\begin{array}{ll}
\nabla^2 u + b = 0 \ in \ (0,1) \\
u(0) = u_D \\
\frac{\partial u}{\partial x}(1) = p
\end{array} \right.
\end{eqnarray}
を満たす $u:(0,1) \rightarrow \mathbb{R}$ を求める.
重み付き残差法
$u(x) = \nabla^2 u(x) + b$ を残差といい、残差にかけて領域 $\Omega=(0,1)$ で積分したら0になるように決めた関数を重み関数という.
\int_{\Omega} w(x)r(x) d\Omega = 0
ここで、重み関数をかけて、残差が0になる $u$ を求める数値解法を重み付き残差法という.
Galerkin法
近似関数の集合
$\phi = (\phi_1,\ldots, \phi_m)^T$ を $\phi_1 (0) = 0, \ldots , \phi_m (0) = 0$ を満たす $m$ 個の 1 次独立な既知関数とする．近似関数の集合を，$\alpha = (\alpha_i)_i \in \mathbb{R}^{m}$ を未定乗数として，
v_h (\alpha) = \left\{ \sum_{i \in \{1,\ldots,m\}}\alpha_i\phi_i = \alpha \cdot \phi \middle| \alpha \in \mathbb{R}^{m} \right\}
とおく．このとき，$\phi$ を基底関数とよぶ.
定義
$u_h (\alpha) − u_D \in U_h$ と $v_h (\beta) \in U_h$ を式 $\int_0^1 \frac{du}{dx} \cdot \frac{dv}{dx} dx = \int_0^1 bv dx + pv(1)$ の $u − u_D$ と $v$ に代入すれば，$\alpha \in \mathbb{R}$ $m$ を未知ベクトルとする連立1次方程式が得られる．
その方程式の解 $\alpha$ を用いて，問題の近似解を $u_h (\alpha) = u_D + \phi \cdot \alpha$ によって求める方法を Galerkin 法という．
Galerkin法の簡単な例
ODE: $ \frac{d^2u}{dx^2} + 1 = 0$
境界条件: $ u(x=0)=0, u(x=1)=0$
を考える
近似解を $$\tilde{u}(x) = \sum_{i=1}^N a_i g_i (x)$$ として,簡単のため $$N=1, g_1(x)=\sin(\pi x)$$ とする.
$$\tilde{u}(x) = a_1 \sin(\pi x)$$
ここで、$a_1$ を求める.
Galerkin法より、重み関数も基底関数も $g_1(x)$ となるため
\begin{align}
a_1 \int^{1}_{0} g_1(x) g''_1(x) dx &= - \int^{1}_{0} g_1(x) dx \\
- a_1 \pi^2 \int^{1}_{0} \sin ^2(\pi x) dx &= - \int^{1}_{0} \sin ^2(\pi x) dx \\
a_1 &= \frac{4}{\pi^3}
\end{align}
よって、
\tilde{u}(x) = \frac{4}{\pi^3} \sin(\pi x)
となる.
必要な関数解析の知識?
関数解析とは?
まず、関数解析とは何かについて、ざっくり話します。
「線型代数→有限次元のベクトル空間」に対して、「関数解析→無限次元のベクトル空間」を扱うもの.
線型代数にはなかった操作（微分とか）を位相空間論で補おうっていう分野です（認識が違っていたらすみません……）
Lax-Milgramの定理
Galerkin法の解の存在と一意性は、以下の定理で保障されます.
Lax-Milgramの定理
$V$ : Hilbert空間
$V'$ : 双対空間
$|| \cdot ||$: Vのnorm
- $ a:V \times V \rightarrow \mathbb{R} $ $a$ は双線形性をもつ
- $ |a(x, y)| \leq c || x || \cdot || y ||$
- $l \in V'$($l:V \rightarrow \mathbb{R}$ は、連続線形汎関数)
このとき, $a(u, v) = l(v)$ の解 $u \in V$ が一意に存在し,$|| u || \leq \frac{1}{\alpha} || l ||_{V'}$ となる.
資料によっては、Galerkin法に先に、Lax-Milgramの定理を入れているものもありました。
Galerkin Finite Element Method
今後の展開
- Galerkin法のさらに詳しい補足記事
- Galerkin法の実装
- Ritz法の数理的なところ
を余裕があったら書こうと思います.