Determination of hydraulic parameters of non-linear consolidation clay layers by type curve method

The consolidation of clay layers is of great significance for groundwater environmental protection, groundwater storage utilization, and land subsidence. In this study, the governing equation for the excess pore water pressure during the non-linear consolidation process of clay layers under load conditions is obtained based on the one-dimensional non-linear consolidation theory. Analytical solutions are then derived for clay layers with single or double drainage caused by the dissipation of the excess pore water pressure. With these analytical solutions, the groundwater dynamics and deformation of the clay layer are analyzed. Correspondingly, a type curve method is proposed to calculate the hydraulic parameters of the clay layer through laboratory experiments, which verifies the reliability of the analytical solutions. The study results show that the deformation of the clay layer predicted by the non-linear consolidation theory is smaller than that predicted by the linear consolidation theory. The deformation of the clay layer increases with the increase in the thickness of the clay layer, the compressive index, and the overburden load, while it decreases with the increase in the initial void ratio and the initial effective stress. The stable time, at which the consolidation of the clay layer is completed, increases with the increase in the compression index and the thickness of the clay layer, while it decreases with the increase in the initial void ratio, the initial effective stress, and the initial hydraulic conductivity. It does not vary with the load pressure. Conclusively, the deformation prediction based on the non-linear consolidation theory is more accurate and applicable to further load pressures.


