寒冷地河川における平面2次元河氷変動計算モデルの開発

北見工業大学 地域未来デザイン工学科

吉

川土木研究所寒地土木研究所

横

山

泰

弘
洋

寒冷地河川における
平面2次元河氷変動計算モデルの開発
Development of 2-dimensional river ice calculation model on ice-covered rivers
吉川 泰弘1・横山洋2Yasuhiro YOSHIKAWA and Hiroshi YOKOYAMA
1
2

北見工業大学地域未来デザイン工学科 准教授
土木研究所寒地土木研究所 主任研究員

要旨

本研究の最終目的は，アイスジャム現象を解明して定式化し，平面2次元河氷変動計算モデルを開発し，iRICソフトウェアのソルバーとして組み込むことである．本研究では，河氷変動現象を，①河氷の形成，②河氷の解氷（融解，破壊），③河氷の流下，④河氷の堆積（アイスジャム），⑤アイスジャムの決壊，⑥河氷の流下，⑦河氷の堆積に区分して研究を進めた．
2020年度は，2019年度に実施した②河氷の解氷の現象に関する現地観測結果の解析結果を踏まえて解氷時期を予測するモデルを構築した．既往ソルバーの
Nays2DHに，2019年度に組み込んだ河氷厚計算式と水温計算式に加えて，河氷
流下堆積式を計算モデルに組み込んだ．本モデルの動作確認を行い，課題を明確にした．
≪キーワード：結氷河川；アイスジャム；解氷；数値計算≫― 21 ―

1. はじめに2018年3月9日の季節外れの大雨と暖気により，北海道内の河川水位は急激に上昇し，各地で被害をもたらした．石狩川水系辺別川では護岸工事中に現場代理人が流され死亡し，沙流川では取水堰において車両3台流出する被害が起きた．浸水被害は，石狩川水系で6箇所，静内川水系で1箇所，網走川水系で2箇所，釧路川水系で2箇所の計11箇所に上る．被害は無いが氾濫危険水位を超過した箇所は，北海道全域で6箇所であった．3月9日に発生した災害は北海道全域に広がった．北海道全域で同時多発的に発生した「2018年3月北海道アイスジャム災害」1）は，前例がなかったため，現場では何が起きたのか把握できなったとの声も聞かれた．今後，同様の気象条件となり，同じ災害を引き起こす可能性がある．
一般に，冬期間の北海道河川は，気温低下により河氷が形成されて結氷する．春先になると，気温上昇により河氷は解氷されて下流へと流下する．流下する河氷が河道内で堆積するとアイスジャムとなり，上流では急激に水位を上昇させて氾濫し，下流では水位を低下させて取水障害を引き起こす．さらにアイスジャムが解消されると，下流へと多量の河氷が流下して被害を拡大させる．これらの現象を定量的に評価する手法として，数値計算モデルが有効である．本研究の目的は，アイスジャム現象を解明して定式化し，平面2次元河氷変動計算モデルを開発し，iRICソフトウェアのソルバーとして組み込むことである．
アイスジャム現象に至るまでの現象を区分すると以下となる．
①河氷の形成②河氷の解氷（融解，破壊）③河氷の流下④河氷の堆積（アイスジャム）⑤アイスジャムの決壊⑥河氷の流下⑦河氷の堆積（アイスジャム）
平成30年度の研究内容は，②の河氷の破壊現象を明らかにするために現地観測を実施した．③④について，著者らが導いた河氷流下堆積式2）の妥当性を現地観測結果との比較により確認した．③④⑤⑥については，これらの現象を再現可能な実験を実施した．なお，①と②の融解については，既往研究により河氷の形成と融解までの現象を再現できる計算式を提案し，現地観測結果の比較から妥当性の確認3）している．
平成31年度（令和元年）の研究内容は，②河氷の解氷の現象を明らかにするために2018年度に実施した現地観測結果を解析した．③～⑥の現象について，2018年度に実施したこれらの現象を再現可能な水理実験の結果を解析した．最終年度のiRICソフトウェアのソルバーとしての組み込みに向けて，既往のソルバーのNays2DHに，河氷厚計算式と水温計算式の計算モデルを組み込んだ．さらに，河氷厚と水温の計算値の妥当性を時系列の観測値から確認した．
令和2年度の研究内容は，昨年度の現地観測の解析結果を踏まえて，気温と降雨から解氷時期を予測するモデルを構築した．既往のソルバーのNays2DHに，河氷流下堆積式の計算モデルを組み込んだ．組
み込んだモデルにおいて，動作確認および水理実験の再現を試みた．本研究の結果，以下の事が明らかとなった．

― 22 ―

