Groundwater and Contaminant Hydrology

1. Introduction

Groundwater and Contaminant Hydrology has a range of research relating to the transport and fate of contaminants in soils and groundwater. The scope of the center includes: 1) the development of new sampling and site characterization techniques; and 2) other improved groundwater remediation techniques.

Contaminant hydrology is the study of processes that affect both ground and surface water pollution. It draws on the principles of hydrology and chemistry. Contaminant hydrology and water quality research seeks to understand the role of soil properties and hydrologic processes on ground and surface water pollution and develop strategies to mitigate their impacts. Research is done at all scales varying from soil pore to basin scale and covers both traditional and emerging contaminants. Groundwater and contaminant hydrology studies include fate and transport of jet fuel leakages from oil depots, producing water injection in shallow wells from the oil and gas exploration field concession areas, veterinary pharmaceuticals from land-applied manure, pathogen losses from manure application, fate and transport of disposal wastes in unlined evaporation ponds from pharmaceutical industries, impacts of tile drainage on sediment and nutrient pollution on Rivers, sediment-turbidity relationships, water quality modeling, and TMDL and paired watershed studies.

Some research institutes address national and international needs for subsurface contaminant characterization and remediation across a spectrum of approaches - laboratory experiments, field tests, and theoretical and numerical groundwater flow and transport investigations. Some of the developing countries most critical subsurface contamination issues, including the chemical evolution of highly alkaline radioactive waste in storage tanks; reduction, re-oxidation, and diffusion of uranium forms in sediments; hydraulic properties of unsaturated gravels; and the natural production of transport-enhancing mobile nanoparticles in the subsurface. Inverse modeling of reactive transport and joint hydrologic and geophysical inversion are investigated to develop new tools and approaches for estimating field-scale reactive transport parameters and characterizing contamination sites.

Contaminants can migrate directly into groundwater from below-ground sources (e.g. storage tanks, pipelines) that lie within the saturated zone. Additionally contaminants can enter the groundwater system from the surface by vertical leakage through the seals around well casings, through wells abandoned without proper procedures, or as a result of contaminant disposal of improperly constructed wells [1].

1.1. Governing processes of contaminant transport

Generally three processes can be distinguished which govern the transport of contaminants in groundwater: advection, dispersion and retardation. Dispersion and density/viscosity differences may accelerate contaminant movement, while retardation processes can slow the rate of movement. Some contamination problems involve two or more fluids. Examples include air, water and organic liquids in the unsaturated zone, or organic liquids and water in an aquifer. Tracers are useful for characterizing water flow in the saturated and unsaturated zone.

• Advection

The term advection refers to the movement caused by the flow of groundwater. Groundwater flow or advection is calculated based on Darcy's law. Particle tracking can be used to calculate advective transport paths [2]. Particle tracking is a numerical method by placing a particle into the flow field and numerically integrating the flow path.

• Dispersion

Dispersive spreading within and transverse to the main flow direction causes a gradual dilution of the contaminant plume. The dispersive spreading of a contaminant plume is due to aquifer heterogeneities. Dispersion on the macroscopic scale is caused by variations in hydraulic conductivity and porosity. Solute transport can be influenced by preferential flow-paths, arising from variations of hydraulic conductivity, at a decimetre scale.

• Retardation

Two major mechanisms that retard contaminant movement are sorption and biodegradation. If the sorptive process is rapid compared with the flow velocity, the solute will reach an equilibrium condition with the sorbed phase and the process can be described by an equilibrium sorption isotherm. The linear sorption isotherm can be described by the equation:

$$C^{*}=K_{d}C$$E1

Where C*= mass of solute sorbed per dry unit weight of solid (mg/kg), C= concentration of solute in solution in equilibrium with the mass of solute sorbed onto the solid (mg/l) and K_d = distribution coefficient (L/kg)

• Non aqueous phase liquids (NAPL)

Organic liquids that have densities greater than water are referred to as DNAPL (dense nonaqueous phase liquids). Nonaqueous phase liquids that have densities less than water are called LNAPLs (light nonaqueous phase liquids). Contamination by LNAPL typically involves spills of fuels like gasoline or jet fuel.

1.2. Groundwater flow and contaminant transport modeling

The preliminary steps in modeling groundwater flow and contaminant transport include development of a conceptual model, selection of a computer code, and developing model design [3]. Defining a numerical groundwater flow model is based on parameters like (a) sources and sinks of water in the field system; (b) the available data on geohydrologic system; (c) the system geometry i.e. types and extent of model layers; (d) the spatial and temporal structure of the hydraulic properties; and (e) boundary condition. The widely used MODFLOW [4] and MT3D solute transport [5] numerical codes use finite differences schemes and are considered very reliable. MODFLOW is a three-dimensional modular finite-difference model of U.S. Geological Survey widely used for the description and prediction of the behavior of groundwater system. The program uses variable grid spacing in x and y directions. Parameter estimation can be approached to find the set of parameter values that provides the best fit of model results to field observations,. At first stage, the model computes drawdown, direction of flow and hydraulic heads on each nodal point using a finite difference grid system. Using the steady-state hydraulic heads calculated by the model as the initial condition, the MT3D model is run to simulate contaminant transport in a groundwater system. Once a model is calibrated, it can be used to make predictions for management or other purposes [6].

