This study deals with a two-dimensional (2D) contaminant transport problem subject to depth varying input source in a finite homogeneous groundwater reservoir. A depth varying input source at the upstream boundary is assumed as the location of disposal site of the pollutant from where the contaminant enters the soil medium and ultimately to the groundwater reservoir. At the extreme boundary of the flow site, the concentration gradient of the contaminant is assumed to be zero. Contaminant dispersion is considered along the horizontal and vertical directions of the groundwater flow. The governing transport equation is the advection–dispersion equation (ADE) associated with linear sorption and first-order biological degradation. The ADE is solved analytically by adopting Laplace transform method. Crank–Nicolson scheme is also adopted for the numerical simulation of the modelled problem. In the graphical comparison of the analytical and numerical solutions, the numerical solution follows very closely with the analytical solution. Also, Root Mean Square (RMS) error and CPU run time are obtained to account for the performance of the numerical solution.

  • Two-dimensional contaminant transport is discussed in a finite homogeneous porous medium.

  • Exponentially decaying depth-dependent source is assumed at the upstream of the domain.

  • The governing equation is solved analytically and numerically.

  • The RMS error is observed very less.

  • Verified the input data used in the present study from the existing literature.

Graphical Abstract

Graphical Abstract
Graphical Abstract

Contamination of the groundwater reservoir is a major concern throughout the world. The unmanaged disposal of garbage from domestic and industrial, sewage, chemicals and pesticides has raised the pollution in the groundwater reservoir (Delay et al. 1997; Rausch et al. 2005). A strong reason for the pollution is the continuous population growth over a few decades. At present, natural resources such as water, air etc. are of great concern regarding their pollution. The input sources like septic tank leakage, pesticide distribution over fields, underground oil tank leakage etc., may contaminate the groundwater reservoir and can be studied by mathematical modelling (Barcelona et al. 1990; Sun & Sun 1996; Todd & Mays 2004). Many studies have been carried out regarding groundwater reservoir contamination. Not only field and experimental studies dealing with contaminant transport problems but also mathematical aspects are beneficial to model this type of problem. Transportation of contaminants in the groundwater system is generally studied mathematically using ADE which is driven by Fick's law of diffusion (Fetter et al. 1999) and law of conservation of mass (Bear 1979; Domenico & Schwartz 1990; Zheng & Bennett 2002). By using this concept, many researchers have proposed various models. Groundwater contamination is of great concern in the field of environmental pollution, engineering geology, chemical, soil physics, petroleum engineering etc. Some studies have been published previously regarding groundwater contamination, in which the impacts of various transport parameters, like advection, dispersion etc., on the pollutant concentration distribution were discussed (Ebach & White 1958; Ogata 1970; Bear 1972; van Genuchten 1981; Kehew 2000). Experimental and physical studies of the chemicals, fertilizers and pesticides transport in the groundwater system have been discussed previously and various groundwater remediation phenomena have also been proposed (Ritter & Chirnside 1984; Kourakos et al. 2012).

The contaminant transport phenomenon in the groundwater system has been mainly analysed by advection and dispersion processes. The combined effect of dispersion and advection processes is described by ADE, which is the governing equation of a solute transport model problem (Batu 2005). Many studies have been reported in the past dealing with the groundwater pollution. Liu et al. (1998) proposed a one-dimensional (1D) analytical solution for the steady flow for multilayer porous media. They discussed the solution for the arbitrary inlet and initial condition by applying a generalized integral transform technique. Kumar (1983a, 1983b) developed analytical solutions for adsorbing and non-adsorbing porous media in which flow velocities were considered unsteady. The variation in seepage velocities with time was discussed under the effect of resistance in flow. Furthermore, Guerrero et al. (2013) established a new analytical solution by applying the Duhamel theorem for the solution of advection–dispersion with time-varying boundary conditions. Moreover, Guerrero & Skaggs (2010) and van Genuchten et al. (2013) developed a series of analytical solutions for the one-, two- and three-dimensional solute transport equations subject to the various types of inlet and outlet boundary conditions. Latinopoulos et al. (1988) obtained analytical solutions for chemical transport in a 2D aquifer under the constant seepage flow, linear sorption and first order decay for various sets of initial and boundary conditions. Simmons et al. (2001) studied different densities, groundwater movement and contaminant distribution in a heterogeneous porous medium and discussed some approaches, resolutions and future challenges on groundwater contamination. Sander & Braddock (2005) overviewed analytical solution for the transportation of solute with instantaneous injection in the stream subject to temporary storage and a range of mathematical explanations to the transient, unsaturated distribution of water and plume through horizontal soil medium. A precise solution for 2D ADE under anisotropic dispersion in a cylindrical geometry subject to some initial conditions was discussed by Massabó et al. (2011), however Zhan et al. (2009) explored the 2D contaminant distribution in an aquifer–aquitard system with first and third type boundary conditions. Murio (2008) and Chen et al. (2008) derived the concept of a 1D linear time-fractional diffusion equation by implicit finite difference approximation and discussed the solution with distance-dependent dispersivity for 2D ADE by a power series solution. Furthermore, Savović & Djordjevich (2012) discussed 1D ADE with the variable coefficient in semi-infinite medium by applying a finite difference method. Later, Djordjevich & Savović (2013) studied longitudinal and transversal pollutant dispersions subject to space- and time-dependent flow in the porous medium.

Input boundary conditions highly influence the fate and transport phenomenon of the pollutant in an aquifer. Previously, many results have been reported for plume migration subject to time-invariant fixed inlet conditions. However, in broad areas of science and engineering, the contaminant transport problem may involve different types of space and temporally dependent inlet boundary conditions. In recent years; exponential input sources at the inlet boundary have been largely used to predict contaminant movement in a radioactive garbage treatment. Chen & Liu (2011) and Chen et al. (2011) examined the generalized mathematical solutions for the ADE in a bounded spatial domain with time-dependent inlet boundary conditions and the authors also incorporated the effect of reaction term rate on contaminant distribution in the aquifer.

