Numerical Analysis of Concrete Gravity Dam Seepage Characteristics Evolution considering theCalciumLeachingEffect

During the long-term service life of hydraulic structures, the calcium compounds in cement-based materials decompose in the aqueous environment, leading to the continuous change of seepage characteristics. To study the influence of calcium leaching on the concrete dam seepage characteristics, we proposed a new mathematic model of the cement-based material calcium leaching model under advection-diffusion-driven leaching. A solid-liquid nonequilibrium model is adopted to model the decomposition of calcium hydroxide (CH) and calcium silicate hydrate gel (C-S-H). To calculate the porosity more accurately, the proposed model takes the effect of different calcium compound decomposition on the porosity increase in consideration, respectively. Shimantan dam is selected for the three-dimensional (3D) calcium leaching analyses. .e 3D finite element model of this dam is analyzed using COMSOL Multiphysics software that is based on the finite element method. Based on the proposed model, seepage characteristics evolutions of the Shimantan dam are studied. Good agreement between the numerical results and the monitored data indicates the accuracy of this simulation..e result shows that after 100 a leaching duration, the uplift pressure increases by 40.8%, and the leakage quantities of the dam body and foundation increase by 48 and 17 times. .e rise of uplift pressure and leakage changes caused by curtain deterioration are the main influences of calcium leaching on the dam seepage..e parameter sensitivity results show that it is necessary to reduce CH content in cement-based materials to obtain better calcium leaching durability. .is model and simulation results can guide the operation of concrete dams under advection-diffusion-driven leaching.


Introduction
Roller-compacted concrete gravity dams have been widely built worldwide due to their shorter construction duration, lower cost, and convenience. Seepage is one of the critical problems in dam operation. Seepage characteristics are closely related to concrete gravity dam safety. Uplift pressure and exit gradient values are the leading causes of different dam failures [1]. Accurate seepage analysis of hydraulic structures is a challenging task [2]. Over their long terms of service, the decomposition of calcium compounds such as that in the form of calcium hydroxide (CH) and calcium silicate hydrate gel (C-S-H) is evitable since the pH of the aqueous environment (pH � 7) is much lower than that of the concrete pore solution (pH � 12.5-13 [3,4]). Under the effect of the concentration gradient, calcium ions leach into the aqueous environment, causing increasing porosity and permeability and degrade material nature. As a result, the seepage characteristic of concrete dams will constantly change in the leaching process.
For example, Fengman concrete gravity dam was built in the year 1937. e impoundment started in the year 1942. After nearly 80 years of operation, Fengman concrete gravity was demolished and rebuilt due to severe leakage and concrete deterioration problems. e aggregate and mortar in the dam concrete are separated, and binding contents almost leached out. Honeycombs and cavities are common in all dam sections. Calcium leaching is the foremost cause of leakage and concrete deterioration of Fengman dam [5]. Gutianxi flat slab buttress dam started to impound in the year 1961. In 1990, there were four dam sections with serious leakage and seven dam sections with calcium leaching. In 2000, there were eight dam sections with serious leakage, twenty dam sections with calcium leaching. e overall strength of concrete slabs is decreased from 49.6 MPa to 37.9 MPa by 23.6% [6].
Calcium leaching has an obvious effect on the diffusivity [7,8], permeability [3,4], and strength [9][10][11][12][13][14] evolution of cement-based materials. e development of the calcium leaching model could guide the operation of hydraulic structures. For modelling research of the calcium leaching process, solid-liquid equilibrium equation has been widely used [15][16][17][18][19]. is equation is a phenomenological chemical equilibrium model relating the calcium concentration of the skeleton and the pore solution on the concept of isotropic damage mechanics [20]. Kamali et al. [21,22] proposed a simplified model to predict the leaching depth considering the initial porosity, water-cement ratio, temperature, and solution composition. For leached sample mechanical properties, Kuhl [20,23] proposed the changing mechanical and transport properties are related to the total porosity and calculated e porosity with average molar volume. Since then, scholars have proposed series of chemomechanical models considering damage mechanics behaviour [24][25][26][27][28][29], the multiscale effect [17,30,31], creep, and damage behaviour [32]. ese models assume a local thermodynamic equilibrium between calcium compounds in the solid skeleton and calcium ions in the pore solution. e leaching process is driven by diffusion. We can calculate the calcium compound decomposition rate with the solid-liquid equilibrium curve. For advection-diffusion-driven leaching, the thermodynamic equilibrium state is broken. erefore, we should not calculate the leaching rate from the solid-liquid curve. Ulm [33] proposed a relation describing the calcium compounds decomposition rate for advection-diffusion-driven leaching in cement-based materials based on chemoporoplasticity theory. When there is a sudden change in capillary Ca 2+ concentration, the decomposition rate is dependent on the actual "distance" from the equilibrium curve [34]. Gawin [34][35][36] adopted Ulm's rate equation and proposed a hydrochemomechanical model of cementitious materials exposed to contact with deionized water. In the calculation of porosity, Gawin adopted Kuhl's calculation method [20]. e porosity in the leaching process is a function of the average molar volume of different calcium compounds and the total calcium content.
In the existing calcium leaching simulation, the key parameters such as diffusivity and permeability are functions of porosity. However, the porosity is generally calculated with average molar volume and ignored the effect of different calcium compounds decomposition on the porosity evolution. Decomposition of CH generates capillary pores, while the decomposition of C-S-H generates gel pores and capillary pores. erefore, to model the leaching process more accurately, especially in the advection-diffusion-driven leaching, we should consider the effect of different calcium compounds decomposition. e long-term service performance of dams has always been one of the most important research fields in hydraulic engineering [37]. e main focuses of these studies are the evolution of dam deformation and stress. However, as a water storage structure, the long-term seepage characteristic evolution of dams under the effect of calcium leaching has not been investigated in the literature. Seepage characteristics are closely related to dam safety. e evolutions of dam seepage characteristics under the effect of advection-diffusion-driven leaching are closely associated with long-term service safety and benefits and need further exploration. e objective of this study is to investigate the concrete dam seepage evolution in the leaching process. In the existing advection-diffusion-driven leaching simulation, the porosity is generally calculated with average molar volume and ignored the effect of different calcium compounds decomposition on the porosity evolution. To improve the simulation accuracy, we proposed a new mathematic model of hydrochemo behaviour of cement-based materials. e decomposition of CH and C-S-H are considered, respectively. Kozeny-Carman (KC) relation is adopted to model the permeability coefficient evolution. Shimantan dam is selected for the three-dimensional (3D) calcium leaching analyses.
e 3D finite element model of this dam is modelled using COMSOL Multiphysics software. Based on the proposed model, seepage characteristics evolutions of Shimantan dam are obtained. e accuracy of this simulation is verified by the water head and leakage monitoring data. e influences of initial CH content and hydraulic conductivity are also discussed. e proposed numerical model can be used to analyze the antiseepage ability deterioration of grout curtains, concrete face slabs, and core walls under the effect of advection-diffusion-driven leaching.
ese simulation results will guide the dam safety monitoring and evaluation. e article is organized as follows. First, the basic theory of the model is presented. en, a numerical example is solved to validate the proposed model. Next, the seepage characteristic evolution of Shimantan concrete gravity dam is analyzed with the presented model. Behind, the influences of CH content and hydraulic conductivity are discussed. Finally, the conclusions are presented.

