Blue Kenue

TELEMAC に使える pre-post GUI を2種試しました。

Blue Kenue

https://nrc.canada.ca/en/research-development/products-services/software-applications/blue-kenuetm-software-tool-hydraulic-modellers

カナダ産 Win 用 GUI です。SHP（ライン） を読み込めるので、CAD データを変換すればセンターや盛土位置などを取り込めます。それを利用しオートメッシュ可。やはりこれが便利（ダイアログ表示にはバグがありました）。

境界条件設定位置を指定できるのですが、そのパラメータ等の設定はできません。

SALOME-HYDRO

https://www.salome-platform.org/contributions/edf_products

Linux の GUI です(「SALOME」 なのでフランスでしょうか？)。

メッシュだけでなく設定ファイルも作成できるようです。が、私の環境では SHP 読み込み段階でエラーがかかり先に進めませんでした。

で、Blue Kenue を利用。

設定ファイル（cas）は example から似たケースをコピーし、テキストエディタを使って修正しました。境界条件の時系列ファイルも同様です。これが速いでしょう。

コマンドプロンプトでバッチファイル替わりに py ファイル、そのあとに設定ファイルを指定し、TELEMAC を実行。で、計算が始まります。

流れのみですが、7万メッシュ、0.1秒間隔、24時間の計算が14.5時間かかりました。parallel だともっと早くなるのでしょうか？

で、再び Blue Kenue 。

結果を3次元で可視化します。が、おかしな結果。地形が正しく読めておらず、すべて0mとなっています。

どうも、メッシュに標高を与える際に「Bottom」の名称ではダメで、「BOTTOM」でないと認識してくれないようです。そういえば、設定ファイル内のキーワードはすべて大文字になっていましたね。ルールなのでしょう。

修正し、計算した結果は良好。

これでデータ作成から不定流、可視化まで一通り試すことができました。

残るは河床変動計算です。

TELEMAC-MASCARET インストール

TELEMAC-MASCARET

http://www.opentelemac.org/

自由表面流れのソルバーです。フランス産。

マニュアルを見ると、2次元の河床変動計算が可能でした。

支配方程式を河川砂防技術基準と比較すると、連続の式にソースタームを加えて発展させた形です。運動方程式の外力項（風なども含めたソースタームと書かれています）にマニングの式を使えます。

河床変形の式も河川砂防技術基準と一致。

掃流砂量式には日本人の名前が載っていませんでしたが、浮遊砂と共に扱えることを確認できました。

OKです。

ソルバーですので CUI です。Fortran ソースが配布されています。

pre-post の GUI は他サイトのフリー版を活用する方針です。

インストールには SVN クライアントや Fortran コンパイラーが必要です（Delft もそうでしたが、海外の河川・海岸屋さんはこれがデフォなのでしょうか？）。 

さらに Python3 ではなく 2.7 を利用します。Win の環境を汚したくないのですが、conda の仮想環境を利用できそうにありません。Ubuntu の Docker も配布されていません。VS ＋ Intel VF で 個別にコンパイルして構成しようかとも考えましたが、それも難しいというか、説明自体がありません。

結局、新たな PC（Win10） に automatic installer でインストールすることにしました。

順序は Python 2.7 と関連ライブラリのインストール、TortoiseSVN のインストール、gfortranなどコンパイラーのインストール、SVN サーバーに接続してソースをダウンロード、コンパイル、コマンドプロンプトの設定までが一連となっています。この一連作業を automatic installer が実施してくれるのですが、Forum を見ると失敗が続出しています。複雑ですからね。

失敗後のフォローがサイトに掲載されています。私の環境ではそれ以前のソースの DL 段階で躓いたので、そこからやり直しました。

http://wiki.opentelemac.org/doku.php?id=instructions_to_finish_installing_the_telemac_system_manually_after_a_failure_of_the_automatic_installation

インストール作業は躓き後のフォロー含め2～３時間程度。

