2次元弾性圧密解析

圧密解析を弾性解析で解くことがあることを、先輩から教えてもらいました。

ソースは以下の通り。 実務的だなあと感じましたね。
土質工学会「地盤工学における数値解析の実務」
29年前の古い図書です。

全応力のE50から有効応力のEへの変換は、どこかで見たことがあるなあと思ったら、コレでした。
http://phreeqc.blogspot.jp/2015/04/blog-post_47.html

圧密試験のmvからEへの変換は以下の通り。
平面ひずみ状態の仮定ですから、3次元にそのまま使えるかどうかはまた別のお話。以前の講習会では気にしなくてよいとも伺ったことがありますが。

これ、知りませんでした（有名な式ならどこかで見てるかもしれませんが）。調べる必要がありますね。

  先輩がおっしゃるには、やはり合わせ込みが必要とのこと。何に合わせるかと言うと、2次元圧密計算。本末転倒の様な気がしますが、応力状態を求めたいので仕方ありません。

しかし、応力が必要であれば、塑性も考えなくてはダメでは？と思います。 弾塑性でしょうね。

ロッキング

多層地盤の支持力をチェックしようと SoilPlus 2012 で計算しました。

モデル自体はDXFを取り込めますのですぐにできます。2次元では平面ひずみ1次要素しか対応していませんので、しかたなくそれでモデル化しました。
が、計算すると、何かおかしい。

よくわからないので、地盤工学会「3.弾塑性有限要素法をつかう」にある簡易なモデルでチェックしてみました。結果、理論解より1オーダー高い値が出ています。変位と支持力の関係をプロットすると、ほとんど弾性の挙動でした。
そんなものなのかと思いながら、3次元のソリッド要素で、境界条件により平面ひずみ状態を作り2次元的な解析を行いました。すると、理論解よりやや高い程度の値が出てきました。

原因がよくわからず、社内のプロやサポートに聞きましたが、満足な回答が得られません。モデルをチェックしてもおかしなところはありません。そんなこんなで3夜ほど考えていました。

まさか？と思い、2次元モデルでせん断成分の減退積分を行ったところ、オーダーが合致しました。

ロッキングでした。

しかし、ソリッド要素の結果に比べ、まだ値が高いので、体積成分に対し減退積分を行いました。結果、ソリッドの結果よりも低くなり、より理論解に近づきました。が、アワーグラスが発生しておりダメ。しかも、時間がかかる。
メッシュを細かく切るしかないのでしょうね。


今回のように、明らかにおかしいと思える場合は良いのですが、ちょっと違うと気付かないこともあるでしょう。そんな心配やチェックの手間を考えるなら、SoilPlusでは2次元モデルは扱わず、ソリッド要素を使うしかありません。平面ひずみならDXFを読み込んで1m押し出すだけなので、手間はそれほど変わりませんし。

ソフトにはクセがあります。今回の場合はクセというか仕様？なので私の理解が足りなかったのが原因でしょうけど。まあ、どんなソフトも中身と自分のレベルをよく知ることからですね。

ロッキングとその回避、復習しておきましょう。

孔内水平載荷試験の再現 その4

以前、孔内水平載荷試験の結果を数値計算で再現すれば、ダイレイタンシー角ψを同定できるのでは？と考えたことがあります。その時は、SoilPlusが関連流れ則しか扱えなかったため、同定作業自体ができませんでした。
http://phreeqc.blogspot.com/2011/08/3.html

今、MIDAS社からGTSをお借りしていますので、そちらで確かめてみることにしました。プリはSoilPlusもMIDAS社から提供を受けているので、操作性は同じです。ただ、こちらは軸対象で水平方向の荷重を扱えますので、2次元で計算できます。3次元より計算時間が早いので、すぐに結果がわかります。

計算結果はこんな感じです。ゾンデ60cm間に1000kN/m2を10分割でかけました。

変位のグラフの例です。

ゾンデの先、余掘り分15～20cmをモデル化し忘れたので、メッシュを切りなおそうとしました。が、同じ領域でメッシュを切りなおした場合、前の節点の情報をクリアしないようです。バグでしょうか？
仕方ないので、そのまま比較することに。

