今回は有限要素法の収束性と誤差について解説します。私たちが直面する多くの構造問題では、厳密解がない、または分かっていない問題が大部分です。しかし、それでも安心して有限要素法によるCAE を使えるのは、要素分割を密にすると、ある一定の解に近づくことが分かっているからです。今回はこの収束性を考えます。
第３回「有限要素法の収束性と誤差」を解説していきます。
第３回もくじ
1. 有限要素法の収束性の条件
有限要素法では分割を密にするに従って、その解は厳密解に近づくことが保証されます。一般に変位法による有限要素法では、分割を増すにつれて、変位の小さい方から厳密解に近づきます。これは仮定した要素内変位場によるひずみエネルギーが、真のひずみエネルギーより必ず小さいからです。有限要素法での収束性を保証する変位関数の条件は、ツェンキーヴィッツの教科書に記載されています（図1）。
ここで変位関数とは、要素内の変位分布を表現するために定義された関数です。変位関数は通常、低次多項式で表され、最も簡単な変位関数は、下に示す三角形1 次要素でのx、y の1 次関数です。これを例に説明します。
ここでa i 、b i は未定係数です。上式を用いてひずみは次式で表されます。
ここで、a 1 、 a 2 、b 1 、b 2 は定数なので、収束性を保証する変位関数の条件の条件2 は満たされます。要素内のひずみを0 と置くと、a 1 =b 2 =a 2 +b 1 =0 からb 1 =-a 2 が得られ、これらを変位関数に代入すれば、次式となります。
ここで、a 0 、b 0 は剛体変位を表し、a 2 =-b 1 は剛体回転角θを表します。また、隣り合う要素間にある辺上でも変位分布は1 次関数なので、明らかに要素間で変位は連続です。つまり、三角形1 次要素は、収束性を保証する変位関数の3 条件を満たしているので、収束性は保証されます。次は簡単な例題で要素分割の収束性を検討してみましょう。
・自由端に集中荷重を受ける片持ちはり（面内変形）
解析モデルの概要を図2 に示します。解析ケースを分割数=1×1、2×2、4×4、8×8 の4 つの構造モデルとします。材料定数は、縦弾性係数E=200GPa、ポアソン比ν=0.3とし、境界条件は左端を完全固定し、右端にFy=-10kN の荷重を作用させます。
せん断変形の効果を含まない単純はり理論でのたわみは、次式で与えられます。
ここでI は断面2 次モーメントです。最大たわみは荷重端（x=0）に生じ、次式で与えられます。
直応力σ x はx、y 方向に線形分布し、せん断応力τ xy はy 方向に2 次分布します。
σ x の最大値は固定端の上下表面で±6,000MPa、荷重端で0となります。τ xy の最大値は、はり中立面上で-150MPa、上下表面で0 となり、長手方向には一定です。
CAE ソフトNastran を用いて、四角形1 次要素で縦横の分割数をそれぞれ1、2、4、8 と増やしたときの、荷重端変位の収束状況を図3 に示します。横軸は各モデルの節点数で、縦軸は荷重端変位です。
Nastran による結果は、8×8 分割で誤差1% と、良好な収束性を示しています。CFEM とあるのが筆者自作プログラムによる結果で、定ひずみ四角形要素を用いたものがCFEM2、アイソパラメトリック1 次要素を用いたものがCFEM3 です。CFEM3r は要素剛性マトリックスの計算に次数低減積分法を適用した結果です。CFEM3r はNastran の結果にほぼ一致しており、Nastran は次数低減積分法を用いていることが分かります。
CFEM2、CFEM3 では収束するものの、8×8 分割でも理論解の各45%、65% 程度で精度が悪く、なかなか真の解に到達しません。これは、1次要素のメッシュでは曲げの変形形状を正確に表現することができないため、剛性を高く評価してしまうためです。この現象をせん断ロッキングといいます。
せん断ロッキングを回避する手法の一つに、次数低減積分法があります。次数低減積分法とは、要素剛性マトリックスの数値積分において、通常の積分点数n（有限要素法で実際に計算する点）よりも1 少ないn-1点を用いる計算方法で、結果の精度が格段に向上することが知られています（参考：O.C.Zienkiewicz, R.L.Talor, J.M.Too, Reduced IntegrationTechniques in General Analysis of Plates and Shells, Intl. J. Numer.Methods Eng., 3, 1971 年, P.275-290）。通常、1 次要素ではn=2 が、2次要素ではn=3 が用いられます。多くのCAE ソフトウェアでは2 次元・3 次元ソリッド要素、シェル要素で次数低減積分法が用いられています。図4 にNastran による結果例（変形＋σ x 濃淡図）を示します。
・自由端に集中荷重を受ける片持ちはり（面外変形）
前項と同じ例題で、板の面外変形の場合を考えます。要素の特性としては、前項で用いた平面応力要素に板曲げ要素（モーメントのみが働く要素）を重ね合わせて、3 次元空間での平板シェル要素にしています。同じく分割数を増していったとき、荷重端での面外変位の収束状況を図5 に示します。
Nastran による結果は、1×1 分割でも誤差3% と良好なまま、理論解に収束しています。CFEM-SHELL とCFEM-SHELLr とあるのは、筆者が作成したプログラムによる結果です。CFEM-SHELL はアイソパラメトリック定式化によるもので、ゆっくりと収束しています。CFEM-SHELLr はそれに次数低減積分法を用いた結果で、良好な収束性を示しています。Nastran は応力- ひずみマトリックスの計算に独自の工夫を凝らしているので、分割数の少ない状態から良い結果が得られているようです。参考として、図6 にNastran による結果例（面外変形＋板上面σ x 濃淡図）を示します。
これまでの検討で、変位法による有限要素法では、矩形要素に分割した片持ちはりモデルでの計算で分割数を増すにつれて、理論解に一様収束することが分かりました。また、同じ形状の要素であっても、要素の定式化（各現象を式に落とし込む工程）によっては精度にかなりの差があることも分かりました。CAE ソフトウェアによる実問題の大型計算の前に、これらを確認することをお勧めします。
また、CAE のプリプロセッサを用いた矩形領域の自動分割では、図7に示すように、長方形でない台形や平行四辺形などの矩形形状に分割されることが一般的です。
その場合の要素定式化の違いによる精度の検討が、米国MSC ソフトウェア社の創始者の一人であるリチャード・マクニール博士によってなされています。表1 に示すように台形や平行四辺形の場合、極端に精度が悪くなります。
2. CAE 解析での誤差
全6回で下記の内容を解説しています。
- 第1回 CAEとは？
- 第2回 有限要素法とは
- 第3回 有限要素法の収束性と誤差
- 第4回 有限要素法と他数値解法との比較
- 第5回 CAE解析の種類
- 第6回 有限要素法のテクニック
