Impact of Pumping Rate on Contaminant Transport in Groundwater—A Numerical Study

: Public supply wells are commonly considered one of the most signiﬁcant sources of freshwater on Earth. Therefore, potential well water contamination can conceivably be regarded as a crucial issue that is closely correlated with both environmental protection and water demand. In the present study, a three-dimensional numerical model is developed to simulate unsteady and spatially varying groundwater ﬂow, along with contaminant migration. Besides, the proposed model is capable of investigating well water quality by the change of the wells’ pumping rates. The developed model uses a ﬁnite-volume time splitting numerical technique to solve governing groundwater ﬂow and soluble contaminant transport equations. Comparison of the numerical simulation results with analytical solutions, as well as experimental and ﬁeld data, clearly demonstrates the satisfactory performance of the present model. The fundamental aim of the study is to evaluate the effect of pumping rate and its variations on pollution migration through saturated porous media. To meet this purpose, contaminant concentrations and contaminants’ travel time were studied under different pump ﬂow rate conditions. The modeling results revealed that choosing an optimum range for the pumping rate increases contaminant travel time and reduces aquifer vulnerability. group wells as did well case. this study, a combined numerical model was presented to make accurate and pre- cise simulations of both groundwater flow and contaminant transport in transient condi- tions through saturated porous media. the numerical model, unsteady governing equations were numerically solved by using the finite volume time splitting scheme. The performance and accuracy of the numerical models were successfully examined by sev- eral test cases, which indicated the current numerical model has enough capability to pre- dict radial groundwater toward a pumping well and contaminant migration through aquifers, including advection, dispersion, and diffusion The effect of pumping rate on solute migration was evaluated by the numerical model, and the following results drawn. The reduction in pumping rate was accom- panied by a decrease in flow velocity, which decreased both the advection and dispersion rates, as well as the growth of contaminant travel time. Discontinuing the pumping brought about a dramatic drop in flow velocity, which induced a remarkable rise in travel time, along with an increase in pollution concentrations downstream.


