【iRIC】Nays2DHやNays2DFloodで計算前後で地形が変わる理由について

# はじめに
河床変動計算をしていないのにiRICでの計算前後で格子点の標高が変わっていて困惑した人は多いのではないでしょうか？
実際フォーラムにもそのような質問はよく来ますし、[よくある質問](https://i-ric.org/help/faq/)に掲載される程度にはよくある話です。

よくある質問での回答には以下のように書かれています。

> 
> Q : 読み込みデータと出力データの形状が違います。
>
> A : 
> nays2dhやnays2dfloodでは以下のように計算・出力されています。
> 1. 初期条件のマッピングでは、格子点に標高がマッピングされる
> 1. 計算時は、セル中央の平均値（周囲の4つの格子点の平均）として計算される
> 1. 結果出力時は、格子点の周囲の4つのセル中央値の平均値として、格子点に出力される
> 
> 詳細は、以下のTipsを参照してください。
> https://i-ric.org/download/iric-notes/

この記事の結論はまさにこの質問の通りなのですが、文字だけだと分かりづらそうなのと、このよくある質問がGoogleの検索結果に引っかからないので少し詳しく記事にしました。

# 図解

理由ですが、上記でも述べているようにNays2DHやNays2DFloodでは計算を行う際に格子点ではなく、セルの中心で計算をしており、計算結果として出力される値は格子点からセルの値に計算された値を更に格子点の値として計算し直した値だからです。

下の図で言うと地形データをマッピングした計算格子が一番左で、計算結果のElevationとして表示されるのが一番右の値です。

![path52-6.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2626673/d5d3fb00-dfa1-16e0-751c-00d72dba93dc.png)

# 実際の例

## 計算前
こんな格子を作ってみました。
長さ10メートル幅5メートル勾配1/10程度の凸凹とした格子です。
格子サイズは1m×1mです。

**平面図**
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2626673/632452ad-3c71-ee46-0f27-9bb23a9221e9.png)

**3D表示**
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2626673/3b69662e-10df-0772-bd14-b4028c762d1d.png)

## 計算後
では河床変動はオフで計算したあとの格子を見てみましょう。
すると以下の図のように凹凸が平滑化されてしまっていることがわかります。

**平面図**
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2626673/438e9854-6ac9-6643-b5d4-c7dc3920b7f0.png)

**3D表示**
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2626673/7bb69514-6e05-829d-3456-4cb6a58307bd.png)

**縦断比較(赤：計算格子、黒:計算結果)**
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2626673/fb4f6dae-eccf-708f-9ce6-70822328a65d.png)

## セルの値との比較

まず、最初の計算格子の格子点の値とセルの値を見てみると、以下のように上下の格子点の平均になっていることがわかります。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2626673/c54133b1-90b2-22dc-45a6-64dbd7095e05.png)

次にセルの値と計算結果の値を比較するとこちらも同様に上下のセルの値の値の平均になっていることがわかります。ただしセル→格子点の場合、上下どちらかにしかセルがない格子点は上下どちらかと同じ値になります。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2626673/a5c13907-da96-1f9a-be3c-92bc0670deac.png)

# 計算にも凹凸を反映するには
上記のように計算格子の凹凸は平滑化されてしまうことがわかりました。
ただし、計算を行う上で凹凸を活かしたまま計算したいこともありますよね。

そんなときは**格子を細かくする**とよいです。

## 例
先ほどの格子の縦断方向の格子数を3倍にしてみました。
斜面の部分にも格子点があることがわかりますね。

**平面図**
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2626673/4c4f370f-0efc-3ec9-af3d-303a2f4d6cd9.png)

**3D表示**
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2626673/4155b8ef-5383-c3dc-3b42-45351059380a.png)

では、計算後の格子をみてみましょう。
このように少しギザギザが滑らかになっていますが、先程よりも計算格子の地形に近いことがわかります。

**平面図**
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2626673/591614ae-9c6a-6570-4176-b2e468c284c3.png)

**3D表示**
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2626673/621c45cc-9995-ee33-16c6-afd5f202a1e1.png)

**縦断比較(赤：計算格子、黒:計算結果)**
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2626673/770b6f79-280c-5d2e-b647-399b605b6469.png)

さらに、格子の数を最初の10倍にすると、以下のようにほぼ凹凸を表現できています。

**縦断比較(赤：計算格子、黒:計算結果)**
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2626673/ee322b1a-95f0-da5d-b613-795afc68ff69.png)

このように計算前後で微地形が平滑化されることで地形が変わり困っている場合は格子サイズを小さくしてみてください！

# まとめ

* Nays2DHやNasy2DFloodで計算を行うと、計算の過程で細かい凹凸は平滑化される
* 凹凸などの微地形を計算でも反映したい場合は、地形変化点だけではなく、斜面にも格子点があるように格子をつくるとよい