J. Jpn. Soc. Soil Phys.
土壌の物理性No. 146, p.3 ∼ 16 （2020）

群落熱収支モデルと HYDRUS-1D の連結プログラムを用いた畑地の土壌水分量 · 地温予測モデルの構築— ダイズ栽培圃場における蒸発 · 蒸散の配分割合の検討 —坂井 勝 1 · 高橋由奈 1 · 丸山篤志 2

Simulation of soil water & heat transport in a vegetated ﬁeld using a coupled programof the canopy energy budget model and HYDRUS-1D— Investigation of partitioning evapotranspiration in a soybean ﬁeld —Masaru SAKAI1 , Yuna TAKAHASHI1 and Atsushi MARUYAMA2

Abstract: Evapotranspiration plays an important role toinvestigate simultaneous movement of water and heat in avegetated ﬁeld. To simulate changes of soil water content
and soil temperature using the Richards equation with root
water uptake, partitioning potential evapotranspiration, Ep ,into potential soil surface evaporation, Egp , and crop transpiration, Ecp , is necessary. The double source model, oneof the canopy energy budget models calculating energybalances at both canopy and soil surface, is a promisingmodel to calculate evaporation rate, Eg , and transpiration
rate Ec , separately. In this study, we developed a numerical model to simulate water movement and heat transport
in the root zone by implementing the double source modelin HYDRUS-1D. Field observation was also carried out ina soybean ﬁeld to test the developed model. Simulated soil
water contents and temperatures using various values ofradiation transmissivity, τ , were compared with observed
data during relatively wet condition. It was indicated that
determining τ based on observed surface coverage, Sc (i.e.,τ = 1 - Sc ) can be a reasonable method to partition Ep to
Egp and Ecp with considering crop growth. Simulation for alonger drying period was also carried out to show changesof Eg , Ec , and energy balance components at both canopy
and soil surface. Eg decreased more rapidly than Ec , resulting in increase of soil temperature amplitude due to increased soil heat ﬂux. Agreements between simulated soil
water contents and temperatures and observed ones improved by accounting for the compensated root water uptake and larger unsaturated hydraulic conductivity for drier
range, due to increasing root water uptakes and evaporation
rates. For proper evaluation of the ratio of Eg and Ec in dry
period, detail investigations for soil hydraulic property androot water uptake property are needed.
1 Graduate School of Bioresources, Mie University, 1577 Kurimamachiya, Tsu, 514-8507, Japan. Corresponding author: 坂井 勝，三重大学大学院生物資源学研究科．
2 Institute for Agro-Environmental Sciences, NARO, 3-1-3 Kannondai,

Tsukuba, 305-8604, Japan.
2020 年 7 月 20 日受稿 2020 年 9 月 7 日受理

Key Words : evapotranspiration, canopy energy budget
model, HYDRUS-1D, radiation transmissivity, numericalsimulation

1. はじめに畑地の土壌水分量と地温は，種子の発芽や根の伸長，根の吸水や養分吸収，有機肥料の無機化など形態変化をもたらす微生物活性等，直接的にも間接的にも作物の生育に大きな影響を与える．畑地で生じる蒸発散は，根圏
土壌中の消費水分量の大部分を占め，潜熱消費をともな
い地温にも影響するため，土壌水分 · 熱動態を明らかにする上で重要な要素である．また，適正な灌漑計画のためにも，蒸発散の定量的な評価が求められている．蒸発散の測定法には，重量ライシメータ法，熱収支ボーエン比法，渦相関法などが挙げられるが，測定機器の設置や維持管理にコストと労力がかかる．そのため，Penman
式を用いて気象データから蒸発散を求め，畑地の消費水分量を算出する方法が提案されている（農林水産省構造改善局, 1997）．さらに近年では，Penman 式等で求まる
可能蒸発散速度を土壌水分移動式であるリチャーズ式に与えて数値計算を行い，土壌水分状態に依存する蒸発散速度の低下を考慮した上で灌漑計画を検討する手法も提案されている（例えば，弓削 · 阿南, 2016; 中村ら, 2016）．

Penman 式や Penman-Monteith 式は，入手が比較的容易な気象データから可能蒸発散速度を計算できるため，非常に有効な方法である．一方で，例えば HYDRUS-1D（Šimůnek et al., 2013）の様なシミュレーションプログラ
ムに適用する場合，得られた可能蒸発散速度 Ep を，土壌面の境界条件に与える可能蒸発速度 Egp と，根の吸水モデルに与える可能蒸散速度 Ecp に配分する必要がある．
実用上は作物と土壌面の放射の配分比とし，ベールの法則により作物の葉面積指数 LAI の関数で与えることが多い（Egp /Ep = exp(−k × LAI); Ritchie, 1972）
．この時，

4

土壌の物理性第 146 号

（2020）

放射の減衰に関わる係数 k の設定が必要である．既往の析を行い，可能蒸発散速度 Ep の Egp と Ecp への配分割

研究では文献値を用いた例が多いが（例えば，宮本ら,合と τ の関係を示した．そして，土壌水分量と地温の観2017；亀山ら, 2018），本来は作物の種類や栽植密度を反測値を再現する τ について検討を行った．さらに，長期映した値を用いる必要がある．放射の配分比を使わずに間の土壌乾燥過程について，蒸発速度と蒸散速度，およEgp と Ecp の比を与える方法として，農林水産省構造改善局（1997）が示す一般値を用いた例（中村ら, 2016），
FAO が提案する Dual crop coefﬁcient を利用した例（弓削 · 阿南, 2016），地表面の被覆率で配分した例（弓削ら,2018）があるが，土壌水分 · 地温の数値シミュレーショ

び熱収支成分変化を示し，土の水分移動特性や根の吸水

ン結果について十分に検討した例は少ない．
また，モデルの簡便化のために地温の日変化を考慮

特性の影響について試算を行った．

2． 方法
2.1 モデルの構築
2.1.1 水分 · 熱移動式

しない Penman 式や Penman-Monteith 式は，蒸発散速度

非等温条件下での土中の水分移動は，水蒸気移動を含

Ep の日変化を求めることはできず，土壌面温度の純放射量への影響による Ep の誤差が指摘されている（近藤,2011）．そのため，畑地の土壌水分 · 熱動態を詳細に解析

むリチャーズ式で定式化され，植物根の吸水を含む鉛

するには，地表面熱収支式を使用する必要がある（斎藤ら, 2007）．Saito et al.（2006）は，裸地条件下の地表面熱収支式と HYDRUS-1D を連結することで，蒸発と土壌
水分 · 熱移動の相互作用を考慮したプログラムを構築した．作物が生育する畑地では，さらに植生の影響も考慮
した群落熱収支モデルを適用する必要がある．Noborio

直 1 次元方向の支配方程式は次式で表される（斎藤ら,2007）．
∂θ
=
∂t
(1)
[ (
)
]
∂
∂h
∂T
∂h
∂T
K
+ 1 + KLT
+ Kvh
+ KvT
− S (z)
∂z
∂z
∂z
∂z
∂z

ここで，θ は体積含水率（m3 m−3 ），t は時間（s），z は上

et al.（2012）は，植被と土壌を合わせて一つの地表面と向き正の位置（m），h は土中水圧力（m），T は地温（◦ C），

して熱収支式を考える単層モデルを，牧草地の蒸発散のK（m s−1 ）と KLT（m2 ◦ C−1 s−1 ）は圧力勾配，温度勾配による液状水の透水係数，Kvh（m s−1 ）と KvT（m2 ◦ C−1s−1 ）は圧力勾配，温度勾配による水蒸気の透過係数である．KLT ，Kvh ，KvT の定式化は，斎藤ら（2007）に詳し
い．S は根の吸水速度（s−1 ）で，Feddes et al.（1974）が提案するモデルを用いて，可能蒸散速度 Ecp（m s−1 ）か

推定に用いている．しかし単層モデルでは，得られた蒸発散速度を蒸発と蒸散に配分する課題が残される．可能蒸発速度 Egp と可能蒸散速度 Ecp の配分は，厳密には作物と土に供給されるエネルギーと輸送抵抗によって決まる（Campbell, 1987）．そこで，群落熱収支モデルの一つである 2 層モデル（Watanabe, 1994）に注目した．2層モデルは群落の放射透過率 τ で土壌面と植被層に放射を配分し，各層の熱収支を考慮するため，蒸発と蒸散の日変化をそれぞれ個別に計算することが可能である．

ら求める．
∫ 0