Introduction
The significance of groundwater management has increased in recent years, due to increasing population growth rates and water needs. Groundwater pumping aids in meeting the ever-increasing demand for water on a sustainable basis although both the quality and quantity of groundwater resources have degraded in recent years [1]. The growing public concern regarding pollutant emissions includes considering groundwater contamination as a crucial subject in hydrogeology; thus, there has been a dramatic interest increase in the assessment and prediction of contaminant migration in groundwater [2][3][4].
The aquifer vulnerability is directly related to the hydraulic conductivity of porous media, the depth to the water table, the groundwater recharge, pollution source, and groundwater flow conditions [5,6]. Travel time and pollutant concentration are the overriding factors for determining the vulnerability of aquifers.
In recent years, various methods have been used by researchers to investigate contaminant transport through porous media. Both analytical and numerical techniques have frequently been used to study groundwater pollution. Bessel function expansion, Green's Function, and Eigen function expansions have been found as the most suitable analytical techniques [7][8][9]. Classical numerical methods, including finite difference, finite element, and finite volume schemes, have been the common approaches for groundwater modeling in recent years [10][11][12][13]. In several investigations, both experimental and numerical methods have been used to study contaminant transport [14]. In recent years, some novel methods have been proposed to overcome the limitations of the classical numerical methods [15][16][17][18][19].
In previous years, detection and remediation of well water contamination have posed novel challenges to groundwater resource managers. As a result, the effective management of groundwater quality has become one of the most complex problems faced by modern groundwater researchers. Halford et al. discussed improvements in pumped water quality from public supply wells in Antelope Valley (California) obtained through grouting of selected intervals to reduce the entry of arsenic-bearing water [20]. Johnson et al. examined the impact of several pumping scenarios on leakage from an unconfined aquifer to a confined aquifer and concluded that a single inactive multi-aquifer well can contribute nearly 10% of total public water supply well flow over a wide range of pumping rates [21].
Several solutions have been developed for addressing the interaction between the pumping rate and water quality of a given well, using analytical and numerical techniques. Park and Aral solved a multi-objective optimization problem using a genetic algorithm to determine pumping rates and well locations for the prevention of saltwater intrusion, while satisfying desired discharge rates in aquifers [22]. Further, Kalwij and Peralta presented a simulation/optimization model for determining optimal pumping strategies, considering the uncertainty of the aquifer's parameters [23]. Both groundwater flow and pollution migration through long-screened wells were simulated by Zinn and Konikow [24]. Their findings indicated that intra-borehole flow had a direct impact on the groundwater flow within a conceptual aquifer system. Additionally, Sharief et al. used a coupled finite element method, which was embedded with the genetic algorithm (GA) to find an optimal pumping pattern for the remediation of aquifers [25]. They concluded that using the pump and treat method can remediate the pollution plume.
Later, the impact of seasonal changes on both the quality and age of water pumped from public supply wells was investigated by Bexfield and Jurgens [26]. They demonstrated that either abatement or cessation of pumping from public supply wells causes contaminant migration through intra-borehole flow. In another study, Cyriac and Rastogi investigated an optimal pumping policy using a coupled finite element-particle swarm simulationoptimization model [27]. Their model considered the collective drawdown of the wells to find the maximum allowable pumping discharge, along with the location of the wells to mitigate pollution. Rodriguez-Pretelin and Nowak developed a management scheme fulfilling the groundwater demand and decreasing the impact of transient flow on the abstraction zone through optimally redistributing pumping rates within well fields [28]. Finally, Dey and Prakash carried out density-dependent simulations to design a pumping management strategy in a coastal aquifer [29]. Their method determined the maximum pumping rate without seawater intrusion.
It is of great importance to accurately predict the migration path and pollution plume during different pumping rates, in the management and protection of groundwater resources. In recent studies, such as studies by Medici et al. [30] and Medici and West [31], groundwater flow and aquifer vulnerability have been studied at steady-state conditions. Steady-state flow conditions are fairly rare, and unsteady conditions can occur as a result of daily and seasonal fluctuations as well as time-dependent hydraulic stresses, such as pumping or recharge. Consequently, transient flow scenarios should be incorporated, and the influence of unsteady-state groundwater flow on contaminant fate should be considered.
The objective of this study is to investigate the effect of different pumping rates on the contaminant plume, by the proposed numerical model, which can simulate the transient groundwater flow and pollution transport through saturated porous media. The model can be used as an effective tool for predicting pollution migration and making decisions in the management of water resource systems.

Governing Groundwater Flow Equation
The equation of transient groundwater flow through non-homogeneous anisotropic porous media, as presented by Darcy's Law, is derived by combining a continuity equation with Darcy's Law. This equation is known as the Boussinesq equation. The groundwater flow in saturated porous media is presented by the following equation [32]: where K LT −1 is the hydraulic conductivity matrix, which includes the hydraulic conductivity of the soil in different directions; h [L] is the hydraulic head; S s L −1 is specific storage; and t [T] represents time.

Governing Contaminant Transport Equation
The transport of solute in porous media, as in groundwater, can be described by three processes: advection, hydrodynamic dispersion, and physical-chemical or biochemical reactions. The transport of contaminants in flowing groundwater is derived from the principle of conservation of mass. A general form of the solute-transport equation is as follows: where C(x, y, z, t) is the contaminant concentration ML −3 ; D is the coefficient of hydrodynamic dispersion L 2 T −1 ; V is the seepage velocity LT −1 ; λ is the contaminant decay rate; and R d is the retardation factor that is used for considering the difference between flow seepage and solute contaminant velocity, as a result of physical, biological, and chemical mechanisms. Contaminant dilution due to hydrodynamic dispersion is represented by the first term on the right side of Equation (2). The second term of the equation constitutes advective transport, describing the movement of solutes in flowing groundwater. Note that in this research, a non-reactive solute contaminant was considered to simulate pollution migration; therefore, the third term on the right side of Equation (2) was neglected. The hydrodynamic dispersion coefficients in different directions are given as where D 1 , D 2 , and D 3 are hydrodynamic dispersion coefficients in three Cartesian coordinate axes. α L , α TH , and α TV are longitudinal, horizontal transverse, and vertical transverse dispersivity, respectively. v 1 , v 2 , and v 3 are velocity vectors in the Cartesian coordinate axes. D * is the coefficient of molecular diffusion, and |V| is the average velocity vector. The longitudinal dispersivity is proportional to the measurement scale. For porous media, the flow direction dispersivity is as follows [33]: where L [L] is the flow distance (path). Both transverse dispersivities are generally smaller than longitudinal, and in most cases, they are 1/10 to 1/100 of the longitudinal dispersivity [34]. In recent years, some research has been conducted to get acute insight into the relationship between longitudinal dispersivity and transverse dispersivities [35].

