最新版を正確に読む場合には pdf ファイル をどうぞ。これは web 検索のための簡易旧版です。

• 5.4.1 柱と梁の剛性方程式は厳密
  • 5.4.1.0.1 柱:
  • 5.4.1.0.2 梁:
  
  
• 5.4.2 要素剛性行列は特異
• 5.4.3 変分原理
  • 5.4.3.1 汎関数
  • 5.4.3.2 エネルギ原理と変分
  
  
• 5.4.4 有限要素法の近似度

5.4 有限要素の特徴

5.4.1 柱と梁の剛性方程式は厳密

文献[126]に示されている方法を用いて，等断面の場合， 剛性方程式(5.18)および式(5.23)が 一般的な支配方程式から直接誘導できることを示しておく。 そもそも柱と梁の影響線（Green関数）が有限要素法の変位関数と 同じであることから， 剛性方程式から得られる節点での解が厳密になると考えればいい。

5.4.1.0.1 柱:

柱のつり合い式と軸力を再度記しておくと

(5.48)

であった。この式(5.48a)を2回積分すると

$$u(x)=c_1 x+c_2-\dfrac{1}{EA}\int_0^x\int_0^\xi p(\eta)\dint\eta\dint\xi$$

となる。$x=0$, $\ell$の変位を$u_1$, $u_2$と記しておくと，この式から

$$u_1=u(0)=c_2, \quad u_2=u(\ell)=c_1 \ell+c_2 -\dfrac{1}{EA}\int_0^\ell\int_0^\xi p(\eta)\dint\eta\dint\xi$$

となることから

$$c_2=u_1, \quad c_1=\dfrac{u_2-u_1}{\ell}+\dfrac{1}{EA\ell} \int_0^\ell\int_0^\xi p(\eta)\dint\eta\dint\xi$$

である。$c_1$の式にある積分に対して無理やり部分積分を実行すると

$$\begin{eqnarray*} c_1&=&\dfrac{u_2-u_1}{\ell}+\dfrac{1}{EA\ell} \int_0^\ell \D*... ...EA} \int_0^\ell p(\xi)\left(1-\dfrac{\xi}{\ell}\right) \dint\xi \end{eqnarray*}$$

と求められるので，変位は

$$u(x)=u_1\left(1-\dfrac{x}{\ell}\right)+u_2 \dfrac{x}{\ell}+... ...nt \xi -\int_0^x \int_0^\xi p(\eta)\dint\eta\dint\xi \right\}$$

と表される。これを式(5.48b)に代入すると，軸力は

$$N(x)=EA \dfrac{u_2-u_1}{\ell}+ \int_0^\ell p(\xi)\left(1-\dfrac{\xi}{\ell}\right)\dint\xi -\int_0^x p(\xi)\dint\xi$$

と求められる。これから$x=0$, $\ell$の軸力を求め， それがそれぞれ外力$-F_1$, $F_2$に一致する（境界条件）ので

(5.49)

となる。 これは，等価節点外力項も含めて式(5.18)の剛性方程式に一致する。

5.4.1.0.2 梁:

梁も同様である。つり合い式と断面力は

(5.50)

である。これも式(5.50a)を積分すると

$$w(x)=c_1 x^3+c_2 x^2+c_3 x+c_4+w_p(x), \quad w_p(x)\equi... ...^\xi \int_0^\eta q(\zeta)\dint\zeta \dint\eta \dint\xi \dint s$$

となる。この特解$w_p$の積分に対しても柱の場合と同様の部分積分を 何度か実行すると

$$w_p(x)=\dfrac{1}{EI}\left\{ \dfrac{x^3}{6}\int_0^x q(s)\din... ...{2} q(s)\dint s -\int_0^x\dfrac{s^3}{6} q(s)\dint s \right\}$$

と求められる。あとは柱の場合と同様に，式(5.50)の 断面力等に代入して，$x=0$, $\ell$の反力を変位を用いて表せば， 最終的に

(5.51)

と求められる。 これも，等価節点外力項も含めて式(5.23)の剛性方程式に一致する。