The waste released by some large industrial plants and buried into the soil sub-surface, and other input sources such as pesticides used in a crop field, medical wastes or garbage dumped into the soil sub-surface etc., may disperse through the soil medium and ultimately reach the groundwater reservoir. Some research studies have also been published in the past on pollutant dispersion with depth-dependent sorption and degradation in finite and semi-infinite porous media. A mathematical study of the 1D ADE with depth-dependent decay coefficients was presented earlier (Flury et al. 1998; Diaw et al. 2001). The analytical solution of a mathematical formulation for simulating the water distribution and plume migration in saturated–unsaturated porous media is derived using discontinuous finite element methods. Later on, the work of Flury et al. (1998) was extended by Gao et al. (2013) and described the contaminant distribution in a porous medium with depth-dependent reaction coefficient and time-dependent inlet boundary conditions. Leij et al. (2016) discussed analytical solutions for colloidal distribution with time- and depth-dependent reactions in porous medium. Farenhorst et al. (2009) studied variation in soil attribute and herbicide sorption coefficients with depth to simulate the pesticide root zone model. Kumar & Yadav (2015) elaborated a 1D solute transport for uniform and varying pulse-type input point sources through the heterogeneous medium, and Xie et al. (2019) derived a new technique for the solution of 1D reactive solute transport with a depth-dependent reaction coefficient in dual permeable media. Chaudhary et al. (2020) discussed 1D pollutant transport using the Dirichlet and Neumann type boundary conditions and demonstrated the pollutant distribution in different porous media. Also, Savović & Djordjevich (2020) used explicit finite difference methods to solve 1D ADE and observed the impact of constant and oscillating type boundary conditions on pollutant transport in the groundwater system. Various types of mathematical techniques, such as the Laplace transform method, Hankel transform method, perturbation approach, and Green's function method have been used previously to predict the desired analytical results of the contaminant migration problem in one-, two- and three-dimensions (Park & Zhan 2001). Cohn & Logan (1995) discussed a mathematical approach for the reactive-diffusive model and the dispersion of chemical tracer with non-linear convection.

In the present study, the Laplace transform technique and Crank–Nicolson finite difference method were adopted to solve the governing transport equation analytically and numerically, respectively. Various mathematical models based on dispersion theory were used mainly to predict the contaminant concentration distribution and remediation in groundwater and soil. When the pollutant leaches through soil under some actions and decomposition such as sorption and degradation, it may affect the contaminant strength in the groundwater system. van Genchutan (1982) extensively discussed the solute transport equation in a most general form and presented a sufficient number of analytical solutions including the zero-order production and/or first order decay subject to the first and third type boundary conditions. Later on, 1D plume migration through porous media with spatial varying retardation, and mathematical explanation for 2D transport model with temporally dependent dispersion coefficient were also studied (Chrysikopoulos et al. 1990; Aral & Liao 1996). The effect of contaminant near the disposal site can be viewed as well as in shallow aquifers or reservoirs mostly affected by the dispersion of contaminants compared with the deeper well (Gilliom 2007). The location of the disposal site also has importance in soil and groundwater pollution and its remediation. It can be seen from the previous literature that the analytical study of the solute transport problem in semi-infinite or infinite porous domain has been the most prefered area due to the lack of progress of analytical solutions in the finite porous domain. Savović & Djordjevich (2013) studied temporally dependent dispersion along uniform flow field and spatially dependent dispersion along non-uniform flow for the 2D ADE in a semi-infinite porous medium. The explicit finite difference method was used to derive the numerical solution.

In this study, 2D ADE is discussed to formulate a model of contaminant transport in a finite homogeneous porous medium. Initially, the reservoir is assumed to have some background contaminants. A depth-decaying input source is considered at the upstream of the finite domain. Such type of depth varying source is fascinating according to hydrological aspects. The flux type boundary conditions are assumed at the downstream of the domain. The longitudinal and transverse dispersion coefficients are taken into consideration, and the seepage flow is considered along the longitudinal direction. Impact of linear sorption and first-order decay are also incorporated in the contaminant transport equation. Laplace transform technique is adopted to solve the equation analytically and for its validation by numerical solution, Crank–Nicolson scheme is used. To our knowledge, no such results have been established previously for 2D ADE with a depth varying input source in a finite porous medium. The main focus of this study was to investigate the 2D contaminant transport behaviour with a depth varying input source along the groundwater flow site in a finite porous medium.

Contamination of the soil medium or groundwater reservoir due to leakage of the septic tank or underground pipelines can be modelled mathematically. The 2D contaminant transport in a homogeneous finite porous medium under sorption condition and first-order decay is mathematically modelled as follows:
(1)
where, is resident plume concentration, is solid-phase concentration, is longitudinal dispersion, is transverse dispersion, is pore-water velocity in x-direction, is solid bulk density of the medium, is time variable, and are spatial variables having domain lengths L and H respectively, is porosity of the porous medium, and is a first-order decay rate coefficient.
Sorption processes take place when a contaminant enters the groundwater reservoir through the soil. The linear sorption isotherm between the liquid and solid-phase contaminant concentrations is defined as Zheng & Bennett (2002):
(2)
where, is the sorption coefficient.
Using Equation (2) in Equation (1) gives
(3)
where, is the retardation factor.
Initially, the porous medium is assumed contaminated with some background concentration before the study of the distribution of input source concentration. A constant initial background contaminant concentration is assigned throughout the porous medium as follows:
(4)
As the contaminants leak from the septic tanks or underground pipelines, an input source concentration decaying exponentially with depth is introduced at the inlet boundary. It can be expressed mathematically as
(5)
where, is a constant input concentration at the inlet boundary and is the depth-decay parameter that is useful for controlling the source concentration.
The contaminant fluxes across the extreme boundaries are assumed to be zero, which is mathematically expressed as follows:
(6)
(7)