Two case studies including i) simulated transport of jet fuel leaking into groundwater, Sindh Pakistan; and ii) deep-seated disposal of hydrocarbon exploration produced water using three-dimensional contaminant transport model, Sindh Pakistan have been discussed to highlight the related issues, implications and concerns.

2. Case study — I

2.1. Introduction

Groundwater is a major source for domestic and industrial uses in many urban settlements of the world. Effluents from industrial areas as well as accidental spills and leaks from surface and underground storage tanks are the main sources of natural groundwater contamination. When such contamination is detected, it becomes essential to estimate the spatial extent of contamination. Conventionally, determination of the extent of contamination is undertaken by taking many samples within time and budgetary constraints from several points, which in general requires the installation of several observation wells. As a result, the cost of such operations can be very high, especially when measurements with higher resolution are required.

A three-dimensional model of the contaminant transport was developed to predict the fate of jet fuel, which leaked from above surface storage tanks in urban site of Karachi, Pakistan. Since the tanks were situated in a sandy layer, the dissolved product entered the groundwater system and started spreading beyond the site. The modelling process comprised of steady-state simulation of the groundwater system, transient simulation of the groundwater system in the period from January 1986 through December 2015, and calibration of jet fuel that was performed in context of different parameters in groundwater system. The fuel was simulated using a modular three-dimensional finite-difference groundwater model (PMWIN) ModFlow and solute transport model (MT3D) in the 1986-2001 periods under a hypothetical scenario. After a realistic distribution of piezometric heads within the aquifer system, calibration was achieved and matched to known conditions; the solute transport component was therefore coupled to the flow. Jet fuel concentration contour maps show the expanding plume over a given time, which become almost prominent in the preceding years.

Two-dimensional (2-D) solute transport models can be used to predict the effects of transverse dispersion of the contaminant plume (spreading). Additionally, 2-D models are appropriate where the contaminant source may lie within or near the radius of influence of a continuously pumping well. While three-dimensional (3-D) numerical models should only be used if extensive data are available regarding vertical and horizontal heterogeneity, and spatial variability in contaminant concentrations. A localized contaminant transport model for groundwater is developed to gain insight into the dynamics of the leakage of jet fuel from above-ground storage tanks in the metropolitan area of Karachi, Pakistan. Jet fuel consists of refined, kerosene-type hydrocarbons, which are mixtures of benzene, toluene, ethyl benzene, and isomers of xylene [7]. Hypothetical monitoring wells were established to estimate the concentrations of jet fuel over a stipulated time period as a result of continuous seepage from the storage depot. Although, no specific data on the history of seepage were available, in view of the results inferred from an electrical resistivity sounding survey (ERSS) regarding the nature of the subsurface lithologies coupled with the findings of previous investigators from Mott MacDonald Pakistan (MMP) [8], it is envisaged that seepage from the storage tanks occurred for more than a decade. ERSS is used to obtain the subsurface resistivity values that are assigned to different geological material. MMP [8] conducted the study on soil and groundwater assessment of environmental damage due to oil pollution and remedial measures were suggested for depots / installations / airfields of Shell Pakistan, scattered throughout the country. Appendix A provides the composition of jet-fuel (Table 1).

2.1.1. Site description

The project site is located between longitudes 67º 07´ 20" and 67º 10´ 30" and latitudes 24º 52´ 20" and 24º 54´ 20". The land surface elevation ranges from approximately 15 to 33 meters above mean sea level (masl). In the far south, the Malir River drains into the Arabian Sea (Figure 1). The typical lithology of the site is silty to sandy clay from 0 to 15 ft (4.6 m) bls, gravelly sand from 15 to 43 ft (4.6 to 13 m) bls, and clayey to silty sand from 43 to 80 ft (13 to 24 m) bls. The region is arid with an average annual rainfall of about 200 mm (7.9 in). Out of this, only 10% [9] is considered to recharge the aquifer system (6.34 x 10-9 m/sec).

The vadose zone is contaminated with up to 1300 ppm of total petroleum hydrocarbons (TPH) within the storage site. In the previous study by MMP [8], soil samples were collected from different locations with 3 feet (1 m) below land surface (bls) and analyzed to estimate the concentration of hydrocarbon compound. Soil samples associated with the storage area have indicated higher TPH concentrations. The contaminant plume follows the hydraulic gradient to the southwest.

2.1.2. Literature review

Kim and Corapcioglu [7] developed two-dimensional model to describe areal spreading and migration of light nonaqueous-phase liquids (LNAPLs) introduced into the subsurface by spills or leaks from underground storage tanks. The nonaqueous-phase liquids (NAPL) transport model was coupled with two-dimensional contaminant transport models to predict contamination of soil gas and groundwater resulting from a LNAPL migrating on the water table. Simulations were performed using the finite-difference method to study LNAPL migration and groundwater contamination. The model was applied to subsurface contamination by jet-fuel. Results indicated that LNAPL migration was affected mostly by volatilization. Further, the spreading and movement of the dissolved plume was affected by the geology of the area and the free-product plume. Most of the spilled mass remained as a free LNAPL phase 20 years after the spill. The migration of LNAPL for such a long period resulted in the contamination of both groundwater and a large volume of soil.

