A Finite-Element Study of Elastic Filtration in Soils With Thin Inclusions

Soil environments are heterogeneous in their nature. This heterogeneity creates significant difficulties both in terms of construction practice and in terms of the mathematical modeling and computer simulation of the physical-chemical processes in these heterogeneous soil arrays. From the standpoint of mathematical modeling, the issue is the discontinuity of functions, which characterize the examined processes, on such inclusions. Moreover, the characteristics of such inclusions may depend on the defining functions of the processes studied (head, temperature, humidity, the concentration of chemicals, and their gradients). And this requires the modification of conjugation conditions and leads to the nonlinear boundary-value problems in heterogeneous areas. That is why this work has examined the impact of the existence of thin inclusions on the conjugation conditions for the defining functions of the filtration and geomigration processes on them. The conjugation condition for heads has also been modified while the mathematical model of an elastic filtration mode in a heterogeneous array of soil, which contains thin weakly permeable inclusions, has been improved. The improvement implies the modification of conjugation conditions for heads on thin inclusions when the filtering factor of the inclusion itself is nonlinearly dependent on the head gradient. The numerical solution to the corresponding nonlinear boundary-value problem has been found using a finite-element method. A series of numerical experiments were conducted and their analysis was carried out. The possibility of a significant impact on the head jump has been shown taking into consideration the dependence of filtration characteristics of an inclusion on head gradients. In particular, the relative difference of head jumps lies between 26 % and 99 % relative to the problem with a stable filtration factor for an inclusion. In other words, when conducting forecast calculations, the influence of such dependences cannot be neglected.


Introduction
Natural soils are heterogeneous at the macro level. This heterogeneity manifests itself in two aspects. The first aspect is that the soil array consists of sub-areas, which host soils that are different in their physical-chemical and mechanical properties. The second aspect is that in the soil array there are thin layers of other soils, which, quite often, differ greatly in their characteristics from the main soil. Such layers create significant difficulties both in terms of construction practice and in terms of mathematical modeling and computer simulation of the physical-chemical processes in these heterogeneous soil arrays. In particular, one of the practical aspects is the intensification of shear processes in which precisely these layers may serve as the sliding surfaces. From the point of view of mathematical modeling, the issue is the discontinuity of functions, which characterize the examined processes, on such inclusions. This requires the use of methods that make it possible to find generalized solutions to the relevant boundary-value problems. Moreover, the characteristics of such inclusions may depend on the defining functions of the processes studied (head, temperature, humidity, the concentration of chemicals, and their gradients). And this requires the modification of conjugation conditions and leads to the nonlinear boundary-value problems in heterogeneous areas. This, in turn, leads to changes in the forecast calculations of filtration processes in such soils and, as a result, affects their safe operation. At the same time, due to the nonlinearity of influences and complex interdependence of processes, it is advisable that studying such processes should involve computer simulation and mathematical modeling.
To solve this problem, we have modified the conjugation conditions, which has made it possible to improve the mathematical model of elastic filtration in the soil with a thin inclusion.

