Soil-Gas Diffusivity-Based Characterization of Variably Saturated Agricultural Topsoils

: Soil-gas diffusivity and its variation with soil moisture plays a fundamental role in diffusion-controlled migration of climate-impact gases from different terrestrial agroecosystems including cultivated soils and managed pasture systems. The wide contrast in soil texture and structure (e.g., density, soil aggregation) in agriculture topsoils (0–10 cm) makes it challenging for soil-gas diffusivity predictive models to make accurate predictions across different moisture conditions. This study characterized gas diffusivity and gas-phase tortuosity in soils sampled from managed pasture and cultivated sites in Sri Lanka at 0–10 cm depth, together with selected soil-gas diffusivity data from the literature. Soil-gas diffusivity was measured using a one-chamber diffusion apparatus using N 2 and O 2 as experimental gases. The measured diffusivity, together with literature data representing both intact and repacked soils, were tested against ﬁve existing widely known gas diffusivity predictive models. The tested models tended to mischaracterize the two-region behavior in some of the aggregated soils, suggesting the need of soil-speciﬁc diffusivity models to better describe gas diffusivity in agricultural soils. We suggested a new parametric two-region model, developed in line with literature-based models, to represent both unimodal and bimodal/two-region behavior of selected soils. The new model statistically outperformed the existing predictive models for both intact and repacked soils and, hence, demonstrated its applicability to better characterize site-speciﬁc greenhouse gas emissions under different soil water regimes.