Numerical Solution
The proposed model uses a finite-volume time splitting method approach to solve the governing equations in porous media. The scalar quantities, such as contaminant concentration and hydraulic head, are defined at the center of the cells, while the velocity

Numerical Solution
The proposed model uses a finite-volume time splitting method approach to solve the governing equations in porous media. The scalar quantities, such as contaminant concentration and hydraulic head, are defined at the center of the cells, while the velocity components are located at the faces, as shown in Figure 1. The flux passing through each cell face can be regarded as the most important issue in FVM.
where the term on the left-hand side demonstrates the mass inside the cell at a new time step, and the first term on the right-hand side represents the mass inside the cell at the old time step. Another technique used in the present model is the time splitting method. The key feature of this method is solving the governing equations in different steps. This method not only reduces computation cost, but also elevates the numerical model's precision. The time splitting algorithm, proposed by Yanenko [36], has widely been employed to numerically solve the governing partial differential equations. The aforementioned methods are employed to solve both groundwater and contaminant transport equations.

Groundwater Equation
Finite volume discretization can be written for the Boussinesq equation, giving where (K ) / and (K ) / are fluxes passing from the downstream and the upstream of the cell faces, respectively. A second-order (in space) estimation for space derivation at the cell face can be written as The calculated flux for each of the cell faces is added to the mass of the downstream cell and subtracted from the upstream cell at the same time, whereby conservativity is guaranteed. The change of mass inside the cell i during a specific time interval ∆t can be calculated as mass n+1 where the term on the left-hand side demonstrates the mass inside the cell i at a new time step, and the first term on the right-hand side represents the mass inside the cell i at the old time step. Another technique used in the present model is the time splitting method. The key feature of this method is solving the governing equations in different steps. This method not only reduces computation cost, but also elevates the numerical model's precision. The time splitting algorithm, proposed by Yanenko [36], has widely been employed to numerically solve the governing partial differential equations. The aforementioned methods are employed to solve both groundwater and contaminant transport equations.

Groundwater Equation
Finite volume discretization can be written for the Boussinesq equation, giving are fluxes passing from the downstream and the upstream of the cell i faces, respectively. A second-order (in space) estimation for space derivation at the cell face can be written as In the next part, the time splitting method is used to solve the discretized groundwater equation. The solution procedure consists of one main step, which has three stages. The three stages include calculating the fluxes in three Cartesian coordinate axes. Full details about using the finite-volume time splitting method, as well as simulation of the free surface by an iterative algorithm, can be found in the study by Ahmadi [17,19].

Contaminant Transport Equation
Employing the time splitting algorithm, the contaminant transport equation can be split into two main steps. The first step, which consists of three stages, solves the advection term, and the second step solves the dispersion term in three stages. In the first step, the contaminant concentration is computed in the time step n + 1/2 by solving the advection term. The second step is to compute the contaminant concentration by solving the dispersion term in the time step n + 1.

Advection Term
Finite volume discretization can be written for the advection term, giving are fluxes passing from downstream and the upstream of the cell i faces, respectively. Employing the FVM for the advection term, flux is presented as This step includes three stages in all directions, and the advection terms in all the stages are solved using the Fromm scheme, as an explicit second-order upwind scheme.

Dispersion Term
Finite volume discretization can be written for the dispersion term, giving are the upstream and downstream fluxes of the cell i faces, respectively.
This term is composed of three stages too, and the flux based on FVM is Considering implicit solution procure, this step is solved in three stages. Note that more details about the techniques to solve advection/dispersion equations can be found in Namin et al. [37].

Result
In this section, some test cases are presented to validate both the groundwater flow and contaminant transport models, along with assessing the possible advantages or limitations of the proposed numerical model. In the first part, the results of the numerical model for simulating radial flow toward a pumping well are investigated. In the next part, the performance of the numerical model for predicting contaminant migration is evaluated.

