Pore Pressure Response to Groundwater Fluctuations in Saturated Double-Layered Soil

Analytical solutions are developed for one-dimensional consolidation of double-layered saturated soil subjected to groundwater fluctuations. The solutions are derived by an explicit mathematical procedure using Duhamel’s theorem in conjunction with a Fourier series, when groundwater fluctuation is described by a general time-dependent function and assumed to be the pore water pressure variations at the upper boundary. Taking as an example the harmonic groundwater fluctuation, the relevant response of the excess pore water pressure is discussed in detail, and the main influencing factors of the excess pore pressure distribution are analyzed. A dimensionless parameter θ has been introduced because it significantly affects the phase and the amplitude of excess pore pressures.The influences of the coefficients of permeability and compressibility of soil on the excess pore pressure distribution are different and cannot be incorporated into the coefficient of consolidation in double-layered soil. The relative permeability ratio of two clayey soils also plays an important role on the curves of the distributions of the excess pore pressures. The effects of the thickness of the soil layer on the excess pore pressure distribution should be considered together with the dimensionless parameter θ and the permeability and compressibility of the double-layered soil system.


Introduction
Excess pore water pressures are generated within saturated clayey soil due to groundwater fluctuations arising from sources such as river or ocean tides, seasonal or climatic water table changes, and human factors (e.g., pumping or recharging of water) [1].The more common analysis models for the generation and propagation of excess pore water pressures can be generalized as double clayey soil layers instead of a single layer of soil with pervious top and impervious bottom boundary conditions (PTIB) (Figure 1(a)) or with both pervious top and bottom boundary conditions (PTPB) (Figure 1(b)) [2].On account of the marginal difference between the natural unit weight and the saturated unit weight of a soil layer with high permeability, the total stress in the clay soils generally varies slightly and, accordingly, can be regarded as constant.Therefore, the groundwater table variation in the top sand layer with high permeability can reasonably be considered as the excess pore water pressure variation at the upper boundary of clayey soil layers.
In addition, the evolution of excess pore pressures with time in clayey soil layers generally do not reflect the changes occurring at the boundary.
Conte and Troncone [1] achieved analytical solutions of the excess pore pressures in single clayey soil caused by pore pressure varying at the boundary, which was conducted based on the calculation procedure for the solution of a sinusoidal pore pressure variation at the boundary [3,4].Many forms of analytical procedures of one-dimensional consolidation for layered systems have been provided by other authors.Schiffman and Stein [5] developed a general solution for the one-dimensional consolidation problem that considered a variety of boundary conditions, initial conditions, loading histories, and types of multilayered soil; however, this solution is quite tedious to use. Lee et al. [2] derived a general analytical solution for the problem of layered soils with a time-dependent load, which was more explicit than the one given by Schiffman and Stein [5].Luan and Qian [6] developed an analytical solution for the case of double-layered soil profile with a depth-dependent ramp load.Xie [7] demonstrated a mathematical procedure for an analytical solution of one-dimensional consolidation subject to the time-dependent load and depth-dependent initially nonlinear pore pressure and subsequently developed a procedure for a multilayered system [8].Xu et al. [9] applied the Laplace transform and the Fourier technique to solve a one-dimensional consolidation problem for layered soil under arbitrary loading.An analytical study involving one-dimensional consolidation of double-layered soil with partially drained boundaries was provided by Xie et al. [10].The recent research studies on analytical solutions for consolidation of porous media to predict the variations of the pore pressures and settlement with time were also proposed by Ho et al. [11,12].Hydrogeologists studied the groundwater table response to tidal fluctuation in a coastal aquifer [13,14].
The studies of tide-induced pore pressure were conducted by van der Kamp and Gale for a homogeneous porous medium [15] and by Wang and Davis for a multilayer medium [16].However, no analytical solution was found for onedimensional consolidation with such a complex upper boundary condition that takes the layered characteristics of soil into account.Analytical solutions are proposed for one-dimensional consolidation of double-layered saturated soils subjected to groundwater fluctuations for the simplest multilayer system: the double-layered system.Both PTIB and PTPB boundary conditions are taken into consideration in this paper.Groundwater table variation is assumed to be the pore water pressure variation at the boundary described by a general time-dependent function.The solution is derived by a more complicated mathematical procedure using Duhamel integration in conjunction with a Fourier series.The response of excess pore pressure to groundwater fluctuation and its main influencing factors were analyzed using the proposed solution, which can be significantly different from the solution of a one-layered system.