Introduction
Terrestrial agroecosystems are predominant anthropogenic sources of methane (CH 4 ) and nitrous oxide (N 2 O), two high-potency greenhouse gases which entail inevitable climate footprints such as global warming and region climate shifts. Both CH 4 and N 2 O are strong greenhouse gases having a global warming potential 25-fold and 298-fold higher than that of CO 2, respectively, over a 100-year time horizon [1]. Among different cultivated systems, paddy ecosystems contribute to approximately 20% of human-induced CH 4 emissions, due mainly to the typical submerged conditions practiced in paddy cultivation which favor methanogenesis in the paddy rhizosphere [2]. Essentially, intrinsic soil properties such as soil organic carbon, total nitrogen, organic matter, texture, etc. significantly affect the soil microbial functions and, thereby, can alter gas emissions to a great extent [3], while the diversity of the soil microbial community also plays a crucial role in gas production [4]. In addition, extensive application of nitrogen (N) fertilizers in paddy landscapes favors enhanced denitrification-driven N 2 O emissions. Though paddy emissions occur mainly by ebullition under submerged conditions, paddy fields typically experience seasonal water fluctuations which enable the produced gases to migrate and emit largely through the soil. In fact, alternate wetting and draining (AWD) is an emerging water-saving strategy practiced in countries for paddy cultivation in which the mid-season draining causes the release of produced greenhouse gases, largely into the atmosphere [5][6][7].
Grazed pasture constitutes another prominent agricultural system which is responsible for nearly 28% of anthropogenic N 2 O in the atmosphere [8]. Pastoral N 2 O is primarily derived from applied N fertilizers and also from the manure of the grazing livestock and is produced via nitrification and denitrification in the presence of anaerobic conditions in pasture topsoil. Scattered patches of ruminant urine and feces are prominent sources of N that potentially create N 2 O hotspots in pasture topsoils. The changes in water regimes in the pasture subsurface upon irrigation or rainfall make a major impact on fate and transport of N 2 O in the subsurface.
Once produced, the gases must be transported through the subsurface via the topsoil and emit across the soil-atmosphere continuum in order to make an impactful climate effect. If the topsoil is sufficiently aerated, N 2 O will be oxidized and emitted into atmosphere as dinitrogen (N 2 ), an environmentally benign gas as compared to N 2 O, with no climate forcing effect. Similarly, the draining of paddy soils may aerate the topsoils and, therefore, result in the emission of CO 2 in lieu of CH 4 . Thus, the agricultural topsoils (0-10 cm) become critically important in regulating the climate in local, regional, and global contexts.
Migration of greenhouse gases in agricultural topsoils and their emission into the atmosphere are predominantly diffusion-controlled, particularly in the absence of windinduced pressure gradients. In addition, soil physical properties such as soil texture (soil type) and structure (e.g., soil density, aggregation) play undeniable roles in agricultural CH 4 and N 2 O emissions. In particular, soil aggregation and associated hierarchical structural development in agricultural topsoils have been frequently discussed in the literature [9][10][11][12]. Extensive application of agricultural machinery, high livestock stocking rates, and frequent tillage of agricultural soils are among major controls affecting soil functional structure in pasture and cultivated soils. Mitigation measures, therefore, are highly dependent on how well the diffusion-controlled gas transport processes in agricultural topsoils are understood and how accurately they can be accounted for in predictive numerical tools.
Gas diffusion in soil is commonly characterized by soil-gas diffusivity, D p /D o , where D p (m 3 soil air m −1 soil s −1 ) and D o (m 2 air s −1 ) are the soil-gas diffusion coefficients in soil and in free air, respectively. Measuring D p /D o is, however, experimentally intensive and instrumentally challenging due to the requirements of specific apparatus and tight control of initial and boundary conditions. The predictive gas diffusivity models offer a convenient alternative to predict soil-gas diffusivity from easily measurable properties such as air-filled porosity (ε) and total porosity (Φ). In complex geomedia such as wellstructured agricultural soils, however, the soil physical heterogeneity makes the application of developed models questionable, since widely used models may have not been validated against soils with distinct structural contrast. Therefore, a careful revisiting of existing predictive models, together with measured gas diffusivity data on a wide range of agricultural topsoils, is an important prerequisite in investigating their applicability and to make reliable emission estimates.
Generally, agricultural topsoils are considered to be well-structured aggregated soils, having both inter-aggregate pores (i.e., pores between the aggregates) and intra-aggregate pores (i.e., pores within the aggregates), resulting in a distinctive bimodal pore structure. The two pore regions are generally considered to be functionally analogous with respect to soil-moisture and soil-gas dynamics and, hence, it is common to describe them in terms of two additive mathematical expressions in modeling soil-moisture retention [13] as well as soil-gas diffusivity [9]. However, compaction-induced changes due to animal treading, mechanical implements, and tilling operations on agricultural fields often alter soil pore structure [14]. Compaction essentially reduces the macropore domains and increases the micropore domains of a soil, thus shifting the bimodal nature of the structured soil. Although the predictive models developed for soil-gas diffusivity in unimodal soils are abound in the literature, these models cannot be used directly for aggregated soils as they tend to mischaracterize the two distinct pore regions. Therefore, several models have been modified, extended, and reintroduced to predict the soil-gas diffusivity in well-structured aggregated soils [14,15] with additional model parameters to numerically differentiate the different pore domains. Consequently, most two-region models are descriptive in nature, and use few a priori measurements to inversely estimate the additional model parameters. The classical numerical models also account for the two structural status of measured soil samples, intact or structurally disturbed, which also require close attention when agricultural soils are investigated.
The main objective of this study is to characterize gas transport behavior in cultivated and pasture topsoils (0-10 cm), including both structurally intact and disturbed soils, and also to investigate the applicability of widely known soil-gas diffusivity predictive models against these structurally diverse topsoils. A series of diffusivity measurements on selected cultivated and pasture soils as well as literature data were used to characterize the cultivated and pasture topsoils (0-10 cm). We used measured data from Sri Lankan cultivated and pasture sites, together with additional supporting data from the literature, to represent a wide geographic origin including soils from Japan (Nishi-Tokyo) and the United Kingdom (Lexington). Five widely applied soil-gas diffusivity models were tested against diffusivity data, along with a multi-parameter two-region descriptive gas diffusivity model to better characterize the two-region behavior of topsoils.

