特許 6159456 水域環境シミュレーション装置、水域環境シミュレーション方法及びプログラム

(19)【発行国】日本国特許庁(JP)

(12)【公報種別】特許公報(B1)

(11)【特許番号】6159456

(24)【登録日】2017年6月16日

(45)【発行日】2017年7月5日

(54)【発明の名称】水域環境シミュレーション装置、水域環境シミュレーション方法及びプログラム

(51)【国際特許分類】

   G06F 19/00 20110101AFI20170626BHJP

【ＦＩ】

   G06F19/00 110

【請求項の数】10

【全頁数】46

(21)【出願番号】特願2016-182373(P2016-182373)

(22)【出願日】2016年9月16日

【審査請求日】2016年11月2日

【早期審査対象出願】

(73)【特許権者】

【識別番号】516281621

【氏名又は名称】株式会社グリーン社会研究所

(74)【代理人】

【識別番号】100064908

【弁理士】

【氏名又は名称】志賀 正武

(74)【代理人】

【識別番号】100162868

【弁理士】

【氏名又は名称】伊藤 英輔

(74)【代理人】

【識別番号】100189348

【弁理士】

【氏名又は名称】古都 智

(72)【発明者】

【氏名】渡邉 正孝

【審査官】 宮地 匡人

(56)【参考文献】

【文献】 特開平０６−３４８７９０（ＪＰ，Ａ）

【文献】 特開２００５−１１５７０１（ＪＰ，Ａ）

【文献】 道路環境影響評価の技術的手法「13.動物、植物、生態系」における環境保全のための取り組みに関する事例集（平成２７年度版），国土技術政策総合研究所資料，国土交通省 国土技術政策総合研究所，２０１６年 ３月，No.906，p.3-9-9，ISSN 1346-7328，ＵＲＬ，http://www.nilim.go.jp/lab/bcg/siryou/tnn/tnn0906.htm

【文献】 中国珠江口海域環境モニタリング整備計画調査 最終報告書（概要版），国土環境株式会社，２００１年 ９月，pp.26-30，ＵＲＬ，http://open_jicareport.jica.go.jp/pdf/11667896_04.pdf

【文献】 北村 立実，霞ヶ浦流域モデルによる排出負荷量の算出と河川水質予測，地理情報システム学会講演論文集，２００８年１０月２３日，Vol.17，pp.411-414

(58)【調査した分野】（Int.Cl.，ＤＢ名）

Ｇ０６Ｆ １９／００

ＪＳＴＰｌｕｓ（ＪＤｒｅａｍＩＩＩ）

ＪＳＴ７５８０（ＪＤｒｅａｍＩＩＩ）

(57)【特許請求の範囲】

【請求項1】

シミュレーション対象の水域に設定された格子毎に、当該格子内における橋脚の投影面積及び体積をパラメータとして橋脚による抵抗の影響を示す項を含む式であって乱流エネルギーが物質の混合に及ぼす影響を示す式を含む水流のモデルを記憶するモデル記憶部と、
前記モデルに基づいて前記シミュレーション対象の水域における水流を模擬するシミュレーション処理部と、
を備える水域環境シミュレーション装置。

【請求項2】

水位、鉛直渦動粘性係数、及び、鉛直渦動拡散係数が他の数値よりも桁数の多い数値で表されたデータを記憶するデータ記憶部を備え、
前記シミュレーション処理部は、前記データ記憶部が記憶するデータを用いて前記シミュレーション対象の水域における水流を算出する、請求項１に記載の水域環境シミュレーション装置。

【請求項3】

前記シミュレーション対象の水域内の位置毎に、当該位置に応じた大きさで設定された格子間距離を記憶する格子間距離記憶部と、
前記格子間距離に基づいて格子点の座標を算出する格子点座標算出部と、
を備え、
前記シミュレーション処理部は、前記シミュレーション対象の水域における水流の模擬を、前記データ記憶部が記憶するデータを用いて、前記格子点座標算出部が算出した格子点毎に行う、
請求項２に記載の水域環境シミュレーション装置。

【請求項4】

水平方向における前記格子の行毎及び列毎に前記格子の数を記憶する格子数記憶部を備え、
前記シミュレーション処理部は、水平方向における前記格子の行及び列に応じた２次元配列を、前記格子数記憶部が記憶する前記格子の数に基づいて１次元配列に置き換えて、ループにて格子毎に当該格子内における水流を算出する、請求項１から３のいずれか１項に記載の水域環境シミュレーション装置。

【請求項5】

水平方向における前記格子の行及び列の演算を同一のループの並列処理にて実行するループ処理部を備える、請求項４に記載の水域環境シミュレーション装置。

【請求項6】

前記シミュレーション処理部は、少なくとも珪藻の炭素含有量及び渦鞭毛藻の炭素含有量の合計量に基づいて透明度を算出する、請求項１から５のいずれか一項に記載の水域環境シミュレーション装置。

【請求項7】

前記シミュレーション処理部は、水中の浮遊物を粒子の大きさで４段階に分類し、各大きさの粒子の濃度を算出する、請求項１から６のいずれか一項に記載の水域環境シミュレーション装置。

【請求項8】

前記シミュレーション処理部は、底泥中における有機物、りん化合物、及び、窒素化合物の分解及び水中への溶出のシミュレーションを行う、請求項１から７のいずれか一項に記載の水域環境シミュレーション装置。

【請求項9】

シミュレーション対象の水域に設定された格子毎に、当該格子内における橋脚の投影面積及び体積をパラメータとして橋脚による抵抗の影響を示す項を含む式であって乱流エネルギーが物質の混合に及ぼす影響を示す式を含む水流のモデルを記憶するモデル記憶部を備える水域環境シミュレーション装置が、前記モデルに基づいて前記シミュレーション対象の水域における水流を算出する水流状態算出ステップを含む、水域環境シミュレーション方法。

【請求項10】

シミュレーション対象の水域に設定された格子毎に、当該格子内における橋脚の投影面積及び体積をパラメータとして橋脚による抵抗の影響を示す項を含む式であって乱流エネルギーが物質の混合に及ぼす影響を示す式を含む水流のモデルを記憶するモデル記憶部を備えるコンピュータに、
前記モデルに基づいて前記シミュレーション対象の水域における水流を算出する水流状態算出ステップを実行させるためのプログラム。

【発明の詳細な説明】

【技術分野】

【0001】

本発明は、水域環境シミュレーション装置、水域環境シミュレーション方法及びプログラムに関する。

【背景技術】

【0002】

海洋環境のシミュレーションに関連して幾つかの技術が提案されている。例えば、特許文献１には、潮流の流動等のシミュレーションを行う際に、ＧＩＳ（Geographic Information System）を用いてメッシュの設定等を行うことが記載されている。

【先行技術文献】

【特許文献】

【0003】

【特許文献1】特開２００３−１５０８８９号公報

【発明の概要】

【発明が解決しようとする課題】

【0004】

例えば東京湾のように水域に橋が設けられている場合がある。この場合、海洋環境のシミュレーションを精度良く行うためには、橋脚が水流に与える影響をシミュレーションに反映させる必要がある。

【0005】

本発明は、海洋などの水域に橋が設けられている場合でもシミュレーションを精度良く行うことができる水域環境シミュレーション装置、水域環境シミュレーション方法及びプログラムを提供する。

【課題を解決するための手段】

【0006】

本発明の第１の態様によれば、水域環境シミュレーション装置は、シミュレーション対象の水域に設定された格子毎に、当該格子内における橋脚の投影面積及び体積をパラメータとして橋脚による抵抗の影響を示す項を含む式であって乱流エネルギーが物質の混合に及ぼす影響を示す式を含む水流のモデルを記憶するモデル記憶部と、前記モデルに基づいて前記シミュレーション対象の水域における水流を算出するシミュレーション処理部と、を備える。

【0007】

前記水域環境シミュレーション装置が、水位、鉛直動粘性係数、及び、鉛直拡散係数が他の数値よりも桁数の多い数値で表されたデータを記憶するデータ記憶部を備え、前記シミュレーション処理部は、前記データ記憶部が記憶するデータを用いて前記シミュレーション対象の水域における水流を算出するようにしてもよい。

【0008】

前記水域環境シミュレーション装置が、前記シミュレーション対象の水域内の位置毎に、当該位置に応じた大きさで設定された格子間距離を記憶する格子間距離記憶部と、前記格子間距離に基づいて格子点の座標を算出する格子点座標算出部と、を備え、前記シミュレーション処理部は、前記シミュレーション対象の水域における水流の模擬を、前記データ記憶部が記憶するデータを用いて、前記格子点座標算出部が算出した格子点毎に行うようにしてもよい。

【0009】

前記水域環境シミュレーション装置が、水平方向における前記格子の行毎及び列毎に前記格子の数を記憶する格子数記憶部を備え、前記シミュレーション処理部は、水平方向における前記格子の行及び列に応じた２次元配列を、前記格子数記憶部が記憶する前記格子の数に基づいて１次元配列に置き換えて、ループにて格子毎に当該格子内における水流を算出するようにしてもよい。

【0010】

前記水域環境シミュレーション装置が、水平方向における前記格子の行及び列の演算を同一のループの並列処理にて実行するループ処理部を備えるようにしてもよい。

【0011】

前記シミュレーション処理部は、少なくとも珪藻の炭素含有量及び渦鞭毛藻の炭素含有量の合計量に基づいて透明度を算出するようにしてもよい。
前記シミュレーション処理部は、水中の浮遊物を粒子の大きさで４段階に分類し、各大きさの粒子の濃度を算出するようにしてもよい。
前記シミュレーション処理部は、底泥中における有機物、りん化合物、及び、窒素の窒素化合物の分解及び水中への溶出のシミュレーションを行うようにしてもよい。

【0012】

本発明の第２の態様によれば、水域環境シミュレーション方法は、シミュレーション対象の水域に設定された格子毎に、当該格子内における橋脚の投影面積及び体積をパラメータとして橋脚による抵抗の影響を示す項を含む式であって乱流エネルギーが物質の混合に及ぼす影響を示す式を含む水流のモデルを記憶するモデル記憶部を備える水域環境シミュレーション装置が、前記モデルに基づいて前記シミュレーション対象の水域における水流を算出する水流状態算出ステップを含む。

【0013】

本発明の第３の態様によれば、プログラムは、シミュレーション対象の水域に設定された格子毎に、当該格子内における橋脚の投影面積及び体積をパラメータとして橋脚による抵抗の影響を示す項を含む式であって乱流エネルギーが物質の混合に及ぼす影響を示す式を含む水流のモデルを記憶するモデル記憶部を備えるコンピュータに、前記モデルに基づいて前記シミュレーション対象の水域における水流を算出する水流状態算出ステップを実行させるためのプログラムである。

【発明の効果】

【0014】

本発明によれば、水域に橋が設けられている場合でもシミュレーションを精度良く行うことができる。

【図面の簡単な説明】

【0015】

【図1】本発明の実施形態に係る水域環境シミュレーション装置の機能構成を示す概略ブロック図である。

【図2】本実施形態に係る水域環境シミュレーション装置に対する入出力データの例を示す説明図である。

【図3】本実施形態に係るモデル記憶部が記憶するシミュレーションモデルの概略構成を示す説明図である。

【図4】本実施形態に係る溶出サブモデルの概略構成の例を示す説明図である。

【図5】本実施形態に係る水質モデルが模擬する、珪藻類に関連する物質の循環の例を示す説明図である。

【図6】本実施形態における可変格子の設定例を示す説明図である。

【図7】本実施形態に係る２次元配列の変数に設定された格子間距離を用いるプログラムの例を示す説明図である。

【図8】高速化処理を行わない場合のプログラムの例を示す説明図である。

【図9】本実施形態で高速化処理を行った場合のプログラムの例を示す説明図である。

【図10】本実施形態で鉛直ループを水平ループの内側に配置したプログラムの例を示す説明図である。

【図11】本実施形態でループを並列化したプログラムの例を示す説明図である。

【発明を実施するための形態】

【0016】

以下、本発明の実施形態を説明するが、以下の実施形態は請求の範囲にかかる発明を限定するものではない。また、実施形態の中で説明されている特徴の組み合わせの全てが発明の解決手段に必須であるとは限らない。
図１は、本発明の実施形態に係る水域環境シミュレーション装置の機能構成を示す概略ブロック図である。同図に示すように、水域環境シミュレーション装置１００は、表示部１１０と、操作入力部１２０と、データ入出力部１３０と、記憶部１８０と、制御部１９０とを備える。記憶部１８０は、モデル記憶部１８１と、格子間距離記憶部１８２と、データ記憶部１８３と、格子数記憶部１８４とを備える。制御部１９０は、シミュレーション処理部１９１を備える。シミュレーション処理部１９１は、格子点座標値算出部１９２と、ループ処理部１９３とを備える。

【0017】

