Next Article in Journal
A Review of the Application of Discrete Element Method in Agricultural Engineering: A Case Study of Soybean
Next Article in Special Issue
Special Issue on Numerical Modeling in Civil and Mining Geotechnical Engineering
Previous Article in Journal
Effects of Water Content and Irrigation of Packing Materials on the Performance of Biofilters and Biotrickling Filters: A Review
Previous Article in Special Issue
Implementation of the Non-Associated Elastoplastic MSDPu Model in FLAC3D and Application for Stress Analysis of Backfilled Stopes
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Unsaturated Hydraulic Conductivity Estimation—A Case Study Modelling the Soil-Atmospheric Boundary Interaction

1
UniSA STEM, University of South Australia, Mawson Lakes Campus, Mawson Lakes, SA 5095, Australia
2
Bothar Consulting Limited, Bryansford Village, Newcastle BT33 0AY, UK
*
Author to whom correspondence should be addressed.
Processes 2022, 10(7), 1306; https://doi.org/10.3390/pr10071306
Submission received: 2 June 2022 / Revised: 29 June 2022 / Accepted: 29 June 2022 / Published: 1 July 2022
(This article belongs to the Special Issue Numerical Modeling in Civil and Mining Geotechnical Engineering)

Abstract

:
Pore water pressure changes due to soil-atmospheric boundary interaction can significantly influence soil behaviour and can negatively affect the safety and stability of geotechnical structures. For example, prolonged rainfall events can lead to increased pore water pressure and lower strength; repeated cycles of pore water pressure changes can lead to degradation of strength. These effects are likely to become more severe in the future due to climate change in many parts of the world. To analyse the behaviour of soil subjected to atmospheric boundary interactions, several parameters are needed, and hydraulic conductivity is one of the more important and is difficult to determine. Hydraulic conductivity deduced from laboratory tests are often different from those from the field tests, sometimes by orders of magnitude. The problem becomes even more complicated when the soil state is unsaturated, where the hydraulic conductivity varies with the soil’s state of saturation. In this paper, a relatively simple alternative approach is presented for the estimation of the hydraulic conductivity of unsaturated soils. The method involved a systematic re-analysis of observed pore water pressure response in the field. Using a finite element software, the soil-atmospheric boundary interaction and related saturated/unsaturated seepage of an instrumented slope have been analysed, and results are compared with field measurements. The numerical model could capture the development of suction, positive pore water pressure and changes in water content with reasonable accuracy and demonstrated the usefulness of the hydraulic conductivity estimation method discussed in this paper.

1. Introduction

In many parts of the world, the frequency of extreme weather events has increased due to climate change. January 2019 was the hottest month on the record in Australia [1]; the period of May–July 2007 was the wettest in 250 years and caused extensive flooding in many parts of the UK [2]. Singapore recorded its wettest (since the record began in 1869) December in 2006 [3]. In 2011, Thailand was struck by monsoon and tropical cyclone rains from July to October, causing extreme flooding in the city of Bangkok. More extreme climate scenarios are predicted in Australia, the UK and many other parts of the world [4,5,6,7,8,9].
Climate change is having an impact on the transportation infrastructure in countries like the UK [10]. More than 160 slope failures were reported during the winter of 2000/2001 by the UK rail and road authority [11,12]. In many cases, meteorologically induced pore water pressure and strain softening resulting from pore water pressure cycling have been responsible for the failure of geotechnical structures [13,14,15].
The changes in pore water pressure in soils are often driven by meteorological parameters, vegetation and groundwater hydrology. For more economical and sustainable design and maintenance of geotechnical structures, it is important to understand the mechanics/processes that may affect their behaviour under the current and possible future climate scenarios.
Several past fields and laboratory studies have investigated the effect of atmospheric boundaries on the behaviour of infrastructure slopes [9,16,17,18,19,20]. Among the numerical studies on soil-slope systems, Briggs et al. [21], using a one-dimensional VADOSE/W [22] model, investigated the generation of pore water pressure in a railway embankment; Tsaparas et al. [23] modelled infiltration and compared with field observations of a sloping site in Singapore; Karthikeyan et al. [24] presented estimation of pore water pressure using software called SEEP/W [25] even though there was lack of agreement between the field measurements and model calculations. Rouainia et al. [26] investigated meteorologically induced pore water pressure generation and its effect on slope stability using a coupling of SHETRAN and FLAC-TP flow and showed the capabilities of sophisticated numerical modelling approaches. However, the complexity involved in their analyses was high and needed assumptions on several modelling parameters, which can be difficult to deduce using an objective approach.
One of the more important input parameters required for analysing the interaction between soil and the atmospheric boundary is the hydraulic conductivity (K) of soil. For a particular void ratio, K for saturated soils (Ksat) can be treated as a constant. However, the same does not apply to the case of K for unsaturated soils (Kunsat). It changes with soil suction even if the void ratio remains the same. The soil suction, on the other hand, is dependent on the soil’s degree of saturation. Direct measurement of K will involve measurement of water flow. However, as soon as some flow has occurred, the water content will change, which will lead to a change in suction and eventually the Kunsat. Due to this, it is extremely difficult to measure Kunsat for soils. Several methods have been proposed in the literature to capture the suction-Kunsat correlation [27,28,29] and commonly, Kunsat is expressed as a function of Ksat and soil suction. For example, Van-Genuchten [27] proposed the following equation,
K u n s a t = K s a t 1 a   Ψ n 1 1 + ( a Ψ ) n m 2 1 + ( a Ψ ) n m 2
where a , n and m are curve fitting parameters and Ψ is the matric suction. In essence, with a reasonable estimation of the soil–water characteristics curve (SWCC) and Ksat, Equation (1) can be used to estimate the Kunsat.
A challenge in using Equation (1), often overlooked, is the estimation of Ksat. Several field tests (e.g., single or double ring infiltrometer, disc permeameter, Guelph permeameter, bailout test) and laboratory tests (e.g., Rowe cell, constant/falling head tests, Oedometer tests, triaxial hydraulic conductivity tests) have been developed and routinely conducted to deduce Ksat. However, Ksat deduced from a different field or laboratory tests can be very different, which is partly due to the underlying assumptions involved in those methods. Discrepancies of several orders of magnitude are not uncommon [9]. The field tests are often believed to yield more representative Ksat values compared to the laboratory tests. However, Ksat values deduced from field tests can also be affected by the presence of hidden cracks or fissures in the soil, local heterogeneity and the presence of higher hydraulic conductivity layers, which may not be detected during the field investigation. Thus, the conventional methods of deducing Ksat can have a low degree of reliability. For example, in the study site involved in this investigation (discussed in the next section), more than 100 times difference in Ksat values from different tests were observed. This makes it very difficult, and in some cases, impossible to objectively decide what value of Ksat to be used in a particular analysis [30]. The engineer is often forced to use either an average value or make a subjective judgement which can be influenced by the engineer’s experience with similar sites or the lack of it. For a seepage analysis involving unsaturated soil, the problem becomes even more complex due to the dependency of Kunsat on the soil’s degree of saturation.
This paper discusses a simple and objective approach for estimating Kunsat using a correlation similar to one presented in Eq. [1]. The method involves a systematic back-analysis of field observations to deduce Kunsat as a function of Ksat. The method of estimation of Ksat was originally proposed by Lo et al. [31] to estimate Ksat for prefabricated vertical drain improved soils. Ksat, in their case, was estimated by systematically analysing field settlement data. It allowed them to avoid explicit modelling of the smear zone and reduced Ksat due to smearing in that area which is often very difficult to determine due to limited time and budget associated with field investigation of a practical project. The effectiveness of the method in different problem domains has been demonstrated in the literature by Karim [32], Karim et al. [33], Manivannan et al. [34], Karim et al. [35] and Lo et al. [36]. The approach has been modified in this paper for analyses in the unsaturated domain, i.e., to estimate Kunsat. The effectiveness of the method has been put to the test by modelling a well-instrumented research site located in Southern England [20]. A selected period of monitoring data (from 1 January 2006 to 31 December 2008—a total of 3 years) have been used. A finite element software SEEP/W (version 2021.4) [37] has been used for this purpose.
In the next section, the research site is discussed in terms of its construction, geotechnical properties, instrumentation and monitoring details. Different aspects of the numerical analyses, including the mesh or geometry chosen, boundary conditions and other input parameters, are discussed in the next section. Subsequent sections discuss the results and conclusions drawn from this paper.

