A Process Model for Leachate Treatment in Adsorbent-Amended Constructed Wetlands

Abstract

This paper introduces a model of constructed wetlands (CWs) for on-site landfill leachate treatment, emphasizing the efficacy of zeolite and biochar as adsorbent media to enhance pollutant removal. While zeolite improves nitrification, biochar enhances denitrification and organic matter removal. The model employed a series of continuously stirred reactors and was calibrated and validated using data from 2 parallel CW mesocosms. Each mesocosm incorporated a 250-L vertical subsurface flow CW (V-CW) followed by a 440-L horizontal subsurface flow CW (H-CW). One was an unamended control with gravel media and the other was amended with zeolite in V-CW and biochar in H-CW. The model considered hourly changes in water storage, temperature, and transformations of nitrogen species, chemical oxygen demand (COD), and dissolved oxygen. Model robustness was demonstrated by low normalized root mean square error values (<0.3) for effluent ammonium nitrogen and COD concentrations. The model, validated and subject to a sensitivity analysis, accurately predicted pollutant removal efficiency (<13.3% error for unamended, <3% for adsorbent-amended), making it a valuable tool for evaluating full-scale CW dynamic performance under variable conditions.

Keywords Process-based numerical model; Landfill leachate treatment; Ammonia nitrogen; Zeolite, Biochar; Sub-surface flow wetland

1. Introduction

Approximately 60 million m³ of landfill leachate are generated in the United States each year (Lang et al. 2017). Landfill leachate contains high concentrations of ammonia (NH₃), recalcitrant organic matter, and color (Yang et al. 2023), and must be treated to prevent contamination of ground and surface waters. Previous studies have demonstrated that constructed wetlands (CWs) are cost-efficient and sustainable methods for on-site landfill leachate treatment (Vymazal and Kröpfelová 2009) due to their low complexity and low energy requirements compared with conventional treatment technologies (Arias and Brown 2009; Kadlec and Wallace 2008).

A review of 85 research papers on performance of CWs found that hybrid CWs that combine vertical subsurface flow CWs (VF-CWs) and horizontal subsurface flow CWs (HF-CWs) were the most effective at removing organic matter and total nitrogen (TN) from landfill leachate (Bakhshoodeh et al. 2020). VF-CWs favor aerobic processes, such as nitrification, because oxygen (O₂) diffuses into the unsaturated media, particularly when operated in batch and/or with pulses of influent. HF-CWs are saturated and develop anoxic conditions that favor denitrification. However, due to the limited adsorption capacity of gravel substrates, the performance of gravel-bed CWs in wastewater treatment is frequently insufficient (Kizito et al. 2017). This demands the use of alternative substrate materials as adsorbents and biofilm carriers.

One of the challenges in treating landfill leachate is that high concentrations of free NH₃ can inhibit nitrifying bacteria (Kim et al. 2006). Recent studies have suggested that addition of an adsorbent, such as zeolite, can prevent nitrification inhibition and increase CWs’ efficiency (Gao et al. 2022). Natural zeolites, such as clinoptilolite, are porous aluminosilicate minerals with a high cation exchange capacity and a strong selectivity for ammonium nitrogen (NH₄⁺-N). The adsorbed NH₄⁺ is oxidized to nitrate (NO₃^-) by nitrifying biofilms attached to the zeolite surface. This bioregenerates the ion exchange capacity of the medium without the need for regenerant brines or addition of fresh zeolite (Rodriguez-Gonzalez et al. 2020).

Biochar is an adsorbent that can enhance TN and chemical oxygen demand (COD) removal in wastewater treatment (Gao et al. 2022). Biochar is a carbon-rich, porous material created from thermal decomposition of organic matter, such as agricultural biomass waste, under low-O₂ conditions. Biochar’s heterogeneous surface and numerous pores allow for enhanced organic matter adsorption, O₂ penetration and microbial colonization, resulting in enhanced COD and TN removal in CWs (Deng et al. 2021). Use of biochar in CWs can also enhance plant growth and pollutant uptake (Huang and Gu 2019).

A mathematical model was developed previously involving the physicochemical and biological processes of ion-exchange-assisted biological treatment (Aponte-Morales et al. 2018). Moreover, Payne (2018) developed a kinetic model, using the Homogenous Surface Diffusion Model (HSDM), to show the adsorption of NH₄⁺ on zeolites. Key assumptions included that the adsorbent is homogeneous and spherical, and that the diffusion process follows Fick’s laws. Recently, a model was developed to predict biochar adsorption of organic pollutants using bench-scale experimental data (Kearns et al. 2020). The study demonstrated 2 different adsorption models: HSDM and the pore surface diffusion model (PSDM). These adsorption models use equilibrium parameters, such as Freundlich isotherm coefficients, as well as coefficients for external and internal mass transfer. Other research has also been published on biochar adsorption models using machine learning (Ukoba and Jen 2022). However, integrating zeolite and biochar in a CW performance model for treating landfill leachate is a novel idea for which there is limited research to date.

Several physical, chemical, and biological processes occur simultaneously and interact with each other in CWs (Kadlec and Wallace 2008). These processes can include N and organic carbon (C) transformations, dissolved O₂ (DO) cycling, and the effects of dispersion, variability, and stagnant areas (Langergraber and Šimůnek 2005). Analyzing these complex relationships is challenging with observational data and laboratory experiments alone, thus process-based numerical models can greatly help with the design of adsorbent-amended CWs. However, accurate predictions can be obtained only when the parameters come from systems that are operating under similar conditions, such as climate, wastewater composition, porous filter media, and plant species (Langergraber 2011). These models include, for instance, the Constructed Wetland 2-Dimensions (CW2D) and FITOVERT models, which have already been tested (Langergraber et al. 2009). Models that simulate key CW processes have been well demonstrated to be more precise in predicting effluent concentrations (Yuan et al. 2020). Currently, numerical process models are acknowledged as powerful tools to better understand CW physicochemical and biological processes and improve their design and performance (Samsó et al. 2014). 

