透水係数と回転マトリックス

3次元の浸透流解析でパラスタをしていました。

なかなか合わない個所があり、kx, ky, kz を変更してみることにしました。それが正しいかどうかは別として。

砂岩・泥岩の互層であり、山を歩いた感じでは一様に傾斜しています。その平均的な走向傾斜は把握しています。
また、Dtransuであれば方向余弦の入力で回転した透水係数(軸)も計算に取り込めます。マニュアルに記載されていますし、2年前に確認しています（なぜ符号をそのようにしたのかは全く覚えていません。GTRANの座標系だったかな？）。
http://phreeqc.blogspot.jp/2011/07/dtransu.html

で、いざ反映しようと思うと、あれ、分かりません。方向余弦を出すにも、2軸で回転させてしまうと、わからない角度が出てくるのです。回転マトリックスを使うにしても、平面上を1回転までは計算できるのですが、次に傾斜方向へもう一度回転となると、その方法が分かりません。どちらかの軸だけで1回回転ならOKですが。

調べてみると、走向傾斜はオイラー角で表現するのが便利だとわかりました。これ、知りませんでした。回転後の軸を使って、もう一度回転させる、そう考えれば回転マトリックスを2回使うだけ、積をとるだけでOKでした（数学、もっとやっておくべきですね）。個人的には、こちらのサイトが分かりやすかったです。
http://www6.ocn.ne.jp/~simuphys/daen1-1.html
wikiのGIFがイメージしやすいです。
http://ja.wikipedia.org/wiki/%E3%82%AA%E3%82%A4%E3%83%A9%E3%83%BC%E8%A7%92

上のサイトでは3回転させてますが、走向傾斜であれば2回転でOKなので、もっと簡単な式になります。

この程度ならEXCELで配列数式があるかと思いましたが、ないですね。ま、それほどややこしいものでもないので、sin、cos と手入力で仕上げました。走向傾斜を入れたら、方向余弦を算出するようにしました。もっと簡単な方法もあると思いますが、ま、今回はこれでOKです。


で、目的だった走向傾斜と透水係数の異方性を考慮した計算をしてみました。

結果、今度はほかの箇所が合わなくなりました。うーん。



*********************************************************
2013.11.02
以下にEXCELデータをUPしました。
https://sites.google.com/site/geochemist001/resources/RMatrix