比較は弾完全塑性の関連流れ則、非関連流れ則、ひずみ軟化で行いました。
パラメータはダイレイタンシー角ψと単体以外、試験での実測値や算出値を使用しています。パラスタはしていません。変形係数は試験でポアソン比を仮定して得られた値、K0は土被り圧とP0より算出、φ'peakは深川の方法、せん断抵抗力(ピーク、残留)はσr-σθで求めた値を使用しました。なお、ψは仮にφ-30°として設定してみました。
ピーク、残留をφで指定できないのは、ソフト上の制限です。実務でGTSのひずみ軟化を使用する場合、砂質土をモデル化するには、深度毎に分割しないといけないようです。

で、結果は以下の通り。

うーん、弾性変形の傾きがやや立っていますね。試験時のポアソン比の設定が良くなかったのでしょう。逆に、試験でGがわかるので、パラスタで適切なポアソン比(=変形係数)が設定できそうです。

この結果では、ひずみ軟化が最も合っていますね。試験で求めた値を使うとばっちり合うのでしょうか？結論を得るにはもっとやってみる必要がありそうです。非関連も、まあ、許容範囲でしょうか？関連は駄目（つまり、SoilPlusでこの土の解析はダメ）ですね。

ちなみに、これは載荷区間の変位の平均値をグラフ化しています。実際の孔内水平載荷試験も、土砂の場合は水の注入量から平均的な変形量を算出していますので、こういった簡易な整理も良いと思います。

載荷区間の下の方は、余堀をモデル化していないミスで変位が抑制されています。試しに下15cmをオミットして平均をとると、以下の通り。

それほど変わらなかったですね。ひずみ軟化が良さそうです。

除荷、再載荷をしてみました。

感動！ちゃんと除荷時に変位の大きい側で低下し、負荷時に小さい側で上昇するループを作っています！実際の土もこういった挙動になるのですが、理由はちゃんとあるんですね。構成則を作った昔の人は偉い。すばらしい。

結果としては当ケースの場合、ひずみ軟化が良いでしょう。ψの同定を目的に始めましたが、思わぬ収穫でした。孔内水平載荷試験の結果より、適切な土質定数（の範囲）を設計側に助言することは可能ですね。

うーん、しかし、この結果はひずみ硬化ですよね？弾完全塑性やひずみ軟化としてモデル化しても、結果はひずみ硬化のような挙動を示しています。硬化と軟化、弾完全塑性の区別が良く分からなくなってきました。

まだ理解していない証拠です。もう一度、勉強し直す必要がありそうです。

Cam-clayで砂と粘土を統一？

ここ最近、Cam-clay model の理論を、現場往復の車内で聞いています。

このモデルを初めて聞いた時、e-logp曲線を弾塑性状態の硬化として解釈していることに目から鱗でした。偉大です。これは思いつくようで思いつきませんね（私だけ？）。土質力学の教科書の弊害でしょうか？圧密は一次元、強度は三軸試験というように、別々の単元で論ぜられています。統一的に解釈できることが実務上良いのかどうかは分かりませんが、確かに、カムクレイの方が頭を整理しやすいですね。
修正CC、関口太田モデルなどは、そこから差が理解できれば到達できますので、最初はオリジナルCam-clayでしょう。まだ完全に理解できていない細かな箇所が残っていますので、もう少し時間をかけてみましょう。

聞いている中で、少し気になったことがありました。

「これ、粘土だけ？」

孔内水平載荷試験を行うと、粘土、砂、岩盤に関係なく、必ず弾塑性的な挙動を示します。しかも硬化。土であれば10%以上の孔壁ひずみを示すことが多く、ほぼ100%硬化を示します。繰り返し載荷をすると良く分かります。応力-ひずみ曲線だけ見れば、やや乱暴ですが圧密試験結果と良く似た形です。

しかし、硬化パラメーターを塑性ひずみで同じように表すことができるのかどうかは分かりません。破壊時の排水状態は結果から推定することは可能ですが、有効応力までは困難でしょう。そのあたりのパラメーター（引数）を何にするかはもう少し理解しないと議論につくことすらできませんが、偉い先生方なら砂と粘土を統一的に説明することができるような気がします。