水域環境シミュレーション装置１００は、海洋における水流及び生物の影響をシミュレーションして、海中に溶融する物質の濃度など海洋環境を示す値を算出する。水域環境シミュレーション装置１００は、コンピュータを用いて構成される。
なお、以下ではシミュレーション対象が海洋である場合を例に説明するが、水域環境シミュレーション装置１００のシミュレーション対象はこれに限らない。例えば、水域環境シミュレーション装置１００が湖など海洋以外の水域のシミュレーションを行うようにしてもよい。

【0018】

表示部１１０は、例えば液晶パネル又はＬＥＤ（Light Emitting Diode、発光ダイオード）ネル等の表示画面を備えて各種画像を表示する。例えば、表示部１１０は、海中に溶融する物質の濃度などシミュレーション結果を表示する。
操作入力部１２０は、例えばキーボード及びマウスなどの入力デバイスを備えてユーザ操作を受ける。

【0019】

データ入出力部１３０は、水域環境シミュレーション装置１００へのデータの入出力を行う。例えば、データ入出力部１３０は、シミュレーション用の各種データを取得し、また、シミュレーション結果を示すデータを出力する。
データ入出力部１３０は、例えば、通信回路を含んで構成され、他の機器と通信を行う。あるいは、データ入出力部が、ＵＳＢ（Universal Serial Bus）ポートなど記憶デバイス用のインタフェースを備え、記憶デバイスに対してデータの入出力を行うようにしてもよい。あるいは、データ入出力部１３０が、通信回路と記憶デバイス用のインタフェースとの両方を含んで構成されていてもよい。
なお、シミュレーション用データの一部又は全部が、操作入力部１２０からのユーザ操作にて入力されてもよい。また、上記のようにシミュレーション結果の一部又は全部が、データ入出力部１３０からのデータの出力に加えて、あるいは代えて、表示部１１０に表示されるようにしてもよい。

【0020】

記憶部１８０は、各種データを記憶する。記憶部１８０が、水域環境シミュレーション装置１００の本体に内蔵されていてもよい。あるいは、記憶部１８０が、水域環境シミュレーション装置１００の本体に対して外付けの記憶デバイスを用いて構成されていてもよい。
モデル記憶部１８１は、水域環境シミュレーション装置１００が行うシミュレーションのモデル（シミュレーションモデル）を記憶する。特に、モデル記憶部１８１は、シミュレーション対象の水域に設定された格子毎に、当該格子内に含まれる橋脚の投影面積及び体積をパラメータとして橋脚が水流に及ぼす影響を算出する支配方程式を含む水流のモデルを記憶する。ここでいう格子は、シミュレーション対象の空間（水域）をメッシュに区分した各セル（Cell）である。モデル記憶部１８１が記憶するシミュレーションモデルの詳細については後述する。

【0021】

格子間距離記憶部１８２は、シミュレーション対象の水域内の位置毎に設定された格子間距離（格子点間の距離）を記憶する。
データ記憶部１８３は、シミュレーション対象の水域におけるデータ測定値（観測値）など、シミュレーションのための各種データを記憶する。特に、データ記憶部１８３は、水位、鉛直渦動粘性係数、鉛直渦動拡散係数が他の数値よりも桁数の多い数値で表されたデータを記憶する。
格子数記憶部１８４は、シミュレーション対象の水域に設定された格子について、水平方向における格子の行毎及び列毎に格子の数を記憶する。

【0022】

制御部１９０は、水域環境シミュレーション装置１００の各部を制御して各種処理を実行する。制御部１９０は、水域環境シミュレーション装置１００が備えるＣＰＵ（Central Processing Unit、中央処理装置）が、記憶部１８０からプログラムを読み出して実行することで構成される。

【0023】

シミュレーション処理部１９１は、モデル記憶部１８１が記憶するシミュレーションモデルを用いて水流及び水質等のシミュレーションを行う。特に、シミュレーション処理部１９１は、モデル記憶部１８１が記憶するシミュレーションモデルに基づいてシミュレーション対象の水域における水流を模擬する。また、シミュレーション処理部１９１は、データ記憶部１８３が記憶するデータを用いてシミュレーション対象の水域における水流を算出する。また、シミュレーション処理部１９１は、ループ処理部１９３が実行するループにて格子毎に当該格子内における水流を算出する。これにより、シミュレーション処理部１９１はシミュレーションにおける繰り返し処理を実行する。また、シミュレーション処理部１９１は、シミュレーション対象の水域内に設定された格子毎に水質及び汚濁を算出し、得られた水質及び汚濁に基づいて、格子点毎に透明度を算出する。

【0024】

格子点座標値算出部１９２は、格子間距離記憶部１８２が記憶している格子間距離に基づいて格子点の座標を算出する。
ループ処理部１９３は、プログラムの実行におけるループ処理を行う。特に、ループ処理部１９３は、シミュレーション処理部１９１が行うシミュレーションにて、格子数記憶部１０４が記憶している格子の数の範囲内で行方向のループ及び列方向のループを実行する。また、ループ処理部１９３は、格子の行及び列のうち少なくともいずれか一方を並列処理にて実行する。

【0025】

図２は、水域環境シミュレーション装置１００に対する入出力データの例を示す説明図である。同図に示すように、データ入出力部１３０が、外界境界データ、各河川日平均データ、各点源日平均データ及び初期設定データの入力を受け、シミュレーション処理部１９１が、これらのデータを用いて前処理を行う。
外界境界データは、シミュレーション対象の水域と外界との境界条件を示すデータである。外界境界データは、例えば、日月合成半日周潮Ｋ２、主太陰半日周潮Ｍ２、主太陽半日周潮Ｓ２、主太陰日周潮Ｏ１、主太陽日周潮Ｐ１、日月合成日周潮Ｋ１、周太陰長円潮Ｎ２等の分潮それぞれの遅角及び振幅等の潮汐データを含む。
また、外界境界データは、例えば、シミュレーション対象の水域と外界との境界における、化学的酸素要求量（Chemical Oxygen Demand；ＣＯＤ）、全窒素（Total Nitrogen；ＴＮ）、全りん（燐）（Total Phosphorus；ＴＰ）、塩分濃度、水温の各々に関する境界条件を含む。

【0026】

各河川日平均データは、シミュレーション対象の水域に流入する各河川との境界条件を、１日当たりの平均値で示すデータである。各河川日平均データは、例えば、シミュレーション対象の水域と河川との境界における、水等の流量、化学的酸素要求量、全窒素、全りんの各々に関する境界条件を含む。
各点源日平均データは、例えば工場からの排水など、河川以外にシミュレーション対象の水域に影響を及ぼす点源における境界条件を、１日当たりの平均値で示すデータである。各河川日平均データは、例えば、シミュレーション対象の水域と点源との境界における水等の流量、化学的酸素要求量、全窒素、全りんの各々に関する境界条件を含む。
初期設定データは、シミュレーション対象の水域内における初期値を示すデータである。初期設定データは、例えば、シミュレーション対象の水域内の所定の点における、化学的酸素要求量、全窒素、全りん、塩分濃度、水温の各々に関する初期値のデータを含む。

【0027】

シミュレーション処理部１９１は、これらのデータに対して、換算処理及び分配処理等の前処理を行う。例えば、シミュレーション処理部１９１は、これらのデータに含まれる化学的酸素要求量を、全有機炭素（Total Organic Carbon；ＴＯＣ）及び生物化学的酸素要求量（Biochemical oxygen demand；ＢＯＤ）に換算する換算処理を行う。また、シミュレーション処理部１９１は、これらのデータに含まれる全窒素（全窒素の濃度）を、懸濁態有機窒素、溶存態有機窒素、アンモニア態窒素、硝酸＋亜硝酸態窒素などの各態様(態様毎の濃度)に分配する分配処理を行う。
シミュレーション処理部１９１は、前処理にて得られたデータを記憶部１８０に記憶させ、シミュレーションの実行（モデル演算処理）に用いる。

【0028】

また、データ入出力部１３０は、さらに、地形データ、海表面風データ、日平均気象データ、計算格子・水深データ、橋脚抵抗データ、流動計算基本条件、水質計算基本条件の入力を受ける。
地形データは、例えば、海（例えば海底地形）、陸（特に海岸の地形）、河川（特に河口の位置）等の地形情報と、防波堤の位置等の防波堤情報と、埋立地の位置及び埋立時期等の埋立情報とを含む。

【0029】

海表面風データは、シミュレーション対象の水域の表面（海表面）における風の向き及び強さを示すデータである。海表面風データは、例えば、シミュレーション対象の水域の海表面における毎時風速（例えば１６方位）及び風速の情報を含む。シミュレーション対象領域が広い場合、海表面風データが、シミュレーション対象の水域を分割したエリア毎の情報を含むようにしてもよい。

【0030】

日平均気象データは、シミュレーション対象の水域における気象情報を１日当たりの平均値で示すデータである。日平均気象データは、例えば、海上における風力、風向（例えば１６方位）、気温、湿度、雲量（例えば、１０分比）、降水量、日射量の各情報を含む。シミュレーション対象領域が広い場合、海表面風データが、シミュレーション対象の水域を分割したエリア毎の情報を含むようにしてもよい。

【0031】

計算格子・水深データは、シミュレーション対象の水域に計算格子を設定するためのデータである。計算格子・水深データは、例えば、計算格子の４辺長（格子間距離）及び重心位置のうちいずれか一方又は両方など各格子の情報と、平均水深（計算格子を設定する深さを示す情報）とを含む。
格子間距離記憶部１８２が計算格子・水深データを記憶し、格子点座標値算出部１９２は、格子間距離記憶部１８２が記憶している計算格子・水深データに基づいて格子点の座標を算出する。
なお、後述するように、場所によって格子間距離が異なっていてもよい。

【0032】

橋脚抵抗データは、シミュレーション対象の水域内に橋脚等の構造物が設置されている場合に、水流に対する構造物の影響を模擬するためのデータである。橋脚抵抗データは、例えば、構造物の中心位置、形状及び大きさなど、構造物の位置を立体的に示す情報と、構造物の流水に対する抵抗係数とを含む。

【0033】

ここで、流動計算では拡散と水の流れを計算し、水質計算ではその流れで運ばれる物質の量を計算する。流動計算を行う上で必要な値を流動計算基本条件という。流動計算基本条件の例として、気象データ、外界潮汐水位データ及び河川流量データなどが挙げられる。一方、水質計算を行う上で必要な値を水質計算基本条件という。水質計算基本条件の例として、流動計算で求めた流れデータ、水温データ及び塩分データの他に、河川流入物質データが挙げられる。

【0034】

シミュレーション処理部１９１は、これらのデータを記憶部１８０に記憶させる。そして、シミュレーション処理部１９１は、データ入出力部１３０が取得したデータ及び前処理で得られたデータを記憶部１８０から読み出してシミュレーションを実行する。
シミュレーション処理部１９１は、シミュレーションの実行中に計算ログを記憶部１８０に記憶させる。ここでいう計算ログは、シミュレーション実行中における水域環境シミュレーション装置１００の状態をユーザが確認するための情報である。

【0035】

シミュレーション処理部１９１は、シミュレーションの実行にて流動三次元分布データ、生態系三次元分布データ、底質二次元分布データ、懸濁態データ及び透明度データ等のシミュレーション結果データを生成する。そして、データ入出力部１３０が、制御部１９０の制御に従ってこれらシミュレーション結果データを出力する。
流動三次元分布データは、シミュレーション対象の水域の各部における海水の状態を示すデータである。例えば、流動三次元分布データは、計算格子の各点における流速を３次元座標の各成分（後述するｕ、ｖ及びｗ）で示すデータを含む。また、例えば流動三次元分布データは、計算格子の各点における塩分濃度及び水温を示すデータを含む。

【0036】

生態系三次元分布データは、海中における物質の分布を示すデータである。例えば、生態系三次元分布データは、海中の格子点の各々について、溶存酸素、化学的酸素要求量、全有機炭素、全有機りん、全有機窒素、全りん、全窒素、懸濁態有機炭素、懸濁態有機りん、懸濁態有機窒素、溶存態有機炭素、溶存態有機りん、溶存態有機窒素、溶存態無機りん、アンモニア態窒素、硝酸＋亜硝酸態窒素、珪藻（炭素）、珪藻（りん）、珪藻（窒素）、珪藻（シリカ）、渦鞭毛藻（炭素）、渦鞭毛藻（りん）、渦鞭毛藻（窒素）、難分解性懸濁態有機炭素（陸源）、難分解性溶存態有機炭素（陸源）、珪藻（単位体積当たりの個体数）、渦鞭毛藻（単位体積当たりの個体数）、溶存態シリカ（Available Silica）、懸濁態シリカ（Particulate Biogenic Silica）、動物プランクトン（炭素）、動物プランクトン（りん）、動物プランクトン（窒素）の各濃度を示すデータを含む。

【0037】