2. Newbury Cutting—The Research Site

The research site discussed in this paper, hereafter referred to as Newbury cutting, is located near A34 Newbury bypass in Southern England. The location of the site is presented in Figure 1. The site was extensively instrumented and monitored. Details of the site instrumentation, monitoring, geology and weather data can be found in Smethurst et al. [20] and Smethurst et al. [19]. For the sake of completeness, a brief description is presented here.
The instrumented section of the cut slope was 8 m high and 28 m in length along the sloping direction, as presented in Figure 2. The slope was excavated in 1997. The site soil consists of stiff grey London clay of about 20 m thickness and its properties have been found to vary in the vertical and horizontal directions. The top 2.5 m of the slope, below the original ground level, was extensively weathered and hereafter referred to as weathered London clay. The presence of several bands of silty clay up to 50 mm thick and large flints was also detected during the site investigation. After the cutting was excavated, a 0.4 m layer of topsoil was placed over the cut surface to facilitate the planting of vegetation. To facilitate quick drainage, a fin drain was also installed near the toe of the slope.
A series of field and laboratory tests were conducted by Smethurst et al. [20] to deduce different soil parameters for the research site. The tests included Atterberg limit tests, triaxial (saturated) hydraulic conductivity test, dry unit weight test and field saturated hydraulic conductivity test using the bail-out method. The deduced Ksat, unit weight, and plasticity indices of the site soil are summarised in Table 1. It is to be noted that a large number of laboratory tests on undisturbed samples along with bailout tests (in the field from hand augured boreholes of up to 3 m deep) were used to deduce the Ksat values. 1 to 3 orders of magnitude differences between deduced Ksat values were reported by Smethurst et al. [20]. It was believed to be due to anisotropy and soil fabric, including silt partings and cracks and fissures.
The relationship between the soil water content and suction (also known as soil water characteristics curve—SWCC) for London clay was reported by Croney [39] and is presented in Figure 3 and can be assumed to be representative of the site soil [9]. It is to be noted that Croney [39] presented SWCC for both drying and wetting stages for London clay. As SEEP/W [37] does not allow the use of multiple SWCCs and also considering the difficulties associated with generating consistent wetting SWCC, the drying curve was used for this analysis. It is expected that if the wetting SWCC was used, the observed pore water pressure response would be different as Kunsat is a function of Ksat and SWCC. Change in one may require adjustment in the other parameter. However, its significance was outside the scope of this study.
Moreover, it should be noted that SEEP/W [37] and most other numerical software use SWCC in terms of volumetric water content and the curves presented by Croney [39] were in terms of gravimetric water content. The gravimetric values were converted into volumetric values using the following relationship (Equation (2)) and the transformed SWCCs are presented in Figure 3.
w v o l = w g r × γ d r y γ w
where wvol is the volumetric water content, wgr is the gravimetric water content, γdry and γw are the dry unit weight of soil and unit weight of water, respectively. The average values for γdry for the two different soil layers were used for the calculations, as presented in Table 1.
The vegetation on the slope was predominantly rough grass and herbs with some small shrubs (0.5 m < tall). The grass and herbs were mowed periodically until October 2002 to help the development of shrubs planted on the slope. Several grown-up Beech, Oak and Silver birch trees were located near the top of the slope.
Instrumentation was done to monitor moisture content (TDR moisture probes), suction (standard water-filled tensiometers and equitensiometer), positive water pressure (piezometers), and water table location (dipping boreholes) at depths between 0.3 and 3.5 m. The instruments were installed in four groups, namely, cluster A (near the crest) to D (near the toe). A weather station was also installed at the site. A 350 mm deep interceptor drain was installed on the slope face (across the slope face) to capture surface runoff and interflow. The site was monitored from October 2002.
Figure 4 and Figure 5 present average daily rainfall and runoff recorded (respectively) for the 3 years duration under consideration in this study. Other recorded meteorological parameters for the duration, i.e., average daily temperature, average daily wind speed (2 m above ground surface) and average daily humidity, are presented in Appendix A. These data will later be used in calculations.

