Long-Term Impact of Cover Crop and Reduced Disturbance Tillage on Soil Pore Size and Soil Water Storage

Using laboratory measurements and numerical simulations, we studied the long-term impact of contrasting tillage and cover cropping systems on soil structure and soil hydraulic properties. Complete water retention and conductivity curves 10 for top (0 – 5 cm) and subsurface (20 – 25 cm) samples were characterized and contrasted. Plot-level properties of water storage and retention were evaluated using numerical simulations in HYDRUS-2D software. Soils under no-till (NT) and cover cropping (CC) systems showed an improved soil structure in terms of pore size distribution (PSD) and the hydraulic conductivity (K) under these systems led to increased infiltration rate and water retention. The conventional measurement of water content at field capacity (water content at -33 kPa suction) and the associated plant available water (PAW) showed that 15 NT and CC plots had lower water content at field capacity and lower PAW compared to standard-till (ST) and plots without cover crop (NO). The numerical simulations, however, showed that NT and CC plots have higher profile-level water storage (albeit marginal in magnitude) and water availability following irrigation. Because the numerical simulations consider retention and conductivity functions simultaneously and dynamically through time, they allow the capture of hydraulic properties that are arguably more relevant to crops. The changes in PSD, water conductivity, and water storage associated with 20 NT and CC systems observed in this study suggest that these systems are beneficial to general soil health and improve water retention at the plot scale. List of Acronyms and Symbols CC, Cover crop; HCF, Hydraulic conductivity function; NT, No-till; PAW, Plant available water content; PSD, Pore size distribution; ST, Standard-till; WRC, Water retention curve 25 https://doi.org/10.5194/soil-2021-41 Preprint. Discussion started: 11 June 2021 c © Author(s) 2021. CC BY 4.0 License.


Introduction
Improving soil health-the vitality of a soil in sustaining the socio-ecological functions of its enfolding land (Janzen et al., 30 2021)-is one of the main challenges of our time as we grapple with the demands of growing population and changing climate.
The tools at our disposal to achieve this goal in agricultural lands are collectively known as conservation agriculture practices.
Conservation agriculture is characterized by a combination of three linked principles: (1) reduced mechanical soil disturbance, (2) preservation of a permanent organic soil cover, and (3) diversification of crop species (Kassam et al., 2019;Li et al., 2018;Mitchell et al., 2019). The adoption of conservation agriculture has been growing worldwide at an increasing rate since the 35 1960s. Between 2008 and 2015, the global area under conservation agriculture grew by 69% to 180 M ha (Blanco-Canqui and Ruis, 2018;Kassam et al., 2019). In California's highly productive Central Valley region the cultivated area under conservation agriculture for tomato and corn production has increased from less than 5,000 ha in 2004 to over 140,000 ha in 2012 (Mitchell et al. (2016a).
Conservation agriculture promises two main categories of benefits to soil health and soil functions. First, conservation 40 agriculture, specifically, reduced tillage, eliminates the negative effects associated with standard (conventional) tillage (ST), including degradation of soil structure, erosion, loss of nutrients, and reduction in soil microbial diversity and soil organic matter (Lal et al., 2007;Zuber and Villamil, 2016). Second, conservation agriculture supports the development of healthy soils.
The soil type at the study site is a Panoche clay loam (fine-loamy, mixed, superactive, thermic Typic Haplocambids) which is representative of much of California's Central Valley. For the first 12 years of the conservation agriculture experiment 85 (between 2000 and 2012), tomato and cotton were grown in rotation, followed by a rotation of sorghum with garbanzo beans since 2012. All plots were irrigated by subsurface drip.
The cover crops were a mix of triticale (Triticosecale Wittm), cereal rye (Secale cereale L.), common vetch (Vicia sativa), radish (Raphanus sativus), and clover (Trifolium incarnatum) seeded in 20 cm rows at 89.2 kg ha-in late October. The cover crops are terminated in late March of the following year using a stalk chopper followed by disk incorporation in the ST system 90 or sprayed with a 2% solution application of glyphosate after chopping and left on the surface as a mulch in the NT systems.
Detailed descriptions of the study site and management have been published in previous works (Mitchell et al., 2017(Mitchell et al., , 2015Veenstra et al., 2006;Mitchell et al., 2016b)

Sampling
Sampling was done in mid-November 2017, months after tillage in the ST treatment plots to avoid the immediate effects of tillage since we were primarily interested in the long-term effects of the treatments. Tillage operations have a transitory effect on porosity and associated soil hydraulic properties as structures collapse, mainly driven by wetting and drying cycles post tillage Mapa et al., 1986). The immediate alterations of tillage on soil porosity and hydraulic properties have 100 been shown to diminish rapidly following only a few wetting and drying cycles (Strudley et al., 2008;Alletto et al., 2015;Green et al., 2003).
Undisturbed soil samples from the top (0 -5 cm) and subsurface (20 -25 cm) layers were collected carefully using a 250 cm 3 volume sampling ring (8 cm diameter by 5 cm height). The depths were chosen to correspond with the depth disturbed by disking to incorporate residue in the ST plots (i.e., 0 -20 cm depth) (Mitchell et al., 2015;Veenstra et al., 2006) and the deeper 105 layer. Samples were collected along the strip ridges within the plots away from the trafficked furrows but slightly off-center to avoid drip irrigation tubes that were buried at the center of ridges. A total of 32 samples were collected by taking one surface, and one subsurface sample from four of the eight treatment replicate plots. This resulted in four replicates of surface and subsurface samples per treatment. The samples were stored at 4 °C before laboratory analysis.