底質二次元分布データは、海底に沈殿した物質及び海底から海中に溶出する物質の分布を示すデータである。例えば、底質二次元分布データは、海底の格子点の各々について、りん酸の溶出量、アンモニア態窒素の溶出量、炭素の沈降量、窒素の沈降量、好気層の炭素含有量、好気層の窒素含有量、嫌気層のりん含有量、嫌気層のりん含有量を、単位面積当たりの量で示すデータを含む。後述するように、好気層及び嫌気層は、海底の泥層を２つに分割した各層である。

【0038】

懸濁態データは、海水中の浮遊物質（Suspended Solid；ＳＳ）の濃度分布を示すデータである。例えば、シミュレーション処理部１９１は、浮遊物質である粒子を大きさ（直径）で４段階に分類し、各大きさの粒子の濃度を海水中の計算格子の各々について算出する。そして、懸濁態データは、各大きさの粒子の濃度を海水中の計算格子の各々について示すデータを含む。
透明度データは、シミュレーション対象の水域の各部における透明度を示すデータである。例えば、透明度データは、シミュレーション対象の水域の海水面に設定された格子点の各々における透明度を示すデータを含む。

【0039】

以下、シミュレーションモデルについて説明する。モデル記憶部１８１が、このシミュレーションモデルを記憶する。そして、シミュレーション処理部１９１が、このシミュレーションモデルを用いてシミュレーションを実行する。以下、シミュレーションモデルを単にモデルとも称する。

【0040】

図３は、モデル記憶部１８１が記憶するシミュレーションモデルの概略構成を示す説明図である。同図に示すように、モデル記憶部１８１が記憶するシミュレーションモデルは、水流モデルと、水質モデルと、溶出サブモデルと、境界条件とを含む。
水流モデルは、海水の流動を模擬するためのモデルである。さらに、水流モデルは、水温及び塩分濃度の算出式（保存式）を含む。

【0041】

水質モデルは、シミュレーション対象の水域の各部における水質を算出するためのモデルである。特に、水質モデルは、海中における有機物及び栄養塩の濃度分布を算出するための式を含む。ここでいう栄養塩は、生物が生育するために用いられる塩類である。栄養塩の例として、窒素、りん、及びシリカのいろいろな態様を挙げることができる。
溶出サブモデルは、底泥と海水との相互作用が有機物及び栄養塩の濃度分布に及ぼす影響を模擬するためのモデルである。
境界条件は、例えば海面、海底及び海岸など、シミュレーション対象の水域と周囲環境との境界及び周囲環境がシミュレーション対象の水域に及ぼす影響を示す。

【0042】

モデル記憶部１８１が記憶するシミュレーションモデルでは、シミュレーション対象の水域を区切った格子毎に（すなわち、メッシュ分割したメッシュ毎に）支配方程式を設けてモデルを表している。シミュレーション処理部１９１は、格子毎の支配方程式を用いて計算を行うことでシミュレーションを実行する。
支配方程式を求めるために、シミュレーション対象の水域にｘ軸、ｙ軸、ｚ軸の３次元直交座標を設定する。ｘ軸は東を増加方向とし、ｙ軸は北を増加方向とし、ｚ軸は上向きを増加方向とする。さらに、シミュレーション対象の水域に３次元のメッシュを設定し、メッシュの各格子点をｘ、ｙ、ｚの３次元座標で示す。関数ηを用いて自由水面をη（ｘ，ｙ，ｔ）と表す。ここで、ｔは時刻を示す。また、関数Ｈを用いて海底を−Ｈ（ｘ，ｙ）と表す。
シミュレーション処理部１９１は、格子点毎の連立方程式（支配方程式）でモデルを記載し、サンプリング時間毎に各格子点について連立方程式を解く。支配方程式においては Hydrostaticの仮定とBoussinesq近似を用いる。
以下、スカラー積を「・」、「×」のいずれでも表記する。また、演算子を表記せずに乗算対象を並列することでもスカラー積を示す。

【0043】

＜水流モデル＞
水流モデルは、海水の流動を算出するための連続方程式及び運動量方程式と、乱流を模擬するための保存式と、水温及び塩分濃度の流動を模擬するための保存式とを含む。
（連続方程式）
連続方程式は式（１）のように示される。

【0044】

【数1】

【0045】

ここで、ｕ、ｖ、ｗは、それぞれｘ、ｙ、ｚ方向に対応したアンサンブル平均流速である。以下では、アンサンブル平均流速を単に流速とも称する。
（運動量方程式）
流速ｕに関する運動方程式は式（２）のように示される。

【0046】

【数2】

【0047】

上記のように、ｔは時刻を示す。また、ρ_０は海水の基礎密度を示す。ここでいう基礎密度とは、基準となる圧力（例えば１気圧）における密度である。Ｐは圧力（水圧）を示す。Ｋ_Ｍは鉛直渦動粘性係数を示す。ｆはコリオリパラメータを示し、ｆ_ｕ、ｆ_ｖ、ｆ_ｗは、それぞれコリオリパラメータｆのｕ成分、ｖ成分、ｗ成分を示す。
また、Ｆ_ｕは、モデルのグリッド（サブグリッドスケール）によって直接表現できないsmall-scaleのプロセスによって生じる運動を示す項であり、式（３）のように示される。

【0048】

【数3】

【0049】

ここで、Ａ_Ｍは、水平渦動粘性係数を示す定数である。Ａ_Ｍの値の決定方法については後述する。
式（３）では、small-scaleのプロセスによって生じる運動を乱流拡散に類似の表現で表している。後述する式（５）、式（８）、式（１０）及び式（１２）についても同様である。
また、流速ｖに関する運動方程式は式（４）のように示される。

【0050】

【数4】

【0051】

ここで、Ｆ_ｖは、モデルのグリッドによって直接表現できないsmall-scaleのプロセスによって生じる運動を示す項であり、式（５）のように示される。

【0052】

【数5】

【0053】

また、流速ｗに関する運動方程式は式（６）のように示される。

【0054】

【数6】

【0055】

ここで、ρは海水の密度を示す。海水の密度ρについては式（１６）を参照して後述する。ｇは重力加速度を示す。
（水温及び塩分の保存式）
水温θの保存式は式（７）のように示される。

【0056】

【数7】

【0057】

ここで、Ｋ_Ｈは、鉛直渦動拡散係数を示す定数である。Ｆ_θは、モデルのグリッドによって直接表現できないsmall-scaleのプロセスによって生じる運動を示す項であり、式（８）のように示される。

【0058】

【数8】

【0059】

ここで、Ａ_Ｈは、水平渦動拡散係数を示す定数である。Ａ_Ｈの値の決定方法については後述する。
また、塩分濃度Ｓの保存式は式（９）のように示される。

【0060】

【数9】

【0061】

ここで、Ｆ_Ｓは、モデルのグリッドによって直接表現できないsmall-scaleのプロセスによって生じる運動を示す項であり、式（１０）のように示される。

【0062】

【数10】

【0063】

また、水質の項目Ｍに関する拡散方程式は、式（１１）のように示される。Ｍは、具体的には海水中に含まれる物質の濃度を示す。

【0064】

【数11】

【0065】

ここで、Ｆ_Ｍは、モデルのグリッドによって直接表現できないsmall-scaleのプロセスによって生じる運動を示す項であり、式（１２）のように示される。

【0066】

【数12】

【0067】

（乱流エネルギー、乱流スケールの保存式）
流動場に植物群落等がある場合、それが流水や物質混合に影響を及ぼし得る。これらの場合にも、全体として取り扱う必要のある空間スケールに比べて、植物群落内の個々の要素・形状・存在位置などは不規則でスケールも小さいため、それらを個々に表現することはせず、抵抗係数などのパラメータを設定することで表す。個々の要素は局所的に流れを変化させ、乱れを生成し、さらにそのエネルギー（乱流エネルギー）を散逸させる。この乱流エネルギーが物質混合に及ぼす影響を以下の２つの方程式で示す。
1つ目の方程式は、式（１３）のように示される。

【0068】

【数13】

【0069】

ここで、ｑは、ｑ^２／２で乱流エネルギーを示す。
また、ベクトルｖは、ｉｕ＋ｊｖを示す。但し、ｉは実軸を示し、ｊは虚軸を示す。従って、ベクトルｖの実数成分はｕであり、ベクトルｖの虚数成分はｖである。ｕ、ｖは、それぞれ上述したアンサンブル平均流速である。なお、ベクトルｖを、「ｖ」の上に矢印（「→」）を付して示す。また、式（１３）では、「・」はベクトルの内積を示す。
また、Ｋ_ｑは乱流エネルギーの鉛直拡散係数（鉛直方向の拡散係数）を示す。
また、Ｂ_１は、乱流クロージャー定数の１つを示す。具体的には、Ｂ_１は、乱流エネルギーの生成と減衰とが均衡している状態で経験的に求められる定数であり、例えば、Ｂ_１＝１６．６とすることができる。
また、ｌは、乱流空間スケール（turbulent length scale）を示す。ｌの値は、流れの状況に応じて実験的にまたは経験的に求められる。
Ｆ_ｑは、乱流エネルギーの水平拡散を示す。
また、２つ目の方程式は、式（１４）のように示される。

【0070】

【数14】

【0071】

ここでｑｌは乱流スケール（乱流の長さのスケール）を示す。
また、Ｅ_１、Ｅ_２は、いずれも、モデルを観測値に整合させるためのパラメータである。Ｅ_１及びＥ_２の値は、流速の観測値に基づいて実験的に決定される。例えば、Ｅ_１＝１．８、Ｅ_２＝１．３３とすることができる。
また、κはカルマン定数を示す。Ｆ_ｑｌは乱流スケールの水平拡散を示す。
また、Ｌは、海表面から海底面までの距離を表す尺度であり、式（１５）のように示される。

【0072】

【数15】

【0073】

ここで、ηは、上述した自由水面を示す関数η（ｘ，ｙ，ｔ）である。すなわち、式（１５）のηは海面のｚ座標値（平均水位（ｚ＝０）からの高さ）を示す。また、式（１５）のｈは、海底のｚ座標値を示す。
式（１３）は、平均流の摩擦抵抗に伴う乱流の発生、乱流エネルギーから熱エネルギーへの転換に伴う乱流エネルギーの減衰、乱流エネルギーの平均流による移流と拡散により乱流エネルギーの時間・空間分布が決まることを示している。式（１３）は、いわば乱流エネルギーの予報方程式である。
また、式（１４）は、乱流空間スケールを予報するための支配方程式である。乱流マクロスケールｌは乱流のうち最大規模のスケールを表すものであり、乱流エネルギーの減衰と拡散係数の決定に重要な量である。
海水の密度ρは状態方程式により式（１６）のように示される。

【0074】

【数16】

【0075】

すなわち、海水の密度ρは、水温θと塩分濃度Ｓとをパラメータとする関数にて示される。
また、鉛直渦動粘性係数Ｋ_Ｍ及び鉛直渦動拡散係数Ｋ_Ｈは、式（１７）に従い逐次計算される。

【0076】

【数17】

【0077】

ここで、ｌは上記のように乱流マクロスケールである。ｑは、上記のようにｑ^２／２で乱流エネルギーを示す。Ｓ_Ｍ、Ｓ_Ｈは安定関数である。Ｓ_Ｍ、Ｓ_Ｈは、例えばブラックス・リチャードソン数Ｒ_ｆの関数で示される。ｌ及びｑは、実験的にまたは経験的に求められる。
内湾や沿岸域等では流動場に橋脚等の構造物が存在する場合がある。この場合、構造物により流れ場に乱流が発生する。そこで、式（１３）及び式（１４）に代えて、橋脚による抵抗の影響を示す項Ｆ_Ｄｑ及びＦ_Ｄｑｌを付け加えた乱流エネルギー式及び乱流スケール式を用いる。
項Ｆ_Ｄｑを付け加えた乱流エネルギー式は、式（１８）のように示される。

【0078】

【数18】

【0079】

項Ｆ_Ｄｑは式（１９）のように示される。

【0080】

【数19】

【0081】

ここで、Ｃ_ｆｑは実験定数（実験により定める定数）である。このように、乱流エネルギー式を示す式（１８）に、橋脚による抵抗の影響を示す項Ｆ_Ｄｑを付け加えることで、構造物が流動場に与える影響を、より正確に表現することができる。
また、Ｆ_Ｄｘは式（２０）のように示される。

【0082】

【数20】

【0083】

ここで、Ｃ_ｄは、杭による流水への抵抗係数を示す。
また、Ａは、１格子内の杭の投影面積を示す。1格子内の杭の投影面積とは、1格子内に入る杭の直径と全長の積である。例えば、太さの異なる主杭及び斜杭からなるジャケット構造の場合、１格子内の杭の投影面積Ａを式（２１）で求める。

【0084】

【数21】

【0085】