Soils and Data
In total, 12 soils were considered, with 6 soils representing each intact and repacked soil group. Each group included three pasture soils and three cultivated soils. Undisturbed soils were sampled from Ambewela pasture site (6.8693 • N, 80.7957 • E), Sri Lanka, at a 0-10 cm depth from five locations, along a downgradient transect, which also showed a natural organic matter gradient. Surface cracks, visible to the naked eye, could be seen on the topsoil of the sampling site and retrieved samples, likely due to the extensive presence of agricultural livestock. In addition, three Sri Lankan cultivated (paddy) soils, two intact soils from Kurunegala (7.5255 • N, 80.4390 • E) and Polonnaruwa (7.9400 • N, 81.0174 • E), and one repacked soil from Kandy (7.2794 • N, 80.5910 • E) were also used in this study. All paddy soils had been mechanically tilled before plantation, and samples were taken from a 0-10 cm depth.
In addition to the unpublished data mentioned above, literature data on pasture and cultivated soils, both in intact and disturbed status, were used in this study. Table 1 shows all 12 considered soils from Sri Lanka (Ambewela, Kurunegala, Polonnaruwa, Kandy, and Peradeniya), Japan (Nishi-Tokyo), and the United Kingdom (Lexington), with the names denoting the sampling location. Textural contrast and important soil physical properties of selected intact and repacked soils are also given in Table 1. For further details on the data from the literature, the reader is referred to the related literature mentioned in Table 1.

Methods
The 100 cm 3 undisturbed soil samples (height = 3.95 cm, diameter = 5.70 cm) were collected from Sri Lankan sites using metallic annular cores, ensuring minimum disturbance to the soil. The retrieved samples, leveled and kneaded at both ends to remove redundant soil, were end-capped and wrapped with polythene to preserve ambient moisture before being transferred to the Geotechnical Engineering laboratory for characterization. Repacked samples, on the other hand, were prepared from bulk soil sampled from a 0-10 cm depth, sieved to obtain 2 mm fraction, and repacked to the same bulk density observed in the corresponding layer.
The samples were first characterized for soil-moisture content, particle density, bulk density, particle size distribution (sieve analysis and hydrometer method), and organic matter content (loss on ignition method). The samples were then saturated for 72 h and sequentially drained to the intended moisture levels by stepwise evaporation. At each moisture level, the samples were kept closed for a sufficient amount of time to reach the hydraulic equilibrium and moisture redistribution before diffusivity measurements.
Soil-gas diffusivity was measured following the one-chamber method introduced by Taylor [18] and developed by Schjønning [19]. The custom-fabricated PVC chamber is 20 cm in height and 3.4 cm in internal diameter and is provisioned with two values as inlet and outlet for priming. The chamber was checked for airtightness before the measurement campaign was initiated. For the measurements, the sample was mounted on top of the diffusion chamber and made it airtight. Then, the chamber was flushed with 99.99% N 2 gas to remove all O 2 inside the chamber. The sample was then opened to the atmosphere by allowing the atmospheric O 2 to diffuse through the sample into the chamber. The increase in O 2 concentration inside the chamber was monitored continually with an O 2 sensor attached to the chamber wall. Calculation of D p /D o was performed following both Taylor (1949) [18] and Currie (1960) [20] methods as follows: The Taylor (1949) method is founded on Fick's first law and can be presented by the following equation: where ∆C t is the change in concentration inside the chamber (gm −3 ) (∆C t = C t=t − C t=0 ), H s is the sample height (m), and H c is the height of the chamber (m). D p can be calculated from the gradient of the graph of ln (∆C t /∆C 0 ) vs. time (t).
The Currie (1960) method, on the other hand, is founded on both Fick's first law and the second law and the basic equation can be expressed as follows: where ∆C t is the change in concentration inside the chamber ( , and H c is the height of the chamber (m). D p (soil-gas diffusion) can be derived from the slope of the plot of ln ( ∆C t ∆C 0 ) versus time (t). Table 2 shows an array of unimodal (one-region) predictive models for soil-gas diffusivity that were tested against the measured diffusivity data mentioned in Table 1. Table 2. Selected predictive soil-gas diffusivity models used in this study.