The geometry of the 2D transport model problem is shown in the Figure 1. The input contaminant source decaying exponentially with depth is shown along the z-axis and the contaminant fluxes are assigned at the extreme boundaries. The groundwater flow is considered along the x-direction.

Figure 1

Geometry of the 2D contaminant transport with depth varying input source in a finite porous medium.

Figure 1

Geometry of the 2D contaminant transport with depth varying input source in a finite porous medium.

Close modal
In this study, the dispersion is assumed as directly proportional to the first power of the seepage velocity, i.e., as described by Freeze & Cherry (1979). Let , where is dispersivity constant of the medium that depends on pore geometry of the porous medium and u is the seepage velocity. The groundwater flow may be transient due to the seasonal variation throughout the year and the activities of human beings. For the transient or unsteady flow we can define , , and , where and are constant initial dispersion coefficients, is constant initial seepage velocity, and is the rate of flow resistance coefficient caused by liquid–solid interaction (Kumar 1983a, 1983b; Singh et al. 2016). Substituting these values in Equation (3) yields
(8)
A transformation in new time variable T is introduced as (Crank 1979)
(9)
Using Equation (9) in Equation (8) yields
(10)
where, .
The initial and boundary conditions are expressed in the new time variable as follows:
(11)
(12)
(13)
(14)
Before solving Equations (10)–(14) analytically, the two-dimensional equation was converted into a one-dimensional equation by introducing a new space variable as follows:
(15)
Previously, Carnahan & Remer (1984) used this type of transformation to reduce two- or three-dimensional ADE into one-dimensional ADE. Applying the transformation defined in Equation (15), Equation (10) reduces to the following one-dimensional form.
(16)
where, and
Also, the initial and boundary conditions defined in Equations (11)–(14) are transformed in one-dimensional space, respectively as follows:
(17)
(18)
(19)
Equation (16) can be reduced into the diffusion equation by applying concentration transformation as follows (Yadav et al. 2012; Singh & Das 2015):
(20)
On applying the concentration transformation defined in Equation (20), Equations (16)–(19) reduce into the following forms:
(21)
(22)
(23)
(24)
where, and .
Now, the above diffusion Equation (21) is solved by the Laplace transform technique and the solution of the governing equation in the Laplacian domain is obtained as follows:
(25)
where, p is a parameter in Laplacian domain, is the space variable defined in the Equation (15), is a parameter varying in the range , and .
Applying the inverse Laplace transform in Equation (25), we obtained as follows:
(26)
where,
(27)
(28)
and
(29)
The other symbols used in Equations (27)–(29) are defined as follows:

Finally, the required analytical solution of the 2D governing ADE in the real time domain can be obtained by substituting from Equation (26) into Equation (20). The obtained analytical solution is illustrated graphically in the Results and Discussion section.

Comparison of numerical solutions of the transport equation under absorption and decay was discussed by van Genuchten (1982). Influence of spatial and temporal sizes for the transport equation was also elaborated. In this study, the Crank–Nicolson finite-difference scheme is used to obtain the numerical solution of the governing ADE. Before applying the Crank–Nicolson scheme to approximate Equation (10), the two-dimensional finite xz-plane is first discretized into a finite number of rectangular meshes by introducing intersecting lines parallel to x- and z-axes. The xz-plane is discretized into P and Q number of sub-intervals with sub-interval lengths and , respectively. Again, the time domain T is discretized into M number of sub-intervals with sub-interval length . Let and denote the nodal point indices along x- and z- directions, respectively and denotes the index for time level along the time domain. Now, on the approximation of Equation (10) by Crank–Nicolson finite difference scheme at the grid point, can be expressed in a simplified form as follows (Smith 1985):
(30)
where, , , , , and .
The initial and boundary conditions (11)–(14) are approximated as follows:
(31)
(32)
(33)
(34)

Equations (30)–(34) result in a penta-diagonal system of simultaneous linear algebraic equations for all the longitudinal and transversal spatial nodal points (i.e., k and j) at each time interval (i.e., n) and is solved graphically by MATLAB coding for a fixed time period. The system of linear algebraic equations is solved by direct matrix method in the MATLAB coding. Since the penta-diagonal system of linear algebraic equations becomes very difficult and tedious to evaluate manually, MATLAB coding used to obtain the required solution at a fixed time period. The obtained numerical solution was validated graphically with the analytical solution and good concurrence between them was found that signified the error-free derivation of the analytical solution. The RMS error of the numerical solution with respect to the analytical solution was also found significantly less and is discussed in the Results and Discussion section.

Groundwater contamination due to septic tanks or underground pipelines is widespread in a dense population area. Pollutants through these sources are harmful to humans as well as other species on the Earth's surface, and contaminate all surrounding wells and other groundwater resources. Solute transport with depth varying sources through porous media along the horizontal groundwater seepage flow is solved analytically using the Laplace transform method. The numerical solution was also obtained using the Crank–Nicolson scheme for the validation of analytical solution. In this study, seepage flow is considered as varying sinusoidally, i.e., , where . This type of sinusoidal variation in a groundwater system occurs generally in tropical areas. In tropical areas, the pore-water flow and water level attain maximum values during the highest reach of the winter season when dependency on the groundwater is lowest and minimum value during the highest reach of the summer season when dependency on the groundwater is highest. Hence the fluctuation of groundwater velocity throughout the year in tropical areas represents a sinusoidal behaviour in nature (Kumar & Kumar 1998). Contaminant transport profiles for the solution of the proposed ADE model were obtained by introducing a set of hydrological input values (Gelhar et al. 1992; Singh & Kumari 2014). The longitudinal and transversal spatial domain lengths were taken as and , respectively. For the spatial discretization, and were taken into consideration along the longitudinal and transversal directions, respectively. Also, was taken as the total number of intervals along the temporal domain . All the rest input data used in this study to obtain graphical solutions are listed in Table 1.

Table 1

Input parameters used to find the graphical solutions

