Listen

Groundwater pollution may occur due to human activities, industrial effluents, cemeteries, mine spoils, etc. This paper deals with one-dimensional mathematical modeling of solute transport in finite aquifers. The governing equation for solute transport by unsteady groundwater flow is solved analytically by the Laplace transform technique. Initially, the aquifer is subjected to the spatially dependent source concentration with zero-order production. One end of the aquifer receives the source concentration and is represented by a mixed-type boundary condition in the splitting time domain. The concentration gradient at the other end of the porous media is assumed to be zero. The temporally dependent velocity and the dispersion coefficients are considered. A numerical solution is obtained by using an explicit finite difference scheme and compared with the analytical result. Accuracy of the solution is discussed by using the root mean square error method. Truncation error is also explored for the parameters like numerical dispersion and velocity terms. The impact of Peclet number is examined. For graphical interpretation, unsteady velocity expressions (i.e., such as exponential, sinusoidal, asymptotic, and algebraic sigmoid) are considered. The work may be used as a preliminary predictive tool for groundwater resource and management.

analytical modeling, aquifer, dispersion, Laplace transform, solute, unsteady
Listen

In India, many aquifers are being contaminated by a host of human activities, such as sewage disposal, refuse dumps, pesticides and chemical fertilizer contamination, industrial effluent discharges, and toxic waste disposal (Rausch et al. 2005; Batu 2006). The traditional advection dispersion equation is based on the conservation of mass and Fick's law of diffusion (Fried & Combarnous 1971; Bear & Verruijt 1987; Chrysikopolous et al. 1990) and constitutes the basis of solute transport models that are used for predicting the movement of contaminants in groundwater systems.

There has been some research on solute transport in groundwater systems. Hunt (1978) analyzed one-, two-, and three-dimensional solutions for instantaneous, continuous, and steady state pollution sources in uniform groundwater flow. Freeze & Cherry (1979) provided a relation between dispersion and groundwater velocity in which the dispersion is proportional to a power of the velocity and experimentally observed that the power ranged between 1 and 2. van Genuchten (1981) explored derivations of analytical solutions using the Laplace transform for the solute transport equation with zero-order production and/or first-order decay subjected to first and third type boundary conditions. Zoppou & Knight (1994) evaluated analytical solutions that are still useful for validating numerical schemes for solving the advection–diffusion equation with spatially variable coefficients. Logan (1996) obtained analytical solutions for a scale-dependent dispersion coefficient increasing exponentially with position up to some constant limiting values. Hantush & Marino (1998) developed analytical solutions using the Laplace and Fourier transform methods and superposition principle for the first-order rate model in an infinite porous medium.

Using Peclet and Courant numbers and a new sink/source dimensionless number, Ataie-Ashtiani et al. (1999) discussed truncation errors associated with finite difference solutions of the advection–dispersion equation (ADE) with first-order reaction. Bedient et al. (1999) presented a mathematical model of the ADE for describing the migration and fate of pollutants in groundwater. Neville et al. (2000) presented semi-analytical solutions for a multi-process non-equilibrium model for describing contaminant concentration distribution patterns. Balla et al. (2002) presented a computational case study using a transport model for pollution of underground water due to damage of the waterproofing system in a waste material depository or sewage sludge composting plant.

Lowry & Li (2005) discussed an improved finite analytical solution method for solving the time-dependent ADE that does not discretize the derivative terms rather solving the equation analytically in the space–time domain. Smedt (2006) presented an analytical solution for solute transport in rivers, including the effect of transient storage and first-order decay. Tkalich (2006) explored the derivation of high-order advection–diffusion schemes by employing the interpolation polynomial method. Hill et al. (2007) proposed upscaling models of solute transport in porous media through genetic programming in heterogeneous porous media. Chen (2007) presented an analytical solution of two-dimensional non-axisymmetrical solute transport in a radially convergent flow tracer test with a diffusion coefficient increasing with travel distance. Yeh & Yeh (2007) derived solutions of the transport equation with a point-source term considered as the point-source solution under the condition that the solute was introduced into the flow system from the boundary that was considered as the boundary-source solution. Kumar et al. (2008) also described transport through a heterogeneous porous medium with a time-dependent dispersivity in solute transport modeling. Zhan et al. (2009) explored two-dimensional solute transport in an aquifer–aquitard system by maintaining the mass conservation at the aquifer–aquitard interface. Gao et al. (2012) explored a mobile–immobile model with an asymptotic dispersivity function of travel distance with the concept of scale-dependent dispersion during solute transport in finite heterogeneous porous media. Rezaei et al. (2013) derived a semi-analytical solution to the two-dimensional conservative solute transport in an aquifer bounded by thin aquitards in the Laplace domain. Singh & Das (2015) explored the analytical and numerical solutions for one-dimensional scale-dependent solute transport in heterogeneous media in which analytical and numerical solutions were compared and found very good agreement among them. The root mean square (RMS) error analysis was made to check the accuracy of the solution.

In the present work, we focus on one-dimensional solute transport modeling using the ADE in a finite aquifer with first-order decay and zero-order production. To simplify ADE, different transformations were applied. Non-dimensional parameters were employed for reducing the number of parameters of the ADE. To predict the pattern of contaminant concentration, different types of unsteady velocities, such as sinusoidal, exponentially decreasing, asymptotic, and algebraic sigmoid, were considered. They helped describe the nature of contaminant concentration in time and space.