El-Kadi [10] investigated the US Navy's bulk fuel storage facility at Red Hill located in the island of Oahu. The facility consisted of 20 buried steel tanks with a capacity of about 12.5 million gallons each. The tanks contain jet-fuel and diesel fuel marine. The bottoms of the tanks are situated about 80 feet above the basal water table. The geology of the area is primarily basaltic lava flows. Investigations found evidence of releases from several tanks. Two borings were drilled to identify and monitor potential migration of contamination to the potable water source. A numerical model of the regional hydrogeology at the Red Hill Fuel Storage Facility (RHFSF) was developed to simulate the fate and transport of potential contamination from the jet-fuel tanks and the effect on the saltwater/freshwater transition zone of various pumping scenarios.

Periago et al. [11] investigated infiltration into soil of contaminants present in cattle slurry. Column experiments were performed in order to characterize the release of contaminants at the slurry-soil interface after surface application of slurry with subsequent rainfall or irrigation. The shape of the release curves suggests that the release of substances from slurry can be modeled by a single-parameter release function. They compared prediction of solute transport (a) with input defined by the release function and (b) assuming rectangular-pulse input.

Eric et al. [12] developed a parameter identification (PI) procedure and implemented with the United States Geological Survey's Method of Characteristics (USGS-MOC) model. The test results showed that the proposed algorithm could identify transmissivity and dispersivity accurately under ideal situations. Because of the improved efficiency in model calibration, extended application to field conditions was effective.

Figure 1.

Location of the study area in Southern Pakistan

Jin et al. [13] investigated hydrocarbon plumes in groundwater through the installation of extensive monitoring wells. Electromagnetic induction survey was carried out as an alternative technique for mapping petroleum contaminants in the subsurface. The surveys were conducted at a coal mining site near Gillette, Wyoming, using the EM34-XL ground conductivity meter. Data from this survey used to validate with the known concentrations of diesel compounds detected in groundwater. Groundwater data correlated perfectly with the electromagnetic survey data, which was used to generate a site model to identify subsurface diesel plumes. Results from this study indicated that this geophysical technique was an effective tool for assessing subsurface petroleum hydrocarbon sources and plumes at contaminated sites.

2.2. Model conceptualization and simulation

The groundwater flow system is treated as a two layer. The upper layer (predominantly silty sand) is bounded above by the water table and is 15 ft (4.6 m) thick, while the lower layer (predominantly clayey sand) is 65 ft (20 m) thick. These are unconfined and recharged from the surface by infiltrating rain, but only over permeable surfaces. A small stream runs along east of model area acts as a drain to groundwater, which flows from the northeast to the southwest. With the exception of the stream, all other boundaries are artificial, that is neither constant head, nor constant flow boundaries. The processes that control the groundwater flow are: (i) recharges from infiltrated rainfall; (ii) flow entering the model across the eastern boundary (also across the northern boundary); (iii) flow reaching the stream; (iv) flow leaving the model across the southern boundary; and (v) pumping from one well near tank no. 9. Values of hydraulic conductivities (K) for the layers are taken from the literature [14]. The K value for depth range 16-31 ft (4.9-9.4m) is taken as 1.7x10-6 ft/sec (5.2 x10^-7 m/sec) and for depth range 31-97 ft (9.4-30 m) as 1.5x10^-5 ft/sec (4.6 x10^-6 m/sec).

The transport and fate of hydrocarbons depend on multi physical and chemical processes, including advection, dispersion, volatilization, dissolution, biodegradation, and sorption. When a solute undergoes chemical reactions, its rate of movement may be substantially less than the average rate of groundwater flow. In this study, retardation of the movement of dissolved hydrocarbons is simulated as a sorption process, which includes both adsorption and partitioning into soil organic matter or organic solvents. The MT3D software was used for simulation [5]. It uses a linear isotherm to simulate partitioning of a contaminant species between the porous media and the fluid phase due to sorption. This sorption process is approximated by the following equilibrium relationship between the dissolved and adsorbed phases:

$$S=K_{d}C$$E2

Where S is the concentration of the adsorbed phase (M/M), C is the concentration of the dissolved phase (M/L3) and K_d is the sorption or distribution coefficient (L³/M). K_d values for organic materials are commonly calculated as the product of the fraction of organic carbon in the soil, foc, and the organic carbon partitioning coefficient, K_oc, or K_d = foc Koc. Koc values are contaminant specific and reported in various sources [15-17]. The foc in the uncontaminated soil was estimated to range from 0.001 to 0.02 based on guidelines by [18]. Assuming the linear isotherm, the retardation factor (R) is expressed as follows:

$$R=1+\left(r_b/n_{eff}\right)K_{d }$$E3

Where r_b is the bulk density of the porous material (M/L3) and n_eff is the effective porosity.

For jet fuel, the distribution coefficient K_d is taken as 0.004415 ft³/kg. With these values [7],

$$R=1+\left[\left({\text{48 kg/ft}}^{\text{3}}\right)\text{/0}\text{.25}\right]x0.004415=1.848$$

2.2.1. Numerical ground water flow modeling