Solution Method
2.1.Governing Equations.With all the assumptions of Terzaghi's one-dimensional consolidation theory [17], except for the pore water pressure at the upper boundary, the equation governing the generation and propagation of excess pore water pressures in clayey soil layers of this analytical model (Figure 1 where  is the spatial coordinate;  is time;  1 (, ) and  2 (, ) are the excess pore water pressures in the upper layer (at 0 <  ≤  1 ) and the lower layer (at  1 <  ≤ ) of the doublelayered clayey soil system, which depend on both  and ;  vs1 and  vs2 are the coefficients of swelling or consolidation of the upper soil layer and the lower soil layer, respectively, which are formulated by the expression  vs =  vs /(   vs ) ( = 1, 2);  vs and  vs are the coefficient of permeability in the vertical direction and the coefficient of volume change of the soil, respectively.In practice,  vs as well as  vs and  vs are generally accepted to be constant during the process [4].
According to the analytical model assumptions, the excess pore water pressure on the upper boundary of the double-layered clayey soil system is first assumed to be harmonic with time, which is expressed as follows: where  and  determine the excess pore pressure amplitudes and  is the angular frequency, whose functional relation with the period  is  = 2/.The bottom boundary of the double-layered clayey soil system may be fully permeable or impermeable; thus, where  is the total thickness of this double-layered soil model.
In addition, the seepage and excess pore water pressures of the double clayey soil layers are continuous: where  1 is the thickness of the upper layer soil of doublelayered soil model.Moreover, the initial condition is Overall, the solution of this analytical model is formulated by (1), which should satisfy the above boundary conditions ( 2) and (3a) or (3b), (4), and the initial conditions in (5); the solution follows the analytical procedure of solving the partial differential equation with complicated inhomogeneous boundary conditions.