Verification of the Flow Model
Field experimental and analytical approach results, reported by Arshad et al. [38], were used as a test case to validate the performance of the numerical model to predict groundwater flow. According to the test case, sandy soil with five boreholes at distinct distances was simulated with the numerical model ( Figure 2). The initial condition was obtained from a steady-state condition without pumping. Then, pumping at the rate of 0.042 m 3 /s was started, and results were obtained at a set of time intervals for the pumping well and all boreholes. The simulated drawdowns by the proposed model are in accordance with field data and an analytical approach, commonly known as the Theis equation [39] ( Figure 3).

Verification of the Flow Model
Field experimental and analytical approach results, reported by Arshad et al. [38], were used as a test case to validate the performance of the numerical model to predict groundwater flow. According to the test case, sandy soil with five boreholes at distinct distances was simulated with the numerical model ( Figure 2). The initial condition was obtained from a steady-state condition without pumping. Then, pumping at the rate of 0.042 m 3 /s was started, and results were obtained at a set of time intervals for the pumping well and all boreholes. The simulated drawdowns by the proposed model are in accordance with field data and an analytical approach, commonly known as the Theis equation [39] (Figure 3).  In order to evaluate the performance of the proposed numerical model for predicting the groundwater flow, mean absolute error (MAE), and root mean square error (RMSE) statistical indices were used in this research:

Verification of the Flow Model
Field experimental and analytical approach results, reported by Arshad et al. [38], were used as a test case to validate the performance of the numerical model to predict groundwater flow. According to the test case, sandy soil with five boreholes at distinct distances was simulated with the numerical model ( Figure 2). The initial condition was obtained from a steady-state condition without pumping. Then, pumping at the rate of 0.042 m 3 /s was started, and results were obtained at a set of time intervals for the pumping well and all boreholes. The simulated drawdowns by the proposed model are in accordance with field data and an analytical approach, commonly known as the Theis equation [39] (Figure 3).  In order to evaluate the performance of the proposed numerical model for predicting the groundwater flow, mean absolute error (MAE), and root mean square error (RMSE) statistical indices were used in this research: In order to evaluate the performance of the proposed numerical model for predicting the groundwater flow, mean absolute error (MAE), and root mean square error (RMSE) statistical indices were used in this research: where H Si j is the jth value of the simulated head, and H Obj is the jth value of the observed head (field data). For the numerical model, MAE and RMSE values were found to be 0.19 m and 0.15, respectively. The comparison of both field experiments and analytical results with numerical simulation outcomes represents the perfect applicability of the numerical model.

Verification of the Transport Model
The performance of the developed finite-volume model to predict the contaminant transport was validated by an experiment conducted by Massabo et al. [40]. The apparatus of the experiment mainly consisted of a chamber filled with glass beads to simulate contaminant migration through a semi-infinite saturated porous medium. The physical model can be regarded as quasi-2-dimensional because of its trivial thickness regarding the other dimensions of the experimental model. The solute was introduced by a pulse-like injection, and it was illuminated by a diffuse backlight UV source ( Figure 4).
where is the jth value of the simulated head, and is the jth value of the observed head (field data). For the numerical model, and values were found to be 0.19 m and 0.15, respectively. The comparison of both field experiments and analytical results with numerical simulation outcomes represents the perfect applicability of the numerical model.