5.4.2 要素剛性行列は特異

要素剛性行列は面白い性質を持っている。 それは，例えば柱の場合の式(5.17)からも明らかなように

$$\det\left\vert\begin{array}{rr} l_{11} & l_{12} \\ l_{21}... ...= \left(\dfrac{EA}{\ell}\right)^2\left\{(1)^2-(-1)^2\right\}=0$$

となり，その行列式が零だということである。 したがって要素剛性行列は特異行列であり，そのままでは逆行列 が存在しない。これは梁の場合も同様である。 剛性行列はバネのようなものだと解説してきた。 もしそうなら，本質的にはバネ剛性が零になっているようにも思えるが， 少し違う。図-5.27のフローチャートや これまでの解析例を思い出して欲しい。 まず境界条件を代入して，剛性方程式の左辺が既知の行をまず解いていた。 つまり，剛性行列から境界条件の分，より具体的には幾何学的境界条件の分を 取り除いてから解いていたのである。 構造はどこかを固定していないと安定した状態にはならない。 その逆の意味が，この剛性行列が特異である性質になっている。 具体的にその性質を見てみよう。

図 5.30: 両端共に軸力で圧縮された柱

柱を例に取って，その物理的な意味を 考えてみよう。図-5.30は図-4.24に 似た系で，1本の棒の両端に相反する方向に軸力が作用している。 もし$P_1=P_2$なら図-4.24になり，それは つり合っている系になる。 これを剛性方程式で解いてみよう。 分布外力は無いから$p_i\equiv 0$ $(i=1,2)$となり

$$\left\{\begin{array}{c} P_1 -P_2 \end{array}\right\}= \df... ...ay}\right)\left\{ \begin{array}{r}u_1 u_2\end{array}\right\}$$

を解けばいいのだが，先に述べたようにこの係数行列の逆行列が存在しない ため，一見したところ解けない。 これは，力の境界条件のみを与えており，どこかを適切に止めておくことを 忘れたための解の不定性であるということは既に述べた通りであるが， この方程式は解けないものだろうか。 実は線形代数には次のような定理が存在する。

Alternative Theorem: 行列$\fat{A}$を$n\times n$の正方行列と し，$\fat{u}$と$\fat{f}$を$n\times 1$の行列とする。そこで

$$\fat{A} \fat{u}=\fat{f} \eqno{(a)}$$

を解きたい。これに対し，次の2式を付随的に考える。

$$\fat{A} \fat{v}=\fat{0}, \qquad\qquad \fat{A}\supersc{t} \fat{w}=\fat{0} \eqno{(b, c)}$$

1. もし$\fat{A}$が正則なら $\left(\det\vert\fat{A}\vert\neq 0\right)$ 式(b) (c)共に 無意味な解 $(\fat{v}=\fat{0}$, $\fat{w}=\fat{0})$のみが解で， 式(a)は唯一な解を持ちそれは， $\fat{u}=\fat{A}^{-1} \fat{f}$である。
2. もし$\fat{A}$が特異であれば，式(a)の解は存在しないか 存在して唯一に定まらないかのいずれかである。解が存在するための 条件は適合条件と呼ばれ，式(c)を満足する独立したすべて の $\fat{w}=\fat{w}_1$, $\fat{w}_2, \cdots \fat{w}_m$（ここ に$m$は$\fat{A}$の大きさ$n$から$\fat{A}$の階数 を引いた値）と$\fat{f}$の内積が
   
   $$(\fat{f}, \fat{w}_j)=0, \qquad j=1,2,\cdots,m$$
   
   

   を満足しなければならない条件で与えられる。そのとき解は
   

   $$\fat{u}=\fat{u}_p+\sum_{i=1}^m c_i \fat{v}_i$$
   
   

   となる。ここに$\fat{u}_p$は式(a)を満足する一つの特解であ り，$\fat{v}_i$は式(b)のすべての解である。$c_i$は 任意定数。

この定理を前の問題に適用してみる。この場合

