Influence of Self-Heating on Landfill Leachate Migration

Abstract

The hydrodynamic processes of landfill leachate migration in the base of a solid waste landfill can have a critical impact on the natural environment. In the case of improper operation of municipal solid waste placement facilities, highly contaminated leachate may penetrate into groundwater and subsequently into surface water. This work addresses fundamental issues of multicomponent fluid propagation in a multilayer porous medium, taking into account temperature inhomogeneities caused by waste decomposition with heat release. The regimes of diffusion and convection of leachate water penetrating into soil layers in the base of municipal solid waste facilities are numerically studied. Archival data from a set of field and laboratory measurements in the area of the operating landfill are used to model the features of pollutant propagation and determine migration parameters. The process of pollutant propagation and migration is described by quantitative values of dry residue content in leachate. Factors that have a significant impact on the migration of leachate are considered. The main ones are convective transfer, diffusion, and properties of the geological composition of the landfill base, which are taken into account in the mathematical formulation of the problem. The calculations show that leachate self-heating can substantially intensify leachate filtration and has to be taken into account in the assessment of leachate migration rate.

1. Introduction

The increase in the volume of municipal solid waste (MSW) is a serious concern for the global economy as it leads to significant environmental and economic challenges that cause the need for development of effective mechanisms to manage and storage MSW. Currently, despite the introduction of advanced waste processing methods, the most common method of the treatment is to dispose of waste in landfills. The most obvious negative impact of obsolete waste disposal facilities on the environment is the effect of the seepage water (also referred to as leachate) on groundwater and surface water, as well as the soil massif in the foundation of MSW landfills. MSW is known to have a complex, multicomponent and heterogeneous composition, depending on the region where a landfill is located [1,2,3,4]. Leachate is characterized as a highly toxic liquid that could include heavy metals, organic pollutants, ammonia nitrogen, and other components [5,6]. Leachate is formed during the interaction of waste with atmospheric precipitation seeping into the body of the landfill, as well as a result of other processes occurring in the body of the landfill for waste placement (squeezing moisture, moisture release as a result of chemical reactions, etc.).

In ref. [7], an analytical solution for unsaturated one-dimensional diffusion of contaminants was derived. It was demonstrated that, when the groundwater is deeper than approximately 3 m, it would be difficult for the contaminants to break through the underlying soil layers with higher desaturation rates below the landfill barriers to reach the groundwater. Thus, thicker soil layers of high desaturation can be an effective diffusion barrier used in the design of a landfill. A landfill barrier system can be considered to have been breached if only one pollutant exceeds the threshold concentration [8]. The pollutant that exceeds the threshold concentration can be regarded as the key pollutant indicator for identifying breakthrough for a landfill site. In ref. [9], batch adsorption tests and geotechnical centrifuge modeling were performed to examine the adsorption intensity of different pollutants transported through the landfill’s compacted clay liner. It is shown that leachate migration through MSW landfill liners can be slowed considerably by adsorption.

The quantitative evaluation of leachate accumulation is critical for selection of drainage systems and reservoirs used for temporary storage of leachate during the construction of waste placement facilities. There is a wide range of studies on the assessment of leachate formation, including experimental and field measurements and mathematical modeling [10,11]. Since the end of the last century, one of the most popular models of dynamic landfill simulation is the Hydrologic Evaluation of Landfill Performance (HELP) model [12]. It uses a quasi-two-dimensional hydrological model of fluid filtration inside and outside the landfill [13,14]. Another solution is the LandSim (Sandfill Performance Simulation) program, which allows modeling the migration of pollutants from a point source at a landfill to soil and groundwater over a long period of time [15,16].

In recent years, a number of studies have shown that the HELP model has limitations and shortcomings [17]. For example, HELP does not take into account waste compaction and changes in their physical and mechanical properties and simulates the transport and formation of leachate under steady-state conditions with constant parameters. In this regard, detailed integrated models are being developed. For example, ref. [18] presents a three-dimensional model of a landfill consisting of layers formed by square cells of horizontal cross-section. It takes into account the effects of surface runoff, changes in hydraulic conductivity of waste with depth, and biodegradation. In ref. [19], an alternative water balance model is proposed for predicting the quantitative assessment of leachate formation at operating and closed landfills. The model takes into account different inflow rates, water losses, and demands, as well as aging and compaction, which allows the estimation of the change in hydraulic and physical properties of the buried waste. The results obtained showed that waste compaction significantly affects the prediction of leachate production during the landfill operation phase, and ignoring these processes can lead to an underestimation of up to one order of magnitude. An improved version of this model, incorporating new and more detailed approaches to landfill volume discretization, different waste disposal methods, a new surface water balance for actual evaporation and runoff, and a weather generator, is presented in ref. [20]. The importance of taking into account water squeezing from waste for quantitative estimation of leachate production was also shown in refs. [21,22].

