Influence of soil suction on runoff during whiplash events

. The influence of soil suction on rainfall-runoff during whiplash events is explored. A series of 600 finite element simulations were completed using PLAXIS LE Groundwater. A model 2D soil profile was developed in the software to simulate a fat clay soil native to Oklahoma. Rainfall intensity, initial soil suction, and surface slope elevation were treated as log normally distributed random variables. Desiccation cracking was explicitly included in the model as a function of initial soil suction. The depth of the desiccation cracking was determined according to the soil’s tension capacity. the results indicate that rainfall-runoff is greatly influenced by soil suction. For a given rainfall intensity, the portion of the rain that becomes runoff varies according to soil suction; however, this relationship is not strongly correlated. The simulation results exhibit significant variability, which emphasize the complexity of this relationship. The occurrence of desiccation cracking reduced but did not eliminate the potential for rainfall-runoff for desiccation crack depths less than 5 m for the soil in this study.


Introduction
Flash flooding can occur rapidly with little or no warning. Accurate prediction of these events is complicated due to the complexity in predicting storm systems. They are further complicated by the variability of ground surfaces on which the rain falls. The amount of runoff that contributes to the flash flood is related to the features of the surfaces in contact with the rainwater.
Soil permeability in an unsaturated state can vary greatly from the saturated state. Soil permeability is a function of soil suction. A decrease in soil permeability near the surface could lead to an increase in runoff that was not previously accounted for. Unsaturated soils act more like an impervious surface if conditions prior to the rainfall event dry the soil.
Whiplash events are events such as severe droughts followed by pluvial, or heavy rainfall, events or vice versa [1 -2]. These events are becoming more common in a changing climate causing an increase in the land area covered by unsaturated soils prior to severe rainfall events.
There are ongoing efforts to better predict flash flooding events considering soil classification and moisture content (e.g., [3 -5]). However, the influence of soil suction is not currently considered in these predictions. It is well established that the relationship between moisture content and soil suction is soil dependent. There is a need to better understand the relationship between soil suction and runoff during rainfall events to predict the likelihood of flash floods more accurately. This paper presents a series of numerical simulations predicting rainfall runoff for unsaturated soils. The * Corresponding author: tbounds@ou.edu rainfall intensity, initial soil suction, and ground slope were varied in the study. Desiccation cracking was explicitly incorporated into the model and directly linked to the initial soil suction.

