物理演算と粒子法

3ds MAX の massfx, mparticles を触ってみました。

触ってみて初めて気づいたのですが、CGの分野ではアニメーションを作るのに、手で（指示で）ものを動かすのですね。そりゃ、物理演算が可能であれば、そちらの方が容易にアニメーションを作成できるのでしょう。やはり、この分野のほうが、進むはずです。

以下の理論も、途中から理解できなくなりました。動画は今年4月ごろに話題になっていたようです。別の動画には GTX 580 レベル（張り付けた動画は680のようです）でリアルタイム演算ともありましたが、すごいですね。
http://mmacklin.com/pbf_sig_preprint.pdf

下側は Rigid Body との組み合わせですが、土の物性を与えた粒子に置き換えておけば、土砂の流れもシミュできそうです。

土石流や深層崩壊の場合、最初は土や岩石で、崩壊後はほぐれた土砂が流体のように移動し、停止すれば土（堆積・天然ダム）になり、水は流れ去る、といった mass としての物性の変化（各相の物性自体は変化していません）があります。私は土石流や天然ダムのプロではないので知らないだけなのかもしれませんが、それらを1相として計算する場合、土砂の崩壊や落石の集合体に近いイメージで取り扱っているのでしょうか？

mass として簡略化するよりは、固相に水を与えて崩壊させるなど忠実にモデル化する方が、結果的には近道なのかもしれません。

ま、私には遠いゴールです。少しづつ、追っていきましょう。

粒子法 その2

最近、また粒子法について、手を出そうかと思い始めました。

2年前、概要は把握しましたが、それ以降は手つかずの状態でした。
http://phreeqc.blogspot.jp/2011/08/blog-post_12.html

で、再度、教科書を読んでいるのですが、やはりそれ以上はダメ。
SPH の基礎理論はわかります。数式も追えます。MPS は概要のみ。数式の一部がまだ理解できていません。離散化はOK。でも、それを実装するのは全くダメ（CPUのみなら、FEMより簡単だと思うのですが）。そんな状態です。要は、理解していないということです。

これ、ゴールの遠いところが、腰の重い原因の一つでしょうね。
計算以外にも、GPGPU やメニーコア対応は当然の分野。ポスト処理としての陰関数曲面の作成や、レイトレーシングなど、数値計算とは異なった知識も必要です。手に入る情報は簡単か難しいの両極端。市販ソフトも売られている現状、それらを一つづつ追う気力が出てきません。3ds MAX などの物理演算で何とかならないか、と、近道を探したくもなります。

ま、いずれも言い訳なのですが。


土木分野で粒子法がメジャーになるのは、幸か不幸かまだ先でしょう。良い手法ですが、変化の遅い分野ですので。
焦らずに、ゆっくり実装していきましょう。

3軸圧縮試験 供試体の直径

一軸や三軸の供試体の直径は35mmか50mmが多いと思います。

地盤工学会基準はφ35～100mmですが、φ35mmやφ50mmなら乱れの少ない試料を成形し、標準的な三軸圧縮試験機で対応できます。

先日、先の技術者と話をしていて知ったのですが、中型三軸圧縮試験機では、φ100mm以外にφ75mm（会社によってはφ65mmも）のペデスタルがあるそうです。基準に、以下の記載があるので、シンウォールやトリプルの内径に合わせたペデスタルが作られているとのことでした。

「側面の成形によって過度な乱れを与えるおそれがある場合は，サンプルチューブから取り出した試料の側面の成形を省略してもよい。」

教えてくださった方は「最近の流行」と言われていましたが、気を使ってくださったのでしょう。φ75mmは結構有名なようですので、昔からあったのでしょうね。そりゃそうです。トリプルで採取した試料を、わざわざリモールドでφ100mmにすることはないですから。

乱れの少ない試料として採取できるのが基本は砂質土まで。特殊なサンプラーを使わない限り、トリプルで礫質土の採取は困難です。したがって、最大粒径10mm程度までは対応可能なφ50mmの標準的な試験機で、通常は事が足りるといった状況になります。そのためか、会社には中型三軸圧縮試験機がないので、φ75mmのペデスタルは使ったことも見たこともありませんでした。