Processing ModFlow for Windows (PMWIN5), a modular 3-D finite-difference groundwater model, is used to configure the flow field [4]. The model consists of 41 columns and 39 rows in each layer (Figure 2). The size of cells is 410 ft x410 ft (125 m x 125 m) outside the fuel storage domain and 205 ft x 205 ft (62.5 m x 62.5 m) within the storage domain. Automatic calibration of the water table was made with algorithm - UCODE and a perfect match obtained with the known condition prior to developing the transport model [6]. Using the steady-state hydraulic heads calculated by PMWIN5 as the initial condition, the solute transport model MT3D was run to simulate the dispersion of the dissolved jet-fuel plume [5]. The parameters adjusted were the retardation factor R for each cell within the finite-difference grid, and the dispersion coefficient. Concentration-time curves have been calculated for ten monitoring wells. PMPATH [19] is used to retrieve the groundwater flow model and simulation result from PMWIN5. A semi-analytical particle-tracking scheme is used to calculate the groundwater flow paths, travel times, and time-related capture zones resulting from pumping a neighboring well at the storage facility [20]. As a preprocessor to modeling and creating input data files, the PMWIN5 utility package was used. Prior to initiating the modeling work, a groundwater information system was established with all data in binary and / or ASCII files that could be exported to other softwares.

2.2.2. Locations of hypothetical wells

The dissolved phase jet-fuel plume was traced using a combination of ten hypothetical monitoring wells (Figure 2) known as MW-1 through MW-10. The wells served to identify lithology, observe water levels, and monitor concentrations of organic compounds. The wells extend to a depth of 80 ft (24 m). In addition, actual well was completed to a depth of 100 ft (30.5 m) near storage tank no.9 in case of emergency need. In the modeling study, this well was used to track the time-related capture zone. The general layout of the storage tanks over the finite-difference grid is shown in Figure 3. The location of the pumping well is marked as a small red square in Figure 2 and Figure 3.

Figure 2.

Model design indicating finite difference grid and locations of hypothetical monitoring wells

Figure 3.

Layout of the storage tanks marked on a finite-difference grid

2.2.3. Model calibration

The model was calibrated for steady-state conditions. Since it was speculated that the seepage of the fuel might have started as early as 1986, simulation of the groundwater flow was begun in January 1971, the opening of storage facility. In the steady-state phase the only input comes from constant head boundaries (along the east and west of the model) and from infiltrated rainfall. All output goes into constant head boundaries (the stream and the southwest boundary of the model). The differences in the known and simulated heads were calibrated to less than 0.30 ft (0.091 m) by making slight adjustments in the K values of both the layers. Transient calibration of groundwater flow was accomplished using the time-variant hydraulic head values. Parameters such as recharge rates during each stress period, hydraulic heads in the stream and along the model boundaries, aquifer storage properties, pumping rates, and time-dependent capture zone were adjusted during the calibration. To be objective and consistent, the recharge from infiltration was made equal to 10% of rainfall in each month. Effective porosity of the aquifer was varied between 15% and 25% until the value of 25% was determined to be the best predictor for the model. Constant pumping rates of 1500 US gallons/hr (1.58 x 10^-3 m³/sec) and 1000 US gallons/hr (1.05 x 10^-3 m³/sec) were used in layer 1 and layer 2 respectively. The period from 1986 through 2000 was divided into seven stress periods, each of 2 years in duration. From 2000 to 2001, one stress period was assigned. The length of a stress period was made equal to the number of days in that month.

2.3. Results and discussion

The effect of the pumping well is clearly visible as a cone of depression (Figure 4). The drawdown was determined to be 14.0 ft (4.3 m) near the storage facility.

Figure 4.

Cone of Depression visible around pumping well developed in layer 1

The model assumes uniform recharge from infiltrated rainfall to every “recharging” cell. Although effective porosity, hydraulic conductivity, and recharge may vary in space and time, the model is expected to have produced a reasonable configuration of the groundwater flow pattern throughout the whole period of simulation. The time-related capture zones produced due to constant pumping are shown in Figure 5 and Figure 6. Water balances, which was calculated for each year of the simulated period, showed a perfect match between the input and output components.

Figure 5.

Capture zone of the pumping well with arrows indicating flow directions

Figure 6.

Days-capture zone calculated by PMPATH

2.3.1. Calibration of plume dispersion

For simulation of the movement of dissolved jet fuel, lateral hydraulic conductivity values equal to 0.864 ft/day (0.263 m/day) for layer 1 and 1.29 ft/day (0.393 m/day) for layer 2 were accepted, while the vertical hydraulic conductivity were taken as 0.0864 ft/day (0.0263 m/day) and 0.129 ft/day (0.0393 m/day) for each layer, respectively [14]. The hydraulic gradient and flow-net were obtained by running the flow component of the model derived from water level information in the previous study.

The United States Environmental Protection Agency (USEPA) and the Georgia Environmental Protection Division (GAEPD) recommend that the value for longitudinal dispersion should be one-tenth of the distance from the place where a contaminant enters the groundwater system to the down-gradient receptor (a well, stream, or other point of compliance). The distance from the storage facility (tank no. 7) to the pumping well is approximately 100 ft (30.48 m). In all calibration runs, as recommended, the value for longitudinal dispersion was set at 10 ft (3 m). USEPA and GAEPD also recommend for a solute transport model that the value for transverse dispersion equal one-third for the longitudinal dispersion. For this model, transverse dispersion would equal 3.3 ft (1.0 m). In the simulation of the fate of jet fuel, the transverse dispersion coefficient was varied within a range of 2.0 ft to 3.3 ft (0.61 m to 1.0 m). In the model a value of 0.001 ft²/day (1 x 10^-5 cm²/sec) was used for molecular diffusion. With a retardation factor of 1.80, dissolved jet fuel takes 1.33 years to travel a lateral distance of about 70 to 80 ft (21 to 24 m) in groundwater beneath tank no. 8. The best value of the microbial decay coefficient for jet fuel is estimated to be 1-10 / day with a microbial yield coefficient for oxygen of 0.52 [14].