A systematic and methodological approach to model subsurface flow CWs assumes plug flow with longitudinal dispersion and variable saturation. However, such a model involves complex mathematical equations regarding time and vertical direction, making it difficult to apply in CW design (Langergraber et al. 2009). On the other hand, the Continuous Stirred Tank Reactor (CSTR) modeling strategy provides a reasonable description of the numerous processes. A prior study demonstrated that non-ideal plug flow simulations can accurately represent large-scale subsurface flow CWs (Pálfy et al. 2017). Currently, CW process models consider the N, C, and O₂ cycles, and temperature and water storage variations (Langergraber 2011). However, no prior studies have developed numerical models for adsorbent-amended CWs for landfill leachate, which can be useful tools for CW design and performance evaluation.

The goal of this study was to develop a numerical process model to predict adsorbent-amended CW performance under varying landfill leachate loads and concentrations. The model considers hourly changes in water storage; N, COD, and DO transformations; temperature; and effect of zeolite and biochar amendment. The model was programmed in Python 3.7 and calibrated and validated using data collected from mesocosm-scale systems operated at a landfill (Yang et al. 2023). An uncertainty analysis was performed to predict system performance under varying hydraulic loading rates (HLRs). To the best of the authors’ knowledge, this study presents the first attempt at integrating adsorptive media mechanisms into a dynamic model of CWs for landfill leachate treatment.

2. Methods

2.1. Mesocosm Scale Systems

Details on design, operation, and performance of the mesocosms can be found in Yang et al. (2023). Briefly, 2 mesocosm-scale hybrid VF-HF CWs were set up at the Southeast Landfill in Hillsborough County, Florida, United States (Figure 1). An unamended hybrid VF-HF CW with a conventional gravel medium was used as a control. The adsorbent-amended VF-CW contained 10% zeolite and 90% gravel by volume (23 kg zeolite with a porosity of 0.3) to facilitate nitrification. The adsorbent amended HF-CW contained 13% biochar and 87% gravel by volume (2.6 kg biochar with a porosity of 0.6) to enhance COD removal, denitrification, and plant growth. The systems were planted with cattails (Typha latifolia) and cordgrass (Spartina alterniflora). A theoretical hydraulic residence time (HRT) of 11 days was maintained by pulsing 24 L of raw landfill leachate per day through each system (15 mins every 2 h). Sensors for water level were installed in the HF-CWs, and the data were used to adjust the water balance. Influent and effluent samples were collected and analyzed weekly over approximately 2 years. Overall, the pilot system demonstrated excellent performance, with zeolite addition increasing nitrification rates by 16% to 93% and biochar addition increasing removal of organic matter by 21% to 44% (Yang et al. 2023). 

Fig. 1 Schematic diagram of mesocosm-scale system. Modified from Yang et al. (2023).

2.2. Model Overview

A numerical process model was developed to predict unamended and adsorbent-amended hybrid VF-HF CW performance under varying wastewater loads and concentrations. This model considers hourly changes in water storage; N, COD, and DO transformations; and effects of temperature and adsorbent addition. Although the measured pH for the reactors ranged from 7.69 to 8.36—the model constant neutral pH—it is well known that pH affects biodegradation rates (Kadlec and Wallace 2008); therefore, this is a limitation of the model. The VF-CW was modeled as an aerobic system and the HF-CW was modeled as an anoxic system to account for the effects of DO on kinetic parameters. The model was developed in Python 3.7 and was informed by data collected from the mesocosm-scale system. Model parameter values were obtained using a combination of model optimization and literature review (see tables S2 and S3 in the Supplementary Material). The steps for developing the numerical process model are synthesized in Figure 2. A more detailed overview of the interactions among processes is presented in in the Supplementary Material, Figure S2.

Fig. 2 Landfill leachate CW model development and evaluation steps.

2.2.1. Water Balance Module

To solve for hourly water volume of both VF-CW and HF-CW, the following mass balance equation was used:

where V is the volume (m³), t is time (h), Q_i is the inflow rate (m³/h), Q_o is the outflow rate (m³/h), P is precipitation (m/h), A_s is the wetland surface area (m²), and ET is evapotranspiration rate (m/h). The hydraulic loading rate (HLR) associated with Q_i was computed from 24 L of water, fed to the system intermittently. The adjusted Thornthwaite method was employed to estimate evapotranspiration (ET) in this study, with the use of hourly temperature data:

where T is the hourly air temperature (°C) and a is derived from a polynomial equation as a function of I, the annual heat index.

To understand the hydraulics in the VF-CW, preliminary computer simulations were performed using Hydrus-1D (Šimůnek et al. 2016). The hydraulic loading rate of the whole system was 0.0039 cm/min and the simulated bottom flux was 0.004 cm/min. Therefore, the VF-CW outflow rate was assumed to be equal to the inflow rate for the same tank. The HF-CW was considered to be a saturated system and the inflow rate of this tank was equal to the outflow rate of the vertical tank. The HF-CW tank’s outflow rate was assumed to be equal to the HF-CW inflow rate.

2.2.2. Effect of Temperature

Temperature affects the biochemical processes regulating the removal efficiency of CWs (Kadlec and Reddy 2001). Although these effects are particularly important in N cycle reactions (mineralization, nitrification, and denitrification), the effects of temperature were incorporated in evaluating the kinetic parameters for all modules by adjusting the rate constants for temperature dependence using the modified Arrhenius equation:

where k₂₀ is the removal rate constant at 20 °C, T is any given temperature (°C), k is the removal rate constant at the given temperature (°C), and θ is the temperature coefficient. The θ values for organic N mineralization values were assumed to be 1.08 based on model calibration. The value of θ for nitrification is close to 1.10 (Kadlec and Reddy 2001). The values for all the rest of the temperature coefficients were assumed to be 1.04 based on model calibration.

2.2.3. Effects of Dissolved Oxygen

Dissolved oxygen is a critical factor in N cycling. In nitrification, aerobic bacteria convert NH₄⁺-N to NO₃^--N in the presence of O₂, while in denitrification, anaerobic bacteria convert NO₃^--N to N gas in the absence of O₂. Since the availability of O₂ affects the rates and efficiency of these processes, the DO mass balance was solved:

where DO is the dissolved O₂ concentration (mg/L), HT_res is the respiration rate by heterotrophs (mg/L/h), AT_res is the respiration rate by autotrophs (mg/L/h), and J_O2 is the O₂ mass flux (mg/L/h), which is calculated using the following equation:

where k_R is the temperature-dependent mass transfer coefficient (h⁻¹), DO_s is the DO concentration when saturation is reached (mg/L), and DO is the average DO concentration measured in the water column (mg/L). The availability of DO determined if the system was aerobic or anaerobic, and thus also determined the values of rate constants to be used in the modules (see the Supplementary Material for parameter values). 