Materials and Methods
e nature of the calcium leaching process of cement-based materials in hydraulic engineering is calcium compounds decomposition in low pH solutions and then precipitation into the surrounding environment under the concentration gradient and hydraulic gradient. e main chemical reactions equations are given in the following equation [20]: e pore solution of cement-based materials is assumed to be in thermodynamic equilibrium with the calcium compound in the solid skeleton. When the calcium ion concentration of pore solution decreased below 22 mol/m 3 , CH starts to decompose. After the complete decomposition 2 Advances in Civil Engineering of CH and the calcium ion concentration decreased to 18 mol/m 3 , C-S-H starts to decay. ere are two modes of calcium ion transport in hydraulic structures: advection and diffusion. Advection is the convection of calcium ions in the pore solution that flows through porous media. e driven force of advection is the hydraulic gradient. Diffusion is the movement of calcium ions in the pore solution from high ion concentration to low ion concentration. e driven force of diffusion is the concentration gradient.
From the concept of thermodynamic, there are two states between calcium ions in the pore solution and calcium compounds in the solid skeleton: equilibrium state and nonequilibrium state. In the equilibrium state, the decrease in calcium ion concentration induces a thermodynamic imbalance, which drives the decalcification of the solid skeleton until thermodynamic equilibrium is reached [33]. In other words, the values of calcium content in the skeleton are determined by calcium ions concentration in the pore solutions. e driven force of calcium leaching in the equilibrium state is diffusion. In the nonequilibrium state, the decalcification of the solid skeleton will not lead to equilibrium. ere is no relation between calcium ion concentration and calcium content in the solid skeleton. e driven force of calcium leaching in the nonequilibrium state is advection. Figure 1 shows the general sketch of the physical model for calcium leaching in concrete dams. e calcium ion concentration in the water environment is much lower than in the pore solution. Under the concentration gradient and hydraulic gradient effect, the calcium ion, which comes from CH and C-S-H decomposition, participates into the surrounding water environment. As a result, the porosity and transport properties of dam materials increase. e dam seepage characteristics constantly change in the process of calcium leaching. Figure 2 shows the calcium leaching phenomenon in corridors of the Shimantan concrete gravity dam overflow sections after 14 years of operation. As we can see from the pictures, the white leached matter cascades down, covering almost the whole inner surface of the corridors. e white leached matter is calcium carbonate, which comes from the carbonization of leached CH. Figure 3 shows all the phenomena of the calcium leaching process under equilibrium (Figure 3(a)) and nonequilibrium (Figure 3(b)) conditions. Equation (1) shows the decomposition reaction equations of CH and C-S-H. In diffusiondriven leaching, the calcium ion concentration in the pore solution is decreasing gradually. When the calcium concentration of the pore solution falls below a critical level, Ca(OH) 2 crystals start to decompose. After all the CH leaches out, C-S-H begins to deteriorate [38,39]. While in advectiondiffusion-driven leaching, calcium ion concentration is lower than the critical concentration of CH and C-S-H. erefore, CH and C-S-H decompose simultaneously. According to the assumption that CH leaching generates capillary pores and C-S-H leaching causes gel pores and capillary pores [3], the diffusivity and permeability increase continuously. e increased permeability coefficient and diffusion coefficient accelerate the leaching process.

Governing Equations.
In the advection-diffusion-driven leaching model, we apply the following six assumptions: (1) Silicon does not leach at all, and the rehydration and leaching of unhydrated cement particles are ignored either, so only CH and C-S-H are considered in the leaching process [39]; (2) CH leaching generates capillary pores, and C-S-H leaching not only causes gel pores but also contributes to capillary pores [3]; (3) Calcium ions in solution are not reacting to form any new compounds [38] (4) Calcium ions in solution are not reacting to form any new compounds (like CaCO 3 , for instance) [38] (5) e materials are and remain saturated over time, and isothermal conditions are preserved [38]; (6) e flow in the capillary pore is laminar, which can be described by Darcy's law e diffusion equation is suitable for leaching in bridge piers, pile foundations, and nuclear waste containers. For calcium leaching in concrete dam face slabs, concrete core walls, and grout curtains, we should add the u∇c to governing equation to express the convection effect of advection flow on calcium ions. u is the Darcy velocity, c is the calcium ion concentration in the pore solution, and ∇c is the concentration gradient. u∇c is the intensity of advection flow, and the physical meaning of u∇c is the quantity of calcium ion pass-through per unit area per unit time under advection flow. e advection flow rate is not invariant; with increasing porosity and permeability, the flow rate changes with redistribution of the water pressure. erefore, the seepage continuity equation is necessary for expressing the evolution of the seepage field. Equation (2) presents the governing equation describing the fluid flow [40,41] and calcium ion dispersion [42] in the permeable porous medium.
In this equation, u is the Darcy velocity (m/s), k is the permeability of the cement-based materials (m 2 ), ρ is the water density (kg/m 3 ), g is the gravitational acceleration (m/ s 2 ), P is the water pressure (Pa), ε p is the porosity, t is the time (s), Q m is a mass source term (kg/(m 3 ·s)), c is the concentration of the species (mol/m 3 ) (in this study, c is the concentration of calcium ions in the pore solution), D denotes the diffusion coefficient (m 2 /s), and R is the calcium compounds decomposition reaction rate (mol/(m 3 ·s)).
In this study, retardation factor R refers to the decomposition rate of calcium compounds in the solid skeleton. e decomposition of calcium compound is assumed    to be an irreversible process under advection-diffusiondriven leaching. In other words, the calcium ion concentration of the solid skeleton will not increase. In practice, calcium ion in the pore solution is nonconservative and can be adsorbed onto the soil particles. e adsorption and chemical reactions that might occur in the soils or rocks are complex. According to Gerard's work [38], calcium ions in solution assume not reacting to other compounds. In this study, considering the enormous scope of the model, the effect of adsorption and chemical reaction in the soils or rocks is assumed to be limited. erefore, to simplify the numerical model, calcium ions in solution are deemed not to react to form any new compounds, like CaCO 3 [38]. e calcium leaching of cement-based material contains two processes. e first stage is the decomposition of calcium compounds, and the second one is the migration of calcium ions in the porous media. Because of the calcium compounds decomposition, the hydraulic conductivity keeps changing in the leaching process. As a result, the flow rate changed. Compared with the existing diffusion-driven leaching, another assumption is needed. Since calcium leaching is a relatively slow process, the water flow in the porous media is still Darcy flow. In this way, Darcy's law determined the flow rate.

