Groundwater contaminant transport modeling using MODFLOW and MT3DMS: a case study in Rajshahi City

Water is a basic necessity for the functioning of all life forms that exist on the Earth. Like all sources of water, groundwater is not safe from pollution. Human activities such as using pesticides on agricultural fields (Yu et al. 2015), unplanned dumping of industrial waste, underground pipeline leakage (Liang He 2009), coal mines (Adhikari & Mal 2021), landfill sites (Nyika 2021), etc. are contributing to groundwater pollution. Among these sources, landfill leachate is one of the crucial sources of groundwater pollution.

A municipal solid waste (MSW) landfill site is considered to be an important source of groundwater contamination due to the leakage of leachate. Leachate is formed mainly due to the percolation of rainfall through the wastes, infiltration through soil pores, and finally joining the groundwater table. Leachate is a complex mixture of pollutants having high chemical oxygen demand, high ammonium nitrogen content, and lasting toxicological characteristics (Han et al. 2016). Many factors affect leachate composition such as the type of waste, the method of exploiting landfill, the availability of oxygen, the hydrogeological condition as well as the age of the landfill (Chofqi et al. 2004). Assessing the severity of groundwater pollution, authorities banned the use of several chemicals. But banning the use of pollutants does not immediately eradicate their effect on groundwater (Beegum et al. 2020). The impact of the contaminant on groundwater stays many more years after initial exposure (Adhikari & Mal 2021).

To minimize the effect of leachate, proper natural and synthetic liners in the landfill site should be provided. But for most landfill sites in Bangladesh, the liner is not provided, which makes the site more prone to groundwater contamination. The unlined landfills are more susceptible to contaminating the groundwater with potentially hazardous chemicals, which makes the groundwater unsafe for drinking purposes (Reyes-López et al. 2008). However, lined landfills also pose a threat to groundwater quality because the liners fail eventually (Sizirici & Tansel 2015). After contaminants travel through the vadose zone (Zhang et al. 2019) and reach the groundwater, the pollutants in leachate can mix with the aquifer system.

The amount of waste produced in Rajshahi City each day is 350 metric tons which increase to 400 metric tons during the summer season. Among them, 12 metric tons or 3.43% of total waste is recognized as hazardous waste. Among the total collected waste, only 210 metric tons is dumped on a landfill site (Ahsan et al. 2014). The remaining 140 metric tons of waste are dumped straight into drains, waterbodies, and open spaces. The total assessment revealed that about 80% of total waste is organic, 15.72% of plastic or nylon, and about 15% of the total waste is metal (Olanrewaju et al. 2009).

The increasing waste generation and disposal resulted in increased groundwater pollution and the unsuitability of the use of soil within the area for agricultural purposes. Areas such as food shops, houses, tube wells, and ponds that are close to the landfill site and along the direction of groundwater flow possess a greater risk of groundwater pollution (Singh & Garg 2016) So it is obvious that continuous monitoring of groundwater pollution with quantitative measurements of different contaminants in groundwater is the prime concern to manage the groundwater pollution.

Researchers have explored different tools to quantify the contaminants in groundwater. The use of deep neural network (DNN) emulators for a generalized groundwater flow and contamination transport model has proven to be a good foundation for risk assessment of groundwater contamination (Yu et al. 2020). Simplified models and Cellular Automata (CA) approach in contaminant transport modeling shows the suitability of this approach in unsteady conditions (Milašinović et al. 2019). HYDRUS 1D, along with MODFLOW and MT3DMS, shows prominent results in evaluating the distribution of contaminants in the unsaturated zone (Zhang et al. 2019). Among all the tools for quantifying groundwater pollution from a source, MODFLOW, coupled with MT3DMS is widely used and also provides satisfactory results (Singh & Garg 2016). It helps to quantify the concentration and migration distance of pollutants and the consecutive affected area. In Bangladesh, MODFLOW coupled with MT3DMS was used for analyzing the effect of salinity intrusion in coastal areas (Kumar Adhikary et al. 2011; Rahman et al. 2019).

This study assesses the probable concentration of heavy metals (Cr, Cu, Pb, Mn, and Zn) in the affected area using MODFLOW (Harbaugh et al. 2017), coupled with MT3DMS-USGS (Bedekar et al. 2016). A groundwater flow model using MODFLOW is developed for the study area and the effect of five major heavy metal exposure on groundwater is evaluated using MT3DMS-USGS. To reduce the complexity of modeling processes, certain geochemical assumptions have been introduced. Like, it was anticipated that each layer and each cell of the model would have similar hydraulic conductivity (K_d). The heavy metals do not interact chemically with one another or with any other substances. Nevertheless, depending on migration distance and duration, water tends to dilute the concentration of heavy metals. The adsorption and buffering capacity of the aquifer has no impact on the migration of contaminants.

The results show the affected area, migration distance, and affected wards of Rajshahi City for respective contaminants and generate an empirical equation of migration distance with respect to time for quantitative measurement that is not considered in previous studies.

