平面トラスの有限要素解析①【要素剛性行列の算出】

こんにちは、トライです。

今回は平面トラスの有限要素解析（有限要素法、以下FEM）の1回目として、FEMの基本的なフローと要素剛性行列の算出方法を解説します。

FEMは機械分野、建築分野の構造解析の基本となる解析方法です。

平面トラスはFEMの基本理論を理解するのに最適な構造です。

本ブログでのプログラミング記事のモットーである

とにかくとっかりやすく、分かりやすい解説

でお送りするので、基本理論を理解したい人はぜひお読み下さい。

平面トラスのFEMのフローについて

上の図のように、節点 i と節点 j から構成されるトラス要素について考えます。

図中に表記されている各記号の意味は以下の通りです。

• ヤング係数 E
• 断面積 A
• 部材長さ L
• 要素座標系 x, y（材軸方向が x 軸, 材軸と直交方向が y 軸）
• 全体座標系 X, Y
• 要素座標系と全体座標系の成す角度 θ

平面トラスのFEMのフローは簡単にまとめると以下のようになります。

すべてを1つの記事にまとめると膨大な量になるので、今回の記事では

1. 各部材の要素座標系における要素剛性行列を算出する
2. 各部材の全体座標系における要素剛性行列を算出する【座標変換】

について解説します。

要素座標系における要素剛性行列

平面（二次元）トラスを考える前に、更に簡単な一次元トラスで考えます。

今回扱うようなピン接合のトラス架構は力学的には、

部材に曲げ応力やせん断応力は発生せず、軸力のみで外力に抵抗する架構

です。

つまり、曲げ剛性やせん断剛性はなく、軸剛性のみを持つ架構ということです。

よって、材軸方向のみの一次元トラス部材の剛性行列 k_bar は以下のように表されます。

これを剛性方程式 { f } = [ k ] { u } に当てはめると以下のようになります。

※ { } 表記はベクトル、[ ] 表記は行列を意味しています。

※ { f } は荷重ベクトル、[ k ] は剛性行列、{ u } は変位ベクトル

これが平面（二次元）トラスになると材軸方向（x 方向）だけでなく、材軸と直交方向（y 方向）成分も出てくるので、要素座標系における剛性行列は以下の式になります。

また、要素座標系における剛性方程式は以下のようになります。

トラス部材自体の変形は材軸方向の伸び縮みのみで、材軸と直交方向の変形（曲げ変形やせん断変形）は生じません。

そのため、u_iy や u_jy の係数となる行列の2行目と4行目、2列目と4列目はすべて0になります。

全体座標系における要素剛性行列【座標変換】

次に、先ほど求めた「要素座標系における要素剛性行列」を「全体座標系における要素剛性行列」に変換します。

どういうことか分かりやすくするために、例として以下の架構を見てみましょう。

この場合、

• 節点１と節点２から成る部材の要素座標軸（材軸）の向きは右向き
• 節点２と節点３から成る部材の要素座標軸の向きは左上向き
• 節点１と節点３から成る部材の要素座標軸の向きは右上向き

となり、部材ごとに座標軸の向きがすべて異なります。

架構全体の剛性行列を算出する際に、各部材の座標軸の向き、つまり各部材の剛性行列の算出ルールがバラバラなまま単純に足し合わせることはできません。

そのため、各部材の剛性行列を架構全体の座標系（＝共通の座標系）で計算し直す必要があります。

この作業を座標変換と言います。

座標変換をするためには、座標変換行列 T を用います。

座標変換行列は、高校数学Cの行列という単元で勉強する回転行列から求められます。

ちなみに高校数学Cで習う回転行列は以下のような行列です。

この θ に -θ を代入すると、一次元トラスの座標変換行列 T が求まります。

平面（二次元）トラスの座標変換行列 T は以下のようになります。

この座標変換行列 T を用いて、要素座標系の剛性方程式を全体座標系に変換します。

詳細説明は割愛します。式だけ以下に示しておきます。

以上より、全体座標系における要素剛性行列は以下のようになります。

となります。お疲れ様でした。

まとめ：要素剛性行列は全体座標系で算出する

今回はFEMの基本的なフローと平面トラスの要素剛性行列の算出方法について解説しました。

この記事では、ある程度の理論を解説しましたが、ほどほどに理解してもらった後は、全体座標系における要素剛性行列の式

↓

これだけ理解すればOKです。

次回は、各部材の要素剛性行列を接合（アセンブル）し、架構の全体剛性行列を作成する方法から、架構の応力（軸力）、変形を求めるまでを解説します。