今（2012年時点）から40年ほど前の有限要素法(FEM)の普及状況といえば-実は筆者がFEMに初めて接した頃である-海外からの輸入ソフトがようやく使われだし、市販ソフトがぼつぼつ出始めた時期だっただろうか。また、大手の企業では社内開発されたFEMソフトが利用されていた。その頃のFEM解析のパイオニアユーザーといえば、造船メーカであり、鋼橋梁メーカであり、あるいはプラントメーカといった、主にメタルを材料とした大型構造物を製作するいわゆる重工業分野の企業でった。

スケールの大きい構造物を対象とするゆえ必然的に使用される要素タイプもほとんどが梁要素や板要素であった。その当時は、ソリッド要素の出番はほとんどなかったと思う。その理由はいくつかある。良性能のソリッド要素がまだ開発されていなかったこともあるが、当時はまだ3次元メッシュを張るプリプロセッサが準備されていなかったことも大きな理由である。そして、何よりも、計算に使うコンピューター環境-当時は汎用機全盛時代-が非常にシビアな状況であったことがある。節点数がどうしても増大化するソリッド要素を使ったFEM解析は、なかなか当時のコンピューター利用では実際的ではなかったのである。自社でFEM解析できなかった企業では、解析を外部委託していたのだが、それに対する見積もり額が、1節点当たりいくらというような、パソコン世代の今の若い人たちにとっては笑い話のようなことがまかり通っていた時代であった。

 

さて今回、FEM要素ではソリッド要素とともに両横綱を張る板要素について話してみたいと思う(本話の訴点は曲面シェル要素にも共通するはず)。

板要素というのは、FEMソフトの開発者にとっても一番難しい要素なのだが、利用者にとってもつい限界を超えた期待をしがちな要素で、使用に当たっては注意を要する要素なのである。板要素の理論背景を知らずして、経験浅いユーザーが気楽に使える要素ではないことを知ってもらいたい、というのが今回のテーマである。

板要素を利用する場合、２種類の注意点がある。一つは、板の面内に関する事柄で、もう一つは、板の面外方向(面の法線方向)に関することだ。前者については、既に本エッセイの第4話で記述してあることで、これらはいわば板要素使用の経験の浅い初級者レベルの注意点である。ここでは、弾性学あるいは材料力学の知識を動員する中級者レベルの後者の話題である。

以下は、一つの簡単な構造モデルを材料にして、板要素の問題点を紹介していきたい。まず導入部は、穏やかな内容から始める。

図1にあるように細長い角棒が単純支持(但し回転はフリーでも、並進は固定)され、上からは等分布荷重が載荷された問題を考える。この構造をFEM解析するのに、梁要素、板要素、ソリッド要素の3タイプの要素を使って解析してみる。なお、ここの話題では、解析結果の絶対量には興味がなく、比較対照できる相対的数値が分かればいいので、あえて単位は示していない。読者の好みの単位を仮想してもらえば結構かと思う。

図2は図1の原構造を3タイプのFEM要素でモデリングした数理モデルで

このときの解析結果のうち、スパン中央位置の縁応力σxを示したのが下の表である。なお、注意していただきたいのは、ソリッド要素モデルでの支持点位置である。図1では、支持記号を構造下端位置に記してあるが、計算に際しては、板要素モデルとの比較するゆえ、その中立面位置である断面中央であることを言っておく(以下同様)。

表を見れば一目瞭然、すべて妥当な結果が得られている。こういう場合のように、構造系に不連続部が無い問題では、実にFEMは有用な解析ツールであり、初心者も安心して利用できるというものである。
なお、板要素のアウトプットは通常、応力を板厚方向に積分して定義した断面力(梁のそれとは違って単位幅当りの量となります)、あるいは一般化応力と呼ばれる物理量となっているはず。曲げモーメントも断面力の一成分で、板表面での応力は、曲げモーメントの定義式を逆算して得られる次の関係式で求まる。