Ec (t) =−LR∫ 0

S (z,t) dz = Ecp (t)−LR

α (h)β (z)dz(2)

Maruyama and Kuwagata（2010）は 2 層モデルを水田に適用し，稲の生長にともなう蒸発速度や蒸散速度の変化

ここで，Ec は実蒸散速度（m s−1 ），LR は根群域の深さ

や，群落温度や熱収支成分の変化などを土壌湿潤条件下（m），α は h に依存する水ストレス応答関数（–），β はにおいて明らかにしている．一方畑地では，土壌乾燥時の不飽和透水係数の低下による蒸発速度の減少や，乾燥
ストレスにともなう根の吸水速度の減少による蒸散速度

吸水強度分布である (m−1 )．

水蒸気移動を含む土中の熱の保存則は次式で表される（斎藤ら, 2007）．

の減少を考慮し，熱収支式を解く必要がある．また前述と同様に，放射透過率 τ は LAI の関数で与えられるが，そこで本研究では，蒸発散をともなう畑地の土壌水分∂ Cp T
∂ θv
+L0
=
∂t
∂t
[
]
∂
∂T
∂ qL T∂ qv T
∂ qv
λ
−Cw
−Cv
− L0
∂z
∂z
∂z
∂z
∂z

· 熱動態を明らかにするために，群落熱収支モデル（2 層モデル）と HYDRUS-1D を連結したプログラムの構築を

ここで，Cp ，Cw ，Cv は土，液状水，水蒸気の体積熱容量

畑地の土壌水分量や地温を再現し，かつ作物の生育状況を反映した簡易的な τ の設定方法の提案が必要である．

行った．モデルを検証するために，ダイズ栽培圃場で土壌水分 · 地温データ，気象データ，および作物生育データの観測を行った．気象データに基づき土壌水分量と地温を予測するためには，計算の境界条件となる可能蒸発速(3)(J m−3 ◦ C−1 )，θv は水蒸気相率（m3 m−3 ），L0 は水単位
体積当たりの蒸発潜熱（J m−3 ），λ は熱伝導率（W m−1◦ C−1 ）
，qL と qv は液状水と水蒸気フラックスである（m−1s ）．

度 Egp ，土中熱フラックス Gs ，および可能蒸散速度 Ecp

2.1.2 群落熱収支モデル

を正しく評価することが重要である．そこで，まず土壌

植生条件下の地表面熱収支を計算する群落熱収支モデ

が湿潤な乾燥初期に焦点を当て，放射透過率 τ の感度解ルとして，Maruyama and Kuwagata（2010）のモデルを

論文：群落熱収支モデルと HYDRUS-1D を連結した畑地の土壌水分・地温予測モデル

5

ルク輸送係数（–）である．湿潤土壌や田面水からの蒸発
潜熱 Lw Eg を計算する場合，土壌面の水蒸気密度 ρg に土壌面温度 Tg で決まる飽和水蒸気密度 ρsat (Tg ) を与える
（Watanabe, 1994）．本研究では，土壌の乾燥にともなう蒸発速度の低下を考慮するため，Tg と土壌面の土中水圧力 hg から決まる水蒸気密度 ρg （= ρsat × Hr ，Hr ：相対湿度）を計算に用いた．各層の顕熱 · 潜熱フラックス成分の和が地表面と大気間の顕熱 · 潜熱フラックスとなる（H = Hg + Hc ，Lw E = Lw Eg + Lw Ec ）．

4 つ の バ ル ク 輸 送 係 数（CHg ，CEg ，CHc ，CEc ）は ，Maruyama et al.（2017）に従って求めた．計算に必要なパラメータを Table 1 に示す．土壌面−大気間のバル
ク輸送係数 CHg と CEg ，および土壌面と植被面を合わせFig. 1 2 層モデルの概略図．
Schematic illustration of the double source model.

た地表面−大気間の顕熱交換と潜熱交換に対するバルク用いた．このモデルでは，土壌面と植被層の熱収支を個別に計算する 2 層モデルに基づき定式化している（Fig.

1）．

算され，得られた CHg ，CEg ，CH ，CE から植被面−大気間のバルク輸送係数 CHc ，CEc が計算される．この時，バ

Rng = Hg + Lw Eg + Gs

(4)