Introduction
With global population growth and economic development, groundwater pollution and geological disasters induced by human activities are becoming increasingly serious, threatening the sustainable development of human society (Gorelick and Zheng, 2015;Jia et al., 2020;Luo et al., 2020;Shu et al., 2022;Zhang et al., 2023). Clayey soils widely cover the surface of basins, plains, seabed, and groundwater aquifer systems, such as the Yangtze Delta in China and the Virginia coastal plain aquifer system in the U.S. (Konikow and Neuzil, 2007;Guo and Li, 2015). Clay and silt are the main components of aquitards, and the hydraulic characteristics of aquitards are closely related to groundwater environmental protection, groundwater storage utilization, and land subsidence (Li et al., 2017;Li et al., 2019;Guo et al., 2022). In addition, clay soil is a significant component of the liner systems installed at hazardous waste disposal sites (e.g., municipal solid waste landfills and nuclear waste disposal sites), thus it has a positive effect on the protection of the water and soil environment around the site, and is an important component of building foundations, dam seepage control systems, mine tailings dams, etc. (Neuzil, 1986;Xue et al., 2005;Zhou et al., 2013;Nicholls et al., 2021;Yan et al., 2021;Zhu et al., 2021). Therefore, it is important to determine the deformation of the clay layer and its hydraulic parameters for groundwater environmental protection and geological disaster prevention.
Many studies have been conducted on the consolidation and drainage of aquitards (Zhuang et al., 2020;Zhuang et al., 2021). Zhou et al. (2013) derived an analytical solution to describe the change of excess pore water pressure in an aquitard based on Terzaghi's one-dimensional consolidation theory. The drawdown was considered to be a constant in an adjacent aquifer and was verified by laboratory tests. Li et al. (2016) proposed a method for estimating the water release from an aquitard during the entire mining history or a limited period, then applied the method to the aquifer system in the Su-Xi-Chang area, but the groundwater level data were insufficient in the study. In addition, Li et al. (2018) proposed a governing equation of drawdown variation to describe the one-dimensional large-strain consolidation behavior of an aquitard without creep effects and derived an analytical solution for the deformation of a large-strain aquitard under a sudden drawdown in an adjacent aquifer.
In addition, many scholars in soil mechanics and hydrogeology have studied the consolidation and deformation of clay layers (Chen et al., 2021;Liu et al., 2021;Wang et al., 2021). Xie and Leo, (2004) derived a completely explicit analytical solution for one-dimensional large-strain consolidation in thick and thin soil layers and compared it with the classical small-strain theory. They also extensively studied the non-linear consolidation properties of clays. Davis and Raymond, (1965) pointed out that Terzaghi's one-dimensional linear consolidation theory cannot accurately predict the dissipation rate of pore pressure. Gibson et al. (1981) proposed a consistent theory for finite-strain consolidation of thick and uniform clay layers under load and found that the traditional consolidation theories may seriously underestimate excess pore water pressure in soft soil layers. Considering the complex non-linearity of materials and the geometric characteristics of the soil consolidation process, Wu et al. (2010) derived an analytical solution and obtained the variation law of soil effective stress and void ratio under different conditions through finite element numerical analysis. Ma et al. (2021) established a settlement calculation model for soft soil foundations, developed a method to determine the deformation modulus of the soil before and after damage under load, and proposed a non-linear settlement calculation method considering the structural characteristics of soft clay. The effect of temperature and chemical properties of pore water on the seepage and compression of cohesive soils was also analyzed (Abdullah et al., 1999;Romero et al., 2001;Deng et al., 2014).
Moreover, research on soil permeability is also in full swing. Zhou et al. (2020) and Zhao et al. (2020) proposed a novel double threshold segmentation algorithm to segment cracks, pores, and grains, and pore-scale variables are defined and extracted from these X-ray CT images to study the geometric characteristics of microstructures of porous geomaterials. Zhao and Zhou, (2020a), Zhao and Zhou, (2020b),  ;; Zhou and Zhao, (2019) successfully predicted the permeability of soil by establishing a model through X-ray CT imaging test, analyzed the internal microstructures of FC and their effects on the compressive strength and permeability, and studied the hydraulic transport properties and fluid flow in single fractures from different projection directions. Zhao et al. (2021) established a digital quantum mechanism-based neural network (DQNN) to study the permeability using digital porosity, coordination number, and pore network size.
The prediction accuracy of consolidation and deformation processes for an aquitard is obviously affected by the uncertainties and errors in the values of hydraulic properties of the clay layer. Based on laboratory testing data, Gregory et al. (2006) proposed several methods to calculate the compressive index and precompression stress from soil compressive test data. An analytical method was proposed to determine the hydrogeological parameters using the data of deformation for a large-strain aquitard (Li et al., 2018). Given that previous analytical solutions were limited by the assumption of a constant head of an overlying unpressurized aquifer, Zhuang et al. (2015) proposed a type curve method for linear drawdown and deformation cases to estimate the hydrogeological parameters of an aquitard in an aquifer system.
Terzaghi's one-dimensional linear consolidation (Terzaghi, 1943) problem is a classic problem in soil mechanics and is the most generally used theory for the consolidation of saturated soils. However, saturated clay soils have obvious non-linear trends during consolidation (Gibson et al., 1967;Gibson et al., 1981;Xie and Leo, 2004). Under the condition that the drawdown of adjacent aquifers is constant, Luo et al. (2020) derived an analytical solution to calculate the drawdown and water release of an aquitard undergoing non-linear consolidation described by Gibson et al. (1967) and proposed a type curve approach to determine the hydraulic parameters of the aquitard. Moreover, many researchers have studied the non-linear consolidation properties of clay by means of analytical solutions or combined analytical solutions with laboratory experiments under the premise that the reduction of permeability is proportional to the reduction of compressibility (C c C k ) (Davis and Raymond, 1965;Xie et al., 2002;Wu et al., 2010;Luo et al., 2020). To the authors' knowledge, no analytical method has been proposed to determine the hydraulic parameters of clay layers undergoing non-linear consolidation Frontiers in Earth Science frontiersin.org 02 under an overburden load condition. Such a new solution would be useful to analyze the non-linear consolidation phenomenon of clay in practical engineering problems, such as soft soil building foundations and landfill clay liners.
Based on the one-dimensional non-linear consolidation theory developed by (Gibson et al., 1967;Gibson et al., 1981), this study first derived the analytical solutions of the excess pore water pressure and the deformation in accordance with Darcy's law under constant overburden load. Next, the deformation of clay layers with single or double drainage, i.e., singly (SDCL) or doubly draining clay layer (DDCL), was analyzed, and the factors for the deformation were explored by analyzing the sensitivity of multiple parameters. In addition, through a laboratory single and double drainage consolidation test, a type curve approach based on the analytical solutions of non-linear consolidation theory was developed to determine the hydraulic parameters of clay layers and to verify the theoretical research presented in this study.