Previous research has supported the use of the Monod relationship for incorporating the effects of DO on determining microbial kinetic parameters (Zakir Hossain et al. 2019):

where k is reaction rate constant (h⁻¹); k_max is the temperature-dependent maximum reaction rate constant (h⁻¹); K_d is the half-saturation constant (mg/L), which is also dependent on temperature; and DO is the DO concentration (mg/L).

2.2.4. Pollutant Mass Balance

Removal of N in subsurface flow CWs relies mostly on microbially mediated transformations (e.g., ammonification, nitrification, denitrification), with plant uptake having a limited impact (Yang et al. 2023). Nitrification is mostly favored by aerobic conditions and depends on O₂ availability and pH, whereas the denitrification process is dependent on low levels of DO and the availability of organic C as electron donor (Langergraber and Šimůnek 2005). 

Nitrification is often modeled as a single-step process in CWs, where NH₄⁺-N is oxidized to NO₃^--N, ignoring nitrite (NO₂^-) production (Mayo and Bigambo 2005). This assumption is valid when NO₂^--oxidizing bacteria are uninhibited and can quickly convert NO₂^- to NO₃^-, eliminating the need to explicitly model the intermediate step. Denitrification is a 4-step process (NO₃^- → NO₂^- → NO → N₂O → N₂); however, it is often assumed that denitrification is a single-step process for model simplification (Barker and Dold 1997).

Mass balances for N species were simulated using a steady-state CSTR in-series model (or Tanks in Series [TIS] model), which is widely used for simulating continuous flow process (Wong et al. 2006; Chapra 2008; Kadlec and Wallace 2008). A TIS model assumes a series of CSTRs, in which the effluent from one tank is the influent of the next tank. The number of TIS is based on the residence time required for the desired treatment process. CWs exhibit deviations from idealized plug flow or CSTR reactors due to short-circuiting and dispersion. Therefore, the TIS model is more suitable for describing hydraulics and pollutant behavior in CWs; it serves as an intermediary between the plug flow and CSTR models and considers the number of TIS as a mathematical representation of the hydraulic performance of the CW, rather than the physical number of cells within the system. In the current model, we assume the number of TIS to be 3, based on model optimization and commonly used values in the literature (Kadlec and Wallace 2008). The model also assumes a constant flow rate in and out, immediate and complete mixing within each tank, uniform temperature throughout the reactor, and no differences in velocity or stagnant areas. The general mass balance for each tank in the TIS is represented as:

where Q_i is the inflow rate, C_i is the influent concentration, Q_o is the outflow rate, C is the effluent concentration, r is the reaction rate, and V is the water volume in each h solved with the water balance.

The mass balances for N species are solved in a similar manner. The following mass balance equation was used for organic N:

where (OrgN)_n−1 is the effluent concentration of organic N (mg/L) in (n − 1)^th tank in TIS model, (OrgN)_n is the concentration of organic N (mg/L) in the n^th tank in TIS model, k_pd is the rate constant for plant decomposition (h⁻¹), k_m is the rate constant for mineralization/ammonification (h⁻¹), V is the tank volume (m³), and Q is the flow rate (m³/h). The following equation was used for calculating effluent NH₄⁺-N concentration:

where (NH₄⁺)_{(n − 1)} is the effluent concentration of NH₄⁺-N (mg/L) in (n − 1)^th  tank in TIS model; (NH₄⁺)_n is the concentration of NH₄⁺-N (mg/L) in n^th tank in TIS model; k_n is the first-order rate constant for nitrification (h⁻¹), which is dependent on availability of DO and temperature; k_pu,ammonia is the first-order rate constant for plant uptake of NH₄⁺ (h⁻¹); V is the tank volume (m³); and Q is the flow rate (m³/h). The following equation was used to calculate effluent nitrate concentration:

where (NO₃^-)_n is the effluent concentration of NO₃^- (mg/L) in n^th tank in the TIS model; (NO₃^-)_n−1 is the concentration of NO₃^- (mg/L) in (n−1)^th tank; k_dn is the first-order rate constant for denitrification (h⁻¹), which is dependent on temperature; k_pu,nitrate is the rate constant for plant uptake of NO₃^- (h⁻¹); V is the tank volume (m³); and Q is the flow rate (m³/h).

Models of microbially facilitated nitrification-denitrification reactions often follow a modified Stover-Kincannon model, which is used to calculate the rate of change in substrate concentration at steady state as a function of loading rate (Nga et al. 2020). The following mass balance equation was used to estimate the effluent COD concentration:

where (COD)_n is the effluent concentration of COD (mg/L) in n^th tank in TIS model, (COD)_n−1 is the concentration of COD (mg/L) in (n − 1)^th tank, u_max is the temperature-dependent specific maximum substrate utilization rate (mg/L/h), K_b is the temperature-dependent saturation constant (mg/L/h), V is the tank volume (m³), and Q is the flow rate (m³/h).

2.2.5. Adsorption Media Module

Adsorption models for zeolite and biochar (Table 1) were carried out using the HSDM, an approach that has shown to be effective in simulating adsorption kinetics (Payne 2018). HSDM assumes the adsorption media to be homogeneous and the particles to be spherical. 

Table 1 Properties of adsorbent media

Following the HSDM approach, the NH₄⁺-N mass balance equation was modified to incorporate adsorption: 

where J_NH4⁺ is the flux of NH₄⁺-N from the bulk liquid to the solid phase (mg/m²/h) and A_adsorbent is the interfacial area of the adsorbent sites (m²/m³). The adsorption flux J_NH4⁺ is calculated as follows:

where ρ is the particle density (g/m³), D_s is the surface diffusivity (m²/h), q_NH4⁺ is the maximum adsorption capacity (mg/g), and r is the radial coordinate. Here, r was assumed to be equal to R, the radius of a zeolite particle (m). Also, since the HSDM assumes particles are spherical, A_adsorbent is calculated as follows:

where M is the mass of zeolite in the system (g) and V is the volume (m³). The COD mass balance was also modified to consider adsorption: 

where COD_i is the initial concentration of COD (mg/L), u_max is the specific maximum substrate utilization rate (mg/L/h), K_b is the saturation constant (mg/L/h), V is the tank volume (m³), Q is the flow rate (m³/h), J_COD is the adsorption flux (mg/m²/h) and A_adsorbent is the surface area of adsorbent (m²/m³). Here J_COD and A_adsorbent are calculated using similar equations to (13) and (14), respectively.