automatic installer では v7p0 まででしたので、手動で最新の v8p1r1 も入れました。後者 の parallel の起動方法がわからなかったのですが、single ではどちらも examples の動作を確認。環境構築完了です。

ひとまず先に進むことに。

統計モデリング

先日、以下の図書をベースにした統計モデリングにかかわる話を聞きました。

樋口「予測に生かす統計モデリングの基本」

以前読んだ「地下水モデル」にデータ同化が記載されており、どのようなものかを知りたいと購入していた図書でした。が、意味が分からず寝かせていました。

その後、逆解析に粒子フィルタを利用された方の話を聞いたこともあり、この機会に引っ張り出して読み直しました。

ま、結局、最後まで読んでも具体的にどのように役立つのか、どのようにシミュレーションに組み込めば良いかは、相変わらず理解できませんでした。

基本コードが手元にあるので、いずれ迫られたら読んでみましょう。

以下、備忘録です。6章まで。

*********************************

加法定理の一般化、同時確率から確率変数を1つ取り除く作業が周辺化

→周辺確率を求める作業

P（A1）が周辺確率。

P（A|B)が「Bが所与のもとでのAの」条件付確率

周辺化

P（A)=∑B P(A,B)

書き下すと（加法定理）

P(A1)=P（A1,B1)+P（A1,B0)

P(A0)=P(A0,B1)+P（A0,B0)

全体：100、A1：30、B1：60、A1∩B1：10の時、

P(A1)=10/100+20/100=30/100

乗法定理

P(A,B)=P(A|B)p(B)=P(B|A)p(A)

書き下すと

P（A1,B1)=P(A1|B1)p(B1)=P（B1,A1)P(A1)

10/100=10/60*60/100=10/30*30/100

ベイズの定理（乗法定理を周辺確率で割る）

P(A|B)=p(B|A)p(A)/P(B)

=p(B|A)p(A)/∑A P(B,A)∵周辺化

=p(B|A)p(A)/∑A(P(B|A)p(A))∵乗法定理

*********************************

状態空間モデル：システムモデルと観測モデルの連立

線形：ガウス状態空間モデル

システムモデル（計算更新時のシステム誤差 vt を含む）

xt=Ftxt-1+Gtvt　 vt～N(0,Qt)

観測モデル（観測誤差wtを含む）

yt=Htxt+wt　 wt～N(0,Rt)

非線形：非線形・非ガウス状態空間モデル（データ同化で利用）

xt=ft(xt-1,vt)　 vt～p(v|θsys)

yt=ht(xt,wt) 　wt～p(w|θobs)

θは未知→最尤法で推定

マルコフ性1：xt～p(xt|xt-1)

マルコフ性2：yt～p(yt|xt)

*********************************

予測分布:p(xt|y1:t-1)

フィルタ分布:p(xt|y1:t) 

平滑化分布:p(xt|y1:T)

情報処理での「フィルタ」は信号処理のフィルタとは別

情報処理での「平滑化」はデータ解析の平滑化とは別。回顧による知識発見が目的。工学系の応用問題は、リアルタイムに適切な処理を施すことが重要で、平滑化分布から得る意味があまりない。

固定ラグ平滑化

固定点平滑化：データ取得毎に初期分布の改善（データ同化で利用）

※固定点平滑化は地下水シミュよりも、力学シミュで役立ちそうです。が、具体的な処理方法が思い浮かぶまで理解できていません。うーん。

Delft3D GRID

RGFGRID でメッシング。

ひとまずチュートリアルを流してみましたが、イマイチ。切れないわけではないのですが、面倒の一言。いえ、フリー版にしてはよくできていると思うのですが、オートメッシュに慣れてしまうとそれ以外を使えなくなります。

計算機能自体は多いので魅力的なのですが、河床変動の式が見当たらず、pre にも不満。猫に小判＋暖衣飽食に慣れすぎ状態です。

結論：残念ですが Delft3D は寝かせ、次のソフトに進みます。

*********************************

20200628追記

