2.25. Methane Model

The representation of processes in the methane biogeochemical model integrated in CLM [CLM4Me; (Riley et al. 2011a)] is based on several previously published models (Cao et al. 1996; Petrescu et al. 2010; Tianet al. 2010; Walter et al. 2001; Wania et al. 2010; Zhang et al. 2002; Zhuang et al. 2004). Although the model has similarities with these precursor models, a number of new process representations and parameterization have been integrated into CLM.

Mechanistically modeling net surface CH4 emissions requires representing a complex and interacting series of processes. We first (section 2.25.1) describe the overall model structure and flow of information in the CH4 model, then describe the methods used to represent: CH4 mass balance; CH4 production; ebullition; aerenchyma transport; CH4 oxidation; reactive transport solution, including boundary conditions, numerical solution, water table interface, etc.; seasonal inundation effects; and impact of seasonal inundation on CH4 production.

2.25.1. Methane Model Structure and Flow

The driver routine for the methane biogeochemistry calculations (ch4, in ch4Mod.F) controls the initialization of boundary conditions, inundation, and impact of redox conditions; calls to routines to calculate CH4 production, oxidation, transport through aerenchyma, ebullition, and the overall mass balance (for unsaturated and saturated soils and, if desired, lakes); resolves changes to CH4 calculations associated with a changing inundated fraction; performs a mass balance check; and calculates the average gridcell CH4 production, oxidation, and exchanges with the atmosphere.

2.25.2. Governing Mass-Balance Relationship

The model (Figure 2.25.1) accounts for CH4 production in the anaerobic fraction of soil (P, mol m-3 s-1), ebullition (E, mol m-3 s-1), aerenchyma transport (A, mol m-3 s-1), aqueous and gaseous diffusion (\({F}_{D}\), mol m-2 s-1), and oxidation (O, mol m-3 s-1) via a transient reaction diffusion equation:

(2.25.1)\[\frac{\partial \left(RC\right)}{\partial t} =\frac{\partial F_{D} }{\partial z} +P\left(z,t\right)-E\left(z,t\right)-A\left(z,t\right)-O\left(z,t\right)\]

Here z (m) represents the vertical dimension, t (s) is time, and R accounts for gas in both the aqueous and gaseous phases:\(R = \epsilon _{a} +K_{H} \epsilon _{w}\), with \(\epsilon _{a}\), \(\epsilon _{w}\), and \(K_{H}\) (-) the air-filled porosity, water-filled porosity, and partitioning coefficient for the species of interest, respectively, and \(C\) represents CH4 or O2 concentration with respect to water volume (mol m-3).

An analogous version of equation (2.25.1) is concurrently solved for O2, but with the following differences relative to CH4: P = E = 0 (i.e., no production or ebullition), and the oxidation sink includes the O2 demanded by methanotrophs, heterotroph decomposers, nitrifiers, and autotrophic root respiration.

As currently implemented, each gridcell contains an inundated and a non-inundated fraction. Therefore, equation (2.25.1) is solved four times for each gridcell and time step: in the inundated and non-inundated fractions, and for CH4 and O2. If desired, the CH4 and O2 mass balance equation is solved again for lakes (Chapter 9). For non-inundated areas, the water table interface is defined at the deepest transition from greater than 95% saturated to less than 95% saturated that occurs above frozen soil layers. The inundated fraction is allowed to change at each time step, and the total soil CH4 quantity is conserved by evolving CH4 to the atmosphere when the inundated fraction decreases, and averaging a portion of the non-inundated concentration into the inundated concentration when the inundated fraction increases.

../../_images/image14.png

Figure 2.25.1 Schematic representation of biological and physical processes integrated in CLM that affect the net CH4 surface flux (Riley et al. 2011a). (left) Fully inundated portion of a CLM gridcell and (right) variably saturated portion of a gridcell.

2.25.3. CH4 Production