Analytical Solution of Harmonic
Groundwater Fluctuations.Equation ( 2) implies that the excess pore water pressure at the upper boundary of the model varies harmonically with time.Furthermore, the solution procedure described as Duhamel's theorem [1,3,18,19] is used to derive (1) with the boundary and initial conditions given by ( 2) and (3a) or (3b), (4), and (5), which is obtained by the following integration: where  1 (0, ) is determined by ( 2) and   (, ) denotes the solution to ( 1) and (3a) or (3b), (4), and ( 5), when the excess pore water pressure at the upper boundary is kept at unity (as the value  1 (0, ) = 1).First, the excess pore water pressure in double clayey soil layers with pervious top and impervious bottom boundary conditions (Figure 1(a)) is systematically obtained in this section.In particular,   (, ) is reduced to a separable partial differential equation: where  is a function of  and  and  is a function only of ;  1 and  1 are for 0 <  ≤  1 and  2 and  2 are for  1 <  ≤ .Subject to all the assumptions identified in a previous section,   and   should satisfy the following equations: 1 (0, ) = 0, ( 14) The function   () is determined by (8), with the boundary conditions in ( 9)-( 12); the result takes the following form: A more complex mathematical procedure is used to derive the solution for (13), with the boundary conditions in ( 14), ( 15), (17), and ( 18) and the initial conditions in (16), in which   (, ) can also be reduced to a separable partial differential equation, expressed as follows: Substituting (20) into differential equation ( 13) and introducing a constant , which is defined as the eigenvalue, we have Deriving (21) yields where ,  1 ,  2 ,  1 , and  2 are undetermined coefficients.
Based on the boundary conditions in ( 14) and ( 15), coefficients  1 and  2 are given in the forms Substituting (24) into the boundary conditions in ( 17) and (18) yields To guarantee nontrivial solutions of  1 and  2 , the determinant of the second-order correlation matrix composed of (25) is identically vanishing.Consequently, the eigenequation is obtained as follows: in which  =  1 +  2 ,  2 is the thickness of the lower layer soil of double-layered soil model;  = √ vs2  vs2 / vs1  vs1 is defined as the interlaminar seepage-compression comprehensive parameter ratio [6].The eigenvalue   of (26), which is proved to be infinite, can be derived by an iteration process in general, except for some special conditions; for example, when  vs1 =  vs2 and  vs1 =  vs2 , (26) yields   = (2 − 1)/2.Considering ( 24) and ( 25), in conjunction with ( 22), (23a), and (23b), the general solution of (20) can also be expressed in series form where   is defined as the modal participating factor; the eigenfunction   is represented as follows: The above eigenfunctions are proved to be orthogonal to the coefficient of the volume change of the soil  vs by the following equation: where and   is defined as follows: In particular, the orthogonal property demonstrated in (29) is familiar and useful in dealing with the mathematical problems of the double-layered system [6,7,10].
Moreover, an inspection of the initial condition in (16) combined with (29) reveals Substituting ( 19) and ( 26) into (31a), after some algebraic manipulation, renders   , such that Consequently, the solution for   (, ) is developed, and the final form for   (, ) is represented as follows: where  1 ,  2 , and   are determined by (28a), (28b), and (31b), respectively.Meanwhile, the parametric coefficient   is determined by (26).Therefore, substituting  for  in (2) and − for  in (7), at the same time performing the partial derivative of   (,  − ) with respect to , (6) is developed into the following simplified calculation formulas: where and  V denotes the time factor of Terzaghi's theory.
Furthermore, the excess pore water pressures in double clayey soil layers with both pervious top and bottom boundary conditions (Figure 1(b)) can be obtained using a similar method, when the corresponding boundary conditions are described in (2) and (3b).After performing the same steps previously described, the solution for   (),   (, ) is where the eigenvalues   ,  1 (),  2 (),   , and   are determined by the equations Based on (35a) and (35b), the final solutions for (6) in this case are given by where and  and  V are defined by (33f) and (33g), respectively.Note that (33a) and (33b) or (36a) and (36b) are applicable only when 0 <  ≤ ; for simplification, ( 2) is used to calculate the values of excess pore water pressure at  = 0 [1].

Analytical Solution of General Groundwater Fluctuations.
Based on the derived equations, an analytical solution to (1) is derived using the Fourier series, when a more general time-dependent excess pore pressure () due to an arbitrary groundwater fluctuation is applied to the top surface of double-layered soil.
First, () should be expanded in harmonic components using the Fourier series theorem, as outlined by Xu et al. [9] and Conte and Troncone [1,18], such that in which the series amplitudes   and   associated with the frequency   = 2/ (with  = 1, 2, . .., and  is the period of ()) are represented as follows, respectively: and  0 can be obtained from (10) by setting  = 0; thus, Second, the corresponding excess pore water pressure for each of the harmonic components in (37) is calculated by (33a) and (33b) or (36a) and (36b) in which   (, ) is indicated as   (, ),  and  are substituted by   and   , respectively, and  is replaced by   .Finally, the actual excess pore water pressure at any time and depth is formulated by using the principle of superposing all the components of the determined   (, ): in which  0 is represented by (40) and   (, ) is represented by (32a) and (32b) or (35a) and (35b).The analytical solutions are so far available for any arbitrary groundwater fluctuation of double-layered system.[1] and Ying et al. [19] assumed that there is only a single clayey soil layer.In this section, these solutions are demonstrated to be special cases of the solution presented in this paper.Taking as an example the PTIB boundary condition, if there is a single clayey soil layer ( vs1 =  vs2 ,  vs1 =  vs2 ), from (26), the expression for   becomes