ParametersValues
Longitudinal dispersion   
Transversal dispersion   
Seepage flow velocity   
Initial concentration   
Input source concentration   
Distribution coefficient   
Depth controlling rate   
Degradation rate coefficient   
Flow resistance coefficient   
Bulk density of soil   
Porosity   
ParametersValues
Longitudinal dispersion   
Transversal dispersion   
Seepage flow velocity   
Initial concentration   
Input source concentration   
Distribution coefficient   
Depth controlling rate   
Degradation rate coefficient   
Flow resistance coefficient   
Bulk density of soil   
Porosity   

Two-dimensional contaminant transport in the xz-plane is shown in Figure 2. It depicts the contaminant distribution profiles with depth varying input sources along the horizontal groundwater flow in a finite gravel medium at time . The input source strength at the inlet boundary, decaying exponentially with depth, is at maximum value. The contaminant concentration decreases with space towards the extreme boundaries due to dispersion, sorption and first-order decay. The effect of longitudinal and transverse dispersions on the contaminant transport profile can be observed in the space. The contaminant distribution along the x-direction advances up to larger distances due to considering the longitudinal dispersion coefficient greater than the transversal dispersion coefficient. The effect of depth varying input source can be observed in the graph throughout the medium.

Figure 2

Analytical solution of the two-dimensional governing transport equation subject to the depth varying input source.

Figure 2

Analytical solution of the two-dimensional governing transport equation subject to the depth varying input source.

Close modal

Figure 3 shows the validation of the analytical and numerical solutions in the gravel medium at . It can be observed that the numerical solution closely converges with the analytical solution. The concentration distribution pattern of the numerical and analytical solutions advanced from the inlet boundary and decreased with space towards the extreme boundaries. The accuracy of the numerical solution is discussed by evaluating the RMS error in the Table 2.

Table 2

Influence of the spatial and temporal step lengths on RMS error of the numerical solution over the analytical solution

Position Analytical solutionNumerical solutionNumerical solutionNumerical solution
Δx = 0.002Δz = 0.00125ΔT = 0.0045Δx = 0.002Δz = 0.00125ΔT = 0.0045Δx = 0.004Δz = 0.0025ΔT = 0.0045Δx = 0.002Δz = 0.00125ΔT = 0.0018
(0.0000, 0.0000) 1.0000 1.0000 1.0000 1.0000 
(0.0040, 0.0025) 0.8882 0.8912 0.8913 0.8916 
(0.0080, 0.0050) 0.7813 0.7854 0.7854 0.7862 
(0.0120, 0.0075) 0.6805 0.6843 0.6843 0.6854 
(0.0160, 0.0100) 0.5867 0.5895 0.5894 0.5908 
(0.0200, 0.0125) 0.5008 0.5021 0.5021 0.5037 
(0.0240, 0.0150) 0.4234 0.4233 0.4233 0.4249 
(0.0280, 0.0175) 0.3546 0.3534 0.3535 0.3551 
(0.0320, 0.0200) 0.2943 0.2928 0.2929 0.2944 
(0.0360, 0.0225) 0.2424 0.2411 0.2413 0.2427 
(0.0400, 0.0250) 0.1984 0.1979 0.1982 0.1994 
(0.0440, 0.0275) 0.1616 0.1626 0.1630 0.1639 
(0.0480, 0.0350) 0.1314 0.1343 0.1346 0.1354 
(0.0520, 0.0325) 0.1068 0.1120 0.1124 0.1129 
(0.0560, 0.0350) 0.0872 0.0948 0.0952 0.0956 
(0.0600, 0.0375) 0.0718 0.0819 0.0823 0.0826 
(0.0640, 0.0400) 0.0599 0.0724 0.0729 0.0730 
(0.0680, 0.0425) 0.0510 0.0658 0.0662 0.0663 
(0.0720, 0.0450) 0.0448 0.0615 0.0619 0.0619 
(0.0760, 0.0475) 0.0411 0.0590 0.0594 0.0594 
(0.0800, 0.0500) 0.0399 0.0582 0.0586 0.0586 
CPU time 0.1516 0.6429 0.8210 3.6228 
RMS error  0.0086 0.0089 0.0090 
Position Analytical solutionNumerical solutionNumerical solutionNumerical solution
Δx = 0.002Δz = 0.00125ΔT = 0.0045Δx = 0.002Δz = 0.00125ΔT = 0.0045Δx = 0.004Δz = 0.0025ΔT = 0.0045Δx = 0.002Δz = 0.00125ΔT = 0.0018
(0.0000, 0.0000) 1.0000 1.0000 1.0000 1.0000 
(0.0040, 0.0025) 0.8882 0.8912 0.8913 0.8916 
(0.0080, 0.0050) 0.7813 0.7854 0.7854 0.7862 
(0.0120, 0.0075) 0.6805 0.6843 0.6843 0.6854 
(0.0160, 0.0100) 0.5867 0.5895 0.5894 0.5908 
(0.0200, 0.0125) 0.5008 0.5021 0.5021 0.5037 
(0.0240, 0.0150) 0.4234 0.4233 0.4233 0.4249 
(0.0280, 0.0175) 0.3546 0.3534 0.3535 0.3551 
(0.0320, 0.0200) 0.2943 0.2928 0.2929 0.2944 
(0.0360, 0.0225) 0.2424 0.2411 0.2413 0.2427 
(0.0400, 0.0250) 0.1984 0.1979 0.1982 0.1994 
(0.0440, 0.0275) 0.1616 0.1626 0.1630 0.1639 
(0.0480, 0.0350) 0.1314 0.1343 0.1346 0.1354 
(0.0520, 0.0325) 0.1068 0.1120 0.1124 0.1129 
(0.0560, 0.0350) 0.0872 0.0948 0.0952 0.0956 
(0.0600, 0.0375) 0.0718 0.0819 0.0823 0.0826 
(0.0640, 0.0400) 0.0599 0.0724 0.0729 0.0730 
(0.0680, 0.0425) 0.0510 0.0658 0.0662 0.0663 
(0.0720, 0.0450) 0.0448 0.0615 0.0619 0.0619 
(0.0760, 0.0475) 0.0411 0.0590 0.0594 0.0594 
(0.0800, 0.0500) 0.0399 0.0582 0.0586 0.0586 
CPU time 0.1516 0.6429 0.8210 3.6228 
RMS error  0.0086 0.0089 0.0090 
Figure 3