3. Numerical Analyses

SEEP/W [37] is a finite element software that allows analysis of groundwater seepage and pore water pressure distribution within porous media such as soil in both saturated and unsaturated states. The types of analyses can range from saturated steady-state problems to complex saturated/unsaturated time-dependent (transient) problems. The pore water pressure distribution from such analyses could be used as an input for other stress or stability analyses to assess the serviceability and ultimate limit state behaviour of the slopes.
The software requires SWCC and soil hydraulic conductivity function (i.e., variation of K with soil suction) as inputs for material definition. The effect of climate variables on the problem can be modelled using an equivalent hydrological boundary condition (explained later). Different modelling aspects of the problem are discussed in the next few sections.

3.1. Model Geometry

The discretised geometry used for the analysis is presented in Figure 6. To minimise the effect of imposed boundary conditions on the analysis results, the geometry was extended by 20 m towards the right and by 15 m towards the left. Several trials were run to ensure the effect of boundary and mesh size was minimum or insignificant to the observed results. Details of the trials are presented in Appendix B. A global mesh size of 0.8 m was used in the analysis. Except for the surface layer, the geometry was discretised using an automatically generated unstructured mesh using triangular and quadrilateral elements. The elements close to the boundary (to a depth greater than the maximum depth of measurement during field investigation) were then refined to half the global element size. The surface layer was discretised using quadrilateral elements. Quadrilateral elements have been shown to perform better in such scenarios as the gradients of primary unknowns are steeper in the direction normal to the surface [37]. The nodal convergence was checked. A total of 4200 nodes and 3965 elements (4-nodded quadrilateral and 3-nodded triangular) were used to discretise the geometry.
The soil in the slope had two distinct subdivisions, namely, weathered London clay and grey London clay layers. They were modelled using different materials. Because of the presence of plant roots and other organic matters, the K of the top surface of the soil slope was higher than the layers below. Separate surface layers (of 0.4 m thickness with element thickness of 0.1 m) using different material models were added to the geometry to account for these differences (i.e., surface layers corresponding to weathered London clay and grey London clay layers). The pore water pressure and other variables usually change rapidly close to the exposed soil surface, i.e., in the surface layer.

3.2. Material Models

SEEP/W [37] allows the soil to be modelled as saturated only or as a combination of saturated and unsaturated soil states. All four soil materials (i.e., grey and weathered London clay and their respective surface layers) were modelled using the saturated/unsaturated option. This way, the soil could stay in a saturated or unsaturated state during an analysis or change state depending on the analyses conditions. The SWCCs for the weathered and grey London clay have been discussed above. The SWCC for the surface layers was kept similar to the layers below. There were no specific K data available for the surface layers. The Ksat values for the soils were assumed as 10 times the values of the layers below [40].

3.3. Calculation for Input Hydrological Boundary Condition

The water balance equation [41] was used for the calculation of input hydrological boundary conditions to represent the effect of the atmospheric boundary (e.g., rainfall, solar radiation, humidity, wind speed) and vegetation as below,
R R O E T S + R E = 0
where R is the rainfall, RO is the runoff, ET is the evapotranspiration, S is the change in stored water within the soil, and RE is the net recharge from surrounding soils. If we consider inflow and outflow for soil mass from the surrounding soil is equal, i.e., RE = 0, we can say S is equal to the magnitude of water percolating into the soil, net surface infiltration (NSF) through the surface boundary. Thus we can rewrite Equation (3) as below,
R R O E T = N S F
The measured daily rainfall and runoff and calculated ET data were used for this purpose. To estimate ET, a reference ET was calculated first using the Penman-Monteith equation as below [42],
E T r = 0.408 Δ R n G + γ 900 T + 273 u 2 e s e a Δ + γ 1 + 0.34 u 2
where ETr is the reference evapotranspiration calculated in a 1-day time step, Δ is the slope of saturation vapour curve, Rn is the net radiation flux in MJ/m2/day, G is the soil heat flux density in MJ/m2/day, γ is the Psychrometric constant in kPa/°C, T is the mean temperature in °C 2 m above the ground level, u2 is the wind speed in m/s at 2 m above the ground, and (esea) is the saturation vapour pressure deficit in kPa.
The Penman–Monteith equation [42] parameters presented in Table 2, along with the weather data presented in Appendix A have been used. Once the daily ETr was calculated, ET was deduced using the following equation,
E T = E T r   for   0     SMD     RAW E T = E T r × T A W S M D T A W R A W   for   SMD     RAW
where SMD is the soil moisture deficit, RAW is readily available water, and TAW is the total available water. Calculated ET values are presented in Figure 7 and calculated daily NSF values are presented in Figure 8. NSF values represent the infiltration of water into the soil due to rainfall and also the removal of water due to evapotranspiration. NSF can be of either positive or negative magnitude. A positive value will indicate water is infiltrating into the ground and the rainfall is dominating the process. A negative value, on the other hand, will indicate water is leaving the soil due to evapotranspiration. The calculated NSF values were used as an input hydrological boundary condition for this analysis and represented the effect of climate and vegetation on the problem. It is to be noted that while using eq [3] to estimate NSF, the value of runoff needs to be known. However, it may not be a readily available measurement. In the absence of measured runoff values, a process outlined in Appendix C can be used to estimate runoff and thus NSF. Similar approaches have been used in the literature [9,43,44].
It is to be noted that the use of a 1-day time step means everything occurring within a particular 24 h period is averaged out over that period and the effect of a shorter duration event may not be captured very accurately. It is possible to conduct a similar analysis with smaller time steps (e.g., hourly) if hourly measurements of all related parameters were available.