2.3. Model Calibration and Validation

The mesocosm-scale CWs were constructed, planted, and put into operation on 28 August 2020. The model was calibrated using 6 months of data from the CWs after the plants were well established (25 November 2021 to 14 May 2022) and validated using another 6-month dataset (15 May 2022 to 13 October 2022). Since samples were collected weekly, influent NH₄⁺-N and COD concentrations were assumed to vary linearly between data points. Auto-calibration was carried out using Platypus, a multi-objective evolutionary algorithm (MOEA) open-source framework for optimization in Python (Hadka 2017). Initial ranges of input parameters came from the scientific literature when available (Martin and Reddy 1997; Welander et al. 1997; Lahav and Green 2000; Alvarez and Becares 2008), then were adjusted based on the sensitivity analysis (described in section 2.4). Hourly rainfall and temperature data collected from a local weather station (Lakeland Linder International Airport Station in Lithia, Florida, United States) were used for precipitation and evapotranspiration (Weather Underground 2024). Water level sensors were used to calibrate the water balance. Pollutant concentrations, measured weekly, were assumed to be constant between collection times. Model fitness was evaluated using the normalized root mean square error (NRMSE) as presented by Ranatunga et al. (2017):

where y_simulated is the value simulated by the model, y_experimental is the concentration value observed by the experimental data, and N is the number of the experimental data points. NRMSE is a widely recognized metric for evaluating the accuracy of prediction tools (Ranatunga et al. 2017). Its advantage lies in its unit independence, making it suitable for comparing models across different scales and applications. A perfect model would achieve an NRMSE value of 0, indicating no prediction error. Previous studies have regarded an NRMSE value below 0.3 as indicative of good model performance (Verhamme et al. 2016; Lei et al. 2023; Sadeghi et al. 2024).

Model validation was performed based on the mean percent removal efficiency. For this, removal efficiency was calculated at each data point for both the simulated and observed scenarios. The mean of these removal efficiencies was then computed for both datasets. Using these mean values, the percent error between the observed and simulated removal efficiencies was calculated. This method provided a quantitative measure of the model’s predictive capability in terms of removal efficiency.

2.4. Parameter Sensitivity Analysis

A sensitivity analysis was carried out to understand how the model’s key outputs (effluent concentrations) reacted to various model parameters. The analysis was done using the Sobol method (Saltelli 2002) on 12 and 11 different parameters for NH₄⁺-N and COD adsorption modules, respectively. The Sobol method is a variance-based sensitivity analysis method that quantifies the contribution of each input parameter to the output variance. It is based on decomposing the output variance into components that can be attributed to individual or combinations of input parameters. This method determines first-order and total-order indices to evaluate each parameter’s individual and cumulative contribution. The first-order sensitivity index measures the contribution of each parameter alone, while the total-order sensitivity index measures the contribution of each parameter and its interactions with other parameters. If the total-order indices are significantly higher than the first-order indices, that means there are higher order interactions among the parameters. Both NH₄⁺-N and COD adsorption modules were analyzed using over 10,000 model runs to evaluate the global variability of the parameters. Indices reported in the results represent the contribution that each parameter makes to the total variance and include all variance caused by its interactions with any other input parameter. Overall, this method provides a powerful tool for understanding the sensitivity of a model to its input parameters and for identifying which parameters are most critical for improving the model’s performance.

3. Results and Discussion

3.1. Parameter Sensitivity Analysis

Results from the parameter sensitivity analysis for the NH₃ adsorption-media processes show that surface diffusivity of zeolite and the maximum NH₃ adsorption capacity were the most sensitive parameters (Table 2). In other words, slight changes in adsorbent properties can have the greatest effect on the predicted NH₄⁺-N concentration. NH₄⁺-N adsorption appeared to be least sensitive to temperature among the parameters evaluated, which is expected as seasonal temperature differences in Florida are not as extreme as in higher latitudes. Furthermore, the sensitivity analysis for COD adsorption indicated that the most sensitive parameters are biochar maximum adsorption capacity and surface diffusivity (Table 3), demonstrating that biochar properties are critical for COD removal from landfill leachate in CWs. 

Table 2 Sensitivity analysis results for NH₄⁺-N adsorption

Table 3 Sensitivity analysis results for COD adsorption

3.2. Unamended CWs Calibration and Validation

Parameters for the unamended CWs were calibrated first to fit the simulation results with the experimental data. The average influent concentrations of NH₄⁺-N and COD were 357 mg/L and 446 mg/L, respectively. Parameters were maintained within the literature range most of the time; however, some model parameters were not found in previous studies and thus, best fit values were chosen both for VF-CW and HF-CW. The calibrated parameters and literature ranges are presented in Table S1 in the Supplementary Material.

The simulated results for NH₄⁺-N for both the VF-CW and HF-CW capture well the NH₄⁺-N reduction trend, as the NRMSE values for the validation period are 0.37 and 0.28 for the VF-CW and HF-CW, respectively (Figure 3). These findings regarding NH₃ removal efficiency indicate a higher efficiency in the VF-CW, with a removal efficiency of 43.1%, compared to the HF-CW, which showed an efficiency of 21.9%. These findings align with observations of low error percentage (7.8% and 8.9% for VF-CW and HF-CW, respectively). The aerobic conditions in the VF-CW favor nitrification, which is a major step in the NH₄⁺-N removal process. On the other hand, anaerobic conditions with limited O₂ supply impede the growth and activity of nitrifying bacteria (Schmidt et al. 2002).

Fig. 3 Time series of observed and simulated NH₄⁺-N concentrations in the unamended system.

The simulated and experimental results for NO₃^--N are shown in Figure 4. The results for NO₃^--N concentration calibration/validation show lower NRMSE values for both VF-CW (NRMSE = 0.23) and HF-CW (NRMSE = 0.29) indicating good model fitness. NO₃^--N accumulation occurs as a function of the imbalance between nitrification and denitrification (Vymazal 2005). The nitrification rate constants were higher in the VF-CW (0.70 h⁻¹) compared to the HF-CW (0.006 h⁻¹). These findings indicate that the VF-CW system, with an average concentration of 63.74 mg/L, had a higher accumulation of NO₃^--N because of effective nitrification and limited carbon availability for denitrification compared to the HF-CW system, which had an average concentration of 44.97 mg/L. 

