有限体積法

有限要素法や差分法はよく聞くのですが、有限体積法は今まで触ったことがありませんでした。

応力など、なめらかな形状が必要な場合にはFEM、地下水などはボクセルの差分で十分だと思っています。FVM はどういった時に使うのでしょうか？周りのソフトは FEM か FDM ですので、昔流行った手法かと思いきや、そうでもないようです。ネットでは、流体解析では FVM というのも見かけます。よくわかりません。

そもそも、有限体積法という名前しか知りませんでした。
調べてみると、以外に容易かもしれません。読んだのは平瀬創也「C#で学ぶ偏微分方程式の数値解法」。
質量保存を考える際、微小体積を仮定しその中の収支を追いますが、FVM もまさにそれでした。格子点周りに微小体積（コントロールボリューム）を想定しすれば、積分形の式を利用しやすくなります。さらにそれを差分形式で近似すればおしまいのようです。
でも、わざわざコントロールボリュームを仮定するのが、まだピンと来ていません。結局差分形式で表現するなら、最初から差分法でも良いような気がします。

非構造格子でも対応できるそうですが、そこまではまだ理解できていません。
そういえば、MODFLOW-USG も驚くような格子になっていましたね。差分法のままなのでしょうか？どうなっているのでしょう。

FVM の解説本は少ないようです。ま、構造格子のみだと上記の本で問題ないので、とりあえずは良しとしましょう。

在宅勤務

今年のはじめだったでしょうか、米Yahoo!の CEO が在宅勤務禁止の意志を示しました。

その後、禁止になったのかどうか知りません。その真意も分かりません。ただ、会社で人と人が顔を突き合わせた方が、ディスカッションが増える、より良いアイデアが生まれるというような意図は理解できます。時代に逆行しているとの論争も出ましたが、個人的には行き過ぎた在宅勤務を引き戻し、モニター越しでなく、面を合わせて話し合うというのは基本的かつ重要なことだとおもます。逆行しているとは思えませんでした。

私の勤めている会社は、在宅勤務を認めていません。

家では仕事をしていないことになっています。

Yahoo!より、進んでいるのかもしれません。

DEM と DTM

国土地理院さんが基盤地図情報を更新されています。

先日のUPで、さらに LP ベースの 5DEM が追加されています。以前、写真測量ベースのDEMしかなかった箇所も整備されていました。いいですね。

そういえば、「DEM と DTM の違いは？」と聞かれたことがります。個人的には、グリッドデータになっているものが DEM という解釈をしていましたが、正確には答えられませんでした。同じ解釈をしている方もいらっしゃいますが、そうでもない解説もあります。
http://gis.stackexchange.com/questions/5701/what-is-the-difference-between-dem-dsm-and-dtm

DTM: bare-earth representation with irregular spaces between points (non-raster).
DEM: gridded raster representation of the DTM.

国土地理院さんの解説では「同義」と解釈してもよさそうです。
http://www1.gsi.go.jp/geowww/Laser_HP/faq07.html

DEMに類似する用語にDTM(Digital Terrain Model)という略語があり、日本語では数値地形モデルといいます。これら2つの言葉は異なる語源からなりますが、同義語として扱われる場合が多く、一般に「DEM」が使われる場合が多いようです。

国土地理院さんの動きを見る限り、取り急ぎ、全国整備が目標なのかと思われます。
今後は時系列で整備していただくと、より便利になりますね。どの地域に何年のLPがあるのか分かれば、崩壊前後のデータがあるのか、すぐに判断できます。
今後に期待です。

河川の流れの基礎方程式

LSFLOW は連続の式（質量保存則）と運動方程式（ナビエ・ストークス）を利用しています。

その詳細は非公開と聞いていたのですが、土研資料にはソースまで載っていますし、砂防学会「地震砂防」という本に、支配方程式の展開や離散化の方法も載っていました。砂防分野では有名なのでしょう。

式の展開を初めて見たときは、独特な定式化だと感じていました。特に、クーロン則などを適用している外力の箇所。解釈は結構、自由なんですね。ま、再現できれば良いわけです。

また、地すべりの運動では鉛直方向の流速が水平方向に比べて無視できるという仮定のもと、式の簡略化が行われています。これを見たときは独特だなあと感じつつ、昔のことですから、計算容量を小さくしたかったのだろうなあ、程度に考えていました。

ところが先日、河村哲也「河川のシミュレーション！」という本を読んでいて、同じような定式化を見ました。河川の流れの基礎方程式というそうです。「鉛直方向は無視」という簡略化は河川砂防分野では常識なんでしょう（私が知らなかっただけのようです）。こういう考え方が河川砂防分野のベースにあったからこそ、LSFLOW も、あのような定式化をしていたのでしょうね。もっといろいろ学ばないといけません。