Listen
The contaminant transport phenomena in a groundwater layered medium or reservoir, i.e., aquifer, can be predicted by analytical solute transport modeling of the ADE. The modeling depends on the distance as well as time. The dominant process of mass transfer in groundwater is the advection and diffusion; in contrast, it refers to the solute transport through the action of random motions in the groundwater reservoir. The one-dimensional ADE for an isotropic homogeneous finite aquifer with zero-order production and first-order decay can be expressed as:
1
Here, is the longitudinal dispersion coefficient (i.e., representing longitudinal dispersion), is the volume averaged dispersing solute concentration in the liquid phase, is the volume averaged dispersing solute concentration in the solid phase, is the unsteady uniform downward pore seepage velocity, is the longitudinal direction of flow, is time, is the zero order production rate coefficients for solute production in the liquid phase, is the first-order decay rate coefficient in the liquid phase, is the bulk density of the porous medium, and n is the porosity of geological formation.
As the contaminant goes from the solid phase into the liquid phase under the linear isotherm condition, we can express mathematically:
2
where is the distribution coefficient, defined as the concentration adsorbed by the solid phase to the liquid phase into the groundwater reservoir, i.e., aquifer.
Equation (1) was solved analytically with appropriate initial and boundary conditions. As the aquifer is not solute free initially, i.e., at ⁠, a linear combination of initial concentration and a function of space-dependent along with zero-order production term is considered. The loss or gain of solute mass occurs due to chemical components within the liquid, radioactive decay and biodegradation, and growth due to bacterial activities. The loss of gain term is usually described through the sources of solute. The solute can grow in liquid phase and the solid phase. Its growth is frequently expressed by a zero-order production. In this study, initially the growth of solute along the space was a linear combination of the initial concentration taken into consideration. This can be written as:
3
Due to the increasing human activity at the earth surface, the solute concentration in groundwater increases in time. Hence, a mixed type boundary condition in the splitting time domain at the source is considered as follows:
4
Due to no mass flow at the other end of the domain, a flux type boundary condition is considered as follows:
5
Using Equation (2), Equation (1) can be written as:
6
where
7
The dispersion theory, proposed by Freeze & Cherry (1979), was employed here. As dispersion is directly proportional to the seepage velocity, i.e., ⁠, and so, ⁠, where, A is constant that depends upon the pore geometry of the aquifer. Let and ⁠, where, and are the initial seepage velocity and initial dispersion coefficient, respectively. Also, let and ⁠, where, and are the initial first-order decay rate coefficient and the initial zero-order production rate coefficient for solute production in the liquid phase, respectively. is the non-dimensional expression where is the flow resistance coefficient. Using this concept, Equation (6) can be written as:
8
A new time variable (Crank 1975):
9
is introduced and Equation (8) can be written as:
10
The following non-dimensional variables:
11
are used to reduce the number of parameters of Equation (10) and one can then write:
12
where
13
The initial and boundary conditions given in Equation (3), Equation (5), and Equation (6) can be written in non-dimensional form, respectively, as follows:
14
15
16
Now, the following transformation:
17
was used to remove the advection term from Equation (12). The initial and boundary conditions accordingly were transformed and the Laplace transform technique was employed. The analytical solution can be obtained as follows (see Appendix, available with the online version of this paper):
18
19
where
20a
20b
20c
20d
20e
21a
21b
21c
21d
21e
21f
21g
21h
21i
21j
21k
21l
21m
22a
22b
22c
22d
22e
22f
22g
22h
22i
22j
22k
22l
23a
23b
23c
23d
23e
23f
23g
23h
23i
23j
23k
24a
24b
24c
24d
24e
24f
24g
24h
24i
24j
24k
25a
25b
25c
25d
25e
25f
25g
25h
25i
25j
25k
25l
25m
26a
26b
26c
26d
26e
26f
26g
26h
26i
27a
27b
27c
27d
27e
27f
27g
27h
28a
28b
28c
28d
28e
28f
28g
29
30
31
Listen
The finite-difference technique has been applied in various works for numerically modeling the solute transport phenomena in groundwater systems (i.e., Molenkamp 1968; van Genuchten & Gray 1978; Richard 1988; Hogarth et al. 1990; Dudley et al. 1991; Moldrup et al. 1996; Zheng & Bennett 2002; Meerschaer & Tadjeran 2004; Ataie-Ashtiani & Hosseini 2005; Wang et al. 2005; Mohebbi & Abbaszadeh 2013). In order to solve numerically the ADE (Equation (12)) together with the initial and boundary condition transformed in domain ⁠, a suitable space transform was used as:
32
Equation (12) together with the initial and boundary condition may be written as:
33
34
35
36
A finite difference technique was derived by using Taylor's series expansion (Mickley et al. 1957; Carnahan et al. 1969). The explicit finite difference scheme is commonly used, even though it may require an extended computing time because of its restrictive stability criterion. Here, we used a general form of the explicit finite difference approximation with forward time and central space forward difference scheme to approximate Equations (33)–(36). We obtained as follows:
37
38
39
40
where superscript j refers to time, subscript i refers to space, is the time increment, and is the space increment.
The space domain Z and time domain T are discretized by a rectangular grid of points with mesh size and ⁠, respectively. So, one can write as follows:
The contaminant concentration at a point in space with the jth subinterval of time T was defined by ⁠.
Listen
The finite difference scheme is convergent if the discretization error approaches zero as the grid spacing and tends to zero (Carnahan et al. 1969). The stability test of the finite difference scheme was proposed by a matrix method (Smith 1978) and this was used by Notodarmojo et al. (1991,). Kwok (1992) investigated the stability properties of the various two-level, six-point finite difference schemes for the approximation of the convection–diffusion equation. The solution was convergent, subject to the satisfaction of the stability criterion. The finite difference scheme of the governing partial differential Equation (37) can be written as follows:
41
where
42
43
and
44
Equation (41) was expressed in matrix form as:
45
where A contains the entire constants.