Mathematical model and analytical solutions 2.1 Governing equation for the consolidation of clay layers
The schematic diagram of the one-dimensional non-linear consolidation of a clay layer under an overlying load condition is shown in Figure 1, in which the thickness of the clay layer is l. Using the Lagrangian coordinate system, a (L) is positive vertically downward, and the origin is located at the top of the clay layer. Gibson et al. (1967); Gibson et al. (1981) gave a governing equation for the non-linear consolidation of large-deformed soil layers (one-dimensional non-linear consolidation theory) governed by the void ratio e, while ignoring the self-weight effect of the soil. z za where e is the void ratio (dimensionless), e 0 is the initial void ratio (dimensionless), k(e) is the hydraulic conductivity (LT -1 ), which is the function of e; σ′ is the effective stress (ML -1 T -2 ), t is the time (T), and γ w is the unit weight of water (ML -2 T -2 ). The one-dimensional non-linear consolidation theory is assumed here : 1) the water flow in the clay layer obeys Darcy's law while one-dimensional consolidation of the saturated clay layer occurs; 2) the compression of the clay layer equals the cumulative water release from the clay layer per unit area, because the compressibility of soil particles and pore water is significantly less than that of soil skeleton, and the soil particles and pore water are assumed to be incompressible; 3) σ′ and k(e) of the clay layer corresponding to the e are in accordance with the compression and permeability models (Eqs 2, 3), respectively (Mitchell, 1993): where C c is the compression index (dimensionless), σ 0 ′ is the initial effective stress (ML -1 T -2 ), C k is the penetration index (dimensionless), and k(e 0 ) is the initial hydraulic conductivity; 4) the clay layer undergoes normal consolidation, and the non-linear variations of the soil compressibility during the consolidation period are assumed to obey Eq. 4 (Wu et al., 2010), and the creep effect is not considered here.
1 + e 0 − C c log 10 σ′ σ 0 where e f is the final void ratio, which can be calculated by Eq. 2 using the ultimate effective stress σ ′ f at the end of clay layer consolidation, and σ ′ f is equal to the sum of σ 0 ′ and additional effective stress caused by the overlying load above the clay layer. Furthermore, the principle of effective stress is observed in the clay layer, where p is the load due to some overlying structure or layers (ML -1 T -2 ), and u(a, t) is the excess pore water pressure at position a in the clay layer at time t (ML -1 T -2 ). It is assumed that the volumetric compressibility and hydraulic conductivity of the clay layer decrease proportionally (i.e., C k C c ) during consolidation and the corresponding consolidation process (Davis and Raymond, 1965;Xie et al., 2002;Wu et al., 2010;Luo et al., 2020). The initial effective stress (σ 0 ′ , geostatic stress) is distributed uniformly through the entire thickness and has little effect on a thin clay layer, and the calculation error of the water released from a thick clay layer will be greater (Gibson et al., 1967;Gibson et al., 1981). Substituting Eqs 2-5 into Eq. 1, the governing equation for the dissipation of excess pore water pressure through a clay layer undergoing one-dimensional non-linear consolidation can be obtained, namely, Eq. 6 uses excess pore water pressure u as the primary variable, and which is closely related to the coefficient of consolidation (c v ) in Terzaghi's one-dimensional linear consolidation theory (L 2 T -1 ). In addition, c w relies on the values of k(e 0 ), σ 0 ′ , C c , e 0 and e f . Frontiers in Earth Science frontiersin.org 03