次に、図3にあるように原構造のスパン中央の位置で段差がある問題を考えてみる。この構造を板要素の数理モデルで解析する場合、ちょっと厄介なことが発生する。スパン中央を境にして左右の板厚が違う点は当然のことだが、中立面の位置もずれることになる。平面の数理モデルで表現する板要素は、その平面は中立面に取られるのが標準的手段である。したがって、通常はFEMメッシュを張った全体の構造面は全要素の中立面がつながった面として扱われる。

隣接要素間で中立面が連続することが前提となっている板要素モデルでは、段差のある構造だとはなはだ都合が悪いのである。そこで、その対応に考えられたのが、オフセット機能というものである(図4参照)。オフセット機能とは、FEMユーザーが任意に決めた-とは言っても、実際の構造物から乖離した位置では非現実的であり、通常は構造物上のどこかの表面か中立面に設定するはず(図4では底部表面に設定)-構造面と各要素が持つ中立面とのズレを考慮した要素剛性を求める手段をいう。この機能を使用した場合、オフセット量がゼロの要素でない限り、変位量は節点位置でアウトプットされるのに対して、曲げモーメント等の断面力は節点位置のものではなく、あくまでも中立面での位置の量であることに注意する必要がある。

図3の構造解析をするには、オフセット機能を考慮した板要素を使ったFEM解析をすればいいのだが、実はこの場合、板には曲げモーメントだけでなく軸力も発生していて、話題を進めるのに少し邪魔になる。そこで、軸力が発生しない、より簡単な問題で代替してみたい。それが、図5にある板厚だけが違った中立面が連続する構造モデルである。この構造でも、ここでの話題の本質がずれることはない。

上図のように、左側の板厚を倍にして、他の諸元は図1のままで解析した結果を表示したのが図6と図7である。グラフは、板の上表面、すなわち圧縮側表面の左端から右端に向けての応力分布を示しているものである(ただし、圧縮応力を示す負号は略)。また、グラフにはソリッド要素を使用した解析結果も併せて表示してある。両図で注目していただきたいのは、左端からの距離が8の位置での板要素、ソリッド要素の結果の違いである。

左右両端からスパン中央に向かって、両者の結果がほぼ一致していたのが、スパン中央付近では、際立った違いを見せている。混乱されてはいけないので、このスパン中央位置での値はどの場所であるか、改めて図8で示めしておく。

A、B点で、どうしてソリッド要素、板要素両者の結果で大きく差異が出るのか？一言でいえば、弾性力学の基本式である構成則の基本量に、前者が歪-応力を採用しているのに対して、後者は曲率-曲げモーメントを採用しているからである。一見すると、なめらかに変化する板要素の結果の方が正しいように錯覚しがちだが、これはとんでもない間違いで、ソリッド要素の結果がこの応力場の正しい姿を暗示しているのである。以下、そのこと解説していくが、A点、B点の両者では、事情が違うため分けて解説する。

■A点での応力σ_x

弾性学の見地から言えば、A点のみならずX=8の位置での右側面の構造表面は外力が作用しない自由表面なので、σ_x=0 しかあり得ない。ソリッド要素の結果はそれを示そうとしているのである。図6では、10弱の値が出ているが、これはメッシュがこの程度のためであって(説明が遅れたが、メッシュ数は図2と同様16メッシュである)、さらに密メッシュにしていくと、極めてゼロに近づくことを知ることになる。

一方、板要素の結果の方はというと、まず板の理論を再確認することが必要となる。と言っても、完全な議論が必要ではなく、ここの話では、ポアソン比がゼロで横方向の影響は関係せず、また軸力も考慮外なので、実質梁理論と同じ簡単な基本式で済む。ここのところは、文章だけで説明するよりは数式を使う方が理解しやすいと思うので、式を少し使用してみる。