2.3.2. Strategy development for release of Jet fuel

The previous integrity test run on the storage tanks containing 10 million liters of jet fuel indicated no loss. The date when the leak initially began is unknown, although inventory records indicated that the leak was not present before tank integrity testing. The product has been detected in several hypothetical-monitoring wells (notably in MW-2 and MW-4) and in many soil samples taken within several tens of feet of the tank. The initial concentration of jet fuel entering the system is not of prime concern for the modeling. The product of the influx (in L³/T) and the concentration (in M/L³) gives the total mass of jet fuel entering the system in a certain time interval. For the purpose of calibrating the jet fuel input, the initial concentration used, based upon field data [8], varied from 0.095 to 0.19 g/ft³ (0.0027 to 0.0054 g/m³). The initial mass of jet fuel, as simulated by the model, was equal to each of four cells “injecting” at a mass rate of 95 to 190 g/ft³ (2.7 to 5.4 g/m³) following the initial period of 15 years during which no groundwater contamination was assumed (Table 1).

Phase	Stress period	Condition
Safe Period	15 years (1971 to 1986)	No leakage
Hazardous Period	10 years (1986 to 1996)	Low to moderate leakage
Risk Assessment	5 years (1996-2001)	Moderate leakage
Future Prediction	14 years (2001 – 2015)	Accretion in leakage

Table 1.

Strategy developed for the plume modeling scenarios

Using steady-state hydraulic heads as initial conditions, the evolution of the plume was modeled over nine stress periods as a result of continuous seepage from cells (18,16,1; tank7), (19,16,1; tank 10), (19,17,1; tank 8), and (20,17,1; tank9) as shown in Table 2.

Stress Period	Time interval(years)	Elapsed Time(sec)	Period
1	2	6.30 x 10^7	1986 – 88
2	2	12.60 x 10^7	1988 – 90
3	2	18.92 x 10^7	1990 – 92
4	2	25.23 x 10^7	1990 – 94
5	2	3.15 x 10^8	1994 – 96
6	2	3.78 x 10^8	1996 - 98
7	2	4.41 x 10^8	1998 –00
8	1	4.73 x 10^8	2000 – 01
9	14	9.14 x 10^8	2001 – 15

Table 2.

Stress period used in time-dependent solute transport modeling of jet fuel

2.3.3. Calibration scenario

Parameters describing various processes are used after calibration with different combination of parameters (Table 3).

Parameters	Value
Longitudinal Dispersion	10 ft
Transverse Dispersion	3.3 ft
Molecular Diffusion	0.001 ft^2/day
Distribution Coefficient	0.004415 ft^3/kg
Retardation Factor (R)	1.80 to 1.848
Decay Coefficient	1 x10^-9 day^--1
Hydraulic Conductivity K (Layer-1)	1x10^-5 ft/sec (0.864 ft/day)
Hydraulic Conductivity K (Layer-2)	1.49x10^-5 ft/sec (1.29 ft/day)
Effective Porosity (Layer-1)	0.25
Effective Porosity (Layer-2)	0.30

Table 3.

Preliminary and final values of parameters used in modeling

The release of the jet fuel is simulated in four cells, all along columns 18 to 20 from row 16 to row 17. The area of injection is equal to 42025 ft² (3,904 m²). The concentration of jet fuel at the source (95 to 190 mg/ft³ [2.7 to 5.4 g/m³]) maintained constant throughout the designated “hazardous period” simulation period (1986-1996). The concentration was increased slightly (about 0.01 %) from 1996 through 2011 and further up to longer time duration of 4 years i.e., up to 2015. The plume simulations are shown in Figure 7.

Figure 7.

Simulated Jet-fuel plumes(1986 to1988 and 1988 to 1990)

The jet-fuel break-through curves for the hypothetical monitoring wells are shown in Figure 8. Conventionally, determination of the extent and level of contamination is undertaken by taking multiple measurements in wells [21-23]. However, higher spatial resolution generally requires installation of monitoring wells, which is costly [24].

Figure 8.

Concentration versus time based on data from 10 monitoring wells

The modeled concentration at MW-1, MW-2, and MW-6 was much higher than the concentrations in the remaining wells. The maximum level of concentrations recorded in MW-6 was 2990 μg/ft³. Figure 9 reflects the plume spreading of year 2001.

The shape of the plume is elliptical, with the major axis in the direction of groundwater flow. This shape results from advection and longitudinal dispersion. The lateral spread of the plume results from transverse dispersion and molecular diffusion. Upgradient spread of the plume results from molecular diffusion [25-26]. The plume travels toward the stream, which is still far away in the west. By the end of 2015 the effect of the plume becomes evident and monitoring wells MW-4, MW-5, and MW-6 indicated increased concentration of jet fuel (Figure 10). The resultant plume appears to be spreading more in the elliptical path but in the direction of groundwater flow.

Figure 9.

Extent of simulated plume in 2001

Figure 10.

