Numerical solution of the Burgers equation associated with the phenomena of longitudinal dispersion depending on time.

In this study, the Burgers equation governing the time-dependent dispersion phenomena is solved numerically using the finite difference scheme and the Runge-Kutta 4 algorithm with appropriate initial and boundary conditions. Two time-dependent dispersion functions have been implemented to analyze the spatio-temporal variation in the domain. For the values of KL and KA < 1.2 years, a significant retention of the mass of solute is observed when the dispersion function is asymptotic. The results obtained show that the concentration profiles are similar when the values of KL and KA ≥ 1.2 years. These results demonstrate the importance of the nature of the dispersion function on the retention capacity of solutes in the porous medium.

Introduction

The fate of properties such as salinity in underground media is of great interest due to increasing concern about the deteriorating environment and human health caused by poor quality groundwater (Ghafari et al., 2020). The transport of dissolved contaminants by heterogeneous hydrogeology is based on expressions of their functional parameters. Several processes acted simultaneously on the chemical constituents transported through the soil. This requires that quantitative descriptions of chemical transport include feasible processes such as the intensity of uptake of contaminants from the subterranean environment. Predicting the concentration of pollutants in the soil matrix is very important to minimize the risk and vulnerability of aquifers (Xie et al., 2019). In addition, the accuracy of the predicted model is crucial to adequately assess and predict the behavior of contaminant transport in the subsoil. Solute transport in groundwater systems is traditionally modeled by the classical non-linear advection-dispersion equation which can add to equilibrium uptake and first order degradation (Ding et al., 2021). The common hypothesis of the dispersion of pollutants in aquifers with constant porosity, regular infiltration flow rate and dispersion coefficient has been approached by many researchers from different angles. Here are some examples. Xue et al. (2020) relied on the new algorithm obtained by the precision of the Crank-Nicolson scheme to solve two-dimensional parabolic equations by alternating the implicit technique of direction. Kumar and Unny (1977) focussed his research on the application of Runge-Kutta methods to solve nonlinear partial differential equations using a specific fluid flow problem. Runge-Kutta of order 4 (RK4), remains up to date a very powerful tool for solving ordinary and partial differential equations. It has the advantage of being simple to program and quite stable for the current functions of physics (Kumar and Unny, 1977; Barletta et al., 2020). Most research has been directed towards improving the precision or flexibility of the classic method from Runge-Kutta methods.

Vázquez-Suñé et al. (2004) considered dispersivity as the most important and difficult parameter to assess in seawater intrusion models, in order to assess the effect of spatial variability of hydraulic conductivity on effective dispersion in seawater intrusion problems. Jaiswal et al. (2011) relied on an analytical solution of the advection-dispersion equation of an input concentration constant along an unstable horizontal flow in a semi-infinite shallow aquifer to determine the level of pollutants in groundwater. Corey et al. (1970), Warrick et al. (1971) have conducted field and experimental studies, and they have suggested that the coefficient of dispersion is not constant but increasing with time.

An unconditionally stable Crank-Nicolson finite difference scheme was used by Ravi et al. (2014) to analyze the constant and longitudinal dispersion profile of contaminants. In his studies, Ravi et al. (2014) did not take into account the time-dependent dispersions coefficients, the decay rate constant and the zero-order production rate coefficient of solute in the liquid phase as done by Das et al. (2017), Guleria et al. (2020) with a linear dispersion advection equation model. Our objective is to exploit the more stable fourth-order Runge-Kutta method (RK4) to evaluate the profile of C (X, T) contaminants through a porous medium and to conduct a comparative study between the profiles of the contaminants obtained from an asymptotic and linear dispersion coefficient, and on the other hand, from the initial and boundary conditions used in the work of Das et al. (2017). The discussion will focus on determining the dispersion parameters related to the concentration profiles using the two time-dependent dispersion coefficients.

Materials and methods

Physical model

Water containing pollutants such as sewage eventually seeps into the aquifer which is a reservoir of drinking water. As water passes through the soil, pollutants are mixed, adsorbed, dispersed and diffused by the flowing stream Fig. 1. Several efforts are being made by the scientific community to develop more precise and economical models capable of predicting the transport and concentration of solutes in unsaturated soils. The transport and mixing of contaminants is governed by the advection and dispersion equation of pollutants in the terrestrial aquifer Xie et al. (2019).