TELEMAC のメッシュ（非構造格子）を読み込めましたが、Delft3D 用に書き出すことはできませんでした。構造格子でないとダメなのかな？

Delft3D GUI

GUI は別途申し込みが必要です。

数日でライセンスファイルと DL リンクが送付されてきました。こちらはインストーラーがついており、サクッと入れることができました。

先にコンパイルしていた exe や dll をどこに入れるのかを探せなかったので、エラーが出るまで進めることにしました。

ざっくり、マニュアル3～5章がtutorial、6章が計算、7章が可視化です。

5章の通り flow のチュートを進め、6章で計算。ここで、エラ-が出ました。x64 フォルダ以下にexeが見当たらないとのことです。

場所がわかったので、GUI インストール先の x64 フォルダに先日コンパイルした x64 フォルダの中身をコピー。

そうすると、計算が完了しました。

可視化もOK。

ひとまず環境構築が終了です。

今回のチュートではメッシュの切り方が含まれていませんでした。これは必須なので、あとで見つけねば。

Delft3D インストール

Delft3D のインストールに、思いのほか時間がかかりました。

Win での流れを補足しておきます（プレコンパイル版の配布はありません）。

https://oss.deltares.nl/web/delft3d/source-code

1. Subversion （ソースコード管理サーバー）から 最新リビジョンを確認して DL （Subversion client を使用）

2. src フォルダ内の prepare_sln.py （Python script） で VS の sln ファイルを作成 （VS、.NET、Intel Parallel Studio XE の Version指定）

3. slnファイルを VS で開いてコンパイル

Subversion を知らなかったのですが、Git よりも古いそうです。

クライアントはサイトに記載のあった TortoiseSVN を使用。適当な場所にフォルダを作成し、右クリックから SVN Checkout。私の場合、必要なユーザー名は Delft3D に登録したメールアドレスの＠より前の部分として登録されていました。

DL、コンパイルに必要なソフトはサイトの Prerequisites にまとめられています。Win では VS ＋ IPS XE の指定なので、無償ソフトのみでは対応できません（本体はオープンソースなのに！）。VS でも必用とされるモジュールを指定し、インストール or 追加する必要があります。

examples が動けば成功だそうです。

今回、trunk 版はダメ、リビジョン 65991 の single は動きましたが、parallel はダメでした。 

次に進みましょう。

Delft3D

HEC-RAS が2次元河床変動計算に対応していないため、Delft3D を試そうかとマニュアルを確認。

https://oss.deltares.nl/web/delft3d

日本では、港湾関係で利用されている例がありました。どうも、波の影響を含めた浜？の形成や河口付近の河床変動計算がこのソフトの主目的のようです。なお、RIC の Ver.2ではソルバーとして含まれていたようですが、現行の Ver.3 では見当たりませんでした。

マニュアルVer.3.15 は690ページ。HEC-RAS より薄いのですが、欲しい記述を探すのが大変です。

まずは流体。NS 方程式がベースです（運動量方程式9.7）。連続の式（式9.3）もそうなのですが、球の座標系？のため複雑に見えます。

掃流砂・浮遊砂の区分も取り扱っています。掃流砂量式には芦田-道上が含まれていました（表11.1、式11.183）。

浮遊砂は移流分散方程式で解いています（式11.1）。濃度扱いですね。

あとは河床形状の計算式だけなのですが、これを探せませんでした。

どこに書いてあるのでしょう？これを扱うプロが周りにいないので聞けません。

コチラの文献には見慣れた形の式がありました（式2.35）。これが採用されていることを期待したのですが。3章でも説明されていないということは私のなにか大きな勘違いなのでしょうか？

https://repository.tudelft.nl/islandora/object/uuid:e6c74ee1-a795-419f-92b2-d023c2c6a050

うーん、ひとまず置いて進めることにしましょう。

風水害のAI活用避難システム

SIP2019

https://sip2019.go.jp/

終わったものと思いきや2期目が始まっていたのですね。知りませんでした。

2期が対象とするのは12課題。IoT や AI といったキーワードが目に留まります。その中に AI を活用した避難システム構築計画がありました。