2. 解氷時期予測手法に関する検討（河氷の解氷②）
冬期の北海道の河川では，気温が低下すると河道内に河氷が形成され，解氷期になると気温上昇や降雨により河氷が融解・破壊されて流下する．これらの流下した河氷が流速の遅い箇所などで堆積しアイスジャムが発生する．アイスジャムは急激な水位上昇や浸水，死亡事故 1），取水口の閉塞4）などの被害を引き起こす．アイスジャムによる被害を抑えるためには，事前に解氷時期を予測することが望ましい．アイスジャムによる取水障害については， 𝐹𝐹� 値を用いて取水障害の危険度を評価5）して図-2.1

準線形貯留型モデルの概念図

いる．アイスジャムによる浸水については，氷板厚による検討6）が実施されているが，この検討に用いられている氷板厚の計算式では，解氷時期を十分には予測できていない．その要因として，気温上昇による融雪や降雨による水深上昇を考慮していないことが指摘されている6）．
本研究は，入手しやすい気象データを用いて，簡易的に解氷時期を予測する手法を開発することを目的とし，気温上昇により融雪が促され，融雪水が河川へと流入して水深が上昇する現象を再現できる解
氷時期予測モデルを構築し，構築したモデルの検証を実施した．
（1）解氷時期予測モデルの構築
解氷時期予測モデルの構築にあたって，気温上昇から融雪量，融雪量と雨量から流量，流量から水深の3つの現象に分けてモデルを構築した．融雪量と雨量による水深上昇は，準線形貯留型モデル7）を用いた．概念図を図-2.1に示し，計算モデルを次に示す．
ௗ௦ௗ௧ൌ ‫ݎ‬௘ Ԣ െ ‫ݍ‬ܵ ൌ ݇Ԣ‫ݍ‬௧݇Ԣ ൌ ೎ '
ଶ

（1）
（2）
（3）

𝑟𝑟� ��mm/h�：有効雨量と有効融雪量の和，‫ݍ‬ሾmm/hሿ：流出高， ܵሾmmሿ：タンクの貯留量， ݇Ԣሾhሿ：貯留係数，𝑡𝑡� ���：洪水到達時間であり，次式に示す角屋・福島の式8）を用いた．

𝑡𝑡� � �𝐴𝐴�.�� 𝑟𝑟� ′��.�� /60（4）

𝐴𝐴�km� � ：流域面積， ‫ܥ‬ሾ無次元ሿ：土地利用形態で決まる定数8）であり，自然丘陵山地として290を与えた． 𝑡𝑡� の単位はminからhourにするために単位変換を行った．また， 𝑡𝑡� の値は，過去のアイスジャムによる氾濫が発生した際の 𝑟𝑟� ′ を用いて，洪水到達時間を算出し，その値を平均した値を与えた．
a）気温上昇から融雪量

気温上昇から融雪量を求めるために，Degree-day法9）を用いた．本研究では1時間毎の融雪量を求めるため，次式を用いた．

ܵ‫ ܯ‬ൌ ܶ ή ‫ܭ‬（5）ܵ‫ܯ‬ሾmmhሿ：1時間毎の融雪量， ܶሾԨሿ：1時間毎の気温， ‫ܭ‬ሾmm/h ή Ԩሿ ：融雪係数であり1.7÷24を与えた10）．融雪時期は，気温が0℃以上の時に融雪するものとした9）．

― 23 ―

b）融雪量と雨量から流量流量を求めるために，まず，次式を用いて有効融雪量と有効雨量を求めた．
𝑆𝑆𝑆𝑆� � 𝑆𝑆𝑆𝑆 � 𝑓𝑓

（6）

𝑟𝑟� � 𝑟𝑟 � 𝑓𝑓

（7）

𝑟𝑟�� � ��� � 𝑟𝑟�

（8）𝑆𝑆𝑆𝑆� �mm/h� ：有効融雪量，ܵ‫ܯ‬ሾmm/hሿ ：式（5）より求めた1時間当たりの融雪量，𝑟𝑟� ���/�� ：有効雨量，‫ݎ‬ሾ݉݉/݄ሿ ：観測雨量， 𝑓𝑓 [無次元]：流出係数であり，本研究では，各河川の流出係数が不明だったため山地河川とし，0.75を与えた11）．

次に準線形貯留型モデルにより流出高を求める．既往研究12）を参考にして，式（1）に式（2）を代入

して1次のテイラー展開で表すと次式となる．
��𝑞𝑞� � �������� ����� ����� � ����
���

��
�

（9）

𝑞𝑞� �mm/h� ：流出高，𝑞𝑞� �mm/h� ：1時間前の流出高，

𝑟𝑟�� ��mm/h� ：有効雨量と有効融雪量の和，𝑟𝑟�� ��mm/h� ：1時間前の有効雨量と有効融雪量の和，∆t�h� ：時間であり，本研究では1を与えた．この式（9）に，式（3）と式（8）で得られた値をそれぞれ代入し，流出高を求める．