と、ここで思いだしたのが、過去に上司から勧められた解析プログラム「GEOASIA」。All Soils All States All Round Geo-analysis Integration だそうです。凄いですね。

http://www.geoasia.jp/about/index.html

当時は興味がなかったので、そのまま放置でした。しかし、砂・粘土など材料に無関係で、有効応力解析が可能、しかも、「SYSカムクレイモデル」を使用しているそうです。SYSが何を意味するのか知りませんが、上記のような話でしょうか？急に興味が沸いてきました。

いずれにしても、もう少し理解しないと議論の場に上がれません。

まずは、基礎力を付けることを考えましょう。

支持力係数

反復回数1000回でも、結果は全く変わりませんでした。

ｃ＝10kN/m2で固定した場合の結果は、以下の通り。

値はOKです。やはり、粘りがないですね。施工時の観測では困るでしょうね。弾性域は200kN/m2程度でしょうか？

これを含めた検証結果より、改変コードは正しく動作していると言えるでしょう。

ところで、支持力係数は道示ではグラフからの読み取りになっており、式が示されていません。古い土研資料が根拠のようですが、持っていないですね。
http://jsce.jp/pro/node/1212
今回は、建築基礎構造設計指針の Nc, Nq 式を使いました。

支持力係数の算定式はいくつかあります。式やグラフは、いろんな参考書に載っています。
技報堂出版「土の力学」は10年以上前の講習会で購入した教科書ですが、そこに係数のグラフが載っています。当時、そこに書き込みをしていました。

「Terzaghi は理論的に少し違う。Caquot・Kerisel が正解。」

極限支持力算出式で、根入れ分の土の荷重は見るが、強度を見ない所が理論的でない（安全側）といったことは覚えています。しかし、支持力係数の「理論的に少し違う」は覚えていないですね。何だったのでしょう？

理工図書「土質力学-全訂新版-」でも式が以下のように書かれています。

「上述の解析は本来、有効応力法であるクーロンの破壊基準に基づき、またプラントルの支持力破壊の手法を踏襲しているので、支持力の解析も有効応力法であるが、第4章で述べたような全応力法による見かけのせん断強さを用いても差し支えないことは、土圧や斜面の安定の場合と同様である。それらの結果に差があるのは現在ではやむを得ない。」

これは1991年の4版なので、ちょうど20年前ですね。20年経ってもあまり変わっていないです。

Terzaghiの式による比較

Ga3d改変コードの検証の続きです。

Terzaghiの式では、Nc, Nq, Nγ の3つの支持力係数を使用します。このうち、理論解が求められるのは前２つですので、γ=0として（Nγの影響を取り除いて）比較しました。

ピーク強度 c=10kN/m2、φ=30°kN/m2、残留強度 c=１kN/m2、φ=25°kN/m2 です。
結果は以下の通り。

値は合っていますが、いきなり破壊するんですね。粘土の方が降伏してからも、文字通り粘りが出るんですね。

ただこれは「ある収束値での修正NRの繰り返しが○○回以上で発散」という人的な設定値によるものも大きいと思います。今は500回なので、1000回ならもう少し伸びるでしょうか？あるいはcを落とさずにφだけ落としたらどうなるでしょう？破壊時におかしな値になっているのも気になります。

もう少し、やってみましょう。

プラントル解との比較

Ga3d にひずみ軟化を取り入れた改変コードの検証を行っています。

最初はFLACの例題と同じモデルで検証しようと思いましたが、モデル作成が難しいので止めました。薄々感じてはいるのですが、この改変コード、プレ処理を工夫しないと実務では使えませんよね。ま、完成してから考えましょう。

簡単なモデルで検証！ということで、支持力問題にて検証しました。ex11ａ.dat での比較です。

まずは c =100kN/m2、φ=0の粘性土をモデルに与え、残留強度として70kN/m2に入れ替えた場合を計算してみました。実現象としてあるかどうかは別として、あくまで粘着力 c 入れ替えの動作確認が目的です。

メッシュ幅は0.5mに変えました。これは、メッシュが細かいほどプラントル解に近づくという理由からです。本にも書かれています。

100mm（1mm×100回）の強制変位を与え、得られた接点力を合計し、極限支持力をチェックしました。結果は以下の通り。

せん断ひずみの形は変化ないですね。