The understanding of leachate migration patterns beyond MSW landfills in the case of leakage is essential for environmental protection. Modern landfills are equipped with a protective bottom layer (screen) that prevents penetration of leachate into groundwater. At the same time, in areas with soils with low permeability, such as clay soils, obsolete landfills could use rocks as a landfill screen [23,24]. Although those landfills are also equipped with a drainage system, which allows leachate to smoothly drain into temporary storage tanks for further treatment [25,26,27], the comprehensive analysis of the possible leachate diffusion beyond the screen is extremely relevant for these objects [16,28,29,30].

In ref. [31], the model of groundwater flow and leachate transfer was developed to identify the risk of groundwater pollution. The results show that this impact largely depends on the hydrogeology of the site, the volume of water entering the aquifer, and the concentration of the pollutant at the source. In ref. [32], leachate migration at an unsealed MSW landfill in a granite region was studied. In particular, two landfills located in the granite zone, where groundwater is heavily polluted, were selected as study objects. Based on field studies, the hydrogeological model of surveyed landfills was developed. It was found that the main sources of leachate migration on the side slope of the granite valley are granite cracks and highly permeable rocks at the mouth of the valley. The leachate plume was observed reaching more than 200 m from the lower boundary of the landfill. An attempt to determine the degree of groundwater contamination in a coastal aquifer was made in ref. [33]. For this purpose, the MODFLOW 2000 and MT3DMS 5.2 software products were used and the applicability of a permanent reactive barrier made of nanoiron around the landfill was assessed as a possible alternative to reduce the concentration of pollutants. Predictive modeling shows that in this case, the concentration of pollutants in groundwater can be reduced. Despite the large number of studies on the process of leachate formation, the process of leachate distribution in the foundation of the landfill base soil rocks in a three-dimensional approach is insufficiently studied.

Heat generation in the landfill body is one of the important factors influencing MSW landfills. Considerable volumes of heat can be produced in the early phases following the disposal of garbage in landfills, mostly as a result of food waste decomposition [34]. The generated heat also leads to a long-term increase in waste temperature, although aerobic and anaerobic reactions can be considered as the main sources of heat in landfills [35]. During field mesuaments, MSW landfills have been found to have temperatures rise up to 80 °C [36,37,38]. Elevated temperature landfills are defined as landfills that have high temperatures across a wide area [39]. Temperature effects on the hydraulic properties of unsaturated rooted soils are considered in ref. [40]. In particular, an experimental study was performed to capture a soil water-retention curve and hydraulic conductivity function of unsaturated rooted soils for different temperatures (25, 45, and 60 °C)

Studies of leachate migration are usually carried out in the framework of the isothermal approach. Field studies aimed at a comprehensive analysis of the thermal aspects of MSW landfills depending on operational and climatic conditions were conducted in ref. [41]. In general, it was demonstrated that a significant amount of heat can be released from landfills, which leads to a considerable increase in leachate temperature. In ref. [42], approaches to numerical modeling for temperature prediction in MSW landfill bodies were developed. Heat release rate functions of varying complexity were obtained empirically, which reflect temperature time trends under heat generation as a result of the biological decomposition of waste. The obtained results are consistent with data obtained from MSW landfills located in Michigan, USA [43].

Despite the variety of waste management techniques, the most common method is still landfill disposal. Landfill sites can be divided into two main categories. The first category is landfills built according to modern environmental standards and regulations. The second category is older operating landfills, which were built in the 1980s and 1990s. These older landfills can be a significant source of environmental pollution due to the potential for leachate to seep into the groundwater and surface water near the location of objects. Solid municipal waste at disposal sites is subject to decomposition processes, and as a result of which, organic and inorganic soluble compounds are transformed into a solution in the presence of precipitation. The leachate formed in the body of the landfill is especially dangerous for the environment, since it is a toxic solution with a high concentration of ammonium ions, chlorine, and other components reaching several tens of grams per liter. The negative impact of the leachate formed in the landfill body on the components of the natural environment is its penetration into groundwater and subsequently into surface water. Drainage systems are organized at landfills to remove the leachate. However, over time, changes in the foundation of the landfill can lead to the formation of depression areas where the leachate accumulates without being able to flow into the drainage system. In the current study, the effect of heating in the body of the landfill on the process of groundwater pollution due to potential seepage of leachate because of its accumulation in the relief depressions of the MSW landfill is numerically investigated. In the framework of analysis, we consider a landfill that has been in operation for more than forty years. The modeling was performed paying special attention to the development of Rayleigh–Taylor instability, which is an instability of the interface between two fluids of different densities. It occurs when a less dense fluid is placed beneath a more dense one and causes the appearance of plumes of the lighter fluid to push upward and the heavier fluid downward, which could potentially increase the spreading rate of the leachate.