Calcium Compounds Decomposition Model.
In the diffusion-driven calcium leaching process, the time required to dissolve CH and C-S-H is shorter than the diffusion time [38]. As a result, there is an equilibrium state between calcium compounds and calcium ions. Buil and Berner first proposed the solid-liquid equilibrium curve in 1992 [43]. Equation (3) is a typical solid-liquid equilibrium curve modified by Nakarai [15].
where s Ca is the actual calcium concentration in the solid skeleton (mol/m 3 ), c is the concentration of calcium ions in the pore solution (mol/m 3 ), x 1 is the calcium ion concentration in the pore solution when sharp decomposition of C-S-H starts (mol/m 3 ), x 2 is the calcium ion concentration in the pore solution when CH is completely dissolved and C-S-H starts to decompose (mol/m 3 ), c satu is the saturation concentration in deionized water at room temperature (mol/ m 3 ), C CSH0 and C CH0 are the nonleached C-S-H gel and CH content, respectively, which can be calculated from the cement chemical composition (mol/m 3 ), x 1 is taken as 2 mol/m 3 by Gerard [44], Wan [39], Nakarai [15], and Jain [45], x 2 is taken as (c satu − 3) mol/m 3 by Wan [39], Nakarai [15], and Phung [3], or (c satu − 1) mol/m 3 by Gerard [44], Jain [45], and c satu is taken as 22 mol/m 3 by Phung [3], Jain [45], and Wan [43] or 20 mol/m 3 by Gerard [45] and Kuhl [20], or 20.5 mol/m 3 by Gawin [34]. In this study, the value of x 1 , x 2 , and c satu are taken as 2, 19, and 22 mol/m 3 . According to Allen's test [18] and Jennings' work [19], the initial Ca/Si ratio of C-S-H is 1.70, and the Ca/Si ratio of partially leached C-S-H is 0.8. e solid-liquid equilibrium curve provides a relationship between calcium ion concentration in pore solution and calcium compounds in the solid skeleton, which is suitable for the locally diffusion-controlled reaction. One of the most critical assumptions of the solid-liquid equilibrium curve is the time required to dissolve the calcium compounds in the solid skeleton shorter than the diffusion time [33,38]. However, in advection-diffusion-driven leaching, the calcium ion concentration in the pore solution may be lower than the equilibrium concentration. e calcium ion in the pore solution is no longer in an equilibrium state with calcium compounds in the solid skeleton. Under the migration effect of advection flow, calcium compounds leaching time in the solid skeleton is no longer faster than the diffusion time. erefore, the solid-liquid equilibrium curve is not suitable for advection-diffusion-driven leaching.
Ulm et al. [33] proposed a nonequilibrium decomposition model of calcium compounds using chemoporoplasticity theory. Gawin [36] neglected the elastic deformation and plastic hardening-softening phenomena term and the decomposition rate presented in the following equation: where s Ca is the actual calcium concentration in the solid skeleton (mol/m 3 ), t is the time (s), η is the microdiffusion of Ca 2+ ions in different compounds, which is dependent on the calcium content (mol/(J·s)), A s is the chemical affinity, being a driving force of chemical reactions (J/m 3 ), R is the gas constant (J/(mol·K)), T is the temperature (K), τ leach is the characteristic time of leaching (s), c Ca is the present calcium concentration in the pore solution (mol/m 3 ), (c eq Ca , s eq Ca ) is a point on the solid-liquid equilibrium curve, and κ(s Ca ) is the equilibrium constant.
A s in equation (4) is not convenient to be used in simulations because it requires knowledge of the pair of values (c eq Ca , s eq Ca ). Gawin [36] transferred this equation into a form which is more suitable for practical applications as in the following equation:

Advances in Civil Engineering
where s 0 Ca is the initial calcium concentration in the solid skeleton (mol/m 3 ), s Ca is the actual calcium concentration in the solid skeleton (mol/m 3 ), s * Ca is the equilibrium value corresponding to current calcium concentration in the pore solution (mol/m 3 ), and s * Ca can be calculated by equation (3). Figure 4(a) shows the solid-liquid equilibrium curve with different CH content and equilibrium e CH content is defined as C CH /(C CH + C CSH ). e total calcium compound in the solid skeleton is 15,000 mol/m 3 . e solid-liquid equilibrium curve is calculated with equation (3). As we can see from Figure 4(a), the CH decomposition stage becomes shorter as CH content decreases. Figure 4(b) shows S Ca S * Ca κ(s)ds as a function of the calcium content, and the temperature T is 25.
Taking the actual "distance" from the equilibrium curve as the measurement of calcium compounds decomposition rate is suitable for advection-diffusion-driven leaching. In advection-diffusion-driven leaching, advection flow carries the calcium ions in the pore solution out. As a result, the calcium ion concentration in the pore solution is always lower than in the equilibrium stage. erefore, we should not calculate the calcium compounds in the solid skeleton from calcium ion concentration in the pore solution with equation (4). Under this circumstance, we can calculate the calcium compound content by integrating the decomposition rate with time.
e driven force of calcium leaching determines the model used in the simulation. For diffusion-driven leaching, we should use the solid-liquid equilibrium model. While for advection-diffusion-driven leaching, we should use the nonequilibrium model. e dimensionless Péclet number (Pe number) represents the relative proportion of advection and diffusion. Pe number is defined as in the following equation [46].
where u is the mean flow velocity in the porous media (m/s), R is the characteristic length, for cement-based material, and R is the average diameter of particles [47]. According to the code for mix design of hydraulic concrete (DL/T 5330-2015), the average diameters of two-grade roller-compacted concrete, three-grade compacted concrete, and regular concrete are 0.02 m, 0.04 m, and 0.02 m. is study takes the average diameter of grout curtain as 0.01 m. D e is the molecular diffusion, m 2 /s. When Pe number is less than 1, diffusion is dominant. When Pe number is more than 100, the effect of diffusion can be ignored [48].