Laboratory measurements 110
To measure the long-term impact of NT and CC practices on soil structure, we measured the soil bulk density ( ), total porosity, pore size distribution (PSD), and soil hydraulic properties (water retention curve, WRC, and hydraulic conductivity function, HCF) of the soil cores.
The saturated hydraulic conductivity ( ) was measured using the falling-head method. For this method, soils were saturated by immersing sample cores in degassed, 0.01 M CaCl2 solution so that the water level was close to the rim. of the saturated 115 soil was then measured by the falling-head method using the KSAT instrument (METER Group, Inc., Munich, Germany) by allowing a 5 cm column of degassed, 0.01 M CaCl2 solution to flow through the soil core. The setup was such so that the flow direction was downward. Following the measurement, soil WRC, and HCF data were determined simultaneously using the evaporation method as developed for the HYPROP instrument (METER Group, Inc., Munich, Germany). The HYPROP simultaneously measures, at high frequency (10 min), suction inside the soil cores at two different depths along with weight 120 loss while saturated soil cores dry. This allows for the calculation of WRC, ( ), and HCF, ( ). Following the HYPROP measurements, soil water retention in the range from 10 3 to 10 6 cm was determined by using the WP4C instrument (Decagon Devices, Inc, Pullman, WA, USA). We use the conventional definition for field capacity ( ) and permanent wilting point ( ) as the volumetric water content with the corresponding volume of water retained in the soil at −33 kPa and −1,500 kPa suction, respectively. and 125 are approximations of water retained after internal drainage has ceased, and the soil water content limit at which plants cannot recover from turgidity, respectively (Hillel, 1998). We calculated plant available water (PAW) as the difference between and (i.e., = − ). In addition to the saturated hydraulic conductivity, we also calculated the unsaturated hydraulic conductivity near field capacity water content at -10 kPa.
Throughout this manuscript, the term water suction, ℎ, is used to represent the soil water matric potential, , such that ℎ = 130 − (cm).