2. Characteristics of the Solid Municipal Waste Landfill Under Consideration

On the territory of the landfill under consideration, field measurements have been regularly carried out within the framework of the program of environmental monitoring of the quality of atmospheric air, soil, ground and surface water, and radiation. The landfill is located on a low slope near a small river. The terrain is uneven, with a general inclination towards the northeast, towards the floodplain of the river, which flows 300 m to the northeast of the site. The adjacent area is a hilly region of the watershed, covered with meadows, shrubs, and small trees. Leachate from the landfill flows along the impermeable base and collects in two depressions. Along the perimeter of the landfill, there is a dam up to 6 m tall. The area of the landfill is approximately $\mathrm{10,000} {\mathrm{m}}^2$. The landfill was constructed in the 1970s using the natural layer of clay loam as insulation. The layer spreads evenly over the landfill base with a thickness between 0.5 and 1.5 m. The drainage system for leachate removal consists of a ring drainage trench and seepage water retention ponds. It is worth noting that the landfill complex is located on a hilly area, and due to the fact that the facility has been operating for over forty years, there is uneven terrain within the area caused by deformations in the geological layers during operation. These changes can lead to the formation of depressions under the landfill material, where waste liquid can accumulate without being able to flow down to the drainage system. These depressions can reach depths of 4–5 m (see Figure 1 and Figure 2). The leachate in these depressions can be accumulated in a relatively short period of time, typically several dozen days, due to precipitation or heavy snow melt, as is typical for this region during the spring season.

Under these conditions, a significant increase in temperature is observed in the waste layer above depressions, which is caused by the aerobic and anaerobic decomposition of waste. This raises the question of the possible diffusion of leachate into the soil at the base of the landfill and its further penetration into groundwater in the presence of elevated temperatures caused by self-heating in the waste mass. The rate of groundwater filtration in this area depends on the slope and permeability of the rocks and can range from 3 to 30 m per day. Thus, if leachate does penetrate into groundwater, it can spread over tens of kilometers in a relatively short time.

The landfill body can be considered as a layer of technogenic deposits represented by bulk heterogeneous non-compacted and compacted substrates. It consists of various types of waste, including household waste, industrial waste, and construction waste. In particular, polyethylene takes 20–25% of the composition, followed by construction waste (crushed stone, pebbles, concrete, and brick) at 15–30%, textile at 10–20%, waste paper at 8–20%, plastic at 9–15%, wood at 8–20%, glass at 7–15%, metal (in the form of cables, wires, and reinforcement) at 7–15%, and mineral soil at 10%. Figure 2 shows the geological section of a landfill site. The geological composition of the site includes three layers of porous media up to the groundwater level. The insulation layer is a clay loam substance having a thickness of 4.3 m. Below the clay loam layer, there is the siltstone bed with a maximum thickness of 1.5 m. The bottom layer is mudstone of approximately 4 m thickness. The permeability and porosity of the media were taken from the archived data of field measurements in 2017 [44,45] and are presented in

Table 1. Filtration properties of soils in the base of MSW landfill.

Medium	Schematic Representation in Figure 2	Filtration Coefficient, ${\mathit{K}}_{\mathit{f}}$, m/day	Porosity	Permeability Coefficient, K, ${\mathbf{m}}^2$
Waste		0.6	0.4	$5.87$
Clay loam		0.0065	0.702	$6.38·{10}^{-8}$
Siltstone		0.481	0.481	$7.85·{10}^{-7}$
Mudstone		0.542	0.542	$2.35·{10}^{-7}$

.