Verification of the Transport Model
The performance of the developed finite-volume model to predict the contaminant transport was validated by an experiment conducted by Massabo et al. [40]. The apparatus of the experiment mainly consisted of a chamber filled with glass beads to simulate contaminant migration through a semi-infinite saturated porous medium. The physical model can be regarded as quasi-2-dimensional because of its trivial thickness regarding the other dimensions of the experimental model. The solute was introduced by a pulselike injection, and it was illuminated by a diffuse backlight UV source (Figure 4). The analytical solution for simulating the contaminant transport through saturated porous media was presented by Ogata and Banks [41]. Comparisons among the numerical model, experimental data, and analytical solution are shown in Figures 5 and 6. The curves in these figures are representative of the plume characteristics in the middle longitudinal cross section (y = 0) and at a fixed transverse cross section in the middle part of the longitudinal domain (x = 0.14 m) at time 400 s. The analytical solution for simulating the contaminant transport through saturated porous media was presented by Ogata and Banks [41]. Comparisons among the numerical model, experimental data, and analytical solution are shown in Figures 5 and 6. The curves in these figures are representative of the plume characteristics in the middle longitudinal cross section (y = 0) and at a fixed transverse cross section in the middle part of the longitudinal domain (x = 0.14 m) at time 400 s.
The accuracy of the simulation model was investigated by the mentioned statistical indices, and MAE and RMSE values were found to be 0.16 m and 0.14 m, respectively. Evaluation of experimental data and numerical simulation outcomes revealed that the numerical results had an acceptable agreement with both experimental and analytical results. The accordance indicates that the proposed model can adequately predict solute migration.  The accuracy of the simulation model was investigated by the mentioned statistical indices, and and values were found to be 0.16 m and 0.14 m, respectively. Evaluation of experimental data and numerical simulation outcomes revealed that the numerical results had an acceptable agreement with both experimental and analytical results. The accordance indicates that the proposed model can adequately predict solute migration.

Numerical Implementation
The proposed numerical model aims to investigate the impact of pumping rates on both groundwater flow and contamination transport through saturated porous media. At the beginning step, the proposed model simulated radial groundwater flow toward a pumping well to find the steady-state condition of groundwater flow. The steady state was obtained provided that groundwater flow parameters, such as piezometric heads, velocity, and pressure, did not change at any single point in the porous media. Then the pumping rate was changed, altering the groundwater flow properties and making the radial flow unsteady. The unsteady-state groundwater flow was simulated again to find a steady-state radial groundwater flow towards the pumping well; furthermore, the duration, through which the unsteady groundwater flow had converted to the steady-state  The accuracy of the simulation model was investigated by the mentioned statistical indices, and and values were found to be 0.16 m and 0.14 m, respectively. Evaluation of experimental data and numerical simulation outcomes revealed that the numerical results had an acceptable agreement with both experimental and analytical results. The accordance indicates that the proposed model can adequately predict solute migration.

Numerical Implementation
The proposed numerical model aims to investigate the impact of pumping rates on both groundwater flow and contamination transport through saturated porous media. At the beginning step, the proposed model simulated radial groundwater flow toward a pumping well to find the steady-state condition of groundwater flow. The steady state was obtained provided that groundwater flow parameters, such as piezometric heads, velocity, and pressure, did not change at any single point in the porous media. Then the pumping rate was changed, altering the groundwater flow properties and making the radial flow unsteady. The unsteady-state groundwater flow was simulated again to find a steady-state radial groundwater flow towards the pumping well; furthermore, the duration, through which the unsteady groundwater flow had converted to the steady-state

Numerical Implementation
The proposed numerical model aims to investigate the impact of pumping rates on both groundwater flow and contamination transport through saturated porous media. At the beginning step, the proposed model simulated radial groundwater flow toward a pumping well to find the steady-state condition of groundwater flow. The steady state was obtained provided that groundwater flow parameters, such as piezometric heads, velocity, and pressure, did not change at any single point in the porous media. Then the pumping rate was changed, altering the groundwater flow properties and making the radial flow unsteady. The unsteady-state groundwater flow was simulated again to find a steady-state radial groundwater flow towards the pumping well; furthermore, the duration, through which the unsteady groundwater flow had converted to the steady-state condition, was simulated by the numerical model. The pollutant migration for each of the time steps, through and after the pumping rate change, was predicted. investigated by considering four different pumping rate stages. These stages include 100%, 70%, 30% of the constant pumping rate, as well as ceasing the pumping. The pumping rate was changed from 100% to another rate as soon as the contaminant was injected into the porous media. The contaminant was injected into the porous media during a 24-h period. It should be noted that the model used a variable grid size to achieve ideal efficiency for the simulation of both groundwater flow and contaminant transport through saturated porous media. The grid size was reduced by half around the well to obtain maximum precision and flexibility. condition, was simulated by the numerical model. The pollutant migration for each of the time steps, through and after the pumping rate change, was predicted. This concept is illustrated by two hypothetical examples with a simplified homogeneous aquifer system, with 250 45 45 m dimensions, including a single well or group of wells under the pumping rate of 0.2 m /s (17,280 m /day). The upstream and downstream conditions for the examples are assumed to be fixed heads, which are 45 and 43 m respectively (Figure 7). The effect of the pumping rate on contaminant transport was investigated by considering four different pumping rate stages. These stages include 100%, 70%, 30% of the constant pumping rate, as well as ceasing the pumping. The pumping rate was changed from 100% to another rate as soon as the contaminant was injected into the porous media. The contaminant was injected into the porous media during a 24h period. It should be noted that the model used a variable grid size to achieve ideal efficiency for the simulation of both groundwater flow and contaminant transport through saturated porous media. The grid size was reduced by half around the well to obtain maximum precision and flexibility. The input variables and aquifer parameters are listed in Table 1. The input variables and aquifer parameters are listed in Table 1.