3.4. Other Boundary Conditions

The left, right and the bottom boundary of the problem was modelled as impermeable (no flow). Part of the ground surface with road surfacing was modelled as impermeable as well. The mesh and geometry sensitivity analyses show that the choice of these boundary conditions was appropriate. The fin drain installed near the toe of the slope was not modelled explicitly. Rather its effect was modelled using a zero pore water pressure boundary. A flux boundary condition was applied at the ground surface, as discussed earlier.

3.5. Initial Condition

The initial pore water pressure condition in the model was established using an initial water table. An educated guess was made by analysing the piezometer pore water pressure data on the date 1 January 2006. This might have some effect on the predicted pore water pressure profile, and to the author’s understanding, a more complicated seepage analysis could be conducted to establish the initial pore water pressure profile. However, as will be shown later, the effect of assuming an initial water table gets diminished as the analysis progresses and somewhat becomes irrelevant beyond a few months of analysis, even when a very crude guess is made.

3.6. Estimation of K

For modelling unsaturated soils, SEEP/W [37] allows the choice of different functions that deduces Kunsat by relating soil suction to Ksat. A user can choose from available functions such as Fredlund et al. [29], Green and Corey [28] and Van-Genuchten [27]. In this investigation, the Van-Genuchten [27] function (Equation (1)) was used. Estimating an appropriate value for Ksat was a challenge as large (more than 100 times) differences were observed for values from different tests. As indicated earlier, in such scenarios, one has a choice of using an average value or using subjective judgment.
An alternative approach was used to estimate a representative Ksat value in this study. The first year of field-measured pore water pressure data was back analysed. Few trial analyses were conducted with Ksat being systematically varied until the best match was found. The field average Ksat as presented in Table 1, was used as the first trial value. SWCC and other parameters in the study remained unchanged (e.g., the ratio between Ksat in surface layers and the underlying layers). The best match Ksat values were then used to analyse the problem further. Piezometer readings from two different depths at two locations, A and C (see Figure 2 and Figure 6), were used.
In principle, this process is not very different from conventional techniques (especially field tests) used to deduce Ksat. Interpretation of many of the field tests also involves several assumptions and back-fitting of field observations. For example, in a Guelph permeameter test, the flow of water into the soil is observed, and collected data are interpreted using saturated-unsaturated flow theories to estimate Ksat.
Figure 9 presents the effect of using different K values on the pore water pressure profile for location C at 2.5 m depth. Three different analyses are presented here, i.e., analysis using the field average Ksat, analysis using the best match Ksat and analysis using three times the best match Ksat. The figures show that Ksat may play a significant role in such pore water pressure analyses. The best match Ksat values for the weathered London clay and grey London clay were 0.015 m/d and 0.0025 m/d, respectively—approximately four and eight times their respective field average values. Plots at other locations are not presented here due to space limitations.

4. Results

The analyses were run from 1 Jan 2006 to 31 December 2008, a total duration of 3 years, including the first year used for the estimation of Ksat. The results are discussed below.
Figure 10a,b show the field measured and calculated pore water pressure from two analyses (i.e., one using matched K and the other using average field value) at 1.5 and 2.5 m depths from the surface at instrumented section A and Figure 11a,b presents the same for location C. At location A, both analyses captured the negative pore-water pressure peaks with reasonable accuracy. Both analyses (using matched and field K) overestimated the magnitude of positive pore water pressure throughout the measurement period. At location C, the qualitative trend (peaks and troughs) was much better captured in the matched K analysis and the use of average field K in the analyses led to a significant underestimation of negative pore water pressure developed during the period of June-July of 2006. The number of peaks and troughs in the field measured values were better captured by using matched K in the analyses. The average K analysis trend was relatively smoother and often did not capture the smaller peaks and troughs in the trend (short duration events). Overall, the use of matched K improves the accuracy of the calculation.
The near-surface soil suctions were recorded by standard water-filled tensiometers at different depths (0.3, 0.6 and 0.9 m from the surface at locations A and C). Figure 12a,b present the measured suction and estimated responses from the numerical analyses (both using field average K and matched K) at 0.6 m depth at instrumented locations A and C, respectively. At location A the suction magnitude was underestimated between June and July 2007 and overestimated at other times by both analyses. The difference between the calculations from the two analyses was very small. At location C, the suction response was significantly overestimated by both analyses. The overestimation was slightly higher for the case of the field K analysis. The qualitative trend was reasonably well captured by both analyses with reasonable accuracy.
Figure 13 presents the measured and calculated volumetric water content variation with time. At location C, both analyses calculated volumetric water content with reasonable accuracy; however, at location A, both analyses significantly underestimated the changes in water content with time. This difference is likely to be due to local heterogeneity of the soil and the SWCC used in the analyses is unable to represent its behaviour. It is to be noted that in a saturated/unsaturated seepage analysis, the hydraulic conductivity and SWCC both affect the changes in water content (as Kunsat is a function of soil suction) and consequently developed pore water pressure. Estimating the Kunsat based on observation of pore water pressure response has allowed prediction of the same with greater accuracy. It is to be noted that the estimation of Kunsat could also be achieved by comparing the predicted and field observed water content data, and in such a case, it is expected that the prediction of water content could have improved accuracy. It is envisaged that with a better estimation of SWCC for different materials used in the simulation, a better estimation of water content could be achieved.
The interaction between soil, vegetation and the atmospheric boundary is complex and can be influenced by a number of variables. This paper presents a simplified process to model the interaction to capture the related changes in water pressure and soil water content and outlines an objective approach to estimate one of the most important parameters in seepage analyses (i.e., hydraulic conductivity). The simplification and assumptions include,
  • the use of an SWCC from literature for a similar soil (not based on soil tested on this particular site)
  • the use of the drying curve for representing both wetting and drying behaviour
  • the assumption of uniform soil properties within different soil layers even though the site investigation revealed the grey London clay to be variable along with the vertical and horizontal directions.