Since depressions in the landfill foundation can have different sizes, a general problem was considered for a part of the depression taken as an excavation of a rectangular section of a part of the waste layer and soil layers of the landfill base (Figure 3). Two cases of leachate distribution in the thickness of the solid municipal waste landfill base were considered: isothermal and with a significant change in temperature in the landfill body under conditions of self-heating during aerobic and anaerobic decomposition of waste. A section of the porous medium, the upper layer (1), in which the leachate is located, was studied. This is a waste layer represented by bulk soils that are heterogeneous, uncompacted and compacted, and loose, the so-called garbage sole. Down the section there is the clay loam layer (2), which acts as a screen preventing leachate from penetrating under the waste disposal massif. The third (3) and fourth (4) layers are natural layers characteristic of the given area, a siltstone layer and a mudstone layer. The mudstone layer contains unconfined groundwater, and its thickness to the impermeable horizon is 3 m. All other layers have a thickness of 1 m. The area considered is 10 m by 3 m in the horizontal direction and 6 m in the vertical direction. In the initial state, the leachate fills the waste layer. In the first version of the calculations, the development of the Rayleigh–Taylor instability in the isothermal case is studied, and in the second version of the calculations, the effect of temperature heterogeneity on the distribution of the leachate in the soil environment is considered.

3. Mathematical Model of Leachate Propagation

Three-dimensional numerical modeling of leachate propagation from an undrained depression of a landfill was carried out using a laminar flow model in a porous medium. The problem was solved within the framework of non-stationary isothermal and non-isothermal approaches. To solve the problem, a uniform spatial grid consisting of four-sided elements was generated.

The equations of motion in the layer in tensor form have the form:

$${\displaystyle \frac{\partial \left(m\rho \right)}{\partial t}}+{\displaystyle \frac{\partial }{\partial x_i}}\left(\rho {\mathrm{v}}_i\right)=0,$$

(1)

$${\displaystyle \frac{\partial }{\partial t}}\left({\displaystyle \frac{1}{m}}\rho {\mathrm{v}}_i\right)+{\displaystyle \frac{\partial }{\partial x_j}}\left({\displaystyle \frac{1}{m^2}}\rho {\mathrm{v}}_i{\mathrm{v}}_j\right)=-{\displaystyle \frac{\partial p}{\partial x_i}}+{\displaystyle \frac{\partial }{\partial x_j}}\left[{\displaystyle \frac{\rho \mu }{m}}\left({\displaystyle \frac{\partial {\mathrm{v}}_i}{\partial x_j}}\right)\right]-\rho {\displaystyle \frac{\mu }{K}}{\mathrm{v}}_i+\rho g_i,$$

(2)

where t is time, m and K are the porosity and permeability of the porous medium, $\rho$ is the density of the liquid, ${\mathrm{v}}_i$ are components of the velocity vector ($i=1,2,3$), $\mu$ is kinematic viscosity of the liquid, $g_i$ is the acceleration due to gravity.

The dependence of the liquid density on the impurity concentration c and temperature was assumed to be linear:

$$\rho ={\rho }_0+{\beta }_{c}C-{\beta }_{T}T,$$

(3)

where ${\rho }_0=999.993 {\mathrm{kg}/\mathrm{m}}^3$, C is impurity concentration, ${\beta }_c$ and ${\beta }_T$ are coefficients of concentration and volume expansion.

The equation of impurity transport is

$${\displaystyle \frac{\partial }{\partial t}}\left(m\rho C\right)+{\displaystyle \frac{\partial }{\partial x_i}}\left(\rho {\mathrm{v}}_{i}C\right)=-{\displaystyle \frac{\partial }{\partial x_i}}{J}_i$$

(4)

Here, $\mathbf{J}$ is the diffusion flux of impurity, determined by the expression

$$\mathbf{J}=-\rho D\nabla C,$$

(5)

where D is the given coefficient of molecular diffusion.

The heat transfer equation is

$$\begin{array}{cc}\hfill \left(\rho C_{p}m+{\rho }_s{C_p}_s\left(1-m\right)\right){\displaystyle \frac{\partial T}{\partial t}}+& \rho C_p{\mathrm{v}}_i{\displaystyle \frac{\partial }{\partial x_i}}T=\hfill \\ \hfill =& \left(\kappa m+{\kappa }_s\left(1-m\right)\right){\displaystyle \frac{\partial }{\partial x_i}}\left({\displaystyle \frac{\partial T}{\partial x_i}}\right).\hfill \end{array}$$

(6)

Here, $C_p$ is the specific heat capacity of the liquid, ${\rho }_s$ is the density of the porous medium, $C_{ps}$ is the specific heat capacity of the porous medium, $\kappa$ is the thermal conductivity coefficient of the liquid, ${\kappa }_s$ is the thermal conductivity coefficient of the porous medium.

The lower boundary of the calculation area is the first aquifer and is assumed to be impermeable and solid, and conditions of adhesion and impermeability for the substance were set on it. The surface of the liquid in the waste layer was assumed to be free of tangential stresses. On the lateral boundaries, conditions of continuation along the flow were set so that the leachate migrates without a jump beyond the limits of the calculation area (conditions of continuity of the liquid flow). Conditions of equality of mass and heat flows were set on the boundaries between the layers. Boundary conditions are listed in