先日、試験を見せていただいたのですが、やはり大きいですね。供試体の高さが直径の2倍以上ですから、φ75mmで最低15cmになります。4供試体であれば、トリプルでの採取長は80cm程度以上必要でしょう。これは厳しいですね。静的3軸＋動的3軸なら、積算上は2本、実際は3本採ることになるでしょうね（積算にφ75mmの中型3軸はありませんが）。

調査計画を立てるためには、各試験の細かいところまで把握しておく必要がありますね。

単6電池

スタイラスペンの電池を買いに行こうと電気屋さんへ。

AAAAの単6電池なのですが、売っていません。
PC屋かと思い行ってみましたが、これも置いていません。

結局、ネットで購入。まだまだ浸透していないんですね。

3軸CUBと限界状態線

先日、砂質土の三軸CUBについて他の技術者と話をしていたのですが、ちょっと話が合わないところがありました。

p'-q軸上であれば、最悪、1供試体でもφ'を求めることができると考えていたのですが、その方は通常の有効応力経路で数点プロットしないとダメと言われる。

うーん？と思い調べてみると、私の間違いでした。
正確には、正規圧密粘土で c' = 0 とみなせるのであれば、1供試体でも可能です（当然、複数やるべきですが）。しかし、それ以外の土は適用できるかどうか不明。限界状態線の傾き q/p'=M の q/p' を、三軸圧縮試験のモール円のσ'1、σ'3とクーロンの破壊線の傾きφ'で表現し、M = 6sinφ'/(3-sinφ')　の関係を導いてますから。Mを使っていますので、砂は論外ということでしょう。

導出の過程（といっても中学レベルです）は地盤工学会「地盤技術者のためのFEMシリーズ2 弾塑性有限要素法がわかる」 p122です。非排水強度 cu と有効応力のφ'の関係も載っています。が、こちらはモール・クーロンとの違いを理解できていません。拡張ミーゼスを知りませんので。ちょっと調べてみましたが、「拡張」の意図する処がわかりませんでした。

しかし、惜しいですよね。感覚的には砂もCUBであればOKですけど。感覚ですから、ま、通用しません。SYSカムクレイもカムクレイをベースに発展させていますので、いずれは何とかなるかもしれませんが。

Civil3D 2014 とタッチ

出張用PCとして、VAIO Duo 11 SVD1122AJB + Win8 Pro を購入し、1か月がたちました。

http://phreeqc.blogspot.jp/2013/08/vaio-duo-11-svd1122ajb.html

これ、購入して正解でした。

一番良いのは、マウスを使用せずに済むこと。マウスは現場や車内で使えません。PCに搭載されているパッドはマウスに比べ操作性が劣ります（操作が遅くなるのでプチストレス）。そこでペンと指の出番。両手で直感的にタッチ操作できるので、案外手早く操作できます。これだけでも、タッチパネル搭載機種を購入した価値があったと思います。

欠点はスタイラスペンの右クリック。持ち直すたびにボタンの位置を探すか、ペン先で画面を長押しする必要があります。指での長押しは、どこを押しているのか見えないことがあり、細かい作業には向いていません。これらはプチストレス。

あと、現段階で例外的にタッチ操作の難しいのが Civil3D。

2014 からタッチに対応しており、拡大・縮小は指でも可能です。が、マウスほど思うように動きません。また、画面移動が全くできません。2本指でスワイプすれば移動できるはずですが、機能しません。結果的にマウスの使用が必須となりました。惜しい。

あと、「タッチ」パネルも as AutoCAD では表示されるのですが、Civil3D では表示されず。サポートに聞けば不具合とのこと。うーん。

しかし、今回使用してみて、CAD とタッチは相性が良いと感じました。いえ、惜しいと感じることが多かったので、改善に伴い処理速度は格段に向上すると思ったのです。オブジェクト選択後、中を1タップで右クリックメニューが出るようにするとか、ダブルタップにコマンドを追加できるようにしるとか。細かい部分をタップすれば、V-nas のように拡大窓が出るようにしても良いでしょう。まだまだこれからだと思います。