Figure 1

Infiltration and dispersive system of pollutants in the porous medium.

Physical description of the time-dependent dispersion coefficient

Dispersion is a physical phenomenon that occurs when moving away from the injection site, the mass of solute dilutes with time to occupy an increasing volume with a correspondingly decreasing concentration see Fig. 2.

Figure 2

Geometry of the time-dependent dispersion problem.

Mathematical model

Over the past decades, a large number of analytical and numerical solutions have been developed to estimate the fate and transport of various constituents in the underground environment. The application of these solutions is generally limited to non-variable subterranean flow fields and relatively simple initial and boundary conditions nevertheless; these solutions play an important role in contaminant transport studies, providing initial or rough estimates of the distribution of solute concentrations in soils and aquifer systems. The convection-dispersion equation has remained the basis of most analytical and numerical studies of solute transport. The dispersion equation which describes the distribution of the concentration of miscible fluids (i.e. water contaminated or salted with fresh water) in heterogeneous underground media remains one-dimensional in the direction of flow fluids along the x axis (ox) and the other oy and oz components become negligible. The transport equation can be written as equation (1) (Kajal et al., 2014; Singh et al., 2015): Where C is the concentration of the contaminated fluid, F is the concentration in the solid phase, ϕ the porosity of the medium, with DL the longitudinal dispersion coefficient on the macroscopic scale, and () the component of the infiltration rate of salt water along the x-axis, and L is the length of dispersion in the direction of flow. Two cases were considered by Lapidus and Amundson (1952) for the solid phase concentration and its derivative, which are equations (2) and (3): Here there is an equilibrium and unbalanced relationship between the concentrations in the two phases. and are the constants generally referred to as the order factor for the distribution of pollutants in aquifers. By combining the equations (3) and (2) in (1), we obtain equation (4): Where R is the delay factor describing the adsorption of solutes in the porous medium, its expression is equation (5): The component of the infiltration rate is related to the dispersal concentration of the pollutant. Ravi et al. (2014) assumes that the rate of infiltration is inversely proportional to the concentration of the pollutant as given in equation (6): This infiltration rate is the cause of the nonlinearity in the advection – dispersion equation. represents the proportionality concentration of salt water dispersion. The new parameters of the independent variables introduced to simplify equation (4) are defined as and .

Therefore, equation (4) reduces to equation (7):

Description of the time-dependent dispersion model

Field and experimental evidence from studies has been suggested by Corey et al. (1970), Warrick et al. (1971), Sauty (1980), Pickens and Grisak (1981), Sposito et al. (1986) that the dispersion coefficient is not constant but apparently increasing as a function of the displacements in time or equivalently with the distance of displacement of the solute. The apparent increase in the dispersion coefficient has been called the scale effect (Fried and Combarnous, 1971). The scale effect is generally attributed to the heterogeneity of the porous formation, in particular in the heterogeneity of the hydraulic conductivity. Stochastic analyzes have shown that the dispersion depends on the transport time and increases until it reaches an asymptotic value (Gelhar et al., 1979). The theoretical deterministic analysis of Bourgeois (1984) also established that the dispersivity in a stratified aquifer is time dependent. Therefore, there is ample evidence that the dependency scale causes dispersion to vary over time. Therefore, a time-dependent dispersion model can be used to provide a rough description of the transport scale in our study by exploiting the Burger equation. In this study, a linear and asymptotic form of the time-dependent dispersion coefficient is integrated into the nonlinear PDE equation (7) to study the phenomena of longitudinal dispersion of solutes in the underground environment; this has been mentioned in various works (Guleria and Swami, 2018; Guleria et al., 2020). This type of dispersion will make it possible to better represent the results for a transport of solutes in a given porous medium. The time-dependent dispersion coefficient is studied in two forms in this work and equation (7) becomes equation (8):

Time dependent linear dispersion function

The time-dependent linear dispersion increases with time without boundary conditions in the porous medium (Gao et al., 2010; Guleria et al., 2020). The linear formula for the time dependent dispersion expression in equation (8) is given by equation (9) below:

Asymptotic dispersion function depend on time

The time-dependent asymptotic distance initially increases with time and eventually approaches an asymptotic value in the porous medium (Guleria et al., 2020). The formula for time-dependent asymptotic dispersion expression in equation (8) is given by equation (10) below: Where is the maximum dispersion coefficient for the asymptotic time-dependent and uniform dispersion function for the linear time-dependent dispersion coefficient; is the effective diffusion coefficient of time; is the time-dependent asymptotic coefficient equivalent to the mean distance travelled by the pollutants in the aquifer; is the linear coefficient depending on time. The values of the parameters of the time-dependent dispersion coefficient were determined by a sensitivity analysis of Barry and Garrison (1989). The same values of the parameters of the linear and asymptotic time-dependent dispersion coefficient exploited in the work of Guleria et al. (2020) which are also exploited in this work, are presented in the Table 1.

Table 1

Parameters of the time-dependent dispersion coefficient.

Parameters	Values	Parameters	Values
D0(L)	5.89 cm2/min	D0(A)	4.51 cm2/min
KL	4500 min	KA	200 min
Dm	0 cm2/min	Dm	0 cm2/min

Fig. 3 illustrates the analysis of the asymptotic limit of the time-dependent dispersion coefficient given in relation 10, which remains important and more practical in porous media. The values of KA depend on the extent of the pre-asymptotic zone. KA equal to zero indicating constant dispersion. The smaller the value of KA, the more the dispersion approaches the asymptotic value. The values of KA other than zero corresponding to the times for which the dispersion reaches half of its asymptotic value.

Figure 3

Asymptotic dispersion coefficient for various KA using data from Table 1.

Initial and boundary conditions

From the source of the pollutants, the concentration disperses more in the porous medium in the direction of flow with a front at the arrival of the boundary see Fig. 4. At the initial moment it is assumed that the aquifer is contaminated, a certain initial background concentration exists in the aquifer and it is represented by a linear combination of an initial concentration and the term of zero order production with rapid infiltration given by equation (11): Where is the initial background concentration, is flow velocity, and γ is the zero order production rate coefficient for liquid phase solute production. A contaminant in radioactive waste decaying exponentially with time is imposed upon entering the aquifer as a linear combination of a source concentration with an initial background concentration at the origin, to describe the transport of solutes in a natural or artificial system as expressed by equation (12) (Das et al., 2017; Singh et al., 2021): λ is the decay rate constant and is the flow velocity of fluids in the porous medium.

Figure 4

Geometry of the problem.

At the other end of the aquifer, solute transport may not be affected and therefore an exit boundary condition is prescribed as a non-flow boundary condition (Das et al., 2017). The mathematical expression of these phenomena is expressed by equation (13) as:

Numerical solution of the mathematical model

The one-dimensional solute transport equation is a parabolic-type partial differential equation in which the finite difference technique is commonly used to obtain the numerical solution (Togue-kamga et al., 2012). An explicit finite difference technique is usually used, and ultimately results in restrictive stability criteria Das et al. (2017).

Fourth-order Runge-Kutta (RK4) method

The numerical technique is used to solve the equation (8) expressed as equation (14) (Togue et al., 2019; Kumar and Unny, 1977): The semidiscrete equation (14) is solved using the RK4 scheme according to the algorithm given by equation (15): With (ΔT) the time step, C(T) the value of the pollutant concentration at time t and C(T + ΔT) the value of the pollutant concentration at time (T + ΔT).

The derivatives of the equation (8) are approximated by the truncated expansion of the numerical finite difference approximation scheme to determine the first and second order spatial derivatives are obtained from equations (16) and (17) (Natarajan et al., 2020): The first order temporal discretization is given by equation (18) written as: The indices (i) and (j) indicate the nodes of discretization along (X) and (T) respectively. ΔX is the spatial step. Thus, equation (8) can be written in a discrete form such as equation (19): Where the function in equation (19) is expressed as equation (20): , , with (i = 0,1,2,...,Nx) and (j = 0,1,2, ...,Nt).

