講 座Lecture
地すべり解析における有限要素法の利用 第３回

Finite element method for landslide analysis No.３

２．
１ 弾塑性FEMの考え方（後編）
Fundamental theories of elasto-plastic FEM
講座委員会委員長

鵜飼恵三／群馬大学工学部

Keizo UGAI／Faculty of Engineering, Gunma University
キーワード：有限要素法，地すべり，弾塑性，基礎理論
Key words：finite element method, landslide, elasto-plastic, fundamental theory２．
１．
５ 二次元弾塑性連続体の基礎方程式

り複雑になる。以下では弾塑性体の簡単な応力ひずみ関

二次元弾塑性連続体（以下，弾塑性体と呼ぶ）が満た係式について考えてみよう。詳細は文献２）を参照されさねばならない３つの基礎方程式を列挙する。力学的にたい。

安定であるには，次のような力とモーメントのつりあい

２．
１．
６ 弾完全塑性モデル

式を満たす必要がある。

地盤の設計では，変形が小さいとき土を弾性体と仮定

      …………………
  
，は，，方向の垂直応力，はせん断応力である。X ，Y は単位面積あたりに作用する，方向の物体力である。方向を水平，方向を鉛直（上向きを正とする）にとると，物体力が自重のみであるときは，X ＝０，Y ＝− となる。 は土の単位体積重量である。
２番目の式として，変位とひずみの関係式が満たされすることが多い。一方，斜面安定のような安定問題では，地盤が破壊状態にあると想定して完全塑性体とする仮定がよく用いられる。このような２つの仮定を合わせ持つ
土のモデルが，図−５に示される弾完全塑性モデルである。このモデルでは，降伏応力 （もしくは降伏ひずみ）まで土は弾性的性質を示し，それ以降は一定の応力を保ちながら大きな塑性変形が進行する。現実の土は

ひずみ硬化（２．
１（前編）の図−１）や軟化を起こしてねばならない。変形が小さいとき次式が成り立つ。

変形が進行するため，厳密にはより複雑な応力ひずみ関
…………………………………
， は， ， 方向の垂直ひずみ， はせん断ひずみ係式が必要になる。しかしながら，弾完全塑性モデルは，

   である。 ，は，方向の変位である。式 とは弾土の弾性係数と強度定数を同時に考慮できるため，簡潔であり理解が容易である。これだけの特徴でも十分な実用性と適用性を有していると考えられる。

塑性体の応力ひずみ関係とは無関係になりたつ。

土は摩擦性材料とみなされ，せん断によって破壊が生

３番目の式は，構成式と呼ばれるものである。これはじると考えられるので，降伏もしくは破壊規準として次材料の物性を表す式で，弾塑性体の場合には，応力ひずのモール・クーロン式がよく用いられる（２．
１（前編）のみ関係がこれに相当する。よく知られているように，弾式性体の応力ひずみ関係はフックの式で表される。これは二次元平面ひずみ状態では次のようになる。


  …………………………………

ここで， 
E はヤング係数，はポアソン比である。式 ，，を参照）。
 !!"#$%!
 ! "#$%!…………………………………

で示される弾性体の基礎方程式をFEMにより離散化し，定式化すると，２．
１（前編）の式

と同じ形の連立一次方程式が得られる。

 ………………………………………
は剛性マトリクスであり，弾性体の剛性に比例する。
に対する節点変位ベ式を解くと，節点力ベクトルが求められる 。
クトル
弾塑性体の場合には，応力ひずみ関係が式よりかな
１）

J. of the Jpn. Landslide Soc., Vol.40, No.2157 (2003)

図−５ 弾完全塑性モデル

５３講 座は土のせん断強度，!はすべり面上の有効垂直応力，
は全応力， は間隙水圧である。
式を主応力で表すと次のようになる。
!&!
!&!'($%!!)*'%!…………………
!
!
と&を，表示すると次のようになる。なお以下においては，有効応力表示のダッシュ記号（!）を省略する。

,, +
&  +…………………

以上に述べた関係式をさらに展開していくと，弾塑性EF
4 7
C7 AD7 AD74 ………………………………………………………

体の応力増分・ひずみ増分関係式が次のように導かれる。

は降伏関数であり，モール・クーロン式の場合には，次式で表される。

図−５や２．
１（前編）の図−１からわかるように，弾&&'($% )*'%…………………
4により，弾性または弾塑性状態から弾性状態へ変化するときR ＝１，弾塑性状態から弾塑性状態へ変化するときR ＝０，弾

塑性体では応力とひずみの間に１対１の関係はなく，弾性状態から弾塑性状態へ変化するとき０＜R ＜１，の値

塑性体がそれまでに受けた応力履歴の影響を大きく受けを取る。詳細は文献２）を参照されたい。

る。このため，弾塑性体の応力と変形状態を表現するた

２．
１．
８ 弾塑性FEMの具体的な計算方法