次に流量を次式により求めた．
஺ή௤

ܳൌ

（10）

ଷ.଺

ܳሾmଷ /sሿ ：流量， 𝐴𝐴�km� � ：流域面積， ‫ݍ‬ሾmm/hሿ ：式（9）より求めた流出高を与えた．

c）流量から水深

水深は，次式のマニングの平均流速公式を用いて，流量から求めた．求める水深は，等流水深となる．
��� � � � 𝐴𝐴 � 𝑅𝑅 � 𝐵𝐵𝐵�𝑖𝑖� �
�

݄ൌ൬

ொ௡

ඥ௜బ ஻

���

య
ఱ

�

൰

�

𝐵𝐵𝐵�

（11）
（12）

ℎ�m� ：水深，ܳሾmଷ /sሿ ：式（10）より求めた流量，��s/𝑚𝑚�/� � ：粗度係数であり，本研究では0.03を与えた， 𝑖𝑖� �無次元� ：河床勾配， ‫ܤ‬ሾmሿ ：川幅である．

以上の手法により，気温上昇により融雪が促され，融雪水および雨水が河川へと流入して，水深を上昇させる現象を計算式により再現した．なお，本検討で用いた各河川の定数を表-2.1に示す．
（2）本計算式の取水障害現象への適応
取水障害が発生する可能性の高い晶氷詰まりの発生過程について，結氷後の発生過程を考えると，

表-2.1

解氷時期予測モデルに用いた各河川の定数
水面幅 流域面積 洪水到達
アメダス
2
（m)
（km ) 時間（h)

河川名

縦断勾配

の増加，次に融解・破壊された河氷の流下と下流

名寄川

1/480

67.96

743

48.05

名寄

側への河氷破壊の進行，さらに流下した河氷の流

網走川

1/232

24

1380

43.45

津別

速低減箇所への堆積（晶氷詰まり）に整理される．

浦幌川

1/1106

92

483

33.9

浦幌

1/697

24

195.5

26.5

旭川

初めに気温上昇による融雪と降雨による河川流量

これらを踏まえ，晶氷詰まり発生の直前では， オサラッペ川― 24 ―

まとまった降雨と気温上昇，それに伴う流量増加（水深上昇）の現象が現れるため，事前に，水深上昇を予測できれば，晶氷詰まりの発生を事前に予測できると考えられる．本モデルの水深変動と晶氷詰まりの現象の比較を行った．
a）調査および現地の概要
取水障害については，過去4年間の名寄川の観測データを用いて構築した計算モデルの妥当性を確認した．観測地点は，天塩川水系名寄川の上流7.05km地点に位置する真勲別頭首工周辺とした．真勲別頭
首工では，2005年1月，2009年2月，2013年1月にアイスジャムによる取水障害が発生している．観測期間は2016年12月13日から2017年3月21日，2017年12月8日から2018年3月9日，2018年12月11日から2019年
3月8日，2019年12月6日から2020年3月12日であり，計算期間は観測期間と同様にした．観測機器は水位計，水温計，ADCP，SWIPを設置している．
なお，ADCP（Acoustic Doppler Current Profiler 1200kHz, Teledyne RD Instruments社製）は，河氷速度と河氷底面の測定 13） のために1台設置し10分毎に測定した．SWIP（Shallow Water Ice Profiler, ASL
Environmental Sciences社製）は，SWIP（周波数546kHz）とSWIP（周波数235kHz）の2台を河氷の存在を確認13）するために設置し10分毎に測定した．
b）晶氷詰まり状況の概要と計算水深
図-2.2に2016年12月13日から2017年3月21日のADCP，SWIPの観測地の経時変化，同時に構築した計算モデルで得られた観測期間中の計算水深を表した．以下に観測した4年間の晶氷詰まり状況の概要と計算水深の変動について記す．
2016年12月13日から2017年3月21日
現地観測では，図-2.2に示すとおり，河氷移動速度が大きく晶氷が詰まり取水障害の危険性が高くなった日時は，2016年12月22日の18時以降，2016年12月28日の15時30分以降，2017年1月13日以降，2017
年1月25日以降，2017年2月17日以降であった．また，図-3より名寄川KP8.4地点のSWIPにて，2016年12月28日の15時30分以降，2017年1月13日以降，2017年2月17日以降で晶氷詰まりによる水位変化が捉えられている．
本モデルの計算結果は，2016年12月21日の18時，2016年12月26日の16時で水深上昇のピークとなっている．このピークからそれぞれ24時間，47時間30分後に晶氷詰まりが発生している．2017年1月25日以降，2017年2月17日以降では，水深上昇のピークよりそれぞれ10時間，9時間早く晶氷詰まりが発生している．
いずれも水深上昇が予測されている最中に発生している．2017年1月13日以降では，晶氷詰まりが発生したが，構築した計算モデルでは水深上昇を予測することができなかった．原因としてはこれらの期間
中に気温が0℃を上回ることがなかったため，気温上昇による融雪や降雨以外の要因が考えられる．