「国家レジリエンス（防災・減災）の強化」

https://sip2019.go.jp/themes/themes-pd09.html

「ＡＩを活用して適時・的確な避難の促進を目指します～実証実験を行うモデル自治体を新たに募集～」

https://www8.cao.go.jp/cstp/stmain/20200319aisystem.html

特徴４つ

・「衛星コンステレーション」地震災害では発災後2時間で大災害か否かの判断

・「逃げ遅れゼロ」風水害ではAIでの避難指示

・「線状降水帯」水蒸気等の観測を革新することで、避難を進める。

・「経済早期復旧」広域経済早期復旧支援システムの開発

AI 活用は風水害を対象とされているのですが、土砂災害は入っていないのでしょうか？

http://www.bousai.go.jp/kyoiku/hokenkyousai/suigai.html

明記されていませんが「線状降水帯」で土砂災害を取り扱われているので省かれているのかもしれません。広島の土砂災害発生時に「線状降水帯」をよく耳にしましたので。

風水害のみの予測であれば、現時点の成果が一定の精度を有していますので、そこに機械学習 やボットを加えるのは容易でしょう。

その他、衛星利用に関し JAXA さんが入っていませんね。これまで、災害と SAR と言えば JAXA さん（とその仲間たち？）でした。今後、省庁内でシステム化し各自治体に配信するのか、それとも自治体でシステムを分散するのかはわかりません。が、人手不要にするのでしょう。小型 SAR 1000機、地震時の広域判定のみであれば、ハードルはグッと低くなります。浸水被害区域の判定を加えることも可能でしょう。

そういえば、国総研さんからも手引きが出されていました。こちらはSIPよりも細かく崩壊地を判定する趣旨ですから、求める精度も要する時間も異なるでしょう。が、将来的には地震・風水害・土砂災害含め統一されたシステムになるでしょうね。

短時間・高精度という点で判定に機械学習を加えたシステムとなることは必須です。現状はまだ十分な制度とは言えませんが、コンペ等で競われているところです。

10～15年と言わず、統合判定システムとして早期に実現してほしいところです。

赤外領域の見方

最近、赤外領域を対象とした分析やカメラをよく目にします。

web 上の製品情報を整理しました。故、汎用性に欠けますが、概ねどの光を対象にしている分析方法や装置か？という程度の理解には使えるでしょう。

区分

可視光：0.38～0.78μm

近赤外光：0.78~2.5μm

中赤外線：2.5~25μm

測定対象

近赤外分光分析：近赤外光

赤外分光分析　：中赤外光

短波赤外カメラ：近赤外光

熱赤外カメラ　：中赤外光（8～14µm）

ハイパースペクトルカメラ：紫外～可視～近赤外光(0.2~2.2μm)

LANDSAT バンド5：近赤外光(0.85 – 0.88μm)

　　　　　バンド9：近赤外光(1.36 – 1.39μm)

　　　　　バンド6：近赤外光(1.56 – 1.66μm)

　　　　　バンド7：近赤外光(2.10 – 2.30μm)

　　　　　バンド10：中赤外光(10.60-11.19μm)