Porosity Evolution.
e main calcium compound content in cement-based materials includes CH and C-S-H. Wan [39] assumed that CH leaching generates capillary pores and that C-S-H leaching generates gel pores. In contrast, Phung [3] proposed that decalcification of C-S-H may also contribute to the increase in capillary porosity due to the spatial connectivity of gel pores. Kuhl [20] proposed a simplified model for calculating porosity with an average molar volume of CH, ettringite, and C-S-H. By integrating the decomposition rate of calcium compounds, we can obtain the total molar amount of dissolved calcium compounds and porosity. e method was adopted by Gawin [36]. e main disadvantage of this method is that it ignores the effect of different calcium compound decomposition on pore growth.
To more accurately calculate the evolution of porosity and permeability coefficient of materials, we should consider the influence of the decomposition of CH and C-S-H on the porosity. e porosity evolution of cement-based materials is shown in the following equations.
where φ is the total porosity, φ 0 is the initial porosity, Δφ CH is the porosity increase due to CH decomposition, Δφ CSH is the porosity increase due to C-S-H decomposition, φ c is the capillary porosity, φ c0 is the initial capillary porosity, and c is the extent of contribution of C-S-H decalcification to capillary porosity; in this study, c � 0.5 [3], V CH and V CSH are the molar volume of CH and C-S-H, V CH � 0.033 × 10 − 6 mol/m 3 , V CSH � 0.04 × 10 − 3 mol/m 3 , η CH and η CSH are the coefficients related to the microdiffusion of Ca 2+ ions in CH and C-S-H, and the values of η CH and η CSH are given in Table 1.
In advection-diffusion-driven leaching, when the calcium ion concentration in the pore solution is between 19 and 22 mol/m 3 , CH decomposes. When the ion concentration is below 19 mol/m 3 , CH and C-S-H decompose simultaneously, combined with equations (4) and (5), and the porosity increase induced by decomposition of CH and C-S-H can be obtained, respectively. According to Morandeau's experimental data [49], the molar volume of C-S-H is linear with the ratio of Ca/Si ratio. In this study, we also considered the dependence of the molar volume of C-S-H on Ca/Si. erefore, the molar volume of C-S-H in the leaching process is defined as a function of Ca/Si. More details are presented in Appendix A.

Permeability Evolution.
During the leaching of cementbased materials, the intrinsic permeability K may increase significantly due to skeleton decalcification and is related to the leaching-induced porosity increase. Saito [50] proposed an exponential formula between porosity and permeability, and the relation presents in equation (8).
6 Advances in Civil Engineering where K(n 0 ) is the initial permeability (m/s), A T is a material parameter that should be determined experimentally, and n is the concrete porosity. is equation works well in accelerated electrochemical tests on mortars. However, to better understand the leaching effect on the increase in cement-based material permeability, it is necessary to explore the evolution of microstructure parameters. e Kozeny-Carman (KC) relation works well for relatively homogeneous porous materials [51][52][53]. e permeability coefficient presents in the following equation.
where K is the permeability coefficient (m 2 ), χ is the microstructure parameter (m 2 ), φ is the capillary porosity, τ is the tortuosity, S a is the specific pore surface (m 2 /m 3 ), and F s is the shape factor. e permeability is defined as a function of porosity in equation (8), neglecting the influence of pore structure evolution on permeability.
e variables related to pore structure are included in the parameter A T .
e Kozeny-Carman relation considered porosity and pore structures such as tortuosity, specific pore surface, and shape factor on the permeability evolution. In this way, the modelling of permeability evolution is more accurate. More details of this equation presented in Appendix B.

Diffusivity Evolution.
In the process of calcium leaching, the diffusivity also increases with the porosity [54][55][56][57]. Garboczi and Bentz [8] proposed the effective diffusion model. We can calculate the effective diffusivity of cement paste from cement hydration and the water-cement ratio. Garboczi and Bentz's model holds for a standard nonleached cement matrix. In the leaching process, the diffusivity increases much more rapidly. Based on Snyder and Clifton's work, van Eijk and Brouwers [58] proposed a revised formula as in the following equation.
where φ c (x, 0) is the initial capillary porosity, H is the Heaviside function, φ c (x, t) is the capillary porosity, and D 0 is the diffusivity of Ca ion in water, D 0 � 4.5 × 10 − 10 m 2 /s.  Advances in Civil Engineering 7 Equation (11) has been widely used in calcium leaching simulations, such as Wan [39]. In this study, we adopt this equation to modelling the evolution of the diffusion coefficient in the leaching process.

Overall Scheme of the Advection-Diffusion-Driven
Leaching Model. In advection-diffusion-driven leaching, the decomposition rate of calcium compounds is a function of calcium ion concentration in the pore solution and current calcium content in the solid skeleton. Under the convection effect of advection and diffusion, the distribution of calcium ion concentration varies with distance and time. erefore, the decomposition rates of different dam structures are not the same. We can calculate the quantity of decomposed calcium compounds from the integral decomposition rate with time. Once the decomposed calcium amount is determined, the porosities, hydraulic conductivities, and diffusivities are obtained. Due to the evolution of hydraulic conductivity and diffusivity, the calcium ion concentration of the pore solution redistributes. In this study, the calcium ion concentration, decomposition rate of calcium compounds, porosity, hydraulic conductivity, and diffusivity varied with distance and time. Figure 5 shows the overall framework of this model. First of all, input the geometric model, material properties, and solution conditions. We can obtain the concentration of calcium ion in pore solution by solving the advection-diffusion equation and seepage continuity equation. Combining with the solid-liquid nonequilibrium model, the decomposition rate of the solid phase is available. By integrating the decomposition rate of calcium compounds with time, the amount of leached and remained calcium compounds are obtained. After that, we can calculate the porosity. e permeability and diffusivity are defined as functions of porosity. erefore, the advection flow rate and calcium ion diffusivity coefficient are obtained. At this moment, the seepage characteristic is changed under the effect of calcium leaching. e advection-diffusion equation and seepage continuity equation need to be solved again. e calcium concentration is redistributed under the impact of changed advection flow rate and calcium ion diffusivity coefficient. In this study, we numerically solve the problem using the COMSOL Multiphysics simulation tool. e finite element method has been adopted.