Table 2. Boundary conditions.

Boundary	Flow Conditions	Concentration Conditions	Thermal
Domain top (landfill body)	Zero shear stresses	$C=C_0$ (300 g/L)	$T=T_h$ (from 30 °C to 70 °C)
Domain bottom (groundwater)	Zero shear stresses	Zero flux	$T=T_c$ (10 °C)
Side boundaries	Zero normal stress	Zero flux	Zero flux
Interlayer boundaries	Equality of ${\mathrm{v}}_{\mathrm{z}}$	Mass balance	Equality of T

.

Three-dimensional numerical modeling of leachate propagation was carried out based on the geometric model shown in Figure 3. It is assumed that the landfill body (the waste layer) is water-saturated due to precipitation throughout the entire calculation time. The behavior of leachate with concentration of impurities $C_0=300$ g/L is considered, which corresponds to the maximum values of dry residue detected during field measurements [44]. The physical parameters of the leachate and porous medium layers under study are presented in

Table 3. Physical parameters of leachate.

Parameter	Value
Thermal conductivity coefficient ($\kappa$), W/(mK)	0.56
Specific heat capacity ($C_p$), J/(kg K)	4180
Thermal expansion coefficient (${\beta }_T$), 1/K	$1.8·{10}^{-4}$
Concentration expansion coefficient (${\beta }_c$)	0.96
Molecular diffusion coefficient (D), ${\mathrm{m}}^2$/s	$1.5·{10}^{-9}$

and

Table 4. Physical parameters of porous media.

Medium	Thermal Conductivity Coefficient, ${\mathit{\kappa }}_{\mathit{s}}$, W/(mK)	Specific Heat Capacity, ${\mathit{C}}_{\mathit{ps}}$, J/(kg K)	Density, ${\mathit{\rho }}_{\mathit{s}}$, ${\mathrm{kg}/\mathbf{m}}^3$
Waste	1.06	1180	$1.25·{10}^3$
Clay loam	1.36	1973	$1.55·{10}^3$
Siltstone	2.1	894	$1.74·{10}^3$
Mudstone	1.69	756	$1.67·{10}^3$

.

The slope of the aquiclude surface in the area of the landfill is assumed to be zero. The initial concentration of the impurity is uniform across all layers of the landfill foundation and equal to the ambient concentration that corresponds to $C=0$ g/L. The field of initial velocity is unperturbed and eventually can be taken equal to zero.

Two cases of leachate distribution in the thickness of the MSW landfill base were considered: isothermal and with the presence of a significant change in temperature in the landfill body under conditions of self-heating of waste. Since the viscosity of the leachate, as an aqueous solution, depends on temperature, then at the initial stage the behavior of the leachate, characterized by different viscosity values corresponding to the limiting values of viscosities that can characterize the leachate, was studied. In this case, viscosity can be determined not only by temperature inhomogeneities but also by the composition of the leachate [46]. At the first stage of our study, we focus on viscosity effect caused by leachate composition in the isothermal approximation. Then, the influence of temperature impact is considered. In this case, following the results of work [46] (Figure 4), we assume that the viscosity decreases with increasing temperature according to the formula:

$$\begin{array}{cc}\hfill \mu =& 1.7533-5.20095·{10}^{-2}T+8.7152·{10}^{-4}{T}^2-\hfill \\ \hfill & -7.57098·{10}^{-6}{T}^3+2.5877·{10}^{-8}{T}^4,\hfill \end{array}$$

(7)

where T is leachate temperature in Celsius. The formula agrees well with experimental data for temperatures between 5 °C and 100 °C (relative error less than 0.1%).

The calculation grid is made up of rectangular elements. Three different grid sizes were considered to determine the optimal spatial step (see Figure 5). The change in the average concentration for the isothermal case shows that the results for a grid with a step of 0.075 m differ from the results obtained on grids with a step of 0.05 m and 0.025 m by more than 5%. On grids with a step of 0.05 m and 0.025, the difference in the maximum concentration value is less than 3%; therefore, a grid with a step of 0.05 m was used to conduct the calculation experiments.

4. Results of Numerical Simulations

4.1. Isothermal Case

In this section, we considered the effect of leachate viscosity in the isothermal case on its distribution in the landfill foundation. As mentioned above, significant formation of leachate over several days can be caused by snowmelt in the spring or heavy precipitation in the summer and autumn. As a result of the accumulation of leachate in the depression of the landfill without the possibility of its movement through the drainage system, a condition arises under which the leachate can accumulate and the level of filtration waters can reach several meters. The case of a maximally saturated aqueous solution (impurity content of 300 g/L) was studied in this work.