Analytical solutions
Under the conditions of SDCL and DDCL, u in the clay layer can be obtained through the governing equation (Eq. 6) with the following initial and boundary conditions (Eqs 8-10). Here, the overlying load is assumed to be constant.
The variations of excess pore water pressure u(a, t) in SDCL and DDCL can be obtained using the variable transformation and characteristic function method, respectively, and expressed as follows: where σ ′ f σ 0 ′ + p, the subscripts 1 and 2 of the symbols correspond to the case of single-sided drainage and double-sided drainage, respectively, and the detailed derivations are provided in Appendix A, B, respectively. The void ratio e(a, t) variation for SDCL and DDCL can be obtained by substituting Eqs 9, 10, respectively, with the assumption proposed by Mitchell (1993).
Similarly, the hydraulic conductivity k(a, t) variation for SDCL and DDCL can be derived, respectively.
Based on Darcy's law of seepage, the dimensionless flow velocity q( a, t) for SDCL and DDCL can be obtained, respectively.
where q(a, t) q( a, t) ln 10σ0 ′ k(e0) is the flow velocity at position a and time t in SDCL and DDCL (L 3 T −1 ), a a/l is the dimensionless form of a, and t c w t/l 2 is the dimensionless time.
According to the water balance principle, the dimensionless rate of water release for the SDCL and DDCL per unit horizontal area q( t) q(1, t) − q(0, t), respectively, where q(0, t) and q(1, t) represent the flow velocity at the upper and lower clay layer interfaces, respectively.
By integrating Eqs 17, 18 over time, the cumulative water fluxes per unit horizontal area of SDCL and DDCL can be obtained, respectively.
Owing to the above assumption about the incompressibility of soil particles, the accumulated released water per unit area of SDCL and DDCL which is equal to the deformation of the clay layer.
σ0 ′ is the deformation of SDCL and DDCL (L). The analytical solutions proposed above are compared with the solution based on onedimensional linear consolidation theory, where the former considers the non-linear variation of compressibility and permeability during the consolidation of the clay layer (Mitchell, 1993). Therefore, the proposed solution is more practical and accurate in estimating deformation and consolidation problems than the conventional linear theory, especially for newly formed soft sedimentary soils. Figure 2A shows that the evolution of q with t. q(0, t) is always equal to 0, which leads to q(1, t) q( t) for SDCL; while the flow values at the upper and lower clay layer interfaces for DDCL are equal, but in the opposite direction. It is obvious that both SDCL and DDCL have the same variation law: 1) In the early stage of the consolidation process, due to the abrupt overlying load at the initial time, u in the clay layer instantaneously reaches the maximum value, and q( t) decreases significantly over time. The value of q( t) of DDCL is greater than that of SDCL.

Dimensionless analysis
2) The value of q( t) gradually decreases and tends to 0 in the second stage of the consolidation process, while u gradually dissipates. When q( t) approaches 0, the consolidation of the clay layer is completed. Figure 2B shows the deformation in SDCL and DDCL. Obviously, the Q vs t curve also has the same variation characteristics as the q vs t curve. The Q( t) of SDCL is equal to Q(1, t), and the Q( t) of DDCL is equal to the sum of Q(0, t) and Frontiers in Earth Science frontiersin.org Q(1, t). The deformation of the clay layer increases significantly in the early stage of the consolidation process, and then gradually approaches a stable value. At the end of consolidation, the dimensionless deformation of DDCL is the same as that of SDCL. The time to complete consolidation for SDCL is four times longer than that for DDCL.