Numerical Example
A numerical example is solved to compare this study and previous studies to validate the proposed model. is example is a fast reaction-advection-diffusion problem solved previously by Gawin et al. [34,35]. It concerns calcium leaching from a cement paste wall (1D problem) exposed to one-side action of deionized water, which is modelled with Dirichlet boundary conditions. e model configuration and boundary conditions are shown in Figure 6. e primary materials' properties are given in Table 2. More details are present in reference [35]. e time history of calcium ion concentration in pore solution and space distributions at different time stations are shown in Figure 7. e simulation results are obtained from the proposed model with and without different calcium compounds decomposition effects. Reference points 1 and 2 are 0.09 cm and 2.88 cm from the upstream side, respectively. e leaching duration is 2 × 10 8 s and 6 × 10 8 s . e boundary conditions and material properties used in this simulation are the same as Gawin's Case 4. e simulation results of this study are consistent with Gawin's results [35], which indicate the rationality of this simulation. However, the time history of calcium ion concentration considering different calcium compound decomposition is larger. is is because the proposed model removed the effect of gel pores produced by C-S-H decomposition on transport properties. Figure 8 shows the time history of hydraulic conductivity and its space distributions at different time stations. As we can see from Figure 8(a), the time history of hydraulic conductivity considering different calcium compounds decomposition is lower. Figure 8(b) shows that within the range of 0.04 m from the upstream side, the hydraulic conductivity without considering different components is one order of magnitude larger than that considered. In the previous model, the pores generated by C-S-H decomposition are all considered capillary pores. In the proposed model, we have taken the effect of different calcium       antiseepage ability of the dam body, two-grade roller-compacted concrete with a 3 m thickness was deployed on the upstream side. On the downstream side, below an altitude of 86.00 m, two-grade roller-compacted concrete with a 1.5 m thickness was deployed. e dam body was made of threegrade roller-compacted concrete. e bed course was made of regular concrete. e geographical location of the Shimantan reservoir is shown in Figure 9. e nonoverflow cross-section of the Shimantan reservoir is shown in Figure 10.
We created a three-dimensional finite element model of Shimantan concrete gravity dam to investigate the seepage characteristic evolution. e model selected the typical nonoverflow dam section. e 3D model of the dam body had six different sections. e concrete face slab is made of two-graded concrete. e maximum thickness of the slab is 4 m, and the minimum thickness is 2 m. e dam body is made of two-graded concrete, and the toothed wall is made of conventional concrete. e thickness of curtain grouting is 2 m, and the depth is 20 m in bedrock. We simplified the cross-section of the drainage pipe in the dam body and dam foundation as square, and the side length is 0.35 m and 0.2 m, respectively. We extended the foundation two times the dam's height at both the upstream and downstream sides of the dam. e foundation's height is taken into account as much as two times the dam's height. We simplified the dam body as a homogeneous porous material to reduce the element number. Totally, 653,914 finite elements are used in the 3D model. e maximum element width is 1.45 m, and the minimum element width is 0.53 m. e finite element meshes of the 3D models are shown in Figure 11.

Calculation Parameters.
Since Shimantan concrete gravity dam has been in operation for more than 20 years, it is not easy to obtain the initial parameters of dam materials. erefore, the parameters can only be obtained through inverse analysis and experience. e initial permeability coefficients of different material zones are obtained by inverse analysis with monitored data, and empirical parameters determine the initial porosity. e initial diffusivity is calculated from equation (7) with initial porosity. e calcium diffusivity in the water is adopted as the diffusivity in the drain pipe. e initial CH and C-S-H contents are 6000 and 9000 mol/m 3 . Bulk density and initial surface area assume Phung's test parameters of S3 [3], the intact and leached bulk density are 30.6/145.8 m 2 /g, and the intact and leached surface areas are 1780/1110 kg/m 3 . e parameters of nonequilibrium solid-liquid dissolution followed Gawin's model. e calculation parameters are given in Tables 1 and 3.

Initial and Boundary Conditions.
e initial conditions and boundary conditions are shown in Figure 12.
For the calcium leaching zones, which include normal concrete, two-grade concrete, three-grade concrete, and grout curtain, the pore solution is assumed to be saturated according to the previous studies [3,39]. e initial calcium ion concentration is assumed to be zero for the nonleaching zones, including the rock layer and drain pipes. e initial calcium concentration is described as in the following equation.
Upstream water pressure boundary applies on the surface of concrete face slab and reservoir basin below 107 m. Downstream water pressure boundary applies on the surface of the dam body and reservoir basin below the elevation of     Figure 9: Geographical location of the Shimantan reservoir.

Advances in Civil Engineering
Calcium ion concentration on the upstream and downstream water pressure boundaries is assumed to be zero. e calcium ion concentration boundary presented in the following equation.

Pe Number.
It is necessary to calculate the Pe number of different dam materials in the initial condition to determine which calcium compounds decomposition model. e Pe numbers in the initial state are given in Table 4.
In the initial condition, the Pe numbers of two-grade concrete, three-grade concrete, normal concrete, and grout curtain are 115. 5, 30.8, 117.0, and 71.7, and greater than 1 indicates that the calcium leaching in Shimantan concrete dam is advection-diffusion-driven leaching. In the calcium leaching process, with the decomposition of calcium compounds, the hydraulic conductivity increases gradually. Since the water level upstream and downstream is unchanged, the flow rate is increasing all the time. Eventhough the diffusion coefficient increases with porosity, judging from equations (10) and (11), the increasing rate of the diffusion coefficient is minor than hydraulic conductivity. erefore, the decomposition model of calcium compounds should adopt the nonequilibrium model. Figure 13 shows the evolution of the overall water head distribution over a hundred years. As we can infer from the figure, when the leaching time is 0 a, the contour lines are around the face slab and the curtain. ere is almost no contour line in the dam body, which indicates the good antiseepage performance of the concrete face slab and the grout curtain. e water head distribution of the dam foundation is between 86 m and 88 m. When the leaching time is 10 a, the overall water head distribution is essentially the same as 0 a. ere is still no contour line inside the dam, indicating that the antiseepage system works well. When the running time is 50 a, the overall water head contour line significantly changes. No contour line is concentrated around the grout curtain, suggesting that the grout curtain has lost antiseepage ability. e water head of the dam foundation rises to 88-90 m. e 86 m and 88 m water head contours are attached to the downstream side of the dam body, which was previously around to the drainage pipe. is finding suggests that the drainage pipe of the dam foundation cannot effectively reduce the uplift pressure of the dam foundation.