めに，ある状態からの応力増分とひずみ増分の関係を求

２．
１（前編）の図−２に示された均質な単純斜面をイめて，時々刻々の状態の変化を解析していくのが普通でメージしながら，斜面の弾塑性変形を計算する方法のうある。

ち代表的な４つのケースについて説明する。最初の３つ

２．
１．
７ 弾塑性体の応力増分・ひずみ増分関係式

弾塑性体の応力増分・ひずみ増分関係式の導出法の概略を以下に示す。詳細は文献２）を参照されたい。
現在の応力とひずみの状態が，各々ベクトルとで与えられる弾塑性体を考える。二次元連続体では，
1-. 1.
.
/.
/ 22

   …………………………
 0.3.
0. 3.
4が加えられ，ひずみ増分
この弾塑性体に応力増分4 が生じたとするとき，4と4 の関係式を導い4 は，次のように弾性成分4 5と塑性成てみよう。
64 に分解される。
分4 4 54 6…………………………………
4と4 5の間にフックの法則が成り立つ。すなわち
474 5 ……………………………………
7は弾性係数マトリクスと呼ばれ，弾性係数8とをこれらは次のように表される。

用いて次のように表される。

:  9 =>;7;<   9 >?…………………
9 9 
4 6と応力との間に次のよう
一方，塑性ひずみ増分

またR は係数であり，応力増分のケースでは，土を弾完全塑性モデルであると仮定する。

 斜面に自重が作用して弾塑性変形が生じるケース

図−２の斜面に自重が作用する場合を考える。自重を加えた後，斜面が崩壊に至る可能性もあるため，少しずつ自重を増やしていくことを想定する。たとえば，

9ずつ増加させてみる。 は土の単位体積重量である。
FEM計算は次の順序で行われる

 自重の一部である 9に対応する等価節点力は斜を求め，式に代入する。剛性マトリクス
i

面形状，メッシュ分割，土の弾性係数などから決定される。与えられた境界条件のもとで式

を解くと，

Gが計算される。
G
 この変位増加量 に対応する弾性応力増加量
Gを計算する。Gがモール・クーロン式またはを越えないなら，その要素については弾性状節点変位の弾性増加量

ii態が確定する。モール・クーロン式を満たすか越える場合には，その要素は弾塑性状態にあると判断される。応力状態がモール・クーロン式を越えることはありえないので，そのような要素については，モール・クーロン式を満たすように応力の補正を行う。

な関係を仮定することが多い。

4 64@ A 4 64@ A 4 64@ A 式をまとめて，ベクトル表示すると
4 64@ A …………………………………
ここで，Aは塑性ポテンシャルと呼ばれ，応力と塑性ひずみの関数である。4@は比例係数である。Aは実際の土では複雑な関数になるが，弾完全塑性モデルでは式にならって次式のようなAを用いることがある。Bは土のダイレイタンシー角である。

A&&'($B……………………………

５４

図−２（再掲） 仮想的な斜面とメッシュ分割例

J. of the Jpn. Landslide Soc., Vol.40, No.2

158 (2003)

鵜飼恵三：地すべり解析における有限要素法の利用 第３回補正された応力に対応する等価節点力を２．
１（前編）生じている。また最小主応力は斜面肩右側で負値（引っの式

張り応力）を示しており，引っ張り亀裂の発生を予想さを使って計算する。通常，この値は，最初に求めた 9に対応する等価節点力と一致しない。こGCと呼ぶ。の不一致の差を残差力 この残差力GCに対応する弾性変位増加量を式から計算し，上と同様な計算を繰り返す。残差力iii

せる。これらの現象は，斜面の破壊直前に見られるものと符号しており，２．
１（前編）の図−４の弾性解析結果からは予想できない。

 自重の作用下で安定している斜面内の地下水位が上昇するケースが所定の許容誤差より小さくなったら，求める計算

図−２の斜面内で，降雨などの浸透により地下水位が結果が得られたと判断する。

 次に 9を追加し，以上と同様な計算を繰り返す。
iv

上昇するケースを想定する。斜面は弾塑性体であると仮に示した方法で自重が作用し，途中で計算が発散するか，もしくは残差力計算を所定する。上述のケース

定の回数以上繰り返しても許容値に達しないときは，安定した状態になっているとする。その状態を出発点（初そこで斜面が破壊したと判断する。計算が収束する期応力状態）として解析を行う。斜面内の間隙水圧分布場合には，加えた自重の和が になるまで同様な計がわかっていれば，間隙水圧ベクトル算を繰り返す。このようにして，弾塑性計算が行わに対応する等価節点力ベクトルを２．
１（前編）の式から計算できるれる。

 以上に示したのは，弾塑性計算アルゴリズムの一ので，それを節点に作用させる。時間と共に間隙水圧が例である。計算の効率性と収束性を追求した多くのに変換させながら，徐々に作用させていく。このやり方アルゴリズムがこれまでに提案されており，非線形は自重の場合と同じである。間隙水圧が大きくなると，v３）計算法と呼ばれる研究の一分野を形成している 。