CAD とタッチ、今後が楽しみです。

LSFLOW と LS-RAPID

準三次元地すべり運動解析プログラムによる地すべり性崩壊の被害範囲の予測, 土木研究所資料, 3057, 1992

LSFLOW に使われている支配方程式が掲載されています。それを順に追っていました。

が、最初から？の状態。
いえ、式が難しいわけではなく、掲載されている式に誤りがあるようです。基本的な連続の式の導出から、中の符号、余計なdxがついていたり、必要なところについていなかったり、NがMになっていたり。あちこちに？な部分がありました。ま、こちらも初心者なので、これはこれで意味があるのかもしれませんが。

基本的には、連続の式とニュートンの運動方程式を使っているだけですね。粘性項と外力項があるので、ナビエ・ストークスに近いと思います（自重項はありません）。あとは外力として内部消費が考慮されていますね。
近似はFDM。今まで、なぜ土石流のような大変形が解けるのか不思議で仕方なかったのですが、基本的には変形を解いているのではなく、流体の移動に伴う高さの変化のみを扱っているようです。ですから、粒子が移動するような可視化はできません。納得です。
計算は準三次元とありますが、高さによって物性の変化はないため、2次元に近いですね。

ソースも公開されていますので、それを読んでいますと、ここにも矛盾が。いえ、input ファイルに書かれているパラメーターよりも、読み込むパラメーターのほうが多くなっています。また、地震波の入力や計算のことは記載されていないのですが、ソースは整備されています。作成した時期にズレがあるのでしょうね。ちなみに、ソースは中央開発さん作のようです。
input ファイルの作成は EXCEL ＋ VBA で、結果の可視化は Surfer などで十分対応できそうです。


よく似たソフトに五大開発さんのLS-RAPIDがあります。
作成は佐々先生のようです。論文になっている改良そりモデルでしょうか？
http://www.godai.co.jp/soft/product/products/LS-RAPID/index.htm

こちらも連続の式と運動方程式を使用していますが、粘性をあつかっておらず、自重、外力（せん断抵抗含む）＋αといったシンプルかつ個人的には理解しやすい構成になっています。一応、間隙水圧も考慮されていますが、Jakyの式は3次元応力状態が対象でなかったはずですし、摩擦係数等のフィッティングという点ではLSFLOWと同じ扱いでしょう。離散化も差分、変形でなく高さの変化みの考慮ですので、本質的には同じつくりですね。ま、LSFLOWはすべった後、LS-RAPIDはすべり発生から扱えるという点が、大きな違いかもしれません。
http://phreeqc.blogspot.jp/2011/07/blog-post_30.html

ポスト処理としての動画が掲載されていますが、波動みたいです。うーん。

両者とも、高さのみの変化なので、可視化に工夫が必要なのでしょう。難しいですが。
やはり粒子法でしょうか？こちらはさっぱり進んでいません。理解しないと。

adsorption isotherm と遅延係数

北海道環境保全技術協会が発行している技術レポートを読んでいました。

自然由来重金属や硝酸性窒素といった問題についてのレポートです。数か月前にまとめて購入していたのですが、パッと見、特に目新しいことは書かれていなかったので後回しになっていました。

「吸着工設計マニュアル」に移流分散と、総量の話が出ていました。洩れても基準値を越えなければOKという思想のようです。国もそうですから、ま、良いのでしょう。
http://phreeqc.blogspot.jp/2011/08/blog-post_26.html

現地発生土でしょうか、火山灰土を吸着層へ利用するレポートもありました。この中で adsorption isotherm に関し、Freundlich と Langmuir に対するコメントがありました。
個人的には、ただのフィッティングなんで、線形含めどれでも OK です。ただし、どの isotherm を仮定するかで誘導される式中の R の定義が異なります。
MT3DMS は Visual MODFLOW の移流分散パッケージに採用されており、非線形も含め3種類の isotherm を使用できます。が、Dtransu は線形しか扱えません。例の部長様のように、Dtransu で解くのに Freundlich でフィッティングする、というようななことをしなければ、個人的には何を使っても問題ないと考えています（基礎をおろそかにしてはいけません）。