Fig. 4 Time series of observed and simulated NO₃^--N concentrations in the unamended system.

The relatively low NRMSE for COD in the unamended system for both VF-CW (0.12) and HF-CW (0.17) suggests that the models fit the COD observations well (Figure 5). The VF-CW showed a slightly improved simulated COD removal efficiency of 17.4% compared to the HF-CW, which had a removal efficiency of 13.4%. This aligns well with the observed data, which show a mean COD removal of 16.5% in the VF-CW and 11.5% in the HF-CW. Typically, aerobic conditions are more favorable for organic matter degradation, as microorganisms utilizing O₂ as a terminal electron acceptor play a key role in the process (Najafi Savadroudbari et al. 2021). COD also plays a crucial role in denitrification, as organic matter in the influent acts as an electron donor and impacts microbial metabolism (Zhao et al. 2021). The low COD removal efficiency can be attributed to the presence of high levels of recalcitrant organic matter in landfill leachate, indicating that further treatment, such as activated C adsorption or advanced oxidation, may be required (Archibald et al. 1998). 

Fig. 5 Time series of observed and simulated COD concentrations in the unamended system.

3.3. Adsorbent-amended CWs Calibration and Validation

The parameters for the adsorbent-amended CWs were calibrated to fit the simulation results with the experimental data. Maximum adsorption capacity and surface diffusivity of both zeolite and biochar were the main calibration parameters. After calibration, parameter values were mostly within the literature range. However, some parameters were not found in previous studies and thus, best fit values were chosen both for VF-CW and HF-CW. As for the surface diffusivity of zeolite, the value was obtained within the literature range (4.77×10⁻¹² m²/s). However, the biochar surface diffusivity value was not found in the literature, so the value 3.37×10⁻¹¹ m²/s was obtained from model optimization. The calibrated parameters and literature ranges are presented in the Supplementary Material, Table S2.

The simulated results for NH₃ for both the amended VF-CW and HF-CW are presented in Figure 6. The NRMSE value for the validation period for VF-CW is 0.36 and that of the HF-CW is 0.23. The values are higher (worse) compared to the unamended CW results, which is expected because of the complex adsorption kinetics involved in these CWs. The simulated removal efficiencies (discussed in more detail in section 3.4), however, matched closely with the observed data; zeolite in the VF-CW helps reduce nitrification inhibition due to high free NH₃ concentrations in the leachate (Aponte-Morales et al. 2018). This likely contributed to the higher NH₄⁺-N removal rate in the VF-CW system compared to the HF-CW system, as biochar is not as effective as zeolite in removing cationic compounds from wastewater. However, the presence of biochar in the HF-CW system can promote plant growth, leading to NH₄⁺-N removal through plant uptake (Rawat et al. 2019).

Fig. 6 Time series of observed and simulated NH₄⁺-N concentrations in the adsorbent-amended system.

The simulation results for COD illustrate there is not much difference in COD removal between the unamended and adsorbent-amended VF-CWs, which is expected as zeolite does not have much adsorption capacity for COD (Figure 7). The NRMSE values are low for both VF-CW (NRMSE = 0.17) and HF-CW (NRMSE = 0.24), indicating good model fitness. The simulation results suggest that the HF-CW system has a much higher COD removal efficiency (38%) compared to the VF-CW system, which only shows a removal efficiency of 16%. This aligns with observed data, as biochar has been found to be effective in removing organic compounds through adsorption. Unlike biochar, zeolite is primarily utilized as a cation exchange material and does not exhibit the same level of organic compound adsorption. Furthermore, biochar can create a favorable environment for organic matter degrading microorganisms, resulting in increased COD removal compared to zeolite, which does not have similar impacts on microbial activity (Dai et al. 2021).

Fig. 7 Time series of observed and simulated COD concentrations in the amended system.

3.4. Removal Efficiency Validation

The mean removal efficiency for each media-amended tank was calculated to compare the experimental and simulated results (Table 4). For the unamended system, the model was able to predict the removal efficiency within a range between 5% and 13%. The removal efficiency for NH₃ was higher for the VF-CW (43.1%) compared to the HF-CW (22%), which matches the observed trends. However, the removal efficiency of COD shows little disparity between the VF-CW and the HF-CW, indicating similar performance between the 2 systems. For the amended system, the removal efficiency for NH₄⁺-N was higher for the VF-CW (67%) compared to the HF-CW (33%), which matches the observed trend. As for COD, the removal efficiency of COD was much higher in the HF-CW (38%) compared to the VF-CW (16%), which is expected due to the addition of biochar. Overall, the results show that the amended CW model can predict the removal efficiency within an error range of 0% to 3%.

Table 4 Removal efficiency comparison

3.5. Effects of Media Amendment on Leachate Pollutant Removal

A direct comparison of the effluent concentrations from the unamended and amended CWs elucidates the model’s capability of simulating the effects of adsorbent amendment (Figure 8). The results for NH₄⁺-N show a 24.3% increase in simulated removal efficiency for NH₃ in the amended VF-CW compared to the unamended system, and a 46.6 % increase in the amended HF-CW compared to the unamended system. As for the COD removal, there was a 24.9% increase in simulated removal efficiency in the amended HF-CW compared to the unamended system. However, simulated COD removal efficiency did not vary much between the amended and unamended VF-CW (16% vs. 17%). This is expected as zeolite does not have much affinity for COD and this has been reflected in the model.

Fig. 8 Effects of adsorbent amendment on effluent concentration.

3.6. Effects of Varying Hydraulic Loading Rate

An analysis with varying HLR was performed to understand the effect of this important design parameter on cumulative pollutant removal efficiency of the amended system (Figure 9). The simulated results were also compared against the available observed data for verification purposes. Note that while the range of HLR simulated covered a range of typical operating conditions for CWs, these simulations covered a larger range than the experiments. As expected, there was a tendency for the removal efficiency to decrease with increasing HLR, up to a value around 17 cm/day, after which the efficiency stabilized. For instance, when increasing the HLR from 2.9 cm/day to 17.2 cm/day, NH₃ removal efficiency decreased from 99% to 67%, but the removal efficiency remained nearly the same at an HLR of 25.7 cm/day. Moreover, COD removal efficiency decreased from 70% to 32% when HLR changed from 2.9cm/day to 17.2 cm/day, but the efficiency was approximately 30% at an HLR of 25.7 cm/day. Increasing HLR can decrease pollutant removal efficiency for several reasons. As the HLR increases, for instance, HRT decreases, leading to shorter contact time between pollutant and treatment system. Also, an increase in HLR can decrease the O₂ supply per unit volume of landfill leachate, leading to reduced pollutant removal efficiency (Zulfikar et al. 2022).