Despite these simplifications, the analyses presented here could capture the qualitative trend and quantitative changes in pore water pressure and water content with reasonable accuracy. The difference between the field averaged value of Ksat and matched Ksat was four to eight times. It is to be noted that standard laboratory or field tests for hydraulic conductivity often produce estimates that can be 100 or even 1000 times different from each other. Furthermore, the calculation here was found to be significantly less sensitive to Ksat compared to seepage or consolidation analyses in saturated soils. Due to this, the differences between the field average K analysis and matched K analysis were close to each other on many occasions, even though the matched K analysis showed a relatively better overall prediction for water pressure and water content. In other situations where more variability exists, the situation may be significantly different, and the use of field average K may produce far inferior predictions.

5. Discussion and Conclusions

A set of numerical analyses has been discussed in this paper. The analysis has been carried out with a finite element software SEEP/W [37]. All the inputs to the model have been determined from laboratory/field test results except for the hydraulic conductivity of the soil, which is likely to be affected by the presence of fissures or cracks in the soil and may not be captured by standard field or laboratory tests. An alternative approach is outlined here for the estimation of the hydraulic conductivity. Surface infiltration was used in the model as an input boundary condition and was calculated using the surface water balance equation. The results show that using the approach, it is possible to capture the seasonal, climate-induced pore water pressure variation in slopes, and the discussed method for estimation of hydraulic conductivity can be a useful tool for seepage analyses in unsaturated soils.

Author Contributions

Conceptualisation, M.R.K. and D.H.; methodology, M.R.K. and M.M.R.; software, M.R.K.; validation, M.R.K.; formal analysis, M.R.K.; investigation, D.H. and M.R.K.; resources, D.H.; data curation, D.H. and M.R.K.; writing—original draft preparation, M.R.K.; writing—review and editing, D.H. and M.M.R.; visualisation, M.R.K. and M.M.R.; supervision, D.H.; project administration, D.H.; funding acquisition, D.H. All authors have read and agreed to the published version of the manuscript.

Funding

Part of this research was carried out when the first author was working as a Research Fellow at the Queens University Belfast and was supported by an EPSRC-funded project, Infrastructure slopes Sustainable Management and Resilience Assessment—iSMART (http://www.ismartproject.org, accessed on 1 June 2022).

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Data related to publication can be accessed via contacting any of the authors in the paper.

Conflicts of Interest

The authors declare no conflict of interest

Appendix A

Figure A1, Figure A2 and Figure A3 showing recorded daily temperature, wind speed and relative humidity at Newbury cutting
Figure A1. Average daily temperature at Newbury cutting.
Figure A1. Average daily temperature at Newbury cutting.
Processes 10 01306 g0a1
Figure A2. Average daily wind speed 2 m above the ground level at Newbury cutting.
Figure A2. Average daily wind speed 2 m above the ground level at Newbury cutting.
Processes 10 01306 g0a2
Figure A3. Average daily humidity at Newbury cutting.
Figure A3. Average daily humidity at Newbury cutting.
Processes 10 01306 g0a3

Appendix B

To investigate the sensitivity of mesh size and the effect of imposed boundary conditions, additional analyses were conducted. Results from analysis based on mesh and geometry presented in Figure 6 were used as a benchmark and deviation from those results due to changes in the mesh and geometry condition was observed. For comparison, the analyses are named GM1—original analyses reported in this paper; GM2—analyses with the same geometry as GM1 but without mesh refinement (element size of 0.8 m except for the surface layers where the element thickness was 0.1 m) and GM3- with reduced depth and length of the geometry and 0.4 m global element size with a surface layer having 0.1 m element thicknesses). The mesh and geometries are shown in Figure A4 and Figure A5, and calculated pore water pressure/suction are compared for location C at three different depths in Figure A6.
As can be seen in Figure A6, there were only very small differences between calculated responses from all the additional analyses. So the mesh size and the boundary could be treated as the reason for the problem being solved and the mesh and geometry combination GM1 was used for all the analyses to remove any concern of numerical ill-conditioning or boundary effect on results.
Figure A4. Geometry and mesh used for sensitivity analysis GM2.
Figure A4. Geometry and mesh used for sensitivity analysis GM2.
Processes 10 01306 g0a4
Figure A5. Geometry and mesh used for sensitivity analysis GM3.
Figure A5. Geometry and mesh used for sensitivity analysis GM3.
Processes 10 01306 g0a5
Figure A6. suction/pore water pressure response at instrumentation line C at (a) 0.3 m depth, (b) 1 m depth and (c) 2 m depth.
Figure A6. suction/pore water pressure response at instrumentation line C at (a) 0.3 m depth, (b) 1 m depth and (c) 2 m depth.
Processes 10 01306 g0a6aProcesses 10 01306 g0a6b

Appendix C

In many cases, the measured runoff data is not available. In such a case, a two-tank soil water storage model (SWSM) [9,43,45] can be used for the estimation of runoff. The model calculates how much water from a particular intensity of rainfall infiltrates the soil based on the available storage capacity and hydraulic conductivity of the soil. The amount of water that cannot infiltrate is treated as RO. That means, even if the soil has enough storage capacity, if the rainfall intensity is more than the soil’s capacity to allow water in (hydraulic conductivity), the excess water will flow as runoff. In this case, both pore water pressure and runoff are a function of K and an iterative process as outlined below will be necessary to estimate K and thus NSF.
  • The first step in the process is to make an initial estimation of Ksat. It is difficult to objectively estimate a value for Ksat as the numbers can be scattered over a large range. A field average can be a good starting point.
  • Following this, RO can be estimated using the two-tank SWSM model and the first estimate of daily NSF can be calculated.
  • The deduced NSF can then be used as an input boundary condition in a first trial numerical analysis. Pore water pressure distribution from the simulation can then be compared with the field measured values.
  • In the next step, Ksat values can be systematically varied and step 1 to 3 repeated until a good match between the observed and calculated pore water pressures can be found. For every chosen Ksat, the calculation for surface RO and NSF needs repeating.