Results and discussion
In order to evaluate the above analytical solutions, the values of u, e, k, and Q of SDCL and DDCL were calculated by case analysis. Then, the values of these parameters based on the one-dimensional linear consolidation theory and the analytical solutions proposed in  Frontiers in Earth Science frontiersin.org 05 this study (based on the one-dimensional non-linear consolidation theory) were compared. According to previous studies (Wu et al., 2010;Luo et al., 2020), the main parameters of the clay layer used in the case study are as follows: C c = C k = 0.0532, σ 0 ′ = 51 kPa, e 0 = 1.101, k(e 0 ) = 2.10 × 10 −5 m/day, p = 100 kPa, and l = 10 m. In the one-dimensional linear consolidation theory, the parameters of the clay layer were considered to be constant, such as e 0 and k(e 0 ), and the coefficient of volume compressibility (m v ) of the clay layer was obtained by m v 0.434Cc (1+e0)σ0 ′ . The values of u at different times (t = 1, 3, 5, 10 years) and positions (a = 1, 3, 5, 7, 9 m) were calculated by the one-dimensional linear consolidation theory, and it was calculated by Eqs 11, 12 for SDCL and DDCL, respectively. The results are shown in Figure 3. Obviously, when t = 0, u is equal to the overlying load p. From Figure 3A and Figure 3C, the value of u in SDCL decreases with the increase in a and t. From Figure 3B and Figure 3D, the value of u in DDCL decreases with the increase in t. Meanwhile, the curve of u is symmetric around a =5 m, which means that the clay layer in the middle is the greatest. It was also found that the value of u predicted by Eqs 11, 12 was greater than that by the linear consolidation theory, which indicated that the dissipation of u predicted by the non-linear consolidation theory was lower. The stable time t s 0.75l 2 /c w predicted by the two theories for the completion of the consolidation of the clay layer is the same. From Figure 3, t s is about 21 years for DDCL, which is close to a quarter of that of SDCL. This is because the length of the drainage path in SDCL is l, which is approximately twice that of DDCL.
The values of e at different times and positions in SDCL and DDCL were calculated, and the results are shown in Figure 4. The results show that in SDCL, e decreases with the increase in a, while in DDCL, it is symmetric around a = 5 m. In both SDCL and DDCL, e eventually decreases to a stable value while the consolidation of the clay layer completes. The values of e f predicted by linear and nonlinear consolidation theories are 1.076 and 1.056, respectively, for the whole section of the clay layer. Moreover, the value of e predicted by non-linear consolidation theory is greater than that predicted by linear consolidation theory. The shapes of the graphs in Figure 4 are similar to those in Figure 3, that is, the variation of e predicted by non-linear consolidation theory is similar to that of u in the consolidation process of the clay layer.
The spatial and temporal distributions of k(e) in SDCL and DDCL were calculated by Eqs 15, 16, respectively. During the consolidation process of the clay layer, in SDCL, k(e) decreases with the increase in a and t ( Figure 5A; Figure 5C); while in DDCL, k(e) has a maximum value in the middle of the clay layer, and the upper and lower parts are symmetrical around a = 5 m ( Figure 5B; Figure 5D). For both conditions, k(e) of the clay layer decreases by 66.2% from the initial value of 7.7 × 10 −3 m/year to 2.6 × 10 −3 m/year when the consolidation of the clay layer completes. That is, the variation of k(e) corresponds to the changes in the void ratio in the consolidation process of the clay layer. The values of Q in SDCL and DDCL were calculated based on linear and non-linear consolidation theories, and the results are shown in Figure 6. During the consolidation process of the clay layer, the deformation gradually increases to a stable value using both theories. In both SDCL and DDCL, the deformation predicted by each theory ends up being equal. Among them, the deformation predicted by Eqs 23, 24 is 0.116 m, which is significantly smaller than the deformation predicted by the one-dimensional linear consolidation theory (0.216 m). The value of t s predicted by both theories are the same, 21 years for SDCL and 84 years for DDCL.

Parametric sensitivity analysis
Under overlying load conditions, the deformation [Q] of a clay layer depends on the hydraulic [k(e 0 )] and soil mechanical parameters [C c , e 0 , σ 0 ′ , l, and p] of the clay layer. Since the Q values were the same for SDCL and DDCL, we analyzed the parametric sensitivity of Q for DDCL by varying C c , e 0 , σ 0 ′ , l, k(e 0 ), and p values taken in the above case study one by one through non-linear and linear consolidation theories. Figure 7 shows that the value of Q predicted by Eq. 24 is smaller than that predicted by the one-dimensional linear consolidation theory. The difference in the final stable deformation calculated by linear and non-linear consolidation theories increases with increasing C c ( Figure 7A), l ( Figure 7D), and p ( Figure 7F) and decreasing e 0 ( Figure 7B) and σ 0 ′ ( Figure 7C), and it does not depend on the value of k(e 0 ) ( Figure 7E). The value of C c has a great effect on the Q value of the clay layer ( Figure 7A). The value of Q predicted by Eq. 24 increases non-linearly with the increase in C c , while that predicted by the one-dimensional linear consolidation theory increases linearly with the increase in C c . When C c = 0.  Predicted Q values in SDCL and DDCL. The impact of e 0 on Q is shown in Figure 7B. From the figure, the Q values of the clay layer predicted by Eq. 24 and the linear consolidation theory decrease non-linearly with the increase in e 0 . When e 0 is 0.5, 1, 1.5, and 2, the calculated Q values by Eq. 24 are 0.159, 0.120, 0.096, and 0.080 m, respectively, and the Q values from the one-dimensional linear consolidation theory are 0.302, 0.226, 0.181, and 0.151 m, respectively. The values of t s increase with the increase in e 0 . When e 0 is 0.5, 1, 1.5, and 2, the obtained t s values are 33. 1, 15.8, 14.4, and 13.2 years in DDCL, respectively.
The effect of σ 0 ′ on Q is shown in Figure 7C. When σ 0 ′ is 10, 20, 50, and 100 kPa, the calculated Q values by Eq. 24 are 0.250, 0.187, Frontiers in Earth Science frontiersin.org 0.115, and 0.073 m, respectively, and the calculated Q values by the one-dimensional linear consolidation theory are 1.10, 0.55, 0.22, and 0.11 m, respectively. The value of Q of the clay layer predicted by Eq.24 decreases non-linearly with the increase in σ 0 ′ , and that predicted by the one-dimensional linear consolidation theory decreases linearly with the increase in σ 0 ′ . The predicted t s decreases with the increase of σ 0 ′ . When σ 0 ′ is 10, 20, 50, and 100 kPa, the obtained t s values are 92.8, 33.1, 13.7, and 6.5 years, respectively. The impact of l on Q is shown in Figure 7D. The Q values predicted by both Eq. 24 and the linear theory increase linearly with the increase in l. When l is 5, 10, 15, and 20 m, the Q values predicted by Eq. 24 are 0.057, 0.114, 0.171, and 0.228 m, respectively, and that calculated by the one-dimensional linear consolidation theory are 0.108, 0.216, 0.323, and 0.431 m, respectively. The predicted value of t s increases with the increase in l. When l is 5, 10, 15, and 20 m, the t s values are 4.2, 22.8, 33.0, and 67.6 years, respectively.
The impact of k(e 0 ) on Q is shown in Figure 7E. The Q values obtained by Eq. 24 and the linear consolidation theory for all cases are 0.114 and 0.216 m, respectively. The values of t s decreases linearly with the increase in k(e 0 ). When k(e 0 ) values are 3.15 × 10 −3 , 6.31× 10 −3 , 1.58 × 10 −2 , and 3.15 × 10 −2 m/year, the t s values are 41.9, 20.8, 8.4, and 4.1 years, respectively.
The impact of p on Q is shown in Figure 7F. When p is 5, 10, 15, and 20 m, the Q values predicted by Eq. 24 are 0.072, 0.114, 0.144, and 0.167 m, respectively, and those by the one-dimensional linear consolidation theory are 0.108, 0.216, 0.323, and 0.431 m, respectively. The predicted value of Q by Eq.24 increases nonlinearly with the increase in p, and that by the one-dimensional linear consolidation theory increase linearly with the increase in p. The value of t s is not affected by p, and has the constant value of about 21 years by both Eq. 24 and the linear consolidation theory.