板の理論というのは、随分昔、本エッセイ第1話の中でも記述した通り、板断面の変形にあたっては、“平面保持の仮定”が導入されることから始まっている。図9を参照しながら、中立軸に垂直に立てた横断線上の変位、歪、応力は次の通りとなる。

式(2)の関係式は、もちろん微小変形が建前なので、変形時の断面回転角が微小であることが前提となっていることは言うまでもない。

次に、応力がゼロである位置(すなわち中立面)回りの曲げモーメントを次の式で定義できる。

式(3)～式(5)に出てくる回転角(たわみ角)の1階微分は実は曲率を意味することを思い出してほしい。ここで、是非確認していただきたいことが二点ある。

1. 板の曲げ理論は、式(5)が意味する曲率-曲げモーメントの構成則をスタートとしていること
2. 式(3)、式(4)が示すように、板の横断面の歪、応力は、平面保持の仮定の導入により中立面の曲率に一意的に依存していること。換言すれば、板の表面応力は、中立面の曲率で決定されてしまうこと。

これで、図8のA点での応力が理屈からはゼロであるにもかかわらず、板要素のFEM解析で、そうはならない理由が分かろうというものだ。すなわち、構造全体が変形すると、A点を含むスパン中央がたわむ。すると、その位置で中立面が曲率を持つ。すると、それに対応する曲げモーメントが発生し、曲げモーメントからは、式(1)により板表面に応力が発生してしまう。これが、ソリッド要素の結果との乖離理由である。ソリッド要素では、平面保持の仮定のような、ある種の拘束というべきものに一切しばられていないため、弾性体の物理を忠実に表現しているのである。

■B点での応力σ_x

結論から言ってしまえば、この位置の応力は正確に求まらない場所なのだ。いわゆる特異点といわれる場所である(本エッセイ第10話参照)。したがって、ソリッド要素の結果も信頼できないものとなる。ちなみに、ソリッド要素モデルで、メッシュをいくら細かくしても、B点の応力は収束するどころか、増加の一途をたどる現象が見て取れる。

勘違いしてはいけないのは、板要素モデルの場合、上でも言ったように中立面の曲率が決まれば、形式的にB点での応力が求まるので、それを正しいと判断してしまうことである。しかも、この問題の場合、採用しているメッシュ数前後で、曲げモーメントも変化しないので、応力も収束していると早計してしまうことである。

この位置の応力を求めることがいかに不安定となるかを理解するには、梁で置き換えてみればいいかもしれない-実際、今扱っている問題は、梁要素モデルでも充分である。曲げモーメントから梁断面の任意位置の応力を求める有名な次の公式で、断面２次モーメントIの値に一体何を採用すればいいのだろう。

B点の左近傍とみなすならば梁高2のIになるし、右近傍とみなすならば梁高1のIとなってしまう。これは一種の不定問題ともいえる。

 

さて、そろそろまとめに入りたい。特異点の問題は別にしても、板要素というのは、梁要素同様、弾性体を取り扱う理論である弾性学の基本量から無視されている物理量がある一方、平面保持の仮定が導入されているという工学的要素なのである。それゆえ、実際の構造物を板要素でモデル化する際には、限界があることを知るべきということだ。

 

隣接要素間で板厚が急激に変化する場合、中立面がずれる場合、あるいは角度をもって板どうしが接合される場合などなど、すべて応力場が不連続となるケースである。上の一例で解説した通り、こういう場所で正確な応力を求めることは、板要素の場合、残念ながらほぼ絶望である。

板構成の構造物では図10にあるラップジョイントがよく使われるが、応力評価する場合、接合部付近から充分離れた位置ですることを推奨する-下図は、ここでの話の材料としている板モデルがスパン中央でラップジョイントされているとし、ソリッド要素モデルで詳細解析した結果である。ラップ部およびその近傍では、とても平面保持した断面回転を前提とした応力場が成り立つとは言えないこと、一目瞭然である。

2012年3月記