Total Water Head Distribution Evolution.
ere is still a dense contour distribution in the face slab, and a sudden drop of the water head appears, which indicates that the concrete face slab still has antiseepage ability after a 50-year leaching duration. When the leaching duration reaches 100 a, the distribution of the water head contour in the concrete face slab is no longer dense. e water head contour arrives inside the dam body, which indicates the concrete face slab has lost antiseepage ability.
In Figure 13, the contour lines 104 and 106 are closer to the upstream surface as shown in Figure 13(a) compared to Figure 13(d), while the contour lines 88 and 90 were closer to the downstream surface in year as shown in Figure 13(d) compared to Figure 13(a). When the leaching duration is 0 a, the grout curtain has excellent antiseepage ability. e contour lines on the upstream side of the grout curtain are denser than those on the downstream side. After 100 years of leaching duration, the grout curtain has lost antiseepage ability. e permeability of dam foundation rock is assumed to be homogeneous. erefore, the contour lines tend to be uniformly distributed in the leaching process. e contour lines 104 and 106 are closer to the upstream surface in year 1, while contour lines 88 and 90 are closer to the downstream surface in year 100. Figure 14 shows the comparison between the simulated water head distribution and monitored data. e upstream water level is 107.0 m, and the downstream water level is 86.0 m; the position of the piezometer tube is shown in Figure 10. For the piezometer tube P 2 9 , the measured value stays at 92 m after 2005 and does not change with the upstream water level. is result may be due to clogging caused by the silt at the bottom of the piezometer tube. e observed value of this point is not representative; as a result, we omitted this value from this analysis.
As we can infer from the figure, the simulated results show good agreement with the monitoring data, indicating the accuracy of this simulation. When the leaching duration is 0 a, the water head decreases to 82.5 m at the position of the drainpipe and 3.5 m lower than the downstream water level, showing a good drainage effect. e antiseepage system reduces the water head by 24.5 m in total. When leached 25 a, the water head of the drainpipe is almost the same as that of the downstream side. is result shows that the drainage pipe can effectively reduce the uplift pressure of the dam foundation. When the leaching duration reaches 50 a, the grout curtain and drainage pipe have essentially lost their seepage prevention ability, and the uplift pressure of the dam foundation has increased significantly. e downstream water level reaches 89.1 m. When the leaching duration comes 100 a, the downstream water level rises to 91.8 m. e water head at the concrete face slab has decreased by 4 m. We can obtain the uplift pressure by integrating the water head along the foundation surface. Compared with the initial time, the uplift pressure increased by 10.9% after 50 a leaching duration. After 100 years of leaching duration, the Advances in Civil Engineering 13 uplift pressure of the Shimantan concrete dam increased by 40.8%.

Leakage Quantity Evolution.
e monitored data shown in Figures 15 and 16 are the average annual leakage from measuring weirs. We calculated the simulated results from the following equation: where q b is the dam body leakage (L/s), v is the flow velocity vertical to the overflow surface (m/s), L 0 is the total length of the dam, 645 m, and L 1 is the width of the numerical model, 6 m. In this study, we assume that leakage quantity is the same in all the dam sections, ignoring the differences between the overflow and nonoverflow dam sections and the influences of the river valley.
As we can see from Figure 15, the leakage quantity of the dam body was 0.884 L/s in 2003 and 1.151 L/s in 2010. e annual leakage of the dam body fluctuated, and overall, there is an increasing trend. From the comparison of numerical results and monitored data, the error of dam foundation leakage is 0.01-0.32 L/s, the error percentage is 3.48-27.46%, and the average error percentage is 0.64%. e simulated results are in good agreement with the monitored data, indicating the accuracy of this simulation. Figure 16 shows the comparion of the simulated dam foundation leakage quantity and monitored data from measuring weirs. e dam foundation leakage is also taken from the average annual amount, and the numerical results are also calculated from equation (15). It is noteworthy that both the monitoring data and the simulated results include the seepage from the right sidewall of the corridor. e leakage mainly comes from the dam foundation; however, there may also be some bypass seepage flow from downstream. As we can infer from the figure, the leakage quantity of the dam foundation was 1.341 L/s in 2003, 0.457 L/s larger than that of the dam body. In the year of 2010, the leakage quantity was 1.400 L/s, 0.249 L/s larger than that of the dam body.
e annual leakage of the dam foundation also fluctuates; overall, there is an increasing trend. From the  comparison of numerical results and monitored data, the error of dam body leakage is 0.01-0.24 L/s, the error percentage is 0.8-21.6%, and the average error percentage is 1.18%. us, the simulated results are in good agreement with the monitored data, indicating the accuracy of this simulation. Figure 17 shows the evolution of the dam body and dam foundation leakage quantity over a century. As we can infer from the figure, in the first ten years, the leakage of dam foundation is more significant than dam body, and then, the leakage of dam body exceeds the dam foundation. e seepage flow goes through different structures and materials for leakage in the dam foundation, including the rock and grout curtains. When the hydraulic conductivity of the curtain is higher than rock, the leakage of the dam foundation should not increase significantly. e water head distribution in Figure 13(c) indicates the hydraulic conductivity of the grout curtain is later than the rock, and the grout curtain has lost antiseepage ability. However, the leakage evolution curve still increases with the leaching duration after 2050. e leakage increase rate becomes faster and faster, probably because of the bypass seepage in the dam body. As we can infer from Figures 13(c) and 13(d), a denser distribution of water head contour in the dam body indicates that the bypass seepage around the dam body is gradually enhanced in the leaching process. For dam body leakage, the seepage flow goes through the concrete face slab. With the increase of the hydraulic conductivity of the slabs, the leakage of the dam body will increase continuously. After 100 years of leaching duration, leakage quantity of dam body and foundation rise 48 and 17 times. Figure 18 shows the evolution of the hydraulic gradient in the concrete face slab central line at 84.3 m. e position of the central line is shown in Figure 10. We calculate the hydraulic gradient from the following equation.