Analytical and numerical solutions of the governing equation.

Figure 3

Analytical and numerical solutions of the governing equation.

Close modal

Figure 4 demonstrates the graphical representation of the analytical solutions obtained by using input data of the present study and from the existing literature (Latinopoulos et al. 1988). The input values of the dispersion coefficient, advection coefficient and others (i.e., , , , , and ) were used from the existing literature to verify the input data used in the present study. From the figure, it is observed that the authors' graphical solution follows the same distribution pattern as the graphical solution obtained using input data from the existing literature.

Figure 4

Analytical solutions of the 2D contaminant transport with depth varying input source obtained by taking input data from this study and Latinopoulos et al. (1988).

Figure 4

Analytical solutions of the 2D contaminant transport with depth varying input source obtained by taking input data from this study and Latinopoulos et al. (1988).

Close modal

The transport of contaminant is examined for the different density values of gravel medium under the sinusoidal varying velocity field at . Figure 5(a) depicts the line plot obtained for and Figure 5(b) and 5(c) are obtained for and , respectively. The larger density of the soil material produced larger retardation in the concentration distribution profile as depicted in Figure 5(a). Pollutant dispersion through the denser porous soil material faces larger drag force due to the solid-liquid interaction. The variation between the contaminant concentration profiles in the different soil densities increases with space. The distribution behaviour of the contaminant concentration with respect to depth and longitudinal directions can be observed from Figure 5(b) and 5(c). This figure depicts the contaminant concentration distribution profiles in the form of contour. The depth varying contaminant source from the inlet boundary (i.e., z-axis) advances in the groundwater flow direction (i.e., x-direction) and attains minimum level at the final boundary. The strength of the contaminant concentration is observed higher along the x-axis at the origin compared with the extreme boundary of the z-axis. Also, the contaminant concentration attained higher value shown in Figure 5(b) compared with the Figure 5(c).

Figure 5

Two-dimensional depth varying input source distribution profile for: (a) and , (b) , and (c) .

Figure 5

Two-dimensional depth varying input source distribution profile for: (a) and , (b) , and (c) .

Close modal

Figure 6 depicts the contaminant concentration distribution profiles in a finite gravel medium at time and , respectively. In Figure 6(a), due to the regular supply of the input source concentration and that partially trapped by soil material, it is observed that the contaminant concentration in the porous medium increases with time but decreases with space and attains a minimum level at the time extreme boundary. It is also observed that, the concentration distribution profiles vary with space and time by imitating each other. Figure 6(b) and 6(c) represent the contour plots for contaminant transport at and , respectively. The contaminant decreased with depth as well as in the longitudinal direction. The contaminant concentration transport profile at the higher time intervals attains greater value along the extreme boundary of the z-axis compared with the extreme boundary of the x-axis as shown in the Figure 6(c).

Figure 6

Contaminant distribution profiles for the two-dimensional depth varying input source at (a) and , (b) , and (c) .

Figure 6

Contaminant distribution profiles for the two-dimensional depth varying input source at (a) and , (b) , and (c) .

Close modal

Figure 7 depicts the contaminant concentration distribution profiles in gravel and clay media at . Larger porosity produces lower retardation in the pollutant concentration distribution. In Figure 7(a), due to the larger porosity of the clay medium, it is observed that the contaminant concentration in the clay medium acquired slightly higher values than the gravel medium. The pattern of the transport profiles reflecting each other in both media. Figure 7(b) and 7(c) show the contour representations of the contaminant distribution in the gravel and clay media, respectively. The depth-decaying contaminant concentration had maximum strength at the inlet boundary, advanced in the longitudinal direction and attained a minimum value at the extreme boundary. The concentration level attained a higher value shown in Figure 7(c) compared with Figure 7(b).

Figure 7

Contaminant distribution profiles for the two-dimensional depth varying input source for: (a) and clay , (b) , and (c) .

Figure 7

Contaminant distribution profiles for the two-dimensional depth varying input source for: (a) and clay , (b) , and (c) .

Close modal

Figure 8 shows the contaminant concentration distribution profiles in the gravel medium for two different values of the distribution coefficient i.e., and . Figure 8(a) depicts the line plot for the contaminant distribution at . The larger distribution coefficient produces larger retardation in the pollutant concentration. It can be inferred from the figure that the contaminant concentration strength attained a lower level for the higher value of the distribution coefficient due to the interaction of the solid-phase concentration with the soil material. Also, the concentration distribution profiles reflecting each other to some extent and then slightly deviated at the extreme boundary. Furthermore, Figure 8(b) and 8(c) show the contour representations of the concentration distribution for the distribution coefficients and , respectively. It can be observed that the contaminant concentration attained a higher value in Figure 8(b) compared with Figure 8(c) due to the increment in the retardation factor on increasing the distribution coefficient.

Figure 8

Contaminant distribution profiles for the two-dimensional depth varying input source for: (a) and , (b) , and (c) .

Figure 8

Contaminant distribution profiles for the two-dimensional depth varying input source for: (a) and , (b) , and (c) .

Close modal
For the accuracy of numerical solution over the analytical solution, RMS error is globally used. The following mathematical formula for RMS error was used in the present study to evaluate the error of the numerical result (Chai & Draxler 2014):
where, and N is sample size of the spatial point concentration values.