　　　　　バンド11：中赤外光(11.50-12.51μm）

おまけ（SAR）

Xバンド：25-37.5mm

Cバンド：37.5-75mm

Lバンド：150-300mm

HEC-RAS 2D unsteady flow

HEC-RAS Ver.5.0.7 で2次元の非定常計算を実施。

参考にした動画はコチラ↓

https://www.youtube.com/watch?v=mUgXQucQTuE

ArcGIS で GeoTIFF を作成し、RAS Mapper に取り込む。

解析範囲をぐるりと囲み、メッシュ生成。

境界条件を設定し、計算条件を整え、実行。

解析結果を RAS Mapper で確認。

以上が大まかな流れです。

1次元に比べ、非常に手軽です。拡散波近似を扱えます。雨を降らせて地表流を生成することもできます。

何より、RAS Mapper が素晴らしい。プリ・ポストとして十分な機能を有しています。これが無償で配布されているのに日本で使われていないのは、勿体ない。日本も負けずにこのようなソフトを作って公開してほしいものです。

残念ながら、当Ver.では ２D河床変動を未サポート。別のソフトで試す必要が出てきました。

プロに紹介いただいていますので、後日試してみましょう。

CNN + KNN

CNNにKNNを組み合わせて精度を上げる手法を散見します。

https://www.kaggle.com/salokr/a-simple-cnn-knn-based-model-to-get-99-score

https://www.ijrte.org/wp-content/uploads/papers/v8i2/B1051078219.pdf

https://www.kaggle.com/lieblos/knn-cnn-ensembling

昨夜、Youtube Live で M5Stack の M5StickV を 画像数枚でAIカメラとして機能させる Brownie のデモ を見ました。アルゴリズムとして、この組み合わせを利用されたそうです。

https://qiita.com/ksasao/items/6d6bcac4c5e92fa692a2

https://github.com/ksasao/brownie

０~１を返す関数なので、その中間値を表示することができます。空のビンを０として１回写し、液体で満たしたビンを100としてもう一度写すと、液体を減らしたビンはその中間の値として表示してくれます。これは面白い。作者によれば量のみならず、距離でも使えそうとのこと。各種センサー代わりに利用可能であり、使い道がありそうです（別の視点で全く受け付けない方もいらっしゃるとは思います）。

ライブカメラでも、同様のアルゴリズムで判定させることは可能でしょう。混雑度判定や水位計の代用品としてに使えそうです。交通、農業、マーケティングなど、様々な分野で需要はあるでしょう。

センサー、カメラ、機会学習。個別技術の進展のみならず、組み合わせや他分野での動向をも追う必要に迫られているようです。いつの時代も大変です。

********************************

20200613追記

Keras で CNN をかけた結果に対し、KNN をかけると精度向上！

と思いきや、imbalanced dataでした。これ、Kを多くすれば意味なく精度は上がります。

imbalanced data で実施する際には、over sampling 等で分布数をそろえたモデルを作成し、予測する必要がありました。

ライブカメラと機械学習

国総研さんの年次報告書では、機械学習を利用した研究報告が散見されます。

道路についているライブカメラにて異常検知を目指したり、交通量調査の省力化を目指したり。特に後者では、以下のようなソフトの作成を目指されているのでしょう。

https://www.youtube.com/watch?v=XwzbFzqhF1Y

https://datafromsky.com/

自動車の場合はそれ自体が発信するデータを使えるため、交通量調査の省力化、高精度化は容易に実現できるでしょう。建設分野と自動車分野との共同で、さらなる飛躍が期待できます。

目的によっては既設のカメラでできること、できないことも出てくるでしょう。新たなカメラに切り替える場合には、カメラや映像、ネットワークの知識のみならず、機械学習の知識が要求されるでしょう（発信側の集計のみで対応できるようになれば不要になりますが）。AWS DeepLens のように機械学習とクラウド側がハードをも一体として要求してしまうかもしれません。あるいは、安価なAIボードを積んだ「AIカメラ」と呼ばれているハードがより進化し、同等以上の働きをするようになるかもしれません。

いろいろと忙しい分野になりそうです。

土石流と浮遊砂

国土技術政策総合研究所 業務 入札公示新着情報に、以下の案内が出ておりました。

「急勾配移動床水理模型実験業務」

河川、砂防の実験は1990年代に盛んだったようです。砂防のプロに聞いても、近年では珍しいのではないか？とのこと。

そういえば、今年の年次報告に河床変動の研究が掲載されていました。このモデルを確立する過程なのでしょうか。

「細粒土砂の挙動を考慮した土砂・洪水氾濫事例の再現」

http://www.nilim.go.jp/lab/bcg/siryou/2020report/index.htm

http://www.nilim.go.jp/lab/bcg/siryou/2020report/ar2020hp024.pdf

高橋モデルに浮遊砂は入っていません（侵食・堆積速度がパラスタなので、実質的には含まれているとも言えますが）。河川で使用されている掃流砂・浮遊砂式を掃流区間の式と置き換えてしまえばこの部分は対応可能でしょう。が、掃流状集合流動区間に組み込むべき式は日本にあるのでしょうか？海外の式が対応できるのかもしれませんが、実験で確かめようという流れなのでしょうか？

実際はこの業務でどのようなことを目指されているのかわかりませんが、来年の今頃にはある程度の結果が報告されているでしょう。

見逃さないようにしましょう。

土石流の Fortran コード

Hyper KANAKO の不具合が多く、代替えコードを使用することになりました。

先輩が調達した土石流の Fortran コードを読んでいたのですが、これがとても読み易い。素晴らしい書き方だと感心しながら読んでいました。

中身は高橋モデルの一部が使用されており、少し手を加えることで河川砂防技術基準に適合させることができそうでした。流れがシンプルかつ明快で、偏微分の離散化もわかりやすく、すぐに読み終わりました。

おそらく、KANAKO の計算部分も同様にシンプルなはずです。なのになぜ不具合が多いのでしょう？

おそらく後付けの GUI 側の不具合がメインなのでしょう。GUI はソフトウェアメーカーさんに任せるのが一番良いのですが、残念ながら FORUM8 さんの開発は止まっています。ユーザーが少なく採算に合わないので Ver.UP されないのだと想像されます。

エンジンと GUI を完全に切り離し、データフォーマットをマニュアルに載せておけば扱いやすくなるでしょう。

ま、Fortran コードも整備しながら、ゆっくり待ちましょう。

「高橋のモデル」平衡土砂濃度

平衡土砂濃度を３区分する水面勾配の閾値（0.03、0.138）と、平衡土砂濃度算出式の解釈に参考となる説明が以下にありました。※式2

高橋「土石流発生のメカニズム (2)」

https://www.jstage.jst.go.jp/article/sabo1973/45/3/45_3_38/_article/-char/ja/

式 (2) では, 勾配が急になればなるほど平衡砂礫濃度が大きくなるが, 密に詰めた充填濃度になれば, 実際上流動は不可能となる。実験によれば, 流動が可能な最高濃度は0.7程度で, 式 (2) がそれ以上の値となる時でも0.7あるいはもっと小さい濃度で流れる1)。一方, 式 (2) は,勾配が小さくなればいくらでも小さい砂礫濃度を与える。しかしながら, 実験によれば, 砂礫濃度は0.2～0.3以下では, 砂礫は全流動層には分散できず, 図-1 (2) に示しているような掃流状集合流動となる。図-2は通常の粗い砂礫を用いた実験によって得られた砂礫濃度と勾配の関係と, 式 (2), 掃流状集合流動に関して著者が理論的に求めた平衡濃度csooを与える式(3)2), および水山の経験式3), さらに, 芦田・高橋・水山の掃流砂量式4)による掃流砂の土砂濃度とを示している。砂礫濃度がいかに勾配に支配されて変化するかがよく分かる。図-2によれば,土石流と掃流状集合流動の境界は流路勾配10. 程度にあり, 十分長い一様な堆積層があれば, 100程度の渓床で土石流が発生できることを示している。なお, 式 (3) は以下のようである。

2013年 Hyper KANAKO 勉強会の第一部で紹介されたスライドによれば、掃流状集合流動と掃流砂の境界：0.03と土石流と掃流状集合流動の境界：0.138は、勾配1.7度、7.9度に相当します。図では3度、10度程度となっていますので、小さめの値が採用されています。その後の実験で決めた閾値だと推察されるのですが、出典がないため追いかけられません。

平衡土砂濃度（掃流状集合流動（式3）、掃流砂）の導出や出典はよくわかりませんでした。理論式だそうですが。

いずれにしても、ある理論式をある範囲で使い分けている、範囲の閾値は実験で決めているというようなことだと推察されます。

ひとまずここまで。