― 25 ―

計算水深(m)

12/28

1/2
1/7
1/12
1/17

1/22

1/27

2/1

2/6

2/11

2/16

2/21

2/26

3/3

3/8

3/13

3/18

晶氷詰まり発生時期

標高(m)― 26 ―

標高(m)

標高

(m)

101.6

102.1

102.6

103.1

103.6

99.0

99.5

100.0

100.5

101.0

101.5

図-2.2

2016年12月13日から2017年3月21日の計算水深，ADCP，SWIP

0

0.4
0.2

99.5

99.0

0.6

100.5

100.0

1

12/23

晶氷詰まり発生（ 2 /17）

計算水深

0.8

12/18

晶氷詰まり発生（12/28 15時30分）

晶氷詰まり発生（ 1 /25）

101.5

12/13

晶氷詰まり発生（ 1 2 / 22 1 8 時 ）晶氷詰まり発生（ 1 /13）

101.0

0

0.2

0.4

0.6

0.8

1

河氷移動速度(m/sec)2017年12月8日から2018年3月9日
構築した計算モデルでは0℃を上回る時に若干の水深上昇が予測され，観測水位も上昇していた．
2017年度は，晶氷詰まりは発生しなかった．
2018年12月11日から2019年3月8日
晶氷詰まりは2019年1月14日から1月15日，2019年2月13日から2月14日に発生している．これはダムの試験放流により流量が増加することによって発生したものと考えられる．構築した計算モデルでは，この晶氷詰まりを計算水深で予測することはできなかった．
2019年12月6日から2020年3月12日
晶氷詰まりは2020年1月11日から1月28日と2020年2月25日から3月3日に発生している．2020年1月11日
から1月28日では，気温が0℃を上回ることがなかったため気温上昇による融雪や降雨による水深上昇は計算されず再現できなかった．気温上昇による融雪や降雨以外の要因が考えられる．2020年2月25日から3月3日では，構築した計算モデルで水深上昇のピークが現れている．この42時間後に晶氷詰まりが発生した．
（3）本計算式のアイスジャム現象への適応
本計算式による水深変動と現地におけるアイスジャム現象の比較を行った．対象としたアイスジャム
現象について記す．2018年3月8日から9日にかけ，温帯低気圧が発達しながら北海道を通過し，道内に季節外れの大雨と暖気をもたらした結果，北海道の多くの河川で比較的規模の大きい出水がみられ，その出水により結氷河川でのアイスジャムが発生し石狩川
水系辺別川では護岸工事中に現場代理人が流されて死亡する事故が起こる14）など，その影響は北海道全域で生じた．2020年3月10日から11日にかけて気温上昇と降雨が重なり，2018年よりは規模は小さいものの道
内の複数河川でアイスジャム発生が確認された．これ
らはいずれも降水量が多く，かつ気温が上昇したことが影響したと推察されている．
本計算式のアイスジャム現象の適応については，
2018年3月に発生した網走川，浦幌川の事例1），2020年（a）網走川津別観測所（津別アメダス）

3月に発生した網走川，浦幌川，オサラッペ川の事
例14）を対象として，構築した計算モデルの妥当性を確認した．
a）アイスジャム発生時の観測水位と計算水深
図-2.3は2018年3月の2河川，図-2.4は2020年3月の3河川の観測水位と計算水深，近傍観測所の毎時の気温と降水量を示し，以下に各河川の発生状況の概要を記
す．なお計算期間は3月1日から3月16日までとした．
網走川（図-2.3（a），図-2.4（a））
アイスジャム発生時刻は，2018年3月9日午後で時刻は不明であり，2020年3月11日9時である．2018年は発生時刻が不明のため計算水深との比較が困難であっ― 27 ―（b）浦幌川 KP0.2（浦幌アメダス）
図-2.3

 河川のアイスジャム発生箇所の水位
各変化量予測値（2018年3月）

た．2020年は計算水深のピークはアイスジャム発生時
よりも3時間ほど早い．
浦幌川（図-2.3（b），図-2.4（b））
アイスジャム発生時刻は，2018年3月9日午後で時刻は不明であり，2020年3月11日9時である．2018年は発生時刻が不明のため計算水深との比較が困難であった．2020年は計算水深のピークはアイスジャム発生時
よりも9時間ほど早い．
オサラッペ川（図-2.4（c））
アイスジャム発生時刻は，2020年3月11日14時である．計算水深のピークはアイスジャム発生時と一致し（a）網走川津別観測所（津別アメダス）ている．
図-2.3，2.4より，構築した計算モデルによる計算
水深は，アイスジャム発生時よりも前にピーク値を示している．計算水深と観測水位のピーク値を比較すると，2018年の網走川では4時間，浦幌川では5時間ほど
計算水深が観測水位よりも先行している．2020年でも網走川では10時間，浦幌川では7時間，オサラッペ川
では1時間と計算水深が観測水位よりも先行している．
一方で，2018年3月1日22時に網走川，浦幌川の水深上昇が予測されたが，実河川では水深は上昇していな
い．今回構築した簡易的な計算モデルでは，ある期間（b）浦幌川 KP0.2（浦幌アメダス）

