Numerical modeling of contaminant advection impact on hydrodynamic diffusion in a deformable medium

Chungung Wng, Shuqing Liu, Xingki Shi, Gunglei Cui, Hongxu Wng, Xing Jin,Kunkun Fn, Songto Hu

a College of Energy and Mining Engineering, Shandong University of Science and Technology, Qingdao, 266590, China

b Key Laboratory of Ministry of Education on Safe Mining of Deep Metal Mines, Northeastern University, Shenyang,110004, China

c Shandong Provincial Geo-Mineral Engineering Co., Ltd., Jinan, 250013, China

Keywords:Consolidation Advection Diffusion Sorption

A B S T R A C T

1. Introduction

Waste landfill can produce large quantities of leachates containing high concentrations of various pollutants(Bou-Zeid and El-Fadel, 2004). Once these pollutants leak from the landfill, both surrounding soil and groundwater can be severely contaminated(Ward et al.,2005;Sang et al.,2006).Various compacted clay liners and composite liners are frequently used to isolate leachate and pollutants from the surrounding soil(Chen et al.,2009).In practice,as the landfill expands, the vertical stress applied on solid waste increases due to the mechanical rolling compaction of the waste or from the mass of the waste itself.These processes can compact the contaminated layer and even destroy the coherence of the liner.Therefore, understanding the mechanism of the consolidationinduced transport of contaminants is critically important for the design of impermeable liners and long-term prediction of landfill leakage.

Contaminant transport models are generally derived from the assumption that the clay or composite liners are of an incompressible material, in which only solute advection or dispersion transport occurs(Ogata,1970;Bear,1972;Freeze and Cherry,1979;Chen et al., 2007; Ciftci, 2017). Relevant solutions of partial differential equations are solved by finite element methods (e.g. Javadi and Al-Najjar, 2007; El-Zein et al., 2012; Ma et al., 2021), boundary element methods(e.g.Tan et al.,1994;ˇSarler and Kuhn,1998a,b;Al-Bayati and Wrobel,2021),spectral methods(e.g.Zhao and Liu,2017; Chen et al., 2018) and meshless methods (e.g. Gingold and Monaghan,1977; Liu et al., 1995; Babuˇska and Melenk, 1997; Liu and Gu, 2002; Xiao, 2004). In particular, (Lin et al., 2020, 2021)developed a backward substitution-based meshless method to solve 1D and 2D linear and nonlinear advection-diffusion-reaction problems. This kind of advection-diffusion problem is applicable for low-compressibility cases (Foose, 2002; Rowe, 2005;Lewis et al., 2009a, b; Shackelford, 2014; Xie et al., 2018). Intermittent waste filling operations can alter the boundary stress conditions on the bottom liner (Zekkos et al., 2006; Yu and Rowe,2018; Rowe and Yu, 2019). In these cases, increasing boundary stress conditions can consolidate the landfill body, and the associated transient flow of contaminant advection may dominate the solute transport process. Attention on the long-term effect of the solid-fluid interaction should be paid to complete characterization.

On the basis of Terzaghi(1925)and Biot’s consolidation theories(Biot, 1941), a series of contaminant transport models has been developed (e.g. Loroy et al.,1996; Smith, 2000; Peters and Smith,2001, 2004; Lee et al., 2009; Lee and Fox, 2009; Lewis et al.,2009a, b; Zhang et al., 2012, 2013, 2021; Pu et al., 2018). These models were developed by coupling small strain based consolidation with solute advection-dispersion transport. Some research groups also pay attention to solute sorption (Smith, 2000; Peters and Smith, 2020; Yan et al., 2021), cation/anion-transport (Van Impe et al., 2002), osmotic consolidation (Barbour and Fredlund,1989), chemical consolidation (Lee et al., 2009; Banaei et al.,2021), and single/double-drained solute transport (Alshawabkeh et al., 2004, 2005). Breakthrough time of effluent contaminant penetrating through large strain clay becomes shorter than those of incompressible or small strain conditions(Peters and Smith,2001).It indicates that assumptions of large strain and nonlinear constitutive behavior are typically required for high compressible materials.