The difference approximation equation was stable if the eigenvalues of A had modulus values less than or equal to unity, i.e., ⁠, where was the eigenvalue of matrix A.

To find the bounds of the eigenvalues of matrix A on applying the Gerschgorin circle method, the stability criterion for the time step was as follows:
46
Listen

Numerical dispersion was first quantified by Lantz (1971). Ataie-Ashtiani et al. (1999) explored the expansion of the Taylor series of solute concentration along the ADE used for determining the truncation error in one dimension. Chaudhari (1971) investigated a second-order error through the examination of the truncated Taylor series approximation with explicit finite difference solution of the one-dimensional ADE. We also explored the truncation error for the various parameters, such as dispersion, seepage velocity, first-order decay, and zero-order production term.

From Equation (37), we obtained as follows:
47
From the Taylor series expansion of each term of Equation (47), we obtained as follows:
48
49
50
where is simply denoted by C.
After imposing the Taylor series expansion on Equation (47), we got the truncation error of the finite difference approximation of order ⁠. The transport parameters were constant within each combination of time and space increments in finite difference calculations. The second-order temporal derivative of C was written in terms of the spatial derivative of C by using the partial differential equation obtained from Equation (47). We obtained as follows:
51
Then the partial differential equation obtained from Equation (47) can be written as:
52

Now, after comparison between Equation (52) and the original partial differential equation, we found different forms of the truncation error, as discussed by Ataie-Ashtiani et al. (1996). These errors can be identified as follows:

  • The second-order truncation error or numerical dispersion was
    53
  • The first-order truncation error or numerical seepage velocity was
    54
  • The zero-order truncation error or numerical first-order decay was
    55

Constant error or numerical zero-order production was
56
After removing the induced numerical errors from the finite difference model, Equation (52) can be written as follows:
57
where
58
59
60
61
Listen
The accuracy of the solution was obtained by comparison of the numerical result with the analytical one. In the numerical solution, the accuracy is the degree of closeness of concentration values of the numerical result obtained with various methods to those of the analytical result. Towler & Yang (1979) adopted a criterion of comparison that was more systematic and consistent: (1) the RMS error and (2) the absolute maximum error between the analytical solution and the numerical solution at all grid points. Roberts & Selim (1984) used the RMS method to calculate the average error at each nodal point of the grid. Singh & Das (2015) explored the accuracy of solution of the solute transport equation in comparison of the analytical result with numerical one. For testing the accuracy of solutions in this paper, we used the RMS error which is the most appropriate method to check the accuracy of the solution (Chai & Draxler 2014). The RMS method was used to calculate the average error at each point which is defined by:
62
where
The difference between the analytical and numerical concentration values is denoted by ⁠, and N is the number of data which were used to evaluate the accuracy of the solution.
Listen
The analytical solution, obtained by Equations (18) and (19), was computed for the following data (Singh & Kumari 2014):
The average porosity of the different geological formation was considered as n = 0.37 (sand), 0.55 (clay) (Freeze & Cherry 1979). We considered four different time-dependent forms of velocity expressions that can be written as follows:
  1. Exponential decreasing form of velocity:
    63
  2. Sinusoidal form of velocity:
    64
  3. Asymptotic form of velocity:
    65
  4. Algebraic form of velocity:
    66
where m is the flow resistance coefficient and is the constant parameter. The first and third ones have been used by Aral & Liao (1996), the second one by Singh et al. (2009), and the last one is based on the properties of the algebraic sigmoid function which include the error function. It starts to progress from a small beginning, accelerates in the rainy season, and then reaches a limit over a period of time.
From Figure 1 it can be seen that the maximum contaminant concentration is observed, i.e., 0.8 in the case of aquitard (i.e., clay) and 0.6 in the case of aquifer (i.e., sand) which are subsequently lowered down to minimum concentration tends to zero at the far end of the domain with respect to the distance, i.e., and ⁠, respectively. The concentration values increases with time at each of the positions in both the media.
Figure 1
Figure 1. Concentration distribution pattern for exponentially decreasing velocity for time with average porosity of the medium.
View largeDownload slide

Concentration distribution pattern for exponentially decreasing velocity for time with average porosity of the medium.

Figure 1
Figure 1. Concentration distribution pattern for exponentially decreasing velocity for time with average porosity of the medium.
View largeDownload slide

Concentration distribution pattern for exponentially decreasing velocity for time with average porosity of the medium.