Methods
The analyses were conducted using PLAXIS LE Groundwater 2D [6]. The model considered was 10 m wide and 9 m in height. The duration of rainfall on the model was set to 300 minutes. Pore-water pressure was constant along the bottom of the model. A unit gradient boundary was used along vertical boundaries. The unit gradient boundary sets the flux out of the model equal to the unsaturated permeability. A unit gradient boundary allows flow out of the model only. The unit gradient boundary condition was found to represent an infinite plane most closely in the horizontal direction during preliminary analyses. The results show no apparent edge effects when the unit gradient boundary condition was used on the vertical boundaries for these analyses. The water table was located along the bottom boundary. The water table did not rise or fall during the analyses due to the bottom boundary conditions and the short rainfall duration.
The software calculates surface runoff parallel to the top model boundary. This means that runoff along the crack walls will contribute to the total runoff despite the water filling the cracks and not actually running off the surface. To avoid this shortcoming, the desiccation cracks were filled with a saturated sandy soil. The sandy soil was extended 15 mm above the model and a smooth surface was created as shown in Figure 1. The sandy soil surface was kept parallel to the clay soil surface. The permeability of the sandy soil was equal to the rainfall intensity. A series of preliminary analyses found the addition of the sand layer impacted the runoff by less than 5% as long as the permeability of the sand was equal to the rainfall intensity. Deviations between the sand permeability and rainfall intensity increased the error in rainfall runoff induced by the sand layer. The preliminary analyses did not consider the influence of desiccation cracking since these cannot be properly modelled in the software without the sand layer. The rainfall was applied to the sandy soil surface. Twentytwo desiccation cracks at approximately equal spacing were included in the model. The spacing of desiccation cracking was based on field observations for similar soils. One of the model layouts is shown in Figure 1.
The soil used in the analysis was collected near Idabel, Oklahoma. The soil has previously been used in other research at the University of Oklahoma. The soil is a highly plastic fat clay (CH) residual soil. The soil experiences desiccation cracking in situ and has significant shrink-swell potential. The liquid limit of the soil is 72% and the plasticity index is 46% as presented in [7]. The Soil Water Characteristic Curve (SWCC) was developed using a chilled mirror hygrometer (WP4) and a pressure plate device [8]. The pressure plate device was used to obtain the low suction range while the WP4 was used for the high suction range. Both wetting and drying SWCC's are shown in Figure 2. Within the software soil hydraulic conductivity is linked to soil suction according to the Modified Campbell Estimation [9]. The Modified Campbell Estimation is an iteration of the Campbell Equation [10] where the permeability levels off at a residual suction. The Modified Campbell Estimation is shown as follows: where k = hydraulic conductivity of the water phase; ks = saturated hydraulic conductivity of the water phase; kmin = calculated minimum hydraulic conductivity; p = parameter used to control the modified Campbell [10] estimation of hydraulic conductivity; af, nf, mf, hr = Fredlund and Xing [11] SWCC fitting parameter, ψ = soil suction. The relationship between soil suction and hydraulic conductivity for the Idabel soil is shown in Figure 3. The sharp hydraulic conductivity decrease is related to changes in the volume fraction and of the water phase which changes rapidly beyond the air entry suction. Rainfall intensity, initial soil suction, and ground slope were treated as random variables. All random variables were assumed to follow a logarithmic distribution so that negative values were not generated. The initial soil suction was incorporated into the model by applying different levels of solar radiation over a period of 10 days. The change in ground slope was incorporated by increasing the height of the left boundary by the randomly generated value. The parameters used to generate the random variables are shown in Table 1. A Cumulative Distribution Function (CDF) for each of the random variables is shown in Figure 4. The near surface suction as estimated by PLAXIS LE is included in Figure 4 for reference. The crack depth was estimated according to the initial soil suction profile as described by [12] and shown below: where d(ua -uw)c = change in suction to produce desiccation cracking, k0 = at rest earth pressure, γzc = overburden pressure, ν = Poisson's ratio, σt = soil tensile strength, E = modulus of elasticity with respect to net normal stress, and H = modulus of elasticity with respect to matric suction. The method relates the crack depth to the soils tension capacity. As soil suction increases desiccation cracks will advance as the soils internal tension capacity is exceeded. The tension for the study soil was determined using a novel soil tension crack measuring device [7]. The crack depth as a function of soil suction is shown graphically in Figure 5. A series of soil suction profiles were developed by varying the amount of solar radiation subjected to the model. The bracketed values indicate the solar radiation (MJ/m 2 /hr) for each soil suction profile. For example, the profile indicated by '+' symbols were generated using a constant solar radiation of 36 MJ/m 2 /hr for 10 days. To use this method a relationship between solar radiation and soil suction was established for the model. This relationship was then used to generate the crack depth according to the initial soil suction profile for the simulation. The maximum crack depth for the study was restricted to 5 m to avoid the influence of the groundwater table.
The relationship between desiccation crack depth and width has been found to vary greatly e.g., [13 -15]. Those researchers found desiccation crack width to range between 1 and 500 % of the desiccation crack depth. It's important to note that many of the cracks included in their studies were shallow. The most common ratio was approximately 3.5% for cracks up to 0.6 m in depth. For this study the crack depth was as much a 5 m, assuming a constant relationship of 3.5% would mean a desiccation crack width of approximately 140 mm. A linear relationship between desiccation crack depth and width seems unlikely when considering greater depths. There is likely an upper limit on desiccation crack width for a given soil. For these reasons desiccation crack width was related to the estimated crack depth using a simple power function for this study. It is assumed that desiccation crack width will increase rapidly for shallow cracks and become limited for greater crack depths. The relationship is presented in Equation 3 and shown graphically in Figure 6.
where Wc = desiccation crack width in mm and Dc = desiccation crack depth in mm. The rainfall was included as a step function meaning that the rainfall started and ended abruptly. The rainfall was applied to the top boundary of the model which coincided with the top of the saturated sandy soil layer. Runoff and infiltration were also recorded along the top boundary.