Because CLM does not currently specifically represent wetland plant functional types or soil biogeochemical processes, we used gridcell-averaged decomposition rates as proxies. Thus, the upland (default) heterotrophic respiration is used to estimate the wetland decomposition rate after first dividing off the O2 limitation. The O2 consumption associated with anaerobic decomposition is then set to the unlimited version so that it will be reduced appropriately during O2 competition. CH4 production at each soil level in the anaerobic portion (i.e., below the water table) of the column is related to the gridcell estimate of heterotrophic respiration from soil and litter (RH; mol C m-2 s-1) corrected for its soil temperature (\({T}_{s}\)) dependence, soil temperature through a \({A}_{10}\) factor (\(f_{T}\)), pH (\(f_{pH}\)), redox potential (\(f_{pE}\)), and a factor accounting for the seasonal inundation fraction (S, described below):

(2.25.2)\[P=R_{H} f_{CH_{4} } f_{T} f_{pH} f_{pE} S.\]

Here, \(f_{CH_{4} }\) is the baseline ratio between CO2 and CH4 production (all parameters values are given in Table 2.25.1). Currently, \(f_{CH_{4} }\) is modified to account for our assumptions that methanogens may have a higher Q\({}_{10}\) than aerobic decomposers; are not N limited; and do not have a low-moisture limitation.

When the single BGC soil level is used in CLM (Chapter 2.21), the temperature factor, \(f_{T}\), is set to 0 for temperatures equal to or below freezing, even though CLM allows heterotrophic respiration below freezing. However, if the vertically resolved BGC soil column is used, CH4 production continues below freezing because liquid water stress limits decomposition. The base temperature for the \({Q}_{10}\) factor, \({T}_{B}\), is 22°C and effectively modified the base \(f_{CH_{4}}\) value.

For the single-layer BGC version, \({R}_{H}\) is distributed among soil levels by assuming that 50% is associated with the roots (using the CLM PFT-specific rooting distribution) and the rest is evenly divided among the top 0.28 m of soil (to be consistent with CLM’s soil decomposition algorithm). For the vertically resolved BGC version, the prognosed distribution of \({R}_{H}\) is used to estimate CH4 production.

The factor \(f_{pH}\) is nominally set to 1, although a static spatial map of pH can be used to determine this factor (Dunfield et al. 1993) by applying:

(2.25.3)\[f_{pH} =10^{-0.2235pH^{2} +2.7727pH-8.6} .\]

The \(f_{pE}\) factor assumes that alternative electron acceptors are reduced with an e-folding time of 30 days after inundation. The default version of the model applies this factor to horizontal changes in inundated area but not to vertical changes in the water table depth in the upland fraction of the gridcell. We consider both \(f_{pH}\) and \(f_{pE}\) to be poorly constrained in the model and identify these controllers as important areas for model improvement.

As a non-default option to account for CH4 production in anoxic microsites above the water table, we apply the Arah and Stephen (1998) estimate of anaerobic fraction:

(2.25.4)\[\varphi =\frac{1}{1+\eta C_{O_{2} } } .\]