当然、変形量も変わりません。

しかし、極限支持力は大きく変わりました。青はピーク強度のみ、赤は残留強度考慮。

プラントル解より少し大きめに出ています。まだメッシュが大きいのでしょう。まあでも、良いところへ行っていると思います。降伏点は変わりませんが（これが必要）、極限支持力ではひずみ軟化モデルのほうが大幅に低下し、それぞれのプラントル解へ近づいています。成功ですね。

次はφのチェックをしましょう。

孔内水平載荷試験結果の解釈と利用

孔内水平載荷試験について引っかかっていたことがありました。

処女載荷の勾配は必ず繰り返し載荷より緩くなります。
ひずみレベルの差として解釈していましたが、前者より同じレベルの区間を取りだして見ても勾配に差があることに違和感がありました。また、繰り返し載荷と同じ勾配が処女載荷に一度も出てこないことも引っかかっていました。
ひょっとするとこれは弾塑性状態にあるのではないか？早々に降伏し、ひずみ硬化の状態にあるのではないか？などと考えたのが数ヶ月前。ちょうど連続体力学や構成則を追っかけていた頃です。

先日の講習会で、引っかかっていた点を講師の方に質問してみました。
回答は、弾塑性状態でOK、押し込んでいるので繰り返し載荷時は勾配が立つとのことでした。ひずみレベルよりは解釈しやすい御回答です。それが実際に起きている現象なのかどうかはまだ分かりません。岩盤の場合はどうなのか？他の問題や扱うスケールとの関連性は？など、まだ完全に理解できていませんので、今後ゆっくり考えて行きましょう。

また、どちらで導いた変形係数を使うかは、扱う問題によって（答えを想定し）使い分けることが重要とのことです。その通りだと思います。基本的に、土砂の孔内水平載荷試験結果の解釈については割り切っています。
http://phreeqc.blogspot.com/2011/01/blog-post_08.html
弾塑性状態から弾性の変形係数を求めているのなら、更なる割り切りが必要になります。割り切った上で「適度な」答えを導けるようになるべきなのかもしれません。

いずれにしても、FEMで変形量を求める場合には、c、φ、ψ、E、νなど、全て工学的（経験上の）判断が重要になるということです。純粋に、精密な試験を行って得られた値を入力すれば、答えの精度が高まるといったものではないのです。
こうなると、経験豊富なベテラン技術者の方が、数値解析に向いているような気がします。
あ、だから、ベテラン技術者はシンプルなツールで答えを出そうとするのか、答えがある程度見えているから。答えが見えない、予測のつかない問題は、シンプルであろうが複雑なツールであろうが、難しいということでしょう。答えを出すのはFEMではなく、技術者ということです。

Ga3dでひずみ軟化

地盤工学会「弾塑性有限要素法をつかう」の中に、Ga3dという弾塑性（静的）解析コードが収録されています。

先日の講習会で空いた時間に、そのフローとソースを見ていました。シンプルであり、書籍中にコメントがあるので、理解しやすいコードです。

ある程度眺めていると、ふと気が付きました。「これ、簡単に残留強度に入れ替えられるんじゃないか？」と。つまり、降伏後に残留強度に落ちるひずみ軟化モデルを容易に反映することができると思えたのです。

講師に確認すると、それでOKとのこと。早速、アルゴリズムを作成し、翌日確認していただきました。合格のようです。後日、FLACの例題と比べてみたり、収束性を見てみたいと思います。

弾塑性計算のアルゴリズムで必ず通るはずの箇所に手を加えていますので、他のソフトでも同じような改変が可能と思います。ただ、市販ソフトは自分で手を加えられないんですよね。プレの使い勝手をとるか、ひずみ軟化の必要性をとるかは問題に応じてといったところでしょう。

弾塑性FEM

地盤工学会の「わかって使うFEM」講習会に参加していました。
弾塑性（静的）やFEMの理論は理解していましたので、復習といった感じです。

今回の目的は静的ではなく、動的、圧密解析でした。これらの理論もパーツとして理解していましたが、まだぼやっとしたところがありました。今回の参加でそれらをつなぐことができ、収穫ありです。それほど難しくはないんですよね、理屈だけなら。まだ完全ではないので、後日、再度数式を追いかけましょう。

