Abstract
Permeability of soil plays an important role in geotechnical engineering and is commonly determined by methods combining measurements with theory. Using the double-scale asymptotic expansion method, the Navier-Stokes equation is numerically solved to calculate the permeability, based on the homogenization method and the assumption that the homogeneous microstructure of the relevant porous media is represented accurately as the Representative Elemental Volume (REV). In this study, the commonly used square model is tested in the calculation of sea clay permeability. The results show large deviations. It is suspected that the square model could not represent the flattened shape of the clay particles and the bound water film wrapping around them. Hence, the Rectangle Particle-Water Film Model (i.e., the R-W model) is proposed. After determining the horizontal and vertical characteristic length of the unit cell using two pairs of initial data, the permeabilities of other different void ratios could be inversely calculated. The results of three types of clay obtained using the R-W model agree well with the experimental data. This shows the efficient feasibility and accuracy of the R-W model by providing a good representation of the clay particles when using the double-scale asymptotic expansion method to calculate clay permeability.
1. Introduction
Permeability of soil is of fundamental importance in geotechnical engineering. Experimental methods, such as the classical Darcy law, are commonly used to obtain the permeability of soil. These methods are convenient to be conducted but always require much time and effort. Moreover, the experimental methods focus only on the external phenomenon from the macroscopic view and ignore the seepage processes inside the porous materials at the microscopic scale. Actually, as a type of porous material, soils can be distinguished at three different scales: the microscopic scale, the Darcy scale, and the macroscopic scale [1]. To produce a more comprehensive and thorough analysis, the seepage process at both the macroscopic and the microscopic scale should be studied.
The homogenization method [2], which is also called the multiscale asymptotic expansion method, is a powerful tool to investigate the characteristics of nonhomogeneous materials. Starting from the physical phenomenon at the microscopic scale, the homogenization method can obtain the characteristic descriptions of the material at the macroscopic scale using asymptotic expansion of the governing equations. The homogenization method has been widely used in the researches of characteristics of composite materials, which generally have a uniform and regular structure. Andreassen and Andreasen [3] used short and self-contained Matlab implementation to provide the calculation method of the elasticity tensor and fluid permeability of composite materials by homogenization. Keip et al. [4] presented a two-scale computational homogenization framework for the simulation of electroactive solids at finite strains. The results of these researches showed that the homogenization method was typically well suited to study composite materials.
In geotechnical engineering, Wang et al. [5] proposed a simplified homogenization algorithm for composite soils to simulate the characteristics of soils via numerical calculations which assumed that the pressure in the microscopic scale was homogeneous. Specific to the seepage in soils, Tang et al. [6] proposed a multiscale method to calculate the permeability coefficient of soils. Based on the homogenization method and the assumption that the homogeneous microstructure of the soils is well represented by the Representative Elemental Volume (REV), the Navier-Stokes equation was numerically solved. Sun et al. [7] used the double-scale asymptotic expansion method to inversely calculate the permeability of clay using a square particle unit cell as the REV. The results were within the allowance but presented significant errors. Considering the calculation of permeability for kaolin clay as an example, the calculated permeability was m2 at the largest void ratio , while the measured permeability m2; also, the calculated value was m2 when the smallest void ratio was determined to be , while the measured value m2.
In this paper, the existing square model, which was used in the inversed calculation by the double-scale asymptotic expansion method, was tested to find the defects of the model. Additionally, a new model would be proposed based on the real features of clay particles aiming to improve the calculated accuracy of clay permeability.
2. Calculation of Clay Permeability by Multiscale Method
2.1. Multiscale Expansion of Navier-Stokes Equations
Previous studies have shown that, for a given type of statistically homogeneous soil, an equivalent REV always exists to represent the structure of the soil. A similar description could be used for nonhomogeneous porous materials [8, 9]. The basic assumption of the calculation performed in this study was that the porous media are homogeneous, and the arrangement of their particles was periodic; thus, an REV should be chosen to represent the porous media, as shown in Figure 1. The pores are fully saturated by an incompressible Newtonian fluid with a small Reynolds number; denotes the solid part of the particle, the space of the fluid in the pore volume, and Γ the interface between the solid and the fluid.
The Stokes flow through the porous media is governed by the Navier-Stokes equation, the continuity equation, and the no-slip boundary condition:where is the velocity vector, is the pressure, and is the viscosity.
Introduce two-scale coordinates as macroscopic variable and as microscopic variable, which are linked by , where is an adequately small ratio between the characteristic lengths of the unit cell and the macroscopic sample. Then, the asymptotic expansions of the velocity and the pressure can be expanded as follows:where and are periodic. Incorporating expansions (2) in (1a)–(1c) and using , (3a)-(3b) could be obtained:From the mathematical view, (3a)–(3c) hold for any value of , which means that all the coefficients should be zero. Hence, and could be obtained, with the following equations:Equations (4a)–(4c), which have a similar form to (1a)–(1c), could be regarded as a special Stokes flow problem in the unit cell. In this special Stokes flow, the macroscopic pressure gradient is the external excitation, while the first-order velocity and the second-order pressure are the corresponding responses. Considering the linearity of (4b), the general solution could be proposed to take the linear function as follows:where is a pore-scale permeability tensor. Considering the term in (3b), it can be obtained that .
Averaging over in the unit cell, (6) could be obtained:where is a volume average velocity in a macroscopic scale and is a macroscopic permeability tensor which must take the form , where is the scalar permeability and is the unit tensor.
2.2. Calculation of the Permeability
Equations (4a)–(4c) could be solved once the unknowns and are periodic, and is a given source term. Then, averaging the resulting velocity over the unit cell, the value of could be determined as follows:The value of could be obtained from (6), and the calculative process could be operated numerically.
2.3. Inverse Calculation
For a given unit cell, as shown in Figure 1, the soil permeability is a function of the void ratio and the characteristic length of the particles present. In a given type of porous media, the characteristic length is a constant; thus, using one pair of data points (i.e., void ratio and permeability), which is measured in the laboratory or in situ, the characteristic length could be calculated numerically. Finally, with the calculated characteristic length, the permeability for various void ratios could be computed by changing the distance between the particles.
The specific procedure is as follows: according to the degrees of freedom of the REV, take one or more pairs of test data (void ratio and the corresponding permeability ) as the initial data. Assign different values to the characteristic length of the REV and calculate the permeability via computer software (COMSOL in this study) until , where is an adequately small value to control the calculating accuracy. Once the characteristic length of the REV is confirmed, the permeability corresponding to other void ratios could be calculated by fixing the characteristic length and changing the distance between the REV particles.
3. Calculation of Clay Permeability
3.1. Calculations of Sea Clay Using Square Model
For the calculation of clay permeability, the REVs are commonly set as regular polygon particles. Sun et al. [7] used a circle particles model and a square particle model as the REVs to inversely calculate the permeability of sand and clay. The calculated values of sand agreed well with the measured values, while the calculated values and the measured values did not agree well for clay. This disagreement showed that there must exist some defects in the REVs of the clay. The defects should be studied and eliminated as far as possible.
Firstly, the REV with square particle, which Sun et al. had used in their study, was tested again to help us to find the defects. As Figure 2 showed, there are 2 parameters in the square model which are and . Since could be calculated when is confirmed along with the void ratio , the REV has 1 degree of freedom. Using 1 pair of measured values (), the characteristic length could be confirmed.
For an equitable comparison, the Champlain sea clay of Louiseville [10] was chosen as the example in our studies, which Sun et al. [7] used in their studies as well. Firstly, the measured value (, m2), which is the largest void ratio and the corresponding permeability tested in the reference, was used as the initial data point; these allowed the characteristic length of the square model to be obtained:Holding the aforementioned characteristic length of the square model as a constant, the permeabilities corresponding to other void ratios could be calculated.
Similarly, the measured value that has the smallest void ratio and the corresponding permeability (, m2) was used to determine the characteristic length of the square model again:Then, the permeabilities corresponding to other void ratios could be calculated in the same way. The 2 groups of calculated values are listed in Table 1, along with the measured values.
Using the measured void ratio as the abscissa and the ratio of the calculated permeability to the measured permeability as the ordinate, the values in Table 1 could be shown in Figure 3.
From the documents published by ASTM [11–13], the values of are generally between 1/3 and 3 times the accurate permeability for one type of clay. As shown in Figure 3, the ratios were within that range, which indicates that the double-scale asymptotic expansion method could be used effectively in the inverse calculation of the clay permeability. However, orderly deviations are also shown in Figure 3.
Theoretically, if the accuracy of the REV is sufficient enough, the results of the characteristic length will not vary observably, irrespective of whichever the pair of data with the largest or smallest void ratio is chosen as the initial value. Comparing (8) and (9), the calculated characteristic length m when the initial and m2 was significantly larger than m when the initial and m2. Additionally, the values of calculated are all larger than the measured permeability when the initial and m2; these deviations increased along with the differences between the used void ratio and the initial void ratio. Similarly, the values of calculated are all less than when the initial and m2; the deviation also increases with the magnitude of the void ratio difference. Fitting the calculated results into straight lines, the slopes of the lines were for m and for m. The slopes were too large to control the calculated accuracy when the difference between the used void ratio and the initial void ratio was large.
The inaccuracies described above should be attributed to the lack of representativeness of the square model. Generally speaking, clay is composed of particles with flattened shape [14]. Hence, the REV should not be in square shape and the model should be modified to represent the clay particles more accurately.
3.2. Calculations for Sea Clay Using the Rectangle Particle-Water Film Model
As Figure 4 shows, clay particles are generally in flattened shape. Hence, the regular polygon shape of REV would reduce the representativeness. To improve the REV, the square particle model was modified into a rectangle one, which could represent the flattened shape of the clay particles.
Additionally, a thin bound water film exists wrapping around the soil particles for both sand and clay at the microscopic scale. The water film acts as a solid due to its low energy level and its inability to move freely and dissolve solutes. During permeation, the water film can be regarded as part of the soil particles. However, during the measurement of the void ratio, the water film would be evaporated and could not be considered in the void ratio. Thus, for the same void ratio, the existence of the water film reduces the passage of permeation and decreases the permeability accordingly.
Synthesizing the flattened shape of the clay particles and the bound water film, the Rectangle Particle-Water Film Model (R-W model) was proposed, as shown in Figure 5. A lamina representing the water film was added to wrap around the rectangle model. The lamina is considered to be a part of the particle when calculating its permeability and is ignored when calculating its void ratio.
Figure 5 shows that there are 6 variable factors to describe the rec. model: the vertical and horizontal characteristic lengths of the particle and , the vertical and horizontal characteristic lengths of the unit cell and , and the vertical and horizontal characteristic thicknesses of the water lamina and . For convenience, the following assumptions were taken into consideration:(1)The vertical and horizontal lengths of the unit cell and maintained the same ratio to the particles, which indicates that(2)Similarly, the vertical and horizontal thicknesses of the water lamina and maintained the same ratio to the particles, which indicates thatIn Section 3.1, it had been found that the characteristic length of the square model m, which was calculated using the initial data with the largest void ratio, would enlarge the results of permeability, while m calculated using the initial data with the smallest void ratio would reduce the calculated results. It could be speculated that the accurate REV should have a characteristic length in between and . Assigning and , then the other 4 parameters could be calculated with void ratio , the hygroscopicity index, and (10) and (11).
From the hygroscopicity index [15], 0.1 cm3 of hygroscopic water is a typical average for every 1 cm3 soil particle for sea clay. This means that, in the progress of permeation, the solid part increases 1.1 times. Converting it to the R-W model, the size of the unit cell combined with the lamina is 1.05 () times the size of the unit cell without the lamina and it could be obtained that . Then, using the modified R-W model, the results of the permeability could be calculated. The results are shown in Table 2 with the associated void ratio and measured permeabilities.
Similarly, using as the abscissa and of the rec. model and the R-W model as the ordinate, Figure 6 could be obtained.
From Figure 6, it is shown that the calculated results for the permeability of the sea clay are improved significantly when using the modified R-W model. Fitting the results into a straight line, the slope was reduced to 0.12, which was much smaller than the square model in Figure 3 to control the calculation to be more accurate.
3.3. Calculations of Kaolin and Illite Clay Using R-W Model
To verify the applicability of the R-W model, the permeabilities of kaolin [16] and illite [17] clay were calculated. Using the hygroscopicity index of kaolin clay and illite clay [15], the size of the unit cell representing the kaolin clay combined with the lamina was 1.05 times the size of the unit cell without the lamina, while that of the illite clay was 1.10 times larger. The calculated results of the kaolin clay and illite clay are listed in Tables 3 and 4, respectively. The calculation of illite clay ignored the significantly inaccurate data pair that occurred when and m2.
Using as the abscissa and of the rec. model and the R-W model as the ordinate for the kaolin clay and illite clay separately, Figure 7 can be obtained.
(a) Kaolin clay
(b) Illite clay
From Figure 7, the ratios of the calculated permeability to the measured permeability () of the kaolin clay and illite clay using the R-W model were all between 1/3 and 3, which demonstrated the feasibility and accuracy of the R-W model in the permeability calculation for clay using the double-scale asymptotic expansion method.
It should be pointed out that the measured permeability of the illite clay was indirectly derived by the consolidation curve of Terzaghi’s Consolidation Theory. This result is shown to have poor stability in the literature [18], which can cause large fluctuations in of illite clay, as shown in Figure 7(b). Conversely, the measured permeability of kaolin clay is obtained using the conventional or modified constant-head test, which produces more accurate results. Thus, the inverse calculation of kaolin clay using the double-scale asymptotic expansion method provides better stability and accuracy, as shown in Figure 7(a).
4. Conclusions
The double-scale asymptotic expansion method used in this study can numerically solve the Navier-Stokes equation and determine the soil permeability based on the assumption that the porous media are homogeneous, and the microstructure of the media is known to be an incompressible Newtonian fluid with a small Reynolds number. Conversely, using several measured values of the void ratio and the corresponding permeabilities of the media, the characteristic length of a reasonable REV can be calculated. Then, the permeability of the other void ratio can be calculated by changing the distance between the unit cells of the REV.
For the calculation of sea clay permeability, the existing square model produced significant deviations. Using the measured value with the largest void ratio as the initial data would lead to a larger characteristic length of the REV and larger calculated results, while using the smallest void ratio as the initial data would lead to smaller results. The straight lines which were fitted by the calculated results had a large slope, which means the deviations would increase more significantly along with the differences between the used void ratio and the initial void ratio. It was speculated that the square model could not represent the flattened shape of the clay particle and ignored the bound water film wrapping around the clay particles.
Thus, the R-W model was proposed, which could represent the flattened shape of the clay particles and the water film around them. Using the measured values with both the largest and the smallest void ratio along with the corresponding permeability and the known hygroscopicity index, the horizontal and vertical characteristic lengths of R-W model can be confirmed. Then, the permeabilities of the other void ratios can be inversely calculated by changing the distance between the unit cells based on the same ratio of the horizontal and vertical characteristic length. The example, where the same sea clay was used, shows higher accuracy and a smoother slope of the fitted straight line using the R-W model. Additionally, the permeabilities of two other types of clay, which were kaolin clay and illite clay, were also calculated. The calculations using R-W model also show satisfactory results, which prove the feasibility and accuracy of the R-W model for the calculation of clay permeabilities using the double-scale asymptotic expansion method.
Competing Interests
The authors declare that they have no competing interests.
Acknowledgments
The authors would like to express sincere gratitude to the National Natural Science Foundation of China (Grant no. 51179168) and the Cultural Relics Protection Project of Zhejiang (Grant no. 2013010) for their financial support of this study.