では実河川の水深を予測できず課題である．しかし，本計算モデルは，既往の手法4）と比べて，アイスジャムの予測精度は向上している．アイスジャム発生の危険性を予測する手法として，本手法と既往の手法6）と併用することで，より広い現象を踏まえた予測手法となる．

(c) オサラッペ川鷹栖橋（旭川アメダス）
図-2.4

-7-

― 28 ―

 河川のアイスジャム発生箇所の水位
各変化量予測値（2020年3月）

3. 平面2次元河氷変動計算モデルの開発
平面2次元河氷変動計算モデルの開発においては，河川の流れ・河床変動解析ソフトウェアiRIC
15）（International River Interface Cooperative）のうち,非定常平面2次元河床変動解析モデル 16）である「Nays2DH」を基本とし，河氷の流れ・河川水温移流拡散・氷板の形成融解現象を再現できる計算モデ
ルを導入して，平面2次元の河氷変動モデルを構築した．
a）流れの基礎方程式への河氷の導入「Nays2DH」の流れの基礎方程式に，河氷の影響17）を考慮した基礎式を以下に示す．

 H ih   hu    hv 
i i



t
x
y
 w t  uh t  vh t  hu 2 x  huv x  huv y  hv 2 
y

（13）


 gh

F
H z  x  Dx  x


x

（14）


 gh

Fy
H z  y  Dy 


y

（15）

ここで，Hi（m）は移動する河氷（氷盤，晶氷）の鉛直方向の厚さ，λiは河氷の空隙率（0.4），ρwは水の密度で1000（kg/m3），ρiは氷の密度で917（kg/m3），Hzは「河床高+水のみが流れる水深+氷板と氷板
下の河氷の喫水深」として，次式により算出した．ここで，移動しない河氷（氷板）の厚さはHis（m）とした．なお，計算モデルでは，これらの式は一般座標系で示した式で計算を行っている．

Hz  Z  h i 1  i  Hi  His 
w

（16）

b）河氷の流れ
河氷の流れの基礎式は，連続の式と河氷流下堆積式として次式17）を用いた．式は，河氷（氷盤，晶氷）
を氷の集合体の氷塊と考え，氷塊に作用する力のモーメントをゼロと仮定して導いた式である．

1  i ui u vi  v H i   H i ui    H i vi 



0
t
x
y

Bi   w  i  gH iBd   w 2

 H i  CL
CD  H i   Cf 
2  Li  Li  2
Bi   w  i  gH iBd   w 2

 H i  CL
CD  H i   Cf 
2  Li  Li  2

（17）

（18）

（19）

H（m）は氷塊の鉛直方向の厚さ，u（m）は氷塊の幅，Bii，v（m/s）はx，y方向の河氷の速さ，Biid（m）
は下流の川幅，g（m/s2 ）は重力加速度で9.8，L（m）は氷塊の河川縦断方向の長さ，Ci
Dは形状抵抗係数，
Cf は摩擦抵抗係数，CLは揚力係数である．

― 29 ―

c）氷板の形成融解
氷板の形成融解の計算は，熱収支式から導出した次式18）を用いた．
 65.2  T
 45.8 
His 
His'   5  a'   2   4 5Tw H 1 3
 10  His  10 

（20）

I B

 0.906  2.770 b
H

（21）



Uw

H2 3

（22）

His（m）：氷板厚，H’is（m）は∆t前の氷板厚，H（m）：水のみが流れる水深である．Ta（℃）：気温（1日の平均値），Tw（℃）：水温（1日の平均値）であり，例えば∆t = 1時間の場合は 24 で割り単位変換した値 Ta /24，Tw /24 が入力値となる．H’i（m）の初期条件は，既往研究18）で妥当性が確認されている
1mm を与えた．α（無次元）は気温に対する氷板形成の程度を表し，大きくなると氷板を増加させる．
β（m1/3/s）は水温と有効水深に対する氷板融解の程度を表し，大きくなると氷板を融解させる係数である．また，係数αは，積雪または晶氷が氷板と一体となる場合は大きく，積雪または晶氷の断熱効果がある場合は小さくなる．係数βは，動水勾配が大きく粗度が小さい場合は大きく，動水勾配が小さく粗度が大きい場合は小さくなる．
係数αは，北海道における河川結氷時の計249回の観測データ解析から得られた式（21）を用いて算出
し，係数βは，式（22）から計算されるUwおよびHを用いて値を得た．なお，流下堆積する河氷が存在する場合には，上式の右辺第三項の融解項を用いて融解をさせた．
d）河川水温移流拡散
一般曲線座標系における河川水温計算式は，横断結氷比Nを用いて大気と河川水，氷板と河川水との熱収支変化を考慮して以下の式で表される．
 1  N  wa  Nw
  wC pTw h    u wC pTw h    u wC pTw h 