$$\fat{A}=\dfrac{EA}{\ell} \left(\begin{array}{cc} 1&-1 -1&... ... \fat{f}=\left\{\begin{array}{c} P_1 -P_2\end{array}\right\}$$

である。$\fat{A}$の階数は，1行目を2行目に足せば 明らかなように1であるから$m=1$となる（$m$は行列$\fat{A}$の 零固有値の数に一致する）。また$\fat{A}$は 対称行列であるから $\fat{A}\supersc{t}=\fat{A}$と なり，$\fat{v}$も$\fat{w}$も独立な解は一つで（この解は 行列$\fat{A}$の零固有値に対応した固有ベクトルに一致する）， いずれも $\fat{v}=\fat{w}=\lfloor 1\;1\rfloor\supersc{t}$でいい。 そこで$\fat{f}$と$\fat{w}$の内積をとってみると

$$(\fat{f}, \fat{w})=P_1\times 1+(-P_2)\times 1=0 \qquad \Rightarrow \qquad P_1=P_2$$

となり，適合条件を満足するためには， 外力同士がつりあっていなければならないことがわかる。 ここでは静止した解を求めようとしているのだから，この結論は当たり前である。

この適合条件が満足している($P_1=P_2=P$)と すると， $\fat{A} \fat{u}=\fat{f}$の特解は 柱が$P$に対して $\slfrac{P\ell}{EA}$縮んだ解でいいから

$$\fat{u}_p=\left\{ \begin{array}{c}\slfrac{P\ell}{EA} 0\end... ...ft\{\begin{array}{c}0 \slfrac{-P\ell}{EA}\end{array}\right\}$$

等でよく，結局，解は上の定理の表現を用いて

$$\fat{u}=\fat{u}_p+\sum_{i=1}^1 c_i \fat{v}_i = \left\{\beg... ...ght\} + c_1 \left\{\begin{array}{c} 1 1 \end{array}\right\}$$

となる。$c_1$が不定だから，両端の 距離が $\slfrac{P\ell}{EA}$だけ縮んだままならば， 柱はいくら剛体的に移動しても境界条件を満足し， つり合いも満足している。

上の例から明らかなように，解が存在するための適合条件の物理的意味は 「与えた両端の外力同士がつり合っていること」である。 つまり両端に異なる外力を与えただけの系では解が存在せず， 静的つり合い位置に留まる可能性があり得ないことを意味している。 数学と物理的直感が完全に一致していると感じられて面白い。 このように，本質的な境界条件が「適切に」与えられて いないと解が唯一には定まらないし，自然境界条件も適切で ないと解が存在しないことがわかる。

5.4.3 変分原理

5.4.3.1 汎関数

図-5.3の例を用いて，少し一般的な説明をしておく。 そのため式(5.10)の$w_a(x)$を$w(x)$とした

$$\int_0^\ell \left\{-EI w''(x)\right\} \left\{-\delta w''(x... ...) \delta w(x)\dint x -C_0 \left\{-\delta w'(\ell)\right\}=0$$ (5.52)

という仮想仕事式を検討してみよう。 ここで（著者は答を知っているので天下り的に） 関数$w(x)$に関するあるスカラー量$\Pi$を

$$\Pi\left[w(x)\right] \equiv \frac12\int_0^\ell EI \left\{-... ...- \int_0^\ell q(x) w(x)\dint x -C_0 \left\{-w'(\ell)\right\}$$ (5.53)

と定義しておこう。 何故$\slfrac12$があるのか等は以下を読み続ければわかるので， ここでは深くは考えないでいい。さて， 式(5.52)の$\delta$という記号をあたかも微分dと同じように 解釈し，関数$w(x)$を少しだけ変動させたときの$\Pi$の 変動分を式(5.53)から算定すると，特に第1項が

$$\delta\left[\frac12\int_0^\ell EI \left\{-w''(x)\right\}^2 ... ...,\left\{-w''(x)\right\} \delta\left\{-w''(x)\right\} \dint x$$

となることがわかれば，それが式(5.52)の 左辺に一致することは容易に確かめられる。 つまり$\Pi$を$w(x)$で全微分したような表現が 式(5.52)の左辺になることは，実際に被積分 関数等に全微分演算の 真似をしてみれば明らかである。すなわち，仮想仕事の原理 はあるスカラー量$\Pi$の変化率が零になる条件であり

$$\delta \Pi\left[w(x)\right] = 0$$ (5.54)

と書いてもよさそうだ。この$\Pi$のような「関数の関数」を汎関数 と呼び，その関数の変動に対する変化率$\delta\Pi$を第1変分 と呼ぶ。仮想仕事の原理は$\Pi$の第1変分が零になる条件であり， 汎関数を基礎とするこのような原理を総じて変分原理 と^5.5呼んでいる。汎関数$\Pi$の第1変分を零にするように$w(x)$を 決定しようとする問題を変分問題 と総称し， その仮想仕事の原理に対応する方程式， 例えば式(5.9a)のような 微分方程式を，この変分問題に対するEuler方程式 と呼んでいる。

こういった汎関数を常に考える必要があるかというと必ずしも そうではない。塑性や粘性等のように何らかのエネルギ的散逸を 現象論的にモデル化した場合等，問題によっては汎関数が 存在しない物理モデルもあり得る。 例えば，梁の動的問題の運動方程式は式(4.88)から

$$m \ddot{w}(x,t)+\eta I \dot{w}''''(x,t)+EI w''''(x,t)-q(x,t)=0$$

と表現でき，これに対する仮想仕事の原理 は形式的に

$$\int_0^\ell m \ddot{w}(x,t) \delta w(x,t)\dint x + \int_0^\ell \left\{-\eta I \dot{w}''(x,t)\right\} \left\{-\delta w''(x,t)\right\}\dint x$$

$$\mbox{} + \int_0^\ell \left\{-EI w''(x,t)\right\} \left\{-\delta... ...^\ell q(x,t) \delta w(x,t)\dint x -C_0(t) \left\{-\delta w'(\ell,t)\right\}=0$$ (5.55)

となるが，これに対する$\Pi$のようなものも定義することはできない。 しかし，この式そのものは後述の有限要素法の基礎式としては重要な 式である。通常，近似解の誤差評価等を行う場合や安定問題を 議論する場合を除いて，汎関数の存在は不要である。また，敢えて この汎関数に物理的な意味や定義を当てはめる必要もほとんど無い。 ただ汎関数が定義できないモデルの場合，後述の有限要素 近似解の精度（収束）には多少問題が残ることが多い[70]らしい。

なお，式(5.53)の汎関数そのものに直接近似解 $w\sim a\psi(x)$を 代入し，$a$に関する微係数を零にすることによっても未知係数$a$を 決定できる。結果的には式(5.54)と全く同じ式になるが， このような汎関数そのものを用いた手法を元来Rayleigh-Ritz法 あるいは単純にRitz法 と呼んでいた。 もちろん結果は式(5.7)を用いて$v=\psi$と 置いたGalerkin法のものと完全に同じになるが，Galerkin法という 名称が古典的な意味で用いられてきたことや， エネルギ原理という言葉からの直接的な関連もあってか，歴史的にはRitz法と いう呼び名がよく使われてきた。しかし，汎関数の定義できない 問題等のことを考えると，Galerkin法の方が汎用性があり， 昨今ではRitz法という呼び名はほとんど聞かなくなっている。

5.4.3.2 エネルギ原理と変分

敢えて式(5.53)を物理的に考察すると，第1項は内部に 蓄えられる内部ひずみエネルギ であり，第2項は分布外力の ポテンシャル，第3項は端モーメントの外力ポテンシャルと呼んでもいい。 ちなみに，この内部ひずみエネルギから外力ポテンシャルを 引いたものを全ポテンシャルエネルギ と呼ぶことがあり，上の変分原理は全ポテンシャルエネルギの停留原理 である。そして仮想仕事の原理 はその停留原理そのものである。 つまり，ある種の関数の集合の中で式(5.53)を停留させる ような関数$w(x)$が求められたとき，その関数$w(x)$が元の境界値問題の 厳密解であることを，この変分原理は述べている。

図 5.31: 汎関数の性質の模式図

代数幾何学的には，内部ひずみエネルギがいわゆる2次形式 であることから，ある関数の集合の上に描かれた全ポテンシャル エネルギのグラフが「曲面」を描き，厳密解のところでその曲面の 接線が零になる原理が仮想仕事の原理であるという解釈もできる。 上では停留と書いたが，もしこの曲面が下に凸で あれば最小原理になる。 曲面が上に凸なのか下に凸なのか調べるには，高校で習った 関数のグラフと微係数との関係から類推されるように，もう一度微分（この 場合は変分）してみてその符号の正負を調べれば いいことが直感的に理解できよう。 仮想仕事式を再度変分したものを第2変分 と称するが，式(5.53)から

$$\delta^2\Pi = \int_0^\ell EI \left\{\delta w''(x)\right\}^2 \dint x > 0$$

のようになり，曲げ剛性$EI$が正である限りこの積分は常に正になる。 つまり，全ポテンシャルエネルギは図-5.31に 模式的に描いたように下に凸の曲面になり， 厳密解あるいはある仮定した関数の集合の中での最適な近似 解に対応する全ポテンシャルエネルギは最小値をとることがわかる。 このように，境界値問題はある制約条件の元で仮想仕事の原理と 等価であり，かつそれは全ポテンシャルエネルギの最小原理 と一致している。何度も書くが，こういった物理的解釈は必ずしも必要ではなく，この 章でも積極的には用いるつもりはないが，章-6では安定^5.6規準の物理的考察に用いる。 そういった意味からも，仮想仕事の原理を単なる数学的表現とは捉えない方が 望ましいかもしれない。

1. 式(5.53)で定義された$\Pi$は，その$w(x)$が 厳密解のときに最小になる。この最小値 $\Pi\subsc{min}$を$w(x)$で表せ。
2. 演習問題5-3の6番で 求めた境界値問題の仮想仕事式を元に， それに対応する汎関数を求めよ。
3. 章-4のp.で 解いた図-4.21の，先端をバネ定数$k$の 線形バネで支持された片持ち梁の仮想仕事式を， つり合い式と 適切な境界条件から求め，それに対する汎関数が，梁の全ポテンシャル エネルギにバネの内部ひずみエネルギ $\slfrac12 k \{w(\ell)\}^2$を 加えたものであることを示せ。

5.4.4 有限要素法の近似度

簡単のために図-5.5に示した柱の例を用いて， 有限要素法の近似について概説しておく。 近似解を$u_a(x)$と記し厳密解を$u_e(x)$とすると， 厳密解はこの場合

$$p_0=-EAu_e'' \quad (0<x<\ell), \qquad u_e(0)=0, \quad EAu_e'(\ell)=P \eqno{(a)}$$

を満足し，近似解は

$$\int_0^\ell EAu_a' \delta u_a' \dint x -\int_0^\ell p_0 \delta u_a \dint x -P \delta u_a(\ell) = 0 \eqno{(b)}$$

を満足している。$\delta u_a$, $\delta u_a'$は近似解の 集合中の任意の二つの$u_a$の差と解釈していい。 式(b)第2項の$p_0$に式(a)の 第1式を代入して1回部分積分したあと，式(a)の 境界条件を代入すると最終的に仮想仕事式は

$$\int_0^\ell EA\left(u_a'-u_e'\right) \delta u_a' \dint x = 0 \eqno{(c)}$$

となる。すなわち，軸力の近似解の誤差 $\left\{N_a(x)-N_e(x)\right\}$に 重み$\delta u_a'$を乗じて全領域で加算したものが零になるように 近似解を求める方法が有限要素法になっていることを示している。 このことは図-5.6を見れば一目瞭然であろう。 そこで用いた近似関数が1次多項式であったため 重みは定数になり，軸力の近似解は厳密解の平均値を与えている。 有限要素法がいわゆる重みつき残差法 として位置付けれられる由縁もここにある。

梁の場合も同様の演算をすると

$$\int_0^\ell EI\left(w_a''-w_e''\right) \delta w_a'' \dint x = 0 \eqno{(d)}$$

という近似になっていることがわかる。ここに$w_a$, $w_e$は それぞれ$w(x)$の近似解と厳密解である。この場合も， 曲げモーメントの誤差 $\left\{M_a(x)-M_e(x)\right\}$に 重み$\delta w_a''$を乗じて全体加算したものを零にするように 近似解が求められていることを示している。図-5.15の 結果からもその実態が理解できると思う。

図 5.32: 近似解の意味

一般に，$0<x<\ell$で定義された二つの関数$f(x)$, $g(x)$の内積 を

$$\left\langle f, g\right\rangle\equiv \int_0^\ell \mbox{w}(x) f^{(n)}(x) g^{(n)}(x)\dint x$$ (5.56)

で定義することがある。上付きの$(n)$は$n$階の微分を 表し，w$(x)$は内積を一般化する重み の関数である。3次元のユークリッド幾何学でのベクトルの内積 がその基底ベクトル方向成分を用いて

$$(\fat{u}, \fat{v})=\fat{u} \fat{\cdot} \fat{v} =\fat{u}\supersc{t} \fat{v}=\sum_{i=1}^3 u_i v_i$$

で表されることを考え，これを無限次元に拡張した上で重みを 乗じて一般化したもの（重みつき内積 ）として式(5.56)を見ると，少しは理解し易いかもしれない。 この内積を零にする関数同士を直交する関数 と呼ぶが，これも幾何学から借りた概念である。 最も代表的な直交関数列としてはFourier級数 でよく用いられる三角関数が挙げられ，それは節-4.7.2で 動的問題の解法に用いた固有関数の性質としても挙げてある。 この「内積」というキーワードを用いて式(c) (d)を眺めると， 「近似解から厳密解への‘ベクトル $\left(EAu_a'-EAu_e'\right)$’と， 近似解を含む‘平面’内の任意の二つの‘ベクトル’の ‘差のベクトル $\left(\delta u_a'\right)$’との内積が零」ということになる。 そして「直交」という概念を使うと， 「近似解を含む‘平面’に厳密解から下ろした‘垂線の足’が， 有限要素解である（図-5.32参照）」ということになる。 つまり，近似解を含む集合の中で，上のように 定義された内積で測ったときに，厳密解に一番近い関数が有限要素解になっている。

また，ベクトルの大きさを測るためのノルムは内積を 用いて $\left\vert\fat{u}\right\vert^2\equiv (\fat{u}, \fat{u})$と 定義されることに対応させて，関数のノルム も上の内積を用いて

$$\left\Vert f \right\Vert^2 \equiv \left\langle f, f\right\rangle$$

と定義してもよさそうである。時折これはエネルギノルム と呼ばれるが，その表現が内部ひずみエネルギに一致しているからである。 このように考えると，上の柱の場合の有限要素解の近似は

$$\left\langle f, g\right\rangle\equiv \int_0^\ell EA f'(x) g'(x)\dint x$$

で定義した内積を用いたノルムによって

$$\left\Vert u_a - u_e \right\Vert \to \mbox{最小}=C(\mbox{要素長})^n \eqno{(e)}$$

にするような近似であることが予想される。 つまり有限要素法とは，このような内積で定義されるノルムという「ものさし」 で測ったときの，関数$u_a(x)$と$u_e(x)$の「距離」を 最小にするような近似解法であり， 厳密解との誤差 が要素長を小さくするに 従って小さくなる（希望としては$n\geq 2$）近似解法だと考えていい。 あるいは，この距離が近似誤差 であることから，上式を一種の最小2乗法 として捉えることもできるだろう。 具体的には，図-5.29に示した弾性床上の梁の例のように， ある1点の解に着目したとしても，変位の微係数で求めた曲げモーメントが 要素長の$n$乗($n\simeq 2$)の速度で厳密解に収束している^5.7ことがわかる。 ただし，内積の定義は解きたい問題毎に異なることには注意する。

最新版を正確に読むためには pdf ファイル をどうぞ。これは web 検索のための簡易旧版です。