Rajshahi is one of the major metropolitan cities in Bangladesh. The Rajshahi City Corporation (RCC) is one of the major self-governing cities in Bangladesh. Geographically, the latitude and longitude of Rajshahi City lie between 24° 20′ N to 24° 24′ N and 88° 32′ E to 88° 40′ E. It is situated within the Barind tract, 23 m above sea level. The approximate area under the RCC is 48 km². According to Bangladesh Population Census 2001 (BBS 2003), the population of Rajshahi City was 388,811. But the population had increased to approximately 449,756 according to the Bangladesh Bureau of Statistics (BBS 2011). The city is located on the bank of the Padma River.

The landfill site located in City Hut, Nowdapara is one of the significant sources of groundwater pollution in Rajshahi City. The daily waste dumping rate in the landfill site is 231.5 tons (210 metric tons) (Islam 2015). The area covered by landfill is approximately 5 acres. On the other hand, the Khulna City Corporation (KCC) produces about 450 tons of MSW per day. KCC and community-based NGOs are taking care of only 42% of the total waste generated, which is about 189 tons (Islam et al. 2019). The area covered by the KCC landfill site is 20 acres (ADB 2014), which is 4 times the area of the landfill site of Rajshahi City. The heavy metal concentrations in leachate extracted from the KCC landfill site have been investigated by Karim et al. (2017). Due to the unavailability of test data for the Rajshahi City landfill site, the heavy metal concentrations were estimated from the KCC landfill site based on area and dumping rate (Table 1).

In this study, a model was developed by integrating the groundwater flow simulation model as MODFLOW and the contaminant transport processes simulation model as MT3DMS-USGS. The outputs from the groundwater flow model (MODFLOW) were used as an input in MT3DMS-USGS for simulating the contaminants’ transport processes. MODFLOW is commonly known as a modular three-dimensional finite-difference groundwater flow model (Anderson et al. 2015). Groundwater contamination transport models like MT3DMS-USGS solve the advection-dispersion-reaction equation in a groundwater flow system under generalized hydrogeologic conditions (Bedekar et al. 2016). ModelMuse is used as the graphical user interface (GUI) for running the MODFLOW and MT3DMS-USGS code. The study area and landfill area were imported into the ModelMuse interface. The grid cells are the basic units to provide all the data as input and obtain corresponding data as output. The study area has been discretized into 100 m × 100 m grids. The water table contour is directly imported into the software for defining the gradient of the groundwater table.

The groundwater table contour lines were assigned as a time-variant specified head boundary with the value of the heads specified in the contour map. At the top active layer, the recharge boundary and evapotranspiration boundary were assigned. The preconditioned conjugate-gradient solver was used to solve the groundwater flow model. Basic transport, advection, dispersion, sink, and source boundaries are used for contaminant transport modeling using MT3DMS-USGS.

The required input parameters for model development, calibration, and validation process are listed in Table 2. Various studies suggested that the hydraulic conductivity of coarse sand and fine sand varies from 10⁻³ to 10⁻⁵m/s (Kumar et al. 2016). Ku (2013) shows that the hydraulic conductivity of saturated fine sand ranges from 10⁻⁴ to 10⁻⁵ m/s. For this study, hydraulic conductivity is considered to be 10⁻⁴ m/s which is also used for calibrating the model.

Though this study focused on the determination of the horizontal migration of the contaminants over time, only two vertical layers are considered (al Mamunul Haque et al. 2012a). The value of specific storage and yield was taken as 10⁻⁵ and 0.2, which is the default value of MODFLOW. The value of longitudinal dispersivity (LD) and the other parameters are specified by validating the model, which is described in Section 2.5. The contaminant transport modeling parameters such as diffusion coefficient, and horizontal and vertical transverse dispersivity ratio were taken as 0.0002, 0.03, and 0.01, respectively (Nyika 2021).

For validation of the model formulated in this study, LD is taken as the main calibration parameter due to its significant impact on contaminants’ transport processes. LD defines the path and travel distance of contaminants in groundwater. Among the two samples, calibration was done corresponding to the laboratory test result of a sample collected from location 1. Calibration was undertaken based on a trial and error process by adjusting the LD value until the model-predicted total dissolved solids (TDS) matches the laboratory test result at location 1. Accordingly, the calibrated LD value was used for the validation task corresponding to the laboratory test result of the sample collected from location 2. The landfill site in Rajshahi was formed in 2002, which makes the landfill more than 20 years old. So, the tentative value of TDS would be less than 1,000 mg/L (Mukherjee et al. 2015).

The water table elevation was found higher in the northwest part and lower in the south part of the study area as shown in Figure 4(a). It depicts that all the groundwater of the study area flows towards the south of the study area and joins with the Padma River.