Single Well
In the first example, a single well in a homogeneous aquifer with a pumping rate of 0.2 m 3 /s was used to study the effects of pumping rate on contaminant transport. The four stages of the pumping rate were simulated by the numerical model, and the pollution plume for each of the pumping rates was obtained. The time taken to reach steady-state groundwater flow conditions, due to the pumping rate change from the initial pumping rate (100%) to 70% and 30%, as well as to cease the pumping were 28.4, 33.2, and 42.3 h, respectively. As mentioned earlier, once the solute contaminant, which had a constant concentration, was injected into the porous media, the pumping rate was changed from the initial pumping rate to another rate. Figures 8 and 9 show the pollution plumes at the plane passing through the pumping well after 24 and 36 h, respectively.

Single Well
In the first example, a single well in a homogeneous aquifer with a pumping rate of 0.2 m /s was used to study the effects of pumping rate on contaminant transport. The four stages of the pumping rate were simulated by the numerical model, and the pollution plume for each of the pumping rates was obtained. The time taken to reach steady-state groundwater flow conditions, due to the pumping rate change from the initial pumping rate (100%) to 70% and 30%, as well as to cease the pumping were 28.4, 33.2, and 42.3 h, respectively. As mentioned earlier, once the solute contaminant, which had a constant concentration, was injected into the porous media, the pumping rate was changed from the initial pumping rate to another rate. Figures 8 and 9 show the pollution plumes at the plane passing through the pumping well after 24 and 36 h, respectively.   Figures 8 and 9 show that contaminant movement, as a result of groundwater flow, was a considerable process through the saturated porous media. This means the advection was the dominant process in the contaminant transport for all pumping rates; on the other hand, the ceased pumping case not only decelerated both advection and dispersion processes, but also made diffusion the predominant process in contaminant migration. Figure  10 shows the migration of the pollution under different pumping rates in the saturated porous media at the section passing through the well at = 10 m over time.  Figures 8 and 9 show that contaminant movement, as a result of groundwater flow, was a considerable process through the saturated porous media. This means the advection was the dominant process in the contaminant transport for all pumping rates; on the other hand, the ceased pumping case not only decelerated both advection and dispersion processes, but also made diffusion the predominant process in contaminant migration. Figure 10 shows the migration of the pollution under different pumping rates in the saturated porous media at the section passing through the well at Z = 10 m over time.
The large velocity gradients attributed to the pumping conditions brought about faster advection and dispersion of the pollution; on the contrary, ceasing pumping considerably affected the groundwater flow, and dramatically reduced the velocity field through porous media, inducing a remarkable decrease in both advection and dispersion rates. The large velocity gradients attributed to the pumping conditions brought about faster advection and dispersion of the pollution; on the contrary, ceasing pumping considerably affected the groundwater flow, and dramatically reduced the velocity field through porous media, inducing a remarkable decrease in both advection and dispersion rates.
In Figure 11, the effect of the pumping rate on both travel time and contaminant concentrations for a single well is presented. Travel time is the time spent for contaminant arrival time at wells, and it directly depends upon the advection. In fact, travel time is a primary measure of the aquifer-solute contact time. The travel time is obtained from numerical simulations of transient flow and transport in the soil media.  In Figure 11, the effect of the pumping rate on both travel time and contaminant concentrations for a single well is presented. Travel time is the time spent for contaminant arrival time at wells, and it directly depends upon the advection. In fact, travel time is a primary measure of the aquifer-solute contact time. The travel time is obtained from numerical simulations of transient flow and transport in the soil media.
In addition, the contaminant concentration at the well point varies with changing pumping rate. Figure 11 illustrates the pumping rate effect on the travel time, along with the maximum concentration at the well point over time. In addition, the contaminant concentration at the well point varies with changing pumping rate. Figure 11 illustrates the pumping rate effect on the travel time, along with the maximum concentration at the well point over time. Figure 11 indicates that as the single-well pumping rate reduces, the travel time negligibly increases, whereas ceasing pumping brings about a great growth in the travel time. Furthermore, by reducing the pumping rate from its initial rate to other rates, the peak of the contaminant concentration at the well point showed a gradual decline to reach its nadir; then it regrew. Decreasing the pumping rate increased the travel time for contaminants, whereby it prepared more time for contaminant dispersion. The difference in the travel time increase due to the reduction of the pumping rate was negligible; therefore, the effect of contaminant concentration reduction was slight. The greater the difference between the initial pumping rate and pumping rate reduction, the more time it took to reach the steady-state condition, with more advection and less dispersion. The final simulation of contaminant transport was related to ceasing pumping when the contaminant concentration at the point of the pumping well was greater than the pumping reduction stages. In this case, the velocity field significantly diminished; therefore, both advection and dispersion rates dwindled. Additionally, ceasing pumping made diffusion the dominant phenomenon in the contaminant transport process and triggered a travel time increase.