We assume that in the initial state, the leachate fills the waste layer so that the filling depth is $h_1=1$ m (Figure 3). In the first series of calculations (C1), the leachate is characterized by the viscosity of $1.31·{10}^{-3}$ Pa·s (water at 10 °C); in the second series of the calculations (C2), the influence of the effect of reducing the viscosity of the leachate on its dynamics is considered taking viscosity equal to $1.04·{10}^{-3}$ Pa·s.

The concentration fields in the vertical sections in the middle of the calculation region are shown in Figure 6. The time dynamics of the spread of the leachate, characterized by a viscosity of $1.31·{10}^{-3}$ Pa·s, is presented as a result of rapid accumulation of the leachate at the site of the landfill depression. Over 12 years of leachate presence in the depression location, the impurity concentration will increase in layer 4; i.e., the leachate will begin to enter the groundwater. A detailed examination of the development of instability for both cases is shown in Figure 7, Figure 8 and Figure 9.

For a detailed examination of the influence of instability development on the process of filtration water distribution in the system of soil layers, the concentration fields of the impurity in the vertical sections in the middle of the calculated section are shown below.

The initial distribution of concentration is almost unchanged during the first days in both cases, for the case of higher viscosity (Figure 7a) and for the less viscous leachate (Figure 7b). However, after 16 days, in the case of the leachate characterized by lower viscosity, the diffusion spreading in the body of the polygon is supplemented by the development of Rayleigh–Taylor instability, the pattern of which is determined by finger-like structures, clearly distinguishable for 20 days from the start of the calculations, while for this time in the more viscous case, the concentration front spreads in a diffusion manner.

The appearance of finger-like structures is also observed in the more viscous case, a month later than the appearance of the wave-like spread of the concentration front in the less viscous case. For both cases, at the initial stage of the Rayleigh–Taylor instability development, finger-like structures are directed into the landfill body since the permeability of the waste layer is two orders of magnitude greater than the permeability of the clay loam layer, which acts as a screen on the landfill. If leachate formation is absent at the landfill, such a process of pure water spreading in the waste layer will lead to dilution of the leachate, eliminating the need for disposal and storage of the leachate. However, as a rule, during landfill operation, the process of leachate formation is continuous throughout the entire period of operation of the MSW landfill. This effect is shown in Figure 7, containing concentration fields for both cases after 4 months and 8 months from the start of the calculations. Finger-like structures develop intensively in the direction of the waste layer, and the period of these structures at the initial stage is determined by the size of the grid; after a while, the distribution pattern for grids of different sizes is similar, and the spatial period is about half a meter. For the isothermal case, finger-like structures lag behind in height for six months. After 8 months, the pattern of concentration torches becomes almost indistinguishable for both cases (Figure 7c,d).

In the presence of a constant source of leachate inflow, an increase in the impurity concentration is also observed in the clay loam layer; while the movement of liquid in the clay loam layer over a year and a half on a macro-scale is not realized, diffusion spreading of the concentration front with dilution occurs. The value of the pollutant concentration is sufficient to achieve the permissible value of the pollutant concentration for the environment (Figure 8). In this case, in the non-isothermal case, after a year and six months, the pattern of the concentration front propagation in the insulating layer acquires a wave character. Rayleigh–Taylor instability develops in the clay loam layer.

In the isothermal case C1, the Rayleigh–Taylor instability does not substantially affect the spread of the leachate over a period of 5 years, as can be seen from Figure 8. However, in the case of a decrease in the viscosity of the leachate formed in the body of the landfill, the process accelerates and the leachate reaches the groundwater layer (mudstone layer–4) over a period of 4 years. The instability development process is of a random spatial nature (Figure 8d), having a kind of “breakdown” place, in which the speed of leachate spread increases almost by an order of magnitude.

The estimated service life of the considered MSW landfill is between fifty and seventy years, depending on the potential capacity of the massif and the rate of its filling. After 18 years of leachate propagation in the case of higher viscosity, Rayleigh–Taylor instability leads to a breakdown of the siltstone layer and intensive entry of pollutants into groundwater (Figure 9). After 30 years, the patterns of leachate propagation in both cases become identical. The spatial scale of instability structures increases by an order of magnitude, such that the maximum amount of pollutants that determines the composition of the filtered water enters the groundwater.