LSFLOW も学位論文等で FVM に改良されています。LSFLOW3、LSV と呼ばれているコードがそれに該当するのようです。
いずれにしても、inputファイルの作成～計算～可視化（MVS）までの実装は完了しましたので、あとは崩壊前の地形の調整を含めたパラスタのみです。地形の調整は時間がかかるでしょうね。
http://phreeqc.blogspot.jp/2013/10/blog-post_9060.html
もう少し考えないといけません。

Visual Fortran と Win8.1

計算用のPC に Visual Fortran Composer XE 2013 SP1 UP1 を入れています。

先週末、そのPCにて 6コアフルの計算させるていたため、他の作業ができませんでした。良い機会なので、VAIO　＋ Win8.1 に VF を入れて作業しようと思い、インストールに取り掛かることに。
最近知ったのですが、VF のシングル・ユーザー・ライセンスは PC 1台につき１つではなく、1人につき1つなんですね。つまり、個人使用の PC が3台あれば、3台すべてにインストールできるわけです（当然、同時使用はナシ）。

VF をインストールしていると、13/30まで進んだ頃にエラーが発生。何度かやり直してみたのですがエラーで進みません。中途半端に入ってしまったので、uninstall しようとしたら、これもエラー。再起動をかけようとすると、ずーっと再起動中の表示で落ちない。電源ボタン長押しで強制終了し、再度起動する羽目に。
その後、uninstall してから、再度インストールするも、全くダメ。失敗です。
結局、システムの復元でさっぱりしました。

13/30 はちょうどVS 2010 shell をインストールする過程でした。せっかくですので、shell を最新版にすることに。
VS 2013 はexpress 版が3種類ありました。for Windows が2種類。DLページに説明は書かれていませんが、たぶん、ストアアプリ用と、デスクトップ用でしょう。気にせず、for Windows Desktop を DLし、インストール。firewallを落とすのを忘れていたのですが、一発でOKでした。

しかし、VFが2013に対応していないのか、統合できません。
HPで確認すると、Win8、VS 2012までのようです。
http://www.xlsoft.com/jp/products/intel/compilers/fcw/

• オフロードと SIMD 拡張を含む多数の OpenMP* 4.0 の機能をサポート SP1 より追加
• Fortran 2003 規格のサポートを拡張 (ユーザー定義の派生型 I/O) SP1 より追加
• Microsoft* Windows* 8 および Microsoft* Visual Studio* 2012 をサポート

仕方ないのでアンインストールしてから2012を入れることに。2012も express 版を入れてみました。
が、これも統合できません。この時点で朝の4時。続きは翌日に持ち越し。


翌日、気を取り直してPCのメンテから始め、再度 VS 2010 shell 付属の英語版を入れてみました。これはあっさり入りました。こうなると、日本語版も入れてみたくなります。ダメもとで再度チャレンジしてみると、13/30でかなり考えながらも、今度は最後まで到達！何が悪かったのか分かりませんが、入りました。

後々、サポートに聞いてみると、「リリースノートを見なさい」とのお返事。
見てみると、ちゃんと書いてありました。概要は以下の通り。

・Windows XP SP3 ～ 8.1までOK

・IA-3、インテル64 対応アプリケーションのビルドに、Microsoft* Visual Studio* 開発環境あるいはコマンドライン・ツールを使用する場合は、次のいずれか:
o Microsoft* Visual Studio* 2012 2010 Professional Edition 以上、2008 Standard Edition 以降、 2010 Shell、2008 Shell

・IA-32 アーキテクチャー・アプリケーションのビルドに、コマンドライン・ツールの
みを使用する場合は、次のいずれか:
o Microsoft* Visual Studio* Express 2012 for Windows Desktop
o Microsoft* Visual C++* 2010 Express Edition [3]
o Microsoft* Visual C++* 2008 Express Edition 

http://www.xlsoft.com/jp/products/intel/compilers/fcw/2013/Release_Notes_sp1_up1.pdf

2012 PROはもっていませんので、2010 shell で BEST だったようです。そこに行きつくまで、1日かかりました。つくづく、周りにプロがいれば、と感じた週末でした。

Tecplot で水面の座標抽出・差分表示（備忘録）

備忘録です。


① Iso-surface で P = 0m の面（水位）を表示させる。

②“Data” - “Extract” - . “Iso-surface” で水位面を抜き出す（Zoneとして新たに追加される）。すべての zone から一気に抜くことも可。

③“Data” - “Alter” - . “Specify equations”にて、変数を2個追加。名前は何でもOK。ここでは、WL1 （ H を代入）、WL2 （ 0 を代入）を作成。
｛WL1｝={H}
｛WL2｝=0
＊この後、“File” - “Write Data File”で、作成した水位の座標出力可能。
＊ただし、吐き出されるXY座標はZone毎に微妙に異なるので、単純に差分は取れない。

④水位を抽出したZoneでは、それぞれジオメトリが異なっている（Z=H（水位）となっている、XYが微妙に異なる）ため、単純に差分は取れない。そこで、Zone同士を比較する場合、比較先のジオメトリにて、比較元のHを補間してから差分をとる手順になる。まずは補完。
“Data” - “Interpolate” - . “Kriging”にて、Destination zone のWL1に、Source zone のWL1をkrigingで補間・代入。ここでは、Zone146のWL1(=H：水位)をZone74のジオメトリ上のWL1に補間・代入している。