Soil porosity determination
Total soil porosity ( ) was calculated as = 1 − / , where is the particle density of soil, taken as 2650 kg m -3 , and is the soil bulk density determined using the standard core method (Grossman and Reinsch, 2002).
The effective pore size distribution (PSD) was estimated from the slope of the WRC using the differential water capacity 135 (Klute, 1986). For this, the WRC-(ℎ)-was first transformed into a curve of effective saturation ( ) as a function of effective pore radius ( ), ( ). was calculated as = ( − )/( − ), where and are the saturated and residual volumetric water contents estimated from a bimodal constrained van Genuchten model fit (Durner, 1994) of measured WRC. The draining pore radius was approximated using the Young-Laplace equation (1): where [µm] is pore radius, ℎ [cm] is the suction, is the surface tension between water and air (0.0729 N m -1 ), is the 140 contact angle (assumed 0), is the density of water (1000 kg m -3 ), and is the acceleration due to gravity (9.81 m s -2 ). The PSD curves were then calculated as (2): where [-] is the density function of effective pore sizes. Prior to calculating PSD, the ( ) curve was fitted with a cubic smoothing spline to remove noise in the measurement data (Kastanek and Nielsen, 2001;Pires et al., 2008). For a deeper 145 insight, we divided pore sizes into four ranges: intra-microaggregates (<0.2 µm), intra-aggregates (0.2 -10 µm), small https://doi.org/10.5194/soil-2021-41 Preprint. Discussion started: 11 June 2021 c Author(s) 2021. CC BY 4.0 License. macropores (10 -50 µm), and large macropores (50 -1000 µm). These range categories allowed us to perform statistical comparisons on the relative abundance of the pore size ranges among the different treatments.

Soil water storage simulations 150
To measure the interactive impact of changes in WRC and HCF on profile water dynamics and storage, we conducted a numerical simulation of field irrigation. The fate of irrigation water applied on the different treatment plots was simulated in HYDRUS-2D software where water flow is modeled using a modified form of the Richards' equation (Equation 3) which incorporates a sink term to account for water uptake by plant roots (Simunek et al., 2012).
where [L 3 L -3 ] is the volumetric water content, [T] is time, [L] are the spatial coordinates, [LT -1 ] is the unsaturated 155 hydraulic conductivity, [-] are the components of a dimensionless anisotropy tensor, ℎ [L] is the pressure head, and [T -1 ] is the sink term representing the rate of water volume removed due to plant water uptake.
The domain was set up as an axisymmetric cylinder of 18 cm radius and 100 cm depth. Figure  This discretization mesh was refined to have more nodes around the emitter (0.5 cm spacing) and soil layer boundaries (1 cm spacing) to capture expected high rates of changes in soil moisture. The material distribution in terms of soil hydraulic properties was such that the top 0 -20 cm and the subsurface 20 -30 cm were that of those measured in this study (Section 3.3). Soil hydraulic properties for the bottom layers, 30 -60 and 60 -100 cm layers, were predicted from soil characteristics using Rosetta-H5 pedotransfer function (Schaap et al., 2001) and van Genuchten-Mualem (1980) hydraulic model. Soil 165 characteristics for these layers were based on soil properties of C1 and C2 soil horizons (41 -58 and 58 -91 cm depths, respectively) for Panoche soils, Pedon ID S1978CA029001 (National Cooperative Soil Survey, n.d.).
Subsurface irrigation emitter was represented with a sphere of 1 cm radius buried 10 cm below the surface. We simulated the fate of 4.8 cm depth equivalent irrigation applied at an emitter discharge rate of 0.61 l h -1 (0.60 cm h -1 equivalent irrigation depth) in each of the 16 sampled plots. 170 The entire domain surface area (1017.9 cm 2 ) was associated with transpiration and the root water uptake was modeled by the HYDRUS-2D default Feddes' parameters for a tomato plant. Root spatial distribution was implemented using Vrugt et al.

185
showing the finite element mesh, related boundary conditions, and potential root water uptake rate distribution.
The starting pressure head of the entire model domain was set to -1000 cm and simulation was initialized by a 14-week spinup period. The model was run recursively with 2.5 cm equivalent depth irrigation applied at the start of every week for 14 weeks after which the final simulation is run with 4.8 cm equivalent depth irrigation (at the rate of 0.6 cm h -1 for 8 h) application ( Figure 3). The amount of water retained in a given soil profile layer following irrigation is calculated as equivalent water 190 depth changes using Equation 5.
where Δ [L] is equivalent water depth retained in the soil profile hours after irrigation application, is the equivalent water depth in the soil profile hours after irrigation, and 0 is equivalent water depth immediately before irrigation application.