所詮、遅延係数は計算上のフィッティング係数であり、物理的根拠は薄いと考えます。根拠を求めるなら、反応計算を実施すべきでしょう。
http://phreeqc.blogspot.jp/2010/10/blog-post_31.html

トリプルサンプラー

トリプルサンプラー（ロータリー式三重管サンプラー）をバラした写真です。

硬さが中位以上の粘性土、締まりの程度が中位以上の砂質土を対象としています。

スイベルヘッド。ロッドにつながる側です。ベアリングとバネが入っています。
8個の孔から水が入り、下に抜けます。下側にも同じ数の穴が開いています。
レンチのかかっている管がインナーチューブ。この中にアクリルパイプが入ります。外側にはアウターチューブを被せ、スイベルヘッドにつなげます。水はアウターチューブとインナーチューブの間を流れていきます。インナーチューブの上に1個あいている孔には出口側に逆止弁が入っているため、流入しません。

シューとメタル。先端だけはずしました。86タイプです。

メタルが回っても、ベアリングのおかげでシューは回転しません。

シューがメタルより飛び出していますが、押せばメタルに近づきます（インナーチューブを介し、スイベルヘッド内のバネを押すようになります）。掘削水は孔から出ますが、シューが出ているため試料に接触しません。

試料はシューに押し込まれ、アクリルパイプの中に入っていきます。

これを見れば、砂礫だとシューにあたって掘進できなさそう、、小さな礫であれば入るかも、細粒分がないと落ちそう、などと想像できるのではないでしょうか？

また、シンウォールはシューの部分だけ、デニソンはシンウォール＋メタル＆アウターチューブと考えると理解しやすいですね。孔径を考えると、イマイチ、デニソンの必要性を感じませんが。

先端を外したところ。アクリルパイプ内に試料が収まっています。

いつも思うのが、シューの中の試料を押し込んでよいのかどうか。試験時、先端は使わないので、簡単に押し込めそうなら押し込みます。オペさんに何も言わなければ、100%押し込んでます。これは通常のコア取り時の癖でしょうね。

アクリルパイプを取り出した後、採取長を測り、水抜き、シール、凍結して（細粒分が多いとそのままで）運搬します。

ただ、以前にも書きましたが、今は道具がよくなっていますので通常のφ66ダブルでも乱れの少ない試料は採取可能です。腕も必要ですが。

http://phreeqc.blogspot.jp/2011/06/blog-post_21.html

http://phreeqc.blogspot.jp/2012/10/blog-post_6021.html

その「腕」のナレッジ部分を基準化するのは困難でしょうね。非常に惜しいですが。

ま、お金のかかる取り方が標準かつお客様の御要望であれば、それはそれで仕方ないですが。

ライフジャケット

台船へ行き来するのに、ゴムボートを使用しています。

手漕ぎボートが破れかけたので、エンジン付きのしっかりしたボートを用意していたところ、船上でボーリング作業をしていたオペさんが迎えに来いと合図。

岸にいた一人が迎えに行こうとし、ライフジャケットを付けて手漕ぎボートに乗り込んだ瞬間・・・ドボン。後転しながらひっくり返ってしまいました。一同、騒然。で、大爆笑。

水際の水深50cmだったので爆笑で済んだわけですが、あれが水深数mだと慌てたでしょうねえ。

まず、「携帯大丈夫か！」
次に、「ライフジャケットって、ちゃんと膨らむんだね。」と妙に感心。

今まで、自動膨張のジャケットを着てきましたが、膨らんだところは見たことがありませんでした。落ちたのも初めて見ましたし。あの瞬間を見ると、手で引っ張って膨らませるっていうのは難しいかもしれません。安くて暑いジャケットは浮かないそうですし。

うーん。恐るべしライフジャケット。
重要性を再認識しました。

Dtransu の流速 その2

飽和時の計算について、マニュアルには何も書かれていませんでした。


仕方ないので、ソースを見てみました。

ここですね。