また、Δｘはｘ方向（ｘ軸方向）の格子間隔を示す。Ｄは層厚（鉛直方向の格子間隔）を示す。ここで、水域環境シミュレーション装置１００では、多層レベルモデル（鉛直方向に複数の格子が設定されるモデル）が用いられており、シミュレーション演算部１９１は、各層を構成する格子毎に支配方程式の演算を行ってシミュレーション結果を算出する。そこで、シミュレーション演算部１９１は、は、杭の投影面積及び体積を、1格子内の杭の投影面積及び体積毎（すなわち、各層の1格子毎）に求める。
Ｃ_Ｍは付加質量係数を示す。ここで、海域に柱を立てた場合、形状抵抗と水の加速度運動に伴う圧力勾配が存在する。そして、この柱には、水の加速度運動に伴う圧力勾配に体積を乗算した大きさに比例する力が働く。この係数を付加質量係数Ｃ_Ｍと称する。例えば、円柱の場合、負荷質量計数Ｃ_Ｍの値として２．０を用いることができる。
また、Ｖは１格子の体積を示す。
また、Ｖ_ｄは１格子内の杭の体積を示す。1格子内の杭の体積とは、1格子内に入る杭の断面積と全長との積である。例えば、太さの異なる主杭及び斜杭からなるジャケット構造の場合、１格子内の杭の体積Ｖ_ｄを式（２２）で求める。

【0086】

【数22】

【0087】

また、式（１９）のＦ_Ｄｙは式（２３）のように示される。

【0088】

【数23】

【0089】

また、項Ｆ_Ｄｑｌを付け加えた乱流スケール式は、式（２４）のように示される。

【0090】

【数24】

【0091】

項Ｆ_Ｄｑｌは式（２５）のように示される。

【0092】

【数25】

【0093】

このように、乱流スケール式を示す式（２４）に、項Ｆ_Ｄｑｌを付け加えることで、構造物が流動場に与える影響を、より正確に表現することができる。
また、水平渦動粘性係数Ａ_Ｍの値は式（２６）に従って求められる。

【0094】

【数26】

【0095】

ここで、Δｙはｙ方向（ｙ軸方向）の格子間隔を示す。
また、水平渦動拡散係数Ａ_Ｈの値は式（２７）に従って求められる。

【0096】

【数27】

【0097】

また、定数Ｃ_Ｍ及びＣ_Ｈの値は格子のサイズに応じて設定する。例えば、１ｋｍ（キロメートル）×１ｋｍのグリッドに対して、Ｃ_Ｍ＝０．０１、Ｃ_Ｈ＝０．０１に設定する。

【0098】

＜境界条件＞
（自由水面での境界条件の式）
自由水面（海面）での境界条件を以下の第１式〜第３式のように設定する。
自由水面での境界条件の第１式は、式（２８）のように示される。

【0099】

【数28】

【0100】

ここで、（τ_０ｘ，τ_０ｙ）は水面での風による応力を示す。（τ_０ｘ，τ_０ｙ）の値は、例えば気象データ及び理論値に基づいて予め設定される。
自由水面での境界条件の第２式は、式（２９）のように示される。

【0101】

【数29】

【0102】

ここで、Ｈ’は水面での純熱フラックスを示す。Ｓ’は水面での塩分フラックスを示す。ここでいうフラックス（Flux）は、流束、すなわち、単位時間かつ単位面積あたりに流れる量（移動量）である。
水面での純熱フラックスＨ’は式（３０）のように示される。

【0103】

【数30】

【0104】

ここで、φ_Ｓは太陽から地表面に到達する太陽放射量（日射量）を示す。φ_ａは、日光が大気中の空気分子、エアロゾル及び雲粒子等により吸収及び散乱を受けた後、地表面に到達する大気放射量を示す。φ_ｂｒは地表面が黒体として放つ赤外放射量を示す。φ_Ｓ、φ_ａ及びφ_ｂｒの値は、例えば気象データ（気象観測データまたは気象推定値）及び理論値に基づいて予め設定される。
また、φ_Ｓｒは地表面で反射される太陽放射量（日射量）を示し、式（３１）のように示される。

【0105】

【数31】

【0106】

ここで、φ_Ｓｃは晴天時の日射量を示す。Ｃは雲量を示す。φ_Ｓｃ及びＣの値は、例えば気象データ（気象観測データまたは気象推定値）及び理論値に基づいて予め設定される。
また、φ_ａｒは地表面で反射される大気放射量を示し、式（３２）のように示される。

【0107】

【数32】

【0108】

ここで、ｅは水蒸気圧を示す。Ｔ_ａは気温を示す。ｅ及びＴ_ａの値は、例えば気象データ（気象観測データまたは気象推定値）及び理論値に基づいて予め設定される。
また、φ_ｅは蒸発散に伴う潜熱輸送量を示し、式（３３）のように示される。

【0109】

【数33】

【0110】

ここで、Ｗ_ｚは水面上高さｚでの風速を示す。ｅ_ｓは水面上の気温に対する飽和水蒸気圧を示す。ｅ_ｚは海面上高さｚにおける水蒸気圧を示す。Ｔ_Ｓは水表面の水温を示す。
また、φ_ｃは対流熱伝達量を示し、式（３４）のように示される。

【0111】

【数34】

【0112】

ここで、Ｔ_ｚは水表面上Ｚメートル（ｍ）の高さでの気温を示す。例えば、気温Ｔ_ｚ及び風速Ｗ_ｚとして、水面上１〜２ｍ程度の高さでの値を用いる。Ｗ_ｚ、ｅ_ｓ、ｅ_ｚ、Ｔ_Ｓ及びＴ_ｚの値は、例えば気象データ（気象観測データまたは気象推定値）及び理論値に基づいて予め設定される。
自由水面での境界条件の第３式は、式（３５）のように示される。

【0113】

【数35】

【0114】

ここで、水面での塩分フラックスＳ’は式（３６）のように示される。

【0115】

【数36】

【0116】

ここで、［Ｅ’−Ｐ’］は、蒸発−降雨による水面での純淡水量フラックスを示す。Ｓ（０）は水面での塩分濃度を示す。［Ｅ’−Ｐ’］の値は、例えば気象データ（気象観測データまたは気象推定値）及び理論値に基づいて予め設定される。Ｓ（０）の値は、例えば実測値に基づいて予め定められる。
（海底での境界条件の式）
また、海底における応力、純熱フラックス、塩分フラックスについて以下のように境界条件を設定する。海底では、ｚ＝−Ｈである。また、海底ではθ及びＳの鉛直勾配はゼロであり、境界を通しての熱の移動及び塩分の移動は無い。
海底における境界条件の第１式は、式（３７）のように示される。

【0117】

【数37】

【0118】

ここで、（τ_ｂｘ，τ_ｂｙ）は海底における応力を示す。この（τ_ｂｘ，τ_ｂｙ）は式（３８）のように対数近似則で与えられる。

【0119】

【数38】

【0120】

ここで、Ｃ_Ｄは抗力係数を示す。抗力係数Ｃ_Ｄの値は０．００２５、又は、式（３９）による値のいずれか大きい方の値を用いる。

【0121】

【数39】

【0122】

ここで、Ｈ（ｘ，ｙ）は海底地形（海底面のｚの値）を表す。Ｈ（ｘ，ｙ）を単にＨとも表記する。また、ｚ_ｂは海底に最も近い計算格子を示す。ｕ_ｂは、計算格子ｚ_ｂでのｘ方向の流速を示す。ｖ_ｂは、計算格子ｚ_ｂでのｙ方向の流速を示す。
また、ｌｎは自然対数を示す。ｚ_０の値は海底の局所的な祖度に依存する。例えば、ｚ_０の値を１センチメートル（ｃｍ）とする。
また、海底における境界条件の第２式は、式（４０）のように示される。

【0123】

【数40】

【0124】

海底における境界条件の第３式は、式（４１）のように示される。

【0125】

【数41】

【0126】

ここで、Ｅｓは、物質粒子の沈降量を示す。
上述した自由水面及び海底面に加え、河川、排水を流す点源等も境界条件にてモデル化することができる。境界条件は、定数又は関数値等の数値で示されていてもよいし、方程式で示されていてもよいし、これらの組み合わせで示されていてもよい。

【0127】

＜座標変換＞
内湾のように沿岸の急激な水深変化を表現するには直交座標（ｘ，ｙ，ｚ，ｔ）は有効ではない。そこで、上述のように設定した各式について鉛直方向に新しい座標（ｘ^＊，ｙ^＊，ｚ^＊，ｔ^＊）を導入し、全ての方程式の座標変換を行う。具体的には、式（４２）のように変換を行う。

【0128】

【数42】

【0129】

ここで、Ｄ＝Ｈ＋ηである。このため、ｚ＝η（すなわち、水面）でσ＝０となる。また、ｚ＝−Ｈ（すなわち、海底）でσ＝−１となる。

【0130】

＜溶出サブモデル＞
生態系モデルに閉鎖性海域や浅海域で重要な過程である底泥と海水との相互作用を溶出サブモデルとして取り込む。例えば溶出サブモデルの要素として、懸濁態有機炭素、懸濁態有機窒素、ＮＨ_４、ＮＯ_３、懸濁態有機りん、ＰＯ_４、珪素（シリカ）、硫化水素（海水域）及び溶存酸素を採用することができるが、これに限らない。
水域環境シミュレーション装置１００に用いる溶出サブモデルとして、例えばＥＰＡ（United States Environmental Protection Agency）で開発されたモデルなど、既存のモデルを用いることができる。

【0131】

モデルでは底泥中を好気層と嫌気層の２層に分割する。上層は直上水の溶存酸素濃度により好気層もしくは嫌気層のいずれかになる。一方、下層は常に嫌気層となる。水中より有機物が底泥に沈降し、底泥中での有機物の無機化（分解）を３種類の速度（易分解性、難分解性、無分解性）で示す。この分解により有機物が栄養塩に回帰し、その過程の中で酸素要求（ＳＯＤ）が起こることを示している。また、回帰した栄養塩は海水中に溶出する。

【0132】

図４は、溶出サブモデルの概略構成の例を示す説明図である。図４に示すように溶出サブモデルでは、底泥と海水との間での有機物、栄養塩及び溶存酸素の移動が示される。また、上述したように、底泥は好気層と嫌気層とに分けられる。溶出サブモデルでは好気層と嫌気層との間の有機物及び栄養塩の移動が示される。また、溶出サブモデルでは、有機物から栄養塩への無機化（分解）、及び、埋積物としての堆積が示される。
以下、数式を参照しながら溶出サブモデルについて説明する。

【0133】

（有機物）
有機物の分解について、３種類の分解速度のクラスに分類する。Ｇ_１は、２０日で半減するクラス、Ｇ_２は、１年で半減するクラス、Ｇ_３は、堆積層へ沈降するまで分解を受けないクラスである。
好気層における粒子態有機炭素の増減は、式（４３）のように示される。

【0134】

【数43】

【0135】

ここで、「（Ｔ−２０）」は、底泥溶出速度試験の標準水温を２０℃にしていることを示している。
好気層における粒子態有機窒素の増減は、式（４４）のように示される。

【0136】

【数44】

【0137】

好気層における粒子態有機りんの増減は、式（４５）のように示される。

【0138】

【数45】

【0139】

ここで、Ｈ_１は好気層の厚さを示す。Ｈ_１の値は例えば０．１センチメートル（ｃｍ）にする。
また、ｉは、１、２、３のうちのいずれかの値を取る。ｉの値１、２、３は、それぞれクラスＧ_１、Ｇ_２、Ｇ_３に対応付けられる。ＰＯＣ_，１、ＰＯＮ_，１、ＰＯＰ_，１は、それぞれ好気層での粒子態有機炭素濃度、窒素濃度、りん濃度（単位：ｍｇ／ｌ）を示す。ＰＯＣ_，０、ＰＯＮ_，０、ＰＯＰ_，０は、それぞれ海水中での粒子態有機炭素濃度、窒素濃度、りん濃度を示す。ＰＯＣ_，２、ＰＯＮ_，２、ＰＯＰ_，２は、それぞれ嫌気層での粒子態有機炭素濃度、粒子態有機窒素濃度、粒子態有機りん濃度を示す。なお、ＰＯＣ、ＰＯＮ、ＰＯＰは、それぞれ粒子態有機炭素濃度、粒子態有機窒素濃度、粒子態有機りん濃度を示す。
また、ｆ_{ＰＯＣ，ｉ}、ｆ_{ＰＯＮ，ｉ}、ｆ_{ＰＯＰ，ｉ}は、それぞれ粒子態有機炭素濃度、粒子態有機窒素濃度、粒子態有機りん濃度のクラスＧ_ｉへの配分係数を示す。濃度に分配係数を乗算することで、クラス毎の濃度を得られる。分配係数は、例えば予め定数で与えられる。粒子の大きさによって沈降速度が異なるため、クラス毎に分けることでより高精度に演算を行える。

