弱形式ラプラス方程式の有限要素法離散化

弱形式で書かれたラプラス方程式は以下の通り

$$\int_V (\nabla \phi) \cdot (\nabla \delta \phi) dV = \int_{S... ...hi\ni \{ \delta\phi \mid \delta\phi=0 \quad (on \quad S_1) \}$$ (14)

$$\phi = g_1 \qquad (\quad on \quad S_1 \quad)$$ (15)

式を成分で書くと

$$\int_V \frac{\partial\delta\phi}{\partial x_k} \frac{\partial\phi}{\partial x_k} dV = \int_{S_2} g_2 \delta \phi dS$$ (16)

関数 $\phi,\delta\phi$を形状関数$N$を用いて離散化する

$$\phi = N^n\phi^n$$ (17)

$$\delta\phi = N^m\delta\phi^m$$ (18)

形状関数で補間された $\phi,\delta\phi$を式に代入して

$$\int_V \frac{\partial N^m}{\partial x_k} \frac{\partial N^n}... ...k} \phi^n \delta\phi^m dV = \int_{S_2} g_2 N^m \delta\phi^m dS$$ (19)

ここで $\phi^n,\delta\phi^m$は節点の値であり，関数ではないので積分の外に出すことができる

$$\delta\phi^m \int_V \frac{\partial N^m}{\partial x_k} \frac{... ...{\partial x_k} dV \phi^n = \delta\phi^m \int_{S_2} g_2 N^m dS$$ (20)

マトリックス方程式は次の様にかける．

$$\delta\bf {\phi} \bf {A} \bf {\phi} = \bf {\delta\phi} \bf {f}$$ (21)

成分表示すると以下のようになる．

$$\delta\phi^m A^{mn} \phi^n = \delta\phi^m f^m$$ (22)

但し，

$$A^{mn} = \int_V \frac{\partial N^m}{\partial x_k} \frac{\partial N^n}{\partial x_k} dV$$ (23)

$$f^{m} = \int_{S_2} g_2 N^m dS$$ (24)

$\bf {A}$は全体係数マトリックス， $f$は全体自由度項ベクトルと呼ばれる， 構造解析の名残から全体剛性マトリックス， 全体外力ベクトルとも呼ばれることもある．

実際に全体係数マトリックスを作る際には $A^{mn}$をすべての$m$,$n$について 計算することは膨大な計算を必要とする． なぜなら，節点数の２乗回だけ形状関数同士の内積 の計算が必要となるからである．

有限要素法ではここの点がうまくできており， 形状関数同士の内積はほとんどの場合０となる． なぜなら有限要素法では解析領域が単純な形をした 要素に分割されており， 形状関数はある節点周りの要素でしか値をもたないからである． 形状関数$N^n$,$N^m$が両方０でない値をとるのは 節点$m$,$n$の両方が属する要素上のみである． このことから要素単位で積分を実行し， その結果を足し合わせることで 計算に無駄がなく全体係数マトリックスを 作成することができる．

要素内で積分を実行した行列$\bf {A}_e$の成分$A^{ab}_e$は次のようにかける． ここで$a$,$b$は要素内の節点番号である． このように要素内で積分を実行した行列は要素係数マトリックス， または要素剛性マトリックスと呼ばれる．

$$A^{ab}_e = \int_{V_e} \frac{\partial N^a}{\partial x_k} \frac{\partial N^b}{\partial x_k} dV_e$$ (25)

要素内で積分を行った要素係数マトリックスから， 全体係数マトリックスを得るためには マージと呼ばれる手続きで要素係数マトリックスを足し合わせる．

以下にマージの関数の擬似コードを示す． ここで$nelem$は全要素数， $tmat$は全体係数マトリックス, $emat$は要素係数マトリックスである． $nnoel$は１要素に属する節点数である． また，要素内節点番号から全体節点番号を得る配列を$lnods$とした．

要は要素係数行列は要素内節点番号がindexになっていて， 全体係数行列は全体節点番号がindexになっているので 要素内節点番号を全体節点番号に 対応づけながら要素係数行列を全体係数行列に足し合わせるということである．

通常，マージは次のように簡単な式で書く

$$\bf {A} = \sum_e \bf {A}_e$$ (26)

ここまで，式の離散化を議論したが， 実際は式のDirichlet条件も組み込んだ マトリックス方程式を解かなければ解が求まらない． 以下ではDirichlet条件を組み込むための マトリックス方程式の変形について述べる．

式から得られたマトリックス方程式は以下の通り

$$\delta\phi^m A^{mn} \phi^n = \delta\phi^m f^m$$ (27)

ここで$\delta\phi^m$は任意のであったので， 両辺から消去することができる． ただし，Derichlet条件が課せられた境界上では $\delta\phi^m=0$であるためこの限りではない． $\hat{m}$をDerichlet境界条件が課された境界に属さない節点とすると $\delta\phi$を削除したマトリックス方程式は次の様にかける．

$$A^{\hat{m}n}\phi^n=f^{\hat{m}}$$ (28)

ここで$\phi^n$はDirichlet境界上で値が既知である． Dirichlet境界上の値が既知の節点を$\bar{n}$ とする．またDirichlet境界上にない， 値が未知数である，節点を$\hat{n}$とすると， マトリックス方程式は次のように変形される．

$$A^{\hat{m}n}\phi^n = A^{\hat{m}\hat{n}}\phi^{\hat{n}} + A^{\... ...\phi^{\hat{n}} + A^{\hat{m}\bar{n}}g_1^{\bar{n}} = f^{\hat{m}}$$ (29)

ここでDirichlet境界上の節点の値を $g_1^\bar{n}$とした． 未知数について整理すると

$$A^{\hat{m}\hat{n}}\phi^{\hat{n}} = f^{\hat{m}} - A^{\hat{m}\bar{n}}g_1^{\bar{n}}$$ (30)

式を解くことで離散解が得られる．