Close modal
Figure 2 depicts the contaminant concentration pattern for the unsteady sinusoidal form of velocity with zero order production parameter. The effect of this parameter is predicted with respect to the aquifer (i.e., sand). The concentration level decreases with the increasing value of the zero-order production parameter, but the contaminant concentration increases with increasing time. The peak of the contaminant concentration is lower at the source due to the increasing zero-order production parameter. The effective parameters, like zero-order production, first-order decay, etc., in tropical regions are more significant for transport modeling for groundwater bodies.
Figure 2
Figure 2. Concentration distribution pattern for sinusoidal velocity for time with zero-order production in the aquifer.
View largeDownload slide

Concentration distribution pattern for sinusoidal velocity for time with zero-order production in the aquifer.

Figure 2
Figure 2. Concentration distribution pattern for sinusoidal velocity for time with zero-order production in the aquifer.
View largeDownload slide

Concentration distribution pattern for sinusoidal velocity for time with zero-order production in the aquifer.

Close modal
Figure 3 compares the contaminant concentration patterns for exponentially decreasing type unsteady velocity expression in the time domain between aquifer (i.e., sand) and aquitard (i.e., clay). This is the case when the source has been removed from the aquifer. From this figure it is observed that after removing the source, some initial concentration exists at the origin. As a result, the effect of contaminant concentration increases with distance that attains the maximum peak in the case of the aquitard as compared to the aquifer and then decreases up to a harmless concentration with respect to the distance. The concentration pattern initially decreases with the increasing value of time and after covering some distance it takes the reverse pattern with respect to time in the aquifer and aquitard and ultimately goes to a minimum harmless concentration with respect to distance. From this figure it is observed that the concentration pattern is high in the aquitard as compared to the aquifer.
Figure 3
Figure 3. Concentration distribution for the exponential decreasing form of unsteady form of velocity for time with average porosity of the medium.
View largeDownload slide

Concentration distribution for the exponential decreasing form of unsteady form of velocity for time with average porosity of the medium.

Figure 3
Figure 3. Concentration distribution for the exponential decreasing form of unsteady form of velocity for time with average porosity of the medium.
View largeDownload slide

Concentration distribution for the exponential decreasing form of unsteady form of velocity for time with average porosity of the medium.

Close modal
Figure 4 shows a comparison of the analytical result against the numerical one with averaging porosity of sand (0.37). The concentration distribution pattern follows its decreasing nature with respect to the distance in both the results. Initially, the concentration values are minimum in the numerical solution, but after covering some distance the numerical result attains slightly higher values in comparison to the analytical one. Both the results attain the minimum level at a particular point, after that it takes the reverse pattern in the case of exponentially decreasing form of the velocity pattern.
Figure 4
Figure 4. Concentration distribution pattern for the exponentially decreasing form of the velocity pattern for the sand medium.
View largeDownload slide

Concentration distribution pattern for the exponentially decreasing form of the velocity pattern for the sand medium.

Figure 4
Figure 4. Concentration distribution pattern for the exponentially decreasing form of the velocity pattern for the sand medium.
View largeDownload slide

Concentration distribution pattern for the exponentially decreasing form of the velocity pattern for the sand medium.

Close modal
The contaminant concentration pattern for the asymptotic form and algebraic sigmoid form of velocity patterns are almost similar and are shown in Figures 5 and 6, respectively. The contaminant concentration increases slowly with respect to time and decreases rapidly with respect to distance. The effect of the Peclet number is shown in Figure 7. The Peclet number physically measures the relative magnitude of advection versus dispersion. The contaminant concentration increases with the increasing Peclet number and decreases with the decreasing Peclet number. Both patterns with different Peclet numbers reach minimum concentration with distance. The effect of the Peclet number on the solute concentration is observed with respect to the distance. For high Peclet number, the concentration level takes a minimum distance to reach its minimum concentration, whereas for low Peclet number it takes more distance for the same value. The above discussion is for time domain for the analytical solution given in Equation (18).
Figure 5
Figure 5. Concentration distribution pattern for asymptotic velocity with average porosity of the medium.
View largeDownload slide

Concentration distribution pattern for asymptotic velocity with average porosity of the medium.

Figure 5
Figure 5. Concentration distribution pattern for asymptotic velocity with average porosity of the medium.
View largeDownload slide

Concentration distribution pattern for asymptotic velocity with average porosity of the medium.

Close modal
Figure 6
Figure 6. Concentration distribution pattern for algebraic sigmoid velocity with average porosity of the medium.
View largeDownload slide

Concentration distribution pattern for algebraic sigmoid velocity with average porosity of the medium.

Figure 6
Figure 6. Concentration distribution pattern for algebraic sigmoid velocity with average porosity of the medium.
View largeDownload slide

Concentration distribution pattern for algebraic sigmoid velocity with average porosity of the medium.

Close modal
Figure 7
Figure 7. Concentration distribution pattern for algebraic sigmoid velocity with different values of Peclet number with the average porosity of the sand medium.
View largeDownload slide

Concentration distribution pattern for algebraic sigmoid velocity with different values of Peclet number with the average porosity of the sand medium.

Figure 7
Figure 7. Concentration distribution pattern for algebraic sigmoid velocity with different values of Peclet number with the average porosity of the sand medium.
View largeDownload slide