Bimodal (Two-Region) D p /D o Models
The two-region soil-gas diffusivity models are sparse in the literature, and typically consist of two mathematical expressions describing D p /D o in two regions. For example, the recent two-region model from Jayaratne et al. [12] takes the form of: For Region 1: For Region 2: where α 1 , β 1 , and α 2 , β 2 (dimensionless) are parameters characterizing the Region 1 (interaggregate) and Region 2 (intraaggregate), respectively, while w is a weighing factor that numerically demarcates Region 1 and Region 2. In an extended formulation, Shanujah et al. [21] presented a two-region model incorporating an additional parameter ε p (cm 3 cm −3 ) , referred to as a gas percolation threshold, an air-filled porosity below which the gas diffusivity is negligibly small due to extremely high moisture-induced tortuosity in the soil-gas phase in Region 1. Thus, the model for Region 1 takes the form: For Region 1: while for the Region 2, it is the same as in Equation (4). In order to make the formulation more physically meaningful in the two regions, we present here a slightly modified version of the above two models as follows: For Region 1: For Region 2: where w o (dimensionless) is the fraction of total porosity that corresponds to the gas percolation threshold; other symbols denote the same as above. The model implicitly assumes that D p /D o in the intraaggregate pore region exhibits (or can be approximated to) a linear increase with air-filled porosity as assumed in most of the previously developed tworegion models. Note that (ε − ω o Φ)/(ωΦ − ω o Φ) in Equation (6) represents the effective air saturation at a given air-filled porosity (ε) within Region 1, whereas (ε − ωΦ)/(Φ − ωΦ) in Equation (7) represents the same at a given air-filled porosity within Region 2. Despite the different formulation, the new model (Equations (6) and (7)) is intended to perform at the same level as the model from the two previous models (Equations (3)- (5)); therefore, we only tested the new model performance in statistical model tests, as presented below. The gas-phase tortuosity (τ, dimensionless) can be described as the roundabout distance a gas molecule will travel in the gaseous phase between two predefined points which are a unit Euclidean distance apart. Based on this definition, τ can be deduced from soil-gas diffusivity as presented by Ball [22] as follows:

Statistical Analysis
Two statistical indices, RMSE and bias, were used to statistically analyze model performance and compare model predictions. The RMSE evaluates the overall fit of a model to the measured data.
The bias, on the other hand, evaluates whether a model overestimates (positive bias) or underestimates (negative bias) the observations.
where d i is the difference between the observed and predicted diffusivity values and n is the number of diffusivity measurements in a data set. The performance of the selected models and the two-region model against the measured gas diffusivity data for intact soils and repacked soils expressed in terms of RMSE and bias is shown in Tables 3 and 4, respectively. Table 3. Statistical analysis on model performance using RMSE (Equation (9)) and Bias (Equation (10)) for intact soils.  Table 4. Statistical analysis on model performance using RMSE (Equation (5)) and Bias (Equation (6)) for repacked soils.  Figure 1 shows the soil-gas diffusivity against air-filled porosity for six intact soils from Sri Lanka and the UK, representing pasture (Figure 1a-c) and cultivated (Figure 1d-f) soils, together with the two-region descriptive model (Equations (6) and (7)).  In selected topsoils, the aggregation seems to have been disrupted due to the frequent treading of livestock ruminants and application of agricultural machinery, leaving a nonaggregated structure, as we see for most of the soils in Figure 1. A distinct exception to this is Ambewela pasture soil, which showed a linear behavior in the D p /D o vs. ε relation, as shown in Figure 1a, likely due to the cracked nature of the soil formed from extensive animal treading and trampling, which has resulted in an altered soil structure. In fact, large cracks could be seen on the soil surface during sampling. Similar linear variations in cracked soils/fractured geomedia in characteristic diffusivity curves have been observed in literature [23], while the gradient of the linear model provides a useful implication of the average fracture angle of the geosystem.