Table 2 represents the contaminant concentration distribution values obtained by MATLAB R2013 software on a personal computer with a CORE i5 processor for the analytical and numerical solutions. The concentration values of the analytical and numerical solutions were obtained at 21 spatial points at time . The performance of the numerical solution was examined by calculating its CPU execution time and RMS error over the analytical solution. The CPU execution time is observed to be very less. The accuracy of the numerical solution was analysed by evaluating RMS errors for two different longitudinal and transversal spatial sub-interval lengths i.e., , and , , respectively at a fixed temporal sub-interval length (i.e., ). The RMS error was also obtained for , , and to examine the performance of numerical solution with respect to the temporal interval length. It is observed from the Table 2 that the RMS error of the numerical solution performed well for , , and (RMS error is 0.0086) in comparison with , , and (RMS error is 0.0089). This implied that the numerical performance improves on decreasing the spatial sub-interval lengths. Conversely, the numerical accuracy slightly decreased (i.e., from 0.0086 to 0.0090) on decreasing the temporal sub-interval length. Overall, RMS errors were observed to be significantly very less to qualify the fruitful application of the numerical solution in this study.

The 2D contaminant transport with depth varying input source is discussed in a finite groundwater reservoir. Mathematical solution of the governing equation is obtained analytically by Laplace transform technique. The contaminant transport, with depth varying input source, along the horizontal groundwater flow increases with time, but decreases with space. Compared with the gravel medium, the contaminant concentration attained a slightly higher value in the clay medium. Contaminant concentration attained a higher level with lower density levels compared with the higher amount of the density of gravel medium. Contaminant distribution profiles are also obtained for the different values of the distribution coefficient that attains a higher level for the lower amount of the distribution coefficient. Numerical simulation of the governing equation was also carried out by the Crank–Nicolson scheme for the validation of the analytical solution and found excellent agreement between them. RMS error of the numerical solution over the analytical solution demonstrated the excellent performance. This type of mathematical model becomes useful for study of groundwater pollution and its remediation in urban areas. Contaminant transport with depth varying input source in a two-dimensional finite homogeneous porous media may further be investigated for its transport behaviour in a multi-dimensional semi-infinite heterogeneous porous media.

Authors are thankful to Indian Institute of Technology (Indian School of Mines), Dhanbad for providing financial assistance through SRF scheme.

All relevant data are included in the paper or its Supplementary Information.