【0140】

例えば、ｆ_{ＰＯＣ，１}＝０．６５，ｆ_{ＰＯＣ，２}＝０．２０，ｆ_{ＰＯＣ，３}＝０．１５とする。
また、例えば、ｆ_{ＰＯＮ，１}＝０．６５，ｆ_{ＰＯＮ，２}＝０．２８，ｆ_{ＰＯＮ，３}＝０．０７とする。
また、例えば、ｆ_{ＰＯＰ，１}＝０．６５，ｆ_{ＰＯＰ，２}＝０．２５５，ｆ_{ＰＯＰ，３}＝０．０９５とする。

【0141】

また、ｋ_{ＰＯＣ，ｉ}、ｋ_{ＰＯＮ，ｉ}、ｋ_{ＰＯＰ，ｉ}は、それぞれＰＯＣ_，ｉ、ＰＯＮ_，ｉ、ＰＯＰ_，ｉの分解率を示す。ｋ_{ＰＯＣ，１}、ｋ_{ＰＯＮ，１}、ｋ_{ＰＯＰ，１}の値は、例えば０．０３５（単位：ｄａｙ^−ｌ）に設定する。
また、ＰＯＣ_，２、ＰＯＮ_，２、ＰＯＰ_，２の分解率を、例えば０．００１８（単位：ｄａｙ^−１）に設定する。
また、ＰＯＣ_，３、ＰＯＮ_，３、ＰＯＰ_，３の分解率を、例えば０．０（単位：ｄａｙ^−ｌ）に設定する。

【0142】

また、θ_{ＰＯＣ，ｉ}、θ_{ＰＯＮ，ｉ}、θ_{ＰＯＰ，ｉ}は、それぞれ分解の温度係数を示す。
例えば、θ_{ＰＯＣ，１}、θ_{ＰＯＮ，１}、θ_{ＰＯＰ，１}の値をいずれも１．１００に設定する。
また、例えば、θ_{ＰＯＣ，２}、θ_{ＰＯＮ，２}、θ_{ＰＯＰ，２}の値をいずれも１．１５０に設定する。

【0143】

また、例えばθ_{ＰＯＣ，３}、θ_{ＰＯＮ，３}、θ_{ＰＯＰ，３}の値をいずれも０（なし）に設定する。
また、Ｗ_２は、堆積速度（単位：ｍ／ｄａｙ（メートル／日））を示す。Ｗ_ｎｅｔは、海水から底泥への沈降速度（単位：ｍ／ｄａｙ）を示す。Ｗ_１２は、混合物質移動係数を示す。堆積速度Ｗ_２及び海水から底泥への沈降速度Ｗ_ｎｅｔは、例えば定数で予め設定される。混合物質移動係数Ｗ_１２の算出については、式（５８）を参照して後述する。

【0144】

式（４３）の右辺の第１項は、分解による粒子態有機炭素濃度の減少を示す。第２項は、好気層から嫌気層への沈降による粒子態有機炭素濃度の減少を示す。第３項は、海水中から好気層への沈降による粒子態有機炭素濃度の増加を示す。第４項は、好気層と嫌気層との混合による粒子態有機炭素濃度の増加を示す。
式（４４）、（４５）の右辺の各項も、粒子態有機炭素濃度をそれぞれ粒子態有機窒素濃度、粒子態有機りん濃度と読み替えて、式（４３）の場合と同様である。

【0145】

また、嫌気層における粒子態有機炭素の増減は、式（４６）のように示される。

【0146】

【数46】

【0147】

また、嫌気層における粒子態有機窒素の増減は、式（４７）のように示される。

【0148】

【数47】

【0149】

また、嫌気層における粒子態有機りんの増減は、式（４８）のように示される。

【0150】

【数48】

【0151】

ここで、Ｈ_２は嫌気層の厚さを示す。Ｈ_２の値は例えば１０センチメートル（ｃｍ）にする。
式（４６）の右辺の第１項は、分解による粒子態有機炭素濃度の減少を示す。第２項は、埋積物としての堆積による粒子態有機炭素濃度の減少を示す。埋積物として堆積した有機物等は、好気層又は海中へは拡散しないと考えられるため、嫌気層における濃度から除外する。式（４６）の第３項は、好気層から嫌気層への沈降による粒子態有機炭素濃度の増加を示す。第４項は、好気層と嫌気層との増加による粒子態有機炭素濃度の減少を示す。
式（４７）、（４８）の右辺の各項も、粒子態有機炭素濃度をそれぞれ粒子態有機窒素濃度、粒子態有機りん濃度と読み替えて、式（４６）の場合と同様である。
また、保存式の第１式は、式（４９）のように示される。

【0152】

【数49】

【0153】

また、保存式の第２式は、式（５０）のように示される。

【0154】

【数50】

【0155】

ここで、Ｃ_Ｔ０、Ｃ_Ｔ１、Ｃ_Ｔ２は、有機物の濃度を示す。添え字の値「０」、「１」、「２」は、それぞれ海水中、好気層、嫌気層を示す。Ｋ_Ｌ１２は、好気層と嫌気層との間の拡散係数を示す。ｆ_ｄ０は、海域の溶存態分画（溶存分画係数）を示す。
κ_１^２については式（７２）を参照して後述する。κ_２は、好気層の反応速度（単位：ｍ／ｄａｙ）を示す。κ_２の値は、例えば０．０に設定される。
また、式（４９）及び式（５０）では、溶存態分画ｆ_ｄと粒子態分画ｆ_ｐとに分けて計算する。好気層の溶存態分画ｆ_ｄ１は式（５１）のように示される。

【0156】

【数51】

【0157】

好気層の粒子態分画ｆ_ｐ１は式（５２）のように示される。

【0158】

【数52】

【0159】

嫌気層の溶存態分画ｆ_ｄ２は式（５３）のように示される。

【0160】

【数53】

【0161】

嫌気層の粒子態分画ｆ_ｐ２は式（５４）のように示される。

【0162】

【数54】

【0163】

ここで、ｍ_１、ｍ_２は、それぞれ好気層、嫌気層における粒子態濃度を示す。π_１、π_２は、それぞれ好気層、嫌気層における分配係数を示す。ｍ_１、ｍ_２、π_１、π_２の値は、例えば定数で予め設定される。
また、海水−底泥間の物質移動係数ｓは、式（５５）のように示される。

【0164】

【数55】

【0165】

ここで、ＳＯＤは、底泥による酸素消費速度を示す。［Ｏ２（Ｏ）］は水中酸素を示す。
また、Ｊ_Ｔ１は式（５６）のように示される。

【0166】

【数56】

【0167】

Ｊ_Ｔ２は式（５７）のように示される。

【0168】

【数57】

【0169】

ここで、ｋ_{ＰＯＭ，ｉ}は、懸濁態有機物の分解速度を示す。ｉの値「１」、「２」は、それぞれクラスＧ_１、Ｇ_２に対応付けられる。θ_{ＰＯＭ，ｉ}^{（Ｔ−２０）}は、分解の温度係数を示す。２０℃が標準温度に設定されている。ＰＯＭ_ｉ１、ＰＯＭ_ｉ２は、それぞれ懸濁態有機物の好気層、嫌気層における濃度を示す。
また、混合物質移動係数Ｗ_１２は式（５８）のように示される。

【0170】

【数58】

【0171】

ここで、Ｄ_ｐは粒子混合拡散係数（単位：ｍ^２／ｄａｙ）を示す。粒子混合拡散係数Ｄ_ｐの値は、例えば１．２Ｅ−０４（浮動小数点表示）に設定される。
また、θ_ＤｐはＤ_ｐの温度係数を示す。温度係数θ_Ｄｐの値は、例えば１．１１７に設定される。
また、Ｇ_{ＰＯＣ，Ｒ}は好気層での水温２０℃におけるＰＯＣ基準濃度（単位：ｍｇＣ／ｇ）を示す。Ｇ_{ＰＯＣ，Ｒ}の値は例えば０．１に設定される。
また、Ｋ_Ｍ，Ｄｐは、酸素に対する粒子混合半飽和定数（単位：ｍｇ／Ｌ）を示す。Ｋ_Ｍ，Ｄｐの値は、例えば４．０に設定される。
好気層と嫌気層との間の溶存態物質の混合はベントスの活動に伴う混合によって分子拡散よりも大きくなる。この効果を式（５９）で組み込む。

【0172】

【数59】

【0173】

ここで、Ｄ_ｄは間隙水拡散係数（単位：ｍ^２／ｄａｙ）を示す。Ｄ_ｄの値は例えば１．０Ｅ−３（浮動小数点表示）に設定される。
また、θ_ＤｄはＤｄの温度係数を示す。θ_Ｄｄの値は例えば１．０８に設定される。

【0174】

（アンモニア）
好気層におけるアンモニアの溶出は式（６０）のように示される。

【0175】

【数60】

【0176】

ここで、ＮＨ_４（０）、ＮＨ_４（１）、ＮＨ_４（２）は、それぞれ海水中、好気層、嫌気層におけるアンモニアの濃度を示している。すなわち、カッコ内の０、１、２は、それぞれ海水中、好気層、嫌気層を示している。上述したように、ｆ_ｄ０は、海域の溶存態分画係数を示す。
嫌気層におけるアンモニアの溶出は式（６１）のように示される。

【0177】

【数61】

【0178】

ここで、Ｊ_Ｎ１は式（６２）のように示される。

【0179】

【数62】

【0180】

Ｊ_Ｎ２は式（６３）のように示される。

【0181】

【数63】

【0182】

また、アンモニアに関しては、ｆ_ｄ１を表す式として式（５１）に代えて式（６４）を用いる。

【0183】

【数64】

【0184】

ここで、ＮＨ_４は、アンモニアの濃度を示す。ｆ_ｐ１については式（５２）を用いる。アンモニアに関しては、式（５２）のｆ_ｄ１の値は、式（６４）に基づく。
また、アンモニアに関しては、ｆ_ｄ２を表す式として式（５３）に代えて式（６５）を用いる。

【0185】

【数65】

【0186】

ｆ_ｐ２については式（５４）を用いる。アンモニアに関しては、式（５４）のｆ_ｄ２の値は、式（６５）に基づく。
また、アンモニアが受ける硝化作用は式（６６）のように示される。

【0187】

【数66】

【0188】

ここで、κ_{ＮＨ４，１}は、硝化反応速度（単位：ｍ／ｄａｙ）を示す。κ_{ＮＨ４，１}の値は例えば０．１３１に設定される。
また、θ_ＮＨ４は、硝化反応の温度係数を示す。θ_ＮＨ４の値は例えば１．１２３に設定される。
また、Ｋ_{Ｍ，ＮＨ４}は、アンモニアの硝化反応半飽和定数（単位：ｍｇＮ／ｍ^３）を示す。Ｋ_{Ｍ，ＮＨ４}の値は、例えば７２８に設定される。
また、θ_{ＫＭ，ＮＨ４}は、硝化反応半飽和定数の温度係数を示す。θ_{ＫＭ，ＮＨ４}の値は、例えば１．１２５に設定される。
また、Ｋ_{Ｍ，ＮＨ４，Ｏ２}は、酸素に対する硝化反応半飽和定数（単位：ｍｇＯ_２／Ｌ）を示す。Ｋ_{Ｍ，ＮＨ４，Ｏ２}の値は、例えば０．３７に設定される。

【0189】

また、π_１ＮＨ_４は、好気層の分画係数（単位：Ｌ／ｋｇ）を示す。π_１ＮＨ_４の値は例えば１．０に設定される。
また、π_２ＮＨ_４は、嫌気層の分画係数（単位：Ｌ／ｋｇ）を示す。π_２ＮＨ_４の値は例えば１．０に設定される。
また、Ｊ_Ｎ１は、好気層の分解窒素フラックス（単位：Ｍｇ・Ｎ／ｍ^２／ｄａｙ）を示す。Ｊ_Ｎ１の値は例えば０．０に設定される。また、Ｊ_Ｎ２は、嫌気層の分解窒素フラックス（単位：Ｍｇ・Ｎ／ｍ^２／ｄａｙ）を示す。Ｊ_Ｎ２をＪ_Ｎとも表記する。海域や計算軌間等の条件に応じてパラメータを設定している。
なお、表記を明確にするために変数を［］でくくる場合がある。例えば、式（６６）で、表記を明確にするために「ＮＨ_４（１）」を［］でくくっている。

【0190】

（硝酸）
好気層における硝酸の溶出は式（６７）のように示される。

【0191】

【数67】

ここで、ＮＯ_３（０）、ＮＯ_３（１）、ＮＯ_３（２）は、それぞれ海水中、好気層、嫌気層における硝酸の濃度を示している。

【0192】

嫌気層における硝酸の溶出は式（６８）のように示される。

【0193】

【数68】

【0194】

ここで、Ｊ_ＮＯ３１は式（６９）のように示される。

