FEMINGWAY 〜有限要素法解析など構造設計にまつわる数理エッセイ〜

第1話 平板(薄板)要素の歴史

固体力学の分野、すなわち数理弾性学の基礎理論は19世紀に活躍した物理学者たち(中には数学者も)によって構築されたものである。当時の物理学者たちは流体力学も含めた連続体の力学を研究の対象とすることが多かった。高校時代の数学、物理の教科書に出てくる、コーシー(Cauchy)、ポアソン(Poisson)、レイリー(Rayleigh)、キルヒホッフ(Kirchhoff)、ケルヴィン(Kelvin)といった天才たちの恩恵を我々は蒙っているのである。この人たちは、もし、20世紀に生まれていれば、間違いなく量子力学の方面で貢献した人たちであろう。

さて、ここでの主役はキルヒホッフである。彼は弾性学の分野では多くの貢献をしており、固体力学を対象とする有限要素法の参考書でも“キルヒホッフ応力”、“キルヒホッフの板理論”という用語が見当たる。このうち、ここでのテーマに大きく関係するのは後者の方である。

1850年頃、板の曲げ解析にキルヒホッフは、それまで棒の曲げ理論で採用されていた“平面保持の仮定(ベルヌイ-ナビエの仮定とも呼ばれる)”を採用することを提示した。

平面保持の仮定とは、棒の断面は変形後も平面を保持し、しかも変形前に中立軸に垂直に立てた垂線は変形後も中立軸と直交性を保つという仮定である(図1-1)。これを換言すれば、せん断変形による断面のゆがみは無視するというものである。また、変形前後の断面内の伸縮もないとする条件を加えるものだから、断面剛の仮定と言ってもいいものである。

図1‒1 断面の回転

図1‒1 断面の回転

まあ簡単に言ってみれば、棒の中立軸を板の中立面に拡張して考え(図1-2)、変形後の板断面は平面を保持するというものである。だから、対象となる板は、せん断変形が無視できるほどの薄い板厚である、ということが暗黙の前提となっている。今風の言葉で言えば、せん断エネルギーを無視した薄板理論がキルヒホッフの板理論である。

 

図1‒2 板の中立面

図1‒2 板の中立面

これだけの話だと大した内容ではないが、薄板理論では境界条件の考察で困った状況に陥ることになる。偏微分方程式の解法から言えば、1辺での境界条件は2つあればよいのに、断面力で考えると曲げモーメント、ねじりモーメント、せん断力の3種類あることになり、境界条件が余剰となってしまうのである(詳細はチモシェンコ[S.Timoshenko]著の“Theory of Plate and Shells”を参照。邦訳本あり)。

キルヒホッフは“有効せん断力”という物理量を導入して、この問題を解決した。もっとも、後にこの力学的解釈を与えたのはケルヴィンと彼の僚友テート(Tait)であったらしい。このことは、有限要素法で薄板要素を使用した場合にも大いに注意する必要があるところである。これについては次の話のテーマであるので、その説明はそちらに譲る。

ともかく、有限要素法が出現する以前、板の数値解析といえば、差分法、級数解法どちらにしても古典的なキルヒホッフの板理論が主流であり、平板理論を解説した数ある書籍もほとんどこの理論であった。

必然的に、有限要素法初期の段階の1960年代に採用された板要素もキルヒホッフの板理論をベースにしたものであった。日本での汎用型有限要素法プログラムの普及草創期に活躍したSAP IV1 が採用していた板要素もこれであった。

いよいよ、有限要素法の薄板要素の話になるが、その前に要素間の連続性ということを話しておく必要がある。棒要素が点で接合されるのに対して、板要素は要素どうしが辺でつながっている。そこで、変形の連続性を考えれば、隣接する板要素に共通する辺に沿ったたわみの連続性のほか、たわみ角の連続性も必要となる。これを数学的に言えば、たわみを表す変数 W の連続性と W の座標系変数に関する1回微分も連続性が必要ということである。W までの連続性を“C0連続”、1回微分までの連続性を“C1連続”といっている。