Extent of plume spreading in year 2015

2.4. Conclusions and recommendation

Based on the modeling, it is concluded that the jet-fuel plume has neither expanded nor moved considerably. It is less than 250 ft (76.2 m) beyond the storage tanks and is oriented northeast to southwest. The level of concentration found in the simulated monitoring wells is significant, but because groundwater is brackish and thus unlikely to be used, no harmful effects are expected. However, with the continuous process of leaking jet-fuel from the storage depot, the level of concentration is expected to increase over the period of time. Regionally, the jet-fuel expansion will have on prominent effect over a longer areal coverage and will be confined to a localized area.

An interdisciplinary investigation of the processes controlling the fate and transport hydrocarbons in the subsurface is needed. Concentrations should be performed in wells within and down-gradient of the plume, as field data would help develop a stronger argument for the fate of jet fuel in groundwater. Periodic observations need to be carried out in wells to have good control on the changes of groundwater chemistry. Defective storage depots need to be mended to stop the release of jet-fuel in future.

2.5. Appendix A

Composition of contaminant (Jet fuel)

Knowledge of the geochemistry of a contaminated aquifer is important to understand the chemical and biological processes controlling the migration of hydrocarbon contaminants in the subsurface. Originally, the jet fuel (kerosene oil) is the name assigned to a material with a biological origin, but now it is used to describe materials most of which contain carbon and hydrogen and which may contain oxygen, nitrogen, the halogens, and lesser amounts of other elements. The simplest of these are the hydrocarbons, molecules of hydrogen and carbon, many of which are the components of natural gas, petroleum, and coal. Petroleum, however, has a very large number of components ranging from methane to the high molecular weight materials asphalt and paraffin. Typical fractions into which crude oil is separated in an oil refinery and some principal molecular species are shown in Table 4.

Fraction from distillation	Boiling range	Product of secondary treatment	Typical molecular components
Gas	Below 20^o C	Gas; Liquefied Pet. Gas (LPG)	CH4 methane, C2H6 ethane; C 3H8 propane,; C 4H10 butane
Naphtha	20^o – 175^o C	Naphtha gasoline	C 11H24 C 18H38
Kerosene	175^o – 400^o C	Kerosene diesel fuel	C 11H24 C 18H38
		Lubricating oil	C 15H32 C 40H82
Residue	above 400^o C	Asphalt	Heavier hydrocarbons

Table 4.

The Fractions and Representative Components obtained from Crude Oil

3. Case study — II

3.1. Introduction

It is important to maintain the existing quality of groundwater because once contamination occurs; it is sometime difficult or rather impossible to clean the aquifer. There is high probability associated with certain landuses like agriculture, industrial/urban land and drainage wells for contaminating the groundwater. The early detection of such contamination can be executed through proper monitoring of the groundwater quality. Many oil & gas companies are disposing off their producing water into the deep Ranikot formation in Bhit oil-field area of southern Pakistan. The producing water contains Total dissolved solid (TDS) within range of 18,000 - 22,000 mg/l besides oil condensate. There were concerns that the producing water is affecting the fresh water aquifers belonging to the overlying Nari and Kirthar formations. This phenomenon has been studied by the utilization of groundwater contaminant transport model. Injection has been monitored at 2100 meters depth in the Pab sandstone formation. A three-dimensional contaminant transport model was developed to simulate and monitor the migration of disposal of hydrocarbon exploration produced water in Injection well at 2000 meters depth in the Upper Cretaceous Pab sandstone in the study area. Framework of regional stratigraphic and structural geology, landform characteristics, climate and hydrogeological setup were used to model the subsurface aquifer. The shallow and deep-seated characteristics of geological formations were obtained from electrical resistivity sounding surveys, geophysical well-logging information and available drilling data. The modeling process comprised of steady-state and transient simulations of the prolific groundwater system and, predictive simulation of contaminants transport after 1-, 10- and 30-year of injection. The contaminant transport was evaluated from the bottom of the injection well and its short and, long-term effects were determined on aquifer system lying in varying hydrogeological and geological conditions.

3.1.1. Description of study area

The study area of Bhit oil field is located about 43 km south of the Manchhar Lake within longitudes 67 ^o 25' - 67 ^o 48' E and latitudes 26 ^o 01' - 26 ^o 30' N in Dadu district in Sindh province of Pakistan (Figure 11). Manchhar Lake, one of the largest lake in Pakistan and in Asia, is formed in a depression in the western side of the Indus River in Sindh province. The total catchment area of the lake is about 97,125 km² [27]. The surface area of the lake fluctuates with the seasons from as little as 350 km² to as much as 520 km² [28]. The lake is fed by two canals, the Aral and the Danister emerging from the river Indus in their eastern side. Due to less rainfall and contamination of surface water, the Manchhar Lake contains brackish water. The elevation ranges between 45 m at the Manchhar Lake and 163 m towards Bhit study area. On regional scale the area is a part of Gigantic Indus river basin composed of alluvium transported by the river and its tributaries. The main surface water sources are Naig Nai stream originating from Bhit Mountain range in the west, Dhanar Dhoro stream passing close to the oil-field plan, besides other small intermittent streams which remain dry in most parts of the year. The discharge data of these streams were not available. The water supply for local communities is maintained by the springs originating from the nearby Limestone Mountains of the Kirthar Range. According to reconnaissance studies conducted in the area, fresh groundwater source is available at different locations in the alluvium. Presently, there is no significant groundwater development in the area. Only few water wells were constructed to fulfill the requirements of local communities. The major sources of potable water to humans and livestock and for irrigation are unconsolidated aquifers. A water supply scheme consisting of four tube wells at Jhangara village is providing water to the villages between Manchhar Lake and Jhangara. The rainfall is scanty. Average annual rainfall is about 200 mm. It is higher in summer months like July and August due to prevalence of monsoon conditions. The aquifer area is located in the alluvial deposits along the Naig Nai stream.