【0195】

【数69】

【0196】

ここで、Ｊ［ＮＨ_４］は、海水への溶出フラックスを示す。
また、硝酸に関しては、κ_１^２は式（７０）のように示される。

【0197】

【数70】

【0198】

また、硝酸に関しては、κ_２は式（７１）のように示される。

【0199】

【数71】

【0200】

ここで、κ_{ＮＯ３，１}は好気層の脱窒反応速度（単位：ｍ／ｄａｙ）を示す。κ_{ＮＯ３，１}の値は例えば０．１０に設定される。
また、κ_{ＮＯ３，２}は嫌気層の脱窒反応速度（単位：ｍ／ｄａｙ）を示す。κ_{ＮＯ３，２}の値は例えば０．２５に設定される。
また、θ_ＮＯ３は脱窒の温度係数を示す。θ_ＮＯ３の値は例えば１．０８に設定される。
また、硝酸に関しては、ｆ_ｄ１を表す式として式（５１）に代えて式（７２）を用いる。

【0201】

【数72】

【0202】

ｆ_ｐ１については式（５２）を用いる。硝酸に関しては、式（５２）のｆ_ｄ１の値は、式（７２）に基づく。
また、硝酸に関しては、ｆ_ｄ２を表す式として式（５３）に代えて式（７３）を用いる。

【0203】

【数73】

【0204】

ｆ_ｐ２については式（５４）を用いる。硝酸に関しては、式（５４）のｆ_ｄ２の値は、式（７３）に基づく。

【0205】

（硫化物）
好気層における硫化物の溶出は式（７４）のように示される。

【0206】

【数74】

【0207】

ここで、Ｓ（０）、Ｓ（１）、Ｓ（２）は、それぞれ海水中、好気層、嫌気層における硫黄の濃度を示している。
嫌気層における硫化物の溶出は式（７５）のように示される。

【0208】

【数75】

【0209】

ここで、硫黄に関しては、ｆ_ｄ１を表す式として式（５１）に代えて式（７６）を用いる。

【0210】

【数76】

【0211】

ｆ_ｐ１については式（５２）を用いる。硫黄に関しては、式（５２）のｆ_ｄ１の値は、式（７６）に基づく。
また、硫黄に関しては、ｆ_ｄ２を表す式として式（５３）に代えて式（７７）を用いる。

【0212】

【数77】

【0213】

ｆ_ｐ２については式（５４）を用いる。硫黄に関しては、式（５４）のｆ_ｄ２の値は、式（７７）に基づく。
ここで、π_１Ｓ及びπ_２Ｓは、分配係数を示す。π_１Ｓ、π_２Ｓの値は、例えば定数で予め設定される。
ここでＪ_ｓ２は、式（７８）のように示される。

【0214】

【数78】

【0215】

また、α_{Ｏ２，ＮＯ３}は、脱窒作用による酸素消費量（単位：ｇ Ｏ_２／ｍ^３）を示す。α_{Ｏ２，ＮＯ３}の値は例えば２．８５７１に設定される。α_Ｏ２，Ｃは、脱窒作用による有機炭素消費量を示す。ＪＣは、海への炭素フラックスを示す。
また、硫黄に関してはκ_１^２は、式（７９）のように示される。

【0216】

【数79】

【0217】

ここで、κ_{Ｈ２Ｓ，ｄ１}は、溶存態硫化物の酸化反応速度（単位：ｍ／ｄａｙ）を示す。κ_{Ｈ２Ｓ，ｄ１}の値は例えば０．２０に設定される。
また、κ_{Ｈ２Ｓ，ｐ１}は、粒子態硫化物の酸化反応速度（単位：ｍ／ｄａｙ）を示す。κ_{Ｈ２Ｓ，ｐ１}の値は例えば０．４０に設定される。
また、θ_Ｈ２Ｓは、硫化物酸化の温度係数を示す。θ_Ｈ２Ｓの値は、例えば１．０８に設定される。

【0218】

また、Ｋ_{Ｍ，Ｈ２Ｓ，Ｏ２}は、硫化物酸化定数（単位：ｍｇ Ｏ２／Ｌ）を示す。Ｋ_{Ｍ，Ｈ２Ｓ，Ｏ２}の値は、例えば４．０に設定される。
また、π_１Ｓは好気層の分画係数（単位：Ｌ／ｋｇ）を示す。π_１Ｓの値は例えば１００に設定される。
また、π_２Ｓは嫌気層の分画係数（単位：Ｌ／ｋｇ）を示す。π_２Ｓの値は例えば１００に設定される。

【0219】

（溶存酸素）
酸素消費量ＳＯＤ（単位：ｇ Ｏ_２／ｍ^３）は、式（８０）のように、メタン酸化による酸素要求量（ＣＳＯＤ、単位：ｇ Ｏ_２／ｇＮ）と、硝化作用に伴う酸素消費量（ＮＳＯＤ）との和で示される。

【0220】

【数80】

【0221】

また、底泥中の有機炭素分解に伴う酸素要求量ＣＳＯＤは、式（８１）のように示される。

【0222】

【数81】

【0223】

また、底泥中の硝化作用に伴う酸素消費量ＮＳＯＤは、式（８２）のように示される。

【0224】

【数82】

【0225】

ここで、［ΣＨ_２Ｓ（１）］は、好気層の総硫化物濃度（単位：ｇ Ｏ２／ｍ３）を示す。
また、α_{Ｏ２，ＮＨ４}は、硝化による酸素消費量（単位：ｇ Ｏ２／ｇＮ）を示す。α_{Ｏ２，ＮＨ４}の値は、例えば４．５７１４に設定される。θ_ＮＨ４^{（Ｔ−２０）}は、水温２０℃時の硝化の温度定数を示す。θ_ＮＨ４^{（Ｔ−２０）}の値として、例えば１．０７６〜１．１２７の範囲内の定数を用いることができる。ｋ_{ＮＨ４，１}^２は、好気層での消化率定数の値を示す。ｋ_{ＮＨ４，１}^２の値として、例えば０．８９４を用いることができる。ｋ_{Ｍ，ＮＨ４}は、硝化の半飽和定数を示す。ｋ_{Ｍ，ＮＨ４}の値として、例えば０．６３〜１．１９の範囲内の定数を用いることができる。θ_{Ｍ，ＮＨ４}^{（Ｔ＝２０）}は、水温２０℃時の温度定数を示す。Ｋ_{Ｏ２，ＮＨ４}は、嫌気層での酸素濃度の半飽和定数を示す。Ｋ_{Ｏ２，ＮＨ４}の値として、例えば０．２５〜２．０の範囲内の定数を用いることができる。

【0226】

（りん）
りんの溶出は、式（８３）のように示される。

【0227】

【数83】

【0228】

ここで、ＰＯ_４（１）は、好気層でのりん酸の濃度を示す。ＰＯ_４（２）は、嫌気層でのりん酸の濃度を示す。
また、りんに関しては、好気層における分配係数π_１は式（８４）のように示される。

【0229】

【数84】

【0230】

ここで、Ｊ_Ｐ２は、嫌気層の分解りんフラックス（単位：ｍｇ Ｐ／ｍ^２／ｄａｙ）を示す。
また、Δπ_{ＰＯ４，１}は、好気層のりん酸の段階的分画係数（単位：Ｌ／ｋｇ）を示す。Δπ_{ＰＯ４，１}の値は、例えば３００に設定される。
また、π_{ＰＯ４，２}は、好気層のりん酸の段階的分画係数（単位：Ｌ／ｋｇ）を示す。Δπ_{ＰＯ４，１}の値は、例えば１００に設定される。
また、［Ｏ２（０）］_{ｃｒｉｔ，ＰＯ４}は、好気層での分画係数を段階的に減少開始をさせる海水中酸素濃度（単位：ｍｇ／ｋｇ）を示す。［Ｏ２（０）］_{ｃｒｉｔ，ＰＯ４}の値は、例えば２．０に設定される。
また、りんに関しては、ｆ_ｄ１、ｆ_ｐ１、ｆ_ｄ２、ｆ_ｐ２を示す式として、それぞれ式（５１）、式（５２）、式（５３）、式（５４）を用いる。
（シリカ）
シリカの溶出は、式（８５）のように示される。

【0231】

【数85】

【0232】

また、Ｈ_２（ｄＰ_Ｓｉ／ｄｔ）は式（８６）のように示される。

【0233】

【数86】

【0234】

ここで、Ｐ_Ｓｉは懸濁態の生物由来シリカの濃度を示す。Ｓ_Ｓｉは溶存態シリカの生成速度を示す。
また、Ｊ_ＰＳｉは式（８７）のように示される。

【0235】

【数87】

【0236】

また、ｋ_２は、式（８８）のように示される。

【0237】

【数88】

【0238】

また、π_１は、式（８９）のように示される。

【0239】

【数89】

【0240】

また、ｋ_Ｓｉは、藻類由来の粒子態シリカの溶解速度定数（単位：ｄ^−１）を示す。ｋ_Ｓｉの値は、例えば０．５０に設定される。
また、θ_Ｓｉ３は、シリカ酸化の温度係数を示す。θ_Ｓｉ３の値は、例えば１．１０に設定される。
また、［Ｓｉ］_Ｓａｔは、間隙水中のシリカの飽和濃度（単位：ｍｇ Ｓｉ／ｍ^３）を示す。［Ｓｉ］_Ｓａｔの値は例えば４００００に設定される。

【0241】

また、Ｋ_{Ｍ，ＰＳｉ}は、粒子態シリカの分解に対する半飽和定数（単位：ｍｇ Ｓｉ／ｍ^３）を示す。Ｋ_{Ｍ，ＰＳｉ}の値は、例えば５．０Ｅ＋０７（５０００００００）に設定される。あるいは、Ｋ_{Ｍ，ＰＳｉ}の単位としてｍｇ Ｓｉ／ｇを用いるようにしてもよい。この場合、Ｋ_{Ｍ，ＰＳｉ}の値は、例えば１００に設定される。
また、Δπ_Ｓｉ，１は、好気層のシリカの段階的分画係数（単位：Ｌ／ｋｇ）を示す。Δπ_Ｓｉ，１の値は例えば１０に設定される。
また、π_Ｓｉ，２は、嫌気層のシリカの分画係数（単位：Ｌ／ｋｇ）を示す。π_Ｓｉ，２の値は例えば１００に設定される。
また、Ｊ_ＰＳｉは海水中から底泥への粒子態シリカのフラックス（単位：ｍｇ Ｓｉ／ｍ^２−ｄ）を示す。
また、シリカに関しては、ｆ_ｄ２を示す式として式（５３）を用いる。

【0242】

＜水質モデル＞
水質モデルでは、生物の活動の水質への影響を中心に、水中に含まれる物質量（濃度）の変化を演算する。以下では、珪藻類に関連する炭素の循環を例に説明する。
図５は、水質モデルが模擬する珪藻類に関連する物質の循環の例を示す説明図である。
図５に示すように、珪藻が珪藻類に関連する物質量は、珪藻が二枚貝及び動物プランクトンに捕食されることで減少する。また、珪藻の死滅によって懸濁態有機物の量が増加し、珪藻の代謝によって溶存態有機物の量が増加する。また、珪藻の光合成に及び呼吸により、溶存酸素の量が変化する。また、珪藻の呼吸、死滅及び光合成によってりん酸態りん、アンモニア態窒素、及び、シリカの量が変化し、りん酸態りんは、硝酸及び亜硝酸に分解される。
珪藻類に関連する炭素量の変化は式（９０）のように示される。

【0243】

【数90】

【0244】

ここで、Ｐｄｃは珪藻類の炭素量（単位：ｍｇＣ／ｌ）を示す。Ｐｍｄは最適環境下での珪藻類の炭素量の増加速度（単位：ｄａｙ^−１）を示す。Ｐｍｄの値は、例えば観測または試験に基づいて予め設定される。Ｐｍｄの値として、例えば２．２５を用いることができる。
ｆｄ（Ｔ）は珪藻類の水温制限関数を示す。水温制限関数ｆｄ（Ｔ）は、珪藻類の増殖速度に対する水温の影響を模擬する関数である。珪藻類など植物プランクトンの増殖は、環境水温が最適水温に達するまで活性化し、最適水温よりも上昇すると活性が落ちる。この成長（活性）と水温変化の関係は、正規分布曲線に類似する。そこで、水温制限関数ｆｄ（Ｔ）を式（９１）のように設定する。

【0245】

【数91】

【0246】

ここで、Ｔは水温（単位：℃）を示す。Ｔｏｐｔは最適水温（単位：℃）を示す。ＫＴｇ１は最適水温以下の成長にかかわる水温係数（単位：℃^−１）を示す。ＫＴｇ２は、最適水温を超えた場合の成長にかかわる水温係数（単位：℃^−１）を示す。Ｔｏｐｔ，ＫＴｇ１、ＫＴｇ２の値は、例えば観測または試験に基づいて予め設定される。Ｔｏｐｔの値として例えば２０（℃）を用いることができる。ＫＴｇ１の値として例えば、０．００４２（／℃）を用いることができる。ＫＴｇ２の値として例えば、０．００６（／℃）を用いることができる。