Comparison with Existing Analytical Solutions. The solutions by Conte and Troncone
Substituting (42) back in to (31b) yields lim The combination of (33a), (33b), (42), and (43) directly leads to which is essentially the same as the equation derived by Conte and Troncone [1] and Ying et al. [19].The same procedure can be conducted for the case with PTPB boundary condition; the resulting equation takes the form  (, ) = lim where It may be meaningful to distinguish the difference between the behavior of the single system and the layered system.In comparison with the expressions for the single system, the ones of the double-layered system are obviously much more complicated.Equations ( 44) and ( 45) indicate that the dimensionless parameter  may be the primary influencing factor of excess pore pressure distribution in single clayey soil layer, which has been verified by Ying et al. [19].The expressions of excess pore pressures in (33a)-(33g) and (36a)-(36e) indicate that four important parameters, the dimensionless parameter , the permeability ratio of the soil, the compressibility ratio of the soil, and the thickness ratio of the soil layers, are involved in the model.Those parameters may have a significant influence on the excess pore pressure distribution of the double clayey soil layers.The discussion above can be the basis of the following further discussion.

Influence of Various Parameters on Excess Pore Water
Pressures.In this section, we pay principal attention to the influence of various parameters on the distribution of excess pore water pressure subject to groundwater fluctuation, for the purpose of understanding the behavior of the layered system and being aware of the potential complications.
For simplicity, let  = 0 in (2), with the assumption that the groundwater fluctuation is a cosinusoidal wave, and define the following dimensionless parameters: Typical values of those four important parameters (, , , and ) were selected to study their influence on the distribution of excess pore water pressures.In addition, a computation program was developed, and the main results are shown in Figures 2 to 6.The distributions of excess pore water pressure with the PTIB boundary condition and with the PTPB boundary condition are obviously different.However, the influence laws of these important parameters on the cases with PTIB and PTPB boundary conditions are alike.For the sake of simplicity, only the cases with PTIB boundary condition are analyzed in Figures 3 to 6.
Figure 2 shows the variation of the excess pore water pressure with dimensionless parameter  for both PTPB and PTIB boundary conditions in clayey soil layers at  =  when other dimensionless parameters , , and  remain constant (i.e.,  =  vs2 / vs1 = 0.1,  =  vs2 / vs1 = 0.1, and  =  2 / 1 = 1).As can be observed from Figure 2, the greater the dimensionless parameter , the greater the excess pore water pressures and the propagation depth of the excess pore water pressures.When the value of  is comparatively low, the excess pore water pressures attenuate rapidly with depth, the propagation depth of excess pore water pressures is relatively small, and the hysteresis of the propagation of excess pore water pressures exists in the deeper depth.When the dimensionless parameter  increases and becomes greater than a certain value (say  = 10 in this case), the distribution of excess pore water pressures better fits the linear distribution.Obviously, the propagation of excess pore water pressures is related to the dimensionless parameter , which is similar to the characteristics of the propagation of excess pore water pressures in single clayey soil [19].However, there exist the inflection points of the isochrones of the excess pore water pressures in a double-layered system other than a single-layered system.
Figures 3 and 4 show the influence of the permeability and the compressibility, respectively, of soil on the propagation of excess pore water pressures for PTIB boundary condition in clayey soil layers at  =  when other dimensionless parameters ,  are fixed (i.e.,  = 0.2 and  =  2 / 1 = 1).As shown in Figure 3, the influence of the dimensionless ratio  on the excess pore water pressures in the upper and the lower layers is opposite when the compressibility of two layers is the same.The greater the ratio , the smaller the excess pore water pressure in upper layer; in contrast, the behavior is opposite for the lower layer.The inflection points of the excess pore water pressure distribution curves exist in the interface between two clayey soil layers transparently, except for  = 1, which implies that the permeability of the two layers is the same.Figure 4 indicates that the inflection points of the excess pore water pressures distribution curves disappear if the permeability of two layers is the same.As expected, the excess pore water pressure isochrone moves gradually to the left, indicating that higher values of the dimensionless parameter  results in higher attenuation of the excess pore water pressures.The data in Figures 3 and 4 implies that the hysteresis of the propagation of excess pore water pressures exists when the second clayey soil layer is with lower permeability or higher compressibility.Figure 5 shows the influence of the dimensionless parameter  on the propagation of the excess pore water pressures for the PTIB boundary condition in clayey soil layers at  = .Obviously, for the case where the permeability in the lower layer is higher than the permeability in the upper layer (e.g.,  = 5 in Figure 5(a)), the thickness of the lower layer with higher permeability is inversely correlated to the excess pore water pressure attenuation in the lower layer and directly correlated to the excess pore water pressure attenuation in the upper layer.When the permeability ratio of the doublelayered system is less than 1 (e.g.,  = 1/5 in Figure 5(c)), the thickness of the lower layer with lower permeability is directly correlated to the excess pore water pressure attenuation in the lower layer and inversely correlated to the excess pore water pressure attenuation in the upper layer.When the lower layer has higher compressibility (e.g.,  = 5 in Figure 5(c)), the thickness of the lower layer is directly correlated to the excess pore water pressure attenuation, while the opposite results occur for the cases of the lower layer with a lower compressibility (e.g.,  = 1/5 in Figure 5(d)).The results in Figure 5 reveal that the thickness of the lower layer with greater permeability or smaller compressibility is directly correlated to the easy of the excess pore pressure propagation in the double-layered system.
Obviously, permeability, which is the property of the medium, significantly affects the diffusion of the excess pore water pressure in two clayey soils.As can be inferred from the results of Figures 2-5, the slope breaks exist only when the permeability of the two clayey layers are different.
setting  = 0 in (2).The results of the analyses show that the influencing factors for excess pore pressure distribution include the dimensionless parameter , the permeability and the compressibility of soil, and the thickness of soil layer.The dimensionless parameter  has a significant influence on the excess pore pressures distribution, and the excess pore water pressure and the influencing depth of excess pore water pressures are inversely correlated to the phase shift and the amplitude attenuation.Other than in single layer soil, the influence of the coefficients of permeability and compressibility of soil on excess pore pressure distribution is different and cannot be incorporated into the coefficient of consolidation of the soil in double-layered soil.The relative permeability ratio of the two clayey soils makes the excess pore pressure/depth curves of double-layered systems fundamentally different from those of uniform soil.Furthermore, the effects of the thickness of soil layer on excess pore pressure distribution depend on the dimensionless parameter  and the coefficients of permeability and compressibility of the double-layered soil system.