Figure 11.

Location of study area in southern Pakistan

The water table is generally in phreatic to semi-unconfined conditions. The observation wells drilled in the area indicated water table depth of about 12 m and hydraulic head value of 148.8 meter above sea level (masl). The groundwater flow is generally from southwest towards northeast direction. The flow direction of groundwater is true replica of the flow direction of Naig Nai stream draining the area. The groundwater level is mainly influenced by seasonal floods, stream flows and tubewell (water well) discharge. The fluctuations are small in the deep water table in the piedmont plain.

3.2. Material and methods

The data of subsoil properties, aquifer characteristics and existing groundwater conditions were collected through reconnaissance level field investigations including geophysical survey. The electrical resistivity soundings survey (ERSS) was conducted primarily to collect required input data for the modeling. The surface drainage and topographic information were extracted from the topo-sheets prepared by Survey of Pakistan on 1:50,000 scale. The climate data i.e. precipitation and temperature, of 1961-1976 period were acquired of nearby meteorological stations like Karachi, Nawabshah, Moenjodaro and Khuzdar from Meteorological Department and Water and Power Development Authority (WAPDA) Pakistan. The published literatures of the region i.e. see [29-32] had been used to firm up the study results.

3.2.1. Geophysical data analysis

The geophysical/drill logs of the injection well field suggest that the subsurface material is composed of layers of sandstone, limestone, dolomite and clay-stone of different formations. One of the interpretative seismic sections of the Bhit concession area (shown in Figure 12) indicates the deeper Chiltan limestone formation beyond the depth of injection well. The hard rock aquifers are mainly composed of partially fractured limestone and sandstone belonging to Nari and Kirthar formations. Limestone, which is the dominant formation, has solution channels due to water action having secondary permeability characteristics. Further, chances of transport contamination could take place much significantly through fracture zone of limestone and dolomitic formations. The hydraulic properties of the underlying overburden and rocks evaluated through geophysical well log data are shown in Table 5.

S.No.	Formation	Depth(m)	K(m/d)	Transmissivity T (m^2/d)
1	Nari (sandstone)	113	5	130; (S=.01)
2	Kirthar (limestone)	616	24	2000; (S=.04)
3	Ghazij (claystone)	762	zero	Regional Seal
4	Laki (limestone & dolomite)	1,259	0.165	4; (S=0.00005)
5	Dunghan (dolomite and claystone)	-	-	-
6	Ranikot (Lakhra+Bara+Khadro)	1,799	0.68	14.97; (S=0.000007)
7	Pab (sandstone)	2,000 to onward	0.138	17.5; (S=0.000005)

Table 5.

Summary of the Aquifer characteristics of hard rock formations

The hydraulic conductivity of the unconsolidated deposits is about 19 m/day and the effective porosity is 0.25. Several hydrocarbon wells are producing gas in the Bhit concession area. The Ghazij (claystone) was found to be a cap rock (regional seal) over Pab sandstone - an enriched reservoir of hydrocarbon. The Ranikot formation is a prolific rock unit having good transmissivity and storativity to accumulate disposal waste of the producing gas Bhit concession.

Figure 12.

Interpretative seismic depth section indicating several geological formations

3.2.2. Model conceptualization and simulation

The groundwater flow system was treated as a multi-layered system. The upper layer aquifer is mainly unconfined and at depths where silt and clay horizons are present, the sand could probably cause partial confinement in some areas. The chances of transport contamination could take place much significantly through fracture zone of limestone and dolomitic formations. Seven aquifer layers were defined on the basis of physical characteristics of lithological formations (Table 1). The disposal of produced water in the injection well is set at 2,100 m depth in the Pab sandstone formation. The main source of recharge to groundwater is rainfall which is highly variable. During rains, different nullahs carry flows and infiltrations through the piedmonts and alluvial fans and cause subsequent lateral movement at depth. The recharge is higher in summer especially in the months of July and August. It may occur during rainy month of March in winter, a mean value for the limestone recharge may be taken as 200 mm/year based on the recharge data of rainfall. The main discharge components are groundwater extraction from water wells/dug wells, evapotranspiration, and spring discharge. The abstraction of groundwater becomes higher during months of little recharge to groundwater i.e. November to January, which may affect the storage of groundwater for a limited period.

The MODFLOW code and the MT3DMS code were used to solve the flow and transport equations. The model domain comprised of 40 x 30 grid network with total area of 1200 sq km (Figure 13). First, the model was run for steady state condition. The model took up groundwater extraction from the oil-field water wells (TW1 & TW2) and the community tube wells (TW3 & TW4). The discharge rate of the water wells was 0.0083 m³/sec and discharge rate of deep seated injection well (DW) as supplied was 0.00152 m³/sec. The injection of the produced wastewater from the injection well was considered during simulation. Once the flow model was completed and run was carried out, the contaminant transport model was set and simulated to evaluate the groundwater contamination and movement of plume.