【0247】

また、ｆ（Ｉ）は光制限関数を示す。光制限関数ｆ（Ｉ）は、生物の増殖速度に対する光の影響を模擬する関数である。植物プランクトンの増殖にかかわる光制限は、光量が最適光量に達するまで活性化し、最適光量よりも増加すると活性が落ちる。そこで、光制限関数ｆ（Ｉ）を式（９２）のように設定する。

【0248】

【数92】

【0249】

ここで、ａｂは式（９３）のように示される。

【0250】

【数93】

【0251】

また、ａｘは式（９４）のように示される。

【0252】

【数94】

【0253】

ここで、ＦＤは日照率を示す。日照率ＦＤは、例えば気象データに基づいて予め設定される。ＦＤの値として例えば０．４〜０．６の範囲内の定数を用いることができる。
また、Ｋｅｓｓは光の全消散係数（単位：ｍ^−１）を示す。光の全消散係数Ｋｅｓｓは、
植物プランクトンの現存量が増加するにつれて、植物プランクトン自らによる光の吸収（細胞による光の吸収）が大きくなり、深いところでの植物プランクトン群集が利用できる光のエネルギーが減少するという関係を示す。光の全消散係数Ｋｅｓｓは、式（９５）のように示される。

【0254】

【数95】

【0255】

ここで、Ｋｅｂは水中での吸収係数（場の吸光係数）（単位：ｍ^−１）を示す。Ｋｅｂの値として０．７を用いることができる。また、Ｋｅｃｈｌは、植物プランクトン（クロロフィルａ量）による吸収係数（単位：ｍ^２／ｍｇｃｈｌ）を示す。Ｋｅｃｈｌの値として、１７を用いることができる。Ｐｃは、大型藻類の炭素量（単位：ｍｇｃ／ｌ）を示す。Ｐｄｃは、上記のように珪藻類の炭素量（単位：ｍｇＣ／ｌ）を示す。Ｐｃｃは、渦鞭毛藻類の炭素量（単位：ｍｇＣ／ｌ）を示す。ＣＣｈｌは、炭素量とクロロフィルａ量との比（単位：ｇ Ｃ／ｍｇ ｃｈｌ）を示す。ＣＣｈｌの値は、例えば観測または試験に基づいて予め設定される。ＣＣｈｌの値として、例えば、５〜１００の範囲内の値を用いることができる。
また、Δｚは、計算格子の厚さ（すなわち、鉛直方向の格子間距離）を示す。ＺＤは、水面から計算格子の上面までの距離（単位：ｍ）を示す。
また、Ｉｓは、最適光量（単位：ｌｙ／ｄａｙ）を示す。植物プランクトンが光の変化に順応するための時間を考慮して、最適光量Ｉｓを式（９６）のように設定する。

【0256】

【数96】

【0257】

ここで、Ｉｏは、水面での日日射量を示す。Ｉ１は、前日の日日射量（単位：ｌｙ／ｄａｙ）を示す。Ｉ２は、２日前の日日射量（単位：ｌｙ／ｄａｙ）を示す。Ｉｏ、Ｉ１及びＩ２の値は、例えば気象データに基づいて設定される。
Ｄｏｐｔは、植物プランクトンの最大増殖水深（単位：ｍ）を示す。植物プランクトンの最大増殖水深Ｄｏｐｔの値として、１を用いることができる。
また、ｆｄ（Ｎ）は栄養塩制限関数を示す。栄養塩制限関数ｆｄ（Ｎ）は、生物の増殖速度に対する栄養塩の濃度の影響を模擬する関数であり、式（９７）のように示される。

【0258】

【数97】

【0259】

ここで、ＫＨｎは、窒素摂取にかかわる半飽和定数（単位：ｇＮ／ｍ^３）を示す。ＫＨｎの値として、例えば０．０１を用いることができる。ＫＨｐは、りん摂取にかかわる半飽和定数（単位：ｇ Ｐ／ｍ^３）を示す。ＫＨｐの値として、例えば０．００１を用いることができる。ＮＨ_４、ＮＯ_３は、それぞれアンモニア態窒素の濃度、硝酸態窒素の濃度（いずれも単位は、ｇ Ｎ／ｍ３）を示す。ＰＯ_４は、オルトりん酸態りん濃度（単位：ｇ Ｐ／ｍ^３）を示す。

【0260】

また、ＢＭｒｄは、最適水温時の呼吸速度（単位：ｄａｙ^−１）を示す。最適水温時の呼吸速度Ｍｒｄの値として、例えば０．０１を用いることができる。
また、ＫＴｂｍは、呼吸にかかわる水温係数（単位：℃^−１）を示す。呼吸にかかわる水温係数ＫＴｂｍの値として、０．０６９を用いることができる。また、Ｔは水温（単位：℃）を示す。水温Ｔの値は、例えば気象データに基づいて設定される。Ｔｏｐｔは、上述したように最適水温を示す。

【0261】

また、ＢＰｒｄは、最適水温時の死滅速度（単位：℃^−１）を示す。最適水温時の死滅速度ＢＰｒｄの値は、例えば観測または試験に基づいて予め設定される。最適水温時の死滅速度ＢＰｒｄの値として０．０８〜０．０２１５の範囲内の定数を用いることができる。
また、ＫＴｂｐは、死滅にかかわる水温係数（単位：℃^−１）を示す。死滅にかかわる水温係数ＫＴｂｐの値は、例えば観測または試験に基づいて予め設定される。死滅にかかわる水温係数ＫＴｂｐの値として、例えば０．０６９を用いることができる。

【0262】

また、Ｓｖｐｄｃは、植物プランクトンの沈降速度（単位：ｍ／ｄａｙ）を示す。植物プランクトンの沈降速度Ｓｖｐｄｃの値は、例えば観測または試験に基づいて予め設定される。植物プランクトンの沈降速度Ｓｖｐｄｃの値は、例えば０．０２〜０．５の範囲内の定数値に設定することができる。

【0263】

また、ＺＰｍは、動物プランクトンの最適環境下での増殖速度（単位：ｄａｙ^−１）を示す。ＺＰｍは、渦鞭毛藻の現存量Ｐｃｃ、珪藻の現存量Ｐｄｃ、及び、粒子態有機物の現存量ＰＯＣに基づいて、式（９８）のように示される。

【0264】

【数98】

【0265】

ここで、ＺＰｍａｘは、動物プランクトンの増殖速度の最大値（単位：ｄａｙ^−１）を示す。ＺＰｍａｘの値は、実験によって求めることができる。具体的には、動物プランクトンの成体雌を採集して実験室にて翌朝までに産出された約千個の卵を１つの集団として準備する。培養した植物プランクトンを餌として潤沢に与えながら、一定温度で保ち孵化から成体に発育するまでに要した日数を計測し、成長速度を算出する。成長速度は水温上昇によって増大するが、餌供給量が潤沢な実験条件下で得られたことから、この水温における最大成長速度と見なすことができる。
Ｐｒｐｆｃは、植物プランクトンの摂取率を示す。植物プランクトンの摂取率Ｐｒｐｆｃの値は、例えば観測または試験に基づいて予め設定される。植物プランクトンの摂取率Ｐｒｐｆｃの値として、例えば０．５を用いることができる。
Ｐｒｐｆｏｃは、懸濁態有機物の摂取率を示す。懸濁態有機物の摂取率Ｐｒｐｆｏｃの値は、例えば観測または試験に基づいて予め設定される。懸濁態有機物の摂取率Ｐｒｐｆｏｃの値として、例えば０．５を用いることができる。
Ｚｏｏｍｉｎは、動物プランクトンの最少現存量（単位：ｇＣ／ｍ^３）を示す。動物プランクトンがゼロになるとその後の計算で増殖できなくなるため最少現存量を設定している。動物プランクトンの最少現存量Ｚｏｏｍｉｎの値として、０．００１を用いることができる。
Ｚｓｐ２は、植物プランクトン摂取にかかわる半飽和係数を示す。植物プランクトン摂取にかかわる半飽和係数Ｚｓｐ２の値は、例えば観測または試験に基づいて予め設定される。植物プランクトン摂取にかかわる半飽和係数Ｚｓｐ２の値として、０．３を用いることができる。
ｆ（Ｔ）は大型藻類の水温制限関数を示す。ｆｄ（Ｔ）は、上述したように珪藻類の水温制限関数を示す。ｆ（Ｔ）として、式（９１）のｆｄ（Ｔ）をｆ（Ｔ）と読み替えて用いる。
Ｚｃは、動物プランクトンの炭素量（単位：ｍｇＣ／ｌ）を示す。
Ｓｈｃは、二枚貝による摂取量（単位：ｇＣ／ｍ^２／ｄａｙ）を示す。

【0266】

式（９０）の右辺の「Ｐｍｄ・ｆｄ（Ｔ）・ｆ（Ｉ）・ｆｄ（Ｎ）」と「Ｐｄｃ」との積は、珪藻の増殖による炭素量の変化（ここでは、珪藻に関する炭素量の変化）を示す。また、「ＢＭｒｄ・ｅ^{ＫＴｂｍ（Ｔ−Ｔｏｐｔ）}」と「Ｐｄｃ」との積は、珪藻の呼吸による炭素量の変化を示す。「ＢＰｒｄ・ｅ^{ＫＴｂｐ（Ｔ−Ｔｏｐｔ）}」と「Ｐｄｃ」との積は、珪藻の死滅による炭素量の変化を示す。「Ｓｖｐｄｃ・（ｅ／ｄｚ）と「Ｐｄｃ」との積は、珪藻の沈降による炭素量の変化を示す。「Ｚｐｍ・ｆ（Ｔ）・Ｚｃ」は、動物プランクトンが珪藻を捕食することによる炭素量の変化を示す。「Ｓｈｃ」は、上記のように二枚貝が珪藻を摂取することによる炭素量の変化を示す。

【0267】

水質モデルは、珪藻の場合と同様に、渦鞭毛藻、動物プランクトン、大型藻類の各々に関連する物質の循環を模擬する。また、水質モデルは、懸濁態有機炭素の沈降、及び、懸濁態有機炭素から溶存態有機炭素への分解など有機炭素の循環を模擬する。また、水質モデルは、陸源負荷に関する有機物の量の変化を模擬する。また、水質モデルは、炭素の循環の場合と同様に、りんの循環、窒素の循環、溶存酸素の循環、及び、シリカの循環を模擬する。

【0268】

＜透明度の算出＞
透明度の算出式は、式（９９）のように示される。

【0269】

【数99】

【0270】

ここで、式（９９）の右辺の分母「０．１３９×”全ＳＳ濃度”＋０．０１９×”クロロフィルａ濃度”＋０．０４」は、光束消散係数（Ｋ_ｄ）を示す。従って、海域毎の定数は、透明度×光束消散係数（Ｋｄ）となる。透明度を測定し、光束消散係数を算出しておくことで、海域毎の定数を求めることができる。例えば東京湾の場合、海域毎の定数の値を１．５とする。また、伊勢湾の場合及び瀬戸内海の場合、海域毎の定数の値を１．６とする。
式（９９）に示されるように、透明度の算出に際し、全ＳＳ（Suspended Solid、懸濁物質）濃度とクロロフィルａ濃度とを求める必要がある。全ＳＳ濃度は無機ＳＳ濃度と有機ＳＳ濃度との二つに分けて、プログラム演算にて算出する。全ＳＳ濃度を求める式は、式（１００）のように示される。

【0271】

【数100】

【0272】

ここで、上述したようにＰＯＣ、ＰＯＮ、ＰＯＰは、それぞれ粒子態有機炭素濃度、粒子態有機窒素濃度、粒子態有機りん濃度を示す。
無機ＳＳとして、特に流入河川からの無機ＳＳを模擬する。具体的には、河川毎に低水時、出水時それぞれのＬ−Ｑ式を予め求めておく。Ｌ−Ｑ式は、河川の流量Ｑと、浮遊土砂輸送量Ｌとの関係を示す式である。
シミュレーション演算部１９１は、このＬ−Ｑ式に河川の流量を適用して浮遊土砂輸送量Ｌを算出し、浮遊土砂輸送量Ｌに基づいて陸起源ＳＳ（河川から流入する無機ＳＳ）の量を算出する。例えば、河川から流入する浮遊土砂のうち８％は浮遊砂であり河口付近に沈降することから、シミュレーション演算部１９１は、残りの９２％がシミュレーション対象の水域に流入すると算出する。

【0273】

