Abstract
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.
HIGHLIGHTS
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
INTRODUCTION
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.
MATHEMATICAL FORMULATION OF THE CONTAMINANT TRANSPORT EQUATION














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.
Geometry of the 2D contaminant transport with depth varying input source in a finite porous medium.
Geometry of the 2D contaminant transport with depth varying input source in a finite porous medium.










ANALYTICAL SOLUTION OF THE GOVERNING EQUATION






NUMERICAL SIMULATION OF THE GOVERNING EQUATION












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.
RESULTS AND DISCUSSION
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.
Input parameters used to find the graphical solutions
Parameters . | Values . |
---|---|
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 ![]() | ![]() |
Parameters . | Values . |
---|---|
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.
Analytical solution of the two-dimensional governing transport equation subject to the depth varying input source.
Analytical solution of the two-dimensional governing transport equation subject to the depth varying input source.
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.
Influence of the spatial and temporal step lengths on RMS error of the numerical solution over the analytical solution
Position ![]() | Analytical solution . | Numerical solution . | Numerical solution . | Numerical 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 solution . | Numerical solution . | Numerical solution . | Numerical 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 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.
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).
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).
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).
Two-dimensional depth varying input source distribution profile for: (a) and
, (b)
, and (c)
.
Two-dimensional depth varying input source distribution profile for: (a) and
, (b)
, and (c)
.
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).
Contaminant distribution profiles for the two-dimensional depth varying input source at (a) and
, (b)
, and (c)
.
Contaminant distribution profiles for the two-dimensional depth varying input source at (a) and
, (b)
, and (c)
.
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).
Contaminant distribution profiles for the two-dimensional depth varying input source for: (a) and clay
, (b)
, and (c)
.
Contaminant distribution profiles for the two-dimensional depth varying input source for: (a) and clay
, (b)
, and (c)
.
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.
Contaminant distribution profiles for the two-dimensional depth varying input source for: (a) and
, (b)
, and (c)
.
Contaminant distribution profiles for the two-dimensional depth varying input source for: (a) and
, (b)
, and (c)
.

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.
CONCLUSIONS
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.
ACKNOWLEDGEMENTS
Authors are thankful to Indian Institute of Technology (Indian School of Mines), Dhanbad for providing financial assistance through SRF scheme.
DATA AVAILABILITY STATEMENT
All relevant data are included in the paper or its Supplementary Information.