Large strain based consolidation theory, involving material coordinate frame,boundary conditions and dependent variables,was first discussed by Gibson et al. (1990). Considering that mass conservation is strictly enforced by following the motions between the solid phase and fluid phase throughout the consolidation process,Smiles (2000) used a variety of material coordinate systems to analyze the movement of solutes in consolidating clays.Peters and Smith (2001) further extended consolidation-induced solute transport in both spatial and material coordinates (Smith, 2000).However, these large strain simulations did not involve with material self-weight, mechanical dispersion and sorption/desorption.

Fox and Berles (1997) used a piecewise-linear approach to stratify saturated layered soil and established a 1D consolidation settlement 2 (CS2) model. The CS2 model is capable of simulating the relative velocity of fluid and solid phases, as well as variable hydraulic conductivity and compressibility during the large strain consolidation process.An enhanced version of CS2 model makes it additional advantageous for incorporating time-dependent loading and an external hydraulic gradient.Fox(2007)extended CS2 model to consolidation and solute transport 1 (CST1) model for coupling consolidation and solute transport in saturated porous medium.Associated solute transports involve advection, dispersion, firstorder decay reactions, and linear equilibrium sorption. In particular,CST1 model adopts two Lagrangian fields of elements to model consolidation-induced motions of fluid and solid phases separately.Fox and Lee (2008) then incorporated additional capabilities into CST1 to CST2 model,which can accommodate ion diffusion during consolidation and nonlinear nonequilibrium sorption. Varying effective diffusion coefficient can significantly affect solute transport, whereas non-equilibrium sorption becomes more important with the on-going consolidation. Using CS2, Fox et al. (2014)developed CS3 model for large strain consolidation of uncontaminated layered soils. The CS3 model was extended by considering layered settlement,time-dependent unload/reloading and variable external hydraulic gradient.Based on CS2-CS3-CST2 models,Pu and Fox(2016)proposed CST3 model which can couple 1D large strain consolidation with solute transport in layered soils.The Lagrangian framework can make it advantageous for tracking multilayer interfaces motion, and satisfying fluid mass balance.

The aforementioned studies demonstrate that there is a strong correlation between the consolidation degree and solute transport.The effluent solute flux is synchronized with temporal-spatial variations in porosity of the consolidating layer. Nevertheless, the excess pore pressure dissipation can delay the duration of the consolidation process, and in turn generate a transient flux of contaminants advection (Alshawabkeh et al., 2004). However,although the consolidation-induced hydraulic gradient driving solute transport has been investigated, little attention has been paid to ion migration caused by concentration gradient variation.

Based on the large strain consolidation theories, numerical solutions for coupled variable external loading and contaminant transport in saturated clay layers were performed in this study.The numerical solution proposed for the clay consolidation process aimed to analyze elastoplastic deformation, nonlinear compressibility-permeability relationships, and its solution for solute transport models advection, mechanical dispersion, linear equilibrium sorption, and porosity-dependent effective diffusion coefficients. The numerical results are validated using the published laboratory data in the study of Lee et al. (2009) and Lee and Fox(2009). This model further investigates the excess pore pressure generation during the consolidation process of kaolinite.It includes the flux of solutions resulting from the overall settling and followup drainage, and the correlation between advection and hydraulic dispersion is discussed based on the spatio-temporal variation of the solute concentration field within kaolinite. Finally, ionic sorption effects on contaminate transport are investigated using parametrical analyses.

2. Mathematical models

Solute transport in clay involves complex physico-chemical processes, such as advection, dispersion, diffusion and sorption.The following assumptions are used for the mathematical model for contaminant transports in the clay liner: (1) soil is isotropic,homogeneous,and saturated;(2)bulk flow velocity of the solution is resultant of Darcy flow and chemiosmosis flow; (3) hydraulic conductivity accounts for tortuosity between clay particles;(4)the modified Cam-clay model is used for describing soil consolidation;and(5)soil consolidation occurs only in the vertical direction.Fig.1 is a schematic diagram of the solute transport manner in a finegrained soil. The orange area represents the soil particles. The blue gradient area denotes solute concentration distribution, in which the dark blue has a high concentration value and light blue has a low concentration value.