Dw 






t 
J
J
J
J
    Dw

 C T h C T h  
2 w p w 
2 w p w 
 E x

 Ex

 


  


（23）

（24）


wa hwa Tw ' Ta w Cwi

uw0.8
hw0.2

Tw ' T f （25）

ここで，ρwは水の密度，Cpは水の比熱（4.2），Twは河川水温，hは水深，uは  および 方向の流速，
hwaは水面の熱交換係数（20），Taは大気温，Cwiは係数（1622），uwは流速，hwは水深，Tf は氷板底面の温度（0），Exは拡散係数0.001である．

― 30 ―

e）横断結氷比
横断結氷比Nは，結氷判定のための最少河氷厚Himax（m）を用いて，次式で算出した．なお，計算上の安定のため，河氷の速度について，Bi /Bd = Nとして算出した．
N

Hi  His
H i max（26）

4. 本計算モデルの動作確認と再現計算
a）河氷の挙動の動作確認（氷板の形成融解）
河氷の挙動の動作確認として，氷板の形成融解現象を確認した．モデル河道の初期条件，境界条件は以下とした．
・初期条件は，水温は0.2℃，外気温は-200℃とした．短時間で結氷が進むように設定した．
・流動しない固定された氷板を初期条件として与えた．流れの阻害による水位上昇を確認する．
なお，氷板は初期水面の上に設定する．このため初期の流動が収まってから流量を50m3/sに増加させる．その後，流入水の水温を100℃に上昇させ，氷板の融解と水位低下を確認する．
図-4.1に初期条件を示し．図-4.2に経時変化を示す．図-4.2のa）の図において，初期の氷板以外にも周辺に氷が形成されているのが分かる．b）の図では，流量が50m3/sとなっている状況であるが，氷板の上流部で水位が上昇しているのが分かる．c）の図では，100℃の温水が流入することにより氷が溶けている．d）の図では，さらに温水が下流へと流入し，固定された氷板を融解しているのが分かる．

図-4.1

河氷の挙動の動作確認（氷板の形成融解）の初期条件― 31 ―

a) 800 秒経過b) 2000 秒経過c) 4000 秒経過d) 5200 秒経過

図-4.2

河氷の挙動の動作確認（氷板の形成融解）の経時変化

b）河氷の挙動の動作確認（河氷の流下堆積）
河氷の挙動の動作確認として，河氷の流下堆積現象を確認した．モデル河道の初期条件，境界条件は以下とした．
・初期条件は，水温は0.2℃，外気温は-200℃とした．短時間で結氷が進むように設定した．
・流動しない固定された氷板を初期条件として与えた．流れの阻害による水位上昇を確認する．
・流動する河氷は，氷板より上流において初期条件として与えた．河氷の流下状況を確認する．
なお，氷板は初期水面の上に設定する．このため初期の流動が収まってから流量を50m3/sに増加させる．その後，流入水の水温を100℃に上昇させ，氷板の融解と水位低下を確認する．計算上の最小河氷
厚は，0.5mmとした．この値を変えることで，河氷流下時の広がり方が変化し，値が大きいと広範囲に広がりやすくなる．
図-4.3に初期条件を示し．図-4.4，図-4.5に経時変化を示す．図-4.4のa）の図では，流動する河氷は形を維持できずに広がっているのが分かる．b）の図では，流量の増加に伴い，河氷が下流へと流下してい,b）る．図-4.5のa）の図では，流量の増加に伴い，さらに河氷が下流へと流下している．c）の図では，100℃の温水が流入して河氷が融解している．

― 32 ―

図-4.3

河氷の挙動の動作確認（河氷の流下堆積）の初期条件a) 400 秒経過
図-4.4b) 1200 秒経過

河氷の挙動の動作確認（河氷の流下堆積）の経時変化c）アイスジャム水理実験の再現
既往の実氷によるアイスジャム実験19）を参考にして，河氷（氷盤）の移動，障害物による阻害等の影響について再現性確認を行った．計算条件を以下に記す．
計算区間は延長13m，河道形状は複断面形状（低水路幅0.3m，高水敷幅0.25m×2），平面形状は流下
方向：521測線（約0.025m間隔）で横断方向：33測線．粗度係数は0.010（木製水路）を基本として実験による5m地点および7m地点の実測水深0.009m程度となるように調整した．流量は0.7L/s（Q=0.0007m3/s），
起算水位は河床勾配（i＝1/120）として等流水位を与えた．橋脚は2本とし計算上は計算安定のため植― 33 ―