⑤差分計算

Destination zone で指定したzone（ここでは74）を選択し、{WL2}={WL1}-{H} でCompute。

WL1 には Zone146 の値が④の過程で補間・代入されているので、選択した Zone74 の水位との差分が、Zone74 の WL2 に代入される。

⑥Time Stands
差分をとり終わったら、表示したい Zone を Time に追加。

⑦コンター図表示で完成



=========================================================

2013.11.23 追記

kriging 時、Driftはナシです。

Tecplot で水位の差分表示

Tecplot で水面を表示できますので、その差分をソフト内で計算・表示できれば、座標書き出しの必要性が小さくなります。

ま、できるでしょうと思いチャレンジしてみたのですが、できませんでした。水面のみ抽出した2つのZoneを指定し、差分を取るだけなので容易にできそうなのですが、ダメ。いえ、差分を取る計算式を入力する箇所はあるのですが、Zoneを指定すると"incompatible zones "というエラーが出るのです。

ググってみると、フォーラムに以下のような質問がありました。

http://www.tecplottalk.com/viewtopic.php?t=339

I have got a problem with the data -> alter -> specify equations function! I have tried for ages to visualize with this tool an anspect of my flow field, but it doesnt work.
description of my problem:
1. I load two files into tec360 with the same geometry datas but other flow fields. This works, and so i can see 10 zones (5 of each loaded files)
2. What i wanna do now is to subtract the pressure level of the first file of the second one. (zones 1-5 are from the first one, 6-10 from the 2nd one)
3. So i go into the specifiy functions tool and write:
{pdelta} = {p0}[6] - {p0}[1],{p0}[7] - {p0}[2],... (and so on for all the corresponding zones from 1 to 10)
-> but now i cannot compute this value pdelta and visualize it...-> if i try to compute it, there are three error messages which appear:
1. if i don't select any zone in the 'zone to alter' list, the error appears: 'invalid zone set requested'
2. if i select all zones, the following error appears: 'incompatible zones referenced in equation'
3. and when i select only a part of the zones, e.g. for the equation{pdelta} = {p0}[6] - {p0}[1] the zone 1 & 6, the following error appears: 'must alter all zones when creating a new variable'
So my question is, how can i compute and visualize the problem i described? what is wrong with my equation? Is it the syntax or has it to do with my zone selection or status? 

驚きです。
私も、まったく同じことをしようとして、同じところで躓いていました。 欲しい機能もわからないところも同じということです。この質問には解決策が載っていませんでしたので、Tecplot にその機能はないということなのでしょう。マクロを組めば、できるかもしれませんが。

代わりと言っては何ですが、1ステップ毎の計算は可能でした。欲張らないことですね。

http://www.tecplottalk.com/viewtopic.php?t=1124

In order for Tecplot to perform a calculation involving more than one zone (like the time average) all the zones involved have to be "compatible". That is, they have to have the same underlying order, and the same number of data points.

The "incompatible zones referenced in equation" message means that at least one zone had a different oder, or different number or points, or both than the rest of the zones.

For these instructions, think of two incompatible zones, each with five variables (X,Y,Z,T,P), where T is Temperature and P Pressure. And suppose that you want to subract P of zone 1 from P of zone 2:

1) Go to Data/Alter/Specify Equations and create two new variables, by highlight both zones in the Zone to Alter box and entering two equations like these: 

{P1} = 0
{P2} = 0
Then click Compute.
Both of your zones will now have two extra variables. 

2) Still in Data/Alter/Specify Equations, highlight only the first zone in the Zones to Alter box and then erase the two equations you entered moments ago and enter this one:
{P1} = {P}
and click Compute. This will copy P of Zone 1 into P1 of zone 1.
Now you can close out of Data/Alter/Specify Equations.  

3) Go to Data/Interpolate and choose a method. Select Zone 1 as the Source zone and Zone 2 as the Destination zone and P1 as the variable and click Compute. This will interpolate P1 of Zone 1 onto P1 of Zone 2.
P of Zone 1 has now been transferred to 

4) Now close out of Data/Interpolate and go back to Data/Alter/Specify Equations. Erase an equations that are present, and select Zone 2 in the Zones to Alter box.
Then enter an equation like this:
{P2} = {P1} - {P}
and click Compute. P2 of Zone 2 will now hold the result of the calculation of P of zone 1 - P of zone 2

なるほどです。確かに、抽出した水位の座標はそれぞれ異なりますので（水面（Z）だけでなく、なぜかXYも微妙に違う箇所がある）、エラーを吐くわけです。
本来は、上側の質問のように複数のタイムステップ同士を比較したかったのですが、贅沢は言えません。後者で実施してみると、うまくいきました。

続く。