Soil-Gas Diffusivity
In the Peradeniya pasture site, on the other hand, the use of livestock or agricultural implements was not extensive over a couple of years prior to sampling and, therefore, the soil has gradually developed into a well-structured aggregated soil. Consequently, Peradeniya soil has exhibited a distinct two-region behavior (Figure 1b). Note the decline in the rate of increase in gas diffusivity when the gases enter the intraaggregate pore space (region 2), likely due to the encounter of a more tortuous pore network within the aggregates. In fact, in strongly aggregated porous media (e.g., manufactured porous substrates), a distinct threshold in gas diffusivity has been observed in the literature [24] due to the presence of highly tortuous pore network inside aggregates. The UK pasture soil also showed a slight two-region behavior, owing to its weakly aggregated nature, but was not so pronounced as the Peradeniya soil. On the other hand, the three cultivated soils, Kurunegala, Peradeniya, and Polonnaruwa (Sri Lanka), which have undergone a regular mechanical tillage before cultivation, showed unaggregated nature and, hence, yielded a typical unimodal behavior in gas diffusivity. Figure 2 shows the measured soil-gas diffusivity against air-filled porosity, together with the descriptive two-region model (Equations (6) and (7)) for six repacked soils: three pasture soils (two Peradeniya soils and the Nishi-Tokyo soil) and three cultivated soils (Peradeniya, Nishi-Tokyo, and Kandy soils). Interestingly, all six repacked soils showed a two-region behavior. Typically, since the sieving and repacking involve de-structuring of soil aggregates in repacked soils, they are likely to show fewer two-region characteristics than intact soils. Note that the repacking did not cause the crushing of aggregates; hence, the aggregated behavior of sieved soils was preserved in repacked samples. In comparison, the three cultivated soils showed a higher gas percolation threshold compared to the three pasture soils. In Peradeniya repacked soils, the less-compacted soils (1.0 g cm −3 , Figure 2a) have higher interaggregate pore volume than the highly compacted soils (1.30 g cm −3 , Figure 2b), implying a decrease in interaggregate pore space due to high compaction, with an overall decrease in D p /D o .   (6) and (7)). PDY; Peradeniya (Sri Lanka), NTK; Nishi-Tokyo (Japan), KDY; Kandy (Sri Lanka). P and C in parentheses denote Pasture and Cultivated, respectively.
To examine the performance of existing widely used predictive models on agricultural topsoils, we tested five selected predictive models, which comprise both classical and recent models, against the measured D p /D o for undisturbed topsoils (Figure 3). Evidently, excluding the Millington (1959) model [25], all other models tended to markedly underpredict diffusivity in undisturbed soils. In fact, predictions were better for cultivated soils than for pasture soils due to the altered structural properties in pasture soils, as discussed above. The Millington and Quirk (1961) model [26], the widely accepted model for predicting soil-gas diffusivity in gas transport models, tended to overpredict the D p /D o in dry soils and markedly underpredict the D p /D o in wet soils, as reported in many previous works in the literature [17,27,28].  The WLR model, developed specifically for sieved and repacked soils, overpredicted the cultivated intact soils as expected, but underpredicted the undisturbed pasture soils due to the enhanced soil structure leading to pronounced diffusion properties. In the SWLR model, in general, a complexity factor smaller (C m = 1.5) than originally proposed (C m = 2.1) was required to better characterize the pasture topsoils. In Figure 3f, the behavior of the two-region descriptive model (Equations (6) and (7)) shows the degree of improved prediction that can be achieved if a two-region model is used in place of the classical unimodal models.
A similar scatterplot comparison was also made for repacked soils, as depicted in Figure 4. As expected, both WLR [27] and SWLR (with C m = 1.0) [28] models performed better due to their intended applicability for repacked soils, and the Buckingham (1904) model [29] also provided promising results. Furthermore, as shown in Figure 4, the behavior of the two-region descriptive model (Equations (6) and (7)) demonstrates the improved predictability when optimized model parameters are used.  Overall, it is clear that the unimodal models apparently mischaracterized the agriculture topsoils with altered soil structure, leading to a significant disparity in measured and modeled diffusivity data for both intact and repacked soils. Application of a two-region descriptive model, with a few a priori diffusivity measurements and calibrated model parameters, can be used to better characterize agricultural topsoils. Such a descriptive model calibrated with measured data specifically for a given soil can be particularly useful when site-specific agricultural management practices are adopted to enhance crop productivity and minimize agricultural emissions. When access to diffusivity measurements is lacking, our analysis shows that the Buckingham model, albeit the oldest and the simplest among the classical models, can be used for reasonably good predictions for agricultural topsoils.
The parameterization of the two-region model (Equations (6) and (7)) is given in Table 5 for both intact and repacked soils. Excluding Kurunegala soil, all other soils have a percolation threshold less than 10% of the total pore space (i.e., w 0 < 0.1). Numerical analysis further shows that the interaggregate porosity in bimodal soils corresponds to 43-73% of the total porosity. Table 5. Numerical parameterization of two-region soil-gas diffusivity model (Equations (6) and (7)) for pasture and cultivated soils.