Laboratory experimental studies
In order to validate the analytical solutions proposed in this study, two laboratory experimental studies were performed for SDCL and DDCL. Besides, a type curve approach was developed to determine the clay layer parameters based on the experimental results.

Determination of clay layer parameters through type curve analysis
According to Eqs 23, 24, the deformation of the clay layer can be expressed as: where Q( t) is the dimensionless deformation of the clay layer. By changing Eq. 25 and dimensionless time ( t = c w t/l 2 ) into logarithmic forms, the following equations can be obtained: log 10 Q t ( ) log 10 Q t ( ) + log 10 2 1 + e 0 ( ) 2 log 10 t log 10 t + log 10 c w l 2 Eq. (26) and (27) are both decomposed into modes, where the corresponding variable and constant terms are added. Therefore, in the double logarithmic coordinate system, the shapes of the curves presented by the deformation Q and the dimensionless deformation Q are analogous. However, there is an offset in the position of the two curves in the coordinate system (c w /l 2 in the horizontal direction; 2(1+e0) 2 Ccl[1+(e0+e f )/2]log 10 (σ ′ f /σ0 ′ ) in the vertical direction), but the order of magnitude of the X and Y-axes corresponds to the same. Therefore, the two curves were translated to enable the data on the two curves to overlap as much as possible, and one overlapping point was selected as the matching point. Substituting the coordinates of the matching point (Q; Q, t,; t) into Eqs 28-30, the following parameters can be determined.

Results of the SDCL experiment
The single-sided drainage experiment carried out by Zhuo et al. (2022) was designed to test the one-dimensional consolidation of the aquitard under load. The thickness of the clay layer was 20 cm, the initial void ratio was 1.28, the load stage was 0-50-100-200 kPa, and the deformations of the clay layer corresponding to the three stages were 6.9, 6.2, and 6.6 mm, respectively.
In this study, the second-level load condition (50-100 kPa) was taken as an example to determine the hydraulic parameters of SDCL with σ 0 ′ = 50 kPa; l = 19.31 cm at this load stage. The deformation of SDCL and the dimensionless deformation curve with time are plotted