Fig.1. Schematic diagram of soil consolidation-induced solute transport manners.The green circle represents reactive ion (K+), and the blue circle represents non-reactive ion (Br-). Green arrow denotes pore water seepage regime, and light blue arrow denotes pore water osmosis regime.

2.1. Flow equation

The flow equation of the solution through porous media can be derived from the mass conservation equation (Kaczmarek, 2001):

where φ is the porosity of the porous media,ρfis the liquid density,υr= υ/φ is the relative specific discharge,υ is the pore fluid seepage velocity, and υsdenotes the velocity of solid phase.

By incorporating the interactions among the pore hydraulic gradient,gravity and leachate hydraulic gradient,the fluid seepage velocity can be expressed by Darcy’s law (Chang et al., 2019):

where υdis the Darcy velocity,Khis the hydraulic conductivity,γw= ρfgis the bulk density of water,Pfis the pore water pressure,His the unit vector in the direction of gravity, andiis the water head gradient.

The chemiosmosis effect induced by the concentration gradient can be expressed by

where υcis the chemiosmosis rate;χ is the chemical permeability coefficient (Witteveen et al., 2013); andPC=RTΣCiis the chemical osmotic pressure caused by concentration difference(Malusis and Shackelford,2002),in whichRis the Mohr coefficient,andTis the temperature.

The overall velocity of fluids can be expressed by the sum of the Darcy velocity and the chemiosmosis velocity(Neuzil, 2000):

Eq. (1) can be expanded as

Eq. (5) can be simplified by using the material derivative equation (Smith, 2000):

where ds()/dtis the material derivative; ∂()/∂tis the local derivative, which represents the rate of change of the fluid physical quantity with time at a fixed point in space; and υs∇() is the convective derivative, which refers to as the change of a physical quantity from one place to another with different physical quantities.

By substituting Eq. (6) into Eq. (5), we have

Eq. (7) can be rewritten as

It can be obtained from Eq. (6) that

Substituting Eq. (9) into Eq. (8), Eq. (10) can be written as

In order to simplify Eq.(10),the compression coefficient αnand βPare introduced, where αn= -dV/(VdPf) = dn/[(1-φ)dPf]represents the change of porous media skeleton with fluid pressure,and βP= dVw/(VwdPf) = dρf/(ρdPf)represents the volume compressibility of fluid.Eq.(11)can be obtained by introducing the compressibility factor into Eq. (10).

2.2. Deformation equation

According to the modified Cam-clay model (Krabbenhoft and Lyamin, 2012), the yield function can be written as

whereq,pare the stress invariants;I1is the first stress invariant;J2is the second deviatoric stress invariant; andMis the slope of the critical state line in P-Q space determined by internal friction angle;pcis the consolidation pressure, which is determined by plastic volume strain;pc0is the initial consolidation pressure;erefis the reference void ratio; κ, λ are the swelling index and the compression index; and εp:volis the volumetric plastic strain.

The deformation process is described by mechanical equilibrium equation (Maltseva et al., 2017).

where σij,iis the normal stress,fiis the volume force inidirection,andfy= [(1 - φ)ρs+ φρf]g.

When contaminants enter the porous media, the density changes of solid and pore fluid as follows:

where ρf0and ρs0are the initial densities of the liquid and solid phases, respectively;mfcandmsccare ion mass of liquid phase and solid phase, respectively;cis the contaminant concentration;F=kdcis the concentration of adsorbed contaminants for per unit mass of the solid phase(Yadav and Singh,2017),in whichkdis the adsorption coefficient.

The effective stress in the solid phase is expressed as

where αBis the Biot coefficient.

Hydrostatic pressure is introduced to describe the nonlinear relationship of the volumetric elastic strain function. The constitutive relation of saturated soil is (Sanei et al., 2020):

whereprefis the reference stress,pmaccis the effective mean stress,Gis the shear modulus,Ktis the tangent bulk modulus(Zhong et al.,2014), andeis the void ratio.

The elastic strain can be obtained by subtracting the strain caused by the concentration change and plastic strain from the total strain of the porous media:

where εijis the total strain and εij= (ui，j+uj，i)/2(Arrue et al.,2016),while εplis the plastic strain and εpl= λ∂Q/∂S, in which the associated flow rule(Chang et al.,2019)Q=Fycan be used;and εcis the strain caused by concentration change(Wallmersperger et al.,2004):

wheregcis the material properties; andgc= δi，ja, in whichais coefficient of chemical expansion.

Combining Eqs.13-18 yields the general Navie-type equation:

where the third and fourth terms on the left side of the equation are respectively the plasticity and concentration; and the first and second terms on the right side are the pore pressure and the body force,respectively.

2.3. Transport equation

2.3.1. Hydrodynamic dispersion

Hydrodynamic dispersion involves molecular diffusion and mechanical dispersion, both of which can be expressed by Fick’s law. Therefore, the hydrodynamic dispersion equation can be expressed as (Ciftci, 2017):

whereJdis the hydrodynamic dispersion flux; andDis the coefficient of hydrodynamic dispersion,which is theoretically the sum of the effective diffusion coefficient,Dc, and the coefficient of mechanical dispersion,DL, so that we haveD=DL+Dc.

whereDLis the longitudinal dispersion and τ is the tortuosity of the porous medium:

where ξ is the shape factor and φ is the porosity (Fu et al., 2020).

2.3.2. Advection

The advection flux of pollutants is

whereJvis the advection flux.

2.3.3. Ionic sorption

According to the linear equilibrium sorption (Smith, 2000), the sorption capacity of a particle is directly proportional to the concentration of pollutants dissolved in the solution:

wherekdis sorption coefficient of the solid phase.

The linear equilibrium isotherm adsorption capacity in porous media is

2.3.4. Governing equation of pollutant migration

Combining Eqs.(24)-(26),the migration process of pollutants is described by the mass conservation equation:

2.4. Hydraulic conductivity equation

Porosity ratio is defined as a function of the effective stress(Cui et al., 2019; Wang et al., 2021):

where φ and φ0are the porosity and the initial porosity, respectively; and εeand εplare the elastic strain and plastic strain respectively.

Hydraulic conductivity is given by

whereKhis the hydraulic conductivity,kis the intrinsic permeability, and μ is the coefficient of dynamic viscosity. The intrinsic permeability model considering tortuosity is given by Dullien(1979) and Fjar et al. (2008):

wheredgis the grain diameter(assuming spherical grains),and κ0is a modifying factor introduced to account for realistic pore shapes.

Combining Eq.(30)into Eq.(29),the hydraulic conductivity can be rewritten as

2.5. Coupling relations

The coupling of the liquid and solid phases involves soil consolidation, solution flow and solute migration. The soil consolidation results from external loading change, during which the change of the void volume among the solid particles can alter the pore water pressure. Accordingly, the bulk flow of the solution is changed. It can alter solute concentration gradient. The coupling relationships are shown in Fig. 2.

3. Model verification

3.1. Model establishment

In order to verify the proposed numerical model,the associated modeling results are compared with the experimental data obtained from Lee and Fox (2009). The geometry of the model is a cylinder with 51 mm in radius and 71.1 mm in height. The cylindrical clay consists of two layers: the upper layer of 52.1 mm and the bottom layer contaminated with KBr of 19 mm.Considering the symmetrical geometry and external loading of the model, the 3D model can be simplified into a 2D axis-symmetric model in COMSOL MULTI-PHYSICS, as shown in Fig. 3. The 2D model adopts quadrilateral mesh including 7272 elements. It is widely accepted that the soil layer is a semi-infinite body, in which the lateral deformation might be neglected, compared with vertical settlement induced by self-weight stress or additional stress.In this case,both the boundary conditions of the side surface and the bottom layer of the 2D model are roller constrained,while the upper layer is confined by a constant force. Relevant parameters of the proposed model are listed in Table 1.

Fig. 2. Schematic defining cross-couplings among physical processes.

Fig.3. (a)Schematic diagram of experimental facility of consolidation-drainage(Lee et al.,2009);and(b)Simulation model for the consolidation-drainage of a composite kaolinite specimen.

3.2. Model validation