Group of Wells
In this section, the effect of the pumping rate for a group of two wells is investigated. The pumping rate of each well was 0.1 m /s, the sum of which equals the pumping rate of the last example. The time taken to reach steady-state groundwater flow, owing to the change of pumping rate from the initial pumping rate (100%) to 70%, 30%, and ceased pumping were 31.2, 38.5, and 54.1 h, respectively. The pollution plume at the section passing through the middle of both pumping wells for the four pumping stages after 24 h is shown in Figure 12. Note that the contaminant concentration in the group of wells case was similar to the single well case; however, the hydraulic conditions of the points were not the same.

Reduction of Pumping Rate (Percent)
Travel Time (hr) C/C0 Figure 11. The pumping rate effect on the travel time for single-well. Figure 11 indicates that as the single-well pumping rate reduces, the travel time negligibly increases, whereas ceasing pumping brings about a great growth in the travel time. Furthermore, by reducing the pumping rate from its initial rate to other rates, the peak of the contaminant concentration at the well point showed a gradual decline to reach its nadir; then it regrew. Decreasing the pumping rate increased the travel time for contaminants, whereby it prepared more time for contaminant dispersion. The difference in the travel time increase due to the reduction of the pumping rate was negligible; therefore, the effect of contaminant concentration reduction was slight. The greater the difference between the initial pumping rate and pumping rate reduction, the more time it took to reach the steady-state condition, with more advection and less dispersion. The final simulation of contaminant transport was related to ceasing pumping when the contaminant concentration at the point of the pumping well was greater than the pumping reduction stages. In this case, the velocity field significantly diminished; therefore, both advection and dispersion rates dwindled. Additionally, ceasing pumping made diffusion the dominant phenomenon in the contaminant transport process and triggered a travel time increase.