板要素が有限要素法の他の連続体要素と際立って違う点は、例えばソリッド要素などはC0連続であればいいのに対して、C1連続まで考える必要があるという点である。C1連続まで満足する板要素を一般に適合要素、そうでないのを非適合要素という。字面から言えば、前者の方がいい要素のように思ってしまいがちだが、そうとは限らないのが有限要素法の厄介なところである。

筆者の以前からの感想だが、有限要素法とは数学的というよりも工学的というか、かなり経験上の知識が支配する手法であると思っている。

意欲ある方が、有限要素法のプログラムを開発しようと思って、巷に出ている有限要素法の参考書を買ってきてプログラミングしたとしよう。結果は骨組要素を除いて、教科書的なプログラムはできても、実務に耐えるものはできないということになる。これは、「有限要素法は経験工学である」という一面と大いに関係がある。

閑話休題。歴史的な経緯として、有限要素法の初期の研究者、開発者たちは、ともかく適合要素の薄板要素を開発しようとした。これがまた、悪戦苦闘の歴史であり、C1連続の要素を生み出すのはなかなか困難であったのである。1要素の節点数を増やす方法や1要素を細分割する方法など、色々試みられた時代が1960年代である。そんな中で、1965年、有名な HCT 要素(開発者Hsieh & Clough & Tocher の3人の名前に由来)が開発されたのである。これが SAP IV に採用されていた薄板要素である。

HCT要素の定式化は三角形要素がベースとなっている。ここで、誘導過程の式を書き下すつもりは毛頭ないが、概略を言っておくと、1つの三角形をさらに重心位置で3個の小三角形に分割し(図1-3)、この各小三角形内で変位関数を仮定する方法である。最後には、重心点に関する剛性マトリックス内の係数を縮約させて元の三角形に対する剛性マトリックスへと導くというわけである。

図1‒3 三角形の細分割

図1‒3 三角形の細分割

いま、有限要素法を勉強する人が仮に10人いたとしよう。この10人に、この剛性マトリックスの誘導をチャレンジしてもらうとする。最後まで、根気が続いて残る根性者は何人いるだろうかと思えるほど、だらだらと長い式の続く定式化が必要な要素なのである。一般に、適合要素の定式化はかなり煩雑な式が延々展開されるのが通常である。

使用頻度の多い四角形要素はどうするのかと言えば、HCT要素では一旦、四角形の重心位置に仮の節点を設定し、その節点を使って、まず4つの三角形に分割する。各三角形での剛性マトリックスは既に誘導されているのだから、これら4つの剛性マトリックスを合成すると、内部節点を持つ5節点四角形要素の剛性マトリックスが出来上がる。仕上げは、内部節点の自由度が隣接要素との関連を持たない内部自由度であることを利用してマトリックスの縮約を図るのである。こうして、4節点四角形薄板要素の剛性マトリックスが求まるのである。

 

HCT 四角形要素の弱点は2つある。

  1. メッシュ分割数に対する収束性が他の薄板要素に対して悪いこと。これは同じ問題を解析するにしても、HCT 要素ではメッシュを細かく切る必要性を物語っている。
  2. 定式化の方法からいっても、応力や曲げモーメント等の出力は要素重心位置での4つの三角形要素から求まった値の平均値であり、やや厳密性に欠けることと、その位置でしか出力されないことである(もっとも、SAP IV では各三角形での値を四角形の各辺に割り当てているが)。ユーザーは応力のアウトプット位置を節点位置で要求することが多い。これは、応力解析においては、一般に、構造物の表面、縁等で高い応力値が出るからであり、要素中心のみの出力応力では外挿補間しても同じ値になり、意味がなくなってしまう。したがって、応力勾配のきつい解析の場合、その付近は相当、メッシュを細かくする必要がある。