Here, \(\varphi\) is the factor by which production is inhibited above the water table (compared to production as calculated in equation (2.25.2), \(C_{O_{2}}\) (mol m-3) is the bulk soil oxygen concentration, and \(\eta\) = 400 mol m-3.

The O2 required to facilitate the vertically resolved heterotrophic decomposition and root respiration is estimated assuming 1 mol O2 is required per mol CO2 produced. The model also calculates the O2 required during nitrification, and the total O2 demand is used in the O2 mass balance solution.

Table 2.25.1 Parameter descriptions and sensitivity analysis ranges applied in the methane model

Mechanism

Parameter

Baseline Value

Range for Sensitivity Analysis

Units

Description

Production

\({Q}_{10}\)

2

1.5 – 4

CH4 production \({Q}_{10}\)

\(f_{pH}\)

1

On, off

Impact of pH on CH4 production

\(f_{pE}\)

1

On, off

Impact of redox potential on CH4 production

S

Varies

NA

Seasonal inundation factor

\(\beta\)

0.2

NA

Effect of anoxia on decomposition rate (used to calculate S only)

\(f_{CH_{4} }\)

0.2

NA

Ratio between CH4 and CO2 production below the water table

Ebullition

\({C}_{e,max}\)

0.15

NA

mol m-3

CH4 concentration to start ebullition

\({C}_{e,min}\)

0.15

NA

CH4 concentration to end ebullition

Diffusion

\(f_{D_{0} }\)

1

1, 10

m2 s-1

Diffusion coefficient multiplier (Table 24.2)

Aerenchyma

p

0.3

NA

Grass aerenchyma porosity

R

2.9\(\times\)10-3 m

NA

m

Aerenchyma radius

\({r}_{L}\)

3

NA

Root length to depth ratio

\({F}_{a}\)

1

0.5 – 1.5

Aerenchyma conductance multiplier

Oxidation

\(K_{CH_{4} }\)

5 x 10-3

5\(\times\)10\({}^{-4}\)\({}_{ }\)- 5\(\times\)10-2

mol m-3

CH4 half-saturation oxidation coefficient (wetlands)

\(K_{O_{2} }\)

2 x 10-2

2\(\times\)10-3 - 2\(\times\)10-1

mol m-3

O2 half-saturation oxidation coefficient

\(R_{o,\max }\)

1.25 x 10\({}^{-5}\)

1.25\(\times\)10\({}^{-6}\) - 1.25\(\times\)10\({}^{-4}\)

mol m-3 s-1

Maximum oxidation rate (wetlands)

2.25.4. Ebullition

Briefly, the simulated aqueous CH4 concentration in each soil level is used to estimate the expected equilibrium gaseous partial pressure (\(C_{e}\) ), as a function of temperature and depth below the water table, by first estimating the Henry’s law partitioning coefficient (\(k_{h}^{C}\) ) by the method described in Wania et al. (2010):

(2.25.5)\[\log \left(\frac{1}{k_{H} } \right)=\log k_{H}^{s} -\frac{1}{C_{H} } \left(\frac{1}{T} -\frac{1}{T^{s} } \right)\]
(2.25.6)\[k_{h}^{C} =Tk_{H} R_{g}\]
(2.25.7)\[C_{e} =\frac{C_{w} R_{g} T}{\theta _{s} k_{H}^{C} p}\]

where \(C_{H}\) is a constant, \(R_{g}\) is the universal gas constant, \(k_{H}^{s}\) is Henry’s law partitioning coefficient at standard temperature (\(T^{s}\) ),\(C_{w}\) is local aqueous CH4 concentration, and p is pressure.

The local pressure is calculated as the sum of the ambient pressure, water pressure down to the local depth, and pressure from surface ponding (if applicable). When the CH4 partial pressure exceeds 15% of the local pressure (Baird et al. 2004; Strack et al. 2006; Wania et al. 2010), bubbling occurs to remove CH4 to below this value, modified by the fraction of CH4 in the bubbles [taken as 57%; (Kellner et al. 2006; Wania et al. 2010)]. Bubbles are immediately added to the surface flux for saturated columns and are placed immediately above the water table interface in unsaturated columns.

2.25.5. Aerenchyma Transport

Aerenchyma transport is modeled in CLM as gaseous diffusion driven by a concentration gradient between the specific soil layer and the atmosphere and, if specified, by vertical advection with the transpiration stream. There is evidence that pressure driven flow can also occur, but we did not include that mechanism in the current model.

The diffusive transport through aerenchyma (A, mol m-2 s-1) from each soil layer is represented in the model as:

(2.25.8)\[A=\frac{C\left(z\right)-C_{a} }{{\raise0.7ex\hbox{$ r_{L} z $}\!\mathord{\left/ {\vphantom {r_{L} z D}} \right.}\!\lower0.7ex\hbox{$ D $}} +r_{a} } pT\rho _{r} ,\]

where D is the free-air gas diffusion coefficient (m2 s-1); C(z) (mol m-3) is the gaseous concentration at depth z (m); \(r_{L}\) is the ratio of root length to depth; p is the porosity (-); T is specific aerenchyma area (m2 m-2); \({r}_{a}\) is the aerodynamic resistance between the surface and the atmospheric reference height (s m-1); and \(\rho _{r}\) is the rooting density as a function of depth (-). The gaseous concentration is calculated with Henry’s law as described in equation (2.25.7).

Based on the ranges reported in Colmer (2003), we have chosen baseline aerenchyma porosity values of 0.3 for grass and crop PFTs and 0.1 for tree and shrub PFTs:

(2.25.9)\[T=\frac{4 f_{N} N_{a}}{0.22} \pi R^{2} .\]

Here \(N_{a}\) is annual net primary production (NPP, mol m-2 s-1); R is the aerenchyma radius (2.9 \(\times\)10-3 m); \({f}_{N}\) is the belowground fraction of annual NPP; and the 0.22 factor represents the amount of C per tiller. O2 can also diffuse in from the atmosphere to the soil layer via the reverse of the same pathway, with the same representation as Equation (2.25.8) but with the gas diffusivity of oxygen.

CLM also simulates the direct emission of CH4 from leaves to the atmosphere via transpiration of dissolved methane. We calculate this flux (\(F_{CH_{4} -T}\); mol m\({}^{-}\)2 s-1) using the simulated soil water methane concentration (\(C_{CH_{4},j}\) (mol m-3)) in each soil layer j and the CLM predicted transpiration (\(F_{T}\) ) for each PFT, assuming that no methane was oxidized inside the plant tissue:

(2.25.10)\[F_{CH_{4} -T} =\sum _{j}\rho _{r,j} F_{T} C_{CH_{4} ,j} .\]

2.25.6. CH4 Oxidation

CLM represents CH4 oxidation with double Michaelis-Menten kinetics (Arah and Stephen 1998; Segers 1998), dependent on both the gaseous CH4 and O2 concentrations:

(2.25.11)\[R_{oxic} =R_{o,\max } \left[\frac{C_{CH_{4} } }{K_{CH_{4} } +C_{CH_{4} } } \right]\left[\frac{C_{O_{2} } }{K_{O_{2} } +C_{O_{2} } } \right]Q_{10} F_{\vartheta }\]

where \(K_{CH_{4} }\) and \(K_{O_{2} }\) are the half saturation coefficients (mol m-3) with respect to CH4 and O2 concentrations, respectively; \(R_{o,\max }\) is the maximum oxidation rate (mol m-3 s-1); and \({Q}_{10}\) specifies the temperature dependence of the reaction with a base temperature set to 12 °C. The soil moisture limitation factor \(F_{\theta }\) is applied above the water table to represent water stress for methanotrophs. Based on the data in Schnell and King (1996), we take \(F_{\theta } = {e}^{-P/{P}_{c}}\), where P is the soil moisture potential and \({P}_{c} = -2.4 \times {10}^{5}\) mm.

2.25.7. Reactive Transport Solution

The solution to equation (2.25.11) is solved in several sequential steps: resolve competition for CH4 and O2 (section 2.25.7.1); add the ebullition flux into the layer directly above the water table or into the atmosphere; calculate the overall CH4 or O2 source term based on production, aerenchyma transport, ebullition, and oxidation; establish boundary conditions, including surface conductance to account for snow, ponding, and turbulent conductances and bottom flux condition (section 2.25.7.2); calculate diffusivity (section 2.25.7.3); and solve the resulting mass balance using a tridiagonal solver (section 2.25.7.5).

2.25.7.1. Competition for CH4 and O2

For each time step, the unlimited CH4 and O2 demands in each model depth interval are computed. If the total demand over a time step for one of the species exceeds the amount available in a particular control volume, the demand from each process associated with the sink is scaled by the fraction required to ensure non-negative concentrations. Since the methanotrophs are limited by both CH4 and O2, the stricter limitation is applied to methanotroph oxidation, and then the limitations are scaled back for the other processes. The competition is designed so that the sinks must not exceed the available concentration over the time step, and if any limitation exists, the sinks must sum to this value. Because the sinks are calculated explicitly while the transport is semi-implicit, negative concentrations can occur after the tridiagonal solution. When this condition occurs for O2, the concentrations are reset to zero; if it occurs for CH4, the surface flux is adjusted and the concentration is set to zero if the adjustment is not too large.

2.25.7.2. CH4 and O2 Source Terms

The overall CH4 net source term consists of production, oxidation at the base of aerenchyma, transport through aerenchyma, methanotrophic oxidation, and ebullition (either to the control volume above the water table if unsaturated or directly to the atmosphere if saturated). For O2 below the top control volume, the net source term consists of O2 losses from methanotrophy, SOM decomposition, and autotrophic respiration, and an O2 source through aerenchyma.

2.25.7.3. Aqueous and Gaseous Diffusion

For gaseous diffusion, we adopted the temperature dependence of molecular free-air diffusion coefficients (\({D}_{0}\) (m2 s-1)) as described by Lerman (1979) and applied by Wania et al. (2010) (Table 2.25.2).

Table 2.25.2 Temperature dependence of aqueous and gaseous diffusion coefficients for CH4 and O2

\({D}_{0}\) (cm2 s-1)

CH4

O2

Aqueous

0.9798 + 0.02986T + 0.0004381T2

1.172+ 0.03443T + 0.0005048T2

Gaseous

0.1875 + 0.0013T

0.1759 + 0.00117T

Gaseous diffusivity in soils also depends on the molecular diffusivity, soil structure, porosity, and organic matter content. Moldrup et al. (2003), using observations across a range of unsaturated mineral soils, showed that the relationship between effective diffusivity (\(D_{e}\) (m2 s-1)) and soil properties can be represented as:

(2.25.12)\[D_{e} =D_{0} \theta _{a}^{2} \left(\frac{\theta _{a} }{\theta _{s} } \right)^{{\raise0.7ex\hbox{$ 3 $}\!\mathord{\left/ {\vphantom {3 b}} \right.}\!\lower0.7ex\hbox{$ b $}} } ,\]

where \(\theta _{a}\) and \(\theta _{s}\) are the air-filled and total (saturated water-filled) porosities (-), respectively, and b is the slope of the water retention curve (-). However, Iiyama and Hasegawa (2005) have shown that the original Millington-Quirk (Millington and Quirk 1961) relationship matched measurements more closely in unsaturated peat soils:

(2.25.13)\[D_{e} =D_{0} \frac{\theta _{a} ^{{\raise0.7ex\hbox{$ 10 $}\!\mathord{\left/ {\vphantom {10 3}} \right.}\!\lower0.7ex\hbox{$ 3 $}} } }{\theta _{s} ^{2} }\]

In CLM, we applied equation (2.25.12) for soils with zero organic matter content and equation (2.25.13) for soils with more than 130 kg m-3 organic matter content. A linear interpolation between these two limits is applied for soils with SOM content below 130 kg m-3. For aqueous diffusion in the saturated part of the soil column, we applied (Moldrup et al. (2003)):

(2.25.14)\[D_{e} =D_{0} \theta _{s} ^{2} .\]

To simplify the solution, we assumed that gaseous diffusion dominates above the water table interface and aqueous diffusion below the water table interface. Descriptions, baseline values, and dimensions for parameters specific to the CH4 model are given in Table 2.25.1. For freezing or frozen soils below the water table, diffusion is limited to the remaining liquid (CLM allows for some freezing point depression), and the diffusion coefficients are scaled by the volume-fraction of liquid. For unsaturated soils, Henry’s law equilibrium is assumed at the interface with the water table.

2.25.7.4. Boundary Conditions

We assume the CH4 and O2 surface fluxes can be calculated from an effective conductance and a gaseous concentration gradient between the atmospheric concentration and either the gaseous concentration in the first soil layer (unsaturated soils) or in equilibrium with the water (saturated soil\(w\left(C_{1}^{n} -C_{a} \right)\) and \(w\left(C_{1}^{n+1} -C_{a} \right)\) for the fully explicit and fully implicit cases, respectively (however, see Tang and Riley (2013) for a more complete representation of this process). Here, w is the surface boundary layer conductance as calculated in the existing CLM surface latent heat calculations. If the top layer is not fully saturated, the \(\frac{D_{m1} }{\Delta x_{m1} }\) term is replaced with a series combination: \(\left[\frac{1}{w} +\frac{\Delta x_{1} }{D_{1} } \right]^{-1}\), and if the top layer is saturated, this term is replaced with \(\left[\frac{K_{H} }{w} +\frac{\frac{1}{2} \Delta x_{1} }{D_{1} } \right]^{-1}\), where \({K}_{H}\) is the Henry’s law equilibrium constant.

When snow is present, a resistance is added to account for diffusion through the snow based on the Millington-Quirk expression (2.25.13) and CLM’s prediction of the liquid water, ice, and air fractions of each snow layer. When the soil is ponded, the diffusivity is assumed to be that of methane in pure water, and the resistance as the ratio of the ponding depth to diffusivity. The overall conductance is taken as the series combination of surface, snow, and ponding resistances. We assume a zero flux gradient at the bottom of the soil column.

2.25.7.5. Crank-Nicholson Solution

Equation (2.25.1) is solved using a Crank-Nicholson solution (Press et al., 1992), which combines fully explicit and implicit representations of the mass balance. The fully explicit decomposition of equation (2.25.1) can be written as

(2.25.15)\[\frac{R_{j}^{n+1} C_{j}^{n+1} -R_{j}^{n} C_{j}^{n} }{\Delta t} =\frac{1}{\Delta x_{j} } \left[\frac{D_{p1}^{n} }{\Delta x_{p1}^{} } \left(C_{j+1}^{n} -C_{j}^{n} \right)-\frac{D_{m1}^{n} }{\Delta x_{m1}^{} } \left(C_{j}^{n} -C_{j-1}^{n} \right)\right]+S_{j}^{n} ,\]

where j refers to the cell in the vertically discretized soil column (increasing downward), n refers to the current time step, \(\Delta\)t is the time step (s), p1 is j+½, m1 is j-½, and \(S_{j}^{n}\) is the net source at time step n and position j, i.e., \(S_{j}^{n} =P\left(j,n\right)-E\left(j,n\right)-A\left(j,n\right)-O\left(j,n\right)\). The diffusivity coefficients are calculated as harmonic means of values from the adjacent cells. Equation (2.25.15) is solved for gaseous and aqueous concentrations above and below the water table, respectively. The R term ensure the total mass balance in both phases is properly accounted for. An analogous relationship can be generated for the fully implicit case by replacing n by n+1 on the C and S terms of equation (2.25.15). Using an average of the fully implicit and fully explicit relationships gives:

(2.25.16)\[\begin{split}\begin{array}{l} {-\frac{1}{2\Delta x_{j} } \frac{D_{m1}^{} }{\Delta x_{m1}^{} } C_{j-1}^{n+1} +\left[\frac{R_{j}^{n+1} }{\Delta t} +\frac{1}{2\Delta x_{j} } \left(\frac{D_{p1}^{} }{\Delta x_{p1}^{} } +\frac{D_{m1}^{} }{\Delta x_{m1}^{} } \right)\right]C_{j}^{n+1} -\frac{1}{2\Delta x_{j} } \frac{D_{p1}^{} }{\Delta x_{p1}^{} } C_{j+1}^{n+1} =} \\ {\frac{R_{j}^{n} }{\Delta t} +\frac{1}{2\Delta x_{j} } \left[\frac{D_{p1}^{} }{\Delta x_{p1}^{} } \left(C_{j+1}^{n} -C_{j}^{n} \right)-\frac{D_{m1}^{} }{\Delta x_{m1}^{} } \left(C_{j}^{n} -C_{j-1}^{n} \right)\right]+\frac{1}{2} \left[S_{j}^{n} +S_{j}^{n+1} \right]} \end{array},\end{split}\]

Equation (2.25.16) is solved with a standard tridiagonal solver, i.e.:

(2.25.17)\[aC_{j-1}^{n+1} +bC_{j}^{n+1} +cC_{j+1}^{n+1} =r,\]

with coefficients specified in equation (2.25.16).

Two methane balance checks are performed at each timestep to insure that the diffusion solution and the time-varying aggregation over inundated and non-inundated areas strictly conserves methane molecules (except for production minus consumption) and carbon atoms.

2.25.7.6. Interface between water table and unsaturated zone

We assume Henry’s Law equilibrium at the interface between the saturated and unsaturated zone and constant flux from the soil element below the interface to the center of the soil element above the interface. In this case, the coefficients are the same as described above, except for the soil element above the interface:

\[\frac{D_{p1} }{\Delta x_{p1} } =\left[K_{H} \frac{\Delta x_{j} }{2D_{j} } +\frac{\Delta x_{j+1} }{2D_{j+1} } \right]^{-1}\]
\[b=\left[\frac{R_{j}^{n+1} }{\Delta t} +\frac{1}{2\Delta x_{j} } \left(K_{H} \frac{D_{p1}^{} }{\Delta x_{p1} } +\frac{D_{m1}^{} }{\Delta x_{m1} } \right)\right]\]
(2.25.18)\[r=\frac{R_{j}^{n} }{\Delta t} C_{j}^{n} +\frac{1}{2\Delta x_{j} } \left[\frac{D_{p1}^{} }{\Delta x_{p1} } \left(C_{j+1}^{n} -K_{H} C_{j}^{n} \right)-\frac{D_{m1}^{} }{\Delta x_{m1} } \left(C_{j}^{n} -C_{j-1}^{n} \right)\right]+\frac{1}{2} \left[S_{j}^{n} +S_{j}^{n+1} \right]\]

and the soil element below the interface:

\[\frac{D_{m1} }{\Delta x_{m1} } =\left[K_{H} \frac{\Delta x_{j-1} }{2D_{j-1} } +\frac{\Delta x_{j} }{2D_{j} } \right]^{-1}\]
\[a=-K_{H} \frac{1}{2\Delta x_{j} } \frac{D_{m1}^{} }{\Delta x_{m1} }\]
(2.25.19)\[r=\frac{R_{j}^{n} }{\Delta t} +C_{j}^{n} +\frac{1}{2\Delta x_{j} } \left[\frac{D_{p1}^{} }{\Delta x_{p1} } \left(C_{j+1}^{n} -C_{j}^{n} \right)-\frac{D_{m1}^{} }{\Delta x_{m1} } \left(C_{j}^{n} -K_{H} C_{j-1}^{n} \right)\right]+\frac{1}{2} \left[S_{j}^{n} +S_{j}^{n+1} \right]\]

2.25.8. Inundated Fraction Prediction

A simplified dynamic representation of spatial inundation based on recent work by Prigent et al. (2007) is used. Prigent et al. (2007) described a multi-satellite approach to estimate the global monthly inundated fraction (\({F}_{i}\)) over an equal area grid (0.25 \(\circ\) \(\times\)0.25\(\circ\) at the equator) from 1993 - 2000. They suggested that the IGBP estimate for inundation could be used as a measure of sensitivity of their detection approach at low inundation. We therefore used the sum of their satellite-derived \({F}_{i}\) and the constant IGBP estimate when it was less than 10% to perform a simple inversion for the inundated fraction for methane production (\({f}_{s}\)). The method optimized two parameters (\({fws}_{slope}\) and \({fws}_{intercept}\)) for each grid cell in a simple model based on simulated total water storage (\({TWS}\)):

(2.25.20)\[f_{s} =fws_{slope} TWS + fws_{intercept} .\]

These parameters were evaluated at the 0.5° resolution, and aggregated for coarser simulations. Ongoing work in the hydrology submodel of CLM may alleviate the need for this crude simplification of inundated fraction in future model versions.

2.25.9. Seasonal Inundation

A simple scaling factor is used to mimic the impact of seasonal inundation on CH4 production (see appendix B in Riley et al. (2011a) for a discussion of this simplified expression):

(2.25.21)\[S=\frac{\beta \left(f-\bar{f}\right)+\bar{f}}{f} ,S\le 1.\]

Here, f is the instantaneous inundated fraction, \(\bar{f}\) is the annual average inundated fraction (evaluated for the previous calendar year) weighted by heterotrophic respiration, and \(\beta\) is the anoxia factor that relates the fully anoxic decomposition rate to the fully oxygen-unlimited decomposition rate, all other conditions being equal.