FIGURE 8
Determination of parameters of SDCL by the type curve approach.
Frontiers in Earth Science frontiersin.org in Figure 8. The type curve of Q( t) was then superposed onto this figure, by keeping the axes of the two graphs parallel to each other. The matching point refers to the intersection point of the type curve Q( t) with the data curve of Q(t). The coordinates of the matching point were Q 6.220 mm; Q 0.980; t 48,303 s,; t 9.105 in both systems. Substituting these values into Eqs 28-30, the values of c w ; C c ,; k(e 0 ) were calculated to be 7.03×10 −6 m 2 /s, 0.25, and 6.89×10 −8 m/s, respectively. Using the one-dimensional linear consolidation theory; c v was the same as c w , volume compressibility (m v ) and hydraulic conductivity were calculated to be 0.66 MPa -1 and 4.62×10 −8 m/s (Toufigh and Ouria, 2009), respectively. By the analytical solutions of Xie and Leo, (2004), the volume compressibility and hydraulic conductivity were calculated to be 0.66 MPa -1 and 4.60×10 −8 m/s, respectively. The calculated m v and k values from the analytical solutions of Xie and Leo, (2004) are close to those obtained by linear consolidation theory, and the k value predicted by Eq. 29 is larger than that determined by the other two methods.
Using the hydraulic parameters determined by the proposed type curve approach, the predicted and measured deformations of SDCL were compared under two load stages (50-100-200 kPa), as shown in Figure 9. From Figure 9, when the load is in the range of 50-100 kPa, the Q(t) values calculated by the analytical solutions in this study, the one-dimensional linear consolidation theory, and the analytical solutions by Xie and Leo, (2004) are consistent. Under this load stage, the calculated Q(t) values by the analytical solutions in this study, the one-dimensional linear consolidation theory, and the analytical solutions of Xie and Leo, (2004) are 6.27×10 −3 m, 6.35×10 −3 m, and 6.22×10 −3 m, respectively. In contrast, the

FIGURE 9
Comparison of predicted and measured deformation of SDCL.  The Q(t) value predicted by the analytical solutions in this study agrees well with the experimental results, while those predicted by the one-dimensional linear consolidation theory and the analytical solutions of Xie and Leo, (2004) are significantly larger than the experimental results. Among them, the value calculated by the analytical solutions of Xie and Leo, (2004) is smaller than that by the one-dimensional linear consolidation theory.

Results of the doubly draining clay layer experiment
The consolidation deformation characteristics of DDCL under overlying load conditions were studied by laboratory tests. The steps of the test were as follows: (a) The inner walls of the consolidation vessel were cleaned. The cohesive soil was air-dried, crushed, screened, milled, and dried. After cooling to room temperature, the soil sample was weighed with an electronic balance and put into a consolidation container (20 mm high). The soil sample was compacted in three layers with the corresponding heights of 6.67mm, 13.33mm, and 20mm, respectively ( Figures 10A-C). (b) The sample ( Figure 10D) was put into the rock vacuum saturated water test device (DP SJA, institute for Beijing Yaou Depeng Technology Co., Ltd., Beijing, China), the saturated water pressure was set to 20 kPa, the saturated pressure time was set to 300 min, the saturated water time was set to 240 min, and the immersion time was set to 2,440 min. After the soil saturation test was complete, the container of soil samples was removed and the soil samples were taken out ( Figure 10E). (c) The double-sided drainage consolidation test was performed on the saturated soil samples. A single-lever arm oedometer (Institute for Nanjing Ningxi Soil Instrument Co., Ltd., Nanjing, China), which met the requirements of the standard consolidation test, was used in the test ( Figure 10F). The load in the consolidation vessel was gradually increased through the following pressure sequence (e.g., 50, 100, and 200 kPa). The dial indicator reading was recorded during the test until the deflection stabilized. The initial thickness of the test soil sample was 2 cm, and the initial void ratio was 0.63. In order to avoid the influence of the stress state in the early stage caused by the difference in the previous sample preparation, the first load stage (0-50 kPa) was used as the pre-load stage in this test. Therefore, the characteristics of the clay layer, l, σ 0 ′ , and e 0 at the end of the consolidation under the load stage of 0-50 kPa are 1.82 cm, 50 kPa, and 0.54, respectively. The deformation of DDCL and the dimensionless deformation curve with time are plotted in Figure 11. The coordinates of the matching point for DDCL are Q = 0.632 mm, Q = 0.980, t = 3,600 s, and t = 11.218 in both systems. According to the analytical solution proposed in this study, c w = 1.03×10 −6 m 2 /s, k(e 0 ) = 1.06×10 −8 m/s, and C c = 0.18. In contrast, according to the onedimensional linear consolidation theory, c v was the same as c w , while m v and k(e 0 ) were calculated to be 0.71 MPa -1 and 7.32×10 −9 m/s, respectively. By the analytical solutions of Xie and Leo, (2004), the volume compressibility and hydraulic conductivity were calculated to be 0.71 MPa -1 and 7.3×10 −9 m/s, respectively. The values of m v and k calculated by the analytical solutions of Xie and Leo, (2004) are close to those obtained by the one-dimensional linear consolidation theory, and the k predicted by the analytical solutions in this study has a greater value than those predicted by the other two methods. Figure 12 compares the deformations of DDCL under two load stages (50-100-200 kPa) by using the hydraulic parameters determined by the type curve approach. Similar to the performance of the SDCL, when the load stage is 50-100 kPa, the calculated Q(t) values by the analytical solutions in this study, the one-dimensional linear consolidation theory, and the analytical solutions of Xie and Leo, (2004) agree with the experimental results. Under this load stage, the calculated Q(t)