C     DETERMINE RELATIVE HYDRAULIC CONDUCTIVITY
C
        LTB = MTBL3D( MAT )
        DO 61 K=1,NN
          I = LM(K)
C
c-------- confined problem ----

          if( kan2.eq.0 ) then

            CR(I) = 1.0
            TXX(K) = 100.0D0
            go to 61
          end if
c
          IF( P(I).LT.0.0 ) THEN
          CALL INTERP(XYP( 1,NMTB(1,2,LTB)+1 ),TX,P(I),NMTB(2,2,LTB),1)
          TXX(K) = TX
          IF( CR(I).LE.0.0 )
     2    CALL INTERP(XYK( 1,NMTB(1,1,LTB)+1 ),TX,CR(I),NMTB(2,1,LTB),0)
          ELSE
          CR(I) = 1.0
          TXX(K) = 100.0D0
          END IF
C
   61   CONTINUE



      SUBROUTINE SETCR( KODE,P,P1,CR,T,CS,POR,XYK,NK,XYP,NP,XYC,NC )
C
      IMPLICIT DOUBLE PRECISION (A-H,O-Z)
C
      COMMON / CTPRM2/ KAN1,KAN2
C
      DIMENSION XYK(2,NK),XYP(2,NP),XYC(2,NC)
C
      IF( CR.GE.0.0D0 ) RETURN
      IF( KAN2.EQ.0 ) GO TO 2000
      IF( P.LT.0.0D0 )  GO TO 1000
C
C---- SATURATED ---------------------------
C
 2000 CR = 1.0D0
      T  = POR
      CS = 0.0D0
C
      RETURN
C
C---- UNSATURATED -------------------------
C
 1000 IF( KODE.LT.1 ) CM = P
      IF( KODE.GE.1 ) CM = (P+P1)*0.5D0
      CALL INTERP( XYP,T,CM,NP,1 )
      CALL INTERP( XYK,T,CR,NK,0 )
      CALL INTERP( XYC,T,CS,NC,0 )
      T = T*POR*0.01D0
      CS= CS*POR
C
      RETURN
C
C
      END




コメントがないので誤っているかもしれませんが、飽和のフラグを立てると、比透水係数1、体積含水率 = 有効間隙率で計算しているのでしょう。単純に、この条件で負圧となれば不飽和時よりも上昇流や水平方向の流れが発生しやすいということでしょうね。


結局、変換ツールのソースを、圧力水頭がマイナス時は、Vi ⇔ { 0, 0, -0.1} として吐き出すように修正しました。不飽和の部分では鉛直下降流が起こっている、という、よくある簡略化です。

そうすると、TECPLOT の streamline が、思い描いたように綺麗になりました。


うーん、いいのか？

Dtransu の流速

Dtransu の計算結果を Tecplot で可視化しています。

といっても、移流分散ではなく、浸透流。
流速3成分から streamline を描こうとしていたのですが、なかなかできません。

Tecplot の設定は問題なさそうでしたので、変換結果を見ました。すると、節点流速が出ていません。そりゃ、描けないはずです。

Dtransu に付属の変換ツールの詳細を見てみますと、流速を読み飛ばす仕様になっていました。それなら、「読んだ後に節点情報とまとめて書き出せばよい」と考え、ソースを修正。
修正したソースをコンパイルし、実行すると、ちゃんと書き出せました。で、それを Tecplot で読み込み。

が、エラー。
なぜ？と思いながら調べてみますと、原因は私のミスでした。
変換ツールが読み飛ばしていた流速は、要素流速。これを節点流速だと思い込んで、節点情報に付加して書き出していました。Dtransu の浸透流では、要素流速しか書き出さない仕様で、節点流速に直してくれないんですね。移流分散は実流速を節点流速として書き出してくれるのですが。

Tecplot では要素流速の取り込み方がわからなかったのと、それを別途書き出して読み込ませるのも手間なので、結局、ダミーの移流分散として回すことにしました。


計算時間は掛かりましたが、節点流速は書き出されました。
濃度を設定していないため、非数 NaN を 0 に変換し読み込むと、ようやく streamline を描くことができました。長かった。