References

  1. BOM. Monthly Summary for Australia; IDCKGC1A00; Bureau of Meteorology: Melbourne, Australia, 2019. [Google Scholar]
  2. Toll, D.G.; Mendes, J.; Augarde, C.E.; Karthikeyan, M.; Phoon, K.K.; Gallipoli, D.; Lin, K.Q. Effects of climate change on slopes for transportation infrastructure. In Proceedings of the International Conference on Transportation Geotechnics, Nottingham, UK, 24–28 August 2008. [Google Scholar]
  3. NEA. Singapore Experienced Its Wettest December since 1968; National Environment Agency: Singapore, 2006. [Google Scholar]
  4. Murphy, J.M.; Sexton, D.M.H.; Jenkins, G.J.; Boorman, P.M.; Booth, B.B.B.; Brown, C.C.; Clark, R.T.; Collins, M.; Harris, G.R.; Kendon, E.J.; et al. UK Climate Projections Science Report: Climate Change Projections; Met Office Hadley Centre: Exeter, UK, 2009. [Google Scholar]
  5. CSIRO. Climate Change in Australia, Projections for Selected Australian Cities; CSIRO: Canberra, Australia, 2015. [Google Scholar]
  6. ClimateChangeAustralia. Climate Change in Australia. Available online: https://www.climatechangeinaustralia.gov.au/en/ (accessed on 2 November 2020).
  7. Karim, M.R.; Rahman, M.M.; Nguyen, K.; Cameron, D.; Iqbal, A.; Ahenkorah, I. Changes in Thornthwaite Moisture Index and reactive soil movements under current and future climate scenarios-a case study. Energies 2021, 14, 6760. [Google Scholar] [CrossRef]
  8. Karim, M.R.; Rahman, M.M.; Nguyen, H.B.K.; Newsome, P.L.; Cameron, D. TMI soil moisture index for South Australia under current and future climate scenarios. In Proceedings of the 20th International Conference on Soil Mechanics and Geotechnical Engineering, a Geotechnical Discovery Downunder, Sydney, Australia, 1–5 May 2022; Rahman, M.M., Jaksa, M., Eds.; Australian Geomechanics Society: Sydney, Australia, 2022; pp. 2233–2236. [Google Scholar]
  9. Karim, M.R.; Hughes, D.; Kelly, R.; Lynch, K. A rational approach for modelling the meteorologically induced pore water pressure in infrastructure slopes AU—Karim, Md Rajibul. Eur. J. Environ. Civ. Eng. 2019, 24, 2361–2382. [Google Scholar] [CrossRef]
  10. Loveridge, F.A.; Spink, T.W.; O’Brien, A.S.; Briggs, K.M.; Butcher, D.J.E. The impact of climate and climate change on infrastructure slopes, with particular reference to southern England. Q. J. Eng. Geol. Hydrogeol. 2010, 43, 461–472. [Google Scholar] [CrossRef]
  11. Ridley, A.; McGinnity, B.; Vaughan, P.R. Role of pore water pressures in embankment stability. Proc. Inst. Civ. Eng. 2004, 157, 193–198. [Google Scholar] [CrossRef]
  12. Turner, S. Climate Change Blamed as Landslip Incidents Treble. Available online: https://www.newcivilengineer.com/climate-change-blamed-as-landslip-incidents-treble/812153.article (accessed on 1 July 2015).
  13. Sivakumar, V.; Hughes, D.; Clarke, G.; Glynn, D.T. A case study: Delayed failure of a deep cutting in lodgement till. Proc. ICE Geotech. Eng. 2007, 160, 193–202. [Google Scholar] [CrossRef]
  14. Ng, P.B.; Poh, K.K.; Tenando, E. Slope failure after prolonged rainfall. In BCA Seminar- Approach to Structural Inspection and Slope Stability; Building and Construction Authority: Singapore, 2007. [Google Scholar]
  15. Hughes, P.N.; Glendinning, S.; Mendes, J.; Parkin, G.; Toll, D.G.; Gallipoli, D.; Miller, P.E. Full-scale testing to assess climate effects on embankments. Proc. Inst. Civ. Eng. Eng. Sustain. 2009, 162, 67–69. [Google Scholar] [CrossRef]
  16. Bolton, M.D.; Take, W.A. Seasonal ratcheting and softening in clay slopes, leading to first-time failure. Géotechnique 2011, 61, 757–769. [Google Scholar]
  17. Briggs, K.M.; Smethurst, J.A.; Powrie, W.; O’Brien, A.S.; Butcher, D.J.E. Managing the extent of tree removal from railway earthwork slopes. Ecol. Eng. 2013, 61, 690–696. [Google Scholar] [CrossRef] [Green Version]
  18. Clarke, D.; Smethurst, J.A. Effects of climate change on cycles of wetting and drying in engineered clay slopes in England. Q. J. Eng. Geol. Hydrogeol. 2010, 43, 473–486. [Google Scholar] [CrossRef]
  19. Smethurst, J.A.; Clarke, D.; Powrie, W. Factors controlling the seasonal variation in soil water content and pore water pressures within a lightly vegetated clay slope. Géotechnique 2012, 62, 429–446. [Google Scholar] [CrossRef]
  20. Smethurst, J.A.; Powrie, W.; Clarke, D. Seasonal changes in pore water pressure in a grass-covered cut slope in London Clay. Géotechnique 2006, 56, 523–537. [Google Scholar] [CrossRef] [Green Version]
  21. Briggs, K.M.; Smethurst, J.A.; Powrie, W.; O’brien, A.S. Wet winter pore pressures in railway embankments. Proc. Inst. Civ. Eng. 2012, 166, 451–465. [Google Scholar] [CrossRef] [Green Version]
  22. Geo-Slope. Vadose Zone Modeling with VADOSE/W An Engineering Methodology; Geo-Slope International Limited: Calgary, AB, Canada, 2013. [Google Scholar]
  23. Tsaparas, I.; Toll, D.G. Numerical analysis of infiltration into unsaturated residual soil slopes. In Proceedings of the 3rd International Conference on Unsaturated Soils, Recife, Brazil, 10–13 March 2002. [Google Scholar]
  24. Karthikeyan, M.; Toll, D.G.; Phoon, K.K. Prediction of changes in pore-water pressure response due to rainfall events. In Unsaturated soils: Advances in Geo-Engineering, Toll, D.G., Ed.; Taylor & Francis Group: London, UK, 2008; pp. 829–834. [Google Scholar]
  25. Geo-Slope. Seepage Modeling with SEEP/W an Engineering Methodology; Geo-Slope International Limited: Calgary, AB, Canada, 2013. [Google Scholar]
  26. Rouainia, M.; Davies, O.; O’Brien, T. Numerical modelling of climate effects on slope stability. Proc. Inst. Civ. Eng. Eng. Sustain. 2009, 162, 81–89. [Google Scholar] [CrossRef]
  27. Van-Genuchten, M.T. A Closed-form Equation for Predicting the Hydraulic Conductivity of Unsaturated Soils. Soil Sci. Soc. Am. J. 1980, 44, 892–898. [Google Scholar] [CrossRef] [Green Version]
  28. Green, R.E.; Corey, J.C. Calculation of hydraulic concuctivity: A further evaluation of some predictive methods. Soil Sci. Soc. Am. J. 1971, 35, 3–8. [Google Scholar] [CrossRef]
  29. Fredlund, D.G.; Xing, A. Equations for the soil-water characteristic curve. Can. Geotech. J. 1994, 31, 521–532. [Google Scholar] [CrossRef]
  30. Karim, M.R.; Lo, S.-C.R. Estimation of the hydraulic conductivity of soils improved with vertical drains. Comput. Geotech. 2015, 63, 299–305. [Google Scholar] [CrossRef]
  31. Lo, S.R.; Mak, J.; Gnanendran, C.T.; Zhang, R.; Manivannan, G. Long-term performance of a wide embankment on soft clay improved with prefabricated vertical drains. Can. Geotech. J. 2008, 45, 1073–1091. [Google Scholar] [CrossRef]
  32. Karim, M.R. Behaviour of piles subjected to passive subsoil movement due to embankment construction—A simplified 3D analysis. Comput. Geotech. 2013, 53, 1–8. [Google Scholar] [CrossRef]
  33. Karim, M.R.; Manivannan, G.; Gnanendran, C.T.; Lo, S.-C.R. Predicting the long-term performance of a geogrid-reinforced embankment on soft soil using two-dimensional finite element analysis. Can. Geotech. J. 2011, 48, 741–753. [Google Scholar] [CrossRef]
  34. Manivannan, G.; Karim, M.R.; Gnanendran, C.T.; Lo, S.-C.R. Calculated and observed long term performance of Leneghans embankment. Geomech. Geoengin. 2011, 6, 195–207. [Google Scholar] [CrossRef]
  35. Karim, M.R.; Gnanendran, C.T.; Lo, S.-C.R.; Mak, J. Predicting the long-term performance of a wide embankment on soft soil using an elastic-visco-plastic model. Can. Geotech. J. 2010, 47, 244–257. [Google Scholar] [CrossRef]
  36. Lo, S.R.; Karim, M.R.; Gnanendran, C.T. Consolidation and Creep Settlement of Embankment on Soft Clay: Prediction Versus Observation. In Geotechnical Predictions and Practice in Dealing with Geohazards; Chu, J., Wardani, S.P.R., Iizuka, A., Eds.; Springer: Dordrecht, The Netherlands, 2013; pp. 77–94. [Google Scholar]
  37. Geo-Slope. Heat and Mass Transfer Modelling with GeoStudio; Seequent Limited (GEOSLOPE International, Ltd.): Calgary, AB, Canada, 2021. [Google Scholar]
  38. GoogleMap. Available online: https://www.google.com/maps/place/Newbury,+UK/@51.3885288,-1.3609492,13z/data=!3m1!4b1!4m5!3m4!1s0x487402002f595ba9:0xc6646baff4a75c50!8m2!3d51.401409!4d-1.3231139 (accessed on 15 April 2019).
  39. Croney, D. The Design and Performance of Road Pavements; Her Magistye’s Stationery Office: London, UK, 1977. [Google Scholar]
  40. Toll, D.G.; Rahim, M.S.; Karthikeyan, M.; Tsaparas, I. Soil atmosphere interactions for analysing slopes in tropical soils. In Proceedings of the 14th International conference of the International Association of Computer Methods and Advances in Geomechanics, Koyoto, Japan, 22–25 September 2014; pp. 1333–1338. [Google Scholar]
  41. Blight, G. The vadose zone soil-water balance and transpiration rates of vegetation. Geotechnique 2003, 53, 55–64. [Google Scholar] [CrossRef]
  42. Zotarelli, L.; Dukes, M.D.; Romero, C.C.; Migliaccio, K.W.; Kelly, T. Step by Step Calculation of the Penman-Monteith Evapotranspiration ( FAO-56 Method). 2013. Available online: https://edis.ifas.ufl.edu/pdf/AE/AE45900.pdf (accessed on 1 June 2022).
  43. McLernon, M. Climate Driven Pore Awter Pressure Dynamics and Slope Stability within Glacial Till Drumlins in Northern Ireland; Queen’s Unviersity: Belfast, UK, 2014. [Google Scholar]
  44. Hughes, D.; Karim, M.R.; Briggs, K.; Glendinning, S.; Toll, D.; Dijkstra, T.; Powrie, W.; Dixon, N. A comparison of numerical modelling techniques to predict the effect of climate on infrastructure slopes. In Proceedings of the Geotechnical Engineering for Infrastructure and Development—Proceedings of the XVI European Conference on Soil Mechanics and Geotechnical Engineering, ECSMGE 2015, Edinburgh, UK, 13–17 September 2015; pp. 3663–3668. [Google Scholar]
  45. Karim, M.R.; Hughes, D.; Rahman, M.M. Estimating hydraulic conductivity assisted with numerical analysis for unsaturated soil—A case study. Geotech. Eng. 2021, 52, 12–19. [Google Scholar]