a) 1600 秒経過b) 2000 秒経過c) 4000 秒経過d) 5200 秒経過

図-4.5

河氷の挙動の動作確認（河氷の流下堆積）の経時変化

生の抵抗として設定した．（実験による6m地点の実測水深0.013m程度となるように調整）．河氷（氷盤）
はB40mm×L40mm×t6mm×500枚の実験に対してB0.2m×L4.0m×t0.006＝0.0048m3として設定．気温：
10℃，水温：1.6℃（1.3℃～1.9℃の実験値を参考）．上流端境界は氷板流入なし（上流端から流入する氷
板面積・体積なし）
本計算モデルの計算結果について，水深の時系列変化を図-4.6に，水位縦断図，河氷（氷盤）厚のコ
ンター図，水深コンター図を図-4.7，図-4.8に示す．図-4.7は，河氷なしの計算結果であり，実験結果の水深と同様であることを確認している．
図-4.6と図-4.8から，流下した河氷（氷盤）は橋脚位置（6m地点）で一時的に閉塞傾向となり，水位が上昇し，その後，徐々に閉塞が解消する傾向は概ね表現できている．一方で，以下の点において，実験結果と異なっており，モデル改良が必要である．
・
河氷を氷盤の集合として流体のように取り扱っているため，橋脚上流側に縦断的には堆積するが，横断的な閉塞が少ない．
・水際や流速がない箇所で河氷が動かなくなっている．
・一時的な閉塞により，水位上昇しているものの，高水敷には冠水していない．
・水深の時系列変化において，閉塞時の水位上昇が実験よりも小さい．

― 34 ―

水深(m)

0.05
0.04

5m地点（実験結果）

6m地点（実験結果）

7m地点（実験結果）

5m地点（計算結果）

6m地点（計算結果）

7m地点（計算結果）

0.03
0.02
0.01
0.00
0

10

20

30

40

50

60

70

80

90

時間(s)
図-4.6

水深の時系列変化（実験と計算の比較）■河氷なしでの水深確認・5m 地点：0.009m
・6m 地点：0.013m

図-4.7

計算水位の確認（河氷なし）― 35 ―

100

110

120

a) 0 秒経過b) 5 秒経過c) 15 秒経過d) 25 秒経過e) 35 秒経過f) 120 秒経過
図-4.8

アイスジャム実験の再現計算― 36 ―

d）本計算モデルの課題
アイスジャム実験の再現計算における課題
アイスジャム実験の再現計算結果から，流下した河氷（氷盤）は，橋脚位置で一時的に閉塞傾向となり水位は上昇し，その後，徐々に閉塞が解消する状況を概ね再現できている．一方で，以下の点において，実験結果と異なっており，モデル改良が必要である．
・
河氷を氷盤の集合として流体のように取り扱っているため，橋脚上流側には縦断的に堆積するが，横断的な閉塞が少ない．
・水際や流速がない箇所で河氷が動かなくなっている．
・
一時的な閉塞により，水位上昇しているものの，高水敷には冠水していない．
基礎式の2次元化に関する課題
・
川幅Bに関するパラメータ，河氷破壊に関する判定など，平面2次元的な影響を考慮して基礎式を改良する必要がある．
・
河氷の挙動について，一般座標化ができていない．矩形断面だけでなく，湾曲などを含む自然河川において検証する必要がある．
境界条件に関する課題
・
水際や最少水深以下の箇所における水温・河氷の取扱について調整する必要がある．流速0の地点では，河氷が動かなくなるが，定常流状態では一度動かなくなると，再移動するきっかけを失い，全く動かなくなる．
・
河氷の上流端条件について，河氷量（m3/s）として入力するためには，計算を不安定化させないように，緩和係数などによる安定化が必要である．河氷を考慮した等流状態の計算手法等が必要．