しかし、これで終わりではありませんでした。
なぜか不飽和の部分で上昇流が発生。今回は飽和の計算しかしていないので、不飽和部分が圧力水頭を持っているとは考えていませんでした。

この続きは、また後日。

斜面災害

8月末から9月頭にかけて、2度の大雨でした。

「大体、2回目の雨で崩れるよね」なんて言ってたら、その日の午後に呼び出し。

国道が崩壊した現場と、市道がすべった現場の2種類でした。出張の合間に、内業を進めようと帰社したら、逆に仕事が増えてしまいました。

国道は完全に地質屋のミス。でも、こういうこともあるのだと、改めて勉強になりました。

市道は早急に観測体制を確立する必要がありました。伸縮計を設置したのですが、通信機や警報機、回転灯などの配線接続は電気工事屋さんみたいな仕事です。もともと、伸縮計の設置は一人ではできなかったのですが、今ではより多くの知識と人数が必要になりました。

当然、その日のうちにネット上で監視できる体制にするため、通信機やサーバーの設定も必要です。今では携帯、スマホで経時グラフを確認できるようになっています。

開発側も大変でしょうが、使う側もそれなりの努力が必要です。

デニソンサンプラーと削孔径

二重管サンプラーの代表として、デニソンがあります。

デニソンサンプラーの中にはシンウォール管が入ります。そのため、掘削径はシンウォールより一回り大きいφ116mmとなります。これは一般的によく知られていますし、正解です。

しかし、二重管サンプラー = φ116mmは誤り。二重管サンプラーもトリプル同様にいくつかの径が出ています。掘削径 86mm で採取できるタイプも販売されています。地盤工学会基準「地盤調査の方法と解説」にもいくつかの径が紹介されています。

しかし、デニソン以外のサンプラーは、あまり見ません。トリプルの86タイプで硬質粘性土も採取できるからでしょうか？

ちなみに、固い粘性土であれば、小径倍圧サンプラーという手もあります。これなら、φ66mmでも採取可能です。ただし、水圧をかけるため、通常のホースでは破裂しますし、場合によっては機械が浮かないようにアンカーを取る必要があります。水周りを耐圧仕様にする必要があるため、やや面倒です。個人的には避けたいサンプラーです。

粘性土の場合、個人的にはアクリル管よりシンウォール管のほうが取り出しやすいので好きです。デニソンで上手に取っていただけるほうがありがたいですね。

トリプルサンプラーと削孔径

ボーリング調査時の乱れの少ない試料採取では、一般的に、以下の削孔径が必要とされています。

シンウォールサンプリング：φ86mm

デニソンサンプリング：φ116mm

トリプルサンプリング：φ116mm

ところがこれ、トリプルのφ116mm は誤りです（積算上はφ116mmしかありませんので正解ですが）。

市販品の中にはφ86mm で採取可能なトリプルサンプラーがあります。昔からよく使われており、普及しています。

私が頻繁に目にするのがコアーパック社の86TRIPLE。

http://www.corepack.co.jp/seihin/packtube_4047/tripletype_4051.html

砂を取るのにφ116mm は必要ありません。φ86mm で十分です。径が小さいと、試料落下のリスクも小さくなります。

オペさんにも何度か話を伺いましたが、サンプラー径が 116mm から 86mm になっても、礫径の影響はほとんど受けないそうです。礫にあたると、どちらも入らない印象が強いようです。それよりも、落下のリスクが小さい、削孔径の小さい方が容易等のメリットを選択される方が多いようです。

メーカーにも話を伺いましたが、一概に「礫が大きければ大きなサンプラーを」、というわけでもなく、「現場に応じて」とのことでした。要は、腕と経験で判断、ということなのでしょう。腕は重要だと思います。

http://phreeqc.blogspot.jp/2011/06/blog-post_21.html

http://phreeqc.blogspot.jp/2012/10/blog-post_6021.html

分業が進み、現場を眺めているだけでは、径の違いにすら気づかないのでしょうね。いまだにトリプル=φ116mm が常識？となっているのはそのせいでしょう。