Figure 13.

Manchhar Lake treated as constant head boundary; Bhit Range as impervious boundary. DW deep injection well, TW tube well

3.2.3. Simulation of contaminant transport of injection well

The steady-state hydraulic heads were used as initial condition in MT3DMS option available in PM5 to simulate the dispersion of plume. The MT3DMS model simulates the processes i.e. advection, mechanical dispersion, retardation, decay and molecular diffusion related to the fate of contaminant. Initial concentration was set to zero in all the layers. In Advection package, 3^rd order TVD scheme [12, 33] was selected. This method is considered as a good compromise between the standard finite difference and particle tracking approaches. In dispersion package, TRPT (Horizontal transverse dispersivity/Longitudinal dispersivity) was set to 0.3 for all the layers except in layer 3, where it was set to 0.1. In dispersion package, TRPV (Vertical transverse dispersivity/Longitudinal dispersivity) was set to 0.3 for all the layers except in layer 3. For this layer, it was set to 0.1. Longitudinal dispersivity was set to 10 m. There was no sorption selected in chemical reaction package. The injection well was set at layer 7 in the sink/source menu. The concentration of the injection well liquid was considered to be 100 ppm. The model was then simulated for 1-, 10-, and 30-year period for studying the behaviors produced wastewater in and around the injection well.

3.3. Results and discussion

The computations were carried out for three cases i.e. in case-I, Injection well was continuously discharged for one-year period, in case-II it well was simulated for 10-year period and in case-III for 30-year period. The hydraulic heads and drawdown were computed in all three cases. The velocity vectors prominent in layer-1 tend to move in the northeast direction towards Manchhar Lake. The groundwater flow had shown decline in confining layers like 4 and 7. The two tube wells of the oil-fields each discharging at the rate of 8.3 x 10^-3 m³/sec along with community wells were used for the study. The results obtained from the 3-D transport model are shown in Table 6 and transport of contaminant plumes in three simulation periods in Figures 14-16.

Layer	1 Year	10 Years	30 Years
1	-	2.21 x 10^-23	1.96 x 10^-20
2	3.16 x 10^-27	3.88 x 10^-21	1.32 x 10^-18
3	3.14 x 10^-21	3.54 x 10^-16	4.01 x 10^-14
4	3.16 x 10^-15	3.18 x 10^-11	1.06 x 10^-09
5	9.58 x 10^-10	9.53 x 10^-07	2.45 x 10^-05
6	8.04 x 10^-05	7.95 x 10^-03	7.21 x 10^-02
7; Pab (Sandstone)	6.28 x 10^-02	0.626	1.861

Table 6.

Maximum concentration observed in different simulation periods (ppm)

After 30 years of simulation period, only traces of contamination were found in Ghazij Formation. Moreover, it is found that after 1-year period of simulation the produced wastewater will reach upward in layer-5 (Ranikot Formation) emerging from layer 7 (Pab sandstone) as shown in Figure 14. In this period, no contamination was found in layer 1 and 2. In 10-year simulation a plume of produced water moved from layer 7 to layer 5 (Figure 15). Only traces of contamination were found in layer 3 (Ghazi Formation). In Figure 6, plume of produced water contamination indicates movement from layer 7 to layer 4 after 30 years simulation. The layer 3 was found to be acting as a regional confining seal. In this layer, only traces of contamination were present. The movement of produced wastewater was found within a radius of 3 km at the bottom of injection well in the Pab sandstone. The upper aquifers in the alluvial deposit, Nari sandstone, and Kirthar limestone was remain safe from the effects of produced wastewater disposal from the deep seated injection well. The community wells tapping in the upper few tens of meters, naturally oozing springs and the Manchhar Lake located about 43 km from the injection well were also found to be safe from the effects of produced water injection even after contaminant transport simulation of 30-year period. The development of plume was significant in layer 7 and upward in the three cases (shown in Figures 14-16).

Figure 14.

Upward movement of the plume from layer 7 to layer 6 (a); from layer 6 to layer 5 (b) and from layer 5 to layer 4 (c) after 1 year simulation period

Figure 15.

Upward movement of the plume from layer 7 to layer 6 (a); from layer 6 to layer 5 (b) and from layer 5 to layer 4 (c) after 10 year simulation period

Figure 16.

Upward movement of the plume from layer 7 to layer 6 (a); from layer 6 to layer 5 (b) and from layer 5 to layer 4 (c) after 30 year simulation period

4. Conclusions

The results of the three-dimensional groundwater modeling study highlighted the hydrogeological characteristics and features of the contaminant transport of the deep injected tube well in the Bhit oil-field area. The groundwater contaminant transport modeling technique has proved to be effective in simulating produced wastewater plume from the deep seated injection well. The study would provide base for evaluating risks of contaminants on long term basis in similar conditions in future. Risk of expansion of plume regionally does not exist as the disposal of wastewater is made in the deeper horizon well below the aquifers and also the quantity is quite limited.

Thorough understanding of surface hydrology, hydrogeological conditions and contaminant behavior in the aquifer system coupled with application of reliable modeling techniques could be helpful in dealing with water management issues related to contaminant hydrology.