Statistical analysis
All quantitative results are expressed as means of four replicates ± standard error unless otherwise indicated. Differences in 200 means were tested by analysis of variance (ANOVA) and pairwise comparison of treatments done using Tukey's honest significant difference (HSD) test at p < 0.15 significance level unless otherwise stated (Least Significant Difference table are provided in Appendix A Table B1). Hydraulic conductivity values were log-transformed before statistical analysis to make their distribution more normal. The normality of the data and the homogeneity of variances were checked using Shapiro-Wilk's and Levene's tests, respectively. All statistical analyses were performed using R statistical software (R Core Team, 205 2019).

Results and Discussion
Water retention and conductivity properties were measured for each soil sample using the KSAT, HYPROP, and WP4C instruments. The saturated hydraulic conductivities were measured using KSAT, the unsaturated hydraulic conductivities and the WRC using HYPROP and the water retention at extreme dry range using WP4C instruments (

Pore size distribution 215
The mean soil PSDs for the different systems are shown in Figure 5 (A). PSD curves for the individual samples are provided in Figure A3. A wider spread of PSD implies a heterogeneous mix of pore sizes and is indicative of soil with a more developed structure. The maximum pore volume density for the top soils occurred between sizes 15 and 20 µm diameter pores with the exception for NT-CC soils which showed a bimodal distribution with maximum pore volume density around 4 and 518 µm (Table 1)  A unique observation is that the topsoil NT-CC has the widest spread of PSD, with statistically more proportion of the smaller and larger diameter pores (i.e., <0.2 µm and 50 -1000 µm, respectively, at p < 0.15) and a bimodal distribution which is not present in the other systems ( Figure 5B). Several studies have similarly observed an increase in the proportion of larger pores in NT systems (Tavares Filho and Tessier, 2009;Pires et al., 2017;Gao et al., 2019Gao et al., , 2017. The reason for the abundance of 230 small and large pores for the NT-CC systems suggests the formation of tightly packed aggregates with smaller pores and larger interaggregate pores between them. This would be consistent with results from a previous study of our site and others that found higher aggregate stability for the NT-CC systems (Mitchell et al., 2017;Gao et al., 2019). Greenland (1977) suggests soil pore size classification based on equivalent diameter into three groups as transmission (50 -500 µm), storage (0.5 -50 µm), and residual pores (< 0.5 µm). Larger transmission pores are important for infiltration, drainage, and aeration while 235 smaller storage pores are important in retaining water. Increased aeration of soil is beneficial for many soil processes including healthy soil organic matter cycling (Lehmann and Kleber, 2015;Janzen, 2015) and other biogeochemical processes (Ekschmitt et al., 2008;Schmidt et al., 2011). All the treatments had higher relative abundance of the larger macropores (50 -1000 µm) compared to ST-NO plots and ST-NO had the highest proportion ( Figure 5B).  For the subsurface soils, the combined effect of NT and CC increased the spread of PSD however, NT without CC showed a narrower PSD with the highest PSD mode and highest abundance of large macropores compared to other treatments. NT-CC plots showed a significantly higher proportion of intra-aggregate size pores and smaller (< 10 µm) at p < 0.15. Plant roots are important actors in soil structure development, they enhance aggregation by compacting soils through growth and exudation 245 of segmenting materials, and also fragmenting aggregates to create larger interaggregate pores (Jarvis, 2007;Angers and Caron, 1998). Given the reduced tillage in the NT plots, it could be that CC play a more critical role in forming more diverse aggregate sizes and wider PSD. The effect of the CC species should also be considered in this interpretation since it has recently been shown that the effect of CC on soil structure and porosity varies significantly with root morphology and architecture of the CC

Bulk density
The mean across all treatments for the top and subsurface layer soils was 1.19 and 1.46 g cm -3 , respectively (which is equivalent to total porosities of 55 and 45 percent). Among the treatments, there was no statistically significant difference in and total porosity at p < 0.05 but at lower confidence levels, the topsoils under NT-NO system had significantly higher compared to ST-NO (p = 0.078) and NT-CC(p = 0.141) ( Figure 6 and Table B1). One of the concerns of NT practice is that it 255 may lead to soil consolidation and an increase in compaction because of the lack of intensive tillage (Blanco-Canqui and Ruis, 2018; Moret and Arrúe, 2007). Compaction reduces soil pore volume and affects soil fertility by reducing water flow and aeration, which negatively affect soil biological activity and redox potential (Vereecken et al., 2016). Our findings show that continued long-term NT led to a slight increase in compaction, but only when practiced without CC. The changes in PSD (discussed in section 0), however, showed that the NT systems increased the PSD in a manner that suggested a better-developed 260 soil structure with primary and secondary structures.

Hydraulic conductivity 265
The CC treatments tended to have a greater impact on than the tillage treatment for the top layer soils (Figure 7). This is consistent with the increase in infiltration reported previously for our soils by Mitchell et al. (2017). They found that CC increased infiltration by 2.8 times. They suggest several possible explanations for this including increased slaking associated with ST, better formation of macropores, and better continuity of soil pores possibly due to better-established soil structure and biology (Pires et al., 2017;Schwen et al., 2011). The NT-CC systems showed higher than NT-NO (p = 0.011) while 270 The ST plots showed midway between the NT-NO and NT-CC plots. The fact that of NT-NO plots is lower even more than ST plots suggests that CC is even more important when NT is practiced to maintain larger transmission pores without tillage. The effect of CC on ST plots was small and not statistically significant. at 100 cm suction, (100 ), is controlled by smaller pores as opposed to . The NT plots had lower (100 ) compared to ST plots (Figure 7), which implies that when soils are unsaturated, the NT plots will lose water to deep drainage at a slower rate than ST plots. This could possibly 275 mitigate the impact of reduced in the NT plots and lead to an increase of water availability to plants this explanation is consistent with our results from the numerical simulation (see section 3.5).

Water retention
The NT treatments had lower compared to ST (Figure 8). and DeLonge, 2017) which showed that CC did not affect total porosity for treatments practiced longer than 7 years or clay contents > 25 % which match the parameters of our study site. Our findings also agree with the findings of Basche and DeLonge (2017) in terms of , they find that while long-term CC tends to increase , it actually tends to decrease it for soils with >25% clay. Our results showed that while this was the case with ST, it was not the case for NT. For the subsurface layer of NT treatments, was significantly lower for the NT-CC compared with NT-NO treatments.

Simulated water storage 305
The simulation results showed that the difference in soil water content between the different treatments is most distinct in the top 40 cm. Figure 9 shows the vertical distribution of soil moisture following the irrigation for selected times. The 2-dimension distribution of soil moisture is shown in Appendix A Figures A4 and A5. Throughout the dry down following irrigation, the CC plots maintain higher volumetric water content in the top 20 cm. In the underlying 20 -30 cm depth layers, however, the NT-CC plots maintain the lowest soil moisture. While the NT-NO plots maintain a moderate soil water content in the top 20 310 cm compared to the other treatments, these plots maintain the highest water content in the 20 -30 cm depth layers.

315
Changes in water storage over time following 4.8 cm equivalent depth irrigation (see Equation 5) are shown in Figure 10. The results show that immediately following the end of irrigation, the top 40 cm layers start to lose water (to evapotranspiration and drainage) while the deepest layer (60 -100 cm) continues to gain water past 5 days after irrigation.

320
The conventional measure of plant-available water storage ( = − ) relies only on the WRC. Since WRC is a description of soil water status at equilibrium, this measure of plant-available water does not account for the dynamic interactions of water retention and hydraulic conductivity (Twarakavi et al., 2009). An alternative measure of field capacity is the "dynamic field capacity" which can be defined as the amount of water maintained in the soil after excess gravitational 325 water is drained and the rate of downward movement is minimal (Veihmeyer and Hendrickson, 1931). This dynamic field https://doi.org/10.5194/soil-2021-41 Preprint. Discussion started: 11 June 2021 c Author(s) 2021. CC BY 4.0 License. capacity is commonly taken as the water content after three (or sometimes even five) days (Twarakavi et al., 2009;Assouline and Or, 2014). In our simulation, the rate of water drainage for the top and middle layers had significantly decreased after three days ( Figure 10).
Comparison of the treatment averages in volumetric water content and amount of water retained three days after irrigation 330 (that is the dynamic field capacity and water storage at time of field capacity) are shown in Figure 11. The magnitude of differences among all treatments were small but generally favored the NT and CC treatments. For the top 20 cm soils, the only statistically significant difference in change in water storage was between NT-CC and ST-NO plots (p = 0.12) with both the ST-CC and NT-NO showing intermediate storage between the two. For the top 20 cm soils, the water content at dynamic field capacity for the CC plots was higher than those for NO plots, with the ST-NO plots showing statistically significant (p < 0.09) 335 lower water content than both CC plots. At the 20 -40 cm depths, the only statistically significant difference is between NT-NO and NT-CC with the NT-NO plots holding the most amount of water while the NT-CC holds the least amount. Both the ST plots, with and without CC show no difference in water content or water storage change three days after irrigation. These findings of water content at field capacity contrast with the the PAW estimated from the conventional steady-state measures (see Figure 8) which showed that the ST plots, in general, had higher water contents at field capacity and higher 340 PAW. The dynamic water content at field capacity for the subsurface layers 20 -30 cm shows similarity with that of the conventional steady-state field capacity for 20 -25 cm soils in that the NT-CC plots have lower water contents compared to NT-NO (p < 0.06). The ST plots 20 -40 cm have water content at dynamic field capacity closer to that of NT-NO but not statistically different from that of NT-CC (at P < 0.15). Only the dynamic water storage and water contents at field capacity capture the interaction between water retention curve and hydraulic conductivity functions, therefore these measures likely 345 capture soil hydrology more accurately.

Conclusion
The long-term NT and CC practices had an impact on soil pore size distribution (PSD). The NT and CC practiced independently led to a moderate increase in PSD range and had small or no effect on the measured soil hydraulic properties and simulated water dynamics. On the plots where NT and CC are practiced together, the changes in soil structure and hydraulic properties were most pronounced. NT with CC led to the development of bimodal pore size distribution in the top (0 -5 cm) soils with 355 https://doi.org/10.5194/soil-2021-41 Preprint. Discussion started: 11 June 2021 c Author(s) 2021. CC BY 4.0 License. the modes of the PSD around 4 and 500 µm diameter sizes which are in the storage and transmission pore size categories, respectively. While ST is done to improve soil structure for crops and overcome the compaction of the topsoil layer, its effect is transitory. Our results suggest that in the longer term, NT and CC increase soil aggregation and the proportion of larger pores while also maintaining total porosity. CC practices increased the saturated hydraulic conductivity ( ), particularly when practiced in conjunction with NT practices. 360 For the top layer soils (0 -5 cm), the of the NT-CC soils was significantly higher (p = 0.01) than that of the NT-NO soils.
The of NT-CC subsurface layer (20 -25 cm) was significantly higher (p < 0.15) than all other systems.
The measured water retention suggested a decrease in soils' ability to store water. The NT-CC practices decreased the calculated plant-available water (PAW) and water content at field capacity ( ). While these steady-state measures of field capacity and PAW indicate soil's ability to store water, the dynamically simulated water storage in soils is the result of the 365 interaction between soil's water retention characteristics and its hydraulic conductivities. Both the water retention and conductivity were accounted for in the HYDRUS-2D irrigation simulation. The results showed that when both retention and conductivity properties are considered together in the simulation, the top layers of NT systems not only do not have a disadvantage but have a marginally increased ability to store water compared to ST plots.
The changes in PSD, water conductivity, and water storage associated with NT and CC systems observed in this study suggest 370 that these systems are beneficial to the general soil health water retention at the plot scale.

Data availability
The soil data used for this study are available at [will be deposited on Figshare.com upon acceptance of the manuscript for 395 publication].

Author contribution
All co-authors contributed to the study design. Sampling was carried out by SA and TA. Soil analysis and modeling were carried out by SA with supervision from TA. SA prepared the manuscript with contributions from all co-authors.

Acknowledgments 400
This work was made possible with support from the Conservation Agriculture Systems Project and the California Department of Water Resources.