Rnc = Hc + Lw Ec
ここで，Rng と Rnc は土壌面と植被層に対する純放射量（W m−2 ），Hg と Hc は土壌面−大気間と植被面−大気間の顕熱フラックス（W m−2 ），Eg と Ec は土壌面蒸発速
度と蒸散速度（kg m−2 s−1 ），Lw は水単位質量当たりの蒸発潜熱（J kg−1 ），Gs は土中熱フラックスである（Wm−2 ）．各層の純放射成分は次式で表される．
()
Rng = 1 − αg τ RS + τ RL + (1 − τ ) σ Tc4 − σ Tg4)(
Rnc = (1 − α ) RS − 1 − αg τ RS()
+ (1 − τ ) RL + σ Tg4 − 2 (1 − τ ) σ Tc4の関数として表され，植生が小さい時は CHg と CEg が大きく，CHc と CEc は小さい．また，植生の生長にともない CHg と CEg が減少し，CHc と CEc が増加する．このように hc と LAI により，裸地条件を含め，植生の生長を考慮して蒸発速度と蒸散速度の変化を表すことができる．特に CE は気孔コンダクタンス gS （m s−1 ）にも依存し，気孔の開閉による蒸散速度の変化を表す．ここで量の制限がない場合の gS の最大値 gSmax （m s−1 ）を入力値とし，日射量の日変化に基づいて gS の日変化を計(5)算するモデルを用いた．

2 層モデルの計算に必要な主な気象観測データは日射
量 RS ，下向き長波放射 RL ，気温 Ta ，相対湿度 RH ，風τ は放射に対する群落の透過率（–），αg と α は土壌面と群落上のアルベド（–），Tg と Tc は土壌面と群落温度（◦ C），σ はステファン−ボルツマン定数（5.67 × 10−8 W−2m ◦ C−4 ）である．τ については，日射と長波放射で異なる値や関数で与えられるが，ここではモデルを簡略化するため，同じ値を用いることとした（近藤, 2011）．
各層の顕熱 · 潜熱フラックス成分は次式で表される．

Hg = caCHgU (Tg − Ta )Hc = caCHcU (Tc − Ta )()
Lw Eg = LwCEgU ρg (Tg , hg ) − ρa

ルク輸送係数は，群落高 hc （m）と葉面積指数 LAI （–）

は，Maruyama and Kuwagata（2008）が提案した，日射

ここで，RS と RL は群落上の日射と長波放射（W m−2 ），ここで，ca

輸送係数 CH と CE が，主に幾何学的な粗度長 h0s から計(6)

速 U であり，植生データは草高 hc ，葉面積指数 LAI ，群落透過率 τ ，気孔コンダクタンスの最大値 gSmax である．

2.1.3 HYDRUS-1D と 2 層モデルの連結
Saito et al.（2006）は，裸地条件下の地表面熱収支式とHYDRUS-1D を連結し，気象データを用いて，土壌中の水分 · 熱移動を計算するプログラムを作成した．このプ
ログラムでは，各計算時間ステップで純放射 Rng ，顕熱
フラックス Hg ，潜熱フラックス Lw Eg を計算し，蒸発フ
ラックス Eg ，および熱収支式から求まる土中熱フラック
ス Gs を，それぞれ水移動式（（1）式）と熱移動式（（3）
式）の土壌面境界条件として与える．この時，Rng ，Hg ，Lw Eg の計算には，土壌面の温度 Tg と土中水圧力 hg が

Lw Ec = LwCEcU (ρsat (Tc ) − ρa )必要である．そのため各計算時間ステップで，土壌水分は空気の体積熱容量（J m−3 ◦ C−1 ），U は風熱移動式の計算を 2 反復行っている．

量と地温の計算結果を更新しながら，熱収支式，水分 ·速（m s−1 ），Ta は気温（◦ C），ρg と ρa は土壌面と大気の本研究では，Saito et al.（2006）のプログラムを拡張する

sat は飽和水蒸気密度である（kg−3m ）．CHg ，CEg はそれぞれ土壌面−大気間の顕熱交換ことで，蒸散にともなう植物根の吸水を考慮した 2 層モと潜熱交換に対するバルク輸送係数（–），CHc ，CEc はそ

1D ver 4.16（Šimůnek et al., 2013）を連結したプログラ
ムを構築した．Fig. 2 に計算の流れを示す．ある計算時

水蒸気密度（kg m−3 ），ρ

れぞれ植被面−大気間の顕熱交換と潜熱交換に対するバデルのプログラム（Maruyama et al., 2017）と HYDRUS-

土壌の物理性

6

第 146 号

（2020）

Table 1 バルク輸送係数の計算に必要なパラメータ．
Parameters for calculating bulk transfer coefﬁcients.

Bulk transfer coefﬁcient
CE = fCE (hc , LAI, za , h0s , cd , ch , gS )CH = fCH (hc , LAI, za , h0s , cd , ch )CEg = fCEg (hc , LAI, za , h0s , cd , ch )CHg = fCHg (hc , LAI, za , h0s , cd , ch )CEc = fCEc (CE ,CEg )CHc = fCHc (CH ,CHg )roughness length of soil surface
h0s = fh0s (hgeo , u∗geo , kv )bulk stomatal conductancegS = gSmax /(1 + SabsH /Sabs )mean absorbed solar radiation per unit leaf area
Sabs = [(1 − α )Rs − (1 − αg )τ Rs ]/LAI

Input parameters
symboldescription

value

hc
LAI
za
cd
ch
hgeo
u∗geo
kv
α
αg
gSmax
SabsH
τ
Rs

canopy heightleaf area index
reference heightleaf transfer coefﬁcient for momentum
leaf transfer coefﬁcient for sensible heatgeometric roughness of bare soil surface
friction velocity for the groundkinematic viscosity of the air at 20 ◦ C
albedo of the entire surface
albedo of the soil surface
maximum bulk stomatal conductanceabsorbed solar radiation when gS = gSmax /2
radiation transmissivity
solar radiation

2m
0.20
0.06
0.025 m0.03 m s−1
1.53 × 10−5 m2 s−1
0.2
0.2
0.02 m s−1
66 W m−2height of ﬁeld observation

observed in this study
observed in this study
observed in this study

間ステップ t j において，まず気象データと植生データを

2 層モデルに入力し，純放射 · 顕熱 · 潜熱の各熱収支成分の計算を行う（
（5）式，
（6）式）
．この時，計算に必要な土壌面の温度 Tg と土中水圧力 hg は，HYDRUS-1D でj−1計算された 1 つ前の時間ステップの結果を用い（Tg ，hgj−1 ），土壌面と植被層の熱収支式が成り立つ群落温度
Tc の数値解を二分法で求める．ここで計算された蒸散潜
熱 Lw Ecp は，可能蒸散速度 Ecp を反映したものである．
得られた蒸発フラックス Eg と土中熱フラックス Gs を土壌面の水分 · 熱境界条件に，また Ecp を根の吸水モデル（（2）式）に与え，根の吸水を含む土壌中の水分 · 熱移動の数値計算を HYDRUS-1D で行う．Tg ，hg ，および実蒸
散速度 Ec の計算結果を用いて再度 2 層モデルの計算を行い，熱収支成分，および熱収支式が成り立つ Tc を求める．この時，蒸散潜熱 Lw Ec は，実蒸散速度 Ec を反映した値で固定する．この反復計算を行うことで，土壌乾燥によって蒸散速度が低下し潜熱項が低下した場合でも，熱収支式が成立する．そして，再度水分 · 熱移動の計算を行うことで，計算時間ステップ t j の結果を得る．構築した本プログラムでは，通常版の HYDRUS-1D の出力
Fig. 2 2 層モデルと HYDRUS-1D の計算の流れ図．
A ﬂow chart for simulation of the double source model and
HYDRUS-1D.

値に加え，植被層と土壌面の各熱収支成分（（4）式）と群落温度 Tc が出力される．また，通常版では入力値であった可能蒸散速度 Ecp も，計算過程で求められる．一論文：群落熱収支モデルと HYDRUS-1D を連結した畑地の土壌水分・地温予測モデル

7方，通常版では可能蒸発速度 Egp を入力し，系依存型境

2.3 数値計算

界条件で実蒸発速度 Eg が計算されるが（斎藤ら, 2006），

2.3.1 水分 · 熱特性と根の吸水特性

ここでは計算過程で Eg が直接計算されるため（（6）式），

土の水分特性曲線は，吸引法や加圧板法等の室内実

Egp は出力されない．

験で測定することが多い．一方，圃場では土の不均一性

2.2 現場観測や，乾燥密度や土の構造の変化，また乾湿のヒステリシ

三重大学附属農場内のダイズ栽培圃場で現場観測をス等の影響を受け，室内実験で得られた水分特性曲線と行った．圃場の土性は，表層の作土層が砂質ロームで，は異なる可能性がある．そこで，特に変化の大きい表層

およそ 25 cm 以深の下層はロームである．特に下層は，土については，圃場に設置した土壌水分センサーと水ポ

シルト質の土塊を含み，乾燥密度 ρb が 1.60 ∼ 1.65 Mg

テンシャルセンサーの測定値の関係から，水分特性曲線

m−3 と非常に高く，排水性が悪いという特徴がある．

を求めた．また土中水圧力 h < −5, 000 cm の乾燥領域の水分特性曲線については，鏡面冷却式湿度計 WP4-T
2018 年 6 月 6 日に各種土壌センサーとウェザーステー（
Decagon devices）を使って室内実験で別途測定した．
ションを設置した．土壌水分センサー（CS655，Campbell
圃場で採取した撹乱土壌を用い，任意の水分状態においScientiﬁc，ロッド長 12 cm）を深さ 0, 10，20 cm 深からて WP4-T で h を測定し，その時の含水比 ω を炉乾法で鉛直下向きに設置し，それぞれ 0 − 10, 10 − 20，20 −測定した．そして，対象土層の乾燥密度 ρb と水の密度30 cm 深の平均体積含水率とした．水ポテンシャルセン
ρw を使って，θ = ρb /ρw × ω の関係により体積含水率 θサー（MPS-6，Decagon Devices（現 TEROS-21，METERに変換することで，h と θ の関係を得た．30 cm 以深のgroup））を 5, 15 cm 深に，熱電対を 1，2，4，8, 10，20，30
cm 深に設置した．気象センサーとして，雨量計（TE525， 下層土の水分特性曲線は，吸引法，加圧板法，WP4-T で測定した．水分特性曲線と不飽和透水係数を表す不飽和
Campbell Scientiﬁc），温湿度計（HMP-155D，クリマテッ
水分移動特性関数には，Durner モデルを用いた（Durner,ク株式会社）と風向風速計（Met One 034B，Campbell
1994）．
Scientiﬁc）を高さ 2 m に，四成分放射収支計（NR01,
Campbell Scientiﬁc）を高さ 1.5 m に設置した．センサーはデータロガー（CR1000，Campbell Scientiﬁc）に接続
θ − θr= Seθし，30 分間隔でデータを保存した．
s − θr(7)−m−m= w1 [1 + (α1 h)n1 ] 1 + w2 [1 + (α2 h)n2 ] 2
ダイズ（フクユタカ）の播種を，2018 年 6 月 12 日に行った．その際，畝立ては行わず，株間 20 cm，条間 70

cm とし, 1 箇所に 3 粒播種した．ダイズの生育期間中，生育データとして草高， LAI ，地表面被覆率を定期的に測定した．地表面被覆率については，まず地上約 5.5 mの高さから，畑の 3 m × 3 m の領域をデジタルカメラで撮影した（有効画素数 1600 万画素）．得られた写真について画像解析を行い，色相に 50 ∼ 125（緑の領域），彩度と明度に 0 ∼ 255（全領域）を規定し，選択される面積
を被覆面積とした．そして，対象領域の面積（9 m2 ）に対する被覆面積の割合を地表面被覆率とした．
土壌水分センサーの体積含水率 θ と比誘電率 εr の関係式は，一般的には撹乱土壌を用いた室内実験で求める（例えば三石 · 溝口, 2015）．一方，室内実験で求めた検K (Se ) =

Ks

[
(
)m1 ]
[
(
)m2 ])2(
1m
1m+ w2 α2 1 − 1 − Se2/ 2(w1 Se1 + w2 Se2 )l w1 α1 1 − 1 − Se1/ 1(w1 α1 + w2 α2 )2(8)
ここで，θs は飽和体積含水率（m3 m−3 ），θr は残留体積
含水率（m3 m−3 ），Se は有効飽和度（–），wi （i = 1, 2）は重み係数（–）で w1 + w2 = 1 である．αi （m−1 ），ni ，mi （= 1 − 1/ni ）は曲線の形状を与えるパラメータ，Ksは飽和透水係数（m s−1 ），l は間隙結合係数である（–）．
飽和透水係数 Ks には, 100 cm3 採土缶で採取した不撹
乱試料により，変水位透水試験で測定した値を用いた．

量式を現場観測に適用した場合，圃場の土壌構造や不均間隙結合係数 l は水分特性曲線とは独立して不飽和透水

一性，またセンサーと土壌の接触具合から，必ずしも土係数関数の形状を与えるパラメータであるが，ここでは壌水分量を正しく測定することができない．そこで，武Mualem（1976）が提案する 0.5 とした．

藤ら（2015）が提案する原位置キャリブレーション法を土壌の熱特性は，熱特性計（KD2Pro，Decagon devices）

用いて，Topp 式の補正を行った．100 cm3 の採土円筒を使って測定した．撹乱土壌を任意の体積含水率に調整を用いて，土壌水分センサーの設置深度から不撹乱土壌

し, 100 cm3 採土缶に充填後，3 cm デュアルニードルセを採土し（2.5 − 7.5, 12.5 − 17.5，22.5 − 27.5 cm 深ンサー（SH-1，Decagon devices）を用いて，熱伝導率おから各深度 3 試料），採土時の θobs を炉乾法で求めた．
採土時の εr と Topp 式で求まる体積含水率 θTopp との比（θobs /θTopp ）を求め，Topp 式に乗じることで補正した検量線とした．本実験では，θobs /θTopp は深さによる違い
は小さく，3 本の水分センサーに対して 0.90 であった．

よび体積熱容量を測定した．得られた熱伝導率と体積含水率の関係に，Chung and Horton（1987）のモデルを適合した．

λ (θ ) = b1 + b2 θ + b3 θ 0.5(9)土壌の物理性

8

第 146 号

ここで，b1 ，b2 ，b3 は実験定数であり，実験値に適合して（2020）

得られた値はそれぞれ 0.183，3.646，0.146 W m−1 ◦ C−1

d，最小時間ステップを 0.00001 d，最大時間ステップを0.001 d とし，比較的小さな時間ステップが選択されるであった．本研究では，土層全体に対して同じ熱特性を設定とした．

与えた．

水分移動および熱移動の上端境界条件には，2 層モデ

根の吸水モデル（（2）式）の水ストレス応答関数 α (h)ルから計算される土壌面蒸発速度 Eg ，土中熱フラックスには線形型モデルを用いた（Feddes et al., 1974）．別途，

Gs をそれぞれフラックス境界条件として与え，下端境界
条件にはそれぞれ自由排水条件（dh/dz = 0），温度勾配ゼロ条件を与えた（dT /dz = 0）．また，2 層モデルで計算される可能蒸散速度 Ecp を根の吸水モデル（式（2））にポットを用いたダイズ栽培実験を行い，土壌乾燥にともなう蒸散速度の低下から α (h) を求めた（坂井ら, 2018）．
飽和近傍の吸水制限は考慮せず，土中水圧力 h3 = −1, 700

cm と h4 = −15, 000 cm（永久しおれ点）の間で，乾燥にともない α が 1 から 0 まで直線的に減少する関数とした．吸水強度分布 β (z) は，可能蒸散速度 Ecp を各深与えた．初期条件は体積含水率と地温の観測値を与え，観測値のない 30 cm 以深の下層については，20 − 30 cm
深の体積含水率と 30 cm 深の地温で一定値を与えた．

さの吸水速度に割り振るための，正規化された吸水分布

まず，期間①と期間②について土壌が十分に湿潤で，（m−1 ）である（根群域の積分値が 1）．実際の吸水分布

可能蒸発速度 Egp および可能蒸散速度 Ecp で蒸発散が生は，根量分布と必ずしも一致するとは限らない．そのたじる乾燥初期に注目した．気象データには，30 分毎のめ β (z) の決定は難しく，既往の研究では，根量分布に

日射量 RS ，下向き長波放射 RL ，気温 Ta ，相対湿度 RH ，基づいて β (z) を与えるケース（加藤 · 西村, 2015; 中村

風速 U の観測値を用いた．2 層モデルの計算に必要な気ら, 2016），根の伸長深さに対して均一な β (z) を与える

孔コンダクタンスの最大値 gSmax には，気孔コンダクタ

ケース（亀山ら, 2018），根圏の上部 20 % は一定で根圏下ンス測定装置（SC-1，Decagon devices）を用いて圃場の端に向けて直線的に減少する様な β (z) 関数型を用いる

ダイズについて測定を行い，両期間とも 0.02 m s−1 としケース（Hoffman and van Genuchten, 1983；弓削 · 阿南,

た．植生データには，期間①（被覆率およそ 30 %）は2016），土壌水分量の測定値に対して数値計算の逆解析で求めるケース（Hupet et al., 2003; Schelle et al., 2012）

草高 hc = 40 cm，葉面積指数 LAI = 1.0 の一定値を，期等，様々な方法がとられている．本実験では，計算対象

間②（被覆率およそ 60 %）は hc = 50 cm，葉面積指数

過去に同一圃場で測定したダイズ根の伸長深さが，発芽LAI = 2.5 の一定値を与えた．それぞれの期間に対して放射透過率 τ = 0.0，0.1，0.2，0.3，0.4，0.5，0.6，0.7，
0.8，0.9, 1.0 の一定値を与えた（2 期間で計 22 通りの

後 1 か月程度で乾燥密度が高く硬い下層土の上まで達し計算）．

期間において根量分布の測定は行っていない．そこで，ていたこと，また本実験の 20 − 30 cm の体積含水率の次に，期間②についてより長期の乾燥過程の計算を

低下具合から判断して，表層から 30 cm 深まで均一な行った．この時，ダイズの生長を考慮するために，hc とβ (z) を与えた．

LAI については観測値を与え，測定値間には内挿値を用

2.3.2 計算条件

いた．土壌乾燥にともなう土壌面蒸発と蒸散の減少，熱

本研究では蒸発散による土壌乾燥過程を計算するため収支成分の変化を示し，体積含水率と地温の観測値の再に，降雨イベント間の無降雨期間に焦点を当てた．ダイ

現性について検討した．

ズの生育期間における（6 月 19 日発芽），植物体の大きさが異なる期間① 7 月 9 日 ∼ 26 日（草高 40 ∼ 45 cm，

LAI 1 ∼ 2，被覆率 30 ∼ 40 %）と期間② 7 月 31 日 ∼ 8
月 12 日（草高 50 ∼ 70 cm，LAI 2.5 ∼ 4，被覆率 60 ∼80 %）の 2 回の無降雨期間を計算対象とした（後述 3.1

3. 結果と考察
3.1 ダイズの生育データと気象データを避けるため，十分な長さとして地表面から深さ 100 cm

Fig. 3 に生育期間中のダイズの草高，葉面積指数 LAI ，
地表面被覆率の変化を示す．ダイズは 6 月 12 日の播種
後，6 月 19 日に発芽し，8 月末には草高 90 cm，LAI 5，地表面被覆率 100 % となった．7 月 20 日頃に草高と LAIまでの 1 次元鉛直土層とした．また，水分特性曲線の違の増加に停滞が見られたが，無降雨期間が続いたことにいに基づいて（後述 3. 2 参照），0 − 10, 10 − 30，30

よる土壌乾燥が原因だと考えられる．Fig. 4 に 2018 年

cm 以下の 3 層の成層土層とした．計算メッシュは，変化の大きい表層 20 cm は 0.5 cm 間隔とし，20 cm 以深は 1 cm 間隔とした．HYDRUS では計算の収束状況に応じて，適切な計算時間ステップが選択される（Rassam etal., 2004）．ところが，熱移動や水蒸気移動を含む計算に

7 月上旬から 8 月中旬の日降水量，気温，相対湿度，純放射量，風速の観測値を示す．7 月 8 日に 50 mm，7 月
27 日 ∼ 30 日に 57 mm，8 月 13 日に 21 mm の降雨が参照）．
計算対象領域は，水分 · 熱移動の下端境界条件の影響

おいては，大きな時間ステップが選択されて計算が収束しないことがある．そこで，初期時間ステップを 0.001

あり，その間は概ね晴天であり土壌乾燥が進行した．また，7 月 29 日の強い風速は台風 12 号による．本研究では，期間①（7 月 9 日 ∼ 26 日）と期間②（7 月 31 日 ∼8 月 12 日）の 2 回の無降雨期間に注目する．

論文：群落熱収支モデルと HYDRUS-1D を連結した畑地の土壌水分・地温予測モデル

9Fig. 3 ダイズの生長にともなう草高，LAI ，地表面被覆率の変化．

Seasonal variations for crop height, LAI, and surface coverage of soybean.

Fig. 5 現場観測と室内実験で測定した水分特性曲線とDurner モデルによる適合曲線．灰色プロットは WP4-T による測定値．

Soil water retention curves obtained from experimental data
and ﬁtted by the Durner model. Gray plots indicate observed
data with WP4-T.

3.2 圃場の水分特性曲線
Fig. 5 に，現場観測，および室内実験から得られた水分特性曲線を示す．期間①と期間②の 2 回の無降雨期間（乾燥期間）について，0 − 10 cm 深の平均 θ と 5 cm
深の h, 10 − 20 cm 深の平均 θ と 15 cm 深の h の観測値の関係から，圃場の水分特性曲線を求めた．この時，

MPS-6 の測定可能領域（−100 cm > h > −10, 000 cm）
に対して，θ と h の関係を示した．MPS-6 の測定範囲の下限値は −2, 000 kPa（約 −20, 000 cm）とされているが（METER Group, 2019），今回の観測では乾燥が進行しても h = −10, 000 cm 程度で停滞した．これは，MPS-6 の多孔質板，または周囲の土壌の不飽和透水係数が低く，多孔質板と土壌の水分状態が十分に平衡していないためだと考えられる．そのため，ここでは h = −10, 000 cm
を下限値とした．また MPS-6 の測定値は，地温の日変化に合わせて大きな温度依存性を示した．そこで，日々の温度変化が小さく，極端に高温にならない午前 7 時の測定データのみを使用した．Fig. 5 の灰色プロットは，
鏡面冷却式湿度計 WP4-T を使った室内実験の測定値である．水分飽和時（Fig. 5 の h = −1 cm）の θ は，降雨直後に観測された θ の最大値とした．Fig. 5 に水分特性曲線の測定値に適合した Durner モデルを示し，Table 2にそのパラメータを示す．ここで，30 cm 以深の 3 層目の水分特性曲線については，他層の様な階段状の形状が確認できなかったため，w2 = 0 としてパラメータの少な
い van Genuchten モデル（van Genuchten, 1980）で適合
した．

3.3 乾燥初期の数値計算結果
Fig. 4 （a）日降水量，（b）気温，（c）相対湿度，（d）純放射量，（e）風速の観測値．
Observed data of (a) daily precipitation, (b) air temperature,(c) relative humidity, (d) net radiation, and (e) wind speed.

3.3.1 可能蒸発と可能蒸散の配分と放射透過率 τ の関係

7 月 8 日，および 7 月 30 日の降雨終了直後は，体積含水率の観測値に短時間で急激な低下が見られた．これ

土壌の物理性

10

第 146 号

（2020）

Table 2 計算に用いた Durner モデルの水分移動特性パラメータ．
Soil water hydraulic parameters of the Durner model used for simulations.
depth
0 − 10 cm
10 − 30 cm
30 − 100 cm

θr
0
0
0

θs
0.371
0.366
0.340

α1
0.030
0.021
0.00088

n1
2.00
1.98
1.25

α2
0.00050.00025
–

w2
0.55
0.67
0

n2
1.30
1.42
–

Ks (cm d−1 )
54.1
23.0
0.04

l
0.5（−0.5∗ ）0.5（−0.5∗ ）0.5∗ l = 0.5 より大きな不飽和透水係数を与える計算に用いた値．

3.3.2 土壌水分 · 地温変化に基づく放射透過率 τ の検討
ダイズが比較的小さい期間①（7 月 10 日 ∼ 7 月 12日）
，および比較的大きい期間②（7 月 31 日 ∼ 8 月 1 日）の，体積含水率 θ および地温の観測値（灰色線）と計算値（τ = 0.1，0.3，0.5，0.7，0.9）の比較を，それぞれFig. 7 と Fig. 8 に示す．θ の計算値は，観測値に合わせて 0 − 10, 10 − 20，20 − 30 cm 深に含まれる各計算
メッシュの平均値である．両期間とも，θ の観測値は蒸発散により各深さでやや減少した．地温の観測値は日振
Fig. 6 計算に用いた放射透過率 τ と Egp /Ep の関係．
Relationship between radiation transmissivity τ used for simulations and Egp /Ep .

は，マクロポアを通る排水や圃場側方への排水と考えら幅を示し，その幅は深い位置ほど小さい（期間①の 4 cm
深で約 10 ◦ C，20 cm 深で約 3 ◦ C の振幅）
．また，ダイズの生長による日射の遮断により，期間②の地温の振幅は，期間①に比べ小さかった（4 cm 深で期間①の約 10 ◦ C に対し，期間②は約 6 ◦ C の振幅）
．期間①の θ について異れ，再現するには飽和近傍の水分移動特性の詳細な検討なる τ による計算値を比較すると（Fig. 7），小さい τ ほど，土壌面蒸発が小さくなるため 0 − 10 cm 深の θ が等（例えば Durner, 1994）が必要である．ここでは期間

日中やや大きな値を示した．しかし，その差は 0.01 m3①と期間②について，それぞれ急激な排水以降である 7
月 10 日 0 時 ∼ 7 月 12 日 0 時と 7 月 31 日 0 時 ∼ 8 月 2

m−3 以下であり，いずれの τ を用いた場合でも観測値の減少をよく再現した．これは，根圏の深さが 30 cm と日 0 時を計算期間とした．期間①の計算初日である 7 月

やや浅いため，土壌面蒸発による水分減少と，根の吸水

10 日，および期間②の計算初日である 7 月 31 日は，30
cm 深以上の根圏土層の体積含水率が十分に高く（計算値で，θ > 0.22 m3 m−3 ，h > −140 cm），計算された蒸発
速度 Eg および蒸散速度 Ec は土壌乾燥により制限されなかった．そこで，この期間の Eg と Ec の計算結果を，可能蒸発速度 Egp と可能蒸散速度 Ecp とし，可能蒸発散速度 Ep（= Egp + Ecp ）を求めた．Fig. 6 に，任意の時間における可能蒸発散速度 Ep に対する Egp の割合（Egp /Ep ）と，計算に用いた放射透過率 τ の関係を示す．全体的に，Egp /Ep と τ はおよそ一致した．これは，Penman 式等で求めた Ep を，ベールの法則に基づいて，放射の透過率で近似的に Egp と Ecp に配分する方法（Ritchie, 1972）

による水分減少によって生じる水分分布の差が小さいためである．地温について異なる τ による計算値を比較す
ると，τ が大きいほど土壌面に到達する日射が大きいた
め，日中の地温上昇が大きかった．一方で，夜間は τ の
影響は小さかった．期間①では，τ = 0.5 ∼ 0.7 で地温の観測値の振幅を良く再現した．期間②の θ の計算値も，期間①と同様に τ の違いによる差は小さく，いずれの τを用いた場合も観測値を良く再現した（Fig. 8）．また，

τ = 0.3 ∼ 0.5 で地温の観測値の振幅を良く再現した．なお，τ = 0.7 や 0.9 といった大きな τ で，8 月 1 日の地温の計算値に急激な上昇が見られた．これは，期間②は期間①に比べ初期水分量が低く，土壌面の乾燥により土壌の妥当性を示している．特に，放射が植被層のみに配分

面蒸発が低下し，潜熱消費が低下したためである．

される τ = 0 の極端な条件では，Egp がほぼ 0 であるた

，7 月
測定した地表面被覆率 Sc に注目すると（Fig. 3）
10 日頃は約 30 %，7 月 31 日頃は約 60 % である．放め，Egp /Ep = 0 となる．放射が土壌面のみに配分される

τ = 1 でも同様に，Ecp がほぼ 0 であるため，Egp /Ep = 1となる．一方 0.4 ∼ 0.8 の中程度の τ で，Egp /Ep は τ よりやや小さい割合を示した．これは，蒸発と蒸散の割合
は放射の配分割合だけでなく，潜熱交換に関わるバルク
輸送係数にも依存するためである（Campbell, 1987）．

射が遮断される割合が被覆率で表されると仮定すると（τ = 1 − Sc ），放射透過率はそれぞれの期間で τ = 0.7，

0.4 であり，土壌水分量と地温の再現性から評価した τとおよそ一致する．これは，Sc の観測値に基づいて簡易的に τ を与えることで，土壌水分 · 地温変化の予測に適論文：群落熱収支モデルと HYDRUS-1D を連結した畑地の土壌水分・地温予測モデル

Fig. 7 7/10 ∼ 7/12 における体積含水率と地温の観測値と τ = 0.1，0.3，0.5，0.7，0.9 を用いた計算値の比較．
Observed soil water contents and soil temperatures, and simulated results using τ = 0.1, 0.3, 0.5, 0.7, and 0.9 from July
10th to July 12th.

Fig. 8 7/31 ∼ 8/1 における体積含水率と地温の観測値と τ = 0.1，0.3，0.5，0.7，0.9 を用いた計算値の比較．
Observed soil water contents and soil temperatures, and simulated results using τ = 0.1, 0.3, 0.5, 0.7, and 0.9 from July
31th to August 1st.

11

土壌の物理性

12

第 146 号

（2020）

Fig. 9 7/31 ∼ 8/11 の体積含水率および地温の観測値と計算値の比較．灰色は観測値，黒実線は l = 0.5 で補填吸水を考慮しない（ωc = 1）場合，黒点線は l = −0.5 で補填吸水を考慮する（ωc = 0.5）場合の計算値．
Observed and simulated soil water contents and soil temperatures from Jury 31th to August 11th. Gray lines indicate
observed data. Black solid and dashed lines indicate simulated results with l = 0.5 and uncompensated root water
uptake (ωc = 1), and with l = −0.5 and compensated root water uptake (ωc = 0.5), respectively.

切な Egp と Ecp ，および土中熱フラックス Gs を計算できることを示している．近年では無人航空機（UAV）等を使用して，より容易に地表面被覆率の測定が可能であることから，植生の生長にともなう τ の変化を簡便に与える方法となり得る．また同様に，可能蒸発速度と可能蒸散速度の配分率としても有用な方法となる可能性がある．
3.4 長期間の乾燥過程の数値計算結果

Fig. 10 に，蒸散速度 Ec と蒸発速度 Eg の日変化の計
算値（黒色実線）を，可能蒸散速度 Ecp と可能蒸発速度
Egp （灰色線）とともに示す．ここで Ecp と Egp は，30cm 深の根圏が十分に湿潤な条件で別途計算した Ec とEg の結果である．具体的には，初期土中水圧力 h に土壌表面 h = −30 cm で下端 h = 70 cm の水理学的平衡条
件を与え，下端水分境界条件を一定圧力境界条件とし，30 cm 以深の 3 層目の飽和透水係数 Ks を 10 cm d−1 と3.4.1 土壌水分量 · 地温および蒸発散速度
大きな値に変更することで，地下水面が 30 cm 深と浅い位置に留まる計算を行った．ダイズの生長による τ のより長期間の土壌乾燥過程について，ダイズの生育を低下にともない，Ecp は徐々に増加し，逆に Egp は減少考慮した計算を行った．期間②について，7 月 31 日 0
した．Ec は 8 月 2 日までは Ecp と一致した．8 月 3 日
時 ∼ 8 月 12 日 0 時を計算期間とし，被覆率 Sc の測定
値（Fig. 3）を用いて放射透過率 τ を τ ＝ 1 − Sc で与え， 以降は乾燥ストレスにより低下し，8 月 11 日の日中では Ecp がおよそ 0.08 cm h−1 であるのに対し，Ec は 0.02
計算期間中の τ の減少を考慮した．草高 hc と葉面積指cm h−1 であった．Eg は 8 月 2 日 11 時頃までは Egp と数 LAI の増加も同様に観測値を与え，生育データの測定一致したが，その後は減少し，8 月 5 日以降は 0.002 cm
値間は内挿値を用いた．Fig. 9 に，体積含水率 θ と地温の観測値（灰色線）と計算値（黒色実線）の比較を示す． h−1 以下と非常に小さい値を示した．表層の不飽和透水係数の低下により土壌面への水分供給が制限され（Jury
7 月 31 日 ∼ 8 月 1 日の土壌乾燥初期については，Fig.
and Horton, 2006），土壌面の土中水圧力 hg の低下にと
8 と同様に計算は θ の低下および地温の日変化を良く再もない，土壌面の水蒸気密度 ρg が低下したためである
現した．一方，8 月 2 日以降の乾燥が進行した期間については，計算値は観測値の変化を良く再現したが, 10 − （式（6））．土壌乾燥による低下は，土壌面の乾燥に依存する Eg で顕著であり，Ec の低下は湿潤な下層からの吸20，20 − 30 cm 深の θ をやや過大評価し，地温も各深水も含むため比較的緩やかであった．土壌が湿潤な条件
さで約 3 ◦ C 過大評価した．

論文：群落熱収支モデルと HYDRUS-1D を連結した畑地の土壌水分・地温予測モデル

Fig. 1113

7/31 ∼ 8/6 の（a）植被層と（b）土壌面における熱収支成分の計算値．

Simulated energy balance components (a) at the canopy and (b)at the soil surface from Jury 31th to August 6th.
7/31 ∼ 8/11 の（a）蒸散速度，（b）蒸発速度の計算
値，および（c）日蒸発量への配分割合 Eg /E ．灰色は可能蒸散 Ecp · 蒸発速度 Egp ，黒実線は l = 0.5 で補填吸水を考慮しない（ωc = 1）場合，黒点線は l = −0.5 で補填吸水を考慮す
る（ωc = 0.5）場合．

Fig. 10

Simulated (a) transpiration rates, (b) evaporation rates, and (c)partition ratio of daily evaporation amount Eg /E from Jury 31th
to August 11th. Gray lines indicate potential transpiration Ecp
and evaporation rate Egp . Black solid and dashed lines indicate
simulated results with l = 0.5 and uncompensated root water
uptake (ωc = 1), and with l = −0.5 and compensated root water
uptake (ωc = 0.5), respectively.

下では，可能蒸発速度 Egp と可能蒸散速度 Ecp の配分
割合の違いが土壌水分量の計算値に与える影響は小さい（Fig.7，8）．一方，土壌乾燥の進行にともなう Eg と Ecの低下は Egp と Ecp に依存するため，Egp と Ecp の配分
割合は土壌乾燥時の土壌水分量の予測に重要であると言える．

Fig. 10（c）に，日蒸発散量（E = Eg + Ec ）に対する日蒸発量の割合 Eg /E を示す（黒色実線）．ここでは，Fig.

10（a），（b）に示す蒸散 · 蒸発速度の日変化の日積算量（cm）を用いて，割合を計算した．可能蒸発散で蒸発 ·蒸散が生じる 7 月 31 日（Eg = 0.26 cm，E = 0.72 cm,

Eg /E = 0.36）と 8 月 1 日（Eg = 0.22 cm，E = 0.66 cm,
Eg /E = 0.33）では，Eg /E は与えた放射透過率 τ とほぼ一致する．その後は，Ec に比べ Eg が大きく低下するた
め，Eg /E は急激に減少し，8 月 6 日では Eg /E = 0.03
（Eg = 0.01 cm，E = 0.38 cm），8 月 9 日以降は Eg /E = 0を示した．なお，8 月 7 日に一時的な Eg /E の増加が見られたが，曇天のため日射量が小さく（Fig. 4d），Ec が低下したためである（Fig. 10a）．

3.4.2 植被面と土壌面の熱収支成分
Fig. 11 に 7 月 31 日 ∼ 8 月 6 日における，植被層と土壌面の熱収支成分の計算値を示す．図中で，純放射は下向きを正，顕熱，潜熱，土中熱フラックス項は上向きを正で示している．この期間において，τ が 0.4（被覆率 60%）から 0.3（被覆率 70 %）に減少するため，植被層の純放射 Rnc は時間とともに増加し，土壌面の純放射 Rngは減少した．可能蒸散速度 Ecp と可能蒸発速度 Egp で蒸発散が進行する 7 月 31 日 ∼ 8 月 1 日では（Fig. 10），
潜熱輸送量（Lw Ec ，Lw Eg ）が大きく，特に蒸散による潜熱交換 Lw Ec は Rnc とほぼ一致した．一方で顕熱輸送量（Hc ，Hg ）は両層で小さい．顕熱輸送量は気温 Ta との温度差から計算されるため（式（6）），土壌が十分に湿潤な期間では，群落温度 Tc および土壌面温度 Tg と，気温 Ta の差が比較的小さいことを示している．土壌乾燥が進行し Ec と Eg が低下すると，潜熱輸送量が低下し，顕熱輸送量が増加した．特に 8 月 4 日以降は Eg が極端に小さいため（Fig. 10），Lw Eg は終日ほぼ 0 W m−2 を示した．土中熱フラックス Gs は，Rng と Lw Eg ，Hg の収支から求められる．土壌を熱するために使われる熱量（Gs の負値，−Gs ）は，Lw Eg の低下とともに一度増加し（日最大値は，7 月 31 日で 100 W m−2 ，8 月 4 日で 150

W m−2 ），それにともない地温の計算値も上昇した（Fig.
9）．乾燥が進行し Lw Eg がほぼ 0 W m−2 となる期間では，τ の低下にともなう Rng の低下により，−Gs はやや減少した．地温の変動は，群落を透過して土壌面に達す
る放射量と土壌面蒸発速度の両方に影響を受けることを示している．

土壌の物理性

14

第 146 号

3.4.3 植物根の吸水特性と土壌水分移動特性の検討
8 月 2 日以降に，体積含水率の計算値が観測値を過大評価する（Fig. 9）原因の一つとして，乾燥ストレスに対（2020）

加した．それにともない，8 月 2 日 ∼ 4 日の地温の再現性が向上した．一方，8 月 5 日以降では l = 0.5 と同様にEg が極端に低下するため，再び地温の観測値を過大評価する根の吸水特性の変化が考えられる．特に，乾燥進行した．乾燥が進行した際の地温の再現性の向上のために時に 10 − 30 cm 深の根圏下部で観測値の低下が計算値

は，乾燥領域の K や熱特性等について，今後より詳細なより大きいことから，根圏上部の乾燥による根の吸水速検討が必要だと言える．Eg /E は，土壌が湿潤な 7 月 31

度低下を湿潤な下部で補う補填吸水が生じた可能性があ日と 8 月 1 日は変化しない．土壌乾燥が進行した 8 月 2る．そこで，HYDRUS-1D に実装されている補填吸水モ

日以降で大きな値を示し，8 月 6 日では Eg /E = 0.06 とデル（Jarvis, 1989）を適用した計算を行った．これは，なった（日積算値 Eg = 0.03 cm，E = 0.41 cm）．土壌乾

可能蒸散速度 Ecp と補填吸水を考慮しない場合の実蒸散

燥時の畑地の蒸発 · 蒸散の配分割合を明らかにするため

速度 Ec の比で表される水ストレス指数 ω （= Ec /Ecp ）には，土壌の水分移動特性や根の吸水特性の把握が重要を用いて，（2）式の Ecp を増大させるモデルである．水であることを示している．

ストレス応答関数 α が小さい乾燥部に比べ，α が大きい
湿潤部の方が根の吸水速度 S が増大する．

S (z,t) =Ecp (t)α (h)β (z)max (ω , ωc )(10)

4. おわりに本研究では，気象データに基づき蒸発散をともなう畑地の土壌水分 · 熱動態を予測するプログラムを構築する
ために，群落の放射透過率 τ で土壌面と植被層に放射を

ここで，ωc は根の適応係数と呼ばれる ω の閾値であり配分し各層の熱収支を考慮する 2 層モデルと HYDRUS-（0 ≤ ωc ≤ 1），ωc =1 の場合は補填吸水を考慮しないモ
デル，ωc = 0 の場合は常に Ecp で蒸散が生じる完全補填

1D の連結を行った．そして，ダイズ栽培圃場で観測した気象データ，生育データ，土壌水分量 · 地温データに

吸水モデルとなる．ここでは ωc = 0.5 として計算を行っ適用し，無降雨期間の土壌乾燥過程において構築したプ

た．また，8 月 2 日以降に地温の計算値が観測値を過大

ログラムの検証を行った．

評価する原因の一つに，蒸発潜熱 Lw Eg の低下にともな

まず蒸発散が気象条件および植生条件にのみ依存すう土中熱フラックス −Gs の増加が考えられる（Fig. 11）．

る，土壌が湿潤な乾燥初期に注目し，τ の感度解析を行っ

そこで，土壌面への土中水分フラックスを増加し，土壌た．可能蒸発散速度 Ep の可能蒸発速度 Egp と可能蒸散

面蒸発速度 Eg の減少を抑制するために，乾燥領域の不速度 Ecp への配分割合と τ の関係を検討したところ，配飽和透水係数 K を大きくした計算を行った．不飽和透水

分割合 Egp /Ep と τ はほぼ一致した．これは，Penman 式係数モデル（式（8））に含まれる間隙結合係数 l は，水等で計算される Ep を，ベールの法則で Egp と Ecp に配分特性曲線とは独立して，土中水圧力 h の低下に対する分する方法の妥当性を示している．続いて，ダイズの大K の低下の傾きを与えるパラメータである（Sakai et al.,2009）．そこで 30 cm 深以上の表層に対して l = −0.5を与えることで（Table 2），乾燥時に l = 0.5 の場合に比べて大きな K を示す不飽和透水係数を用いた計算を行った（例えば h = −10, 000 cm において，l = 0.5 の時
K = 8.90 × 10−6 cm d−1 ，l = −0.5 の時 K = 2.68 × 10−5
cm d−1 ）．

きさが異なる 2 期間について，体積含水率の低下と地温

Fig. 9 の黒色点線に体積含水率と地温の計算値を，Fig.
10 の黒色点線に蒸散速度 Ec と蒸発速度 Eg の計算値，および日蒸発散量に対する日蒸発量の割合 Eg /E の変化を示す．補填吸水を考慮することで 8 月 3 日から 8 月 5 日までの Ec が増加し，補填吸水を考慮しない場合に比べ，根圏からより多くの水が吸水 · 蒸散によって消費された．
この時, 10 − 20 cm 深の体積含水率の計算値は観測値とよく一致し，20 − 30 cm 深については依然過大評価した

フラックス Gs を計算できることが示された．

の日変化の現場観測値と計算値の比較を行った．乾燥初期においては，τ の違いによる体積含水率の計算値の違いは小さく，地温の日振幅に大きな影響を与えた．また
現場観測値の再現性から，地表面被覆率 Sc の観測値に基づいて簡易的に τ を与えることで（τ = 1 − Sc ），土壌水
分 · 地温変化の予測に適切な Egp と Ecp ，および土中熱
さらに，被覆率の観測値に基づいて τ の変化を与え，より長期間の乾燥過程の計算を行った．土壌乾燥による蒸発速度 Eg の低下は，蒸散速度 Ec の低下に比べ急激であり，蒸発散 E（ = Eg + Ec ）に対する割合 Eg /E も時間とともに低下した．また，Eg の低下にともない土壌面の蒸発潜熱 Lw Eg が小さくなるため，地温が上昇した．土壌乾燥の進行にともなう Eg と Ec の低下は，Egp と Ecp

ものの，観測値の再現性は向上した．なお補填吸水を考に依存するため，Egp と Ecp の配分割合は土壌乾燥時の慮することで，土壌乾燥の進行が促進されるため，8 月 8

土壌水分量 · 地温の予測に重要であると言える．表層土

日以降では Ec が低下した．また，乾燥領域の K を大き壌の低土中水圧力領域の不飽和透水係数を大きくするこ

くすることで（l = −0.5），土壌面の水蒸気密度の低下がとで Eg の低下を抑制し，また根の吸水モデルに補填吸抑制され（式（6）），8 月 2 日 ∼ 4 日で Eg の計算値が増

水を考慮して根圏下層の吸水量を増加させることで，体論文：群落熱収支モデルと HYDRUS-1D を連結した畑地の土壌水分・地温予測モデル
積含水率と地温の観測値の再現性が向上した．今後，根15

Meteorol., 148: 1161–1173.

圏土壌の水分移動特性や根の吸水特性をより詳細に検Maruyama, A. and Kuwagata, T. (2010): Coupling land sur-討し，プログラムに反映することで，土壌乾燥時においface and crop growth models to estimate the effects ofても Eg と Ec の配分割合を正しく評価できると考えら

changes in the growing season on energy balance and wa-

れる．

ter use of rice paddies. Agric. Forest Meteorol., 150: 919–
930.

謝辞

Maruyama, A., Nemoto, M., Hamasaki, T., Ishida, S. andKuwagata, T. (2017): A water temperature simulation

本研究は JSPS 科研費 17K08002 の助成を受けたもの

model for rice paddies with variable water depths. Water

です．また，モデルの構築に有益なご助言を頂いた農研

Resour. Res., 53: doi.org/10.1002/2017WR021019.

機構の桑形恒男上席研究員，および栽培管理等でご尽力

METER Group (2019): TEROS21 operator’s manual.

頂きました三重大学大学院生物資源学研究科の長菅輝義
教授に心より感謝申し上げます．

引用文献

三石正一, 溝口 勝 (2015): 一次補正式を用いた静電容量型

GS3 土壌水分センサーが測定した体積含水率の補正. 土壌の物理性, 130: 19–25.
宮本輝仁, 岩田幸良. 亀山幸司 (2017): HYDRUS を用いた湿潤地域におけるマルチ栽培下の土壌水分動態解明. 土壌の物理性, 136: 3–14.

Campbell, G.S.（中野政詩, 東山勇監訳）(1987): パソコンで学ぶ土壌物理学, p. 131. 鹿島出版会, 東京.

Chung, S.O. and Horton, R. (1987): Soil heat and waterﬂow with a partial surface mulch. Water Resour. Res., 23:
2175–2186.
Durner, W. (1994): Hydraulic conductivity estimation for soils
with heterogeneous pore structure. Water Resour. Res., 30:
211–223.
Feddes, R.A., Bresler, E. and Neuman, S.P. (1974): Field testof a modiﬁed numerical model for water uptake by rootsystems. Water Resour. Res., 10: 1199–1206.

Mualem, Y. (1976): A new model for predicting the hydraulicconductivity of unsaturated porous media. Water Resour.
Res., 12: 513–522.
武藤由子, 渡辺晋生, 山本清仁, 倉島栄一 (2015): 5TE セン
サーで測定した土壌水分量の簡易補正と電気伝導率の検証. 農業農村工学会論文集, 296: I 9–I 17.
中村公人, 大串祥子, 池浦康広, 田中宣多 (2016): 畑地用水計
画のための HYDRUS-1D を用いた土壌水分移動解析における土壌水分特性パラメータの推定例. 土壌の物理性,

134: 25–40.
Noborio, K., Ono, K., Mano, J., Miura, T., Satta, N., Koga, K.,

Hoffman, G.J. and van Genuchten, M.Th. (1983): Soil prop-and Mukaida, Y. (2012): Evaluation of energy-balance-

erties and effect water use: Water management for salinity

based evapotranspiration in a grass ﬁeld. J. Jpn. Soc. Soilcontrol, In Taylor, H.M., Jordan, W.R. and Sinclair, T.R.

Phys., 122: 15–21.

(eds.) Limitation and efﬁcient water use in crop production, pp. 73–85. Am. Soc. of Agron., Madison, WI.
Hupet, F., Lambot, S., Feddes, R.A., Van Dam, J.C. and Van-

農林水産省構造改善局 (1997): 土地改良事業計画設計基準
計画「農業用水（畑）」基準書

技術書.（社）農業土木

学会.

clooster, M. (2003): Estimation of root water uptake pa-

Rassam, D., Šimůnek, J. and van Genuchten, M. Th. (取出rameters by inverse modeling with soil water content data.

伸夫 · 井上光弘監訳) (2004): HYDRUS-2D による土中Water Resour. Res., 39: doi.org/10.1029/2003WR002046.

の不飽和流れの計算. 農業土木学会土壌物理研究部会

Jarvis, N.J. (1989): A simple empirical model of root water
uptake. J. Hydrol., 107: 57–72.
Jury, W.A. and Horton, R.（取出伸夫監訳, 井上光弘, 長 裕幸,
西村 拓, 諸泉利嗣, 渡辺晋生訳）(2006): 土壌物理学 —土中の水 · 熱 · ガス · 化学物質移動の基礎と応用 —, pp.

154–157. 築地書店, 東京.
亀 山 幸 司, 岩 田 幸 良, 宮 本 輝 仁, 佐 々 木 康 一 (2018):

HYDRUS-1D を用いた砂丘未熟土畑における灌水操作方法の検討. 土壌の物理性, 139: 3–12.
加藤千尋, 西村 拓 (2015): 農地土壌水分状態予測に向けた

GCM 予測値の時間ダウンスケーリング手法の検討. 農業農村工学会論文集, 295: 11–19.
近藤純正 (2011): 水環境の気象学, pp. 222–234. 朝倉書店,
東京.

Maruyama, A. and Kuwagata, T. (2008): Diurnal and seasonalvariation in bulk stomatal conductance of the rice canopy
and its dependence on developmental stage. Agric. Forest

HYDRUS グループ, 東京.
Ritchie, J.T. (1972): Model for predicting evaporation froma row crop with incomplete cover. Water Resour. Res., 8:
1204–1213.
斎藤広隆, 坂井 勝, Šimůnek, J., 取出伸夫 (2006): 不飽和土中の水分移動モデルにおける境界条件. 土壌の物理性, 104:

63–73.
Saito, H., Šimůnek, J. and Mohanty, B.P. (2006): Numericalanalysis of coupled water, vapor, and heat transport in thevadose zone. Vadose Zone J., 5: 784–800.
斎藤広隆, Šimůnek, J., 取出伸夫 (2007): 裸地土壌の水分 · 温度長期変動予測 — 近似した気象データによる表面境界
条件の設定 —. 土壌の物理性, 107: 79–96.
坂井 勝, 大西一平, 取出伸夫 (2018): ポット栽培実験における蒸散速度の日変化に基づく水ストレス応答関数の推定. 農業農村工学会全国大会講演要旨集, pp. 330–331.

土壌の物理性

16

第 146 号

（2020）

Sakai, M., Toride, N. and Šimůnek, J. (2009): Water and vapor

van Genuchten, M.Th. (1980): A closed-form equation formovement with condensation and evaporation in a sandy

predicting the hydraulic conductivity of unsaturated soils.

column. Soil Sci. Soc. of Am. J., 73: 707–717.

Soil Sci. Soc. Am. J., 44: 892–898

Schelle, H., Iden, S.C., Fank, J. and Durner, W. (2012):

Watanabe, T. (1994): Bulk parameterization for a vegetated

Inverse estimation of soil hydraulic and root distribu-

surface and its application to a simulation of nocturnaltion parameters from lysimeter data. Vadose Zone J., 11:drainage ﬂow. Boundary-Layer Meteorol., 70: 13–35.
弓削こずえ, 阿南光政 (2016): 畑地における作物根の分布と

doi.org/10.2136/vzj2011.0169.
Šimůnek, J., Šejna, M., Saito, H., Sakai, M. and van
Genuchten, M.Th. (2013): The Hydrus-1D software pack-

土壌面蒸発の変動を考慮した土壌水分動態解析と消費水分量の定量化. 土壌の物理性, 134: 17–24.

age for simulating the movement of water, heat, and multi-

弓削こずえ, 濵田耕佑, 阿南光政, 濱上邦彦, van Dam, J.C.,

ple solutes in variably saturated media, version 4.17, HY-

Kroes, J.G. (2018): 干拓農地における水分ストレスの評

DRUS software series 3, p. 342. Department of Environ-

価と灌漑スケジューリングの最適化. 農業農村工学会論

mental Sciences, University of California Riverside, River-

文集, 306: I 55–I 62.

side, CA.

要旨

蒸発散は，畑地の土壌水分 · 熱動態を明らかにする上で重要な要素である．本研究では，群落熱収支モ
デルの 1 つである土壌面と植被層の熱収支式を表す 2 層モデルと HYDRUS-1D を連結し，蒸発と蒸散をともなう土壌水分量 · 地温変化の予測プログラムを構築した．ダイズ栽培圃場の観測値に対し，モデ
ルの検証，放射透過率 τ の検討を行った．可能蒸発 Egp · 蒸散 Ecp が生じる乾燥初期の土壌水分量と地温の再現性から，地表面被覆率に基づいて簡易的に τ を与えることで，適切な Egp と Ecp の配分割合を計算できることが示された．また，より長期間の乾燥過程の数値計算を行い，土壌乾燥による蒸発 · 蒸散
速度の低下，熱収支成分の変化を示した．土壌水分量と地温の再現性から，補填吸水の考慮による根圏
下層の吸水促進の必要性や，大きな不飽和透水係数を与え蒸発速度低下を抑制する必要性が示された．
キーワード：蒸発散，群落熱収支モデル，HYDRUS-1D，放射透過率，数値計算