一次元の有限要素法

ここ最近の勉強会では,http://www.u-aizu.ac.jp/~niki/
http://www.u-aizu.ac.jp/~niki/feminstr/introfem/introfem.html
の第一章に関して,一次要素を使って実際にプログラムをし,理論解と照らし合わせてあれこれやっています。

今回のブログは前回の勉強会で発表したことを簡単にまとめました。


http://www.u-aizu.ac.jp/~niki/さんは実際の理論解自体を http://www.u-aizu.ac.jp/~niki/feminstr/introfem/introfem.html
では書いてはいませんが,与えられた2階微分の方程式と初期条件,境界条件より理論解は
u(x) = - \frac{b}{2a}x^2 + \frac{R+2Lb}{a}x .... (1)
となります。

有限要素法での解法ですが,自作のプログラムで行います。
基本的には節点,要素,材質,条件等をテキストから読み取り方程式を作成してガウスの消去法で解きます。

入力ファイルは
http://d.hatena.ne.jp/malibu-bulldog/20090124/1232766033;title
勉強会宿題 - malibu-bulldogの日記

ガウスの消去法は
ガウスの消去法 - malibu-bulldogの日記
に書かれています。


今回は微分方程式の各定数をa=1.0 b=1.0,全体の長さを2L=2.0 (つまりL=1.0),初期条件をR=1.0として(1)に代入したものと,節点数が,要素数が,のデータを数値計算でといたものを比較します。


まず,理論解は上記の数値を代入すると
u(x)=-0.5x^2+3x
となります。

一方,要素数を10,節点を0.0から2.0まで0.2刻みで作った入力データを使い,一次元の有限要素法で解いた結果は次のようになりました。

一応理論解と合致していますね。