From 79cc419896b71a3ef79bef3ea1386e85a11cce3d Mon Sep 17 00:00:00 2001 From: Yilin Fang Date: Thu, 5 Jun 2025 16:26:37 -0700 Subject: [PATCH 1/9] Add documentation for hydrology and plant hydraulics, converted from CLM5 technical note. --- components/elm/docs/tech-guide/hydrology.md | 961 ++++++++++++++++++++ components/elm/docs/tech-guide/index.md | 2 + components/elm/docs/tech-guide/phs.md | 679 ++++++++++++++ 3 files changed, 1642 insertions(+) create mode 100644 components/elm/docs/tech-guide/hydrology.md create mode 100644 components/elm/docs/tech-guide/phs.md diff --git a/components/elm/docs/tech-guide/hydrology.md b/components/elm/docs/tech-guide/hydrology.md new file mode 100644 index 000000000000..a8044e365e0c --- /dev/null +++ b/components/elm/docs/tech-guide/hydrology.md @@ -0,0 +1,961 @@ + +<<<<<<< HEAD +======= +The model parameterizes interception, throughfall, canopy drip, snow +accumulation and melt, water transfer between snow layers, infiltration, +evaporation, surface runoff, sub-surface drainage, redistribution within +the soil column, and groundwater discharge and recharge to simulate +changes in canopy water $\Delta W_{can,\,liq}$, canopy snow water +$\Delta W_{can,\,sno}$ surface water $\Delta W_{sfc}$, snow water +$\Delta W_{sno}$, soil water $\Delta w_{liq,\, i}$, and soil ice +$\Delta w_{ice,\, i}$, and water in the unconfined aquifer +$\Delta W_{a}$ (all in kg m^-2^ or mm of H~2~O) +(`Figure Hydrologic processes`{.interpreted-text role="numref"}). + +The total water balance of the system is + +$$\begin{aligned} +\begin{array}{l} {\Delta W_{can,\,liq} +\Delta W_{can,\,sno} +\Delta W_{sfc} +\Delta W_{sno} +} \\ {\sum _{i=1}^{N_{levsoi} }\left(\Delta w_{liq,\, i} +\Delta w_{ice,\, i} \right)+\Delta W_{a} =\left(\begin{array}{l} {q_{rain} +q_{sno} -E_{v} -E_{g} -q_{over} } \\ {-q_{h2osfc} -q_{drai} -q_{rgwl} -q_{snwcp,\, ice} } \end{array}\right) \Delta t} \end{array} +\end{aligned}$$ + +where $q_{rain}$ is the liquid part of precipitation, $q_{sno}$ is the +solid part of precipitation, $E_{v}$ is ET from vegetation (Chapter +`rst_Momentum, Sensible Heat, and Latent Heat Fluxes`{.interpreted-text +role="numref"}), $E_{g}$ is ground evaporation (Chapter +`rst_Momentum, Sensible Heat, and Latent Heat Fluxes`{.interpreted-text +role="numref"}), $q_{over}$ is surface runoff (section +`Surface Runoff`{.interpreted-text role="numref"}), $q_{h2osfc}$ is +runoff from surface water storage (section +`Surface Runoff`{.interpreted-text role="numref"}), $q_{drai}$ is +sub-surface drainage (section +`Lateral Sub-surface Runoff`{.interpreted-text role="numref"}), +$q_{rgwl}$ and $q_{snwcp,ice}$ are liquid and solid runoff from glaciers +and lakes, and runoff from other surface types due to snow capping +(section +`Runoff from glaciers and snow-capped surfaces`{.interpreted-text +role="numref"}) (all in kg m^-2^ s^-1^), $N_{levsoi}$ is the number of +soil layers (note that hydrology calculations are only done over soil +layers 1 to $N_{levsoi}$; ground levels $N_{levsoi} +1$ to $N_{levgrnd}$ +are currently hydrologically inactive; +`(Lawrence et al. 2008) `{.interpreted-text +role="ref"} and $\Delta t$ is the time step (s). + +::: {#Figure Hydrologic processes} +![Hydrologic processes represented in CLM.](hydrologic.processes.png) +::: + +## Canopy Water {#Canopy Water} + +Liquid precipitation is either intercepted by the canopy, falls directly +to the snow/soil surface (throughfall), or drips off the vegetation +(canopy drip). Solid precipitation is treated similarly, with the +addition of unloading of previously intercepted snow. Interception by +vegetation is divided between liquid and solid phases $q_{intr,\,liq}$ +and $q_{intr,\,ice}$ (kg m^-2^ s^-1^) + +$$q_{intr,\,liq} = f_{pi,\,liq} \ q_{rain}$$ + +$$q_{intr,\,ice} = f_{pi,\,ice} \ q_{sno}$$ + +where $f_{pi,\,liq}$ and $f_{pi,\,ice}$ are the fractions of intercepted +precipitation of rain and snow, respectively + +$$f_{pi,\,liq} = \alpha_{liq} \ tanh \left(L+S\right)$$ + +$$f_{pi,\,ice} =\alpha_{sno} \ \left\{1-\exp \left[-0.5\left(L+S\right)\right]\right\} \ ,$$ + +and $L$ and $S$ are the exposed leaf and stem area index, respectively +(section `Phenology and vegetation burial by snow`{.interpreted-text +role="numref"}), and the $\alpha$\'s scale the fractional area of a leaf +that collects water +(`Lawrence et al. 2007 `{.interpreted-text +role="ref"}). Default values of $\alpha_{liq}$ and $\alpha_{sno}$ are +set to 1. Throughfall (kg m^-2^ s^-1^) is also divided into liquid and +solid phases, reaching the ground (soil or snow surface) as + +$$q_{thru,\, liq} = q_{rain} \left(1 - f_{pi,\,liq}\right)$$ + +$$q_{thru,\, ice} = q_{sno} \left(1 - f_{pi,\,ice}\right)$$ + +Similarly, the liquid and solid canopy drip fluxes are + +$$q_{drip,\, liq} =\frac{W_{can,\,liq}^{intr} -W_{can,\,liq}^{max } }{\Delta t} \ge 0$$ + +$$q_{drip,\, ice} =\frac{W_{can,\,sno}^{intr} -W_{can,\,sno}^{max } }{\Delta t} \ge 0$$ + +where + +$$W_{can,liq}^{intr} =W_{can,liq}^{n} +q_{intr,\, liq} \Delta t\ge 0$$ + +and + +$$W_{can,sno}^{intr} =W_{can,sno}^{n} +q_{intr,\, ice} \Delta t\ge 0$$ + +are the the canopy liquid water and snow water equivalent after +accounting for interception, $W_{can,\,liq}^{n}$ and $W_{can,\,sno}^{n}$ +are the canopy liquid and snow water from the previous time step, and +$W_{can,\,liq}^{max }$ and $W_{can,\,snow}^{max }$ (kg m^-2^ or mm of +H~2~O) are the maximum amounts of liquid water and snow the canopy can +hold. They are defined by + +$$W_{can,\,liq}^{max } =p_{liq}\left(L+S\right)$$ + +$$W_{can,\,sno}^{max } =p_{sno}\left(L+S\right).$$ + +The maximum storage of liquid water is $p_{liq}=0.1$ kg m^-2^ +(`Dickinson et al. 1993 `{.interpreted-text +role="ref"}), and that of snow is $p_{sno}=6$ kg m^-2^, consistent with +reported field measurements +(`Pomeroy et al. 1998 `{.interpreted-text role="ref"}). + +Canopy snow unloading from wind speed $u$ and above-freezing +temperatures are modeled from linear fluxes and e-folding times similar +to `Roesch et al. (2001) `{.interpreted-text role="ref"} + +$$q_{unl,\, wind} =\frac{u W_{can,sno}}{1.56\times 10^5 \text{ m}}$$ + +$$q_{unl,\, temp} =\frac{W_{can,sno}(T-270 \textrm{ K})}{1.87\times 10^5 \text{ K s}} > 0$$ + +$$q_{unl,\, tot} =\min \left( q_{unl,\, wind} +q_{unl,\, temp} ,W_{can,\, sno} \right)$$ + +The canopy liquid water and snow water equivalent are updated as + +$$W_{can,\, liq}^{n+1} =W_{can,liq}^{n} + q_{intr,\, liq} - q_{drip,\, liq} \Delta t - E_{v}^{liq} \Delta t \ge 0$$ + +and + +$$W_{can,\, sno}^{n+1} =W_{can,sno}^{n} + q_{intr,\, ice} - \left(q_{drip,\, ice}+q_{unl,\, tot} \right)\Delta t + - E_{v}^{ice} \Delta t \ge 0$$ + +where $E_{v}^{liq}$ and $E_{v}^{ice}$ are partitioned from the stem and +leaf surface evaporation $E_{v}^{w}$ (Chapter +`rst_Momentum, Sensible Heat, and Latent Heat Fluxes`{.interpreted-text +role="numref"}) based on the vegetation temperature $T_{v}$ (K) (Chapter +`rst_Momentum, Sensible Heat, and Latent Heat Fluxes`{.interpreted-text +role="numref"}) and its relation to the freezing temperature of water +$T_{f}$ (K) (`Table Physical Constants`{.interpreted-text +role="numref"}) + +$$\begin{aligned} +E_{v}^{liq} = +\left\{\begin{array}{lr} +E_{v}^{w} & T_v > T_{f} \\ +0 & T_v \le T_f +\end{array}\right\} +\end{aligned}$$ + +$$\begin{aligned} +E_{v}^{ice} = +\left\{\begin{array}{lr} +0 & T_v > T_f \\ +E_{v}^{w} & T_v \le T_f +\end{array}\right\}. +\end{aligned}$$ + +The total rate of liquid and solid precipitation reaching the ground is +then + +$$q_{grnd,liq} =q_{thru,\, liq} +q_{drip,\, liq}$$ + +$$q_{grnd,ice} =q_{thru,\, ice} +q_{drip,\, ice} +q_{unl,\, tot} .$$ + +Solid precipitation reaching the soil or snow surface, +$q_{grnd,\, ice} \Delta t$, is added immediately to the snow pack +(Chapter `rst_Snow Hydrology`{.interpreted-text role="numref"}). The +liquid part, $q_{grnd,\, liq} \Delta t$ is added after surface fluxes +(Chapter +`rst_Momentum, Sensible Heat, and Latent Heat Fluxes`{.interpreted-text +role="numref"}) and snow/soil temperatures (Chapter +`rst_Soil and Snow Temperatures`{.interpreted-text role="numref"}) have +been determined. + +The wetted fraction of the canopy (stems plus leaves), which is required +for surface flux (Chapter +`rst_Momentum, Sensible Heat, and Latent Heat Fluxes`{.interpreted-text +role="numref"}) calculations, is +(`Dickinson et al.1993 `{.interpreted-text +role="ref"}) + +$$\begin{aligned} +f_{wet} = +\left\{\begin{array}{lr} +\left[\frac{W_{can} }{p_{liq}\left(L+S\right)} \right]^{{2\mathord{\left/ {\vphantom {2 3}} \right.} 3} } \le 1 & \qquad L+S > 0 \\ +0 &\qquad L+S = 0 +\end{array}\right\} +\end{aligned}$$ + +while the fraction of the canopy that is dry and transpiring is + +$$\begin{aligned} +f_{dry} = +\left\{\begin{array}{lr} +\frac{\left(1-f_{wet} \right)L}{L+S} & \qquad L+S > 0 \\ +0 &\qquad L+S = 0 +\end{array}\right\}. +\end{aligned}$$ + +Similarly, the snow-covered fraction of the canopy is used for surface +alebdo when intercepted snow is present (Chapter +`rst_Surface Albedos`{.interpreted-text role="numref"}) + +$$\begin{aligned} +f_{can,\, sno} = +\left\{\begin{array}{lr} +\left[\frac{W_{can,\, sno} }{p_{sno}\left(L+S\right)} \right]^{{3\mathord{\left/ {\vphantom {3 20}} \right.} 20} } \le 1 & \qquad L+S > 0 \\ +0 &\qquad L+S = 0 +\end{array}\right\}. +\end{aligned}$$ + +## Surface Runoff, Surface Water Storage, and Infiltration {#Surface Runoff, Surface Water Storage, and Infiltration} + +The moisture input at the grid cell surface,$q_{liq,\, 0}$, is the sum +of liquid precipitation reaching the ground and melt water from snow (kg +m^-2^ s^-1^). The moisture flux is then partitioned between surface +runoff, surface water storage, and infiltration into the soil. + +### Surface Runoff {#Surface Runoff} + +The simple TOPMODEL-based +(`Beven and Kirkby 1979 `{.interpreted-text +role="ref"}) runoff model (SIMTOP) described by +`Niu et al. (2005) `{.interpreted-text role="ref"} is +implemented to parameterize runoff. A key concept underlying this +approach is that of fractional saturated area $f_{sat}$, which is +determined by the topographic characteristics and soil moisture state of +a grid cell. The saturated portion of a grid cell contributes to surface +runoff, $q_{over}$, by the saturation excess mechanism (Dunne runoff) + +$$q_{over} =f_{sat} \ q_{liq,\, 0}$$ + +The fractional saturated area is a function of soil moisture + +$$f_{sat} =f_{\max } \ \exp \left(-0.5f_{over} z_{\nabla } \right)$$ + +where $f_{\max }$ is the potential or maximum value of $f_{sat}$, +$f_{over}$ is a decay factor (m^-1^), and $z_{\nabla}$ is the water +table depth (m) (section `Lateral Sub-surface Runoff`{.interpreted-text +role="numref"}). The maximum saturated fraction, $f_{\max }$, is defined +as the value of the discrete cumulative distribution function (CDF) of +the topographic index when the grid cell mean water table depth is zero. +Thus, $f_{\max }$ is the percent of pixels in a grid cell whose +topographic index is larger than or equal to the grid cell mean +topographic index. It should be calculated explicitly from the CDF at +each grid cell at the resolution that the model is run. However, because +this is a computationally intensive task for global applications, +$f_{\max }$ is calculated once at 0.125° resolution using the 1-km +compound topographic indices (CTIs) based on the HYDRO1K dataset +(`Verdin and Greenlee 1996 `{.interpreted-text +role="ref"}) from USGS following the algorithm in +`Niu et al. (2005) `{.interpreted-text role="ref"} and then +area-averaged to the desired model resolution (section +`Surface Data`{.interpreted-text role="numref"}). Pixels with CTIs +exceeding the 95 percentile threshold in each 0.125° grid cell are +excluded from the calculation to eliminate biased estimation of +statistics due to large CTI values at pixels on stream networks. For +grid cells over regions without CTIs such as Australia, the global mean +$f_{\max }$ is used to fill the gaps. See +`Li et al. (2013b) `{.interpreted-text role="ref"} for +additional details. The decay factor $f_{over}$ for global simulations +was determined through sensitivity analysis and comparison with observed +runoff to be 0.5 m^-1^. + +### Surface Water Storage {#Surface Water Storage} + +A surface water store has been added to the model to represent wetlands +and small, sub-grid scale water bodies. As a result, the wetland land +unit has been removed as of CLM4.5. The state variables for surface +water are the mass of water $W_{sfc}$ (kg m^-2^) and temperature +$T_{h2osfc}$ (Chapter `rst_Soil and Snow Temperatures`{.interpreted-text +role="numref"}). Surface water storage and outflow are functions of fine +spatial scale elevation variations called microtopography. The +microtopography is assumed to be distributed normally around the grid +cell mean elevation. Given the standard deviation of the +microtopographic distribution, $\sigma _{micro}$ (m), the fractional +area of the grid cell that is inundated can be calculated. Surface water +storage, $Wsfc$, is related to the height (relative to the grid cell +mean elevation) of the surface water, $d$, by + +$$W_{sfc} =\frac{d}{2} \left(1+erf\left(\frac{d}{\sigma _{micro} \sqrt{2} } \right)\right)+\frac{\sigma _{micro} }{\sqrt{2\pi } } e^{\frac{-d^{2} }{2\sigma _{micro} ^{2} } }$$ + +where $erf$ is the error function. For a given value of $W_{sfc}$, +`7.66`{.interpreted-text role="eq"} can be solved for $d$ using the +Newton-Raphson method. Once $d$ is known, one can determine the fraction +of the area that is inundated as + +$$f_{h2osfc} =\frac{1}{2} \left(1+erf\left(\frac{d}{\sigma _{micro} \sqrt{2} } \right)\right)$$ + +No global datasets exist for microtopography, so the default +parameterization is a simple function of slope + +$$\sigma _{micro} =\left(\beta +\beta _{0} \right)^{\eta }$$ + +where $\beta$ is the topographic slope, +$\beta_{0} =\left(\sigma_{\max } \right)^{\frac{1}{\eta } }$ determines +the maximum value of $\sigma_{micro}$, and $\eta$ is an adjustable +parameter. Default values in the model are $\sigma_{\max } =0.4$ and +$\eta =-3$. + +If the spatial scale of the microtopography is small relative to that of +the grid cell, one can assume that the inundated areas are distributed +randomly within the grid cell. With this assumption, a result from +percolation theory can be used to quantify the fraction of the inundated +portion of the grid cell that is interconnected + +$$\begin{aligned} +\begin{array}{lr} f_{connected} =\left(f_{h2osfc} -f_{c} \right)^{\mu } & \qquad f_{h2osfc} >f_{c} \\ f_{connected} =0 &\qquad f_{h2osfc} \le f_{c} \end{array} +\end{aligned}$$ + +where $f_{c}$ is a threshold below which no single connected inundated +area spans the grid cell and $\mu$ is a scaling exponent. Default values +of $f_{c}$ and $\mu$ are 0.4 and 0.14, respectively. When the inundated +fraction of the grid cell surpasses $f_{c}$, the surface water store +acts as a linear reservoir + +$$q_{out,h2osfc}=k_{h2osfc} \ f_{connected} \ (Wsfc-Wc)\frac{1}{\Delta t}$$ + +where $q_{out,h2osfc}$ is the surface water runoff, $k_{h2osfc}$ is a +constant, $Wc$ is the amount of surface water present when +$f_{h2osfc} =f_{c}$, and $\Delta t$ is the model time step. The linear +storage coefficent $k_{h2osfc} = \sin \left(\beta \right)$ is a function +of grid cell mean topographic slope where $\beta$ is the slope in +radians. + +### Infiltration {#Infiltration} + +The surface moisture flux remaining after surface runoff has been +removed, + +$$q_{in,surface} = (1-f_{sat}) \ q_{liq,\, 0}$$ + +is divided into inputs to surface water ($q_{in,\, h2osfc}$ ) and the +soil $q_{in,soil}$. If $q_{in,soil}$ exceeds the maximum soil +infiltration capacity (kg m^-2^ s^-1^), + +$$q_{infl,\, \max } =(1-f_{sat}) \ \Theta_{ice} k_{sat}$$ + +where $\Theta_{ice}$ is an ice impedance factor (section +`Hydraulic Properties`{.interpreted-text role="numref"}), infiltration +excess (Hortonian) runoff is generated + +$$q_{infl,\, excess} =\max \left(q_{in,soil} -\left(1-f_{h2osfc} \right)q_{\inf l,\max } ,0\right)$$ + +and transferred from $q_{in,soil}$ to $q_{in,h2osfc}$. After evaporative +losses have been removed, these moisture fluxes are + +$$q_{in,\, h2osfc} = f_{h2osfc} q_{in,surface} + q_{infl,excess} - q_{evap,h2osfc}$$ + +and + +$$q_{in,soil} = (1-f_{h2osfc} ) \ q_{in,surface} - q_{\inf l,excess} - (1 - f_{sno} - f_{h2osfc} ) \ q_{evap,soil}.$$ + +The balance of surface water is then calculated as + +$$\Delta W_{sfc} =\left(q_{in,h2osfc} - q_{out,h2osfc} - q_{drain,h2osfc} \right) \ \Delta t.$$ + +Bottom drainage from the surface water store + +$$q_{drain,h2osfc} = \min \left(f_{h2osfc} q_{\inf l,\max } ,\frac{W_{sfc} }{\Delta t} \right)$$ + +is then added to $q_{in,soil}$ giving the total infiltration into the +surface soil layer + +$$q_{infl} = q_{in,soil} + q_{drain,h2osfc}$$ + +Infiltration $q_{infl}$ and explicit surface runoff $q_{over}$ are not +allowed for glaciers. + +## Soil Water {#Soil Water} + +Soil water is predicted from a multi-layer model, in which the vertical +soil moisture transport is governed by infiltration, surface and +sub-surface runoff, gradient diffusion, gravity, and canopy +transpiration through root extraction +(`Figure Hydrologic processes`{.interpreted-text role="numref"}). + +For one-dimensional vertical water flow in soils, the conservation of +mass is stated as + +$$\frac{\partial \theta }{\partial t} =-\frac{\partial q}{\partial z} - e$$ + +where $\theta$ is the volumetric soil water content (mm^3^ of water / +mm^-3^ of soil), $t$ is time (s), $z$ is height above some datum in the +soil column (mm) (positive upwards), $q$ is soil water flux (kg m^-2^ +s^-1^ or mm s^-1^) (positive upwards), and $e$ is a soil moisture sink +term (mm of water mm^-1^ of soil s^-1^) (ET loss). This equation is +solved numerically by dividing the soil column into multiple layers in +the vertical and integrating downward over each layer with an upper +boundary condition of the infiltration flux into the top soil layer +$q_{infl}$ and a zero-flux lower boundary condition at the bottom of the +soil column (sub-surface runoff is removed later in the timestep, +section `Lateral Sub-surface Runoff`{.interpreted-text role="numref"}). + +The soil water flux $q$ in equation `7.79`{.interpreted-text role="eq"} +can be described by Darcy\'s law +`(Dingman 2002) `{.interpreted-text role="ref"} + +$$q = -k \frac{\partial \psi _{h} }{\partial z}$$ + +where $k$ is the hydraulic conductivity (mm s^-1^), and $\psi _{h}$ is +the hydraulic potential (mm). The hydraulic potential is + +$$\psi _{h} =\psi _{m} +\psi _{z}$$ + +where $\psi _{m}$ is the soil matric potential (mm) (which is related to +the adsorptive and capillary forces within the soil matrix), and +$\psi _{z}$ is the gravitational potential (mm) (the vertical distance +from an arbitrary reference elevation to a point in the soil). If the +reference elevation is the soil surface, then $\psi _{z} =z$. Letting +$\psi =\psi _{m}$, Darcy\'s law becomes + +$$q = -k \left[\frac{\partial \left(\psi +z\right)}{\partial z} \right].$$ + +Equation `7.82`{.interpreted-text role="eq"} can be further manipulated +to yield + +$$q = -k \left[\frac{\partial \left(\psi +z\right)}{\partial z} \right] += -k \left(\frac{\partial \psi }{\partial z} + 1 \right) \ .$$ + +Substitution of this equation into equation `7.79`{.interpreted-text +role="eq"}, with $e = 0$, yields the Richards equation +`(Dingman 2002) `{.interpreted-text role="ref"} + +$$\frac{\partial \theta }{\partial t} = +\frac{\partial }{\partial z} \left[k\left(\frac{\partial \psi }{\partial z} + 1 +\right)\right].$$ + +In practice (Section `Numerical Solution Hydrology`{.interpreted-text +role="numref"}), changes in soil water content are predicted from +`7.79`{.interpreted-text role="eq"} using finite-difference +approximations for `7.84`{.interpreted-text role="eq"}. + +### Hydraulic Properties {#Hydraulic Properties} + +The hydraulic conductivity $k_{i}$ (mm s^-1^) and the soil matric +potential $\psi _{i}$ (mm) for layer $i$ vary with volumetric soil water +$\theta _{i}$ and soil texture. As with the soil thermal properties +(section `Soil And Snow Thermal Properties`{.interpreted-text +role="numref"}) the hydraulic properties of the soil are assumed to be a +weighted combination of the mineral properties, which are determined +according to sand and clay contents based on work by +`Clapp and Hornberger (1978) `{.interpreted-text +role="ref"} and `Cosby et al. (1984) `{.interpreted-text +role="ref"}, and organic properties of the soil +(`Lawrence and Slater 2008 `{.interpreted-text +role="ref"}). + +The hydraulic conductivity is defined at the depth of the interface of +two adjacent layers $z_{h,\, i}$ +(`Figure Water flux schematic`{.interpreted-text role="numref"}) and is +a function of the saturated hydraulic conductivity +$k_{sat} \left[z_{h,\, i} \right]$, the liquid volumetric soil moisture +of the two layers $\theta _{i}$ and $\theta_{i+1}$ and an ice impedance +factor $\Theta_{ice}$ + +$$\begin{aligned} +k\left[z_{h,\, i} \right] = +\left\{\begin{array}{lr} +\Theta_{ice} k_{sat} \left[z_{h,\, i} \right]\left[\frac{0.5\left(\theta_{\, i} +\theta_{\, i+1} \right)}{0.5\left(\theta_{sat,\, i} +\theta_{sat,\, i+1} \right)} \right]^{2B_{i} +3} & \qquad 1 \le i \le N_{levsoi} - 1 \\ +\Theta_{ice} k_{sat} \left[z_{h,\, i} \right]\left(\frac{\theta_{\, i} }{\theta_{sat,\, i} } \right)^{2B_{i} +3} & \qquad i = N_{levsoi} +\end{array}\right\}. +\end{aligned}$$ + +The ice impedance factor is a function of ice content, and is meant to +quantify the increased tortuosity of the water flow when part of the +pore space is filled with ice. +`Swenson et al. (2012) `{.interpreted-text role="ref"} +used a power law form + +$$\Theta_{ice} = 10^{-\Omega F_{ice} }$$ + +where $\Omega = 6$and $F_{ice} = \frac{\theta_{ice} }{\theta_{sat} }$ is +the ice-filled fraction of the pore space. + +Because the hydraulic properties of mineral and organic soil may differ +significantly, the bulk hydraulic properties of each soil layer are +computed as weighted averages of the properties of the mineral and +organic components. The water content at saturation (i.e. porosity) is + +$$\theta_{sat,i} =(1-f_{om,i} )\theta_{sat,\min ,i} +f_{om,i} \theta_{sat,om}$$ + +where $f_{om,i}$ is the soil organic matter fraction, $\theta_{sat,om}$ +is the porosity of organic matter, and the porosity of the mineral soil +$\theta_{sat,\min,i}$ is + +$$\theta_{sat,\min ,i} = 0.489 - 0.00126(\% sand)_{i} .$$ + +The exponent $B_{i}$ is + +$$B_{i} =(1-f_{om,i} )B_{\min ,i} +f_{om,i} B_{om}$$ + +where $B_{om}$ is for organic matter and + +$$B_{\min ,i} =2.91+0.159(\% clay)_{i} .$$ + +The soil matric potential (mm) is defined at the node depth $z_{i}$ of +each layer $i$ (`Figure Water flux schematic`{.interpreted-text +role="numref"}) + +$$\psi _{i} =\psi _{sat,\, i} \left(\frac{\theta_{\, i} }{\theta_{sat,\, i} } \right)^{-B_{i} } \ge -1\times 10^{8} \qquad 0.01\le \frac{\theta_{i} }{\theta_{sat,\, i} } \le 1$$ + +where the saturated soil matric potential (mm) is + +$$\psi _{sat,i} =(1-f_{om,i} )\psi _{sat,\min ,i} +f_{om,i} \psi _{sat,om}$$ + +where $\psi _{sat,om}$ is the saturated organic matter matric potential +and the saturated mineral soil matric potential $\psi _{sat,\min,i}$ is + +$$\psi _{sat,\, \min ,\, i} =-10.0\times 10^{1.88-0.0131(\% sand)_{i} } .$$ + +The saturated hydraulic conductivity, $k_{sat} \left[z_{h,\, i} \right]$ +(mm s^-1^), for organic soils ($k_{sat,\, om}$ ) may be two to three +orders of magnitude larger than that of mineral soils +($k_{sat,\, \min }$ ). Bulk soil layer values of $k_{sat}$ calculated as +weighted averages based on $f_{om}$ may therefore be determined +primarily by the organic soil properties even for values of $f_{om}$ as +low as 1 %. To better represent the influence of organic soil material +on the grid cell average saturated hydraulic conductivity, the soil +organic matter fraction is further subdivided into \"connected\" and +\"unconnected\" fractions using a result from percolation theory +(`Stauffer and Aharony 1994 `{.interpreted-text +role="ref"}, +`Berkowitz and Balberg 1992 `{.interpreted-text +role="ref"}). Assuming that the organic and mineral fractions are +randomly distributed throughout a soil layer, percolation theory +predicts that above a threshold value $f_{om} = f_{threshold}$, +connected flow pathways consisting of organic material only exist and +span the soil space. Flow through these pathways interacts only with +organic material, and thus can be described by $k_{sat,\, om}$. This +fraction of the grid cell is given by + +$$\begin{aligned} +\begin{array}{lr} +f_{perc} =\; N_{perc} \left(f_{om} {\rm \; }-f_{threshold} \right)^{\beta_{perc} } f_{om} {\rm \; } & \qquad f_{om} \ge f_{threshold} \\ +f_{perc} = 0 & \qquad f_{om} `{.interpreted-text +role="ref"}) as + +$$k_{sat,\, \min } \left[z_{h,\, i} \right]=0.0070556\times 10^{-0.884+0.0153\left(\% sand\right)_{i} } .$$ + +The bulk soil layer saturated hydraulic conductivity is then computed as + +$$k_{sat} \left[z_{h,\, i} \right]=f_{uncon,\, i} k_{sat,\, uncon} \left[z_{h,\, i} \right]+(1-f_{uncon,\, i} )k_{sat,\, om} \left[z_{h,\, i} \right].$$ + +The soil organic matter properties implicitly account for the standard +observed profile of organic matter properties as + +$$\theta_{sat,om} = max(0.93 - 0.1\times z_{i} / zsapric, 0.83).$$ + +$$B_{om} = min(2.7 + 9.3\times z_{i} / zsapric, 12.0).$$ + +$$\psi_{sat,om} = min(10.3 - 0.2\times z_{i} / zsapric, 10.1).$$ + +$$k_{sat,om} = max(0.28 - 0.2799\times z_{i} / zsapric, k_{sat,\, \min } \left[z_{h,\, i} \right]).$$ + +where $zsapric =0.5$ m is the depth that organic matter takes on the +characteristics of sapric peat. + +### Numerical Solution {#Numerical Solution Hydrology} + +With reference to `Figure Water flux schematic`{.interpreted-text +role="numref"}, the equation for conservation of mass (equation +`7.79`{.interpreted-text role="eq"}) can be integrated over each layer +as + +$$\int _{-z_{h,\, i} }^{-z_{h,\, i-1} }\frac{\partial \theta }{\partial t} \, dz=-\int _{-z_{h,\, i} }^{-z_{h,\, i-1} }\frac{\partial q}{\partial z} \, dz-\int _{-z_{h,\, i} }^{-z_{h,\, i-1} } e\, dz .$$ + +Note that the integration limits are negative since $z$ is defined as +positive upward from the soil surface. This equation can be written as + +$$\Delta z_{i} \frac{\partial \theta_{liq,\, i} }{\partial t} =-q_{i-1} +q_{i} -e_{i}$$ + +where $q_{i}$ is the flux of water across interface $z_{h,\, i}$, +$q_{i-1}$ is the flux of water across interface $z_{h,\, i-1}$, and +$e_{i}$ is a layer-averaged soil moisture sink term (ET loss) defined as +positive for flow out of the layer (mm s^-1^). Taking the finite +difference with time and evaluating the fluxes implicitly at time $n+1$ +yields + +$$\frac{\Delta z_{i} \Delta \theta_{liq,\, i} }{\Delta t} =-q_{i-1}^{n+1} +q_{i}^{n+1} -e_{i}$$ + +where +$\Delta \theta_{liq,\, i} =\theta_{liq,\, i}^{n+1} -\theta_{liq,\, i}^{n}$ +is the change in volumetric soil liquid water of layer $i$ in time +$\Delta t$and $\Delta z_{i}$ is the thickness of layer $i$ (mm). + +The water removed by transpiration in each layer $e_{i}$ is a function +of the total transpiration $E_{v}^{t}$ (Chapter +`rst_Momentum, Sensible Heat, and Latent Heat Fluxes`{.interpreted-text +role="numref"}) and the effective root fraction $r_{e,\, i}$ + +$$e_{i} =r_{e,\, i} E_{v}^{t} .$$ + +::: {#Figure Water flux schematic} +![Schematic diagram of numerical scheme used to solve for soil water +fluxes.](image2.png) +::: + +Shown are three soil layers, $i-1$, $i$, and $i+1$. The soil matric +potential $\psi$ and volumetric soil water $\theta_{liq}$ are defined at +the layer node depth $z$. The hydraulic conductivity +$k\left[z_{h} \right]$ is defined at the interface of two layers +$z_{h}$. The layer thickness is $\Delta z$. The soil water fluxes +$q_{i-1}$ and $q_{i}$ are defined as positive upwards. The soil moisture +sink term $e$ (ET loss) is defined as positive for flow out of the +layer. + +Note that because more than one plant functional type (PFT) may share a +soil column, the transpiration $E_{v}^{t}$ is a weighted sum of +transpiration from all PFTs whose weighting depends on PFT area as + +$$E_{v}^{t} =\sum _{j=1}^{npft}\left(E_{v}^{t} \right)_{j} \left(wt\right)_{j}$$ + +where $npft$ is the number of PFTs sharing a soil column, +$\left(E_{v}^{t} \right)_{j}$ is the transpiration from the $j^{th}$ PFT +on the column, and $\left(wt\right)_{j}$ is the relative area of the +$j^{th}$ PFT with respect to the column. The effective root fraction +$r_{e,\, i}$ is also a column-level quantity that is a weighted sum over +all PFTs. The weighting depends on the per unit area transpiration of +each PFT and its relative area as + +$$r_{e,\, i} =\frac{\sum _{j=1}^{npft}\left(r_{e,\, i} \right)_{j} \left(E_{v}^{t} \right)_{j} \left(wt\right)_{j} }{\sum _{j=1}^{npft}\left(E_{v}^{t} \right)_{j} \left(wt\right)_{j} }$$ + +where $\left(r_{e,\, i} \right)_{j}$ is the effective root fraction for +the $j^{th}$ PFT + +$$\begin{aligned} +\begin{array}{lr} +\left(r_{e,\, i} \right)_{j} =\frac{\left(r_{i} \right)_{j} \left(w_{i} \right)_{j} }{\left(\beta _{t} \right)_{j} } & \qquad \left(\beta _{t} \right)_{j} >0 \\ +\left(r_{e,\, i} \right)_{j} =0 & \qquad \left(\beta _{t} \right)_{j} =0 +\end{array} +\end{aligned}$$ + +and $\left(r_{i} \right)_{j}$ is the fraction of roots in layer $i$ +(Chapter `rst_Stomatal Resistance and Photosynthesis`{.interpreted-text +role="numref"}), $\left(w_{i} \right)_{j}$ is a soil dryness or plant +wilting factor for layer $i$ (Chapter +`rst_Stomatal Resistance and Photosynthesis`{.interpreted-text +role="numref"}), and $\left(\beta_{t} \right)_{j}$ is a wetness factor +for the total soil column for the $j^{th}$ PFT (Chapter +`rst_Stomatal Resistance and Photosynthesis`{.interpreted-text +role="numref"}). + +The soil water fluxes in `7.103`{.interpreted-text role="eq"},, which +are a function of $\theta_{liq,\, i}$ and $\theta_{liq,\, i+1}$ because +of their dependence on hydraulic conductivity and soil matric potential, +can be linearized about $\theta$ using a Taylor series expansion as + +$$q_{i}^{n+1} =q_{i}^{n} +\frac{\partial q_{i} }{\partial \theta_{liq,\, i} } \Delta \theta_{liq,\, i} +\frac{\partial q_{i} }{\partial \theta_{liq,\, i+1} } \Delta \theta_{liq,\, i+1}$$ + +$$q_{i-1}^{n+1} =q_{i-1}^{n} +\frac{\partial q_{i-1} }{\partial \theta_{liq,\, i-1} } \Delta \theta_{liq,\, i-1} +\frac{\partial q_{i-1} }{\partial \theta_{liq,\, i} } \Delta \theta_{liq,\, i} .$$ + +Substitution of these expressions for $q_{i}^{n+1}$ and $q_{i-1}^{n+1}$ +into `7.103`{.interpreted-text role="eq"} results in a general +tridiagonal equation set of the form + +$$r_{i} =a_{i} \Delta \theta_{liq,\, i-1} +b_{i} \Delta \theta_{liq,\, i} +c_{i} \Delta \theta_{liq,\, i+1}$$ + +where + +$$a_{i} =-\frac{\partial q_{i-1} }{\partial \theta_{liq,\, i-1} }$$ + +$$b_{i} =\frac{\partial q_{i} }{\partial \theta_{liq,\, i} } -\frac{\partial q_{i-1} }{\partial \theta_{liq,\, i} } -\frac{\Delta z_{i} }{\Delta t}$$ + +$$c_{i} =\frac{\partial q_{i} }{\partial \theta_{liq,\, i+1} }$$ + +$$r_{i} =q_{i-1}^{n} -q_{i}^{n} +e_{i} .$$ + +The tridiagonal equation set is solved over $i=1,\ldots,N_{levsoi}$. + +The finite-difference forms of the fluxes and partial derivatives in +equations `7.111`{.interpreted-text role="eq"} - +`7.114`{.interpreted-text role="eq"} can be obtained from equation +`7.82`{.interpreted-text role="eq"} as + +$$q_{i-1}^{n} =-k\left[z_{h,\, i-1} \right]\left[\frac{\left(\psi _{i-1} -\psi _{i} \right)+\left(z_{i} - z_{i-1} \right)}{z_{i} -z_{i-1} } \right]$$ + +$$q_{i}^{n} =-k\left[z_{h,\, i} \right]\left[\frac{\left(\psi _{i} -\psi _{i+1} \right)+\left(z_{i+1} - z_{i} \right)}{z_{i+1} -z_{i} } \right]$$ + +$$\frac{\partial q_{i-1} }{\partial \theta _{liq,\, i-1} } =-\left[\frac{k\left[z_{h,\, i-1} \right]}{z_{i} -z_{i-1} } \frac{\partial \psi _{i-1} }{\partial \theta _{liq,\, i-1} } \right]-\frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta _{liq,\, i-1} } \left[\frac{\left(\psi _{i-1} -\psi _{i} \right)+\left(z_{i} - z_{i-1} \right)}{z_{i} - z_{i-1} } \right]$$ + +$$\frac{\partial q_{i-1} }{\partial \theta _{liq,\, i} } =\left[\frac{k\left[z_{h,\, i-1} \right]}{z_{i} -z_{i-1} } \frac{\partial \psi _{i} }{\partial \theta _{liq,\, i} } \right]-\frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta _{liq,\, i} } \left[\frac{\left(\psi _{i-1} -\psi _{i} \right)+\left(z_{i} - z_{i-1} \right)}{z_{i} - z_{i-1} } \right]$$ + +$$\frac{\partial q_{i} }{\partial \theta _{liq,\, i} } =-\left[\frac{k\left[z_{h,\, i} \right]}{z_{i+1} -z_{i} } \frac{\partial \psi _{i} }{\partial \theta _{liq,\, i} } \right]-\frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta _{liq,\, i} } \left[\frac{\left(\psi _{i} -\psi _{i+1} \right)+\left(z_{i+1} - z_{i} \right)}{z_{i+1} - z_{i} } \right]$$ + +$$\frac{\partial q_{i} }{\partial \theta _{liq,\, i+1} } =\left[\frac{k\left[z_{h,\, i} \right]}{z_{i+1} -z_{i} } \frac{\partial \psi _{i+1} }{\partial \theta _{liq,\, i+1} } \right]-\frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta _{liq,\, i+1} } \left[\frac{\left(\psi _{i} -\psi _{i+1} \right)+\left(z_{i+1} - z_{i} \right)}{z_{i+1} - z_{i} } \right].$$ + +The derivatives of the soil matric potential at the node depth are +derived from `7.94`{.interpreted-text role="eq"} + +$$\frac{\partial \psi _{i-1} }{\partial \theta_{liq,\, \, i-1} } =-B_{i-1} \frac{\psi _{i-1} }{\theta_{\, \, i-1} }$$ + +$$\frac{\partial \psi _{i} }{\partial \theta_{\, liq,\, i} } =-B_{i} \frac{\psi _{i} }{\theta_{i} }$$ + +$$\frac{\partial \psi _{i+1} }{\partial \theta_{liq,\, i+1} } =-B_{i+1} \frac{\psi _{i+1} }{\theta_{\, i+1} }$$ + +with the constraint +$0.01\, \theta_{sat,\, i} \le \theta_{\, i} \le \theta_{sat,\, i}$. + +The derivatives of the hydraulic conductivity at the layer interface are +derived from `7.85`{.interpreted-text role="eq"} + +$$\begin{array}{l} +{\frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta _{liq,\, i-1} } += \frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta _{liq,\, i} } += \left(2B_{i-1} +3\right) \ \overline{\Theta}_{ice} \ k_{sat} \left[z_{h,\, i-1} \right] \ \left[\frac{\overline{\theta}_{liq}}{\overline{\theta}_{sat}} \right]^{2B_{i-1} +2} \left(\frac{0.5}{\overline{\theta}_{sat}} \right)} \end{array}$$ + +where $\overline{\Theta}_{ice} = \Theta(\overline{\theta}_{ice})$ +`7.86`{.interpreted-text role="eq"}, +$\overline{\theta}_{ice} = 0.5\left(\theta_{ice\, i-1} +\theta_{ice\, i} \right)$, +$\overline{\theta}_{liq} = 0.5\left(\theta_{liq\, i-1} +\theta_{liq\, i} \right)$, +and +$\overline{\theta}_{sat} = 0.5\left(\theta_{sat,\, i-1} +\theta_{sat,\, i} \right)$ + +and + +$$\begin{array}{l} +{\frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta _{liq,\, i} } += \frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta _{liq,\, i+1} } += \left(2B_{i} +3\right) \ \overline{\Theta}_{ice} \ k_{sat} \left[z_{h,\, i} \right] \ \left[\frac{\overline{\theta}_{liq}}{\overline{\theta}_{sat}} \right]^{2B_{i} +2} \left(\frac{0.5}{\overline{\theta}_{sat}} \right)} \end{array}.$$ + +where +$\overline{\theta}_{liq} = 0.5\left(\theta_{\, i} +\theta_{\, i+1} \right)$, +$\overline{\theta}_{sat} = 0.5\left(\theta_{sat,\, i} +\theta_{sat,\, i+1} \right)$. + +#### Equation set for layer $i=1$ + +For the top soil layer ($i=1$), the boundary condition is the +infiltration rate (section `Surface Runoff`{.interpreted-text +role="numref"}), $q_{i-1}^{n+1} =-q_{infl}^{n+1}$, and the water balance +equation is + +$$\frac{\Delta z_{i} \Delta \theta_{liq,\, i} }{\Delta t} =q_{infl}^{n+1} +q_{i}^{n+1} -e_{i} .$$ + +After grouping like terms, the coefficients of the tridiagonal set of +equations for $i=1$ are + +$$a_{i} =0$$ + +$$b_{i} =\frac{\partial q_{i} }{\partial \theta_{liq,\, i} } -\frac{\Delta z_{i} }{\Delta t}$$ + +$$c_{i} =\frac{\partial q_{i} }{\partial \theta_{liq,\, i+1} }$$ + +$$r_{i} =q_{infl}^{n+1} -q_{i}^{n} +e_{i} .$$ + +#### Equation set for layers $i=2,\ldots ,N_{levsoi} -1$ + +The coefficients of the tridiagonal set of equations for +$i=2,\ldots,N_{levsoi} -1$ are + +$$a_{i} =-\frac{\partial q_{i-1} }{\partial \theta_{liq,\, i-1} }$$ + +$$b_{i} =\frac{\partial q_{i} }{\partial \theta_{liq,\, i} } -\frac{\partial q_{i-1} }{\partial \theta_{liq,\, i} } -\frac{\Delta z_{i} }{\Delta t}$$ + +$$c_{i} =\frac{\partial q_{i} }{\partial \theta_{liq,\, i+1} }$$ + +$$r_{i} =q_{i-1}^{n} -q_{i}^{n} +e_{i} .$$ + +#### Equation set for layer $i=N_{levsoi}$ + +For the lowest soil layer ($i=N_{levsoi}$ ), a zero-flux bottom boundary +condition is applied ($q_{i}^{n} =0$) and the coefficients of the +tridiagonal set of equations for $i=N_{levsoi}$ are + +$$a_{i} =-\frac{\partial q_{i-1} }{\partial \theta_{liq,\, i-1} }$$ + +$$b_{i} =\frac{\partial q_{i} }{\partial \theta_{liq,\, i} } -\frac{\partial q_{i-1} }{\partial \theta_{liq,\, i} } -\frac{\Delta z_{i} }{\Delta t}$$ + +$$c_{i} =0$$ + +$$r_{i} =q_{i-1}^{n} +e_{i} .$$ + +#### Adaptive Time Stepping + +The length of the time step is adjusted in order to improve the accuracy +and stability of the numerical solutions. The difference between two +numerical approximations is used to estimate the temporal truncation +error, and then the step size $\Delta t_{sub}$ is adjusted to meet a +user-prescribed error tolerance +`[Kavetski et al., 2002]`{.interpreted-text +role="ref"}. The temporal truncation error is estimated by comparing the +flux obtained from the first-order Taylor series expansion +($q_{i-1}^{n+1}$ and $q_{i}^{n+1}$, equations `7.108`{.interpreted-text +role="eq"} and `7.109`{.interpreted-text role="eq"}) against the flux at +the start of the time step ($q_{i-1}^{n}$ and $q_{i}^{n}$). Since the +tridiagonal solution already provides an estimate of +$\Delta \theta_{liq,i}$, it is convenient to compute the error for each +of the $i$ layers from equation `7.103`{.interpreted-text role="eq"} as + +$$\epsilon_{i} = \left[ \frac{\Delta \theta_{liq,\, i} \Delta z_{i}}{\Delta t_{sub}} - +\left( q_{i-1}^{n} - q_{i}^{n} + e_{i}\right) \right] \ \frac{\Delta t_{sub}}{2}$$ + +and the maximum absolute error across all layers as + +$$\begin{array}{lr} +\epsilon_{crit} = {\rm max} \left( \left| \epsilon_{i} \right| \right) & \qquad 1 \le i \le nlevsoi +\end{array} \ .$$ + +The adaptive step size selection is based on specified upper and lower +error tolerances, $\tau_{U}$ and $\tau_{L}$. The solution is accepted if +$\epsilon_{crit} \le \tau_{U}$ and the procedure repeats until the +adaptive sub-stepping spans the full model time step (the sub-steps are +doubled if $\epsilon_{crit} \le \tau_{L}$, i.e., if the solution is very +accurate). Conversely, the solution is rejected if +$\epsilon_{crit} > \tau_{U}$. In this case the length of the sub-steps +is halved and a new solution is obtained. The halving of substeps +continues until either $\epsilon_{crit} \le \tau_{U}$ or the specified +minimum time step length is reached. + +Upon solution of the tridiagonal equation set, the liquid water contents +are updated as follows + +$$w_{liq,\, i}^{n+1} =w_{liq,\, i}^{n} +\Delta \theta_{liq,\, i} \Delta z_{i} \qquad i=1,\ldots ,N_{levsoi} .$$ + +The volumetric water content is + +$$\theta_{i} =\frac{w_{liq,\, i} }{\Delta z_{i} \rho _{liq} } +\frac{w_{ice,\, i} }{\Delta z_{i} \rho _{ice} } .$$ + +## Frozen Soils and Perched Water Table {#Frozen Soils and Perched Water Table} + +When soils freeze, the power-law form of the ice impedance factor +(section `Hydraulic Properties`{.interpreted-text role="numref"}) can +greatly decrease the hydraulic conductivity of the soil, leading to +nearly impermeable soil layers. When unfrozen soil layers are present +above relatively ice-rich frozen layers, the possibility exists for +perched saturated zones. Lateral drainage from perched saturated regions +is parameterized as a function of the thickness of the saturated zone + +$$q_{drai,perch} =k_{drai,\, perch} \left(z_{frost} -z_{\nabla ,perch} \right)$$ + +where $k_{drai,\, perch}$ depends on topographic slope and soil +hydraulic conductivity, + +$$k_{drai,\, perch} =10^{-5} \sin (\beta )\left(\frac{\sum _{i=N_{perch} }^{i=N_{frost} }\Theta_{ice,i} k_{sat} \left[z_{i} \right]\Delta z_{i} }{\sum _{i=N_{perch} }^{i=N_{frost} }\Delta z_{i} } \right)$$ + +where $\Theta_{ice}$ is an ice impedance factor, $\beta$ is the mean +grid cell topographic slope in radians, $z_{frost}$ is the depth to the +frost table, and $z_{\nabla,perch}$ is the depth to the perched +saturated zone. The frost table $z_{frost}$ is defined as the shallowest +frozen layer having an unfrozen layer above it, while the perched water +table $z_{\nabla,perch}$ is defined as the depth at which the volumetric +water content drops below a specified threshold. The default threshold +is set to 0.9. Drainage from the perched saturated zone $q_{drai,perch}$ +is removed from layers $N_{perch}$ through $N_{frost}$, which are the +layers containing $z_{\nabla,perch}$ and, $z_{frost}$ respectively. + +## Lateral Sub-surface Runoff {#Lateral Sub-surface Runoff} + +Lateral sub-surface runoff occurs when saturated soil moisture +conditions exist within the soil column. Sub-surface runoff is + +$$q_{drai} = \Theta_{ice} K_{baseflow} tan \left( \beta \right) +\Delta z_{sat}^{N_{baseflow}} \ ,$$ + +where $K_{baseflow}$ is a calibration parameter, $\beta$ is the +topographic slope, the exponent $N_{baseflow}$ = 1, and $\Delta z_{sat}$ +is the thickness of the saturated portion of the soil column. + +The saturated thickness is + +$$\Delta z_{sat} = z_{bedrock} - z_{\nabla},$$ + +where the water table $z_{\nabla}$ is determined by finding the first +soil layer above the bedrock depth (section +`Depth to Bedrock`{.interpreted-text role="numref"}) in which the +volumetric water content drops below a specified threshold. The default +threshold is set to 0.9. + +The specific yield, $S_{y}$, which depends on the soil properties and +the water table location, is derived by taking the difference between +two equilibrium soil moisture profiles whose water tables differ by an +infinitesimal amount + +$$S_{y} =\theta_{sat} \left(1-\left(1+\frac{z_{\nabla } }{\Psi _{sat} } \right)^{\frac{-1}{B} } \right)$$ + +where B is the Clapp-Hornberger exponent. Because $S_{y}$ is a function +of the soil properties, it results in water table dynamics that are +consistent with the soil water fluxes described in section +`Soil Water`{.interpreted-text role="numref"}. + +After the above calculations, two numerical adjustments are implemented +to keep the liquid water content of each soil layer ($w_{liq,\, i}$ ) +within physical constraints of +$w_{liq}^{\min } \le w_{liq,\, i} \le \left(\theta_{sat,\, i} -\theta_{ice,\, i} \right)\Delta z_{i}$ +where $w_{liq}^{\min } =0.01$ (mm). First, beginning with the bottom +soil layer $i=N_{levsoi}$, any excess liquid water in each soil layer +($w_{liq,\, i}^{excess} =w_{liq,\, i} -\left(\theta_{sat,\, i} -\theta_{ice,\, i} \right)\Delta z_{i} \ge 0$) +is successively added to the layer above. Any excess liquid water that +remains after saturating the entire soil column is added to drainage +$q_{drai}$. Second, to prevent negative $w_{liq,\, i}$, each layer is +successively brought up to $w_{liq,\, i} =w_{liq}^{\min }$ by taking the +required amount of water from the layer below. If this results in +$w_{liq,\, N_{levsoi} } >>>>>> Converted from CLM5 diff --git a/components/elm/docs/tech-guide/index.md b/components/elm/docs/tech-guide/index.md index 6ce636548204..feef6d521697 100644 --- a/components/elm/docs/tech-guide/index.md +++ b/components/elm/docs/tech-guide/index.md @@ -7,3 +7,5 @@ Shortwave radiation model - [TOP Parameterization](top_solar_parameterization.md): Parameterization of sub-grid topographical effects on solar radiation. - [Longwave Radiation](longwave_radiation.md): Longwave radiation model +- [Hydrology](hydrology.md): Hydrology model +- [Plant Hydraulics](phs.md): Plant Hydraulics model diff --git a/components/elm/docs/tech-guide/phs.md b/components/elm/docs/tech-guide/phs.md new file mode 100644 index 000000000000..1f2fc254998f --- /dev/null +++ b/components/elm/docs/tech-guide/phs.md @@ -0,0 +1,679 @@ +<<<<<<< HEAD + +======= +# Plant Hydraulics + +## Roots {#Roots} + +### Vertical Root Distribution {#Vertical Root Distribution} + +The root fraction $r_{i}$ in each soil layer depends on the plant +functional type + +$$r_{i} = +\begin{array}{lr} +\left(\beta^{z_{h,\, i-1} \cdot 100} - \beta^{z_{h,\, i} \cdot 100} \right) & \qquad {\rm for\; }1 \le i \le N_{levsoi} +\end{array}$$ + +where $z_{h,\, i}$ (m) is the depth from the soil surface to the +interface between layers $i$ and $i+1$ ($z_{h,\, 0}$, the soil surface) +(section `Vertical Discretization`{.interpreted-text role="numref"}), +the factor of 100 converts from m to cm, and $\beta$ is a +plant-dependent root distribution parameter adopted from +`Jackson et al. (1996)`{.interpreted-text role="ref"} +(`Table Plant functional type root distribution parameters`{.interpreted-text +role="numref"}). + +**Table 1: Plant functional type root distribution parameters ** + +| Plant Functional Type | $\beta$ | +|----------------------|---------| +| NET Temperate | 0.976 | +| NET Boreal | 0.943 | +| NDT Boreal | 0.943 | +| BET Tropical | 0.993 | +| BET temperate | 0.966 | +| BDT tropical | 0.993 | +| BDT temperate | 0.966 | +| BDT boreal | 0.943 | +| BES temperate | 0.964 | +| BDS temperate | 0.964 | +| BDS boreal | 0.914 | +| C₃ grass arctic | 0.914 | +| C₃ grass | 0.943 | +| C₄ grass | 0.943 | +| Crop R | 0.943 | +| Crop I | 0.943 | +| Corn R | 0.943 | +| Corn I | 0.943 | +| Temp Cereal R | 0.943 | +| Temp Cereal I | 0.943 | +| Winter Cereal R | 0.943 | +| Winter Cereal I | 0.943 | +| Soybean R | 0.943 | +| Soybean I | 0.943 | +| Miscanthus R | 0.943 | +| Miscanthus I | 0.943 | +| Switchgrass R | 0.943 | +| Switchgrass I | 0.943 | + +::: {#Table Plant functional type root distribution parameters} + ----------------------------------------------------- + Plant Functional Type $\beta$ + ---------------------------------- ------------------ + NET Temperate 0.976 + + NET Boreal 0.943 + + NDT Boreal 0.943 + + BET Tropical 0.993 + + BET temperate 0.966 + + BDT tropical 0.993 + + BDT temperate 0.966 + + BDT boreal 0.943 + + BES temperate 0.964 + + BDS temperate 0.964 + + BDS boreal 0.914 + + C~3~ grass arctic 0.914 + + C~3~ grass 0.943 + + C~4~ grass 0.943 + + Crop R 0.943 + + Crop I 0.943 + + Corn R 0.943 + + Corn I 0.943 + + Temp Cereal R 0.943 + + Temp Cereal I 0.943 + + Winter Cereal R 0.943 + + Winter Cereal I 0.943 + + Soybean R 0.943 + + Soybean I 0.943 + + Miscanthus R 0.943 + + Miscanthus I 0.943 + + Switchgrass R 0.943 + + Switchgrass I 0.943 + ----------------------------------------------------- + + : Plant functional type root distribution parameters +::: + +### Root Spacing {#Root Spacing} + +To determine the conductance along the soil to root pathway (section +`Soil-to-root`{.interpreted-text role="numref"}) an estimate of the +spacing between the roots within a soil layer is required. The distance +between roots $dx_{root,i}$ (m) is calculated by assuming that roots are +distributed uniformly throughout the soil +(`Gardner 1960`{.interpreted-text role="ref"}) + +$$dx_{root,i} = \left(\pi \cdot L_i\right)^{- \frac{1}{2}}$$ + +where $L_{i}$ is the root length density (m m ^-3^) + +$$L_{i} = \frac{B_{root,i}}{\rho_{root} {CA}_{root}} \ ,$$ + +$B_{root,i}$ is the root biomass density (kg m ^-3^) + +$$B_{root,i} = \frac{c\_to\_b \cdot C_{fineroot} \cdot r_{i}}{dz_{i}}$$ + +where $c\_to\_b = 2$ (kg biomass kg carbon ^-1^) and $C_{fineroot}$ is +the amount of fine root carbon (kg m ^-2^). + +$\rho_{root}$ is the root density (kg m ^-3^), and ${CA}_{root}$ is the +fine root cross sectional area (m ^2^) + +$$CA_{root} = \pi r_{root}^{2}$$ + +where $r_{root}$ is the root radius (m). + +## Plant Hydraulic Stress {#Plant Hydraulic Stress} + +The Plant Hydraulic Stress (PHS) routine explicitly models water +transport through the vegetation according to a simple hydraulic +framework following Darcy\'s Law for porous media flow equations +influenced by `Bonan et al. (2014) `{.interpreted-text +role="ref"}, `Chuang et al. (2006) `{.interpreted-text +role="ref"}, `Sperry et al. (1998) `{.interpreted-text +role="ref"}, +`Sperry and Love (2015) `{.interpreted-text +role="ref"}, +`Williams et al (1996) `{.interpreted-text +role="ref"}. + +PHS solves for the vegetation water potential that matches water supply +with transpiration demand. Water supply is modeled according to the +circuit analog in `Figure Plant hydraulic circuit`{.interpreted-text +role="numref"}. Transpiration demand is modeled relative to maximum +transpiration by a transpiration loss function dependent on leaf water +potential. + +::: {#Figure Plant hydraulic circuit} +![Circuit diagram of plant hydraulics scheme](circuit.jpg) +::: + +### Plant Water Supply {#Plant Water Supply} + +The supply equations are used to solve for vegetation water potential +forced by transpiration demand and the set of layer-by-layer soil water +potentials. The water supply is discretized into segments: soil-to-root, +root-to-stem, and stem-to-leaf. There are typically several (1-49) +soil-to-root flows operating in parallel, one per soil layer. There are +two stem-to-leaf flows operating in parallel, corresponding to the +sunlit and shaded \"leaves\". + +In general the water fluxes (e.g. soil-to-root, root-to-stem, etc.) are +modeled according to Darcy\'s Law for porous media flow as: + +$$q = kA\left( \psi_1 - \psi_2 \right)$$ + +$q$ is the flux of water (mmH~2~O/s) spanning the segment between +$\psi_1$ and $\psi_2$ + +$k$ is the hydraulic conductance (s^-1^) + +$A$ is the area basis (m^2^/m^2^) relating the conducting area basis to +ground area $\psi_1 - \psi_2$ is the gradient in water potential +(mmH~2~O) across the segment The segments in +`Figure Plant hydraulic circuit`{.interpreted-text role="numref"} have +variable resistance, as the water potentials become lower, hydraulic +conductance decreases. This is captured by multiplying the maximum +segment conductance by a sigmoidal function capturing the percent loss +of conductivity. The function uses two parameters to fit experimental +vulnerability curves: the water potential at 50% loss of conductivity +($p50$) and a shape fitting parameter ($c_k$). + +$$k=k_{max}\cdot 2^{-\left(\dfrac{\psi_1}{p50}\right)^{c_k}}$$ + +$k_{max}$ is the maximum segment conductance (s^-1^) + +$p50$ is the water potential at 50% loss of conductivity (mmH~2~O) + +$\psi_1$ is the water potential of the lower segment terminus (mmH~2~O) + +#### Stem-to-leaf {#Stem-to-leaf} + +The area basis and conductance parameterization varies by segment. There +are two stem-to-leaf fluxes in parallel, from stem to sunlit leaf and +from stem to shaded leaf ($q_{1a}$ and $q_{1a}$). The water flux from +stem-to-leaf is the product of the segment conductance, the conducting +area basis, and the water potential gradient from stem to leaf. +Stem-to-leaf conductance is defined as the maximum conductance +multiplied by the percent of maximum conductance, as calculated by the +sigmoidal vulnerability curve. The maximum conductance is a PFT +parameter representing the maximum conductance of water from stem to +leaf per unit leaf area. This parameter can be defined separately for +sunlit and shaded segments and should already include the appropriate +length scaling (in other words this is a conductance, not conductivity). +The water potential gradient is the difference between leaf water +potential and stem water potential. There is no gravity term, assuming a +negligible difference in height across the segment. The area basis is +the leaf area index (either sunlit or shaded). + +$$q_{1a}=k_{1a}\cdot\mbox{LAI}_{sun}\cdot\left(\psi_{stem}-\psi_{sunleaf} \right)$$ + +$$q_{1b}=k_{1b}\cdot\mbox{LAI}_{shade}\cdot\left(\psi_{stem}-\psi_{shadeleaf} \right)$$ + +$$k_{1a}=k_{1a,max}\cdot 2^{-\left(\dfrac{\psi_{stem}}{p50_1}\right)^{c_k}}$$ + +$$k_{1b}=k_{1b,max}\cdot 2^{-\left(\dfrac{\psi_{stem}}{p50_1}\right)^{c_k}}$$ + +Variables: + +$q_{1a}$ = flux of water (mmH~2~O/s) from stem to sunlit leaf + +$q_{1b}$ = flux of water (mmH~2~O/s) from stem to shaded leaf + +$LAI_{sun}$ = sunlit leaf area index (m2/m2) + +$LAI_{shade}$ = shaded leaf area index (m2/m2) + +$\psi_{stem}$ = stem water potential (mmH~2~O) + +$\psi_{sunleaf}$ = sunlit leaf water potential (mmH~2~O) + +$\psi_{shadeleaf}$ = shaded leaf water potential (mmH~2~O) + +Parameters: + +$k_{1a,max}$ = maximum leaf conductance (s^-1^) + +$k_{1b,max}$ = maximum leaf conductance (s^-1^) + +$p50_{1}$ = water potential at 50% loss of conductance (mmH~2~O) + +$c_{k}$ = vulnerability curve shape-fitting parameter (-) + +#### Root-to-stem {#Root-to-stem} + +There is one root-to-stem flux. This represents a flux from the root +collar to the upper branch reaches. The water flux from root-to-stem is +the product of the segment conductance, the conducting area basis, and +the water potential gradient from root to stem. Root-to-stem conductance +is defined as the maximum conductance multiplied by the percent of +maximum conductance, as calculated by the sigmoidal vulnerability curve +(two parameters). The maximum conductance is defined as the maximum +root-to-stem conductivity per unit stem area (PFT parameter) divided by +the length of the conducting path, which is taken to be the vegetation +height. The area basis is the stem area index. The gradient in water +potential is the difference between the root water potential and the +stem water potential less the difference in gravitational potential. + +$$q_2=k_2 \cdot SAI \cdot \left( \psi_{root} - \psi_{stem} - \Delta \psi_z \right)$$ + +$$k_2=\dfrac{k_{2,max}}{z_2} \cdot 2^{-\left(\dfrac{\psi_{root}}{p50_2}\right)^{c_k}}$$ + +Variables: + +$q_2$ = flux of water (mmH~2~O/s) from root to stem + +$SAI$ = stem area index (m2/m2) + +$\Delta\psi_z$ = gravitational potential (mmH~2~O) + +$\psi_{root}$ = root water potential (mmH~2~O) + +$\psi_{stem}$ = stem water potential (mmH~2~O) + +Parameters: + +$k_{2,max}$ = maximum stem conductivity (m/s) + +$p50_2$ = water potential at 50% loss of conductivity (mmH~2~O) + +$z_2$ = vegetation height (m) + +#### Soil-to-root {#Soil-to-root} + +There are several soil-to-root fluxes operating in parallel (one for +each root-containing soil layer). Each represents a flux from the given +soil layer to the root collar. The water flux from soil-to-root is the +product of the segment conductance, the conducting area basis, and the +water potential gradient from soil to root. The area basis is a proxy +for root area index, defined as the summed leaf and stem area index +multiplied by the root-to-shoot ratio (PFT parameter) multiplied by the +layer root fraction. The root fraction comes from an empirical root +profile (section `Vertical Root Distribution`{.interpreted-text +role="numref"}). + +The gradient in water potential is the difference between the soil water +potential and the root water potential less the difference in +gravitational potential. There is only one root water potential to which +all soil layers are connected in parallel. A soil-to-root flux can be +either positive (vegetation water uptake) or negative (water +deposition), depending on the relative values of the root and soil water +potentials. This allows for the occurrence of hydraulic redistribution +where water moves through vegetation tissue from one soil layer to +another. + +Soil-to-root conductance is the result of two resistances in series, +first across the soil-root interface and then through the root tissue. +The root tissue conductance is defined as the maximum conductance +multiplied by the percent of maximum conductance, as calculated by the +sigmoidal vulnerability curve. The maximum conductance is defined as the +maximum root-tissue conductivity (PFT parameter) divided by the length +of the conducting path, which is taken to be the soil layer depth plus +lateral root length. + +The soil-root interface conductance is defined as the soil conductivity +divided by the conducting length from soil to root. The soil +conductivity varies by soil layer and is calculated based on soil +potential and soil properties, via the Brooks-Corey theory. The +conducting length is determined from the characteristic root spacing +(section `Root Spacing`{.interpreted-text role="numref"}). + +$$q_{3,i}=k_{3,i} \cdot \left(\psi_{soil,i}-\psi_{root} + \Delta\psi_{z,i} \right)$$ + +$$k_{3,i}=\dfrac{k_{r,i} \cdot k_{s,i}}{k_{r,i}+k_{s,i}}$$ + +$$k_{r,i}=\dfrac{k_{3,max}}{z_{3,i}} \cdot RAI \cdot 2^{-\left(\dfrac{\psi_{soil,i}}{p50_3}\right)^{c_k}}$$ + +$$RAI=\left(LAI+SAI \right) \cdot r_i \cdot f_{root-leaf}$$ + +$$k_{s,i} = \dfrac{k_{soil,i}}{dx_{root,i}}$$ + +Variables: + +$q_{3,i}$ = flux of water (mmH~2~O/s) from soil layer $i$ to root + +$\Delta\psi_{z,i}$ = change in gravitational potential from soil layer +$i$ to surface (mmH~2~O) + +$LAI$ = total leaf area index (m2/m2) + +$SAI$ = stem area index (m2/m2) + +$\psi_{soil,i}$ = water potential in soil layer $i$ (mmH~2~O) + +$\psi_{root}$ = root water potential (mmH~2~O) + +$z_{3,i}$ = length of root tissue conducting path = soil layer depth + +root lateral length (m) + +$r_i$ = root fraction in soil layer $i$ (-) + +$k_{soil,i}$ = Brooks-Corey soil conductivity in soil layer $i$ (m/s) + +Parameters: + +$f_{root-leaf}$ = root-to-shoot ratio (-) + +$p50_3$ = water potential at 50% loss of root tissue conductance +(mmH~2~O) + +$ck$ = shape-fitting parameter for vulnerability curve (-) + +### Plant Water Demand {#Plant Water Demand} + +Plant water demand depends on stomatal conductance, which is described +in section `Stomatal resistance`{.interpreted-text role="numref"}. Here +we describe the influence of PHS and the coupling of vegetation water +demand and supply. PHS models vegetation water demand as transpiration +attenuated by a transpiration loss function based on leaf water +potential. Sunlit leaf transpiration is modeled as the maximum sunlit +leaf transpiration multiplied by the percent of maximum transpiration as +modeled by the sigmoidal loss function. The same follows for shaded leaf +transpiration. Maximum stomatal conductance is calculated from the +Medlyn model `(Medlyn et al. 2011) `{.interpreted-text +role="ref"} absent water stress and used to calculate the maximum +transpiration (see section +`Sensible and Latent Heat Fluxes and Temperature for Vegetated Surfaces`{.interpreted-text +role="numref"}). Water stress is calculated as the ratio of attenuated +stomatal conductance to maximum stomatal conductance. Water stress is +calculated with distinct values for sunlit and shaded leaves. Vegetation +water stress is calculated based on leaf water potential and is used to +attenuate photosynthesis (see section `Photosynthesis`{.interpreted-text +role="numref"}) + +$$E_{sun} = E_{sun,max} \cdot 2^{-\left(\dfrac{\psi_{sunleaf}}{p50_e}\right)^{c_k}}$$ + +$$E_{shade} = E_{shade,max} \cdot 2^{-\left(\dfrac{\psi_{shadeleaf}}{p50_e}\right)^{c_k}}$$ + +$$\beta_{t,sun} = \dfrac{g_{s,sun}}{g_{s,sun,\beta_t=1}}$$ + +$$\beta_{t,shade} = \dfrac{g_{s,shade}}{g_{s,shade,\beta_t=1}}$$ + +$E_{sun}$ = sunlit leaf transpiration (mm/s) + +$E_{shade}$ = shaded leaf transpiration (mm/s) + +$E_{sun,max}$ = sunlit leaf transpiration absent water stress (mm/s) + +$E_{shade,max}$ = shaded leaf transpiration absent water stress (mm/s) + +$\psi_{sunleaf}$ = sunlit leaf water potential (mmH~2~O) + +$\psi_{shadeleaf}$ = shaded leaf water potential (mmH~2~O) + +$\beta_{t,sun}$ = sunlit transpiration water stress (-) + +$\beta_{t,shade}$ = shaded transpiration water stress (-) + +$g_{s,sun}$ = stomatal conductance of water corresponding to $E_{sun}$ + +$g_{s,shade}$ = stomatal conductance of water corresponding to +$E_{shade}$ + +$g_{s,sun,max}$ = stomatal conductance of water corresponding to +$E_{sun,max}$ + +$g_{s,shade,max}$ = stomatal conductance of water corresponding to +$E_{shade,max}$ + +### Vegetation Water Potential {#Vegetation Water Potential} + +Both plant water supply and demand are functions of vegetation water +potential. PHS explicitly models root, stem, shaded leaf, and sunlit +leaf water potential at each timestep. PHS iterates to find the +vegetation water potential $\psi$ (vector) that satisfies continuity +between the non-linear vegetation water supply and demand (equations +`11.103`{.interpreted-text role="eq"}, `11.104`{.interpreted-text +role="eq"}, `11.107`{.interpreted-text role="eq"}, +`11.109`{.interpreted-text role="eq"}, `11.201`{.interpreted-text +role="eq"}, `11.202`{.interpreted-text role="eq"}). + +$$\psi=\left[\psi_{sunleaf},\psi_{shadeleaf},\psi_{stem},\psi_{root}\right]$$ + +$$\begin{aligned} +\begin{aligned} +E_{sun}&=q_{1a}\\ +E_{shade}&=q_{1b}\\ +E_{sun}+E_{shade}&=q_{1a}+q_{1b}\\ +&=q_2\\ +&=\sum_{i=1}^{nlevsoi}{q_{3,i}} +\end{aligned} +\end{aligned}$$ + +PHS finds the water potentials that match supply and demand. In the +plant water transport equations `11.302`{.interpreted-text role="eq"}, +the demand terms (left-hand side) are decreasing functions of absolute +leaf water potential. As absolute leaf water potential becomes larger, +water stress increases, causing a decrease in transpiration demand. The +supply terms (right-hand side) are increasing functions of absolute leaf +water potential. As absolute leaf water potential becomes larger, the +gradients in water potential increase, causing an increase in vegetation +water supply. PHS takes a Newton\'s method approach to iteratively solve +for the vegetation water potentials that satisfy continuity +`11.302`{.interpreted-text role="eq"}. + +### Numerical Implementation {#PHS Numerical Implementation} + +The four plant water potential nodes are ( $\psi_{root}$, +$\psi_{xylem}$, $\psi_{shadeleaf}$, $\psi_{sunleaf}$). The fluxes +between each pair of nodes are labeled in Figure 1. $E_{sun}$ and +$E_{sha}$ are the transpiration from sunlit and shaded leaves, +respectively. We use the circuit-analog model to calculate the +vegetation water potential ( $\psi$) for the four plant nodes, forced by +soil matric potential and unstressed transpiration. The unstressed +transpiration is acquired by running the photosynthesis model with +$\beta_t=1$. The unstressed transpiration flux is attenuated based on +the leaf-level vegetation water potential. Using the attenuated +transpiration, we solve for $g_{s,stressed}$ and output +$\beta_t=\dfrac{g_{s,stressed}}{g_{s,unstressed}}$. + +The continuity of water flow through the system yields four equations + +$$\begin{aligned} +\begin{aligned} +E_{sun}&=q_{1a}\\ +E_{shade}&=q_{1b}\\ +q_{1a}+q_{1b}&=q_2\\ +q_2&=\sum_{i=1}^{nlevsoi}{q_{3,i}} +\end{aligned} +\end{aligned}$$ + +We seek the set of vegetation water potential values, + +$$\psi=\left[ \begin {array}{c} +\psi_{sunleaf}\cr\psi_{shadeleaf}\cr\psi_{stem}\cr\psi_{root} +\end {array} \right]$$ + +that satisfies these equations, as forced by the soil moisture and +atmospheric state. Each flux on the schematic can be represented in +terms of the relevant water potentials. Defining the transpiration +fluxes: + +$$\begin{aligned} +\begin{aligned} +E_{sun} &= E_{sun,max} \cdot 2^{-\left(\dfrac{\psi_{sunleaf}}{p50_e}\right)^{c_k}} \\ +E_{shade} &= E_{shade,max} \cdot 2^{-\left(\dfrac{\psi_{shadeleaf}}{p50_e}\right)^{c_k}} +\end{aligned} +\end{aligned}$$ + +Defining the water supply fluxes: + +$$\begin{aligned} +\begin{aligned} +q_{1a}&=k_{1a,max}\cdot 2^{-\left(\dfrac{\psi_{stem}}{p50_1}\right)^{c_k}} \cdot\mbox{LAI}_{sun}\cdot\left(\psi_{stem}-\psi_{sunleaf} \right) \\ +q_{1b}&=k_{1b,max}\cdot 2^{-\left(\dfrac{\psi_{stem}}{p50_1}\right)^{c_k}}\cdot\mbox{LAI}_{shade}\cdot\left(\psi_{stem}-\psi_{shadeleaf} \right) \\ +q_2&=\dfrac{k_{2,max}}{z_2} \cdot 2^{-\left(\dfrac{\psi_{root}}{p50_2}\right)^{c_k}} \cdot SAI \cdot \left( \psi_{root} - \psi_{stem} - \Delta \psi_z \right) \\ +q_{soil}&=\sum_{i=1}^{nlevsoi}{q_{3,i}}=\sum_{i=1}^{nlevsoi}{k_{3,i}\cdot RAI\cdot\left(\psi_{soil,i}-\psi_{root} + \Delta\psi_{z,i} \right)} +\end{aligned} +\end{aligned}$$ + +We\'re looking to find the vector $\psi$ that fits with soil and +atmospheric forcings while satisfying water flow continuity. Due to the +model non-linearity, we use a linearized explicit approach, iterating +with Newton\'s method. The initial guess is the solution for $\psi$ +(vector) from the previous time step. The general framework, from +iteration [m]{.title-ref} to [m+1]{.title-ref} is: + +$$\begin{aligned} +q^{m+1}=q^m+\dfrac{\delta q}{\delta\psi}\Delta\psi \\ +\psi^{m+1}=\psi^{m}+\Delta\psi +\end{aligned}$$ + +So for our first flux balance equation, at iteration [m+1]{.title-ref}, +we have: + +$$E_{sun}^{m+1}=q_{1a}^{m+1}$$ + +Which can be linearized to: + +$$E_{sun}^{m}+\dfrac{\delta E_{sun}}{\delta\psi}\Delta\psi=q_{1a}^{m}+\dfrac{\delta q_{1a}}{\delta\psi}\Delta\psi$$ + +And rearranged to be: + +$$\dfrac{\delta q_{1a}}{\delta\psi}\Delta\psi-\dfrac{\delta E_{sun}}{\delta\psi}\Delta\psi=E_{sun}^{m}-q_{1a}^{m}$$ + +And for the other 3 flux balance equations: + +$$\begin{aligned} +\begin{aligned} +\dfrac{\delta q_{1b}}{\delta\psi}\Delta\psi-\dfrac{\delta E_{sha}}{\delta\psi}\Delta\psi&=E_{sha}^{m}-q_{1b}^{m} \\ +\dfrac{\delta q_2}{\delta\psi}\Delta\psi-\dfrac{\delta q_{1a}}{\delta\psi}\Delta\psi-\dfrac{\delta q_{1b}}{\delta\psi}\Delta\psi&=q_{1a}^{m}+q_{1b}^{m}-q_2^{m} \\ +\dfrac{\delta q_{soil}}{\delta\psi}\Delta\psi-\dfrac{\delta q_2}{\delta\psi}\Delta\psi&=q_2^{m}-q_{soil}^{m} +\end{aligned} +\end{aligned}$$ + +Putting all four together in matrix form: + +$$\left[ \begin {array}{c} +\dfrac{\delta q_{1a}}{\delta\psi}-\dfrac{\delta E_{sun}}{\delta\psi} \cr +\dfrac{\delta q_{1b}}{\delta\psi}-\dfrac{\delta E_{sha}}{\delta\psi} \cr +\dfrac{\delta q_2}{\delta\psi}-\dfrac{\delta q_{1a}}{\delta\psi}-\dfrac{\delta q_{1b}}{\delta\psi} \cr +\dfrac{\delta q_{soil}}{\delta\psi}-\dfrac{\delta q_2}{\delta\psi} +\end {array} \right] +\Delta\psi= +\left[ \begin {array}{c} +E_{sun}^{m}-q_{1a}^{m} \cr +E_{sha}^{m}-q_{1b}^{m} \cr +q_{1a}^{m}+q_{1b}^{m}-q_2^{m} \cr +q_2^{m}-q_{soil}^{m} +\end {array} \right]$$ + +Now to expand the left-hand side, from generic $\psi$ to all four plant +water potential nodes, noting that many derivatives are zero (e.g. +$\dfrac{\delta E_{sun}}{\delta\psi_{sha}}=0$) + +Introducing the notation: $A\Delta\psi=b$ + +$$\Delta\psi=\left[ \begin {array}{c} +\Delta\psi_{sunleaf} \cr +\Delta\psi_{shadeleaf} \cr +\Delta\psi_{stem} \cr +\Delta\psi_{root} +\end {array} \right]$$ + +$$A= +\left[ \begin {array}{cccc} +\dfrac{\delta q_{1a}}{\delta \psi_{sun}}-\dfrac{\delta E_{sun}}{\delta \psi_{sun}}&0&\dfrac{\delta q_{1a}}{\delta \psi_{stem}}&0\cr +0&\dfrac{\delta q_{1b}}{\delta \psi_{sha}}-\dfrac{\delta E_{sha}}{\delta \psi_{sha}}&\dfrac{\delta q_{1b}}{\delta \psi_{stem}}&0\cr +-\dfrac{\delta q_{1a}}{\delta \psi_{sun}}& +-\dfrac{\delta q_{1b}}{\delta \psi_{sha}}& +\dfrac{\delta q_2}{\delta \psi_{stem}}-\dfrac{\delta q_{1a}}{\delta \psi_{stem}}-\dfrac{\delta q_{1b}}{\delta \psi_{stem}}& +\dfrac{\delta q_2}{\delta \psi_{root}}\cr +0&0&-\dfrac{\delta q_2}{\delta \psi_{stem}}&\dfrac{\delta q_{soil}}{\delta \psi_{root}}-\dfrac{\delta q_2}{\delta \psi_{root}} +\end {array} \right]$$ + +$$b= +\left[ \begin {array}{c} +E_{sun}^{m}-q_{b1}^{m} \cr +E_{sha}^{m}-q_{b2}^{m} \cr +q_{b1}^{m}+q_{b2}^{m}-q_{stem}^{m} \cr +q_{stem}^{m}-q_{soil}^{m} +\end {array} \right]$$ + +Now we compute all the entries for $A$ and $b$ based on the soil +moisture and maximum transpiration forcings and can solve to find: + +$$\Delta\psi=A^{-1}b$$ + +$$\psi_{m+1}=\psi_m+\Delta\psi$$ + +We iterate until $b\to 0$, signifying water flux balance through the +system. The result is a final set of water potentials ( $\psi_{root}$, +$\psi_{xylem}$, $\psi_{shadeleaf}$, $\psi_{sunleaf}$) satisfying +non-divergent water flux through the system. The magnitude of the water +flux is driven by soil matric potential and unstressed ( $\beta_t=1$) +transpiration. + +We use the transpiration solution (corresponding to the final solution +for $\psi$) to compute stomatal conductance. The stomatal conductance is +then used to compute $\beta_t$. + +$$\beta_{t,sun} = \dfrac{g_{s,sun}}{g_{s,sun,\beta_t=1}}$$ + +$$\beta_{t,shade} = \dfrac{g_{s,shade}}{g_{s,shade,\beta_t=1}}$$ + +The $\beta_t$ values are used in the Photosynthesis module (see section +`Photosynthesis`{.interpreted-text role="numref"}) to apply water +stress. The solution for $\psi$ is saved as a new variable (vegetation +water potential) and is indicative of plant water status. The +soil-to-root fluxes $\left( q_{3,1},q_{3,2},\mbox{...},q_{3,n}\right)$ +are used as the soil transpiration sink in the Richards\' equation +subsurface flow equations (see section `Soil Water`{.interpreted-text +role="numref"}). + +### Flow Diagram of Leaf Flux Calculations: {#Flow Diagram of Leaf Flux Calculations} + +PHS runs nested in the loop that solves for sensible and latent heat +fluxes and temperature for vegetated surfaces (see section +`Sensible and Latent Heat Fluxes and Temperature for Vegetated Surfaces`{.interpreted-text +role="numref"}). The scheme iterates for convergence of leaf temperature +($T_l$), transpiration water stress ($\beta_t$), and intercellular CO2 +concentration ($c_i$). PHS is forced by maximum transpiration (absent +water stress, $\beta_t=1$), whereby we first solve for assimilation, +stomatal conductance, and intercellular CO2 with $\beta_{t,sun}$ and +$\beta_{t,shade}$ both set to 1. This involves iterating to convergence +of $c_i$ (see section `Photosynthesis`{.interpreted-text +role="numref"}). + +Next, using the solutions for $E_{sun,max}$ and $E_{shade,max}$, PHS +solves for $\psi$, $\beta_{t,sun}$, and $\beta_{t,shade}$. The values +for $\beta_{t,sun}$, and $\beta_{t,shade}$ are inputs to the +photosynthesis routine, which now solves for attenuated photosynthesis +and stomatal conductance (reflecting water stress). Again this involves +iterating to convergence of $c_i$. Non-linearities between $\beta_t$ and +transpiration require also iterating to convergence of $\beta_t$. The +outermost level of iteration works towards convergence of leaf +temperature, reflecting leaf surface energy balance. + +::: {#Figure PHS Flow Diagram} +![Flow diagram of leaf flux calculations](phs_iteration_schematic.*) +::: + +>>>>>>> Converted from CLM5 From ffbf84ecabdb765a6a2f454f3cf36370c83e027c Mon Sep 17 00:00:00 2001 From: Yilin Fang Date: Thu, 5 Jun 2025 16:32:38 -0700 Subject: [PATCH 2/9] Modifications --- components/elm/docs/tech-guide/hydrology.md | 4 +-- components/elm/docs/tech-guide/phs.md | 39 +-------------------- 2 files changed, 2 insertions(+), 41 deletions(-) diff --git a/components/elm/docs/tech-guide/hydrology.md b/components/elm/docs/tech-guide/hydrology.md index a8044e365e0c..be160e2046f4 100644 --- a/components/elm/docs/tech-guide/hydrology.md +++ b/components/elm/docs/tech-guide/hydrology.md @@ -1,6 +1,5 @@ +# Hydrology {#rst_Hydrology} -<<<<<<< HEAD -======= The model parameterizes interception, throughfall, canopy drip, snow accumulation and melt, water transfer between snow layers, infiltration, evaporation, surface runoff, sub-surface drainage, redistribution within @@ -958,4 +957,3 @@ and lakes, which reduces the total amount of runoff available to the river routing model (Chapter `rst_MOSART`{.interpreted-text role="numref"}). ->>>>>>> Converted from CLM5 diff --git a/components/elm/docs/tech-guide/phs.md b/components/elm/docs/tech-guide/phs.md index 1f2fc254998f..539ae1722341 100644 --- a/components/elm/docs/tech-guide/phs.md +++ b/components/elm/docs/tech-guide/phs.md @@ -1,7 +1,4 @@ -<<<<<<< HEAD - -======= -# Plant Hydraulics +# Plant Hydraulics {#rst_Plant Hydraulics} ## Roots {#Roots} @@ -24,39 +21,6 @@ plant-dependent root distribution parameter adopted from (`Table Plant functional type root distribution parameters`{.interpreted-text role="numref"}). -**Table 1: Plant functional type root distribution parameters ** - -| Plant Functional Type | $\beta$ | -|----------------------|---------| -| NET Temperate | 0.976 | -| NET Boreal | 0.943 | -| NDT Boreal | 0.943 | -| BET Tropical | 0.993 | -| BET temperate | 0.966 | -| BDT tropical | 0.993 | -| BDT temperate | 0.966 | -| BDT boreal | 0.943 | -| BES temperate | 0.964 | -| BDS temperate | 0.964 | -| BDS boreal | 0.914 | -| C₃ grass arctic | 0.914 | -| C₃ grass | 0.943 | -| C₄ grass | 0.943 | -| Crop R | 0.943 | -| Crop I | 0.943 | -| Corn R | 0.943 | -| Corn I | 0.943 | -| Temp Cereal R | 0.943 | -| Temp Cereal I | 0.943 | -| Winter Cereal R | 0.943 | -| Winter Cereal I | 0.943 | -| Soybean R | 0.943 | -| Soybean I | 0.943 | -| Miscanthus R | 0.943 | -| Miscanthus I | 0.943 | -| Switchgrass R | 0.943 | -| Switchgrass I | 0.943 | - ::: {#Table Plant functional type root distribution parameters} ----------------------------------------------------- Plant Functional Type $\beta$ @@ -676,4 +640,3 @@ temperature, reflecting leaf surface energy balance. ![Flow diagram of leaf flux calculations](phs_iteration_schematic.*) ::: ->>>>>>> Converted from CLM5 From 9c201ea92311d37845a721767fda4250cbc4f213 Mon Sep 17 00:00:00 2001 From: Yilin Fang Date: Thu, 5 Jun 2025 16:34:40 -0700 Subject: [PATCH 3/9] Modifications --- components/elm/docs/tech-guide/hydrology.md | 2 +- components/elm/docs/tech-guide/phs.md | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/components/elm/docs/tech-guide/hydrology.md b/components/elm/docs/tech-guide/hydrology.md index be160e2046f4..60c63c93991d 100644 --- a/components/elm/docs/tech-guide/hydrology.md +++ b/components/elm/docs/tech-guide/hydrology.md @@ -1,4 +1,4 @@ -# Hydrology {#rst_Hydrology} +# Hydrology The model parameterizes interception, throughfall, canopy drip, snow accumulation and melt, water transfer between snow layers, infiltration, diff --git a/components/elm/docs/tech-guide/phs.md b/components/elm/docs/tech-guide/phs.md index 539ae1722341..3c6db2af327e 100644 --- a/components/elm/docs/tech-guide/phs.md +++ b/components/elm/docs/tech-guide/phs.md @@ -1,4 +1,4 @@ -# Plant Hydraulics {#rst_Plant Hydraulics} +# Plant Hydraulics ## Roots {#Roots} From b449406f9051ea2e8b83546c78ea958470985be2 Mon Sep 17 00:00:00 2001 From: Yilin Fang Date: Fri, 6 Jun 2025 15:20:05 -0700 Subject: [PATCH 4/9] correction --- components/elm/docs/tech-guide/hydrology.md | 176 ++++++++++---------- 1 file changed, 89 insertions(+), 87 deletions(-) diff --git a/components/elm/docs/tech-guide/hydrology.md b/components/elm/docs/tech-guide/hydrology.md index 60c63c93991d..7fe785da36db 100644 --- a/components/elm/docs/tech-guide/hydrology.md +++ b/components/elm/docs/tech-guide/hydrology.md @@ -8,14 +8,15 @@ changes in canopy water $\Delta W_{can,\,liq}$, canopy snow water $\Delta W_{can,\,sno}$ surface water $\Delta W_{sfc}$, snow water $\Delta W_{sno}$, soil water $\Delta w_{liq,\, i}$, and soil ice $\Delta w_{ice,\, i}$, and water in the unconfined aquifer -$\Delta W_{a}$ (all in kg m^-2^ or mm of H~2~O) +$\Delta W_{a}$ (all in $kg\ m^{-2}$ or $mm\ of\ H_2O$) (`Figure Hydrologic processes`{.interpreted-text role="numref"}). The total water balance of the system is $$\begin{aligned} -\begin{array}{l} {\Delta W_{can,\,liq} +\Delta W_{can,\,sno} +\Delta W_{sfc} +\Delta W_{sno} +} \\ {\sum _{i=1}^{N_{levsoi} }\left(\Delta w_{liq,\, i} +\Delta w_{ice,\, i} \right)+\Delta W_{a} =\left(\begin{array}{l} {q_{rain} +q_{sno} -E_{v} -E_{g} -q_{over} } \\ {-q_{h2osfc} -q_{drai} -q_{rgwl} -q_{snwcp,\, ice} } \end{array}\right) \Delta t} \end{array} +{\Delta W_{can,\,liq} +\Delta W_{can,\,sno} +\Delta W_{sfc} +\Delta W_{sno} +} \\ {\sum _{i=1}^{N_{levsoi} }\left(\Delta w_{liq,\, i} +\Delta w_{ice,\, i} \right)+\Delta W_{a} =\left(\begin{array}{l} {q_{rain} +q_{sno} -E_{v} -E_{g} -q_{over} } \\ {-q_{h2osfc} -q_{drai} -q_{rgwl} -q_{snwcp,\, ice} } \end{array}\right) \Delta t} \end{aligned}$$ +

(1)

where $q_{rain}$ is the liquid part of precipitation, $q_{sno}$ is the solid part of precipitation, $E_{v}$ is ET from vegetation (Chapter @@ -32,7 +33,7 @@ $q_{rgwl}$ and $q_{snwcp,ice}$ are liquid and solid runoff from glaciers and lakes, and runoff from other surface types due to snow capping (section `Runoff from glaciers and snow-capped surfaces`{.interpreted-text -role="numref"}) (all in kg m^-2^ s^-1^), $N_{levsoi}$ is the number of +role="numref"}) (all in $kg\ m^{-2} s^{-1}$), $N_{levsoi}$ is the number of soil layers (note that hydrology calculations are only done over soil layers 1 to $N_{levsoi}$; ground levels $N_{levsoi} +1$ to $N_{levgrnd}$ are currently hydrologically inactive; @@ -43,14 +44,14 @@ role="ref"} and $\Delta t$ is the time step (s). ![Hydrologic processes represented in CLM.](hydrologic.processes.png) ::: -## Canopy Water {#Canopy Water} +## Canopy Water Liquid precipitation is either intercepted by the canopy, falls directly to the snow/soil surface (throughfall), or drips off the vegetation (canopy drip). Solid precipitation is treated similarly, with the addition of unloading of previously intercepted snow. Interception by vegetation is divided between liquid and solid phases $q_{intr,\,liq}$ -and $q_{intr,\,ice}$ (kg m^-2^ s^-1^) +and $q_{intr,\,ice}$ ($kg\ m^{-2} s^{-1}$) $$q_{intr,\,liq} = f_{pi,\,liq} \ q_{rain}$$ @@ -61,7 +62,7 @@ precipitation of rain and snow, respectively $$f_{pi,\,liq} = \alpha_{liq} \ tanh \left(L+S\right)$$ -$$f_{pi,\,ice} =\alpha_{sno} \ \left\{1-\exp \left[-0.5\left(L+S\right)\right]\right\} \ ,$$ +$$f_{pi,\,ice} =\alpha_{sno} \ \left\(1-\exp \left[-0.5\left(L+S\right)\right]\right\) $$ and $L$ and $S$ are the exposed leaf and stem area index, respectively (section `Phenology and vegetation burial by snow`{.interpreted-text @@ -69,7 +70,7 @@ role="numref"}), and the $\alpha$\'s scale the fractional area of a leaf that collects water (`Lawrence et al. 2007 `{.interpreted-text role="ref"}). Default values of $\alpha_{liq}$ and $\alpha_{sno}$ are -set to 1. Throughfall (kg m^-2^ s^-1^) is also divided into liquid and +set to 1. Throughfall ($kg\ m^{-2} s^{-1}$) is also divided into liquid and solid phases, reaching the ground (soil or snow surface) as $$q_{thru,\, liq} = q_{rain} \left(1 - f_{pi,\,liq}\right)$$ @@ -93,17 +94,17 @@ $$W_{can,sno}^{intr} =W_{can,sno}^{n} +q_{intr,\, ice} \Delta t\ge 0$$ are the the canopy liquid water and snow water equivalent after accounting for interception, $W_{can,\,liq}^{n}$ and $W_{can,\,sno}^{n}$ are the canopy liquid and snow water from the previous time step, and -$W_{can,\,liq}^{max }$ and $W_{can,\,snow}^{max }$ (kg m^-2^ or mm of -H~2~O) are the maximum amounts of liquid water and snow the canopy can +$W_{can,\,liq}^{max }$ and $W_{can,\,snow}^{max }$ ($kg\ m^{-2}$ or $mm\ of\ +H_2O) are the maximum amounts of liquid water and snow the canopy can hold. They are defined by $$W_{can,\,liq}^{max } =p_{liq}\left(L+S\right)$$ $$W_{can,\,sno}^{max } =p_{sno}\left(L+S\right).$$ -The maximum storage of liquid water is $p_{liq}=0.1$ kg m^-2^ +The maximum storage of liquid water is $p_{liq}=0.1$ $kg\ m^{-2}$ (`Dickinson et al. 1993 `{.interpreted-text -role="ref"}), and that of snow is $p_{sno}=6$ kg m^-2^, consistent with +role="ref"}), and that of snow is $p_{sno}=6$ $kg\ m^{-2}$, consistent with reported field measurements (`Pomeroy et al. 1998 `{.interpreted-text role="ref"}). @@ -137,18 +138,18 @@ role="numref"}) $$\begin{aligned} E_{v}^{liq} = -\left\{\begin{array}{lr} +\left\(\begin{array}{lr} E_{v}^{w} & T_v > T_{f} \\ 0 & T_v \le T_f -\end{array}\right\} +\end{array}\right\) \end{aligned}$$ $$\begin{aligned} E_{v}^{ice} = -\left\{\begin{array}{lr} +\left\(\begin{array}{lr} 0 & T_v > T_f \\ E_{v}^{w} & T_v \le T_f -\end{array}\right\}. +\end{array}\right\). \end{aligned}$$ The total rate of liquid and solid precipitation reaching the ground is @@ -175,44 +176,45 @@ role="numref"}) calculations, is (`Dickinson et al.1993 `{.interpreted-text role="ref"}) -$$\begin{aligned} -f_{wet} = -\left\{\begin{array}{lr} -\left[\frac{W_{can} }{p_{liq}\left(L+S\right)} \right]^{{2\mathord{\left/ {\vphantom {2 3}} \right.} 3} } \le 1 & \qquad L+S > 0 \\ -0 &\qquad L+S = 0 -\end{array}\right\} -\end{aligned}$$ +$$ +f_{\text{wet}} = +\begin{cases} +\left[\dfrac{W_{\text{can}}}{\rho_{\text{liq}}(L + S)} \right]^{2/3}, & \text{if } L + S > 0 \text{ and } \left[\dfrac{W_{\text{can}}}{\rho_{\text{liq}}(L + S)} \right]^{2/3} \leq 1 \\ +1, & \text{if } L + S > 0 \text{ and } \left[\dfrac{W_{\text{can}}}{\rho_{\text{liq}}(L + S)} \right]^{2/3} > 1 \\ +0, & \text{if } L + S = 0 +\end{cases} +$$ while the fraction of the canopy that is dry and transpiring is -$$\begin{aligned} -f_{dry} = -\left\{\begin{array}{lr} -\frac{\left(1-f_{wet} \right)L}{L+S} & \qquad L+S > 0 \\ -0 &\qquad L+S = 0 -\end{array}\right\}. -\end{aligned}$$ +$$ +f_{\text{dry}} = +\begin{cases} +\dfrac{(1 - f_{\text{wet}}) L}{L + S}, & \text{if } L + S > 0 \\ +0, & \text{if } L + S = 0 +\end{cases} +$$ Similarly, the snow-covered fraction of the canopy is used for surface alebdo when intercepted snow is present (Chapter `rst_Surface Albedos`{.interpreted-text role="numref"}) -$$\begin{aligned} -f_{can,\, sno} = -\left\{\begin{array}{lr} -\left[\frac{W_{can,\, sno} }{p_{sno}\left(L+S\right)} \right]^{{3\mathord{\left/ {\vphantom {3 20}} \right.} 20} } \le 1 & \qquad L+S > 0 \\ -0 &\qquad L+S = 0 -\end{array}\right\}. -\end{aligned}$$ +$$ +f_{\text{can,\,sno}} = +\begin{cases} +\min\left( \left[\dfrac{W_{\text{can,\,sno}}}{\rho_{\text{sno}}(L + S)} \right]^{3/20},\ 1 \right), & \text{if } L + S > 0 \\ +0, & \text{if } L + S = 0 +\end{cases} +$$ -## Surface Runoff, Surface Water Storage, and Infiltration {#Surface Runoff, Surface Water Storage, and Infiltration} +## Surface Runoff, Surface Water Storage, and Infiltration The moisture input at the grid cell surface,$q_{liq,\, 0}$, is the sum -of liquid precipitation reaching the ground and melt water from snow (kg -m^-2^ s^-1^). The moisture flux is then partitioned between surface +of liquid precipitation reaching the ground and melt water from snow ($kg\ +m^{-2} s^{-1}$). The moisture flux is then partitioned between surface runoff, surface water storage, and infiltration into the soil. -### Surface Runoff {#Surface Runoff} +### Surface Runoff The simple TOPMODEL-based (`Beven and Kirkby 1979 `{.interpreted-text @@ -231,7 +233,7 @@ The fractional saturated area is a function of soil moisture $$f_{sat} =f_{\max } \ \exp \left(-0.5f_{over} z_{\nabla } \right)$$ where $f_{\max }$ is the potential or maximum value of $f_{sat}$, -$f_{over}$ is a decay factor (m^-1^), and $z_{\nabla}$ is the water +$f_{over}$ is a decay factor ($m^{-1}$), and $z_{\nabla}$ is the water table depth (m) (section `Lateral Sub-surface Runoff`{.interpreted-text role="numref"}). The maximum saturated fraction, $f_{\max }$, is defined as the value of the discrete cumulative distribution function (CDF) of @@ -256,14 +258,14 @@ $f_{\max }$ is used to fill the gaps. See `Li et al. (2013b) `{.interpreted-text role="ref"} for additional details. The decay factor $f_{over}$ for global simulations was determined through sensitivity analysis and comparison with observed -runoff to be 0.5 m^-1^. +runoff to be 0.5 $m^{-1}$. -### Surface Water Storage {#Surface Water Storage} +### Surface Water Storage A surface water store has been added to the model to represent wetlands and small, sub-grid scale water bodies. As a result, the wetland land unit has been removed as of CLM4.5. The state variables for surface -water are the mass of water $W_{sfc}$ (kg m^-2^) and temperature +water are the mass of water $W_{sfc}$ ($kg\ m^{-2}$) and temperature $T_{h2osfc}$ (Chapter `rst_Soil and Snow Temperatures`{.interpreted-text role="numref"}). Surface water storage and outflow are functions of fine spatial scale elevation variations called microtopography. The @@ -319,7 +321,7 @@ storage coefficent $k_{h2osfc} = \sin \left(\beta \right)$ is a function of grid cell mean topographic slope where $\beta$ is the slope in radians. -### Infiltration {#Infiltration} +### Infiltration The surface moisture flux remaining after surface runoff has been removed, @@ -328,7 +330,7 @@ $$q_{in,surface} = (1-f_{sat}) \ q_{liq,\, 0}$$ is divided into inputs to surface water ($q_{in,\, h2osfc}$ ) and the soil $q_{in,soil}$. If $q_{in,soil}$ exceeds the maximum soil -infiltration capacity (kg m^-2^ s^-1^), +infiltration capacity ($kg\ m^{-2} s^{-1}$), $$q_{infl,\, \max } =(1-f_{sat}) \ \Theta_{ice} k_{sat}$$ @@ -363,7 +365,7 @@ $$q_{infl} = q_{in,soil} + q_{drain,h2osfc}$$ Infiltration $q_{infl}$ and explicit surface runoff $q_{over}$ are not allowed for glaciers. -## Soil Water {#Soil Water} +## Soil Water Soil water is predicted from a multi-layer model, in which the vertical soil moisture transport is governed by infiltration, surface and @@ -378,9 +380,9 @@ $$\frac{\partial \theta }{\partial t} =-\frac{\partial q}{\partial z} - e$$ where $\theta$ is the volumetric soil water content (mm^3^ of water / mm^-3^ of soil), $t$ is time (s), $z$ is height above some datum in the -soil column (mm) (positive upwards), $q$ is soil water flux (kg m^-2^ -s^-1^ or mm s^-1^) (positive upwards), and $e$ is a soil moisture sink -term (mm of water mm^-1^ of soil s^-1^) (ET loss). This equation is +soil column (mm) (positive upwards), $q$ is soil water flux ($kg\ m^{-2} +s^{-1}$ or $mm\ s^{-1}$) (positive upwards), and $e$ is a soil moisture sink +term ($mm\ of\ water mm^{-1}\ of\ soil\ s^{-1}$) (ET loss). This equation is solved numerically by dividing the soil column into multiple layers in the vertical and integrating downward over each layer with an upper boundary condition of the infiltration flux into the top soil layer @@ -392,19 +394,19 @@ The soil water flux $q$ in equation `7.79`{.interpreted-text role="eq"} can be described by Darcy\'s law `(Dingman 2002) `{.interpreted-text role="ref"} -$$q = -k \frac{\partial \psi _{h} }{\partial z}$$ +$$q = -k \frac{\partial \psi_{h} }{\partial z}$$ -where $k$ is the hydraulic conductivity (mm s^-1^), and $\psi _{h}$ is +where $k$ is the hydraulic conductivity ($mm\ s^{-1}$), and $\psi_{h}$ is the hydraulic potential (mm). The hydraulic potential is -$$\psi _{h} =\psi _{m} +\psi _{z}$$ +$$\psi_{h} =\psi_{m} +\psi_{z}$$ -where $\psi _{m}$ is the soil matric potential (mm) (which is related to +where $\psi_{m}$ is the soil matric potential (mm) (which is related to the adsorptive and capillary forces within the soil matrix), and -$\psi _{z}$ is the gravitational potential (mm) (the vertical distance +$\psi_{z}$ is the gravitational potential (mm) (the vertical distance from an arbitrary reference elevation to a point in the soil). If the -reference elevation is the soil surface, then $\psi _{z} =z$. Letting -$\psi =\psi _{m}$, Darcy\'s law becomes +reference elevation is the soil surface, then $\psi_{z} =z$. Letting +$\psi =\psi_{m}$, Darcy\'s law becomes $$q = -k \left[\frac{\partial \left(\psi +z\right)}{\partial z} \right].$$ @@ -427,10 +429,10 @@ role="numref"}), changes in soil water content are predicted from `7.79`{.interpreted-text role="eq"} using finite-difference approximations for `7.84`{.interpreted-text role="eq"}. -### Hydraulic Properties {#Hydraulic Properties} +### Hydraulic Properties -The hydraulic conductivity $k_{i}$ (mm s^-1^) and the soil matric -potential $\psi _{i}$ (mm) for layer $i$ vary with volumetric soil water +The hydraulic conductivity $k_{i}$ ($mm\ s^{-1}$) and the soil matric +potential $\psi_{i}$ (mm) for layer $i$ vary with volumetric soil water $\theta _{i}$ and soil texture. As with the soil thermal properties (section `Soil And Snow Thermal Properties`{.interpreted-text role="numref"}) the hydraulic properties of the soil are assumed to be a @@ -450,13 +452,13 @@ $k_{sat} \left[z_{h,\, i} \right]$, the liquid volumetric soil moisture of the two layers $\theta _{i}$ and $\theta_{i+1}$ and an ice impedance factor $\Theta_{ice}$ -$$\begin{aligned} +$$ k\left[z_{h,\, i} \right] = -\left\{\begin{array}{lr} -\Theta_{ice} k_{sat} \left[z_{h,\, i} \right]\left[\frac{0.5\left(\theta_{\, i} +\theta_{\, i+1} \right)}{0.5\left(\theta_{sat,\, i} +\theta_{sat,\, i+1} \right)} \right]^{2B_{i} +3} & \qquad 1 \le i \le N_{levsoi} - 1 \\ -\Theta_{ice} k_{sat} \left[z_{h,\, i} \right]\left(\frac{\theta_{\, i} }{\theta_{sat,\, i} } \right)^{2B_{i} +3} & \qquad i = N_{levsoi} -\end{array}\right\}. -\end{aligned}$$ +\begin{cases} +\Theta_{\text{ice}} \, k_{\text{sat}}\left[z_{h,\, i} \right] \left[\dfrac{0.5\left(\theta_{i} + \theta_{i+1} \right)}{0.5\left(\theta_{\text{sat},\, i} + \theta_{\text{sat},\, i+1} \right)} \right]^{2B_{i} + 3}, & \text{for } 1 \le i \le N_{\text{levsoi}} - 1 \\ +\Theta_{\text{ice}} \, k_{\text{sat}}\left[z_{h,\, i} \right] \left(\dfrac{\theta_{i}}{\theta_{\text{sat},\, i}} \right)^{2B_{i} + 3}, & \text{for } i = N_{\text{levsoi}} +\end{cases} +$$ The ice impedance factor is a function of ice content, and is meant to quantify the increased tortuosity of the water flow when part of the @@ -494,19 +496,19 @@ The soil matric potential (mm) is defined at the node depth $z_{i}$ of each layer $i$ (`Figure Water flux schematic`{.interpreted-text role="numref"}) -$$\psi _{i} =\psi _{sat,\, i} \left(\frac{\theta_{\, i} }{\theta_{sat,\, i} } \right)^{-B_{i} } \ge -1\times 10^{8} \qquad 0.01\le \frac{\theta_{i} }{\theta_{sat,\, i} } \le 1$$ +$$\psi_{i} =\psi_{sat,\, i} \left(\frac{\theta_{\, i} }{\theta_{sat,\, i} } \right)^{-B_{i} } \ge -1\times 10^{8} \qquad 0.01\le \frac{\theta_{i} }{\theta_{sat,\, i}} \le 1$$ where the saturated soil matric potential (mm) is -$$\psi _{sat,i} =(1-f_{om,i} )\psi _{sat,\min ,i} +f_{om,i} \psi _{sat,om}$$ +$$\psi_{sat,i} =(1-f_{om,i} )\psi_{sat,\min ,i} +f_{om,i} \psi_{sat,om}$$ -where $\psi _{sat,om}$ is the saturated organic matter matric potential -and the saturated mineral soil matric potential $\psi _{sat,\min,i}$ is +where $\psi_{sat,om}$ is the saturated organic matter matric potential +and the saturated mineral soil matric potential $\psi_{sat,\min,i}$ is -$$\psi _{sat,\, \min ,\, i} =-10.0\times 10^{1.88-0.0131(\% sand)_{i} } .$$ +$$\psi_{sat,\, \min ,\, i} =-10.0\times 10^{1.88-0.0131(\% sand)_{i}} .$$ The saturated hydraulic conductivity, $k_{sat} \left[z_{h,\, i} \right]$ -(mm s^-1^), for organic soils ($k_{sat,\, om}$ ) may be two to three +($mm\ s^{-1}$), for organic soils ($k_{sat,\, om}$ ) may be two to three orders of magnitude larger than that of mineral soils ($k_{sat,\, \min }$ ). Bulk soil layer values of $k_{sat}$ calculated as weighted averages based on $f_{om}$ may therefore be determined @@ -566,7 +568,7 @@ $$k_{sat,om} = max(0.28 - 0.2799\times z_{i} / zsapric, k_{sat,\, \min } \left[z where $zsapric =0.5$ m is the depth that organic matter takes on the characteristics of sapric peat. -### Numerical Solution {#Numerical Solution Hydrology} +### Numerical Solution With reference to `Figure Water flux schematic`{.interpreted-text role="numref"}, the equation for conservation of mass (equation @@ -583,7 +585,7 @@ $$\Delta z_{i} \frac{\partial \theta_{liq,\, i} }{\partial t} =-q_{i-1} +q_{i} - where $q_{i}$ is the flux of water across interface $z_{h,\, i}$, $q_{i-1}$ is the flux of water across interface $z_{h,\, i-1}$, and $e_{i}$ is a layer-averaged soil moisture sink term (ET loss) defined as -positive for flow out of the layer (mm s^-1^). Taking the finite +positive for flow out of the layer ($mm\ s^{-1}$). Taking the finite difference with time and evaluating the fluxes implicitly at time $n+1$ yields @@ -683,26 +685,26 @@ equations `7.111`{.interpreted-text role="eq"} - `7.114`{.interpreted-text role="eq"} can be obtained from equation `7.82`{.interpreted-text role="eq"} as -$$q_{i-1}^{n} =-k\left[z_{h,\, i-1} \right]\left[\frac{\left(\psi _{i-1} -\psi _{i} \right)+\left(z_{i} - z_{i-1} \right)}{z_{i} -z_{i-1} } \right]$$ +$$q_{i-1}^{n} =-k\left[z_{h,\, i-1} \right]\left[\frac{\left(\psi_{i-1} -\psi_{i} \right)+\left(z_{i} - z_{i-1} \right)}{z_{i} -z_{i-1} } \right]$$ -$$q_{i}^{n} =-k\left[z_{h,\, i} \right]\left[\frac{\left(\psi _{i} -\psi _{i+1} \right)+\left(z_{i+1} - z_{i} \right)}{z_{i+1} -z_{i} } \right]$$ +$$q_{i}^{n} =-k\left[z_{h,\, i} \right]\left[\frac{\left(\psi_{i} -\psi_{i+1} \right)+\left(z_{i+1} - z_{i} \right)}{z_{i+1} -z_{i} } \right]$$ -$$\frac{\partial q_{i-1} }{\partial \theta _{liq,\, i-1} } =-\left[\frac{k\left[z_{h,\, i-1} \right]}{z_{i} -z_{i-1} } \frac{\partial \psi _{i-1} }{\partial \theta _{liq,\, i-1} } \right]-\frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta _{liq,\, i-1} } \left[\frac{\left(\psi _{i-1} -\psi _{i} \right)+\left(z_{i} - z_{i-1} \right)}{z_{i} - z_{i-1} } \right]$$ +$$\frac{\partial q_{i-1} }{\partial \theta _{liq,\, i-1} } =-\left[\frac{k\left[z_{h,\, i-1} \right]}{z_{i} -z_{i-1} } \frac{\partial \psi_{i-1} }{\partial \theta _{liq,\, i-1} } \right]-\frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta _{liq,\, i-1} } \left[\frac{\left(\psi_{i-1} -\psi_{i} \right)+\left(z_{i} - z_{i-1} \right)}{z_{i} - z_{i-1} } \right]$$ -$$\frac{\partial q_{i-1} }{\partial \theta _{liq,\, i} } =\left[\frac{k\left[z_{h,\, i-1} \right]}{z_{i} -z_{i-1} } \frac{\partial \psi _{i} }{\partial \theta _{liq,\, i} } \right]-\frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta _{liq,\, i} } \left[\frac{\left(\psi _{i-1} -\psi _{i} \right)+\left(z_{i} - z_{i-1} \right)}{z_{i} - z_{i-1} } \right]$$ +$$\frac{\partial q_{i-1} }{\partial \theta _{liq,\, i} } =\left[\frac{k\left[z_{h,\, i-1} \right]}{z_{i} -z_{i-1} } \frac{\partial \psi_{i} }{\partial \theta _{liq,\, i} } \right]-\frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta _{liq,\, i} } \left[\frac{\left(\psi_{i-1} -\psi_{i} \right)+\left(z_{i} - z_{i-1} \right)}{z_{i} - z_{i-1} } \right]$$ -$$\frac{\partial q_{i} }{\partial \theta _{liq,\, i} } =-\left[\frac{k\left[z_{h,\, i} \right]}{z_{i+1} -z_{i} } \frac{\partial \psi _{i} }{\partial \theta _{liq,\, i} } \right]-\frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta _{liq,\, i} } \left[\frac{\left(\psi _{i} -\psi _{i+1} \right)+\left(z_{i+1} - z_{i} \right)}{z_{i+1} - z_{i} } \right]$$ +$$\frac{\partial q_{i} }{\partial \theta _{liq,\, i} } =-\left[\frac{k\left[z_{h,\, i} \right]}{z_{i+1} -z_{i} } \frac{\partial \psi_{i} }{\partial \theta _{liq,\, i} } \right]-\frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta _{liq,\, i} } \left[\frac{\left(\psi_{i} -\psi_{i+1} \right)+\left(z_{i+1} - z_{i} \right)}{z_{i+1} - z_{i} } \right]$$ -$$\frac{\partial q_{i} }{\partial \theta _{liq,\, i+1} } =\left[\frac{k\left[z_{h,\, i} \right]}{z_{i+1} -z_{i} } \frac{\partial \psi _{i+1} }{\partial \theta _{liq,\, i+1} } \right]-\frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta _{liq,\, i+1} } \left[\frac{\left(\psi _{i} -\psi _{i+1} \right)+\left(z_{i+1} - z_{i} \right)}{z_{i+1} - z_{i} } \right].$$ +$$\frac{\partial q_{i} }{\partial \theta _{liq,\, i+1} } =\left[\frac{k\left[z_{h,\, i} \right]}{z_{i+1} -z_{i} } \frac{\partial \psi_{i+1} }{\partial \theta _{liq,\, i+1} } \right]-\frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta _{liq,\, i+1} } \left[\frac{\left(\psi_{i} -\psi_{i+1} \right)+\left(z_{i+1} - z_{i} \right)}{z_{i+1} - z_{i} } \right].$$ The derivatives of the soil matric potential at the node depth are derived from `7.94`{.interpreted-text role="eq"} -$$\frac{\partial \psi _{i-1} }{\partial \theta_{liq,\, \, i-1} } =-B_{i-1} \frac{\psi _{i-1} }{\theta_{\, \, i-1} }$$ +$$\frac{\partial \psi_{i-1} }{\partial \theta_{liq,\, \, i-1} } =-B_{i-1} \frac{\psi_{i-1} }{\theta_{\, \, i-1} }$$ -$$\frac{\partial \psi _{i} }{\partial \theta_{\, liq,\, i} } =-B_{i} \frac{\psi _{i} }{\theta_{i} }$$ +$$\frac{\partial \psi_{i} }{\partial \theta_{\, liq,\, i} } =-B_{i} \frac{\psi_{i} }{\theta_{i} }$$ -$$\frac{\partial \psi _{i+1} }{\partial \theta_{liq,\, i+1} } =-B_{i+1} \frac{\psi _{i+1} }{\theta_{\, i+1} }$$ +$$\frac{\partial \psi_{i+1} }{\partial \theta_{liq,\, i+1} } =-B_{i+1} \frac{\psi_{i+1} }{\theta_{\, i+1} }$$ with the constraint $0.01\, \theta_{sat,\, i} \le \theta_{\, i} \le \theta_{sat,\, i}$. @@ -826,7 +828,7 @@ The volumetric water content is $$\theta_{i} =\frac{w_{liq,\, i} }{\Delta z_{i} \rho _{liq} } +\frac{w_{ice,\, i} }{\Delta z_{i} \rho _{ice} } .$$ -## Frozen Soils and Perched Water Table {#Frozen Soils and Perched Water Table} +## Frozen Soils and Perched Water Table When soils freeze, the power-law form of the ice impedance factor (section `Hydraulic Properties`{.interpreted-text role="numref"}) can @@ -854,7 +856,7 @@ is set to 0.9. Drainage from the perched saturated zone $q_{drai,perch}$ is removed from layers $N_{perch}$ through $N_{frost}$, which are the layers containing $z_{\nabla,perch}$ and, $z_{frost}$ respectively. -## Lateral Sub-surface Runoff {#Lateral Sub-surface Runoff} +## Lateral Sub-surface Runoff Lateral sub-surface runoff occurs when saturated soil moisture conditions exist within the soil column. Sub-surface runoff is @@ -921,10 +923,10 @@ $$w_{ice,\, 1}^{n+1} =w_{ice,\, 1}^{n} -q_{subl} \Delta t.$$ Sublimation of ice is limited to the amount of ice available. -## Runoff from glaciers and snow-capped surfaces {#Runoff from glaciers and snow-capped surfaces} +## Runoff from glaciers and snow-capped surfaces All surfaces are constrained to have a snow water equivalent -$W_{sno} \le W_{cap} = 10,000$ kg m^-2^. For snow-capped columns, any +$W_{sno} \le W_{cap} = 10,000$ $kg\ m^{-2}$. For snow-capped columns, any addition of mass at the top (precipitation, dew/riping) is balanced by an equally large mass flux at the bottom of the snow column. This so-called capping flux is separated into solid $q_{snwcp,ice}$ and From ec41a40ec4b4a620b636c2ab8909cd448f5a6043 Mon Sep 17 00:00:00 2001 From: Yilin Fang Date: Fri, 6 Jun 2025 16:07:44 -0700 Subject: [PATCH 5/9] correction --- components/elm/docs/tech-guide/hydrology.md | 79 +++++++++++---------- 1 file changed, 43 insertions(+), 36 deletions(-) diff --git a/components/elm/docs/tech-guide/hydrology.md b/components/elm/docs/tech-guide/hydrology.md index 7fe785da36db..dcb49c3d343d 100644 --- a/components/elm/docs/tech-guide/hydrology.md +++ b/components/elm/docs/tech-guide/hydrology.md @@ -13,9 +13,13 @@ $\Delta W_{a}$ (all in $kg\ m^{-2}$ or $mm\ of\ H_2O$) The total water balance of the system is -$$\begin{aligned} -{\Delta W_{can,\,liq} +\Delta W_{can,\,sno} +\Delta W_{sfc} +\Delta W_{sno} +} \\ {\sum _{i=1}^{N_{levsoi} }\left(\Delta w_{liq,\, i} +\Delta w_{ice,\, i} \right)+\Delta W_{a} =\left(\begin{array}{l} {q_{rain} +q_{sno} -E_{v} -E_{g} -q_{over} } \\ {-q_{h2osfc} -q_{drai} -q_{rgwl} -q_{snwcp,\, ice} } \end{array}\right) \Delta t} -\end{aligned}$$ +$$ +\begin{aligned} +\Delta W_{can,\,liq} + \Delta W_{can,\,sno} + \Delta W_{sfc} + \Delta W_{sno} + +\sum_{i=1}^{N_{levsoi}} (\Delta w_{liq,\, i} + \Delta w_{ice,\, i}) + \Delta W_{a} \\ += \left( q_{rain} + q_{sno} - E_{v} - E_{g} - q_{over} - q_{h2osfc} - q_{drai} - q_{rgwl} - q_{snwcp,\, ice} \right) \Delta t +\end{aligned} +$$

(1)

where $q_{rain}$ is the liquid part of precipitation, $q_{sno}$ is the @@ -271,24 +275,24 @@ role="numref"}). Surface water storage and outflow are functions of fine spatial scale elevation variations called microtopography. The microtopography is assumed to be distributed normally around the grid cell mean elevation. Given the standard deviation of the -microtopographic distribution, $\sigma _{micro}$ (m), the fractional +microtopographic distribution, $\sigma_{micro}$ (m), the fractional area of the grid cell that is inundated can be calculated. Surface water storage, $Wsfc$, is related to the height (relative to the grid cell mean elevation) of the surface water, $d$, by -$$W_{sfc} =\frac{d}{2} \left(1+erf\left(\frac{d}{\sigma _{micro} \sqrt{2} } \right)\right)+\frac{\sigma _{micro} }{\sqrt{2\pi } } e^{\frac{-d^{2} }{2\sigma _{micro} ^{2} } }$$ +$$W_{sfc} =\frac{d}{2} \left(1+erf\left(\frac{d}{\sigma_{micro} \sqrt{2} } \right)\right)+\frac{\sigma_{micro} }{\sqrt{2\pi } } e^{\frac{-d^{2} }{2\sigma_{micro} ^{2} } }$$ where $erf$ is the error function. For a given value of $W_{sfc}$, `7.66`{.interpreted-text role="eq"} can be solved for $d$ using the Newton-Raphson method. Once $d$ is known, one can determine the fraction of the area that is inundated as -$$f_{h2osfc} =\frac{1}{2} \left(1+erf\left(\frac{d}{\sigma _{micro} \sqrt{2} } \right)\right)$$ +$$f_{h2osfc} =\frac{1}{2} \left(1+erf\left(\frac{d}{\sigma_{micro} \sqrt{2} } \right)\right)$$ No global datasets exist for microtopography, so the default parameterization is a simple function of slope -$$\sigma _{micro} =\left(\beta +\beta _{0} \right)^{\eta }$$ +$$\sigma_{micro} =\left(\beta +\beta_{0} \right)^{\eta }$$ where $\beta$ is the topographic slope, $\beta_{0} =\left(\sigma_{\max } \right)^{\frac{1}{\eta } }$ determines @@ -433,7 +437,7 @@ approximations for `7.84`{.interpreted-text role="eq"}. The hydraulic conductivity $k_{i}$ ($mm\ s^{-1}$) and the soil matric potential $\psi_{i}$ (mm) for layer $i$ vary with volumetric soil water -$\theta _{i}$ and soil texture. As with the soil thermal properties +$\theta_{i}$ and soil texture. As with the soil thermal properties (section `Soil And Snow Thermal Properties`{.interpreted-text role="numref"}) the hydraulic properties of the soil are assumed to be a weighted combination of the mineral properties, which are determined @@ -449,7 +453,7 @@ two adjacent layers $z_{h,\, i}$ (`Figure Water flux schematic`{.interpreted-text role="numref"}) and is a function of the saturated hydraulic conductivity $k_{sat} \left[z_{h,\, i} \right]$, the liquid volumetric soil moisture -of the two layers $\theta _{i}$ and $\theta_{i+1}$ and an ice impedance +of the two layers $\theta_{i}$ and $\theta_{i+1}$ and an ice impedance factor $\Theta_{ice}$ $$ @@ -482,7 +486,7 @@ where $f_{om,i}$ is the soil organic matter fraction, $\theta_{sat,om}$ is the porosity of organic matter, and the porosity of the mineral soil $\theta_{sat,\min,i}$ is -$$\theta_{sat,\min ,i} = 0.489 - 0.00126(\% sand)_{i} .$$ +$$\theta_{sat, min, i} = 0.489 - 0.00126(\text{% sand})_{i}$$ The exponent $B_{i}$ is @@ -490,7 +494,7 @@ $$B_{i} =(1-f_{om,i} )B_{\min ,i} +f_{om,i} B_{om}$$ where $B_{om}$ is for organic matter and -$$B_{\min ,i} =2.91+0.159(\% clay)_{i} .$$ +$$B_{\min ,i} =2.91+0.159(\text{% clay})_{i} $$ The soil matric potential (mm) is defined at the node depth $z_{i}$ of each layer $i$ (`Figure Water flux schematic`{.interpreted-text @@ -505,7 +509,7 @@ $$\psi_{sat,i} =(1-f_{om,i} )\psi_{sat,\min ,i} +f_{om,i} \psi_{sat,om}$$ where $\psi_{sat,om}$ is the saturated organic matter matric potential and the saturated mineral soil matric potential $\psi_{sat,\min,i}$ is -$$\psi_{sat,\, \min ,\, i} =-10.0\times 10^{1.88-0.0131(\% sand)_{i}} .$$ +$$\psi_{sat,\, \min ,\, i} =-10.0\times 10^{1.88-0.0131(\text{% sand})_{i}} $$ The saturated hydraulic conductivity, $k_{sat} \left[z_{h,\, i} \right]$ ($mm\ s^{-1}$), for organic soils ($k_{sat,\, om}$ ) may be two to three @@ -528,12 +532,13 @@ span the soil space. Flow through these pathways interacts only with organic material, and thus can be described by $k_{sat,\, om}$. This fraction of the grid cell is given by -$$\begin{aligned} -\begin{array}{lr} -f_{perc} =\; N_{perc} \left(f_{om} {\rm \; }-f_{threshold} \right)^{\beta_{perc} } f_{om} {\rm \; } & \qquad f_{om} \ge f_{threshold} \\ -f_{perc} = 0 & \qquad f_{om} `{.interpreted-text role="ref"}) as -$$k_{sat,\, \min } \left[z_{h,\, i} \right]=0.0070556\times 10^{-0.884+0.0153\left(\% sand\right)_{i} } .$$ +$$k_{sat,\, \min } \left[z_{h,\, i} \right]=0.0070556\times 10^{-0.884+0.0153\left(\text{% sand}\right)_{i} } .$$ The bulk soil layer saturated hydraulic conductivity is then computed as @@ -575,7 +580,7 @@ role="numref"}, the equation for conservation of mass (equation `7.79`{.interpreted-text role="eq"}) can be integrated over each layer as -$$\int _{-z_{h,\, i} }^{-z_{h,\, i-1} }\frac{\partial \theta }{\partial t} \, dz=-\int _{-z_{h,\, i} }^{-z_{h,\, i-1} }\frac{\partial q}{\partial z} \, dz-\int _{-z_{h,\, i} }^{-z_{h,\, i-1} } e\, dz .$$ +$$\int_{-z_{h,\, i} }^{-z_{h,\, i-1} }\frac{\partial \theta }{\partial t} \, dz=-\int_{-z_{h,\, i} }^{-z_{h,\, i-1} }\frac{\partial q}{\partial z} \, dz-\int_{-z_{h,\, i} }^{-z_{h,\, i-1} } e\, dz .$$ Note that the integration limits are negative since $z$ is defined as positive upward from the soil surface. This equation can be written as @@ -621,7 +626,9 @@ Note that because more than one plant functional type (PFT) may share a soil column, the transpiration $E_{v}^{t}$ is a weighted sum of transpiration from all PFTs whose weighting depends on PFT area as -$$E_{v}^{t} =\sum _{j=1}^{npft}\left(E_{v}^{t} \right)_{j} \left(wt\right)_{j}$$ +$$ +E_{v}^{t} = \sum_{j=1}^{n_{\text{pft}}} \left( E_{v,j}^{t} \cdot \text{wt}_j \right) +$$ where $npft$ is the number of PFTs sharing a soil column, $\left(E_{v}^{t} \right)_{j}$ is the transpiration from the $j^{th}$ PFT @@ -631,15 +638,15 @@ $r_{e,\, i}$ is also a column-level quantity that is a weighted sum over all PFTs. The weighting depends on the per unit area transpiration of each PFT and its relative area as -$$r_{e,\, i} =\frac{\sum _{j=1}^{npft}\left(r_{e,\, i} \right)_{j} \left(E_{v}^{t} \right)_{j} \left(wt\right)_{j} }{\sum _{j=1}^{npft}\left(E_{v}^{t} \right)_{j} \left(wt\right)_{j} }$$ +$$r_{e,\, i} =\frac{\sum_{j=1}^{npft}\left(r_{e,\, i} \right)_{j} \left(E_{v}^{t} \right)_{j} \left(wt\right)_{j} }{\sum_{j=1}^{npft}\left(E_{v}^{t} \right)_{j} \left(wt\right)_{j} }$$ where $\left(r_{e,\, i} \right)_{j}$ is the effective root fraction for the $j^{th}$ PFT $$\begin{aligned} \begin{array}{lr} -\left(r_{e,\, i} \right)_{j} =\frac{\left(r_{i} \right)_{j} \left(w_{i} \right)_{j} }{\left(\beta _{t} \right)_{j} } & \qquad \left(\beta _{t} \right)_{j} >0 \\ -\left(r_{e,\, i} \right)_{j} =0 & \qquad \left(\beta _{t} \right)_{j} =0 +\left(r_{e,\, i} \right)_{j} =\frac{\left(r_{i} \right)_{j} \left(w_{i} \right)_{j} }{\left(\beta_{t} \right)_{j} } & \qquad \left(\beta_{t} \right)_{j} >0 \\ +\left(r_{e,\, i} \right)_{j} =0 & \qquad \left(\beta_{t} \right)_{j} =0 \end{array} \end{aligned}$$ @@ -689,13 +696,13 @@ $$q_{i-1}^{n} =-k\left[z_{h,\, i-1} \right]\left[\frac{\left(\psi_{i-1} -\psi_{i $$q_{i}^{n} =-k\left[z_{h,\, i} \right]\left[\frac{\left(\psi_{i} -\psi_{i+1} \right)+\left(z_{i+1} - z_{i} \right)}{z_{i+1} -z_{i} } \right]$$ -$$\frac{\partial q_{i-1} }{\partial \theta _{liq,\, i-1} } =-\left[\frac{k\left[z_{h,\, i-1} \right]}{z_{i} -z_{i-1} } \frac{\partial \psi_{i-1} }{\partial \theta _{liq,\, i-1} } \right]-\frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta _{liq,\, i-1} } \left[\frac{\left(\psi_{i-1} -\psi_{i} \right)+\left(z_{i} - z_{i-1} \right)}{z_{i} - z_{i-1} } \right]$$ +$$\frac{\partial q_{i-1} }{\partial \theta_{liq,\, i-1} } =-\left[\frac{k\left[z_{h,\, i-1} \right]}{z_{i} -z_{i-1} } \frac{\partial \psi_{i-1} }{\partial \theta_{liq,\, i-1} } \right]-\frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta_{liq,\, i-1} } \left[\frac{\left(\psi_{i-1} -\psi_{i} \right)+\left(z_{i} - z_{i-1} \right)}{z_{i} - z_{i-1} } \right]$$ -$$\frac{\partial q_{i-1} }{\partial \theta _{liq,\, i} } =\left[\frac{k\left[z_{h,\, i-1} \right]}{z_{i} -z_{i-1} } \frac{\partial \psi_{i} }{\partial \theta _{liq,\, i} } \right]-\frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta _{liq,\, i} } \left[\frac{\left(\psi_{i-1} -\psi_{i} \right)+\left(z_{i} - z_{i-1} \right)}{z_{i} - z_{i-1} } \right]$$ +$$\frac{\partial q_{i-1} }{\partial \theta_{liq,\, i} } =\left[\frac{k\left[z_{h,\, i-1} \right]}{z_{i} -z_{i-1} } \frac{\partial \psi_{i} }{\partial \theta_{liq,\, i} } \right]-\frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta_{liq,\, i} } \left[\frac{\left(\psi_{i-1} -\psi_{i} \right)+\left(z_{i} - z_{i-1} \right)}{z_{i} - z_{i-1} } \right]$$ -$$\frac{\partial q_{i} }{\partial \theta _{liq,\, i} } =-\left[\frac{k\left[z_{h,\, i} \right]}{z_{i+1} -z_{i} } \frac{\partial \psi_{i} }{\partial \theta _{liq,\, i} } \right]-\frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta _{liq,\, i} } \left[\frac{\left(\psi_{i} -\psi_{i+1} \right)+\left(z_{i+1} - z_{i} \right)}{z_{i+1} - z_{i} } \right]$$ +$$\frac{\partial q_{i} }{\partial \theta_{liq,\, i} } =-\left[\frac{k\left[z_{h,\, i} \right]}{z_{i+1} -z_{i} } \frac{\partial \psi_{i} }{\partial \theta_{liq,\, i} } \right]-\frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta_{liq,\, i} } \left[\frac{\left(\psi_{i} -\psi_{i+1} \right)+\left(z_{i+1} - z_{i} \right)}{z_{i+1} - z_{i} } \right]$$ -$$\frac{\partial q_{i} }{\partial \theta _{liq,\, i+1} } =\left[\frac{k\left[z_{h,\, i} \right]}{z_{i+1} -z_{i} } \frac{\partial \psi_{i+1} }{\partial \theta _{liq,\, i+1} } \right]-\frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta _{liq,\, i+1} } \left[\frac{\left(\psi_{i} -\psi_{i+1} \right)+\left(z_{i+1} - z_{i} \right)}{z_{i+1} - z_{i} } \right].$$ +$$\frac{\partial q_{i} }{\partial \theta_{liq,\, i+1} } =\left[\frac{k\left[z_{h,\, i} \right]}{z_{i+1} -z_{i} } \frac{\partial \psi_{i+1} }{\partial \theta_{liq,\, i+1} } \right]-\frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta_{liq,\, i+1} } \left[\frac{\left(\psi_{i} -\psi_{i+1} \right)+\left(z_{i+1} - z_{i} \right)}{z_{i+1} - z_{i} } \right].$$ The derivatives of the soil matric potential at the node depth are derived from `7.94`{.interpreted-text role="eq"} @@ -713,8 +720,8 @@ The derivatives of the hydraulic conductivity at the layer interface are derived from `7.85`{.interpreted-text role="eq"} $$\begin{array}{l} -{\frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta _{liq,\, i-1} } -= \frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta _{liq,\, i} } +{\frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta_{liq,\, i-1} } += \frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta_{liq,\, i} } = \left(2B_{i-1} +3\right) \ \overline{\Theta}_{ice} \ k_{sat} \left[z_{h,\, i-1} \right] \ \left[\frac{\overline{\theta}_{liq}}{\overline{\theta}_{sat}} \right]^{2B_{i-1} +2} \left(\frac{0.5}{\overline{\theta}_{sat}} \right)} \end{array}$$ where $\overline{\Theta}_{ice} = \Theta(\overline{\theta}_{ice})$ @@ -727,8 +734,8 @@ $\overline{\theta}_{sat} = 0.5\left(\theta_{sat,\, i-1} +\theta_{sat,\, i} \righ and $$\begin{array}{l} -{\frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta _{liq,\, i} } -= \frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta _{liq,\, i+1} } +{\frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta_{liq,\, i} } += \frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta_{liq,\, i+1} } = \left(2B_{i} +3\right) \ \overline{\Theta}_{ice} \ k_{sat} \left[z_{h,\, i} \right] \ \left[\frac{\overline{\theta}_{liq}}{\overline{\theta}_{sat}} \right]^{2B_{i} +2} \left(\frac{0.5}{\overline{\theta}_{sat}} \right)} \end{array}.$$ where @@ -826,7 +833,7 @@ $$w_{liq,\, i}^{n+1} =w_{liq,\, i}^{n} +\Delta \theta_{liq,\, i} \Delta z_{i} \q The volumetric water content is -$$\theta_{i} =\frac{w_{liq,\, i} }{\Delta z_{i} \rho _{liq} } +\frac{w_{ice,\, i} }{\Delta z_{i} \rho _{ice} } .$$ +$$\theta_{i} =\frac{w_{liq,\, i} }{\Delta z_{i} \rho_{liq} } +\frac{w_{ice,\, i} }{\Delta z_{i} \rho_{ice} } .$$ ## Frozen Soils and Perched Water Table @@ -843,7 +850,7 @@ $$q_{drai,perch} =k_{drai,\, perch} \left(z_{frost} -z_{\nabla ,perch} \right)$$ where $k_{drai,\, perch}$ depends on topographic slope and soil hydraulic conductivity, -$$k_{drai,\, perch} =10^{-5} \sin (\beta )\left(\frac{\sum _{i=N_{perch} }^{i=N_{frost} }\Theta_{ice,i} k_{sat} \left[z_{i} \right]\Delta z_{i} }{\sum _{i=N_{perch} }^{i=N_{frost} }\Delta z_{i} } \right)$$ +$$k_{drai,\, perch} =10^{-5} \sin (\beta )\left(\frac{\sum_{i=N_{perch} }^{i=N_{frost} }\Theta_{ice,i} k_{sat} \left[z_{i} \right]\Delta z_{i} }{\sum_{i=N_{perch} }^{i=N_{frost} }\Delta z_{i} } \right)$$ where $\Theta_{ice}$ is an ice impedance factor, $\beta$ is the mean grid cell topographic slope in radians, $z_{frost}$ is the depth to the @@ -883,7 +890,7 @@ the water table location, is derived by taking the difference between two equilibrium soil moisture profiles whose water tables differ by an infinitesimal amount -$$S_{y} =\theta_{sat} \left(1-\left(1+\frac{z_{\nabla } }{\Psi _{sat} } \right)^{\frac{-1}{B} } \right)$$ +$$S_{y} =\theta_{sat} \left(1-\left(1+\frac{z_{\nabla } }{\Psi_{sat} } \right)^{\frac{-1}{B} } \right)$$ where B is the Clapp-Hornberger exponent. Because $S_{y}$ is a function of the soil properties, it results in water table dynamics that are @@ -949,7 +956,7 @@ $$q_{rgwl} =q_{grnd,ice} +q_{grnd,liq} -E_{g} -E_{v} -\frac{\left(W_{b}^{n+1} -W where $W_{b}^{n}$ and $W_{b}^{n+1}$ are the water balances at the beginning and ending of the time step defined as -$$W_{b} =W_{can} +W_{sno} +\sum _{i=1}^{N}\left(w_{ice,i} +w_{liq,i} \right) .$$ +$$W_{b} =W_{can} +W_{sno} +\sum_{i=1}^{N}\left(w_{ice,i} +w_{liq,i} \right) .$$ Currently, glaciers are non-vegetated and $E_{v} =W_{can} =0$. The contribution of lake runoff to $q_{rgwl}$ is described in section From 6d2c9f80b1c946a635ac34237e4c88f6ec8f28ee Mon Sep 17 00:00:00 2001 From: Yilin Fang Date: Mon, 9 Jun 2025 11:51:27 -0700 Subject: [PATCH 6/9] Remove spaces. --- components/elm/docs/tech-guide/hydrology.md | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/components/elm/docs/tech-guide/hydrology.md b/components/elm/docs/tech-guide/hydrology.md index dcb49c3d343d..9253d3ab2ebd 100644 --- a/components/elm/docs/tech-guide/hydrology.md +++ b/components/elm/docs/tech-guide/hydrology.md @@ -15,7 +15,7 @@ The total water balance of the system is $$ \begin{aligned} -\Delta W_{can,\,liq} + \Delta W_{can,\,sno} + \Delta W_{sfc} + \Delta W_{sno} + +\Delta W_{can,\,liq} + \Delta W_{can,\,sno} + \Delta W_{sfc} + \Delta W_{sno} + \sum_{i=1}^{N_{levsoi}} (\Delta w_{liq,\, i} + \Delta w_{ice,\, i}) + \Delta W_{a} \\ = \left( q_{rain} + q_{sno} - E_{v} - E_{g} - q_{over} - q_{h2osfc} - q_{drai} - q_{rgwl} - q_{snwcp,\, ice} \right) \Delta t \end{aligned} @@ -965,4 +965,3 @@ role="numref"}. The runoff term $q_{rgwl}$ may be negative for glaciers and lakes, which reduces the total amount of runoff available to the river routing model (Chapter `rst_MOSART`{.interpreted-text role="numref"}). - From e94dbacc4dbc89f66a91e9af1c5932335023e771 Mon Sep 17 00:00:00 2001 From: Yilin Fang Date: Mon, 9 Jun 2025 12:07:02 -0700 Subject: [PATCH 7/9] correction --- components/elm/docs/tech-guide/hydrology.md | 1 - components/elm/docs/tech-guide/phs.md | 1 - 2 files changed, 2 deletions(-) diff --git a/components/elm/docs/tech-guide/hydrology.md b/components/elm/docs/tech-guide/hydrology.md index 9253d3ab2ebd..3fd2e3c95516 100644 --- a/components/elm/docs/tech-guide/hydrology.md +++ b/components/elm/docs/tech-guide/hydrology.md @@ -20,7 +20,6 @@ $$ = \left( q_{rain} + q_{sno} - E_{v} - E_{g} - q_{over} - q_{h2osfc} - q_{drai} - q_{rgwl} - q_{snwcp,\, ice} \right) \Delta t \end{aligned} $$ -

(1)

where $q_{rain}$ is the liquid part of precipitation, $q_{sno}$ is the solid part of precipitation, $E_{v}$ is ET from vegetation (Chapter diff --git a/components/elm/docs/tech-guide/phs.md b/components/elm/docs/tech-guide/phs.md index 3c6db2af327e..60d3ca50b912 100644 --- a/components/elm/docs/tech-guide/phs.md +++ b/components/elm/docs/tech-guide/phs.md @@ -639,4 +639,3 @@ temperature, reflecting leaf surface energy balance. ::: {#Figure PHS Flow Diagram} ![Flow diagram of leaf flux calculations](phs_iteration_schematic.*) ::: - From 4aaf5850301a39c21c1e00190ba9ef2f4da02d17 Mon Sep 17 00:00:00 2001 From: Yilin Fang Date: Mon, 16 Jun 2025 13:54:25 -0700 Subject: [PATCH 8/9] correction --- components/elm/docs/tech-guide/hydrology.md | 33 +++++++++++++-------- 1 file changed, 20 insertions(+), 13 deletions(-) diff --git a/components/elm/docs/tech-guide/hydrology.md b/components/elm/docs/tech-guide/hydrology.md index 3fd2e3c95516..5ee275504a71 100644 --- a/components/elm/docs/tech-guide/hydrology.md +++ b/components/elm/docs/tech-guide/hydrology.md @@ -471,7 +471,7 @@ used a power law form $$\Theta_{ice} = 10^{-\Omega F_{ice} }$$ -where $\Omega = 6$and $F_{ice} = \frac{\theta_{ice} }{\theta_{sat} }$ is +where $\Omega = 6$ and $F_{ice} = \frac{\theta_{ice} }{\theta_{sat} }$ is the ice-filled fraction of the pore space. Because the hydraulic properties of mineral and organic soil may differ @@ -630,24 +630,26 @@ E_{v}^{t} = \sum_{j=1}^{n_{\text{pft}}} \left( E_{v,j}^{t} \cdot \text{wt}_j \ri $$ where $npft$ is the number of PFTs sharing a soil column, -$\left(E_{v}^{t} \right)_{j}$ is the transpiration from the $j^{th}$ PFT -on the column, and $\left(wt\right)_{j}$ is the relative area of the +$\left(E_{v}^{t} \right)\_{j}$ is the transpiration from the $j^{th}$ PFT +on the column, and $\left(wt\right)\_{j}$ is the relative area of the $j^{th}$ PFT with respect to the column. The effective root fraction $r_{e,\, i}$ is also a column-level quantity that is a weighted sum over all PFTs. The weighting depends on the per unit area transpiration of each PFT and its relative area as -$$r_{e,\, i} =\frac{\sum_{j=1}^{npft}\left(r_{e,\, i} \right)_{j} \left(E_{v}^{t} \right)_{j} \left(wt\right)_{j} }{\sum_{j=1}^{npft}\left(E_{v}^{t} \right)_{j} \left(wt\right)_{j} }$$ +$$r_{e,i} = \frac{\sum_{j=1}^{npft} (r_{e,i})\_j (E_v^t)\_j (wt)\_j}{\sum_{j=1}^{npft} (E_v^t)\_j (wt)\_j}$$ -where $\left(r_{e,\, i} \right)_{j}$ is the effective root fraction for +where $\left(r_{e,i} \right)_{j}$ is the effective root fraction for the $j^{th}$ PFT -$$\begin{aligned} +$$ +\begin{aligned} \begin{array}{lr} -\left(r_{e,\, i} \right)_{j} =\frac{\left(r_{i} \right)_{j} \left(w_{i} \right)_{j} }{\left(\beta_{t} \right)_{j} } & \qquad \left(\beta_{t} \right)_{j} >0 \\ -\left(r_{e,\, i} \right)_{j} =0 & \qquad \left(\beta_{t} \right)_{j} =0 +\left(r_{e,i} \right)\_{j} = \frac{\left(r_{i} \right)\_{j} \left(w_{i} \right)\_{j}}{\left(\beta_{t} \right)\_{j}} & \qquad \left(\beta_{t} \right)\_{j} > 0 \\ +\left(r_{e,i} \right)\_{j} = 0 & \qquad \left(\beta_{t} \right)_{j} = 0 \end{array} -\end{aligned}$$ +\end{aligned} +$$ and $\left(r_{i} \right)_{j}$ is the fraction of roots in layer $i$ (Chapter `rst_Stomatal Resistance and Photosynthesis`{.interpreted-text @@ -718,10 +720,15 @@ $0.01\, \theta_{sat,\, i} \le \theta_{\, i} \le \theta_{sat,\, i}$. The derivatives of the hydraulic conductivity at the layer interface are derived from `7.85`{.interpreted-text role="eq"} -$$\begin{array}{l} -{\frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta_{liq,\, i-1} } -= \frac{\partial k\left[z_{h,\, i-1} \right]}{\partial \theta_{liq,\, i} } -= \left(2B_{i-1} +3\right) \ \overline{\Theta}_{ice} \ k_{sat} \left[z_{h,\, i-1} \right] \ \left[\frac{\overline{\theta}_{liq}}{\overline{\theta}_{sat}} \right]^{2B_{i-1} +2} \left(\frac{0.5}{\overline{\theta}_{sat}} \right)} \end{array}$$ +$$ +\begin{array}{l} +\frac{\partial k\left[z_{h,i-1} \right]}{\partial \theta_{liq,i-1}} += \frac{\partial k\left[z_{h,i-1} \right]}{\partial \theta_{liq,i}} += \left(2B_{i-1} + 3\right) \, \overline{\Theta}\_{ice} \, k_{sat}\left[z_{h,i-1} \right] +\left( \frac{\overline{\theta}\_{liq}}{\overline{\theta}\_{sat}} \right)^{2B_{i-1} + 2} +\left( \frac{0.5}{\overline{\theta}_{sat}} \right) +\end{array} +$$ where $\overline{\Theta}_{ice} = \Theta(\overline{\theta}_{ice})$ `7.86`{.interpreted-text role="eq"}, From d7390ee0a72dc82921d1147aab3eb277fa9daf16 Mon Sep 17 00:00:00 2001 From: Yilin Fang Date: Mon, 16 Jun 2025 14:45:51 -0700 Subject: [PATCH 9/9] correction --- components/elm/docs/tech-guide/hydrology.md | 22 +++++++++++---------- 1 file changed, 12 insertions(+), 10 deletions(-) diff --git a/components/elm/docs/tech-guide/hydrology.md b/components/elm/docs/tech-guide/hydrology.md index 5ee275504a71..231c87d20042 100644 --- a/components/elm/docs/tech-guide/hydrology.md +++ b/components/elm/docs/tech-guide/hydrology.md @@ -730,23 +730,25 @@ $$ \end{array} $$ -where $\overline{\Theta}_{ice} = \Theta(\overline{\theta}_{ice})$ +where $\overline{\Theta}\_{ice} = \Theta(\overline{\theta}\_{ice})$ `7.86`{.interpreted-text role="eq"}, -$\overline{\theta}_{ice} = 0.5\left(\theta_{ice\, i-1} +\theta_{ice\, i} \right)$, -$\overline{\theta}_{liq} = 0.5\left(\theta_{liq\, i-1} +\theta_{liq\, i} \right)$, +$\overline{\theta}\_{ice} = 0.5\left(\theta_{ice\, i-1} +\theta_{ice\, i} \right)$, +$\overline{\theta}\_{liq} = 0.5\left(\theta_{liq\, i-1} +\theta_{liq\, i} \right)$, and -$\overline{\theta}_{sat} = 0.5\left(\theta_{sat,\, i-1} +\theta_{sat,\, i} \right)$ +$\overline{\theta}\_{sat} = 0.5\left(\theta_{sat,\, i-1} +\theta_{sat,\, i} \right)$ and -$$\begin{array}{l} -{\frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta_{liq,\, i} } -= \frac{\partial k\left[z_{h,\, i} \right]}{\partial \theta_{liq,\, i+1} } -= \left(2B_{i} +3\right) \ \overline{\Theta}_{ice} \ k_{sat} \left[z_{h,\, i} \right] \ \left[\frac{\overline{\theta}_{liq}}{\overline{\theta}_{sat}} \right]^{2B_{i} +2} \left(\frac{0.5}{\overline{\theta}_{sat}} \right)} \end{array}.$$ +$$ +\begin{array}{l} +{\frac{\partial k\left[z_{h,i} \right]}{\partial \theta_{liq,i}} += \frac{\partial k\left[z_{h,i} \right]}{\partial \theta_{liq,i+1}} += \left(2B_{i} +3\right) \, \overline{\Theta}\_{ice} \ k_{sat} \left[z_{h,i} \right] +\ \left(\frac{\overline{\theta}\_{liq}}{\overline{\theta}\_{sat}} \right)^{2B_{i} +2} \left(\frac{0.5}{\overline{\theta}_{sat}} \right)} \end{array}.$$ where -$\overline{\theta}_{liq} = 0.5\left(\theta_{\, i} +\theta_{\, i+1} \right)$, -$\overline{\theta}_{sat} = 0.5\left(\theta_{sat,\, i} +\theta_{sat,\, i+1} \right)$. +$\overline{\theta}\_{liq} = 0.5\left(\theta_{\, i} +\theta_{\, i+1} \right)$, +$\overline{\theta}\_{sat} = 0.5\left(\theta_{sat,\, i} +\theta_{sat,\, i+1} \right)$. #### Equation set for layer $i=1$