The dynamics of the leachate concentration along the vertical axis is shown in Figure 10. The concentration profiles in the figure are plotted for coordinates, where the beginning of instability development occurs, and the maximum deviation of the concentration front is observed (so-called direction of maximum growth). At the initial stage of instability development, the disturbance of concentration fronts is symmetrical. However, the perturbations lose their symmetry with respect to this direction over time, and a local minimum appears in Figure 10b, located at the base of the finger-shaped convective structure. This also explains the sharp change in the 9-year curve in the range of 0.5 to 2.5 m.

Changes in the maximum and average concentration values in the layers over time are shown in Figure 11. Figure 11a corresponds to the case of higher viscosity; Figure 11b corresponds to the case of less viscous leachate. Continuous distribution of the leachate under the most unfavorable conditions reaches saturation after 30 years. This time is obtained based on the maximum water saturation. The amount of precipitation determines the amount of leachate, including the amount of leachate that can accumulate in a depression in the landfill area. Within 100 years, clay loam will perform the screening function by 20% since the maximum value will not be reached. An increase in the amount of impurity, as a rule, leads to an increase in the viscosity of the solution, which should lead to a slower spread of the leachate in the soil layer.

4.2. The Influence of Heat Emission in the Body of the Polygon

The temperature change in the landfill body depends on the waste content, the period of waste residence in the landfill, the method of waste conservation, and climatic conditions [41]. Thus, in the work [41], it is shown that the maximum temperature of waste newly deposited in the landfill body varies from 40 to 80 °C. The landfill we are studying is an example of the processes that occur in functioning landfills. Since the area of the landfill object is more than 1 ha and the storage zones are constantly replenished, the processes of elevated temperature are observed throughout the entire service life of the landfill. Cases are considered when, due to heat release in the landfill body, the temperature at the upper boundary at the level of the leachate occurrence is equal to the heating temperature $T_h$, and at the lower boundary the temperature is also constant and equal to $T_c$. In this case, the viscosity depends on the temperature in accordance with the data of the work [47].

We studied heating in the body of the landfill, leading to temperatures from 40 °C to 70 °C. The groundwater layer obviously can be assumed isothermal with temperature equal to 10 °C. Thus, the gradient difference for each case consisted of the temperature gradient. Since the disturbances develop randomly, for the case of increasing the temperature in the body of the landfill to 70 °C, the concentration fields of the impurity are shown in the section for which the greatest increase in disturbances is observed (Figure 12). For the limiting case of maximum temperature, after 8 years, the leachate begins to penetrate into the groundwater.

Changes in the concentration fields in the vertical section in the middle of the calculation region for the isothermal case and for cases of increasing temperature to 40 °C and 70 °C are shown in Figure 13 and Figure 14. At the initial stage of instability development, finger-shaped structures appear, the spatial period of which is infinitely small (see Figure 14). Over time, these structures become larger so that the spatial period becomes half the thickness of the considered layer system. The rate of leachate flow increases with an increase in temperature in the landfill body, almost doubling. Thus, considering heat release processes is essential for forecasting and monitoring the state of a solid municipal waste landfill, especially in cases of depressions caused by geological changes at the landfill site. As can be seen from Figure 13a,b, the development of instability determines the rate of leachate propagation.Thus, when the temperature is heated to 40 °C and 70 °C, the concentration profile remains the same over a period of 5 years (shown by the black lines in Figure 13a,b). For 9 years, if the temperature is increased to 70 °C, the impurity concentration in groundwater reaches 180 g/L, whereas if the temperature remains at 40 °C, the impurity concentration is 20 g/L. The curves for the gray and blue lines are influenced by the asymmetry of the finger-like structure, and in the center of the groundwater layer, there are areas with lower impurity concentrations.

In the waste layer, there is an equivalent change in concentration, with the onset of instability. A sharp decrease in impurity concentration is observed. Over time, the concentration tends to accumulate again. After that, ten years later, the dilution rate of concentration also increases, indicating an oscillatory process. In the groundwater layer, an increase in the amount of infiltration water has been observed with each change in temperature in the landfill body for 25 years (see Figure 15). At a temperature of $T_h=30$ °C, the amount of impurity penetrates 1.9 times more in 20 years compared to the isothermal case. When the temperature rises to 70 °C for 20 years, the amount of impurity increases by five times (see dashed line in Figure 15) Thus, heat release in the landfill body significantly affects the spread of pollution.

5. Discussion