5. まとめ
昨年度の現地観測の解析結果を踏まえて，気温と降雨から解氷時期を予測するモデルを構築した．既
往のソルバーのNays2DHに，河氷流下堆積式の計算モデルを組み込んだ．組み込んだモデルにおいて，動作確認および水理実験の再現を試みた．本研究により，以下の結果が得られた．
●
入手しやすい気象データを入力値とする簡易的な解氷時期予測モデルを構築した．過去に発生した取水障害およびアイスジャムの事例を用いて，本計算モデルから得られた計算結果と比較し，計算値の妥当性を確認した．本研究の解氷時期予測モデルは，既往の手法4）と比べて，取水障害およびアイスジャムの予測精度は向上した．
●
平面2次元河氷変動計算モデルの開発においては，河川の流れ・河床変動解析ソフトウェアiRICのうち「Nays2DH」を基本として，河氷の流れ・河川水温移流拡散・氷板の形成融解現象を再現でき
る計算モデルを導入して，平面2次元の河氷変動モデルを構築した．本計算モデルの動作確認を行った．また，アイスジャム実験の再現を試みた．河氷（氷盤）は橋脚位置（6m地点）で一時的に閉塞傾向となり，水位が上昇し，その後，徐々に閉塞が解消する傾向は概ね表現できた．一方で，実験結果と異なる結果が出ている箇所もあり，今後，モデル改良が必要である．
謝辞
本研究は，（株）開発工営社の濱木道大氏よりiRICソフトウェアへの取り込み作業のご協力，北見工業大学（2021年3月卒業）の寺崎海斗氏より解氷時期に関する解析作業等の取り纏めのご協力を頂いた．
記して謝意を表します．

― 37 ―

参考文献
1）横山洋，吉川泰弘，伊波友生，前田俊一，矢部浩規：2018年3月北海道アイスジャム洪水による被害の現象解明と今後の課題，土木学会論文集B1（水工学），Vol.74，No.5，pp.I_1039-I_1044，
2018．
2）吉川泰弘，朴昊澤，大島和裕，横山洋：寒冷地河川における簡易的なアイスジャム計算モデルに関する検討，土木学会論文集B1（水工学），Vol.74，No.4，pp.I_351-I_356，2018．
3）吉川泰弘，渡邊康玄，早川博，平井康幸：結氷河川における解氷現象と実用的な氷板厚計算式の開発，土木学会論文集B1（水工学），Vol.68，No.1，pp.21-34，2012．
4）入交泰文，吉川泰弘，黒田保孝：天塩川水系名寄川真勲別頭首工における晶氷の発生について，国土交通省北海道開発局技術研究発表会,，第58回，AA-6, 安21（治），2015．
5）諸田匡紀，吉川泰弘，横山洋，森田共胤，村山文保，山本新吾：冬期取水障害における晶氷変動量, 75, 2, pp.I_931-I_936，2019．
の勾配を用いた危険度評価手法の検討，土木学会論文集B1（水工学）6）横山洋，吉川泰弘，伊波友生，矢部浩紀：結氷河川における解氷期の河氷挙動と解氷の簡易予測手法の検討，土木学会論文集B1（水工学）Vol. 75，No2，p. I_931-I_936，2019．
7）橋本健，森田正：土地利用変化を評価する洪水流出モデルに関する研究，土木学会論文報告書，第325号，pp45-50，1982．
8）角屋睦, 福島晟：中小河川の洪水到達時間, 京大防災 研究所年報, 第19号B, 昭和51年4月, pp. 143152．1976．
9）河島克久，飯倉茂弘，杉山友康，遠藤徹，藤井俊茂：鉄道防災に適用可能な日融雪量の簡易推定方法，日本雪氷学会誌，64巻6号，pp.605 – 615，2002．
10）宮田俊介，朝岡良浩，風間聡：全国のAMeDAS観測点におけるDegree-day法および融雪係数の検証，土木学会論文集B1（水工学）Vol.68，No4，I_348，2012．
11） 国土交通省河川砂防技術基準調査編：国土交通省，水管理・国土保全局，第3章，第2節，p.11，平成26年4月．
12）岡本芳美：流出計算マルチ・タンク・モデルに基づく分割法について（Ⅰ）,水利科学，55巻3号，pp.1-21，2011．
13）橋場雅弘，吉川泰弘：天塩川における河川解氷時の河氷の挙動に関する現地観測，土木学会，河川技術論文集，第17巻，pp.365-370，2011．
14）横山洋，吉川泰弘，伊波友生，矢部浩規：2020年3月に発生したアイスジャムの現地調査および発生リスクの評価，土木学会論文集B1（水工学）Vol.76, No.2, I_157-I_162, 2020．
15）http://i-ric.org/index.html
16）清水康行：河道平面形状の形成における河床・河岸の変動特性の相互関係について，
水工学論文集，第47巻，PP643-648，2003．
17）吉川泰弘ほか：寒冷地河川における簡易的なアイスジャム計算モデルに関する検討，土木学会論文
集B1（水工学），Vol.74，No.4，I_1351-1356，2018．
18）吉川泰弘ほか：結氷河川における解氷現象と実用的な氷板厚計算式の開発，土木学会論文集B1（水工学），Vol.68，No.1，pp.21-34，2012．
19）小池太郎，吉川泰弘，横山洋：橋脚部における実氷を用いたアイスジャム氾濫実験，土木学会論文
集B1（水工学）Vol.75，No.2，I_1411-I_1416，2019．

― 38 ―