ただ、理論を理解するのと、答えが評価できるのとは大きな差があります。その辺はプロの領域ですね。地質屋として求められるのは、ある問題にある理論で望みたいと言われたときに、どの程度の設定値（幅）が必要だから、どのような調査を提案すれば良いかを判断するといったところでしょうか？

今回の講習ではパラメーターの扱いについて詳細なお話がありました。前日に質問していましたので、講師の配慮もあったともいます。ありがたいことです。今後は数値解析について一歩引いた目線で扱えそうです。以前、FEMのベテラン技術者に言われたことがありましたが、ようやくたどりついた感じです。

また、講師との雑談や質問に対する回答もかなり有意義でした。
そのあたりは、また後日。

The Verification Problems of FLAC3D

FLAC3Dのマニュアルを読み始めました。
ソフトの特徴がわかりやすいのは、やはり実際の問題を解いている章ですね。

Verification Problems の章では、理論との検証を行っています。まず引き止められたのが、最初の検証 Cylindrical Hole in an Infinite Mohr-Coulomb Material です。マニュには以下のような概要が記載されています。

Stresses and displacements are determined numerically for the case of a cylindrical hole in an infinite elasto-plastic material subjected to in-situ stresses. The material is assumed to be linearly elastic, perfectly plastic, with a failure surface defined by the Mohr-Coulomb criterion, and both associated (dilatancy = friction angle) and non-associated (dilatancy = 0) flow rules are used. The results of the simulation are compared with an analytic solution.

これ、理論解があったんですね。知りませんでした。
孔内水平載荷試験（岩盤）では弾性域で厚肉円筒理論が使われていますが、ここで紹介されている式は塑性領域まで含まれており、導出が異なっているようです。単純な問題ですが、手法が複数あるというのは面白いですね。

また、ここでは理論解との検証が目的で計算を行っています。結果は、当然理論解と近い値を示していますが、関連・非関連流れ則で応力-ひずみ曲線が異なっているようです（グラフしか載っていませんので、正確な値は不明）。以前書いたように、孔内水平載荷試験結果を再現し、ψを同定するというのも、目的は異なりますが同じ手法で良いでしょうし、FLAC3Dで可能ということでしょう。

ただ、理論解があるなら、ψの同定はそれで十分ですね。
文献、取り寄せてみましょう・・・と思いきや、仏語。
どうしましょう。
ま、数式を追いかけるのがメインなので、なんとかなるでしょう。

ダイレイタンシーと斜面安定度

注文していた論文が届きました。