Literature review and problem statement
Study [1] noted that Darcy's law may be violated at high and relatively low filtration rates. The criterion in the first case is the Reynolds number, in the second -the viscosity of a liquid (attributing it to the class of non-Newtonian), and the initial head gradient.
The author of [2] suggested the following dependence between the filtering coefficient k and the head gradient I: where k 0 , k u is the initial and boundary (at I→∞) filtration coefficient; I k is the critical head gradient; k is the empirical constant, which, as noted by the author, is similar to the constant of semi-saturation known in biology and physical chemistry. Using the above dependence, the author derived the following important partial cases The nonlinear function of the dependence of the filtration coefficient on the head gradient used in work [3] was the power dependence in the following form where I is the head gradient. At β=1, we obtain a linear dependence; at β=0, we obtain k=k 0 . This dependence is similar to the function proposed in [4]. Dependences of the filtration coefficient on the head gradient are also given in [5]. The authors of [6] used Hansbo's flow mathematical equations for the filtration rate in saturated clays: c and k are the filtration coefficients at small and large head gradients; i 1 is the initial head gradient. The following dependence is used to show the dependence of the filtration coefficient on a change in the porosity coefficient where k 0 is the initial filtration coefficient corresponding to the initial porosity coefficient e 0 ; c k is the permeability index. In addition, the above work examined specific samples of clay soil both experimentally and using the constructed mathematical model based on the above dependences, which made it possible to find appropriate numerical solutions. There was a good agreement between the resulting approximate solutions and the data from field experiments. Study [7] reported an experimental study into water filtration in low-permeable porous environments at small head gradients [7]. It was established that in a separate case the movement of liquid is described by the nonlinear law. The phenomenon is elucidated by the presence of layers of bound water in the medium's pores. The explanation was substantiated and confirmed by field experiments.
It is noted in work [8] that, since the thickness of pores in clay soils varies from 0.01 to 0.1 μm, the interaction of "solid parts-liquid" significantly affects the flow in these pores. The result was the assumption, put forward and substantiated, about the nonlinearity of the filtration law in saturated clay soils. In particular, the law in the following form was proposed 2 1 , 1 where a 1 , a 2 , b are the parameters that are established experimentally; u is the filtration rate, h is the head. The authors of the work used the proposed law in the theory of consolidation. Experimental verification of that law was also carried out.
Study [9] numerically examined the mathematical model of filtration in a cracked-porous environment. In this case, it is believed that Darcy's law holds in a porous environment, while in the cracks, where the rate of filtration is greater, the Darcy-Forechheimer law is obeyed.
The movement of liquids in thin capillaries and the formation of a near-boundary layer of the bound fluid was experimentally investigated in [10]. Depending on the magnitude of the head gradient in the liquid, the thickness of such a near-boundary layer can be up to 30 % of the diameter of the capillary. As a result, at low pressure, the nonlinear effects start to manifest themselves. At the same time, the dependences for the velocity of fluid movement in capillaries differ from those consistent with Darcy's law. The authors generalized experimental research for porous environments and conducted numerical experiments for relevant tasks in the oil industry.
Several effects were taken into consideration in work [11]. First, nonlinear Darcy's law was used as a filtration law (non-Darcy Hansbo's flow). Second, the permeability coefficient was considered variable and dependent on the coefficient of porosity. These dependences are taken into consideration in the general theory of soil consolidation (general theory of consolidation of the theory of Bio-Bio).
Under Hansbo's nonlinear filtration law, where u is the filtration rate, k is the filtration coefficient at relatively large gradients of head, when Darcy's law holds; k is the filtration parameter used in the power part of the nonlinear law; m is the indicator of power in the filtration law at small gradients of head, i 0 is the initial head gradient; i l is the gradient that makes it possible to overcome the maximum energy of the bound water and set it into motion.
Taking into consideration the condition of continuity and differentiation, at i=i l , we obtain 1 , when .
As regards the filtration coefficient, the following dependence is proposed α is in 0.5 ⁓2, and the total value in 0.5 ⁓1.0 x y z σ σ σ ′ ′ ′ are the effective stresses; p′ is the average effective stress; 0 p′ is the initial value of the average effective stress.
Paper [12] considered the issue of filtration in poorly permeable soils [12]. The mathematical filtration model takes into consideration the presence of the initial head gradient and the nonlinear dependence of porous fluid density on head. As a result, the first factor required the modification of Darcy's law while the second factor led to the inclusion of a filtration term in the equation, which depends on the square of the head gradient in the porous liquid. Moreover, the existence of the initial head gradient predetermined the need to study the problem in the area with a moving boundary. The reported numerical experiments showed a significant difference in predictive calculations once the quadratic term of the head gradient is neglected in the presence of the initial head gradient.
Special attention should be paid to the issue related to thin soil inclusions and conjugation conditions for the defining functions of the filtration and geomigration processes on them.
Thus, the soil array can contain thin inclusions from natural soils. Also, such inclusions are created artificially using geotextile [13]. As shown in [14], the presence of such inclusions significantly changes the forecast calculations of filtration processes, heat-salt transfer processes, the stressed-strained state in natural and artificial soil arrays. In terms of mathematical modeling and numerical methods, the existence of such inclusions requires the modification of approaches to their (processes) mathematical notation, the statement of proper boundary-value problems, their numerical solution.
At present, the scientific literature provides a large body of research into the properties of argillaceous inclusions. Quite often, such inclusions are created as artificial anti-filtration barriers in the design and construction of production waste storage facilities in order to minimize the propagation of pollutants into groundwater. Work [15] describes the properties of such geosynthetic argillaceous inclusions. The authors noted their low hydraulic conductivity and the ability to detect semi-permeable membrane properties, which limit the migration of solutions through bentonite.
Paper [16] examines bentonites that make it possible to limit the migration of pollutants. The so-called Fukakusa clay with different amounts of dry bentonite (5 %, 10 %, 15 % and 20 %) was investigated.
Article [17] explored the interconnected processes of clay consolidation and the migration of contaminated substances through an argillaceous inclusion with semi-permeable properties. The result of numerical modeling showed that depending on the conditions, consolidation can exert a significant impact on the breakout time of pollutants. It also has an impact on the distribution of concentration in argillaceous inclusions not only during the consolidation process but also long after its completion. The effect of consolidation increases with an increase in the thickness of an inclusion, an increase in the applied load, a decrease in loading time, an increase in the values of effective diffusion coefficient, etc.
Work [13] reviewed different types of diffusion coefficients regarding the migration of radionuclides in geoprotective barriers. Three potentially significant complex problems have been identified: a barrier system geochemistry, the effect of a surface and/or interlayer diffusion, and the presence of semi-permeable membrane properties as a result of the inclusion's anions behavior.
Paper [18] addresses argillaceous barriers to waste storage facilities, behaving like semi-permeable membranes. The equations of fluid movement and migration of saline solutions through the so-called argillaceous membrane barriers are given. It was established that the flow of the solution (pollutant) decreases compared to the case when the inclusion has no semi-permeable properties of the membrane. The issues relating to diffusion through argillaceous membranes were considered in [19].
When constructing mathematical models of geomigration processes in environments with inclusions, the corresponding boundary-value problems should be correct. Therefore, such models should be supplemented with the conjugation conditions for unknown functions on such inclusions.
Papers [20,21] examine parabolic equations that describe the process of heat-and-mass transfer in an area with a thin low-permeable inclusion under the conjugation conditions of the "concentrated natural source" type. A new statement of such problems was proposed, in which the principal parabolic equation changes to a system of differential equations of the first order with coefficients from classes of the generalized functions. The authors proved theorems of the existence and uniformity of a generalized solution.
Work [22] investigates a problem about the diffusion of an impurity in a band containing an inner layer of other material. It is believed that at the contact boundaries of layers, the concentration of impurity particles changes in a jump-like manner, although their flow through these boundaries is continuous.
The issues of mathematical modeling of the filtration and geomigration processes in heterogeneous environments are systematically studied in [14]. In particular, the conjugation condition for heads in the filtration problem is derived from the assumption that the head h is linearly changed from hto h + along the normal n to the cross-section of the thin inclusion γ of thickness l. That is, Here h + and hare the head values on an inclusion at ξ=l and ξ=0, respectively. The coordinate system Oξ is tied to the normal to the inclusion. Then the limit conjugation condition for a non-ideal contact takes the following form [14] ( ) In the above condition, a constant quantity is used as the filtration coefficient k 0 of a thin inclusion. However, the soil filtration coefficient depends on the head gradients. And this requires a modification of the conjugation conditions.