The discretization of boundary and initial conditions is necessary to apply this method. The discrete equation of the initial condition and boundary conditions associated with equation (11), (12) and (13) is expressed as equation (21):

Stability condition

Equation (8) can be written in the discrete form as equation (22) below: The coefficients α and β defined in equation (22) are expressed in equation (23) and (24): This implies that the condition of stability of equation (22) is expressed as equation (25): For the conditions of stability equation (25), the coefficients α and β defined in equation (22) must be positive as expressed in equations (26) and (27), thus From equation (25), we obtain the stability condition which is expressed as equation (28): These inequalities fix a strict maximum limit to the size of the time step and represent a serious limitation for the centered finite difference diagram.

Results and discussion

The spatial and temporal variation in the concentration

We analyzed the numerical solutions using the following model parameters: = 0.01, λ = 0.001 = 0.01 km/year, (Das et al., 2017; Guleria et al., 2020; Natarajan, 2016). Figs. 5a and 5b illustrate the spatial and temporal variation in the concentration of the contaminant for different values of and .

Figure 5

Spatial representation of the concentration C (X, T) for de T = 0.3 (a) and T = 0.6 (b) according to the different values of KL and KA.

Fig. 5a illustrates the spatial distribution of the concentration of pollutants by considering the different values of and . The spatial concentration profiles decrease exponentially in the aquifer for the values of and (0.2 years, 0.7 years, and 1.2 years). This profile decreases rapidly when the dispersion function is asymptotic, reflecting greater retention of the mass of solute in the medium. These results are similar to those obtained by Guleria et al. (2020). These authors have studied time moments (concentration per unit time) to interpret the behavior of the solute plume in porous media such as porous media laminated hydraulically with a dispersion as a function of time. The different concentrations all decrease together from the entry of the aquifer independently of the dispersion coefficient, to reach horizontal asymptotes at the different positions X = 0.029, 0.054, 0.058 respectively for three values of (0.2 years, 0.7 years, 1.2 years). The percentage of residual concentration associated with these positions is 0.76; 0.30; 0.24. On the other hand, the different horizontal asymptotes are rather obtained at the different positions X = 0.04, 0.055, 0.058 respectively for the values of ( = 0.2 years, 0.7 years, 1.2 years). The associated residual concentration percentage is 0.55; 0.28; 0.24. These positions grow with the increase of and . The percentage of residual solute concentration remains low when the dispersion function is asymptotic. These results confirm that the increase in assigned to asymptotic time-dependent dispersion coefficients significantly reduces the mass of solute in the subterranean medium compared to associate with linear dispersion coefficients. The concentration profiles in the porous medium get closer when the characteristic and of the medium increase. For the values of and years, a significant retention of the mass of solute is observed when the dispersion function is asymptotic. A similarity of the concentration profiles is observed for the value of and years. This similarity of the concentration profile is due to the temporal variation of the effective dispersion coefficient without upper limit in the underground medium (Guleria et al., 2020).

Fig. 5b shows that the retention of the mass of solute gradually decreases compared to that observed in Fig. 5a. Likewise, for the values of and years, a significant retention of the mass of solute is observed when the dispersion function is asymptotic. Also, an analogy of the concentration profiles is observed for the value of and 1.2 years.

Conclusion

In this study, the Burgers equation associated with a linear and asymptotic time-dependent dispersion function was numerically simulated to determine the spatiotemporal variation in concentrations. Analysis of the results shows that the concentration profiles decrease rapidly when the dispersion function is asymptotic, reflecting greater retention of the mass of solute in the medium. Values of KL and KA for which the concentration profiles are similar were determined. These results demonstrate the importance of the nature of the dispersion function on the retention capacity of solutes in the porous medium.

Declarations

Author contribution statement

Calvia Yonti Madie: Performed the experiments; Contributed reagents, materials, analysis tools or data; Wrote the paper.

Fulbert Kamga Togue, Paul Woafo: Conceived and designed the experiments; Analyzed and interpreted the data.

Funding statement

This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.

Data availability statement

Data will be made available on request.

Declaration of interests statement

The authors declare no conflict of interest.

Additional information

No additional information is available for this paper.