Concentration distribution pattern for algebraic sigmoid velocity with different values of Peclet number with the average porosity of the sand medium.

Close modal
A similar type of concentration pattern was found for the case of the sinusoidal form of unsteady velocity pattern, as shown in Figure 8. The concentration is high in the aquitard as compared to the aquifer and ultimately goes to a minimum harmless concentration with respect to distance and time.
Figure 8
Figure 8. Concentration distribution for sinusoidal velocity for the time with average porosity of the medium.
View largeDownload slide

Concentration distribution for sinusoidal velocity for the time with average porosity of the medium.

Figure 8
Figure 8. Concentration distribution for sinusoidal velocity for the time with average porosity of the medium.
View largeDownload slide

Concentration distribution for sinusoidal velocity for the time with average porosity of the medium.

Close modal
Figures 9 and 10 compare the concentration patterns for the asymptotic form and the algebraic sigmoid form of unsteady velocity patterns, respectively. From these figures it is observed that the concentration pattern attains a maximum peak value at a certain distance and thereafter it begins to decrease with respect to distance. The concentration pattern initially decreases with the increasing value of time and after covering some distance of approximately 0.005 for the case of aquifer (i.e., sand) and 0.012 for the case of aquitard (i.e., clay) the concentration takes a reverse pattern with respect to time and goes to a minimum harmless concentration. The concentration pattern takes on a lesser peak value in the case of asymptotic form and algebraic sigmoid form than for the exponentially decreasing and sinusoidal form of the velocity pattern for time domain ⁠.
Figure 9
Figure 9. Concentration distribution for asymptotic velocity for time with average porosity of the medium.
View largeDownload slide

Concentration distribution for asymptotic velocity for time with average porosity of the medium.

Figure 9
Figure 9. Concentration distribution for asymptotic velocity for time with average porosity of the medium.
View largeDownload slide

Concentration distribution for asymptotic velocity for time with average porosity of the medium.

Close modal
Figure 10
Figure 10. Concentration distribution for the algebraic sigmoid velocity for time with average porosity of the medium.
View largeDownload slide

Concentration distribution for the algebraic sigmoid velocity for time with average porosity of the medium.

Figure 10
Figure 10. Concentration distribution for the algebraic sigmoid velocity for time with average porosity of the medium.
View largeDownload slide

Concentration distribution for the algebraic sigmoid velocity for time with average porosity of the medium.

Close modal
Figure 11 exhibits contaminant concentration for different values of the Peclet number with the algebraic sigmoid form of unsteady velocity pattern for time domain ⁠. The solute concentration value increases with the increasing value of the Peclet number and the peak of the concentration is reduced with the increasing time.
Figure 11
Figure 11. Concentration distribution pattern for algebraic sigmoid velocity for time with different values of Peclet numbers with average porosity of the sand medium.
View largeDownload slide

Concentration distribution pattern for algebraic sigmoid velocity for time with different values of Peclet numbers with average porosity of the sand medium.

Figure 11
Figure 11. Concentration distribution pattern for algebraic sigmoid velocity for time with different values of Peclet numbers with average porosity of the sand medium.
View largeDownload slide

Concentration distribution pattern for algebraic sigmoid velocity for time with different values of Peclet numbers with average porosity of the sand medium.

Close modal
Figure 12 predicts for the sinusoidal form of the velocity pattern for the averaging porosity of sand. It also follows the same type of nature as in Figure 4. The concentration value attains its minimum level in a short distance for the case of sinusoidal form of the velocity pattern in comparison to exponentially decreasing velocity pattern. From Figures 4 and 12, we observed that the decreasing nature of the concentration distribution patterns is faster for the sinusoidal form of velocity pattern (which shows the nature of groundwater contamination in the tropical region in which the fluctuation behavior of groundwater recharge is shown) in comparison to the exponentially decreasing form of velocity (Singh & Singh 2001; Jain et al. 2007) in which contaminant concentration follows a decreasing nature.
Figure 12
Figure 12. Concentration distribution pattern for the sinusoidal form of the velocity pattern for the sand medium.
View largeDownload slide

Concentration distribution pattern for the sinusoidal form of the velocity pattern for the sand medium.

Figure 12
Figure 12. Concentration distribution pattern for the sinusoidal form of the velocity pattern for the sand medium.
View largeDownload slide

Concentration distribution pattern for the sinusoidal form of the velocity pattern for the sand medium.

Close modal
Figures 13 and 14 predict for the averaging porosity of the clay. After the rainy season the contamination fluctuation in groundwater is shown from the asymptotic and algebraic sigmoid forms of velocity patterns (Aral & Liao 1996; Singh & Kumari 2014). The concentration distribution pattern for the asymptotic form of the velocity pattern for the clay is shown in Figure 13. The concentration values are minimum for the numerical result in comparison to the analytical one and beyond some distance it takes a reverse pattern. Both the analytical and numerical concentration values attain their minimum level at a particular position. A similar type of nature of the contaminant distribution pattern was observed for the algebraic sigmoid form of the velocity pattern depicted in Figure 14.
Figure 13
Figure 13. Concentration distribution pattern for the asymptotic form of the velocity pattern for the clay medium.
View largeDownload slide

Concentration distribution pattern for the asymptotic form of the velocity pattern for the clay medium.