Soil
Sample Status Vegetation Type Characterization of Two-Region Soil-Gas Diffusivity Model (Equations (6) and (7)

Gas-Phase Tortuosity
The calculated tortuosity of the functional soil gaseous phase from the measured soil-gas diffusivity and air-filled porosity data (Equation (8)) for the intact and repacked agricultural topsoils are shown in Figure 5, together with calculated tortuosity from the two-region model (Equations (6) and (7)). Note that the tortuosity is presented in decreasing order, which denotes the increasing air-filled porosity. Clearly, the tortuosity is high in the wet region due to the pronounced moisture-induced effects on pore network with constrained gas diffusion, which gradually decrease as the soil drains and becomes minimum in dry soil. The observed minimum tortuosity values (around 1.1 for both intact and repacked soils) are generally lower than those of typical mineral soils (around 1.5-2.0) due to the enhanced pore connectivity in aggregated soils. The uncertainty in tortuosity measurements and predictions is also high in the wet region due to the random connection and disconnection of water bridges across soil particles, yielding a large scatter in tortuosity values. As expected, the scatter in tortuosity is typically higher for intact soils than for the repacked soils due to the intrinsic soil structural effects in intact soils. Notably, the two-region model demonstrated its usefulness in describing the tortuosity in agricultural topsoils.
This study has provided useful experimental and numerical insight into characterization of agricultural topsoils based on measured and modeled soil-gas diffusivity data. It is worthy to be noted that, due to the wide structural diversity of agricultural topsoils, further measurements on both intact and repacked soils will enable a more detailed characterization of soils. It should be noted that all measurements involved have been carried out in laboratory-controlled environments where natural environmental complexities (e.g., dynamics in temperature, evaporation, wind, and humidity) were eliminated. Such additional environmental factors were out of the scope of this study but must be accounted for when making more realistic conclusions. Results, therefore, need to be compared against field-measured data with caution. Figure 5. Pore tortuosity analysis using the two-region descriptive model (Equations (6) and (7)

Conclusions
This study investigated new measurements and literature data on soil-gas diffusivity to characterize agricultural topsoils (0-10 cm) representing pasture and cultivated soil systems. The diffusivity data showed a wide structural diversity in both intact and structurally disturbed samples, with some soils showing distinct two-region behavior. A tworegion descriptive soil-gas diffusivity model, with model parameters calibrated against a few a priori measured data, seemed to be useful in more accurately characterizing the structure-dependent behavior the of considered agricultural topsoils when compared to the widely used predictive gas diffusivity models. The gas-phase pore tortuosity analysis showed pronounced moisture-induced pore tortuosity with a large scatter in wet soils, while the two-region model demonstrated an adequate description of the calculated tortuosity. The overall numerical analysis, thus, highlighted the importance of site-specific characterization of gas diffusivity and pore tortuosity in order to make better predictions of soil-gas emissions as a prerequisite to introduce agricultural management strategies in relation to mitigation of greenhouse gases. The laboratory-controlled measurements of this study must be compared cautiously with the results of field-based studies, which provide more realistic observations in the presence of other environmental factors.

Conclusions
This study investigated new measurements and literature data on soil-gas diffusivity to characterize agricultural topsoils (0-10 cm) representing pasture and cultivated soil systems. The diffusivity data showed a wide structural diversity in both intact and structurally disturbed samples, with some soils showing distinct two-region behavior. A two-region descriptive soil-gas diffusivity model, with model parameters calibrated against a few a priori measured data, seemed to be useful in more accurately characterizing the structuredependent behavior the of considered agricultural topsoils when compared to the widely used predictive gas diffusivity models. The gas-phase pore tortuosity analysis showed pronounced moisture-induced pore tortuosity with a large scatter in wet soils, while the tworegion model demonstrated an adequate description of the calculated tortuosity. The overall numerical analysis, thus, highlighted the importance of site-specific characterization of gas diffusivity and pore tortuosity in order to make better predictions of soil-gas emissions as a prerequisite to introduce agricultural management strategies in relation to mitigation of greenhouse gases. The laboratory-controlled measurements of this study must be compared cautiously with the results of field-based studies, which provide more realistic observations in the presence of other environmental factors.