Notation
: Fourier cosine series amplitude associated with the frequency    0 : Fourier cosine series amplitude when   = 0 , : Amplitudes of the excess pore pressure when it varies harmonically with time   : Fourier sine series amplitude associated with the frequency    vs : Coefficient of swelling/consolidation of the soil (): Time-dependent function describing the pore pressure variations at the boundary : Soil layer thickness  vs1 : Coefficient of permeability in the vertical direction of the upper layer soil  vs2 : Coefficient of permeability in the vertical direction of the lower layer soil  vs1 : Coefficient of volume change of the upper layer soil  vs2 : Coefficient of volume change of the lower layer soil : Period of the time-dependent function describing the pore pressure variations at the boundary  V : Time factor : T ime : Excess pore water pressure  1 (, ): Excess pore water pressure at depth  and time  in the upper layer soil  2 (, ): Excess pore water pressure at depth  and time  in the lower layer soil  1 (, ): Excess pore water pressure at depth  and time  when the excess pore water pressure at boundary is kept at unity in the upper layer soil  2 (, ): Excess pore water pressure at depth  and time  when the excess pore water pressure at boundary is kept at unity in the lower layer soil  1 (, ): th harmonic component of the excess pore water pressure of the upper layer soil  2 (, ): th harmonic component of the excess pore water pressure of the lower layer soil : S p a t i a lc o o r d i n a t e   : Unit weight of water : Cir cula rfr eq uen cy   : Circularfrequencyoftheth harmonic component of Fourier series   ,   : Eigenvalue.

2 MathematicalFigure 1 :
Figure 1: One-dimensional consolidation model of two clayey soil layers due to groundwater fluctuations for (a) PTPB and (b) PTIB boundary conditions.

1 Figure 4 :
Figure 4: Influence of the compressibility on the excess pore water pressures (PTIB).

Figure 5 :
Figure 5: Influence of the dimensionless parameter  on the excess pore water pressures (PTIB).