The modeling results, such as pore pressure, settlement thicknesses, concentrations of effluent solutes and elastic modulus of the kaolinite layer,are plotted in Fig.4.It shows that the numerical solutions are in good agreement with those obtained by the laboratory measurements (Lee et al., 2009; Lee and Fox, 2009), which helps validate the proposed numerical model.Fig.4a represents the pore pressure at the base of the experimental apparatus and shows a pulsed variation for each loading increment.This is because each loading can consolidate the kaolinite layer, consequently generating a transient excess pore pressure. Since the top boundary of the kaolinite is drained, the excess pore pressure can significantly discharge the pore water. This on-going drainage process can attenuate the excess pore pressure and it gradually returns to the initial condition. The results also show that a larger loading increment can decrease the time required to dissipate the excess pore pressure. The results also indicate that the outflow velocity of the pore water increases first,and then decreases as the consolidation of the kaolinite is gradually completed.The consolidation process is consistent with the increase in the elastic modulus observed in Fig. 4c. It appears that each loading increment can stiffen the kaolinite layer subjected to 1D consolidation. Consequently, the void volume is further diminished, and the pore water outflow is retarded.

After the last load increment is completed,the concentrations of K+and Br-in the pore water are shown in Fig.5a.It is clear that the effluent K+concentration are significantly less than Br-concentration due to lower initial values and sorption in the uncontaminated kaolinite layers.Fig.5b shows the concentration distribution of adsorbed K+along the kaolinite layer. The sorption concentrations of K+in the contaminated layer(at the base)are significantly higher than that in the top layer.These modeling results are also in good agreement with the laboratory data. This indicates that the assumption of linear equilibrium sorption (Eq. (25)) is valid when the increase in the K+concentration is primarily due to the consolidation of the clay.

4. Discussion

In the numerical model,the initial concentrations of Br-and K+are 1672 mg/L and 259 mg/L,respectively(Lee et al.,2009), which makes it challenging to evaluate the consolidation and the sorption effect on solute transport.In the following numerical experiments,the same initial concentrations for K+and Br-are used (both concentrations of 10 mol/m3),while the other parameters in Table 1 are unchanged. The Br-concentration variation can represent advection, diffusion and dispersion induced by consolidation,whilst the variation of K+concentration is used to evaluate the sorption effect on solute concentration.

Table 1 Parameters of experimental and numerical modeling.

4.1. Consolidation effects on solute transport

Fig. 6 compares the kaolinite porosity and effluent K+and Brconcentrations with the numerical results of the incompressible media. In the consolidation case, the porosity profile is a stepwise decline with the application of the loading increments, while the porosity is constant for the incompressible cases. Both the breakthrough curves for K+and Br-in the consolidation case are higher than those of the incompressible case.The results demonstrate that the consolidation effect can accelerate the solute transport.Another surprising observation is that the effluent K+concentration is less than the effluent Br-concentration,which may possibly be the result of sorption effects.

Fig. 7 compares temporospatial variations of Br-concentration under the consolidation and incompressible conditions.From 0.7 d to 3.8 d, the Br-concentration profiles in the consolidation condition are generally consistent with those in the incompressible condition,possibly because the void ratio of the kaolinite layer has little change during the first loading. Nevertheless, during the subsequent loading steps, the consolidation-induced Br-concentrations gradually deviate from those in the incompressible case,as shown in Fig. 6. This could be attributed to consolidation-induced bulk drainage of the pore water, which can accelerate Br-migration from the contaminated layer into the uncontaminated layer.

The K+concentration profiles for the liquid and solid phases are shown in Fig.8a and b,respectively.It is assumed that the sorption coefficient of the incompressible media is 0.003, and the sorption coefficient of the compressible media ranges from 0.003 to 0.01.The results show that the profiles of K+concentrations for the incompressible condition are similar to those in the compressible cases. This phenomenon is different from the consolidationinduced K+concentration enhancement. It seems that sorptioninduced K+concentration decline can predominate over the consolidation effect.

In order to further investigate the sorption impact on K+concentration, the sorption coefficient of the clay particle was varied while all other parameters in the consolidation condition model kept constant.With the sorption coefficient being set to 0.01,the K+concentrations of both the solid and liquid phases near the top of the column are significantly lower than those at the bottom.As the sorption coefficients decrease, both concentration profiles become progressively more uniform with the application of each loading increment.

