六面体一次要素による離散化 | |||||||
支配方程式の有限要素式が得られましたので,六面体一次アイソパラメトリック要素を使用した離散化式を求めることにします。 図 1 に直交座標系での六面体一次要素を,図 2 に正規化された局所座標系における六面体一次要素を示します。 一般に,要素の節点数が増えると直交座標系で数値計算を進めることは計算の手数が増えるため得策ではありません。 そこで,要素ごとの積分計算は正規化された局所座標系において行ない,その結果を全体座標系に写像変換するという方針が採られます。 ここでは熱伝導マトリックス [K] を取り扱いその過程を追ってみます。 |
|||||||
図 1 全体座標系での要素 図 2 局所座標系での要素 |
|||||||
(1)
|
|||||||
熱伝導マトリックス [ K ] において,次のことを理解しておきます。 |
|||||||
三次元解析では,要素の積分計算を全体座標系で行なうのではなく,局所座標系で行なうことは理解できたでしょうか。 |
|||||||
(2)
|
|||||||
ここで {T}e は要素の節点温度です。 |
|||||||
(3)
|
|||||||
です。 |
|||||||
六面体一次要素の補間関数 {N} の各成分は次のようになります。 | |||||||
(4)
|
|||||||
式(1)の熱伝導マトリックスの積分計算では,
補間関数 {N} の各成分の局所座標系における偏微分が必要になります。 そこでこれらを求めると以下のようになります。 |
|||||||
|
(5)
|
||||||
六面体1次アイソパラメトリック要素では,全体座標と局所座標は補間関数 {N} と要素の節点座標
( xyz 座標系) {x}e {y}e {z}e を用いて次のように補間されます。 |
|||||||
ただし, |
(6)
|
||||||
したがって,直交座標系での偏微分と局所座標系での偏微分は次の式で関係付けられます。 |
|||||||
(7)
|
|||||||
これをマトリックス表示すると次のようになります。 | |||||||
(8)
|
|||||||
なお[J ]はヤコビ行列であり,ヤコビ行列は次のようなマトリックスです。 | |||||||
(9)
|
|||||||
離散化過程で要求されているのは,式(8)の右辺にある全体座標系での補間関数の偏微分です。 |
|||||||
(10)
|
|||||||
式(10)より,全体座標系での補間関数の偏微分を求めるためにはヤコビ行列の逆行列を求めておく必要があります。 |
|||||||
ここで です。 |
(11)
|
||||||
式(11)を式(10)に代入すると, |
|||||||
(12)
|
|||||||
となります。 したがって,全体座標系における偏微分は次の3式で評価されることになります。 |
|||||||
(13)
|
|||||||
また積分の体積変化係数は次の式で表されます。 これは,局所座標系で求められた物理量を全体座標系での値に変換するために必要になります。 |
|||||||
(14)
|
|||||||
境界条件の影響を与える項は面積積分となり,積分係数の変換は次のように表されます。 *** 六面体要素では,境界面は六面から構成されています。 局所座標において p , q , r のどれかが一定となるように定義することができます。 |
|||||||
(p=±1のとき) |
(15)
|
||||||
(q=±1のとき) |
(16)
|
||||||
(r=±1のとき) |
(17)
|
||||||
以上より,熱伝導マトリックス[K],熱容量マトリックス[C],および熱流速ベクトル{F}の変換は次のようになります。 | |||||||
(18)
|
|||||||
(19)
|
|||||||
(20)
|
|||||||
ここで,面積積分における積分係数の Adpdq は,境界面に応じて前述のように 3 つに場合分けされるものとします。 以上より各マトリックス及びベクトルの被積分関数の変換式が示されたので,数値積分を行なうことができるようになります。 ここでの数値積分は, 1: 被積分関数をf(p,q,r) 2: 局所座標系における分点座標を( pi,qj,rk) 3: 分点における重みを( Wi,Wj,Wk) 4: 分点の数を n とすると,体積積分及び面積積分は次のようになります。 |
|||||||
体積積分 |
(21)
|
||||||
面積積分 |
(22)
|
||||||
分点数は 2〜3ほどとれば精度的には十分です。 |
|||||||
(23)
|
|||||||
温度の肩の添え字は時間指標で,n は現在のタイムステップを表すものとします。 時刻 n における節点温度 {T}n まで求まっているものとすれば,この項は既知量となります。 したがって,未知温度と既知量を分離すると次のような式が得られます。 |
|||||||
(24)
|
|||||||
後退差分近似は完全陰解法とも呼ばれ,数値的に安定しており,時間増分値 dt を大きくしても,解が振動を起して発散してしまうことはありません。 ただし, 時間方向の解像度である dt を大きくすると得られる解の精度は当然低下します。 有限要素法では,空間部分の離散化はメッシュにより行い,メッシュサイズが小さいほうが精度が良くなることは一般的に知られています。 同様にして,時間方向の離散化も dt を小さくするほうが精度が良いことが理解できると思います。 |
|||||||