ひとつはハズレ、もう一つはまあまあ、といったところです。
後者はダイレイタンシーと斜面安定度（いわゆる安全率ではなく、荷重が増えた場合にどこまで持ちこたえられるか、といった定義）の話です。同じ摩擦角であれば、ダイレイタンシー角が大きいほど斜面安定度が増すといった内容です。 非関連流れ則よりも、関連流れ則(ψ=φ')の方が安定度が高くなります。
Majid T. Manzari et al., Significance of Soil Dilatancy in Slope Stability Analysis, J. Geotech. Geoenviron. Eng. 126, 75 (2000)

過去の地すべり学会誌はJ-Stageで公開されています。さらに古い18年前の鵜飼 恵三ほか「モデル化された切土斜面の3次元安定解析と留意点」では、そのメカニズムも解説されています。そのもとになっているのが例の地盤工学会「土の強さと地盤の破壊入門です。」モールのひずみ円～流れ則のページですね。ただ、φmobと応力傾角μが区別されていない点と、不動層の節点を固定している点が気になりますが。
http://www.journalarchive.jst.go.jp/jnlpdf.php?cdjournal=jls1964&cdvol=29&noissue=4&startpage=18&lang=ja&from=jnlabstract
この報告、関連流れ則でないとSSRによる2次元安全率がLEMと同じ値にならない点を（モールのひずみ円と応力円を念頭に）指摘しています。一方、関連流れ則では大きな膨張が生じるため、変形が現実的でないとも述べられています。非常に納得できる解説で、頭の中がすっきりしました。その後の論文も探して変更点がないか洗ってみましょう。

これらと似たような話は地盤工学会「地盤工学における数値解析入門」にも出ています。こちらは支持力の話です。FEMで地すべりや支持力など塑性化を取り扱う場合には、ダイレイタンシー角について注意すべきなのでしょう。

ダイレイタンシー角-変形量（塑性化領域）-安全率の3つを考慮する必要があること、適切なダイレイタンシー角(実務的には問題に応じて)やソフトを選定する必要があることがわかりました。

これで調査方針も決まりですね。

ひずみ軟化が難しい

いえ、低レベルな話です。

地下空洞の安定性や地すべりを扱う場合、破壊後にピーク強度から残留強度に落とす、ひずみ軟化の解析例をよく見ます。特に前者では当たり前のように使用されています。
先のダイレイタンシー同様、変形量やせん断帯の形成に効いてくるようですので、調査時にどの程度気にすべきかテストケースで理解しようとしています。FLACでは残留強度への切り替えを扱えることを、以前教えていただいたのですが、そのレベルに行くまでに使い方が全く分かりません。とりあえず手元にあるSoilPlusでどうなるかを試してみようと考えました。

が、SoilPlusではφ'rの入力ができません。ひずみと硬化・軟化係数、もしくは相当応力との関数としての入力だけです。イメージではわかるのですが、正確な表現手法がわかりません。止まってしまいました。変形のプロなら簡単なのでしょうね。昔、プロがたくさんいらっしゃった環境で、なぜ勉強しなかったのだろうと、今になって反省です。

今回、ひずみ軟化について過去の文献をあさっていると、剛性を落とすもの、cを落とすもの、cφ両方を落とすものなど、いろいろな提案があることに気付きました。また、残留応力状態への応力再配分の方法も、いくつかあるようです。何が主流で、どのように結果が異なっていくのか詳細がわかりません。まだまだ理論も含め理解の幅を広げる必要がありそうです。
地質や観測データから、こういう変形過程を再現したいので、こういう解析が必要（という方針を解析屋さんと話せて）、だからこういう調査が必要！というところまで、できるようになればBESTですが、まだまだです。

ここを乗り越えると、目処がつきます。
年末までには何とかしましょう。

孔内水平載荷試験の再現 その3

CTCさんに変形量が変化しない点について問い合わせてみました。
すぐに回答がありました。

ダイレイタンシー角については大規模オプションのみ反映される項目で、ノーマルのソルバーでは値を入力しても反映されませんとのこと。

・・・マニュアルに書いておいてほしいですよね。

しかし、ここで疑問が。ダイレイタンシー角を入れないと、塑性ひずみが分からない筈です。どうやって決めているのでしょうか？これも回答がありました。

内部的にψ=φとして処理しているとのこと。
関連流れ則を採用しているようです。

MC-MC、DP-DPモデルは可能ですが、MC-DPモデルはダメです。

結局、現段階ではソフト上確認できないことが分かりました。
まあ、こんなものでしょう。他のソフトではどうなんでしょうね。
方針がダメだったわけではないので、いずれ確認しましょう。

孔内水平載荷試験の再現 その2

昨日の3次元モデルをダイレイタンシー角だけ変更して、試算してみました。

40分ほどで結果が出てましたので確認すると、変更前と全く同じ変形量です。設定を間違えたかなと思い、確認しましたが、ダイレイタンシー角はきちんと変更されています。
おかしい、もう一度と思い、S+を立ち上げなおし、再度計算。しかし、結果は変わらずでした。

以前、ダイレイタンシー角を変えることで変形量に差が出た2次元モデルをもう一度計算すると、確かに変形量が変わります。なんでだ？と思ってパラメーターを確認すると、材料値も変えていました。情けない。

それを修正して計算し直すと、3次元モデル同様、変形量に差がでなくなりました。適当なダイレイタンシー角を入れても、まったく変形量が変わりません。私の設定がおかしいのかもしれませんね。

まあ、効かないなら効かないで、お金をかけて調べる必要もないことが分かったわけですが、どうも納得できません。構成則、ちゃんと理解しないといけないですね。

盆明け、プロに確認してみましょう。

孔内水平載荷試験の再現

昨日思いついた孔内水平載荷試験の再現シミュをやってみました。

目的は、ダイレイタンシー角ψの同定です。
試験結果より得られたE,φ',K0,ν,Gを固定し、ψを変えて加圧に応じた変形量を求める方法で実施しました。この試験では、主に平面ひずみ状態での上記ψ以外の値を求めることが可能です。使用したソフトはSoilPlus2011です。

始めは、軸対称条件で、内圧をかけるモデルを作成し、実施しました。掘削過程やケーシングによる変形の制御も施工段階として取り入れました。ただ、押し出してきた地盤をケース設置時に削り取る方法がわかりませんでした。とりあえず載荷箇所ではないので0クリアで済ませましたが、プロはどう対処しているのでしょうか？
モデル自体は2時間ほどで組めたのですが、さあ、計算！となったとき、エラーが発生しました。
SoilPlusの軸対称モデルでは、鉛直方向の荷重しか扱えないようです。残念ですが、やり直しです。

気を取り直して、3次元でモデルを組み直しました。といっても、回転押し出しで済ませばよいので、それほど手間はかかりません。メッシュを切ると1万節点強のモデルになってしまいましたが、テストなので、まあいいでしょうと思い計算させてみました。が、古いPC（Centｒino Duo T7500 2.2GHz）でやっているため、結果が出るまで40分強かかりました。モデルは↓です。

最初はψ=φ'-30=2.9°で実施してみました。
最終ステップの結果をみると、載荷箇所の周辺地盤が塑性化していることから、この点では合格です（実際の試験でも破壊させています）。局所安全率が1を大きく切っている箇所があることはどうなんでしょう？弾完全塑性で0.8とかありなんでしょうか。計算として不合格？これもプロに聞いてみましょう。
問題の変形量は樽型になっているので直接比較はできません。が、壁面付近は全体的に6mm強なので少ないですね。その倍は欲しいところです。（色はデフォルトのままですが、↓のような感じ。）

試験では均等に変形していると仮定して変形係数を算出しますが、理想的にはこのような樽形で算出すべきなのでしょう。キャリブレーション時のゾンデの膨れ具合をみると、もっと複雑な変形をしていることは容易に想像できますが。
この結果を検証するには体積変化しかないなあと思いつつ、それは後回しにしました。

先に、気になった未掘削部の対称軸付近の変な塑性化がno-tension解析で解消するかを試してみました。これは計算として不合格ですから。しかし、これがドツボの入り口でした。
計算が収束しません。no-tension解析だから収束しないというのはおかしいなあ、と思いつつ、それしか原因がないので元に戻し、ψ=32.9°で回してみました。結果、収束せず。その後、値を変えても、ソルバーを変更しても駄目でした。すべて載荷時で引っかかります。結局、何が悪いのか分からずに元の値に戻して計算しましたが、これも駄目。最後はソルバーすら立ち上がらなくなりました。

気を取り直してS+関連の生きているexeをタスクマネージャーから殺し、再度立ち上げなおすと、すんなり動くようになりました。
載荷時のステップを細かくとったのが良かったのか、収束計算対象を変位に変更したのが良かったのか計算も収束します。

今日1日は何だったのかと思いつつ、まあ、方針は悪くなさそうなので、明日も懲りずに試してみることにしましょう。

ダイレイタンシー角

お盆休みの間に、変形解析の実務上の（調査に関わる）問題点を拾い出しておこうと思い、先日よりSoilPlusで2次元変形の試算を行っていました。

弾塑性でも完全塑性、ひずみ軟化、硬化が選択できます。完全塑性でも、通常扱わないパラメーターが必要です。難しいと感じたのは、ダイレイタンシー角、硬化係数でしょうか？

ダイレイタンシー角ψは降伏基準としてMCを選択した場合に必要です。前に紹介した本には砂質土としてψ=φ'-30°の関係が掲載されていますし、一般的にはψ=φ'として計算されることも多いようです。粘性土ではψ=0が多いでしょうか？
ダイレイタンシー角ψは三軸試験でも得られますが、この試験ではε2≠0なので物理的意味は明確ではありません。
http://phreeqc.blogspot.com/2011/08/blog-post_05.html
この値、結果に効いてきます。試算で使ったモデルでは、局所安全率が大きく変わってしまいました。(2011/8/15追記：材料値も変えていました。)
平面ひずみ圧縮試験でポアソン比とともに求める方がよいのでしょうか？もう少し、理解する必要がありそうです。

試算しながら思ったのですが、変形係数を求める際に実施する、孔内水平載荷試験の結果（繰り返し載荷）を再現してやれば、これらの値は同定できる可能性が高いのではないでしょうか？軸対象平面ひずみ状態ですが、上記の一般値よりは、ましな値が得られるのではないでしょうか？

既にどなたかされているかもしれませんね。文献、探してみましょう。

連続体力学

地盤工学会の「地盤の連続体力学入門講習会」を受講しました。

3人の先生による連続体力学の基礎理論、構成式、微小変形論の適用性でした。理解の弱い変形理解の基礎となればと思い受講たのですが、ちょうど良いレベルでした。久しぶりに頭を使いましたね。

基礎理論はOKです。これは特に問題はありませんでした。

構成式は理解度50%程度というところでしょうか？時間の都合上、式の展開はほぼ省略されましたので後日追う必要があります。私にはちょっと難しかったですね。弾塑性のみでしたが、あとは参考書で理解を深める必要があります。講義の最後にセメント量や飽和度を指標とした強度特性の反映も紹介されましたが、構成式の変更は割と単純でしたね。ただ、コードを自作するしかないようで、実務ではOUTです。もちろん、技術論ではありません。有名な先生のお墨付きか、数多くの実績があるコードでないとお客様が嫌がりますから使えないのです。仕方ありません。

微少変形理論は分からないところの方が多かったですね。難しかった。ただ、疑問に思っていた箇所、たとえば、「体積ひずみεv=ε1＋ε2＋ε3を定義した場合、定体積せん断でも体積が減少（過小評価）する理由」などに対し分かりやすく御説明頂きました。
他にも、「テルツァギーの圧密理論の適用範囲」や「段階載荷と一括載荷での自重解析時の変形モードの違い」などを連続体力学を用いて示されました。数式は全く追えませんでしたが、面白い内容でした。

今回は自分に足りない箇所が具体的に理解できました。整理すべき箇所が山ほどあって時間が足りないですね。
こういった基礎力は講習会だけで身に付く物ではありません。一つ一つ式を追っていって、仮定なり適用範囲なりを理解する必要があります。時間をかける必要があり、王道はありません。提示頂いた参考書を購入し、足りない箇所をこつこつ補っていきましょう。

トンネル掘削と透水係数の変化

COMSOLのConference CD 2011 Ed. の中に、トンネル掘削による影響によって、周辺の透水係数の変化を考慮した浸透流計算の実施例がありました。2010年 Bostonでの発表のようです。

原理は簡単で、弾塑性でひずみを計算し、それに応じた透水係数の変化をあらかじめ設定していた関数から導いて、浸透流計算をするといった内容でした。関数の出典も載っていましたので、取り寄せましょう。

昔、ポストグラウトで変位を出してしまった事例を経験しましたが、その時は水の解析は水理担当者、変形の解析は力学担当者として別々で行っていました。本来はこういった連成解析が必要だったんでしょうね。

COMSOLは柔軟ですね。

数値解析による地すべりとトンネル影響評価手法

標記の資料が土研から（案）として出ています。
http://www.pwri.go.jp/team/landslide/outcome.htm

今日はこの資料を読んでいました。

実はこれに似た内容をH18の技術士の体験論文で書いています。（見事不合格）
ずいぶん前からある考え方ですがオーソライズされたものは見ませんね。

内容は、地すべりブロックとトンネルの安全な離隔距離を、2次元弾塑性解析の変位量で判断したというようなもの。パラメーターの設定法といった趣旨ではなかったですね。弾塑性は市販ソフトで良いのがありますし、2次元ならモデル化・解析ともに簡単なのでマニュアル化するほどのことではないとも思いますが、確かに各機関にオーソライズされたマニュアル類はあまり見ませんね。

こういった地すべりとトンネルといった問題は、弾完全塑性ではなくひずみ軟化、あるいは弾塑性＋塑性領域のみ残留強度としたLEMなどが良い手法だと思います。土研さんが弾塑性をオーソライズした形になりましたので、今後どうなっていくか楽しみです。