4.2. Solute transport mechanisms

Fig.9 shows the flow patterns of the two ions breaking through the top layer of the kaolinite.Each loading increment can generate a pulsed flow of the two ions. For the first loading increment, the effluent solution flux increases rapidly and then gradually declines.However,the advection fluxes of both the K+and Br-remain zero at this stage (Fig. 9b). A possible explanation of this effect may be that the excess pore pressure has been dissipated before the two ions are migrated to the top of the kaolinite layer.This explanation also matches the diffusion flux of the effluent ion. As shown in Fig.9c,only K+and Br-begin to discharge from the kaolinite layer after the first loading increment run for 1 d. The ionic diffusion fluxes reach peak at approaching 3 d(4 mol/m2/s for K+and 5 mol/m2/s for Br-). It suggests that the distribution of the ionic concentrations is always variable, even after dissipation of the excess pore pressure.

Fig. 4. Comparisons of the modeling results and measured results (a) Pore pressure;(b) Soil settlement and effluent K+ and Br-concentrations; and (c) Elastic modulus.

Fig. 5. Distributions of ion concentrations in the kaolinite body. (a) Free-state K+ and Br- concentration distributions; and (b) Adsorbed-state K+ concentration.

Fig. 6. Comparison of the soil deformation/solute transport for consolidation and incompressible cases. (a) Porosity of the kaolinite body; and (b) Breakthrough concentrations of K+ and Br-.

Interestingly,the K+mechanical dispersion flux still remains at zero during the first loading process, as shown in Fig. 9d. This phenomenon can be explained by the fact that K+advection has been greatly diminished as a result of the excess pore pressure dissipation.The results suggest that the K+transport manner at the beginning of the consolidation process is dominated by ionic diffusion. The subsequent loading events also trigger a series of pulsed advection and diffusion/dispersion processes.Fig.9b clearly shows that the advection amplitudes of K+and Br-flux are related to incremental loading. This correlation is interesting because the ever-increasing consolidation degree can accelerate the drainage process. Another important effect is that once the consolidation process is generated, the ionic diffusion will decay spontaneously with elapsed time due to the solute concentration distribution becoming homogeneous over the kaolinite. In reality, advection occurs simultaneously with mechanical dispersion. Similar to the ionic diffusion variations, both mechanical dispersion fluxes of K+and Br-also decline with the follow-up loadings. However, the reduction of the mechanical dispersion flux of Br-is obviously greater than that of K+.The different dispersion behaviors between K+and Br-result from K+being easier to adsorb on the clay particles. Free K+near the particles can be readily absorbed by the particles,which can increase the K+concentration gradient around the particles and, as a result, activate K+dispersion.

Fig. 7. Comparisons of Br- concentrations in the kaolinite for the consolidation and incompressible conditions (C0,Br- = 10 mol/m3 =799 mg/L, kp = 3 ×10-3 m3/kg).

In order to investigate the spatial distribution of K+associated with ion sorption, Fig. 8b plots free versus adsorbed K+distributions in the kaolinite layer with different sorption coefficients.Increasing the sorption coefficient reduces the free K+concentration in the slurry and increases the adsorbed K+concentration on the particles at the same depth.The higher the sorption coefficient of the particle is, the lower the K+concentration at the top of the kaolinite layer is.

Fig. 8. Comparison of distribution of K+ concentration for compressible case and incompressible case. Squares are concentration data using incompressible case, liners are concentration data using compressible media: (a) Free-state K+ and Br- concentration; and (b) Adsorbed-state K+ concentration (C0,K+ = 10 mol/m3 = 391 mg/L).

In order to quantify the solute mass of various flow patterns induced by the consolidation effect, the consolidation of the kaolinite was modeled by assuming different compression coefficients(ranging from 0.27 to 0.4).Comparison of the advection flux and the hydrodynamic diffusion flux is evaluated by

whereEis the proportion of the advection mass flux,cais ionic advection mass flux, andcnis total mass flux.