Hydraulic Gradient Evolution.
where i is the hydraulic gradient, and H x , H y , and H z are the gradients of the water head H. e maximum value of hydraulic gradient at the beginning is 11.93, which is less than the allowable value of 150.
ere is a sudden drop in the area 3 m from the upstream side, and the lowest hydraulic gradient is 0.30. is result occurs because of the existence of a dam body drain pipe. With the increase of leaching duration, the hydraulic  gradient of the upstream side decreases, while it rises slightly in the area 1-2 m away from the upstream side. is is due to the increase of hydraulic conductivity in the upstream concrete face slab. When leached 50 a, the upstream side hydraulic gradient decreased to 3.36, while the region behind the drainpipe increased to 2.15.
is result is due to the decomposition of calcium compounds on the upstream side. e hydraulic conductivity of the upstream side is larger than the downstream side. When leached 100 a, the hydraulic gradient of the upstream side decreased to 1.75, and the downstream side increased to 3.41. e increase of hydraulic gradient in the downstream side also reflects bypass seepage in the dam body. Figure 19 shows the evolution of the hydraulic gradient in the grout curtain central line at 61.00 m. e grout curtain central line is shown in Figure 10. e hydraulic gradient is 2.6 at the initial moment. As the leaching duration increases, the hydraulic gradient decreases. When the leaching duration is 100 a, the hydraulic gradient is 0.11, 95.8% lower than that initially. As we can see from the evolution of the hydraulic gradient curves, the decrease of the upstream side is much faster than that of the downstream side. is result occurs because the seepage flow controls the calcium ion migration direction. Calcium ions are carried from the upstream side to the downstream side. e calcium ion concentration on the upstream side is lower than that in the other parts. As a result, the decomposition of calcium compounds is much faster, and the permeability coefficient is more significant.

Sensitivity Analysis of CH Content and Initial Hydraulic
Conductivity. To better understand the calcium leaching process in cement-based materials, we carry out parameter sensitivity analysis of CH content and initial hydraulic conductivity. e total water head distribution shows that the curtain is the most vulnerable structure in the leaching process. erefore, the evolution of CH concentration in the solid skeleton, hydraulic conductivity, hydraulic conductivity increase times of grout curtain with different CH content, and initial hydraulic conductivity are studied.
In the sensitivity analysis of CH content, the total calcium compound is 15,000 mol/m 3 , and the CH content is 0.4, 0.2, 0.05, and 0.02. ese four groups of CH contents are named grout curtain A, B, C, and D. In the sensitivity analysis of initial hydraulic conductivity, the initial hydraulic conductivity is 0.01K, 0.1K, 1.0K, 10K, K � 1.0 × 10 − 8 m/s. e other parameters are the same as the initial simulation. Figure 20 shows the evolution curves of the residual CH in the curtain midline within 100 years. As we can see from the figure, in the whole leaching process, the CH concentration in the upstream side is lower than the downstream side, indicating the upstream side is more vulnerable to calcium leaching. When the leaching time is 50 a, the leached CH amount of these four grout curtains at x � 0 m point is 652, 651, 639, and 300 mol/m 3 , and the decomposition rate of CH in grout curtain A is slightly faster than B and C. For grout curtain D, the CH is totally leached within 0.7 m. When the leaching duration is 100 a, the leached CH amount of these four grout curtains at x � 0 m point is 1287, 1284, 750, and 300 mol/m 3 . CH in grout curtain C is completely leached within 0.4 m, and CH in grout curtain C is totally leached. Figure 21 shows the evolution curves of hydraulic conductivity at the grout curtain central line (elevation 75 m) in 100 a. e initial hydraulic conductivity of all grout curtains is 1.5 × 10 − 8 m/s. As we can see from Figure 21, when the leaching duration is 50 a, the hydraulic conductivity of grout curtains A, B, and C is almost the same. However, when the leaching duration is 100 a, the hydraulic conductivity curves of grout curtains A, B, and C are different. Within 0.4 m to the upstream side, the hydraulic conductivity of grout curtains A and B varied from 1.0 × 10 − 6 m/s to 3.1 × 10 − 7 m/s, and for grout curtain C, the hydraulic conductivity is from 3.3 × 10 − 7 m/s to 2.8 × 10 − 7 m/s. Combined with the CH distribution, shown in Figure 20, the hydraulic conductivity increases rapidly in the CH decomposition stage. After the complete decomposition of CH, the rise of the curtain hydraulic conductivity will slow down. For grout curtain C, because of the low initial CH content, after 100 a of calcium leaching, the hydraulic conductivity order of magnitude is still 10 − 8 m/s. ese results may provide a way to improve the calcium leaching durability of grout curtains. We can use additives such as fly ash or silica fumes to lower CH content in grout curtains [59].
Furthermore, a small number of admixtures cannot effectively slow down the increase of the hydraulic gradient of the grout curtain. It is necessary to reduce the CH content by adding a large number of admixtures. e less CH content, the better the calcium leaching durability of the grout curtains. Figure 22 shows hydraulic conductivity evolution curves of the grout curtain center line with different initial hydraulic conductivities, and the initial hydraulic conductivities are 0.01K, 0.1K, 1.0K, 10K, K � 1.0 × 10 − 8 m/s is the actual hydraulic conductivity of the grout curtain. e upstream and downstream water levels are 107 m and 86 m. e corresponding initial Pe are 5.1, 50.6, 299.4, and 544.5, and the leaching duration is 100 a. e initial CH and C-S-H content remain the same. As we can see from the figure, the smaller the initial hydraulic conductivity, the fewer increase times in the same time. After 100 years of leaching duration, the maximum increase is on the upstream side. For grout curtains with initial hydraulic conductivity 10 K and 1.0 K, the increase times are 67.8 and 66.9, respectively. e increased times of hydraulic conductivity decrease with the distance to the upstream side. However, for grout curtains with initial hydraulic conductivity 0.1 K and 0.01 K, the hydraulic conductivity increase times at the downstream side are higher than the middle part. When Pe number is more than 100, we can ignore the effect of diffusion. For these two curtains, the Pe numbers are 50.6 and 5.6. With the decrease of Pe number, the diffusion effect becomes more and more evident in the leaching process. If the hydraulic gradient is 0, the increased hydraulic conductivity times on the upstream and downstream sides should be symmetrically distributed. Eventhough the hydraulic conductivity evolution is significantly affected by the initial hydraulic conductivity, the final hydraulic conductivity is determined by calcium composition.    4.6. Discussion. In the existing advection-diffusion-driven calcium leaching model, the porosity is generally calculated with average molar volume, leading to inaccurate calculation of hydraulic conductivity. Decomposition of CH generates capillary pores, and decomposition of C-S-H causes gel pores and capillary pores. Furthermore, the influence of microstructure changes on hydraulic conductivity is also ignored. In this study, we proposed a new advection-diffusion-driven calcium leaching model of cement-based materials and considered that the effects of CH and C-S-H decomposition on pore growth are considered, respectively. In this way, on the one hand, we can use the model to study the influence of cement components on calcium leaching. On the other hand, the porosity of materials evolution can be determined more accurately. Moreover, we adopted the Kozeny-Carman (KC) relation to model the permeability coefficient evolution. In this way, the influence of microstructure changes on the permeability coefficient is taken into consideration.
KC relation requires evolution curves of bulk density and specific pore surface area in the leaching process. For this issue, Phung measured bulk density and specific pore surface area of intact and leached cement paste samples and assumed these two parameters are linearly varying with the decrease of Ca/Si [3]. We adopted this assumption in this study. However, for more accurate simulation, more data of bulk density and specific pore surface are needed.
In the previous literature, the transport properties of cement-based materials are assumed invariant in the dam seepage characteristic analysis. However, due to the calcium leaching phenomenon, the hydraulic conductivity and diffusivity of cement-based materials are increasing, so the seepage characteristics of the dam are not constant. With the increase of leaching degree, the seepage characteristics change constantly. To better understand concrete dam service performance in the whole life cycle, it is necessary to analyze the seepage characteristics in the leaching process. According to the seepage characteristics evolution of the Shimantan concrete gravity dam, the following engineering measures are suggested to increase the efficiency of antiseepage systems of concrete dams and prevent calcium leaching.
e main antiseepage structure of the concrete dam is the concrete face slab and the grout curtain. For the calcium leaching in the concrete face slab, we can spray the antiseepage coating such as polyurea elastomer on the upstream side of the face slab to prevent the reservoir water from penetrating. For grout curtains, it is impossible to avoid calcium leaching by spray coating materials. However, we can use admixtures such as fly ash and silica powder to reduce CH content. With a lower CH content, fewer capillary pores are produced by CH decomposition, and permeability coefficient increase in the calcium leaching process will significantly slow down.
In this study, the material properties of cement-based materials in concrete dams are porosity, diffusivity, hydraulic conductivity, and calcium content in the solid skeleton. e focus of this study is the calcium leaching effect on the dam seepage characteristics. In the actual operation of reservoir dams, many other factors affect the calcium leaching process, such as the stress state, water quality, and temperature. For example, the stress state of materials will affect the porosity and transport properties. Suppose the mechanical properties are considered in the analysis. It is necessary to provide the evolution of the constitutive model of the cement-based materials in the leaching process and the hydromechanical-chemical (HMC) coupling mechanism, which makes the model very complex. In this study, we ignored these factors to simplify the calculation. To solve the calcium leaching issues in hydraulic structures more accurately, the influence of these factors needs to be further studied.