Fig. 9 Removal efficiency for varying HLRs.

3.7. Limitations

Several limitations in the model could be addressed in future studies:

• There were not sufficient data available to calibrate the model for organic N, which may have affected the accuracy of the predictions.

• The model does not consider the bio-regeneration of biochar and zeolite, which can lead to higher errors in the prediction of effluent concentrations. This is particularly relevant for biochar, as the experimental data showed a decline in adsorption capacity after 160 days, which was not captured by the model (Figure 7).

• As already stated, the model assumes constant pH. The model also does not consider the limiting effects of C availability on the denitrification process and the HF-CW system has been modeled as fully anoxic to account for the effects of low levels of DO on the denitrification process for model simplification.

• The model could be improved by performing a tracer study to make more accurate assumptions about the number of tanks in the mixing model. The reaction rates utilized in the model were acquired through a combination of literature review and model optimization. However, more accurate results could have been achieved if the reaction rates were obtained through batch reactor experimental results.

• Model fitness metrices other than NRMSE need to be explored, as well, to compare and evaluate the model results.

4. Conclusions and Recommendations

A numerical process model was developed to predict the performance of CWs in response to varying wastewater loads, concentrations, design parameters, and adsorbent-media characteristics. This model produces simulation results for both unamended and adsorbent-amended hybrid VF-HF CWs using a TIS model for NH₄⁺-N and COD removal. This model considers a range of factors including water storage; N, COD, and DO transformations; temperature effects; and the enhancement potential of zeolite and biochar amendments. The model was calibrated and validated using data collected from a mesocosm-scale system installed at an operating municipal landfill. A parameter sensitivity analysis was carried out prior to model calibration. Additionally, an uncertainty analysis was conducted to assess the model’s ability to predict system performance under different HLRs. The resulting model is scalable and could provide performance estimates for large-scale CWs designed to treat landfill leachate.

The results of the study indicate that the developed model can successfully simulate the concentrations of NH₄⁺-N and COD in CW effluent under varying conditions. The model’s performance was evaluated using the NRMSE, which was found to be <0.3, demonstrating its efficacy in predicting pollutant removal efficiency. The error in predicting removal efficiency was less than 13% for the unamended system and 3% for the adsorbent-amended system. The results further show that the addition of zeolite and biochar to the CWs significantly improved the removal of pollutants, as there was a 24% and 13% increase, respectively, in simulated removal efficiency of NH₄⁺-N in the adsorbent amended VF-CW and HF-CW compared to the unamended system. Additionally, there was a 25% increase in simulated removal efficiency of COD in the amended HF-CW compared to the unamended system. However, no major difference was observed between the results for the unamended and amended VF-CW for COD removal. A sensitivity analysis demonstrated that slight changes in adsorbent properties can have the greatest effects on the predicted NH₄⁺-N concentration. The analysis for varying HLR in the adsorbent-amended system revealed that by increasing the HLR from 2.9 cm/day to 20 cm/day, NH₄⁺-N removal efficiency decreased from 99% to 63% and COD removal efficiency decreased from 70% to 31%. Further studies are needed to assess the long-term performance of the adsorbent-amended CWs and the economic feasibility of using adsorbents in a full-scale CW.

Supplementary Material

The online version of this article contains a link to supplementary material that includes: Figure S1 Hydraulic Simulation Results for VF-CW; Figure S2 Overview of Module Interactions; Table S1 Van Genuchten Parameters for Hydraulic Simulation; Table S2 Calibrated Values for Adsorbent-Amended CW Model Parameters; Table S3 Calibrated Values for Unamended CW Model Parameters.

Acknowledgments

This work was supported by the Hinkley Center for Solid and Hazardous Waste Management (Grant No. SUB00001960). The authors would also like to thank Larry Ruiz (Landfill Manager), Joseph Kinsey (Landfill Supervisor), Kelvin Melton (Landfill Operator), and Frank Barrella (Landfill Mechanic) at Southeast Hillsborough County Landfill, as well as USF students Thanh (Misty) Lam, Megan Kramer, and Joshelyn Guimaray. Thanks also to Dr. Jay Martin and the anonymous reviewers for their valuable feedback during the peer review process.

Author Contributions Statement

Conceptualization: MEA, SJE; methodology: MEA, LM, IN, SJE; data analysis: LM, IN, XY; laboratory analyses: XY; writing original draft: LM, IN; review/editing original draft: MEA, SJE, XY; resources: MEA, SJE; supervision: MEA, SJE; project administration: MEA, SJE; funding acquisition: MEA, SJE. All authors have read and agreed to the published version of the manuscript.

Conflict of Interest Statement

No conflicts to disclose.

Data and Code Availability Statement

The modeling code is available at
https://github.com/IshfaqunNisa/CW-Model.

Related Publication Statement

Previous drafts of this manuscript were prepared as part of the master’s theses of the first author, Ishfaqun Nisa, and co-author Lillian Mulligan. The current manuscript, however, is a substantial synthesis of those theses. 

ORCID iDs

Ishfaqun Nisa https://orcid.org/0009-0005-9514-2187

Mauricio E. Arias https://orcid.org/0000-0002-8805-6353

Lillian Mulligan https://orcid.org/0009-0000-1918-5971

Xia Yang https://orcid.org/0000-0002-6009-0592

Sarina J. Ergas https://orcid.org/0000-0001-6820-7196

References

Alvarez J, Becares E. 2008. The role of harvest and plant decomposition in constructed wetlands. In: Russo R, editor. Wetlands: ecology, conservation and restoration. Hauppauge (NY): Nova Science Publishers. p. 213–245.

Aponte-Morales VE, Payne KA, Cunningham JA, Ergas SJ. 2018. Bioregeneration of chabazite during nitrification of centrate from anaerobically digested livestock waste: experimental and modeling studies. Environ Sci Tech. 52(7):4090–4098. https://doi.org/10.1021/acs.est.7b06248.