Aral
M. M.
Liao
B.
1996
Analytical solutions for two-dimensional transport equation with time-dependent dispersion coefficients
.
Journal of Hydrological Engineering
1
(
1
),
20
32
. DOI: 10.1155/2011/493014.
Barcelona
M.
Wehrmann
A.
Keely
J.
Pettyjohn
W. A.
1990
Contamination of Ground Water
.
US Department of Energy Office of Scientific and Technical Information
.
Batu
V.
2005
Applied Flow and Solute Transport Modeling in Aquifers: Fundamental Principles and Analytical and Numerical Methods
.
CRC Press
,
Boca Raton, FL
.
Bear
J.
1972
Dynamics of Fluids in Porous Media
.
Dover
,
Mineola, New York
.
Bear
J.
1979
Hydraulics of Groundwater
.
McGraw-Hill Inc.
,
New York
.
Carnahan
C. L.
Remer
J. S.
1984
Non-equilibrium and equilibrium sorption with a linear sorption isotherm during mass transport through an infinite porous medium: some analytical solutions
.
Journal of Hydrology
73
(
3–4
),
227
258
.
doi:10.1016/0022-1694(84)90002-7
.
Chai
T.
Draxler
R. R.
2014
Root mean square error (RMSE) or mean absolute error (MAE)?–Arguments against avoiding RMSE in the literature
.
Geoscientific Model Development
7
(
3
),
1247
1250
.
https://doi.org/10.5194/gmd-7-1247-2014
.
Chaudhary
M.
Thakur
C. K.
Singh
M. K.
2020
Analysis of 1-D pollutant transport in semi-infinite groundwater reservoir
.
Environmental Earth Sciences
79
(
1
),
24
.
https://doi.org/10.1007/s12665-019-8748-4
.
Chen
J. S.
Chen
J. T.
Liu
C. W.
Liang
C. P.
Lin
C. W.
2011
Analytical solutions to two-dimensional advection–dispersion equation in cylindrical coordinates in finite domain subject to first-and third-type inlet boundary conditions
.
Journal of Hydrology
405
(
3–4
),
522
531
.
doi:10.5194/hess-2018-462
.
Chen
J. S.
Liu
C. W.
2011
Generalized analytical solution for advection–dispersion equation in finite spatial domain with arbitrary time-dependent inlet boundary condition
.
Hydrology and Earth System Sciences
15
(
8
),
2471
2479
.
doi:10.5194/hess-15-2471-2011
.
Chen
J. S.
Ni
C. F.
Liang
C. P.
2008
Analytical power series solutions to the two-dimensional advection–dispersion equation with distance-dependent dispersivities
.
Hydrological Processes: An International Journal
22
(
24
),
4670
4678
.
doi:10.1002/hyp.7067
.
Chrysikopoulos
C. V.
Kitanidis
P. K.
Roberts
P. V.
1990
Analysis of one-dimensional solute transport through porous media with spatially variable retardation factor
.
Water Resource Research
26
(
3
),
437
446
.
doi:10.1029/WR026i006p01189
.
Cohn
S.
Logan
J. D.
1995
Mathematical analysis of a reactive-diffusive model of the dispersal of a chemical tracer with nonlinear convection
.
Mathematical Models and Methods in Applied Sciences
5
(
01
),
29
46
.
doi:10.1142/S0218202595000036
.
Crank
J.
1979
The Mathematics of Diffusion
.
Oxford University Press
,
New York
.
doi:10.1002/hyp.7067
.
Delay
F.
Porel
G.
de Marsily
G.
1997
Predicting solute transport in heterogeneous media from results obtained in homogeneous ones: an experimental approach
.
Journal of Contaminant Hydrology
25
(
1–2
),
63
84
.
https://doi.org/10.1016/S0169-7722(96)00020-4
.
Diaw
E. B.
Lehmann
F.
Ackerer
P.
2001
One-dimensional simulation of solute transfer in saturated-unsaturated porous media using the discontinuous finite element method
.
Journal of Contaminant Hydrology
51
(
3–4
),
197
213
.
doi:10.1016/S0169-7722(01)00129-2
.
Djordjevich
A.
Savović
S.
2013
Solute transport with longitudinal and transverse diffusion in temporally and spatially dependent flow from a pulse type source
.
International Journal of Heat and Mass Transfer
65
,
321
326
.
doi:10.1016/j.ijheatmasstransfer.2013.06.002
.
Domenico
P. A.
Schwartz
F. W.
1990
Physical and Chemical Hydrogeology
.
Wiley Publication
,
New York
.
Ebach
E. A.
White
R. R.
1958
The mixing of fluids flowing through beds of packed solids
.
AIChE Journal
4
(
2
).
http://dx.doi.org/10.1002/aic.690040209
.
Farenhorst
A.
Mc Queen
D. A. R.
Saiyed
I.
Hilderbrand
C.
Li
S.
Lobb
D. A.
Lindstrom
M. J.
2009
Variations in soil properties and herbicide sorption coefficients with depth in relation to PRZM (pesticide root zone model) calculations
.
Geoderma
150
(
3–4
),
267
277
.
doi:10.1016/j.geoderma.2009.02.002
.
Fetter
C. W.
Boving
T. B.
Kreamer
D. K.
1999
Contaminant Hydrogeology
, vol.
500
.
Prentice Hall
,
Upper Saddle River, NJ
.
Flury
M.
Wu
Q. J.
Wu
L.
Xu
L.
1998
Analytical solution for solute transport with depth-dependent transformation or sorption coefficients
.
Water Resource Research
34
(
11
),
2931
2937
.
doi:10.1029/98WR02299
.
Freeze
R. A.
Cherry
J. A.
1979
Groundwater
, vol.
7632
.
Prentice-Hall Inc.
,
Englewood Cliffs
, p.
604
.
Gao
G.
Fu
B.
Zhan
H.
Ma
Y.
2013
Contaminant transport in soil with depth-dependent reaction coefficients and time-dependent boundary conditions
.
Water Research
47
(
7
),
2507
2522
.
doi:10.1016/j.watres.2013.02.021
.
Gelhar
L. W.
Welty
C.
Rehfeldt
K. R.
1992
A critical review of data on field-scale dispersion in aquifer
.
Water Resources Research
28
(
7
),
1955
1974
.
doi:10.1029/92WR00607
.
Gilliom
R. J.
2007
Pesticides in US streams and groundwater
.
Environmental Science & Technology
41
(
10
),
3408
3414
.
doi:10.1021/es072531u
.
Guerrero
J. S. P.
Pontedeiro
E. M.
van Genuchten
M. T.
Skaggs
T. H.
2013
Analytical solutions of the one dimensional advection–dispersion solute transport equation subject to time-dependent boundary condition
.
Chemical Engineering Journal
221
,
487
491
.
https://doi.org/10.1016/j.cej.2013.01.095
.
Guerrero
J. P.
Skaggs
T. H.
2010
Analytical solution for one-dimensional advection–dispersion transport equation with distance-dependent coefficients
.
Journal of Hydrology
390
(
1–2
),
57
65
.
https://doi.org/10.1016/j.jhydrol.2010.06.030
.
Kehew
A. E.
2000
Applied Chemical Hydrogeology
.
Prentice Hall
,
Englewood Cliffs
. .
Kumar
N.
1983a
Dispersion of pollutant in semi-infinite porous media with unsteady velocity distribution
.
Hydrology Research
14
(
3
),
167
178
.
https://doi.org/10.2166/nh.1983.0014
.
Kumar
N.
1983b
Unsteady flow against dispersion in finite porous media
.
Journal of Hydrology
63
(
3–4
),
345
358
.
https://doi.org/10.1016/0022-1694(83)90050-1
.
Kumar
N.
Kumar
M.
1998
Solute dispersion along unsteady groundwater flow in a semi-infinite aquifer
.
Hydrology and Earth System Sciences
2
(
1
),
93
100
.
https://doi.org/10.5194/hess-2-93-1998
.
Kumar
A.
Yadav
R. R.
2015
One-dimensional solute transport for uniform and varying pulse type input point source through heterogeneous medium
.
Environmental Technology
36
(
4
),
487
495
.
doi:10.1080/09593330.2014.952675
.
Latinopoulos
P.
Tolikas
D.
Mylopoulos
Y.
1988
Analytical solutions for two-dimensional chemical transport in aquifers
.
Journal of Hydrology
98
(
1–2
),
11
19
.
https://doi.org/10.1016/0022-1694(88)90202-8
.
Leij
F. J.
Bradford
S. A.
Sciortino
A.
2016
Analytic solutions for colloid transport with time- and depth-dependent retention in porous media
.
Journal of Contaminant Hydrology
195
,
40
51
.
doi:10.1016/j.jconhyd.2016.10.006
.
Liu
C.
Ball
W. P.
Ellis
J. H.
1998
An analytical solution to the one-dimensional solute advection–dispersion equation in multi-layer porous media
.
Transport in Porous Media
30
(
1
),
25
43
.
https://doi.org/10.1023/A:1006596904771
.
Massabó
M.
Cianci
R.
Paladino
O.
2011
An analytical solution of the advection dispersion equation in a bounded domain and its application to laboratory experiments
.
Journal of Applied Mathematics
2011
.
doi:10.1155/2011/493014
.
Murio
D. A.
2008
Implicit finite difference approximation for time fractional diffusion equations
.
Computers & Mathematics with Applications
56
(
4
),
1138
1145
.
doi:10.1016/j.camwa.2011.02.051
.
Ogata
A.
1970
Theory of Dispersion in A Granular Medium
.
US Government Printing Office
,
Washington
.
https://doi.org/10.3133/pp411I
.
Park
E.
Zhan
H.
2001
Analytical solutions of contaminant transport from finite one-, two-, and three-dimensional sources in a finite-thickness aquifer
.
Journal of Contaminant Hydrology
53
(
1–2
),
41
61
.
doi:10.1016/S0169-7722(01)00136-X
.
Rausch
R.
Schäfer
W.
Therrien
R.
Wagner
C.
2005
Solute Transport Modelling: An Introduction to Models and Solution Strategies
.
Gebrüder Borntraeger Verlagsbuchhandlung
,
Germany
.
Ritter
W. F.
Chirnside
A. E.
1984
Impact of land use on groundwater quality in Southern Delaware
.
Groundwater
22
(
1
),
38
47
.
doi:10.1111/j.1745-6584.1984.tb01474.x
.
Sander
G. C.
Braddock
R. D.
2005
Analytical solutions to the transient, unsaturated transport of water and contaminants through horizontal porous media
.
Advanced in Water Resources
28
(
10
),
1102
1111
.
doi:10.1016/j.advwatres.2004.10.010
.
Savović
S.
Djordjevich
A.
2012
Finite difference solution of the one-dimensional advection–diffusion equation with variable coefficients in semi-infinite media
.
International Journal of Heat and Mass Transfer
55
(
15–16
),
4291
4294
.
doi:10.1016/j.ijheatmasstransfer.2012.03.073
.
Savović
S.
Djordjevich
A.
2013
Numerical solution for temporally and spatially dependent solute dispersion of pulse type input concentration in semi-infinite media
.
International Journal of Heat and Mass Transfer
60
,
291
295
.
doi:10.1016/j.ijheatmasstransfer.2013.01.027
.
Savović
S. M.
Djordjevich
A.
2020
Explicit finite difference solution for contaminant transport problems with constant and oscillating boundary conditions
.
Thermal Science
24
(
3 Part B
),
2225
2231
.
doi:10.2298/TSCI190722422S
.
Simmons
C. T.
Fenstemaker
T. R.
Sharp
J. M.
Jr.
2001
Variable-density groundwater flow and solute transport in heterogeneous porous media: approaches, resolutions and future challenges
.
Journal of Contaminant Hydrology
52
(
1–4
),
245
275
.
doi:10.1007/s10040-004-0408-3
.
Singh
M. K.
Das
P.
2015
Scale dependent solute dispersion with linear isotherm in heterogeneous medium
.
Journal of Hydrology
520
,
289
299
.
doi:10.1016/j.jhydrol.2014.11.061
.
Singh
M. K.
Kumari
P.
2014
Contaminant concentration prediction along unsteady groundwater flow
. In:
Modelling and Simulation of Diffusive Processes. Simulation Foundations, Methods and Applications
(
Basu
S.
Kumar
N.
, eds).
Springer
,
Cham
, pp.
257
275
.
doi:10.1007/978-3-319-05657-9_12
.
Singh
M. K.
Singh
V. P.
Das
P.
2016
Mathematical modeling for solute transport in aquifer
.
Journal of Hydroinformatics
18
(
3
),
481
499
.
https://doi.org/10.2166/hydro.2015.034
.
Smith
G. D.
1985
Numerical Solution of Partial Differential Equations: Finite Difference Methods
.
Oxford University press
,
New York
.
Sun
N. Z.
Sun
A.
1996
Mathematical Modeling of Groundwater Pollution
.
Springer, Verlag New York, Inc
.
doi:10.1007/978-1-4757-2558-2
.
Todd
D. K.
Mays
L. W.
2004
Groundwater Hydrology
.
John Wiley & Sons
,
New York
.
van Genuchten
M. T.
1981
Analytical solutions for chemical transport with simultaneous adsorption, zero-order production and first-order decay
.
Journal of Hydrology
49
(
3–4
),
213
233
.
https://doi.org/10.1016/0022-1694(81)90214-6
.
van Genuchten
M. T.
1982
A comparison of numerical solutions of the one-dimensional unsaturated-saturated flow and mass transport equations
.
Advanced in Water Resources
5
(
1
),
47
55
.
doi:10.1016/S0309-1708(96)00059-0
.
van Genuchten
M. T.
Leij
F. J.
Skaggs
T. H.
Toride
N.
Bradford
S. A.
Pontedeiro
E. M.
2013
Exact analytical solutions for contaminant transport in rivers 1. The equilibrium advection–dispersion equation
.
Journal of Hydrology and Hydromechanics
61
(
2
),
146
160
.
https://doi.org/10.2478/johh-2013-0020
.
Xie
S.
Wen
Z.
Jakada
H.
2019
A new model approach for reactive solute transport in dual-permeability media with depth-dependent reaction coefficients
.
Journal of Hydrology
577
,
123946
.
doi:10.1016/j.jhydrol.2019.123946
.
Yadav
S. K.
Kumar
A.
Kumar
N.
2012
Solute transport from a pulse type source along temporally and spatially dependent flow: analytical solution
.
Journal of Hydrology
412
,
193
199
.
doi:10.1016/j.jhydrol.2011.02.024
.
Zhan
H.
Wen
Z.
Huang
G.
Sun
D.
2009
Analytical solution of two-dimensional solute transport in an aquifer-aquitard system
.
Journal of Contaminant Hydrology
107
(
3–4
),
162
174
.
doi:10.1016/j.jconhyd.2009.04.010
.
Zheng
C.
Bennett
G. D.
2002
Applied Contaminant Transport Modelling
, vol.
2
.
Wiley-Interscience
,
New York
,
353
pp.
This is an Open Access article distributed under the terms of the Creative Commons Attribution Licence (CC BY-NC-ND 4.0), which permits copying and redistribution for non-commercial purposes with no derivatives, provided the original work is properly cited (http://creativecommons.org/licenses/by-nc-nd/4.0/).