第44話 有限要素法は経験工学か
要素剛性マトリックスの定式化方法といい、個別問題での要素選択といい、「有限要素法は経験工学か」と思いたくなる点が、この手法には所々ある。数学的一貫性(CONSISTENCY)のなさも有限要素法には結構あるというのが今回のテーマである。
まず、話の展開の中で利用するため、線形解析における代表的平衡方程式を並べてみる。
ここに、K:要素剛性マトリックス=∫BtDBdV、M:質量マトリックス=∫NtρNdV、C:減衰マトリックス、KG:初期応力マトリックス=∫Btσ0dV、U:変位ベクトル、F:荷重ベクトル、表面荷重では ∫NtqdS、P(t):振動外力、Σ:要素数について集計することを意味する。
さらに、N:変位-節点変位関係マトリックス、B:歪み-節点変位関係マトリックス、D:応力-歪み関係マトリックス、σ0:初期応力、q:面荷重強度分布関数、ρ:質量密度とする。
有限要素法というのは、一番の基本となる要素剛性マトリックス K については実に慎重に求めている。解の収束性や要素間の連続性に気を配った要素開発がされてきた。これを構成する B マトリックスの開発には長い歴史がある。特に要素間で微係数を連続させる必要がある板とかシェル要素の開発には苦難の足跡がある。それでも数学的には妥協せざるを得ない状況もあり、その産物の結果が非適合要素といったものである。ここには、「結果が良好であれば、いいではないか」という経験工学が見られる。それでも、K マトリックスは総体に慎重に定式化されてきたといえる。それに比べると、M マトリックスや KGマトリックスの定式化には、いささか荒っぽいところがある。
振動解析において、慣性力を形成する質量マトッリクス M は数学的に一貫性を保った定式化をとるならば、K マトッリクス作成時に使用した N マトリックスを使用すべきである。実際、そうした M マトッリクスがある。これを整合質量マトリックス(Consistent Mass Matrix)という。
一方、N マトリックスなど使用せずに、要素全体の質量を適当な方法で各節点に配分してしまう集中質量マトリックス(Lumped Mass Matrix)というのもある。この方法では、ビーム要素や板要素の場合に、回転自由度に対応する回転慣性は特別な配慮をしない限り導入できないことになるが、実用的にはかえってそれらを考慮しなくてもいい場合が多い。
M マトッリクスの行列としての姿を見ると、前者が非対角項にも非ゼロの値を持つのに対して、後者は対角行列となる(図45-1)。プログラミングの立場では圧倒的に後者が有利である。
ところで、集中質量マトリックスの考えは、有限要素法の先輩格にあたる変位法の範疇にある、質点系の動力学が土木・建築分野において広く使用されてきたことが背景にある。
数学的に一貫性を保つ整合質量マトリックスがあれば、何も集中質量マトリックスなど採用せずともいいのでは、という意見も出るかもしれない。しかし、事はそう単純ではない。かえって集中質量マトリックスの方が良好な結果を出す場合も多いのである。実用的には、これで充分であるという長い歴史があり、それならプログラミングも楽な集中質量マトリックスというわけで、今も存在価値が大きいのである。ここにも、経験知識が生かされている有限要素法が見られる。
次は、数学的には固有値解析になる線形座屈解析で登場してくる初期応力マトリックスKG である。質量マトリックス M に比べると、こちらはやや深刻な状況になることもある。コンシステントな初期応力マトリックス KG を作成しようとする場合、要素タイプによっては一貫性を保つことができない状況に陥ることになる。これは N マトリックスと B マトリックスの違いが背景にある。M マトリックスの作成に N マトリックスを必要とするのに対し、KG マトリックスには B マトリックスを必要とする。N マトリックスには比較的簡単な補間多項式が採用されているのに対し、B マトリックスは、良好な剛性マトリックス開発という目的のために、各種のテクニックが要素タイプごと使用された経緯がある。
B マトリックスをそのまま KG マトリックスの被積分関数内に採り入れると、数式処理できないこともあるのである。そこで、B マトリックスの作成時に使用された形状関数とは別の、もっと簡単な形状関数の採用ということになる。一貫性はなくなるが、結果はこれまた、実用的には充分なのである。
最後の例は、荷重項の話である。
よく、等価節点力という言葉を聞かないだろうか。この言葉の意味には少し注意が必要なのだが、ここでは要素表面に載荷している荷重を要素構成節点に振り分けたときの節点力という意味で使用する。
梁要素の場合の等価節点力は、たわみ角法などの変位法に馴染んでおられる年配の方にはすぐ理解できるだろう。梁の中間に分布荷重なり集中荷重が載った場合の梁両端での節点力、節点モーメント、これが有限要素法では等価節点力と呼ばれている。当然、等価節点力には曲げモーメント力が含まれている。
ところが、梁要素と同じく回転自由度を持つ板要素では、慣習的に等価節点力の計算で曲げモーメント力を無視してしまう。実際、これでうまくいっているのである。
荷重ベクトルの計算に一貫性を持たせるため、剛性マトリックス K を求める際に使用した形状関数のマトリックス N を使用すると、必然的に回転自由度に対応する曲げモーメント力の項が生じてくる。しかし、一般にはそうしなくて、回転自由度の対応項をもたない別の補間関数(形状関数)を採用するのである。プログラム開発者でもない限り、こんな裏話は知る由もないであろう。
有限要素法の教科書をひも解くと、変分学、エネルギー法に始まってマトリックス代数が延々と続いているものだから、いかにも有限要素法は徹頭徹尾、数学的理論かと思われるが――実際、大筋にはそうなのだが――個別の詳細部では結構、経験学が入っているのである。
そんなことで、標準的な教科書を読破した方が、よし、これから有限要素法のプログラムを開発してやろうと殊勝な志を持っても、実用的なプログラムができるかといえば、なかなか難しいところがある。有限要素法は経験工学だという面も知っておく必要がある。
2007年11月記