Group of Wells
In this section, the effect of the pumping rate for a group of two wells is investigated. The pumping rate of each well was 0.1 m 3 /s, the sum of which equals the pumping rate of the last example. The time taken to reach steady-state groundwater flow, owing to the change of pumping rate from the initial pumping rate (100%) to 70%, 30%, and ceased pumping were 31.2, 38.5, and 54.1 h, respectively. The pollution plume at the section passing through the middle of both pumping wells for the four pumping stages after 24 h is shown in Figure 12. Note that the contaminant concentration in the group of wells case was similar to the single well case; however, the hydraulic conditions of the points were not the same.  Figure 13 indicates that advection was the dominant process in contaminant transport for all pumping stages; nevertheless, in the ceased pumping stage, both advection and dispersion rates were reduced. The pollution migration at the section through the middle of both wells and at = 10 m under different pumping rates over time, is shown in Figure 13.  Figure 13 indicates that advection was the dominant process in contaminant transport for all pumping stages; nevertheless, in the ceased pumping stage, both advection and dispersion rates were reduced. The pollution migration at the section through the middle of both wells and at Z = 10 m under different pumping rates over time, is shown in Figure 13. Figure 13 indicates that advection was the predominant phenomenon in pumping stages, and diffusion was the dominant process in the ceased pumping stage. Figure 14 shows the travel time of the group of wells, along with the contaminant concentration at the point located between the pumping wells.
As the figure shows, the travel time and the peak concentration at the point between two pumping wells were similar to the trend of the single well peak concentration because of the aforementioned reasons.
The time taken to reach steady-state groundwater flow conditions under different pumping rates is presented in Table 2.   Figure 13 indicates that advection was the predominant phenomenon in pumping stages, and diffusion was the dominant process in the ceased pumping stage. Figure 14 shows the travel time of the group of wells, along with the contaminant concentration at the point located between the pumping wells.  As the figure shows, the travel time and the peak concentration at the point between two pumping wells were similar to the trend of the single well peak concentration because of the aforementioned reasons.
The time taken to reach steady-state groundwater flow conditions under different pumping rates is presented in Table 2. The case where a group of wells took a longer time to reach steady-state conditions, rather than the single well, was because of the interaction between the wells. A comparison between the pumping rate impact on both travel time and contaminant concentration for the single well and the group of wells is presented in Figure 15. The case where a group of wells took a longer time to reach steady-state conditions, rather than the single well, was because of the interaction between the wells. A comparison between the pumping rate impact on both travel time and contaminant concentration for the single well and the group of wells is presented in Figure 15. Based on Figure 15, the effect of the pumping rate change on the travel time, for both a single well and a group of wells, is the same. Therefore, the number of wells is not an overriding factor for determining contaminant travel time under various pumping rates. Conversely, it indicates that the peak of the contaminant concentration depends to a great extent on the number of wells, and it is intensified when the number of wells increases. The reason that the contaminant concentration in the group of wells condition increases is that, in this case, the pollution plume takes a longer time to reach the steady-state condition in comparison to the single well; consequently, the pollutant not only further conveys, but also disperses less. Unsteady duration is a period during which the hydraulic Based on Figure 15, the effect of the pumping rate change on the travel time, for both a single well and a group of wells, is the same. Therefore, the number of wells is not an overriding factor for determining contaminant travel time under various pumping rates.
Conversely, it indicates that the peak of the contaminant concentration depends to a great extent on the number of wells, and it is intensified when the number of wells increases. The reason that the contaminant concentration in the group of wells condition increases is that, in this case, the pollution plume takes a longer time to reach the steady-state condition in comparison to the single well; consequently, the pollutant not only further conveys, but also disperses less. Unsteady duration is a period during which the hydraulic gradient and direction of groundwater flow dramatically change, which yields a high rate of advection. In the case of the group of wells, the advection rate was much higher than the advection rate in the single well case. Therefore, contaminant concentrations for the group of wells did not decrease as much as they did for the single well case.

Conclusions
In this study, a combined numerical model was presented to make accurate and precise simulations of both groundwater flow and contaminant transport in transient conditions through saturated porous media. Within the numerical model, unsteady governing equations were numerically solved by using the finite volume time splitting scheme. The performance and accuracy of the numerical models were successfully examined by several test cases, which indicated the current numerical model has enough capability to predict radial groundwater flow toward a pumping well and contaminant migration through aquifers, including advection, dispersion, and diffusion terms.
The effect of pumping rate on solute migration was evaluated by the numerical model, and the following results were drawn. The reduction in pumping rate was accompanied by a decrease in flow velocity, which decreased both the advection and dispersion rates, as well as the growth of contaminant travel time. Discontinuing the pumping brought about a dramatic drop in flow velocity, which induced a remarkable rise in travel time, along with an increase in pollution concentrations downstream.
Comparing the results of the single well with the group of wells revealed that the single well required less time to reach the steady-state groundwater flow condition, with the same inverse correlation between the travel time and the pumping rate in contaminant transport. By contrast, the effect of the pumping rate on the contaminant concentration was different. In fact, the peak contaminant concentration, as a result of the pumping rate reduction in the group of wells, was bigger than the single well case. The results illustrate that managing the pumping rate can yield reduced pollution at the well point. Thus, in order to minimize contamination, the optimum range for the pumping rate, which is attributed to the lowest well pollution, can be determined by numerical simulations.