Figure 13
Figure 13. Concentration distribution pattern for the asymptotic form of the velocity pattern for the clay medium.
View largeDownload slide

Concentration distribution pattern for the asymptotic form of the velocity pattern for the clay medium.

Close modal
Figure 14
Figure 14. Concentration distribution pattern for the algebraic sigmoid form of the velocity pattern for the clay medium.
View largeDownload slide

Concentration distribution pattern for the algebraic sigmoid form of the velocity pattern for the clay medium.

Figure 14
Figure 14. Concentration distribution pattern for the algebraic sigmoid form of the velocity pattern for the clay medium.
View largeDownload slide

Concentration distribution pattern for the algebraic sigmoid form of the velocity pattern for the clay medium.

Close modal

In this paper, the RMS error was used to check the validity of numerical solution against the analytical one, as shown in Tables 1 and 2. The two parameters, and ⁠, play an important role to investigate the performance of the numerical solution. In the explicit finite difference scheme is restricted under the stability condition. Thus, in this present study, the accuracy was investigated by selecting different mesh sizes. The RMS error was investigated for ⁠, 0.05, and 0.07 for the particular time period 20 years in the time domain for the exponential decreasing and the asymptotic form of the velocity patterns, and 30 years with within the sinusoidal and algebraic sigmoid form of the velocity patterns in sand and clay media, respectively.

Table 1

The RMS values for the sand medium at particular 20 year in 0 < T ≤ TP

DistanceAnalytical resultNumerical result
ΔZ = 0.02ΔZ = 0.05ΔZ = 0.07
Case i: For exponential decreasing form of the velocity pattern 
 0.0019 0.5125 0.0119 0.0149 0.0169 
 0.0055 0.3188 0.0159 0.0248 0.0308 
 0.0093 0.1736 0.0198 0.0348 0.0448 
 0.0130 0.0837 0.0238 0.0448 0.0588 
 0.0167 0.0376 0.0278 0.0548 0.0727 
 RMS error  0.2719 0.2665 0.2634 
Case ii: For asymptotic form of the velocity pattern 
 0.0019 0.3613 0.0120 0.0150 0.0170 
 0.0055 0.0687 0.0160 0.0250 0.0310 
 0.0093 0.0130 0.0200 0.0350 0.0449 
 0.0130 0.0094 0.0240 0.0449 0.0589 
 0.0167 0.0091 0.0280 0.0549 0.0729 
 RMS error  0.1583 0.1585 0.1596 
DistanceAnalytical resultNumerical result
ΔZ = 0.02ΔZ = 0.05ΔZ = 0.07
Case i: For exponential decreasing form of the velocity pattern 
 0.0019 0.5125 0.0119 0.0149 0.0169 
 0.0055 0.3188 0.0159 0.0248 0.0308 
 0.0093 0.1736 0.0198 0.0348 0.0448 
 0.0130 0.0837 0.0238 0.0448 0.0588 
 0.0167 0.0376 0.0278 0.0548 0.0727 
 RMS error  0.2719 0.2665 0.2634 
Case ii: For asymptotic form of the velocity pattern 
 0.0019 0.3613 0.0120 0.0150 0.0170 
 0.0055 0.0687 0.0160 0.0250 0.0310 
 0.0093 0.0130 0.0200 0.0350 0.0449 
 0.0130 0.0094 0.0240 0.0449 0.0589 
 0.0167 0.0091 0.0280 0.0549 0.0729 
 RMS error  0.1583 0.1585 0.1596 
View Large
Table 2

RMS values for the clay medium with their averaging porosity at particular 30 years for T > Tp

DistanceAnalytical resultNumerical result
ΔZ = 0.02ΔZ = 0.05ΔZ = 0.07
Case i: For the sinusoidal form of the velocity pattern 
 0.0028 0.1094 0.0118 0.0148 0.0168 
 0.0083 0.2425 0.0158 0.0248 0.0307 
 0.0138 0.2524 0.0198 0.0347 0.0447 
 0.0193 0.1869 0.0238 0.0447 0.0586 
 0.0248 0.1135 0.0278 0.0546 0.0725 
 RMS error  0.1726 0.1596 0.1514 
Case ii: For the algebraic sigmoid form of the velocity pattern 
 0.0028 0.1001 0.0119 0.0149 0.0169 
 0.0083 0.1608 0.0159 0.0249 0.0309 
 0.0138 0.0868 0.0199 0.0349 0.0449 
 0.0193 0.0301 0.0239 0.0449 0.0589 
 0.0248 0.0126 0.0279 0.0549 0.0728 
 RMS error  0.0818 0.0780 0.0774 
DistanceAnalytical resultNumerical result
ΔZ = 0.02ΔZ = 0.05ΔZ = 0.07
Case i: For the sinusoidal form of the velocity pattern 
 0.0028 0.1094 0.0118 0.0148 0.0168 
 0.0083 0.2425 0.0158 0.0248 0.0307 
 0.0138 0.2524 0.0198 0.0347 0.0447 
 0.0193 0.1869 0.0238 0.0447 0.0586 
 0.0248 0.1135 0.0278 0.0546 0.0725 
 RMS error  0.1726 0.1596 0.1514 