The aim and objectives of the study
The aim of this study is to quantify the impact of the existence of thin inclusions and the modified conjugation conditions for the defining functions of the filtration and geomigration processes on them. This would make it possible to draw conclusions about the degree of significance or insignificance of changes in the forecast calculations of filtration characteristics in a heterogeneous array of soil, which contains thin weakly permeable inclusions.
To accomplish the aim, the following tasks have been set: -to build an improved mathematical model of an elastic filtration mode in a heterogeneous array of soil, which contains thin weakly permeable inclusions; -to find numerical finite-element solutions to the corresponding nonlinear boundary-value problem with the modified conjugation conditions, which describes the constructed mathematical model; -to conduct a series of numerical experiments and analyze them for the quantitative impact of the presence of thin inclusions and modified conjugation conditions on head jumps on these inclusions.

Mathematical model of elastic filtration in the soil with a thin inclusion
Consider the process of elastic filtration in a heterogeneous array of soil with a thin inclusion in a one-dimensional case, whose mathematical model is described by the following boundary-value problem [14]: , ; , , , ; , , , , , .
Here, Ω 1 =(0; ξ), Ω 2 =(ξ; l), 0<ξ<l; h l (t), h 0 (x), Q(h, t) are the known functions, h(X, t) is the head in a porous liquid; k is the filtration coefficient; f(X, t) is the function that sets the intensity of the internal sources (drains) of a liquid. The coefficient η is called the coefficient of rock elastic capacity. According to [23,24], at short depths (10- In this mathematical model, formula (11) is a modified conjugation condition for heads, derived on the basis of the following considerations.
By analogy with [25], suppose that due to the subtlety of an inclusion the filtration processes in the cross-section of a given inclusion are stationary (or at least quasi-stationary). Thus, consider the following filtration problem for an inclusion: Here, k γ (h, ∇h) is the filtration coefficient of a thin inclusion that is nonlinearly dependent on heads and their gradients. From equation (6) (8) where h 2 =const is also a not yet known constant. Considering (8) and boundary conditions (7), we obtain a system of linear algebraic equations (SLAE) According to [14], the conjugation condition will be derived on the basis of the law of preservation of the flow of liquid through a single area of the transverse cross-section of the surface of an inclusion along the normal over a small period of time Δt. Since the flow (10) and (9) produce the final modified conjugation condition for a non-ideal contact for heads on an inclusion, whose filtration coefficient depends on the heads themselves and their (head) gradients (11) produces a classic conjugation condition for a non-ideal contact [5].

A finite-element scheme to solve the problem of elastic filtration in the area with an inclusion
Similarly to [14], we shall introduce a series of definitions and spaces. Let H 0 be the space of the functions ( ), x φ which, at each interval (0; ξ), (ξ; l), belong to the Sobolev space ( ) x H ∀ϕ ∈ satisfies integral relations (12), (13), is termed the generalized solution to the boundary-value problem (1) to (5).
The approximate generalized solution to the boundary-value problem (1) to (5) is found in the form .
Solving problem (1) to (5) by a finite-element method from a weak statement(12), (13) of the problem, taking into consideration (14), we obtain a Cauchy problem After time sampling involving a completely implicit linearized difference scheme [26], system (15) produces where τ is the time step, A (p) =A(t p ), t p =pτ. One can calculate the integral in the l ij definition by using quadrature formulae and head values and their gradients' values from the previous time step.

Results of numerical experiments and their analysis
The software implementation of the algorithm for solving the above boundary-value problem employed the modern Python visual programming technology.
In the simulation problem, the soil layer l=30 m thick was considered. The depth of an inclusion level is ξ=15 m. A step in the variable x is 0.1 m. A time step is=0.5 day. The initial distribution of heads is h 0 (x)=1 m. At the upper boundary of the soil, we set the boundary condition of the first kind with a head value of 10 m. The following soil filtration and inclusion parameters were set: k 0 =0.01 m/day, 0 0.0001m day, k γ = 4 1 5 10 . m − η = ⋅ Results of the numerical experiment given in Table 1 outline the following trends. As the thickness of an inclusion decreases, the head jump value decreases. However, the relative differences in head jumps for the cases of a variable and stable filtration coefficient increase from 26 % to 43 %. The "-" sign before relative differences means that they are reduced for the case of a variable filtration coefficient. The results of solving this problem, although it is a model, show that with such relative differences one must not neglect the dependence of the filtration coefficient on the head gradient for thin inclusions. Table 2 gives the values of heads and their jumps for a model problem (case 2). Table 2 The values of heads and their jumps in a numerical experiment for the case k(I)=k 0 I β , β=1 The head values at k γ =const correspond to the use of the classic conjugation condition, and, at k γ =k γ (∇h), to the modified. Parentheses show a relative increase or decrease (as a percentage) in the head jump when applying the classic conjugation condition and the modified conjugation condition. In Table 2, such relative changes range from −58 % to −76 %. At the same time, with an increase in the thickness of an inclusion, the amplitude of differences increases while the relative difference falls. The results of solving the model problem show that one must not neglect the dependence of the filtration coefficient of thin weakly permeable inclusions on the head gradients. Table 3 gives the values of heads and their jumps for a model problem (case 3). Table 3 The values of heads and their jumps in a numerical experiment for the case k(I)=k 0 I β , β=2  Table 3 show that at the quadratic dependence of the filtration coefficient on the head gradient, the head jumps decrease while the relative differences, compared to the problem with a stable filtration coefficient, increase.

Discussion of results of studying the problem of elastic filtration in a heterogeneous soil array
By comparing the results given in Table 1 and Table 2, one can note that the application of the classic (k γ =const) and modified (k γ =k γ (∇h)) conjugation conditions produces different indicators of the head jump value. This difference ranges from −58 % to 76 %. In addition, all tables show that increasing the thickness of a weakly permeable inclusion in the soil array increases the amplitude of differences in the head value from above (h + ) and the head value from below (h -). Although, in this case, the value of the relative difference falls. In addition, the results of our study provided information on the form of setting a dependence of the filtration coefficient on the head gradient. According to data in Table 2 and Table 3, with the quadratic dependence of the filtration coefficient on the head gradient, the head jumps are smaller compared to its linear dependence. However, there is an increase in the relative difference compared to the problem with a stable filtration coefficient.
The application of the dependence of the filtration coefficient on the head gradient in the classic conjugation conditions on a thin inclusion causes a series of contradictions. In particular, to the left and right of the inclusion, the head gradients are different. And there is currently no justification for what kind of dependence to use (to the left or right of an inclusion). In this case, the use of the modified integral conjugation condition has allowed us to solve the issue associated with the nonlinear filtration parameters of thin argillaceous inclusions regarding the head gradients. The proposed conjugation condition takes into consideration the results of field experiments and is free from contradictions of the classical conjugation condition for the case of nonlinear dependences. The justification of the need to modify the conjugation condition, as well as the results in Tables 1-3 for model problems, indicate the prospects for further research in this area. This is especially true of argillaceous soils as materials of thin inclusions whose nonlinearities are strongly expressed in filtration processes.
Our study could be advanced in the following ways: -to study the impact of other nonlinear dependences of the filtration coefficient on the gradient head on the magnitude of head jumps; -to study the impact of the derived dependences on other physical-chemical processes in soils (for example, heat transfer); -to theoretically investigate the accuracy of the derived finite-element solutions to boundary-value problems with the modified conjugation conditions.

Conclusions
1. We have built an improved mathematical model of the elastic filtration mode in a heterogeneous soil array containing thin weakly permeable inclusions. The improvement implies the modification of conjugation conditions for heads on thin inclusions when the filtration coefficient of the inclusion itself is nonlinearly dependent on the head gradient. After all, the use of classic conjugation conditions for such nonlinear dependences faces the ambiguity of choice whether a head gradient to the left or right from a thin inclusion should be considered. The modified conjugation condition has been derived as a solution to the nonlinear boundary-val-ue problem in a thin inclusion and thus is free from such contradictions.
2. The numerical solutions to the corresponding nonlinear boundary-value problem have been derived by using a finite-element method (FEM). The schematic algorithm for the finite-element solution to the resulting problem is given. The possibility of discontinuous solutions requires the use of numerical methods that make it possible to find the approximate generalized solutions to the relevant boundary-value problems. This is what justifies the application of a finite-element method in our study. 3. Numerical modeling has shown that taking into consideration the dependence of filtration characteristics of an inclusion on the head gradients has a significant impact on the magnitude of head jumps. In particular, the relative difference in head jumps lies between 26 % and 99 % relative to the problem with a stable filtration coefficient for an inclusion.