##############################################################
有限要素法( FEM )の基礎
##############################################################

=========================================================
有限要素法でやりたいこと
=========================================================

次の偏微分方程式を解くことを考える．

.. math::
   L(u) = f \ \ \ \ ( e.g. \dfrac{d^2 u}{dx^2} + \dfrac{du}{dx} + u = f )
   
有限要素法では，任意の空間で偏微分方程式の形によらず一般的な解を与えるために次の操作をとる．

(1) 有界な空間をいくつかの領域に分割する．
(2) 領域内の解を近似解で表す．
   
簡潔にいうと，(1)ちぎって(有限要素化)，(2)近似解をペタペタ貼り付ける（フィッティング）である．
近似解は試行関数 :math:`\Psi_i` の重ね合わせで表し，ある区間の近似解 :math:`u_M` は

.. math::
   u_M = \sum ^M_{i=1} a_i \Psi_i

として表す．これが全ての区間で真値と一致する( :math:`u=u_M`  )のであれば，近似解は厳密解である．
厳密解とまでいかなくても，十分な精度の近似解を得ることができれば十分であるから，どうすれば精度の高い近似解を得ることができるかを考える．精度の良い近似解を得るためには，(1)ちぎる際に細かくちぎって解をできるだけ簡単な形でフィットできるようにする．これは折れ線で近似しても十分滑らかに見えるくらいまで細かく分割していけば精度が良くなることであり，つまりは課金（計算資源(お金)の投入）である．(2)次に，関数のフィッティングをうまく頑張ってぴったり当てはまるような関数でフィットしよう．実際，どんな関数を試行関数として選べば良いか，が，実質上重要となる．



=========================================================
重み付き残差法
=========================================================
厳密解に近い近似解を得よう．
ある区間における近似解と厳密解の誤差を評価する．

.. math::
   R = L(u_M) - f

至る所で R=0 であれば，勿論それは厳密解である．
厳密解に運良くたどり着けることができない場合でも，残差(R)を平均的にゼロにする解を見つけることができるなら，十分幸せといって良いだろう．至る所で厳密というのは難しいので，平均的には厳密解と一致する近似解を求める．
誤差の単純平均ではなく，誤差の重みづけ平均を考える．

.. math::
   \int_V w(u_M) R(u_M) dV = 0

例えばw(u_M) = 1 であった場合， :math:`\int_V R(u_M) dV = 0` となり，これは誤差の単純平均がゼロである．上式は，重み関数wによる重みづけ平均誤差であり，重み関数wの分だけ自由度が増える．このように，有限要素法では，

.. note::
   『全体の複雑な関数系を求める操作』(偏微分方程式の解法)が『区分的に基底関数をフィットする操作』へ分解された．




=========================================================
重み付き残渣法(MWR)の例
=========================================================
発熱源(ソース/シンク)を有する物体が片方で熱浴に接している様子を考える．
この系の定常状態を表す方程式は

.. math::
   \lambda \Delta T + Q = 0

ここで， :math:`\lambda` は熱伝導度， :math:`Q` は熱ソース/シンク， :math:`T` は温度である．
ここに，試行関数の重ね合わせで表現した近似解 :math:`T=\sum^n_{j=1} a_j \Psi_j`  を代入する．

.. math::
   R(a_j,x,y) &= \lambda \sum_j ( \dfrac{\partial}{\partial x^2} + \dfrac{\partial }{\partial y^2} ) a_j \Psi_j + Q \\
   &= \lambda \sum_j a_j ( \dfrac{\partial \Psi_j}{\partial x^2} + \dfrac{\partial \Psi_j}{\partial y^2} ) + Q


n個の重み関数（任意）を用いて， n個の式を作り，:math:`\int_V w_i R dV=0` をゼロにすることを考える．勿論，n個の重み関数を使うのはn個の :math:`a_j` に対して式を閉じさせるためである． :math:`i=1,2,\cdots,N` のテキトーな :math:`w_i` に対して，

.. math::
   \sum_j a_j [ \int_V \lambda ( \dfrac{\partial}{\partial x^2} + \dfrac{\partial }{\partial y^2} ) \Psi_j w_i dV ] + \int_V w_i Q dV = 0

この式のうち， :math:`b_{ij} = [ \int_V \lambda ( \dfrac{\partial}{\partial x^2} + \dfrac{\partial }{\partial y^2} ) \Psi_j w_i dV ]` の部分は， :math:`\Psi_j, w_i` を決めれば，予め計算できる．つまりここは，係数行列となるので，上式全体が行列として書ける．

.. math::
   [B] \{a\} = \{Q\}

試行関数 :math:`\Psi_j` に対する重ね合わせ係数列 :math:`a_j` は連立方程式の形で書けることがわかった．勿論，係数行列Bから逆行列を求めれば，係数列aを得ることができ，どれだけ各試行関数を重ね合わせれば近似解として良いものとなるかがわかる．