Case ii: For the algebraic sigmoid form of the velocity pattern 
 0.0028 0.1001 0.0119 0.0149 0.0169 
 0.0083 0.1608 0.0159 0.0249 0.0309 
 0.0138 0.0868 0.0199 0.0349 0.0449 
 0.0193 0.0301 0.0239 0.0449 0.0589 
 0.0248 0.0126 0.0279 0.0549 0.0728 
 RMS error  0.0818 0.0780 0.0774 
View Large

In both the tables, was fixed. Tables 1 and 2 were tabulated for the RMS error in the aquifer (i.e., sand) and aquitard (i.e., clay) for four different types of the velocity patterns, respectively. The RMS error decreases with the increasing grid space for the exponential decreasing form of the velocity pattern for the sand medium, which was observed from Case i of Table 1. In the asymptotic form of the velocity pattern the RMS error increases with respect to the increasing grid space in Case ii of Table 1. The RMS error attains its minimum value with the increasing mesh size for Cases i and ii in Table 2 for the clay medium. The RMS error was evaluated for the accuracy of solution for the sinusoidal form of the velocity pattern tabulated in Case i in Table 2, and for the algebraic sigmoid form of the velocity pattern tabulated in Case ii of Table 2. In both the velocity patterns the result is more accurate for the maximum value of the mesh size, except in the case of the asymptotic form of the velocity pattern where the result is more accurate in the case of minimum value of the mesh size.

Listen
One-dimensional ADE in multilayer porous media was solved analytically using generalized integral transform technique by Liu et al. (1998), where the analytical solution was derived under arbitrary initial and boundary conditions. In this present paper, the authors have shown the validation of the model equation with the existing research work done by Liu et al. (1998). The analytical solution obtained in Equations (18) and (19) was computed for the same set of input data, except some parameters which had an effect on the solute transport modeling, which have been compared with the input values taken by Liu et al. (1998) as follows:
The following inputs are taken for validation purposes in this present paper. The concentration distribution pattern for Input (i) and (ii) are predicted for the different geological formations with their averaging porosity and shown in Figures 15 and 16 for different time domains. The concentration distribution pattern for a particular time, 10 years, for sand and clay medium are predicted and shown in Figure 15, and for 25 years is also predicted and shown in Figure 16. The concentration values in each of the positions in Input (i) are higher as compared to Input (ii) in both the media as observed in Figure 15. The clay medium attains maximum concentration values as compared to the sand medium, but both patterns ultimately follow minimum concentration values with their respective distance. After removing the source, contaminant concentration increases with distance up to a certain limit, then decreases with respect to the distance in both the medium for Input (i) and (ii), as seen in Figure 16. The concentration pattern attains its maximum peak in the case of the aquitard (i.e., clay) as compared to the aquifer (i.e., sand) in both Input (i) and (ii). Initially, minimum concentration values are attained for Input (i) as compared with Input (ii), but after covering a certain distance the reverse pattern is observed. At the end of the position, minimum concentration values are attained by Input (ii), which have been taken in the present paper for illustration.
Figure 15
Figure 15. Concentration distribution pattern with different input for time with average porosity of the medium.
View largeDownload slide

Concentration distribution pattern with different input for time with average porosity of the medium.

Figure 15
Figure 15. Concentration distribution pattern with different input for time with average porosity of the medium.
View largeDownload slide

Concentration distribution pattern with different input for time with average porosity of the medium.

Close modal
Figure 16
Figure 16. Concentration distribution with the different input for time with average porosity of the medium.
View largeDownload slide

Concentration distribution with the different input for time with average porosity of the medium.

Figure 16
Figure 16. Concentration distribution with the different input for time with average porosity of the medium.
View largeDownload slide

Concentration distribution with the different input for time with average porosity of the medium.

Close modal
Listen

Employing the concept of linear isotherms, analytical solutions for the ADE with respect to the solid and liquid phase are derived. The mixed type boundary condition is employed at the source in the splitting time domain. The contaminant concentration patterns for different types of velocity patterns are evaluated. The following conclusions can be drawn from this study:

  1. The impact of contaminant concentration for linear isotherms with the distribution coefficient is significantly observed in the splitting time domain for different velocity patterns, such as exponentially decreasing, sinusoidally varying, algebraic, and asymptotic forms.

  2. The contaminant concentration values depend upon the decreasing or increasing values of the zero-order production term and first-order decay rate coefficient.

  3. The contaminant concentration distribution behavior is predicted for different geological formations in two time domains, i.e., and ⁠.

  4. Comparison of the analytical result with the numerical result is taken into account. Accuracy of the solution is significantly observed by using RMS error. Truncation error of various parameters is also explored, which causes the inconsistency among analytical and numerical results.

  5. The validation of the model is made with the result of an existing solution given by Liu et al. (1998) and the same trend for contaminant concentration was found.

Listen

The authors are grateful to the Indian School of Mines, Dhanbad for providing financial support to PhD candidate under the ISMJRF scheme. The authors are also grateful to the editor and reviewers for the comments which helped improve the quality of the paper.

Ataie-Ashtiani
B.
Lockington
D. A.
Volker
R. E.
1996
Numerical correction for finite-difference solution of the advection-dispersion equation with reaction
.
J. Contam. Hydrol.
23
,
149
–
156
.
Google Scholar
Crossref
Search ADS
 
