GTSのSRM

年末より GTS でせん断強度低減法のテストをしています。

3月までお借りしているので、その間に実務への適用性を見極めようと考えています。特に3次元ＳＲＭの動作確認が目的です。

今のところ、いくつかのエラーで計算がまわっていません。サポートに伺うと、原因がまだ分からないとのこと。やはり、ソフトによってはクセや回すための細かなノウハウが必要となる場合もあるので、購入前に実際に手を動かしてテストしないといけないのでしょう。また、事前にサポートの反応を見ることも必要ですね。

それでも、いくつかの特徴はつかめました。
ソルバーによっては 32bit 版のものがあります。逆に、64bit 版や、スレッド数の（段階的な）指定が可能な並列化対応ソルバーでは、大きなモデルを扱えます。
SRM では、Fs = 1.0以上を初期値として指定し、そこから収束し難い方へ計算を向かわせる手法が取られているようです。Ga3d とは異なります。なぜ1.0以上なのか説明書には書かれていません。安全率の小さい側から少しずつ壊していかないとなかなか収束しないのでは？と思いサポートへ質問したところ、サポート側にも説明書は残っておらず、詳細は良く分からないようでした。ただ、他社のソフトと同様のアルゴであり、計算もそちらの方が早いとのことです。
ということは、説明の書かれている他のソフトを使用すれば良いということでしょうか？そういえば、FLAC も11月に Ver.7 が出ていましたね。

とりあえずエラーの件、アルゴの件で MIDAS 社の反応はつかめました。
まだ、時間はありますので、ゆっくり検討してみましょう。

不飽和斜面の安全率

久しぶりの休日です。

学会HPのチェックをしていると、地盤工学ジャーナルで面白い論文が出ているのを見かけました。

川﨑 元, 西垣 誠, 実用的な不飽和土用三軸圧縮試験装置の開発とそれを用いて測定した不飽和土のせん断強度について, 地盤工学ジャーナル,Vol. 6 (2011) , No. 1
http://www.jstage.jst.go.jp/article/jgs/6/1/39/_pdf/-char/ja/

内容は2本立て。
1つ目は通常の三軸圧縮試験機を改良することで実務的に不飽和三軸試験が行えることの報告。
2つ目は、締固めたマサ土について、Öberg と Sällfors が提案している（飽和時の試験結果から導く）不飽和時のせん断強度推定値と、上記試験機で得られたサクション一定排水せん断試験結果の比較、盛土安全率への影響度などの議論です。2つ目の結果は概要に以下のように記されています。

1) Öberg と Sällfors が提案している従来の予測法では，本装置で測定したせん断強度よりせん断強度を低めに見積もることになる。そしてその差はサクション変化によるダイレイタンシーのせん断強度への寄与分の差と，湿潤過程におけるサクション 0 の不飽和ケースの粘着力と飽和ケースの粘着力の差で構成されている。

2)  従来の方法で推定したせん断強度と本装置で測定したせん断強度との差の影響は湿潤過程にある盛土の安定解析上，無視できない

面白いですね。つまり、少なくともサクション0の試験をして予測式に使用する必要がありますよ、と。
それで出るなら簡単ですね。試験機については本当に簡単に改良できるのか、プロに聞いてみましょう。

サクションによる安全率の変化について、粘着力の寄与分、φの寄与分を述べられています。あまりに明確な結果であり、面白いですね。斜面では降雨浸透により見かけの粘着力が減少し崩壊が起きると言われますが、その通りなのでしょう。

拘束圧に応じてポアソン比は不変ですが、変形係数は異なります（=せん断剛性が異なります）。では、ダイレイタンシー角はどうなんでしょう？結果は変わっていますね。低拘束圧の方が正のダイレイタンシーが出やすいというのは感覚として納得できますが、深く考えたことはなかったですね。注意すべき点です。

安全率の評価ではGa3dのSSR法が採用されています。ピーク強度とダイレイタンシー寄与分考慮後の強度の2つで分けて計算されています。拘束圧依存性も反映できるよう、モデルを2層としています。この程度の反映は実務でも必須なのでしょう。

前者のケースでは実際に比べて過大な安全率が出ますし、後者のケースでは過小な安全率が得られものと推定されます。降伏後にピーク強度からダイレイタンシー寄与分考慮後の強度（≒残留強度）に落とせば、それらの安全率は中間になると思います。試験結果では緩やかに落ちていますので、私が岩盤を意識して改良したコードは使えませんね。FLACなどひずみに応じてc・φを設定するタイプのひずみ軟化を扱えるコードを用いて計算した方が良いのでしょう。

全般的には実務に対して2段階程度上を行く内容です。
一つは、安定計算では飽和度を無視していること。2次元・3次元、あるいは順解析・逆解析においても、すべり計算において不飽和帯の強度を考慮することはありません。港湾や軟弱粘土では、ほぼ飽和として扱うので問題ありませんが、盛土や地すべりなどでは無視です。つまり、試験値を使ったとしてもそれは飽和の試験結果であり、計算上は不飽和も含めた平均強度として扱っています。不飽和斜面の崩壊も実務ではオーソライズされたものがありませんので、なかなか取り込めないのが現状でしょう。
もうひとつはSSR法の使用です。斜面安定に関する数値計算は土研（トンネルと地すべり）や盛土工指針（浸透流）でオーソライズされた感はありましたが、SSR法はさらに上を行く話でしょう。

いろいろ考えさせられる論文でした。
時間のあるときに出会えてラッキーでしたね。

支持力係数

反復回数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回）の強制変位を与え、得られた接点力を合計し、極限支持力をチェックしました。結果は以下の通り。

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

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

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

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

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

コンパイラーによる計算結果の差

Ga3dにひずみ軟化を取り入れ、コンパイルしてみました。

最初は、Win7 64bit ＋ Intel Visual Fortran Composer XE 2011 で行いました。
ところが、計算が進みません。おかしい。

originalのソースをコンパイルしても同様。バンド幅の計算で誤った答えを出して止まります。

色々試しましたが、うまくいかないのでWinXP 32bit + Intel VF9 でコンパイルしてみました。
結果、配布されている original の exe （コンパイラーは Compaq）と答えは一致しましたが、outファイルの0に近い値で10^-13程度の誤差が生じます。これは仕様上、仕方ないところでしょう。

結局、このソースでは以下のような動作確認結果となりました。

 OS     exe   結果
32bit  32bit   ◎（答えは一致）
64bit  32bit   ×
64bit  64bit   ○（答えが32bitと微妙に違う箇所あり。影響のない程度）

コンパイラーによって、結果に影響のない程度の差が出てくるのは仕方のないことかもしれませんが、同じコンパイラーなのに64bitと32bitで差が出てくるのはいやらしいですね。まあ、PCを使用した数値計算の限界や収束設定値を考えて、誤差が影響のない程度であれば問題ないわけですが。

最適化オプションによっても微妙に答えが変わってくると思いますので、こういったコードの配布時にはプロジェクトのプロパティ―（セッティング）等の情報も参考程度に公開していただきたいものです。

Ga3dでひずみ軟化

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

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

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

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

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