Figure 4(b) shows the predicted groundwater table elevation map considering the effect of precipitation and evapotranspiration after 50 years. In this map, the north of the study area shows a high elevation of the groundwater table and flows toward the river, which confirms that the aquifer is feeding the river (al Mamunul Haque et al. 2012b). The groundwater table may fluctuate with the seasonal demand, which further makes the model more complicated, hence, the modeling is performed without considering the seasonal effects.

The contaminant transport modeling is done for five heavy metals that can be found in groundwater and have significant long-term effects on human health if not treated properly. For each contaminant, the developed model is applied for predicting the travel distance after 1, 3, 5, 10, 15, 20, and 50 years. The areas where the concentration of the contaminants is more than their standard limit are considered the alarming zones.

Table 3 showed the variation in TDS concentration between the model projection and laboratory testing. As sample location 1 is closer to the landfill site, its groundwater quality is less affected by other contamination sources. Sample location 2 is far from the landfill site. As a result, other sources of groundwater pollution surrounding the site also affect the quality thereby increasing TDS concentration in laboratory tests compared to model prediction.

All the model projection results show that contaminant migrates toward the direction of groundwater flow (Mescher 2018). The intensity of pollution over time decreases as pollution migration increases (Han et al. 2016). The areas that are close to the source and along the direction of groundwater flow possess a greater risk of groundwater pollution (Singh & Garg 2016).

Table 4 shows the summary of empirical equations for each contaminant. All the equations follow a logarithmic trend where the year (x) must be counted from 1 and above. The coefficient of determination (R²) for each equation is more than 0.9 which supports the acceptability of these equations.

As the landfill site is located at Ward 17, it becomes more susceptible to groundwater contamination (Figure 12). Hence, all the selected contaminants are adequately available resulting in a higher variation in concentration from their standard. As seen in Figures 12–15, the concentration of contaminants increases rapidly with time in the initial years. After 20 years, the concentration variation becomes somewhat constant due to the dilution of groundwater after traveling long distances.

The Cr, Pb, and Mn enter the groundwater of Ward 18 after 5 years while Cu and Zn enter after 10 years. Pb and Mn enter the groundwater in Ward 16 after 5 years, Cu and Cr after 10 years, and Zn after 15 years.

The variation in the concentration from the standard for each contaminant shows an increasing trend. For each ward, the Cr, Pb, and Mn vary more than 90% from their standard compared to other contaminants. This suggests that the groundwater is more affected by these three pollutants. The concentration of Cu varies between 50 and 80% from their standards. For Zn, the concentration does not exceed 40% of its standards. Helal 2013 exhibits that, in most of the wards during 2009 and 2010, the concentration of Pb and Mn were higher than the advised limit for drinking water, although Cu and Zn concentrations were below the recommended range.

It is obvious that groundwater flow is a complex issue and depends on a lot of factors. For instance, excessive groundwater use (due to overpopulation, logging, industrial development, etc) would also affect groundwater fluctuations and influence the reported sea level rise (Koutsoyiannis 2020). Precipitation, groundwater table, evapotranspiration, and other input parameters of the model have a long-term impact on one another, which is known as the Hurst phenomenon (Koutsoyiannis 2004). This phenomenon appears to have an impact on the hydrological processes which is more apparent in longer time series data. In order to improve the model's accuracy, it should be analyzed.

In this study, a groundwater flow and contaminant transport model were developed for the RCC for five selected pollutants including Cr, Cu, Pb, Zn, and Mn for different projection years, e.g. 1, 3, 5, 10, 15, 20, and 50.

• The outcomes from groundwater flow modeling using MODFLOW confirm that the groundwater moves from the northwest part to the south of the study region, finally meeting the Padma River stream.

• The findings from the groundwater contamination transport modeling using MT3DMS indicate that the affected area due to contaminants increases with time.

• An empirical equation of logarithmic trend was developed for each pollutant to project the maximum travel distance of the contaminant within the desired time. It is evident that during the initial years, the migration distance of contaminants varies rapidly with time. With the increasing time, the contaminant starts to dilute causing the curve to become shallow.

• Ward numbers 14, 16, 17, and 18 of the RCC are more susceptible to groundwater pollution by heavy metals.

• Among the contaminants, the concentration of Cr, Pb, and Mn varies more than 90% from their standards. It suggests that these three contaminants have a greater impact on groundwater pollution. The concentration of Cu varies between 50 and 80% from its standards. For Zn, the concentration does not exceed 40% of its standard.

Before using a piece of land as a landfill, an adequate liner should be constructed to reduce groundwater contamination caused by leachate. The everyday usage of various chemicals should be decreased, and we should switch to using natural alternatives. Especially the products which contain high levels of Mn (nonferrous alloy), Cr (stainless steel), and Pb (ceramic tiles, paints, and lipstick) should be avoided. Hence, concerned authorities should be aware of it and precautions must be taken like regular monitoring of the groundwater and managing the landfill site to protect groundwater from being polluted.

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

The authors declare there is no conflict.