Results
The results from the simulation were collected and compared to better understand the relationship between soil suction and rainfall runoff for the study soil. A fair amount of variability exists in the results emphasizing the complexity of this problem. Many of the following comparisons consider the runoff as a portion of the total amount of rainfall introduced into the model. That is, a runoff percentage of 50 in the following charts would indicate that half of the water introduced as rainfall into the model became runoff and the other half infiltrated the soil.
The relationship between rainfall runoff and rainfall intensity is well established. However, for the simulations in this study there is a fair amount of scatter likely due to the varying permeability of the study soil. The inclusion of desiccation cracking also had an effect on the amount of rainfall runoff. A comparison between the runoff percentage and rainfall intensity is shown in Figure 7. The portion of rainfall that ultimately became runoff increases asymptotically towards total rainfall amount (i.e., 100%) as the rainfall intensity increases. The soil absorbed a portion of the rainfall for all simulations included in this study. It is noteworthy that for a given rainfall intensity the portion of rain that becomes runoff can vary greatly. The variation is partly due to changes in permeability because of soil suction variations. Many of the outliers in the data represent simulations with desiccation cracking. The focus of this study is the effect of soil suction on rainfall runoff. When comparing the two, there is a wide range of rainfall runoffs for a given suction. This relationship holds true even when the percentage of runoff is considered, see Figure 8. A line indicating the suction at which desiccation cracking starts is also included in the figure. The suction shown in the graph is represents the soil suction 0.75 m below the surface at the middle of the model. Prior to the formation of desiccation cracking the majority of the rainfall becomes runoff. There is a slight decrease in this relationship as the soil approaches desiccation cracking. Once desiccation cracking begins the runoff tends to decrease rapidly. The runoff becomes zero or near zero after a suction of approximately 450 to 500 kPa. This suction corresponds to a crack depth of 5 m since the desiccation crack depth was restricted. The inclusion of desiccation cracking allows for more surface area for water to absorb into the soil. The desiccation cracks also provide a pathway to higher permeability layers. The runoff can still be quite large at moderate to high suctions as shown in the figure. In these analyses desiccation cracks were uniformly distributed and systematically developed. Desiccation cracks do not uniformly develop in a systematic manner. Preferential path desiccation crack development is a more realistic approach. Desiccation cracks assuming this type of development are horizontal as well. Horizontal desiccation cracking would likely not have much of an impact on the results. Once desiccation crack development began the runoff decreased dramatically.
The results for runoff seem a little counterintuitive. It would be expected that runoff would increase as suction decreases prior to the formation of desiccation cracks. During this time the soil permeability is decreasing yet the soil surface is intact. This would result in a less permeable surface. However, for these simulations there is a general decrease in the rainfall runoff as suction increases. Some of this could be related to the soils water storage capacity. The total water storage capacity of each model was the same since the dimensions of the model and the water table were consistent. The actual water in the model decreased as solar radiation was applied to the model and the suction increased. Meaning there was more water storage capacity available in the model. The change in water storage from a model time of zero to right before rainfall began and the rainfall runoff are shown in Figure 9. Δ water storage increases as the soil suction increases. In other words, there is more space available for water in the model as the soil dries. The data in the figure agrees well with what is observed in Figure 8. Intuitively the soil permeability would control the infiltration rate. For this study the soil permeability does not change drastically for the suctions of interest (see Fig. 3). This could explain why the soils water storage capacity is governing the infiltration behaviour. In future research this relationship will be explored further. The surface slope was also treated as a random variable in these analyses. The elevation changes across the top of the model varied from 0 to 2.32 m. This correlates to a change in slope up to approximately 13º. A comparison between rainfall runoff and surface slope is shown in Figure 10. For the shallower slopes there doesn't appear to be any relationship and there is a lot of scatter in the data. For slopes less than approximately 6º the rainfall runoff is more controlled by the soil. However, beyond 6º the surface slope does appear to be impacting the results. Notice that the minimum rainfall runoff increases as the surface slope increases. Based on the data there appears to be a minimum rainfall runoff possible for surface slopes greater than 6º. However, there are not many data points for this region due to the sampling methods adopted in this study. In future research the sampling intervals will be expanded to include more simulations in this region.

Conclusions
This paper presented a numerical study investigating the effects of soil suction on rainfall runoff. PLAXIS LE Groundwater was used for the simulations. A fat clay with high shrink swell potential was used as the model soil in this study. Pertinent soil parameters for the soil were determined in the laboratory. The change in permeability as a function of soil suction was estimated using the Modified Campbell Estimation [9] within the software. Rainfall intensity, initial soil suction, and surface ground slope were treated as log normally distributed random variables. The initial soil suction was estimated by the software considering the variations in solar radiation. A larger constant solar radiation created a more unsaturated soil profile prior to the rainfall event. Desiccation cracking was explicitly incorporated into the model according to a relationship proposed by [12]. Using the randomly generated variables, 600 simulations were performed.
From this study if was found that the potential for rainfall runoff is greatly reduced once desiccation cracking begins. the desiccation cracking reduced but did not eliminate rainfall runoff until a crack depth of 5 m. Desiccation cracking is difficult to model and many unknowns remain. In future research the spatial distribution of desiccation cracking will need to be studied and incorporated. A simple relationship between desiccation crack depth and width was adopted for this research. A more robust relationship for desiccation crack width also needs to be developed. It was concluded that rainfall runoff varies greatly with soil suction. That is, the same soil with the same rainfall intensity and duration can have very different amounts of runoff due to varying soil suction. The water storage capacity was also found to have a profound effect on the amount of rainfall runoff. The surface slope appears to limit the amount of rainfall infiltration as the slope increases beyond 6º. The relationship between rainfall runoff and soil suction is complex and more research needs to be done to better understand the effects of the various variables involved. From this initial study it is evident that changes in permeability as a function of soil suction cannot fully explain the simulation results.