Archibald F, Roy-Arcand L, Méthot M, Valeanu L. 1998. Recalcitrant organic compounds (chemical oxygen demand sources) in biologically treated pulp and paper mill effluents: their fate and environmental impact in receiving waters. Water Environ Res. 70(7):1314–1326. https://doi.org/10.2175/106143098X123705.

Arias ME, Brown MT. 2009. Feasibility of using constructed treatment wetlands for municipal wastewater treatment in the Bogotá Savannah, Colombia. Ecol Eng. 35(7):1070–1078. https://doi.org/10.1016/j.ecoleng.2009.03.017.

Bakhshoodeh R, Alavi N, Oldham C, Santos RM, Babaei AA, Vymazal J, Paydary P. 2020. Constructed wetlands for landfill leachate treatment: A review. Ecol Eng. 146:105725. https://doi.org/10.1016/j.ecoleng.2020.105725.

Barker PS, Dold PL. 1997. General model for biological nutrient removal activated-sludge systems: model application. Water Env Res. 69(5):985–991.
https://doi.org/10.2175/106143097X125678.

Chapra SC. 2008. Surface water-quality modeling. Long Grove (IL): Waveland press.

Dai Z, Xiong X, Zhu H, Xu H, Leng P, Li J, Tang C, Xu J. 2021. Association of biochar properties with changes in soil bacterial, fungal and fauna communities and nutrient cycling processes. Biochar. 3(3):239–254. https://doi.org/10.1007/s42773-021-00099-x.

Deng S, Chen J, Chang J. 2021. Application of biochar as an innovative substrate in constructed wetlands/biofilters for wastewater treatment: performance and ecological benefits. J Clean Prod. 293:126156. https://doi.org/10.1016/j.jclepro.2021.126156.

Gao B, Yang X, Dasi EA, Lam T, Arias ME, Ergas SJ. 2022. Enhanced landfill leachate treatment in sequencing batch biofilm reactors (SBBRs) amended with zeolite and biochar. J Chem Technol Biot. 97(3):759–770. https://doi.org/10.1002/jctb.6964.

Hadka D. 2017. Platypus: A free and open source python library for multiobjective optimization. Fort Collins (CO); [accessed 0000-Mo-Da]. https://github.com/Project-Platypus/Platypus.

Huang L, Gu M. 2019. Effects of biochar on container substrate properties and growth of plants—a review. Horticulturae. 5(1):14. https://doi.org/10.3390/horticulturae5010014.

Kadlec RH, Reddy KR. 2001. Temperature effects in treatment wetlands. Water environment research: a research publication of the Water Environment Federation. 73(5):543–557. https://doi.org/10.2175/106143001x139614.

Kadlec RH, Wallace S. 2008. Treatment wetlands. 2nd ed. Boca Raton (FL): CRC Press.

Kearns J, Dickenson E, Knappe D. 2020. Enabling organic micropollutant removal from water by full-scale biochar and activated carbon adsorbers using predictions from bench-scale column data. Environ Eng Sci. 37(7):459–471. https://doi.org/10.1089/ees.2019.0471.

Kim D-J, Lee D-I, Keller J. 2006. Effect of temperature and free ammonia on nitrification and nitrite accumulation in landfill leachate and analysis of its nitrifying bacterial community by FISH. Bioresource Technol. 97(3):459–468. https://doi.org/10.1016/j.biortech.2005.03.032.

Kizito S, Lv T, Wu S, Ajmal Z, Luo H, Dong R. 2017. Treatment of anaerobic digested effluent in biochar-packed vertical flow constructed wetland columns: Role of media and tidal operation. Sci Total Environ. 592:197–205. https://doi.org/10.1016/j.scitotenv.2017.03.125.

Lahav O, Green M. 2000. Bioregenerated ion-exchange process: The effect of the biofilm on ion-exchange capacity and kinetics. Water SA. 26. https://www.wrc.org.za/wp-content/uploads/mdocs/WaterSA_2000_01_jan00_p51.pdf.

Lang JR, Allred BM, Field JA, Levis JW, Barlaz MA. 2017. National estimate of per- and polyfluoroalkyl substance (PFAS) release to U.S. municipal landfill leachate. Environ Sci Technol. 51(4):2197–2205. https://doi.org/10.1021/acs.est.6b05005.

Langergraber G. 2011. Numerical modelling: a tool for better constructed wetland design? Water Sci Technol. 64(1):14–21. https://doi.org/10.2166/wst.2011.520.

Langergraber G, Giraldi D, Mena J, Meyer D, Peña M, Toscano A, Brovelli A, Korkusuz EA. 2009. Recent developments in numerical modelling of subsurface flow constructed wetlands. Sci Total Environ. 407(13):3931–3943. https://doi.org/10.1016/j.scitotenv.2008.07.057.

Langergraber G, Šimůnek J. 2005. Modeling variably saturated water flow and multicomponent reactive transport in constructed wetlands. Vadose Zone J. 4(4):924–938. https://doi.org/10.2136/vzj2004.0166.

Lei G, Zeng W, Yu J, Huang J. 2023. A comparison of physical-based and machine learning modeling for soil salt dynamics in crop fields. Agr Water Manage. 277:108115. https://doi.org/10.1016/j.agwat.2022.108115.

Martin JF, Reddy KR. 1997. Interaction and spatial distribution of wetland nitrogen processes. Ecol Model. 105(1):1–21. https://doi.org/10.1016/S0304-3800(97)00122-1.

Mayo AW, Bigambo T. 2005. Nitrogen transformation in horizontal subsurface flow constructed wetlands I: Model development. Physics and Chemistry of the Earth, Parts A/B/C. 30(11):658–667. doi:10.1016/j.pce.2005.08.005.

Najafi Savadroudbari E, Fallah N, Davarpanah L, Nasernejad B. 2021. A comparison between anaerobic and aerobic biological treatment for real wastewater containing high concentration of dimethylamine: case study of wastewater from artificial leather production. Desalin Water Treat. 220:93–100. https://doi.org/10.5004/dwt.2021.26936.

Nga DT, Hiep NT, Hung NTQ. 2020. Kinetic modeling of organic and nitrogen removal from domestic wastewater in a down-flow hanging sponge bioreactor. Environ Eng Res. 25(2):243–250. https://doi.org/10.4491/eer.2018.390.