また、無機ＳＳを粒径に応じて１μｍ未満、１μｍ以上〜４μｍ未満、４μｍ以上〜１６μｍ未満、１６μｍ以上の４つに区分する。例えば、１μｍ未満、１μｍ以上〜４μｍ未満、４μｍ以上〜１６μｍ未満の各区分の割合をそれぞれ２％、８％、５０％、４０％に設定しておく。シミュレーション演算部１９１は、Ｌ−Ｑ式に基づいて算出した無機ＳＳ量に割合を乗算して区分毎の無機ＳＳ量を算出する。そして、シミュレーション演算部１９１は、１μｍ未満、１μｍ以上〜４μｍ未満、４μｍ以上〜１６μｍ未満、１６μｍ以上の各区分の沈降速度（単位：ｍ／ｄａｙ）を、それぞれなし（０）、０．２８、４．５３、１８．０７としてシミュレーションを行う。

【0274】

有機ＳＳ濃度はＰＯＣ、ＰＯＮ及びＰＯＰから求めることができる。ＰＯＣ、ＰＯＮ、ＰＯＰは、それぞれ水質モデルにて層毎、かつ、メッシュ毎（格子毎）に求められる。層毎、かつ、メッシュ毎に求めたＰＯＣ、ＰＯＮ及びＰＯＰを集計することで、有機ＳＳ濃度を求めることができる。
また、全ての植物プランクトンがクロロフィルａ色素を含むことから、クロロフィルａ濃度は、海水にある植物プランクトンの量を指す。例えば、植物プランクトンのうち代表的な珪藻及び渦鞭毛藻の量（珪藻の炭素量Ｐｄｃ、及び、渦鞭毛藻の炭素量Ｐｃｃ）を算出し、その合計量をクロロフィルａ濃度として用いることができる。

【0275】

次に、シミュレーションを実行するためのプログラムについて説明する。
＜可変格子の採用＞
図６は、可変格子の設定例を示す説明図である。同図に示すように、詳細に求めたい箇所の格子を細かくし、重要性が薄い所は大きい格子を設定する。これにより、格子を細かく設定した箇所については詳細な演算結果を得ることができる。また、大きい格子を設定した箇所については、演算量を低減させることができる。

【0276】

ここで、格子間の距離（Δx及びΔy）の値を２次元配列の変数に設定する。また、プログラム内で計算に用いる格子間距離を適切に扱うことで、可変格子の計算を可能にしている。
図７は、２次元配列の変数に設定された格子間距離を用いるプログラムの例を示す説明図である。同図は、フラックスを求めるときのループの例を示している。変数Ｈ１はｘ方向の格子間距離（Δｘ）を示す。変数Ｈ２は、ｙ方向の格子間距離（Δｙ）を示す。
但し、２次元配列の変数に設定された格子間距離を用いる箇所は図７に示す箇所に限らず複数個所存在する。

【0277】

＜安定計算の改良＞
コンピュータの計算では丸め誤差が発生する。倍精度浮動小数点数を用いれば、丸め誤差を少なくなり、精度を高めることができる。しかし、すべての項目を倍精度浮動小数点数で計算すると、単精度浮動小数点数で計算する場合と比べて１０倍以上の時間かがかるなど演算時間を要し、また、コンピュータに負荷がかかる。
そこで、水位、垂直渦動粘性係数、及び、垂直渦動拡散係数を倍精度浮動小数点数で表現し、他のデータについては単精度浮動小数点で表現して演算を行う。これにより、精度良く演算を行うことができ、かつ、全てのデータを倍精度浮動小数点で表現する場合よりも演算時間及びコンピュータの負荷の増加を抑制することができる。

【0278】

＜演算の高速化（１）＞
本モデルでは、計算時間の大部分がI,J,Kのループ処理で時間が掛かっている。ループの改良を行うことで、処理速度を改良した。改良したループは生態系を含み100を超えており、同様な処理を行っている。下記に例として、元のプログラムを示し、その後の改良方法を説明する。

【0279】

図８は、高速化処理を行わない場合のプログラムの例を示す説明図である。同図のＩ、Ｊ、Ｋのループは、それぞれｘ方向、ｙ方向、ｚ方向の座標を示す。図８のプログラムでは、Ｉ及びＪの全範囲に対してループを実行している。
図９は、高速化処理を行った場合のプログラムの例を示す説明図である。図９のプログラムでは、Ｉ、Ｊの最大値を予めＩＮＵＭ（ｎ）、ＪＮＵＭ（ｎ）に格納しておく。そして、格納した最大値の範囲内でループを実行する。このように、水域の範囲外の部分に対してループの実行を抑制することで、例えば東京湾のシミュレーションでは演算時間を３分の１に短縮することができた。

【0280】

＜演算の高速化（２）＞
鉛直ループを水平ループの内側に配置して、ＩＦ文の実行回数を減少させた。ここでいう鉛直ループは、格子点を鉛直方向（ｚ軸方向）に順に処理するループである。水平ループは、格子点を水平方向（ｘ軸方向、ｙ軸方向それぞれ）に順に処理するループである。
図１０は、鉛直ループを水平ループの内側に配置したプログラムの例を示す説明図である。同図に示すプログラムでは、変数Ｉ、Ｊのループが水平方向のループの例に該当し、変数Ｋのループが鉛直方向のループの例に該当する。変数Ｉ，Ｊのループは、行Ｌ１１で始まり、行Ｌ１７で終了する変数ｎのループによって実行される。変数Ｉ、Ｊの値は、それぞれ行Ｌ１２、Ｌ１３で設定される。また、変数Ｋのループは、行Ｌ１５で始まり、行Ｌ１６で終了する。
行Ｌ１４のＩＦ文は、変数Ｉ，Ｊの値が決まっている状態で行う必要があり、変数Ｉ、Ｊのループの内側で行う必要がある。一方、行Ｌ１２のＩＦ文には変数Ｋは出現していない。このため、行Ｌ１１のＩＦ文が変数Ｋのループの内側にある必要は無い。

【0281】

水平ループ（変数Ｉ、Ｊのループ）が鉛直ループ（変数Ｋのループ）の内側に位置する場合、ＩＦ文も鉛直ループの内側に位置することになり、鉛直ループを実行する毎にＩＦ文を繰り返し実行する必要がある。
それに対し、図１０の例では、行Ｌ１２のＩＦ文が、行Ｌ１５から始まる鉛直ループの外側に位置している。これにより、図１０の例では、鉛直ループを実行する毎にＩＦ文を繰り返し実行する必要がない。

【0282】

このように、ＩＦ文が鉛直ループの外に位置することで、鉛直ループの並列化を比較的容易に行うことができる。また、ＩＦ文を実行回数が減る点でＣＰＵ負荷を軽減させることができる。
また、図１０の例では、配列の宣言についても鉛直ループの方が内側に位置するようにしている。これにより、ベクトル計算が行えるように最適化している。具体的には、図９の例では、プログラム上部での変数宣言で、ＡＶＤ（ＭＡＸＮＵＭ，Ｋ）としている。これに対し、図１０の例では、プログラム上部での変数宣言で、ＡＶＤ（Ｋ，ＭＡＸＮＵＭ）としている。
図１０のように鉛直ループを水平ループの内側に位置させ、また、配列の宣言についても鉛直ループの方を内側に位置させることで、実験では、鉛直ループが水平ループの外側に位置し、配列の宣言についても鉛直ループの方が外側に位置している場合よりも処理時間を５分の２に短縮することができた。

【0283】

＜演算の高速化（３）＞
ループの並列化により処理時間を短縮させる。並列には、例えばＯｐｅｎＭＰ（Open Multi-Processing）の技術を用いる。
図１１はループを並列化したプログラムの例を示す説明図である。同図の例では、変数Ｉ、Ｊのループを変数ｎのループにて実行している。具体的には、変数ｎのループは、行Ｌ２２で開始し、行Ｌ２５で終了する。このループ内の行Ｌ２３、Ｌ２４で、それぞれ変数Ｉ、Ｊの値を設定している。
また、行Ｌ２１で並列化を宣言しており、変数ｎのループを並列実行する。これにより、実験では、１２コアのＣＰＵを用いて１１倍程度の高速化を行うことができた。

【0284】

以上のように、モデル記憶部１８１は、シミュレーション対象の水域に設定された格子毎に設定され、格子内に含まれる橋脚の投影面積及び体積をパラメータとして橋脚が水流に及ぼす影響を算出する支配方程式を含む水流のモデルを記憶する。また、シミュレーション処理部１９１は、モデル記憶部１８１が記憶しているモデルに基づいて、シミュレーション対象の水域における水流を算出する。
これにより、水域環境シミュレーション装置１００では、水域に橋が設けられている場合でもシミュレーションを精度良く行うことができる。特に、水域環境シミュレーション装置１００では、橋脚の形状にかかわらず投影面積及び体積をパラメータとして用いることができ、従って、橋脚の形状毎に支配方程式をカスタマイズする必要がない。この点で、ユーザの負担を軽減することができる。

【0285】

また、格子間距離記憶部１８２は、シミュレーション対象の水域内の位置毎に設定された格子間距離を記憶する。格子点座標算出部１９２は、格子間距離記憶部１８２が記憶している格子間距離に基づいて、格子点の座標を算出する。
これにより、水域環境シミュレーション装置１００では、シミュレーション対象の水域の形状に応じて柔軟に格子を設定することができる。
なお、実験の結果、隣り合う格子間距離の変化率を±８パーセント（％）以内にすることで、計算を安定させることができた。また、水域環境シミュレーション装置１００またはユーザが、座標を変数（パラメータ）とする２次関数を用いて格子間距離を算出し、格子を設定するようにしてもよい。これにより、スムーズな水流を再現することができる。

【0286】

また、データ記憶部１８３は、水位、鉛直渦動粘性係数、及び、鉛直渦動拡散係数が他の数値よりも桁数の多い数値で表されたデータを記憶する。例えば、データ記憶部１８３は、水位、鉛直渦動粘性係数、及び、鉛直渦動拡散係数を倍精度小数で記憶し、他のデータを単精度小数で記憶する。そして、シミュレーション処理部１９１は、データ記憶部１８３が記憶するデータを用いてシミュレーション対象の水域における水流を算出する。
このように、データの種類に応じて精度を使い分けることで、水域環境シミュレーション装置１００では、必要な記憶容量の増大を抑制し、かつ、高精度にシミュレーションを行うことができる。

【0287】

また、格子数記憶部１８４は、水平方向における格子の行毎及び列毎に格子の数を記憶する。そして、ループ処理部１９３は、格子数記憶部１８４が記憶する格子の数の範囲内で行方向のループ及び列方向のループを実行する。シミュレーション処理部１９１は、ループ処理部１９３が実行するループにて格子毎に当該格子内における水流を算出する。
これにより、水域環境シミュレーション装置１００では、格子数にかかわらずループを一定回数実行する場合よりもループの実行回数を低減させることができる。これにより、水域環境シミュレーション装置１００では、水域環境シミュレーション装置１００の処理負荷を低減させ、かつ、シミュレーションの実行時間を短縮することができる。

【0288】

また、ループ処理部は、水平方向における格子の行及び列の演算を同一のループの並列処理にて実行する。
これにより、水域環境シミュレーション装置１００では、鉛直方向のループを水平方向のループの内側に配置することができ、ループの上限を判定するＩＦ文の実行回数を減少させることができる。これにより、水域環境シミュレーション装置１００では、ループの並列化率を高め、水域環境シミュレーション装置１００の処理負荷を低減させ、かつ、シミュレーションの実行時間を短縮することができる。

【0289】

また、シミュレーション処理部１９１は、シミュレーション対象の水域内に設定された格子毎に水質及び汚濁を算出し、得られた水質及び汚濁に基づいて、格子点毎に透明度を算出する。
このように、水域環境シミュレーション装置１００によれば、シミュレーション対象の水域における透明度を算出し、この透明度に基づいて水域の水質を評価することができる。

【符号の説明】

【0290】

１００ 水域環境シミュレーション装置
１１０ 表示部
１２０ 操作入力部
１３０ データ入出力部
１８０ 記憶部
１８１ モデル記憶部
１８２ 格子間距離記憶部
１８３ データ記憶部
１８４ 格子数記憶部
１９０ 制御部
１９１ シミュレーション処理部
１９２ 格子点座標値算出部
１９３ ループ処理部

【要約】

【課題】水域に橋が設けられている場合でもシミュレーションを精度良く行えるようにする。
【解決手段】水域環境シミュレーション装置が、シミュレーション対象の水域に設定された格子毎に設定され、格子内に含まれる橋脚の投影面積及び体積をパラメータとして橋脚が水流に及ぼす影響を算出する支配方程式を含む水流のモデルを記憶するモデル記憶部と、前記モデルに基づいて前記シミュレーション対象の水域における水流を算出するシミュレーション処理部と、を備える。
【選択図】図１

【図1】

【図2】

【図3】

【図4】

【図5】

【図6】

【図7】

【図8】

【図9】

【図10】

【図11】