上昇していく場合には，間隙水圧の増加量を等価節点力式

からわかるようにせん断強度が小さくなるので，斜本ケースの一例として，図−２の斜面に自重を作用さ面の弾塑性変形が徐々に進行し，斜面崩壊へ近づいていせ，弾塑性解析をした結果を図−６に示す。土質定数は，%#$%を１．１でkPa，B＝０°とした。図−６の結果は，cと"く。このようにして，間隙水圧の上昇による斜面破壊のE ＝２０００００kPa， ＝０．
２５， ＝２０kN／m３， ＝２０°，c＝１２

シミュレーションが実施されることになる。

割り低減させたときの結果で，斜面はほぼ破壊直前の状斜面の土は膨張する。またダイレイタンシーが生じる場態にある。図−６より傾斜部付近の変位は左斜め下方に合も土は膨張（もしくは収縮）する。土の膨張により水

ところで間隙水圧の上昇は有効応力を減少させるため，圧も変化するはずだが，通常の斜面安定解析ではこの影響を考慮することは少ない。このような解析は土の排水状態を仮定するのと同じであり，非連成解析と呼ばれることがある。一方，膨張（土の体積変化）と間隙水圧の相互作用を考慮する方法は連成解析と呼ばれる。連成解
析は軟弱地盤の変形解析でよく用いられる。

 自重の作用下で安定している地盤の一部が掘削されるケース
ケース

の前半と同様に，図−２の斜面が自重の作用下で安定した初期応力状態になっているとする。今，図−２の傾斜部CDに沿った要素が表面から１つの要素の厚さだけ掘削される場合の解析を行う。掘削される部分と残り（掘削後）の斜面を仮想的に分離して考えると，掘削後の斜面には掘削される部分から作用していた力が外力として作用している。掘削後の形状の斜面に対してこの外力をゼロになるまで徐々に減少させると，最終的に掘削後の斜面の変形と応力状態が得られる。このときの計算法はケース

と同じである（外力の方向と大きさは異なるが）。
文献４）に示したように，斜面内の初期応力状態と掘削施工過程は，いずれも掘削に伴う斜面の弾塑性変形に大きな影響を与えるため，注意が必要である。

 斜面が動的な地震力を受けるケース

５，６）

図−６ 図ー２の斜面の弾塑性解析（破壊直前の状態）

J. of the Jpn. Landslide Soc., Vol.40, No.2159 (2003)

地震時の斜面安定問題では，力とモーメントのつりあい式は式の右辺を次のように置き換えたものになる。

５５

講 座

  H JH 
I ……………
I
すなわち土の慣性力が式 に加わる。Hは土の密度で，H K（Kは重力加速度）の関係がある。U ，V は，，方向の入力地震変位であり，t は時間を表す。
地震時においても変位とひずみの関係式は式と同じである。地震時には，応力が正負方向に繰返し作用するため，土の応力ひずみ関係はこれを反映して大変複雑になる。地盤振動の増幅作用が地震力の大きさに強く影響するため，静的な問題では重宝される弾完全塑性モデルを，地震時の土の応力ひずみ関係として安易に利用することは慎まねばならない。
動的な基礎方程式のFEMによる離散化は，静的な場合と同様に行われるが，時間軸方向の離散化に関しては差分法が利用される。
以上の４ケースはいずれも実際上重要なものであるが，斜面安定問題に限れば，ケース５６

算例の報告がこれまで比較的少ないように思われる。実際の斜面崩壊や模型実験結果を対象にした弾塑性FEMによる解析例が，今後数多く報告されることを大いに期待したい。
参考文献
１）http : //geotech.ce.gunma-u.ac.jp/ga/の FEMのいろは
２）田中忠次他（２０００）：地盤工学における数値解析入門，pp．
１２９−１４４，地盤工学会．
３）久田俊明，野口裕久（１９９５）：非線形有限要素法の基礎と応用，丸善．
４）山上拓男，鵜飼恵三（２００１）：斜面の安定と変形解析総説（LEMとFEMの応用），地すべり，Vol．
３８，No．
３，pp．
９−１９．
５）鵜飼恵三，井田寿朗，若井明彦（１９９５）：動的弾塑性FEMによる地震時斜面のすべり解析，地すべり，Vol．
３２，No．
１，pp．
８−１１．
６）若井明彦，鵜飼恵三（２００３）：地震被害予測のための動的弾塑
性有限要素法とその発展，土と基礎，Vol．
５
１，No．
２，pp．
１３
−１
５．
（原稿受付２００３年２月２７日，原稿受理２００３年３月５日）

とケースに関する計J. of the Jpn. Landslide Soc., Vol.40, No.2160 (2003)