Ataie-Ashtiani
B.
Lockington
D. A.
Volker
R. E.
1999
Truncation errors in finite difference models for solute transport equation with first-order reaction B
.
J. Contam. Hydrol.
35
,
409
–
428
.
Google Scholar
Crossref
Search ADS
 
Balla
K.
Keri
G.
Rapcsak
T.
2002
Pollution of underground water – a computational case study using a transport model
.
J. Hydroinform.
4
(
4
),
255
–
263
.
Google Scholar
Crossref
Search ADS
 
Batu
V.
2006
Applied Flow and Solute Transport Modeling in Aquifers: Fundamental Principles and Analytical and Numerical Methods
.
CRC Press
,
Boca Raton, FL
,
USA
.
Google Scholar
 
Bear
J.
Verruijt
A.
1987
Modeling Groundwater Flow and Pollution
.
Reidel
,
Dordrecht
,
The Netherlands
.
Google Scholar
 
Bedient
P. B.
Rifai
H. S.
Newell
C. J.
1999
Ground Water Contamination: Transport and Remediation
, 2nd edn.
Prentice-Hall PTR
,
Upper Saddle River, NJ
,
USA
.
Google Scholar
 
Carnahan
B.
Luther
H. A.
Wilkes
J. O.
1969
Applied Numerical Methods
.
Wiley
,
New York
,
USA
, pp.
429
–
530
.
Google Scholar
 
Chrysikopolous
C. V.
Roberts
P. V.
Kitanidis
P. K.
1990
One-dimensional solute transport in porous media with partial well-to-well recirulation: application to field experiment
.
Water Resour. Res.
26
(
6
),
1189
–
1195
.
Google Scholar
 
Crank
J.
1975
The Mathematics of Diffusion
, 2nd edn.
Oxford University Press
,
London
,
UK
.
Google Scholar
 
Freeze
R. A.
Cherry
J. A.
1979
Groundwater
.
Prentice-Hall
,
Upper Saddle River, NJ
,
USA
.
Google Scholar
 
Fried
J. J.
Combarnous
M. A.
1971
Dispersion in porous media
.
Adv Hydrosci.
7
,
169
–
281
.
Google Scholar
Crossref
Search ADS
 
Gao
G.
Zhan
H.
Feng
S.
Fu
B.
Huang
G.
2012
A mobile-immobile model with an asymptotic scale-dependent dispersion function
.
J. Hydrol.
424–425
,
172
–
183
.
Google Scholar
Crossref
Search ADS
 
Hill
D. J.
Minsker
B. S.
Valocchi
A. J.
Babovic
V.
Keijzer
M.
2007
Upscaling models of solute transport in porous media through genetic programming
.
J. Hydroinform.
9
(
4
),
251
–
266
.
Google Scholar
Crossref
Search ADS
 
Hunt
B.
1978
Dispersion sources in uniform groundwater flow
.
J. Hydraul. Div.
HY1
,
74
–
85
.
Google Scholar
 
Jain
S. K.
Agarwal
P. K.
Singh
V. P.
2007
Hydrology and Water Resources of India
.
Springer
,
Dordrecht
,
The Netherlands
, p.
87
.
Google Scholar
 
Mickley
H. S.
Sherwood
T. K.
Reed
C. E.
1957
Applied Mathematics in Chemical Engineering.
McGraw-Hill
,
New York
,
USA
.
Google Scholar
 
Neville
C. J.
Ibaraki
M.
Sudicky
E. A.
2000
Solute transport with multiprocess nonequilibrium: a semi analytical solution approach
.
J. Cont. Hydrol.
44
(
2
),
141
–
159
.
Google Scholar
Crossref
Search ADS
 
Notodarmojo
S.
Ho
G. E.
Scott
W. D.
Davis
G. B.
1991
Modelling phosphorus transport in soils and groundwater with two-consecutive reactions
.
Water Res.
25
(
10
),
1205
–
1216
.
Google Scholar
Crossref
Search ADS
 
Rausch
R.
Schafer
W.
Therrien
R.
Wagner
C.
2005
Solute Transport Modelling: An Introduction to Models and Solution Strategies
.
Gebruder Borntrager
,
Berlin
,
Germany
.
Google Scholar
 
Richard
A. S.
Jr
1988
Small grid testing of finite difference transport schemes
.
J. Hydraul. Eng.
114
(
3
),
329
–
336
.
Google Scholar
Crossref
Search ADS
 
Singh
M. K.
Kumari
P.
2014
Contaminant concentration prediction along unsteady groundwater flow
.
Modell. Simul. Diff. Process.
XII
,
257
–
276
.
Google Scholar
Crossref
Search ADS
 
Singh
P.
Singh
V. P.
2001
Snow and Glacier Hydrology
.
Kluwer Academic Publishers
,
Amsterdam
,
The Netherlands
, p.
78
.
Google Scholar
 
Smith
G. D.
1978
Numerical Solution of Partial Differential Equations: Finite Difference Methods
, 2nd edn.
Oxford University Press
,
Oxford
,
UK
.
Google Scholar
 
Zheng
C.
Bennett
G. D.
2002
Applied Contaminant Transport Modeling
, 2nd edn.
Wiley
,
New York
,
USA
.
Google Scholar
 
© IWA Publishing 2016

Supplementary data

APPENDIX- pdf file