FIGURE 11
Determination of parameters of DDCL by the type curve approach.

FIGURE 12
Comparison of predicted and measured deformation of DDCL.
Frontiers in Earth Science frontiersin.org by the analytical solutions in this study, the one-dimensional linear consolidation theory, and the analytical solutions of Xie and Leo, (2004) are 6.03 × 10 −3 m, 6.45 × 10 −3 m, and 6.32 × 10 −3 m, respectively, while the measured value of the experiment is 6.10 × 10 −3 m. During the load stage of 100-200 kPa, the Q(t) value predicted by the analytical solutions in this study agrees well with the experimental results, while the Q(t) values predicted by the one-dimensional linear consolidation theory and the analytical solutions of Xie and Leo, (2004) are significantly larger than that predicted by the analytical solutions in this study. Under the load stage of 100-200 kPa, the calculated Q(t) values by the analytical solutions in this study, the one-dimensional linear consolidation theory, and the analytical solutions of Xie and Leo, (2004) are 6.23 × 10 −3 m, 1.29 × 10 −2 m, and 1.20 × 10 −2 m, respectively, while the measured values obtained in the laboratory test is 6.25 × 10 −4 m. Therefore, based on these experimental results, the analytical solutions in this study predict the deformation of the clay layer more accurately than the one-dimensional linear consolidation theory and the analytical solutions of Xie and Leo, (2004).

Conclusion
In this study, analytical solutions (ignoring creep effects) for SDCL and DDCL undergoing non-linear consolidation, excess pore water pressure, and deformation under overlying load conditions were derived. Combined with the case analysis based on the analytical solution presented in this paper, it is found that under the condition of an overlying load, the deformation in SDCL is equal to that in DDCL. There is a lag in t s for SDCL compared to DDCL. The drainage path of the consolidated drainage of SDCL is twice as large as that of DDCL, resulting in a t s lag that is about four times that of DDCL. In addition, there is a difference in u under the two theories, i.e., the u value by the non-linear consolidation theory is larger than that by the linear consolidation theory. The deformation predicted by the non-linear consolidation theory is smaller than that by the linear consolidation theory. Larger C c , l, and p and smaller e 0 and σ 0 ′ lead to a larger difference in the amount of deformation under the two theories. The t s predicted by the two theories are the same. Through the sensitivity analysis of the parameters, it is obvious that the deformation in SDCL and DDCL is controlled by C c , e 0 , σ 0 ′ , l, and p. However, the variation of k(e 0 ) does not affect the deformation in SDCL and DDCL, and t s does not depend on the overlying load of the clay layer. Besides, a type curve approach is proposed in this study based on the non-linear consolidation theory, and it agrees well with the measured data from the test. Compared with the linear consolidation theory and the analytical solution of Xie and Leo, (2004), the analytical solution based on the non-linear consolidation theory in this study is more accurate and suitable for predicting the consolidation and deformation of the clay layer under further load conditions.

Data availability statement
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material.