Similar to leachate spreading at the base of a landfill caused by the destruction of the natural seal, in geothermal systems, a transient boundary layer is formed by the diffusion of dissolved substances. For instance, fertilizers or pollutants may dissolve in water at the ground surface, increasing its density. When convection starts, the solute is quickly transported downward. A similar process can occur in groundwater beneath salt lakes, where evaporation from the surface increases the density of brine. Brine migration due to convection has implications for the possible use of salt lakes as disposal sites for pumped saline groundwater [48]. Another application has arisen from proposals for large-scale underground storage of carbon dioxide to reduce atmospheric emissions and thus limit the impact of hydrocarbon use on the global climate. Typical storage sites are located at depths greater than 1000 m, and in such subsurface settings the density of the carbon dioxide-rich phase is approximately one-half to two-thirds that of the formation water. Following injection into permeable rock beneath a suitable low-permeability seal, some carbon dioxide will rise due to buoyancy and accumulate beneath the seal. In the process, the carbon dioxide dissolves in the formation water (typical solubility is 2–5% by weight, depending on salinity). Unusually for a gas, the dissolved carbon dioxide increases the density of the fluid, and thus the system becomes unstably stratified [49].

The onset of convection significantly accelerates further dissolution of the carbon dioxide and is important for assessing the safety of storage over hundreds to thousands of years. A complication in natural porous environments is that in many cases the solute can react with minerals in the rock, changing both the permeability and the density of the fluid. These changes can either reduce or enhance convection, depending on whether the geochemical reactions lead to pure precipitation or pure dissolution. Such a relationship goes far beyond the standard problem considered in the present work. Such assessments are necessary, and they are carried out within the framework of in situ measurements [50]. In refs. [48,51], a sudden process of the emergence of a complex plume formed near the leading edge of a suddenly saline semi-infinite surface [48,51] or near the leading edge of a suddenly heated semi-infinite surface [51] is observed. In our work, a process in which highly concentrated brine accumulates in a polygon depression is investigated. The process of leachate accumulation is determined by times of the order of half a year to a year, whereas the process of instability development is determined by scales of more than five years. Therefore, in the present process the suddenness described in refs. [48,51] is not observed. However, the nature of the instability itself is the same. Thus, if the accumulation process is not sudden, instability may develop over time, which leads to an acceleration of the spread of the concentration front.

6. Conclusions

The process of leachate migration is considered, taking into account temperature heterogeneity caused by aerobic and anaerobic decomposition of waste. Calculations show that at the initial stage, leachate spreads slowly by diffusion only. This behavior is observed for periods of several months in both isothermal conditions and in the presence of temperature gradients caused by self-heating effects within the landfill. After six months, an increase in temperature leads to a finger-like pattern of lower concentration moving upwards through the waste layer of the landfill, while the bottom concentration front continues spreading downward by diffusion (see Figure 15). The upward motion can be treated as a positive effect because it decreases the mean impurity concentration in the waste layer. In all calculations, the downward motion results in a picture typical to Rayleigh–Taylor instability after several years from the leakage starting. However, there is a substantial impact of leachate self-heating causing reducing the liquid viscosity on the migration rate. Calculations show that if there is a temperature difference of 70 °C, the instability arises two times earlier and develops more intensely over the time, which causes critical growth (up to five times) of mean concentration at groundwater level (see concentration after 20 years in Figure 14).

The leachate’s self-heating affects the rate of pollution spreading by two mechanisms. Firstly, an increase in temperature is observed in the body of the landfill, which is heating from the above effect for foundation soils. As a result, the heating suppresses natural convection in the porous medium under the landfill and reduces leachate migration. On the other hand, leachate viscosity decreases with temperature growth, which reduces the threshold of Rayleigh–Taylor instability and intensifies the convection flows associated with it. The calculations show that for the parameters of natural soils, the second mechanism prevails, and, therefore, internal heating intensifies leachate migration. Thus, the self-heating effect has to be taken into account in the assessment of the leachate migration rate.

It is worth noting that the ambient temperature can have an impact on the processes occurring within geological layers. In particular, daily fluctuations in air temperature affect temperature changes at a depth of 1–1.5 m. This is due to the transfer of solar heat flux through the molecular thermal conductivity of rocks and the convection of air and water vapor, including infiltrating precipitation and groundwater. Annual fluctuations can cause temperature changes up to 20–40 m deep. At these depths, heat is transferred mainly through molecular thermal conductivity and the motion of groundwater. Below 20–40 m, there is a neutral zone, where the temperature remains stable. The landfill we are considering has a waste thickness of more than 30 m above the geological layers of the landfill foundation. Therefore, seasonal temperature variations do not affect the temperatures in the porous media we are concerned with, and, thus, they are not taken into account in the modeling.