According to the proportion of advection flux from Fig.10, the stepwise loading process can promote the advection and retard the hydrodynamic diffusion. In addition, the solute advection flux is positively related to the consolidation degree.This finding suggests that intermittent pile-loading can gradually homogenize the solute concentration distribution, and during that time the solute transport mechanism undergoes a transition from concentration gradient-controlled diffusion to excess pore pressure-controlled advection.

5. Conclusions

A consolidation-induced contaminant transport model has been developed in this study. This model accounts for (i) multiple physical processes involving soil consolidation, pore water flow,ionic sorption, solute advection and diffusion; (ii) the soil consolidation governed by the modified Cam-clay model, which can describe elastic and plastic deformations of the compacted clay layers; (iii) the pore water flow controlled by the pore pressure evolution induced by the consolidation effect, while the solute transport is influenced by solute concentration distribution; and(iv)The soil permeability controlled by voids between clay particles and by tortuosity of clay.

The results reveal that waste emplacement operations can not only consolidate the waste layer,but also influence the drainage of the slurry and contaminants.For a saturated medium subjected to 1D loading process, the consolidation can evoke solute advection and alter its hydrodynamic dispersion. The advection process largely dependents on the excess pore pressure. Although the excess pore pressure dissipation inhibits solution drainage, the solute distribution is expanded. The long-term effects of the external stress/loading can inhibit the solute diffusion, which is controlled by the concentration gradient, and promote the pore pressure-controlled advection.

The sorption capacity of the kaolinite particles also affects the mechanical dispersion of K+. As the adsorption capacity increases,the solute concentration gradient adjacent to the solid particles increases simultaneously due to the increase in the number of free ions in proximity to the particles. This also has an effect on the acceleration of the mechanical dispersion of the ions.

Fig.9. Effect of episodic loading on the flow patterns of the effluent K+and Br-at the top boundary of the kaolinite: (a) Loading, (b) Advection, (c) Molecular diffusion, (d)Mechanical dispersion, and (e) Total flux.

Fig.10. Proportion of advection mass flux to total mass flux with increasing vertical stress.

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgments

We would like to thank the editor and anonymous reviewers for giving valuable improvements to the paper.The financial support of the National Natural Science Foundation of China (Grant No.41772154), Natural Science Foundation of Shandong Province(ZR2019MA009), Science and Technology Project of Qingdao West Coast New Area(2019-47) is gratefully acknowledged.

List of symbols

φ Porosity of the porous media

ρfLiquid density

υrRelative specific discharge

υ Pore fluid seepage

υsVelocity of solid phase

υdDarcy velocity

KhHydraulic conductivity

γwBulk density of water

PfPore water pressure

υcChemiosmosis velocity

PcChemical osmotic pressure caused by concentration difference

RMohr coefficient

TTemperature

χ Chemical permeability coefficient

αnChange of porous media skeleton with fluid pressure

βPVolume compressibility of a fluid

P,qStress invariants

I1First stress invariant

J2Second deviatoric stress invariant

MSlope of the critical state line in P-Q space determined by internal friction angle

pcThe consolidation pressure

λ Compression index

κ Swelling index

pc0Initial consolidation pressure

erefThe reference void ratio

εp.volVolumetric plastic strain

σij,iNormal stress

fiVolume force in i direction

ρf0Initial densities of the liquid

ρs0Initial densities of solid phases

mfcIon mass of liquid phase

mscIon mass of solid phase

FConcentration of adsorbed contaminants in per unit mass of the solid phase

cContaminant concentration

kdThe adsorption coefficient

prefReference stress

αBBiot coefficient

pmaccEffective mean stress

GShear modulus

KtTangent bulk modulus

εijTotal strain

εplPlastic strain

εcStrain caused by concentration change

gcMaterial properties

aCoefficient of chemical expansion

JdHydrodynamic dispersion flux

DCoefficient of hydrodynamic dispersion

DcEffective diffusion coefficient

DLCoefficient of mechanical dispersion

τ Tortuosity of the porous medium

ξ Shape factor

JvAdvection flux

εeElastic strain

kPermeability

dgGrain diameter

κ0A modifying factor introduced to account for realistic pore shapes