Pálfy TG, Gourdon R, Meyer D, Troesch S, Molle P. 2017. Model-based optimization of constructed wetlands treating combined sewer overflow. Ecol Eng. 101:261–267. https://doi.org/10.1016/j.ecoleng.2017.01.020.

Payne KA. 2018. Mathematical and numerical modeling of hybrid adsorption and biological treatment systems for enhanced nitrogen removal (dissertation). Tampa (FL): University of South Florida. 137 p. https://digitalcommons.usf.edu/etd/7702.

Ranatunga T, Tong STY, Yang YJ. 2017. An approach to measure parameter sensitivity in watershed hydrological modelling. Hydrolog Sci J. 62(1):76–92.
https://doi.org/10.1080/02626667.2016.1174335.

Rawat J, Saxena J, Sanwal P. 2019. Biochar: a sustainable approach for improving plant growth and soil properties. In: Abrol V, Sharma P, editors. Biochar - an imperative amendment for soil and the environment. London (UK): IntechOpen; [accessed 2023 Feb 12]. https://doi.org/10.5772/intechopen.82151.

Rodriguez-Gonzalez L, Miriyala A, Rice M, Delgado D, Marshall J, Henderson M, Ghebremichael K, Mihelcic J, Ergas S. 2020. A pilot-scale hybrid adsorption–biological treatment system for nitrogen removal in onsite wastewater treatment. J Sustain Water in Built Environ. 6:04019014. https://doi.org/10.1061/JSWBAY.0000898.

Sadeghi M, Shayannejad M, Sadraddini AA, Salehi S. 2024. Estimation of soil infiltration coefficient in the furrow irrigation using the combination of the hydrodynamics model and Richard’s equation. Phys Chem Earth, Parts A/B/C. 135:103649. https://doi.org/10.1016/j.pce.2024.103649.

Saltelli A. 2002. Sensitivity analysis for importance assessment. Risk Anal. 22(3):579–590. https://doi.org/10.1111/0272-4332.00040.

Samsó R, Meyer D, García J. 2014. Subsurface Flow Constructed Wetland Models: Review and Prospects. In: The Role of Natural and Constructed Wetlands in Nutrient Cycling and Retention on the Landscape. Springer. p. 149.

Schmidt I, Sliekers O, Schmid M, Cirpus I, Strous M, Bock E, Kuenen JG, Jetten MSM. 2002. Aerobic and anaerobic ammonia oxidizing bacteria – competitors or natural partners? FEMS Microbiol Ecol. 39(3):175–181.
https://doi.org/10.1111/j.1574-6941.2002.tb00920.x.

Šimůnek J, van Genuchten MTh, Šejna M. 2016. Recent Developments and Applications of the HYDRUS Computer Software Packages. Vadose Zone Journal. 15(7):vzj2016.04.0033. https://doi.org/10.2136/vzj2016.04.0033.

Ukoba K, Jen T-C. 2022. Biochar and application of machine learning: a review. In: Bartoli M, Giorcelli M, Tagliaferro A, editors. Biochar - productive technologies, properties and applications. London (UK): IntechOpen; [accessed 2023 Jan 24]. https://doi.org/10.5772/intechopen.108024.

Verhamme EM, Redder TM, Schlea DA, Grush J, Bratton JF, DePinto JV. 2016. Development of the Western Lake Erie Ecosystem Model (WLEEM): Application to connect phosphorus loads to cyanobacteria biomass. J Great Lakes Res. 42(6):1193–1205. https://doi.org/10.1016/j.jglr.2016.09.006.

Vymazal J. 2005. Horizontal sub-surface flow and hybrid constructed wetlands systems for wastewater treatment. Ecol Eng. 25(5):478–490. https://doi.org/10.1016/j.ecoleng.2005.07.010.

Vymazal J, Kröpfelová L. 2009. Removal of nitrogen in constructed wetlands with horizontal sub-surface flow: a review. Wetlands. 29(4):1114–1124. https://doi.org/10.1672/08-216.1.

Weather Underground. 2024. Lithia, FL Weather History. TWC Product and Technology LLC [accessed 2024 Jul 19]. https://www.wunderground.com/history/daily/us/fl/lithia/KFLLITHI4/date/2021-2-2.

Welander U, Henrysson T, Welander T. 1997. Nitrification of landfill leachate using suspended-carrier biofilm technology. Water Res. 31(9):2351–2355.
https://doi.org/10.1016/S0043-1354(97)00080-8.

Wong THF, Fletcher TD, Duncan HP, Jenkins GA. 2006. Modelling urban stormwater treatment—A unified approach. Ecol Eng. 27(1):58–70. https://doi.org/10.1016/j.ecoleng.2005.10.014.

Yang X, Arias ME, Ergas SJ. 2023. Hybrid constructed wetlands amended with zeolite/biochar for enhanced landfill leachate treatment. Ecol Eng. 192:106990. https://doi.org/10.1016/j.ecoleng.2023.106990.

Yuan C, Huang T, Zhao X, Zhao Y. 2020. Numerical models of subsurface flow constructed wetlands: review and future development. Sustainability. 12(8):3498. https://doi.org/10.3390/su12083498.

Zakir Hossain SM, Al-Bastaki N, Alnoaimi AMA, Ezuber H, Razzak SA, Hossain MM. 2019. Mathematical modeling of temperature effect on algal growth for biodiesel application. Renew Energ Environ Sustainability. 4:8. https://doi.org/10.1051/rees/2019005.

Zhao Y, Zhang L, Zhang M, Wu J, Li S, Ran D, Sun L, Zhu G. 2021. Effect of chemical oxygen demand concentration on nutrient removal in simultaneous nitrification, denitrification and phosphorus removal system in high-altitude areas. Water. 13(19):2656. https://doi.org/10.3390/w13192656.

Zulfikar Z, Nasrullah N, Kartini K, Aditama W. 2022. Effect of hydraulic retention time on the levels of biochemical oxygen demand and total suspended solid with simple integrated treatment as an alternative to meet the household needs for clean water. Open Access Maced J Med Sci. 10(E):6–11. https://doi.org/10.3889/oamjms.2022.7828.

Peer Review File

Received February 14, 2024
Accepted September 23, 2024
Published December 16, 2024

Editors
Sara W. McMillan, Editor in Chief
James Bays, Associate Editor

Banner image by Xia Yang.