Summary and Conclusions
We proposed a new mathematic model of hydrochemo behaviour of cement-based materials to investigate the concrete dam seepage characteristic evolution in the calcium leaching process. e proposed model had considered the effect of different calcium compounds decomposition on the porosity evolution. A three-dimensional finite element model of the Shimantan concrete gravity dam has been built. e water head distribution, leakage, and hydraulic gradient evolution in the leaching process were obtained. e main conclusions are as follows: (1) e simulation results of the proposed model are consistent with previous studies. Compared with the existing model, the effect of different calcium compound decomposition is taken into consideration. e proposed model removed the effect of gel pores produced by C-S-H decomposition on transport properties. erefore, our simulation of hydraulic conductivity and diffusion coefficient evolution is more accurate. (2) According to the seepage evolution of dam areas, the most important influence of calcium leaching on concrete gravity dams is the increasing uplift pressure along the dam support caused by the grout curtain's malfunction. e simulated leakage of the dam body and foundation showed good agreement with the monitoring data, indicating the rationality of this numerical modelling. (3) After 100 years of leaching duration, leakage quantity of dam body and foundation increase 48 and 17 times. e uplift pressure of the Shimantan concrete dam increased by 40.8%. e hydraulic gradient of concrete face slabs and grout curtains continued to decrease during the leaching process. (4) rough the sensitivity analysis, the initial hydraulic conductivity significantly affected the calcium leaching rate. Calcium composition determines the final hydraulic conductivity. To obtain better calcium leaching durability, we should reduce the CH content as much as possible.

A. Molar Volume of C-S-H in the Leaching Process
e molar volume of C-S-H in the leaching process is defined as where V CSH is the current molar volume of C-S-H, m 3 /mol, V CSH0 is the initial C-S-H molar volume, V CSH0 � 0.04 × 10 − 3 m 3 /mol, y 0 is the initial Ca/Si ratio, and y 0 � 1.7, y is the current Ca/Si ratio. In diffusion-driven leaching, y is a function of calcium concentration in pore solution [38]. However, in advection-diffusion-driven leaching, the calcium concentration in pore solution is generally lower than the equilibrium state. erefore, y is no longer a function of calcium concentration. In this study, we assume y is a function of current C-S-H concentration in the solid skeleton as y � 0.9 · C CSH x 1 /x 2 C CSH0 , 0.9 + 0.8 where C CSH is the current C-S-H concentration in the solid skeleton (mol/m 3 ), and C CSH0 , x 1 , and x 2 are same as equation (3).

B. Pore Structure Evolution in KC Relation
e specific surface area of pores can be obtained from the Brunauer, Emmett, and Teller (BET) nitrogen gas adsorption test, and the bulk density is obtained from mercury intrusion porosimetry (MIP) measurements. Phung [3] assumed that the specific surface area of pores and the bulk density linearly vary with a decrease in the Ca/Si ratio in a portlandite/C-S-H mixture. He also introduced the lumped term Ω to simplify the calculation of the shape factor F s and tortuosity τ. e lumped term present in the following equation.
where Ω 0 is the lumped term for sound materials, Ω l is the lumped term for leached materials, and d l is the degradation degree. e degradation degree is defined as where d l is the leaching degree, C 0 CH is the initial concentration of CH, (mol/m 3 ), and C CH is the current concentration of CH (mol/m 3 ).

Data Availability
All data, models, and code generated or used during the study are included within the article.

Conflicts of Interest
e authors declare that there are no conflicts of interest.