今では古典的要素ともいえる適合要素(説明が前後してしまったが、本書で対象としている有限要素法は変位を仮定する変位型有限要素法である。市販の汎用コードはほとんどがこれである)の開発とは別に、C1連続にとらわれない非適合要素の開発もなされていたのが薄板要素開発の実状であった。O. C. Zienkiewicz の参考書などで紹介されているアイソパラメトリック型の薄板要素もその1つである。

非適合薄板要素の実例として、ここでは DKT要素を紹介しておこう。この要素は1980年代に開発され、その有効性は広く認められているものである。数個の離散的でしか、キルヒホッフの仮定を満足しない所からその名がついている。“Discrete Kirchhoff Theory”がその由来である。ただし、この要素も元々は三角形要素の開発から出発しており、“Discrete Kirchhoff Triangle”だという説もある。後者の説をとるなら、四角形要素の場合、DKQ(“Discrete Kirchhoff Quadlaterial”)というべきか。どちらにしても、HCT要素と違って三角形要素も四角形要素も同じ考えで定式化できる方法である。

HCT要素がたわみ W のみを独立変数にとって変位関数を仮定した方法に対して、DKT要素では W だけでなく2つのたわみ角βx、βyも変位関数に含める方法である。DKT要素の定式化のスタートは、三角形要素では6節点、四角形要素では8節点の2次要素でスタートする。すなわち、各辺の中間位置に節点をもつ仮想要素をまず考える。

適合要素の誘導式に比較すると、ずいぶん煩雑さの取れた数式となるが、それでも結構、面倒な定式化の式が続き、この過程の中で中間節点ではキルヒホッフの仮定が満足されるという条件を導入すると、中間節点での自由度が消去されて、3節点三角形要素と4節点四角形要素のDKT要素ができあがる(図1-4)

図1‒4 DKT三角形要素

図1‒4 DKT三角形要素

DKT要素はかなり実績がある要素であり、今では薄板要素での代表的な要素となっている。

 

【追記】
ここで紹介した板要素は多くの市販のコードで採用されていた変位仮定型の古典的要素の一部である。その他、有限要素法の進化の過程で応力仮定型要素や、変位/応力の両方を仮定するハイブリッド型要素も存在していることを一言いっておく。

2000年1月 記

  1. SAP IV とは汎用型有限要素法プログラムの源流の1つであった SAP(Structural Analysis Program)の4番目のバージョンのソフトウェアである。SAP は1970年、カリフォルニア大学バークレー校にいた Wilson 教授が中心となって、後に非線形有限要素法プログラム ADINA の開発で有名になった Bathe が開発したソフトウェアである。SAP IV は1974年にリリースされており、ソースプログラムで提供されたものだから、日本の工学分野の教育機関や企業でも多く導入され、その後の日本での有限要素法の普及に少なからず貢献したソフトウェアである。 []

Advertisement

読者からの寄せられたコメント

  1. 今村純也 より:

    DKT要素の定式化のスタートは、三角形要素では6節点、四角形要素では8節点の2次要素でスタートする。
    との認識は正しいと思います。ただ、”ベクトルθ(ゴチ)要素の”が重要です。大前提は平面保持仮定ですから、面内変位・ひずみを表すために変数θからスタートとの認識は正しいとするものです。ただし、θには(curlθ=0)が付帯します。DKT,DKQ要素には付帯させていない点は問題です。
    すなわち、ベクトル場θは非回転成分であり、∇φで表されるので、w=φであり、ひずみエネルギーはdivθに比例して表される点が、上述の根拠でもあります。
    C1連続要素のくだりで”工学”で逃げている印象を受けます。工学的とは”判断”とも言います。最近、C1連続要素が重調和関数式の力学問題に適用できない理由がやっと分かってきました。第3境界値問題にしか適用できないからです。
    “2008年、還暦を前に”ではまだお若い。手を動かし続けていれば、まだまだプログラミングはできる、というのが私の持論です。絵描きや、楽器奏者のように!
    ちょっと書き過ぎたようで、ごめんなさい。(まだ、あるが。)

    追伸:原田さんは計算工学会ではお見かけしないようですが?

今村純也 へ返信する コメントをキャンセル

ページ上部へ