Figure 1. Location map of Newbury cutting with an arrow pointing to the site location [38].
Figure 1. Location map of Newbury cutting with an arrow pointing to the site location [38].
Processes 10 01306 g001
Figure 2. A cross-section of the Newbury cutting (dimensions after Smethurst et al. [20]).
Figure 2. A cross-section of the Newbury cutting (dimensions after Smethurst et al. [20]).
Processes 10 01306 g002
Figure 3. Transformed SWCCs for the Newbury cutting soil.
Figure 3. Transformed SWCCs for the Newbury cutting soil.
Processes 10 01306 g003
Figure 4. Measured daily rainfall at Newbury cutting.
Figure 4. Measured daily rainfall at Newbury cutting.
Processes 10 01306 g004
Figure 5. Measured daily runoff at Newbury cutting.
Figure 5. Measured daily runoff at Newbury cutting.
Processes 10 01306 g005
Figure 6. Discretised geometry was used for the SEEP/W analysis (dimensions are in m).
Figure 6. Discretised geometry was used for the SEEP/W analysis (dimensions are in m).
Processes 10 01306 g006
Figure 7. Calculated ET for Newbury cutting.
Figure 7. Calculated ET for Newbury cutting.
Processes 10 01306 g007
Figure 8. Net surface flux data for Newbury cutting.
Figure 8. Net surface flux data for Newbury cutting.
Processes 10 01306 g008
Figure 9. Different trial Ksat values and related pore water pressure at Section C 2.5 m depth from the surface.
Figure 9. Different trial Ksat values and related pore water pressure at Section C 2.5 m depth from the surface.
Processes 10 01306 g009
Figure 10. Field measured and calculated pore water pressure at (a) 1.5 m and (b) 2.5 m depths in section A.
Figure 10. Field measured and calculated pore water pressure at (a) 1.5 m and (b) 2.5 m depths in section A.
Processes 10 01306 g010
Figure 11. Field measured and calculated pore water pressure at (a) 1.5 m depth and (b) 2.5 m depth at section C.
Figure 11. Field measured and calculated pore water pressure at (a) 1.5 m depth and (b) 2.5 m depth at section C.
Processes 10 01306 g011
Figure 12. Measured and calculated soil suction at 0.3 m depth (a) at section A and (b) at section C.
Figure 12. Measured and calculated soil suction at 0.3 m depth (a) at section A and (b) at section C.
Processes 10 01306 g012
Figure 13. Measured and calculated volumetric water content 0.6 m depth (a) at location A and (b) at location C.
Figure 13. Measured and calculated volumetric water content 0.6 m depth (a) at location A and (b) at location C.
Processes 10 01306 g013
Table 1. Hydraulic conductivity, unit weight and plasticity index of grey and weathered London Clay at the Newbury test site (after Smethurst et al. [20]).
Table 1. Hydraulic conductivity, unit weight and plasticity index of grey and weathered London Clay at the Newbury test site (after Smethurst et al. [20]).
PropertyGrey London ClayWeathered London Clay
RangeAverageRangeAverage
Ksat from triaxial hydraulic conductivity test (m/day)3.37 × 10−6–5.7 × 10−51.99 × 10−54.32 × 10−5–1.38 × 10−47.52 × 10−5
Ksat from borehole bailout tests (m/day)1.99 × 10−4–3.8 × 10−43.2 × 10−43.11 × 10−3–4.32 × 10−3 3.72 × 10−3
Dry unit weight (kN/m3) 13.2–15.214.613.2–16.216.0
Plasticity index (%)32.5–36.434.831.7 *31.7 *
* one out of five of weathered London Clay samples showed plasticity.
Table 2. Parameters used for the calculation of evapotranspiration.
Table 2. Parameters used for the calculation of evapotranspiration.
ParameterValue
Psychrometric constant (kPa/°C)0.000665 × atm. pressure
Solar constant0.082
Latitude (rad)0.895877505
Albido or canopy reflection coefficient0.23
Stefan–Boltzman constant (Mj/K4/M2/day)4.903 × 10−9
Elevation above sea level (m)105
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Karim, M.R.; Hughes, D.; Rahman, M.M. Unsaturated Hydraulic Conductivity Estimation—A Case Study Modelling the Soil-Atmospheric Boundary Interaction. Processes 2022, 10, 1306. https://doi.org/10.3390/pr10071306

AMA Style

Karim MR, Hughes D, Rahman MM. Unsaturated Hydraulic Conductivity Estimation—A Case Study Modelling the Soil-Atmospheric Boundary Interaction. Processes. 2022; 10(7):1306. https://doi.org/10.3390/pr10071306

Chicago/Turabian Style

Karim, Md Rajibul, David Hughes, and Md Mizanur Rahman. 2022. "Unsaturated Hydraulic Conductivity Estimation—A Case Study Modelling the Soil-Atmospheric Boundary Interaction" Processes 10, no. 7: 1306. https://doi.org/10.3390/pr10071306

APA Style

Karim, M. R., Hughes, D., & Rahman, M. M. (2022). Unsaturated Hydraulic Conductivity Estimation—A Case Study Modelling the Soil-Atmospheric Boundary Interaction. Processes, 10(7), 1306. https://doi.org/10.3390/pr10071306

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop