Development of Field Pollutant Load Estimation Module and Linkage of QUAL2E with Watershed-Scale L-THIA ACN Model

Abstract

The Long Term Hydrologic Impact Assessment (L-THIA) model was previously improved by incorporating direct runoff lag time and baseflow. However, the improved model, called the L-THIA asymptotic curve number (ACN) model cannot simulate pollutant loads from a watershed or instream water quality. In this study, a module for calculating pollutant loads from fields and through stream networks was developed, and the L-THIA ACN model was combined with the QUAL2E model (The enhanced stream water quality model) to predict instream water quality at a watershed scale. The new model (L-THIA ACN-WQ) was applied to two watersheds within the Korean total maximum daily loads management system. To evaluate the model, simulated results of total nitrogen (TN) and total phosphorus (TP) were compared with observed water quality data collected at eight-day intervals. Between simulated and observed data for TN pollutant loads in Dalcheon A watershed, the R² and Nash–Sutcliffe efficiency (NSE) were 0.81 and 0.79, respectively, and those for TP were 0.79 and 0.78, respectively. In the Pyungchang A watershed, the R² and NSE were 0.66 and 0.64, respectively, for TN and both statistics were 0.66 for TP, indicating that model performed satisfactorily for both watersheds. Thus, the L-THIA ACN-WQ model can accurately simulate streamflow, instream pollutant loads, and water quality.

1. Introduction

Estimation of pollutant loads in watersheds is very important for long-term watershed management plans such as non-point source (NPS) pollution management or total maximum daily loads (TMDLs) [1]. To estimate pollutant loads at a watershed, long-term monitoring of flow and water quality is essential. However, due to limitations in monitoring flow and pollutants in terms of labor and budget [2], various computer models have been developed and used worldwide [3,4,5,6,7]. Diverse input data are needed to operate complex watershed-scale hydrology or water quality models such as SWAT (Soil and Water Assessment Tool) [8], HSPF (Hydrologic Simulation Program–Fortran) [9], and SWMM (Storm Water Management Model) [10,11]. For this reason, simpler computer models, such as PLOAD (Pollutant Loading Application) [12], STEPL (Spreadsheet Tool for Estimating Pollutant Load [13,14], and L-THIA (Long Term Hydrologic Impact Assessment) [15,16,17,18], have been used at the screening level. Among these models, L-THIA has been used in various studies and long-term watershed management plans because only daily rainfall data, in conjunction with a land use or soil dataset, are needed for daily direct runoff simulation and pollutant loads using the Natural Resources Conservation Services curve number (NRCS-CN) method and event mean concentration (EMC) databases for representative land uses [17,19,20,21]. The ArcView/ArcGIS/Web/Web GIS interface to the L-THIA model has been developed and used for almost two decades [15,16,18,20,22,23]. The Low Impact Development (LID) module was developed and integrated into L-THIA for environmentally friendly urban land use planning [24].

However, the current L-THIA model has limitations in estimating direct runoff during low-flow season because it does not compute direct runoff depending on the amount of rainfall [17,25,26]. In addition, the baseflow component for each hydrologic response unit (HRU) is not simulated [17]. When applying the L-THIA model to a watershed, lag time of the direct runoff component (to account for time needed to reach streams), the baseflow component contribution from unconfined aquifers, and flow routing through channel networks should all be considered. These limitations were addressed in the watershed-scale L-THIA asymptotic curve number (ACN) model developed by Ryu et al. [27].

With these capabilities enhanced, the watershed-scale L-THIA ACN model could be used at ungauged watersheds for streamflow estimation. However, it does not simulate pollutant loads from HRUs, or evaluate the effects of watershed routing on water quality changes through stream networks.

Thus, the objectives of this study are to (1) develop a module for estimation of pollutant loads from fields based on EMC databases for various land uses and amounts of rainfall; (2) to add the QUAL2E model (The enhanced stream water quality model) for stream water quality modeling to the watershed-based L-THIA model; and (3) to evaluate the accuracy of the watershed-scale L-THIA ACN model with water quality (watershed-scale L-THIA ACN-WQ model) by comparing estimated water quality with observed water quality.

2. Materials and Methods

2.1. Development of the Water Quality Simulation Module of the Watershed-Scale L-THIA ACN Model

The watershed-scale L-THIA ACN model, developed by Ryu et al. [27], does not include a water quality component. Therefore, a field pollutant load estimation module and a stream water quality module were added to the L-THIA model (Figure 1).

In the L-THIA model, pollutant loads for eight representative land uses were computed by multiplying EMC data collected in the USA by direct runoff volume [17,28].

2.1.1. Development of the Field Pollutant Load Estimation Module

To develop the field pollutant load estimation module, it is assumed that pollutant loads through direct runoff and baseflow should be estimated from each HRU, and move to the stream of each subbasin. For this, the EMC data for direct runoff and baseflow components for various land cover types are needed. However, EMC data on infiltration from various land cover types are not available in South Korea, as is the case in many other countries. Thus, pollutant load from the baseflow component was estimated by multiplying the baseflow by adjusted EMC data. The adjusted EMC data is calculated by multiplying the EMC concentration by a user-defined coefficient to explain adsorption or degradation through soil and the loss of some portions to deep aquifers [29]. The EMC data were derived from monitoring data collected at various locations; thus, regional characteristics of EMC should be considered when these EMC data are applied. The EMC adjustment coefficients were used in pollutant estimation from each HRU in this study (Equation (1)), which allows ease of model use, providing reasonable accuracy using a limited number of input parameters, similar to the original L-THIA model [15,16,17,18].

$$\begin{array}{l}EM{C}_{DR}=EMC\times Adj\_EM{C}_{DR}\\ EM{C}_{BF}=EMC\times Adj\_EM{C}_{BF}\end{array}$$

(1)

where EMC_DR is the pollutant concentration in direct runoff (mg/L), EMC_BF is the pollutant concentration in baseflow (mg/L), EMC is the initial pollutant concentration of both direct runoff and baseflow (mg/L), Adj_EMC_DR is the adjustment coefficient for EMC in the direct runoff, and Adj_EMC_BF is the adjustment coefficient for EMC in the baseflow.

The pollutant loads through either direct runoff or baseflow from each HRU (Equation (2)) are summed together for each subbasin (Equation (3)) and then these loads, later converted to concentration of pollutants, are used as input data to the instream water quality simulation (Equation (4)):

$$\begin{array}{l}{L}_{DR,hru}=\left[\frac{EM{C}_{DR}\times Q_{DR,HRU}}{1000}\right]\\ L_{BF,hru}=\left[\frac{EM{C}_{BF}\times Q_{BF,HRU}}{1000}\right]\end{array}$$

(2)

where L_DR,hru is the pollutant load of direct runoff (kg), L_BF,hru is the pollutant load of baseflow (kg), EMC_DR is the pollutant concentration in direct runoff (mg/L), EMC_BF is the pollution concentration in baseflow (mg/L), Q_DR,HRU is the amount of direct runoff (m³) discharged to the main channel, and Q_BF,HRU the amount of baseflow (m³) into the main channel:

$$\begin{array}{l}{L}_{DR,sub}={\displaystyle \sum _{HRU}^n{L}_{DR,hru}}\\ L_{BF,sub}={\displaystyle \sum _{HRU}^n{L}_{BF,hru}}\end{array}$$

(3)

$$\begin{array}{l}{C}_{DR,sub}=\left[\frac{L_{DR,sub}\times 1000}{Q_{DR,sub}}\right]\\ C_{BF,sub}=\left[\frac{L_{BF,sub}\times 1000}{Q_{BF,sub}}\right]\end{array}$$

(4)

where L_DR,sub is the pollutant load of direct runoff in subbasin (kg), L_BF,sub is the pollutant load of baseflow in subbasin (kg), C_DR is the pollutant concentration of direct runoff in the subbasin (mg/L), C_BF is the pollutant concentration of baseflow in the subbasin (mg/L), Q_DR,sub is the amount of direct runoff discharged to the streams in the subbasin (m³), and Q_BF,sub is the amount of baseflow into the main channel in the subbasin (m³).

The EMC data from the “Project of the long-term monitoring for the nonpoint source (NPS) pollution”, monitored by the Ministry of Environment (MOE) in South Korea, were used (

Table 1. Monitored event mean concentration (mg/L) of individual land cover for various rainfall magnitudes.

Land Cover	BOD5 (Biocheminal Oxygen Demand)				TN (Total Nitrogen)				TP (Total Phosphorous)
	<10 mm	10–30	31–50	>50 mm	<10 mm	10–30	31–50	>50 mm	<10 mm	10–30	31–50	>50 mm
Residential area	6.48	6.66	4.66	3.08	5.44	5.71	3.17	4.06	0.373	0.289	0.210	0.282
Manufacturing area	24.96	22.36	10.72	4.97	3.11	3.49	7.09	2.94	0.342	0.423	0.338	0.358
Regional public facility area	32.21	49.52	22.70	14.52	7.14	8.62	6.08	2.47	0.823	0.757	0.425	0.287
Recreational facility area	18.53	10.76	10.21	6.38	6.28	3.69	4.89	1.67	0.845	0.352	0.460	0.269
Road	10.61	7.17	8.89	3.31	6.18	2.68	2.76	2.01	0.221	0.231	0.217	0.176
Commercial area	6.52	8.12	4.48	4.12	5.81	5.12	2.59	5.68	0.277	0.331	0.212	0.453
Upland	2.26	3.40	3.57	3.50	2.38	2.83	2.33	2.41	0.215	0.386	0.301	0.436
Orchard	0.00	0.68	5.11	11.37	0.00	1.93	6.75	6.89	0.000	0.692	1.358	2.233
Green house	4.77	7.13	6.34	5.78	1.88	3.04	5.18	3.12	0.585	1.155	2.948	2.595
Paddy	0.00	1.10	2.09	2.63	0.00	1.41	3.88	5.54	0.000	0.154	0.423	0.962
Pasture	1.18	1.16	1.22	1.64	1.33	2.26	2.48	2.81	0.043	0.059	0.046	0.061
Forest	3.80	5.01	3.46	3.23	1.69	2.69	2.30	3.33	0.193	0.325	0.232	0.372
Bare land	6.52	8.12	4.48	4.12	5.81	5.12	2.59	5.68	0.277	0.331	0.212	0.453

) because the EMC data for various land cover types and rainfall ranges are available in these datasets. The NPS pollution loads have been measured over seven years (2008–2015) in a small subbasin for thirteen representative land covers. The data compiled through this long-term monitoring project have information on NPS pollution, such as the amount and ratio of NPS pollution in runoff, and EMC data [30,31,32,33]. In recent years, the MOE of Korea has determined the EMC and unit loads for each land cover using the 7-year data from the project. The data were validated through government hearings [34].

2.1.2. Incorporation of Simplified QUAL2E Instream Water Quality Model

In this study, the simplified QUAL2E model [35], which has been widely used in one-dimensional water quality simulation, was added into the channel routing module of the watershed-scale L-THIA ACN model. The mechanism of instream water quality changes in the QUAL2E model is given in Brown and Barnwell [35].

To simulate water quality through stream networks using the simplified QUAL2E, flow and pollutant concentration parameters from each subbasin must be prepared. Pollutant concentrations calculated in Section 2.1.1 were classified as total nitrogen (TN) and total phosphorus (TP), without being further divided into specific categories such as nitrate nitrogen (NO₃-N), nitrite nitrogen (NO₂-N), ammonia nitrogen (NH₃-N), organic nitrogen (organic-N), organic phosphorus (organic-P), and inorganic-phosphorus (inorganic-P). However, for simulation of instream water quality changes, division into these categories is required.

Thus, pollutant concentrations of TN and TP in both direct runoff and baseflow were subdivided by partitioning coefficient parameters, which must be prepared by model users based on water quality data collected in the study watershed or nearby areas. After partitioning TN and TP, the six pollutant concentrations (NO₃-N, NO₂-N, NH₃-N, organic-N, organic-P and inorganic-P) of the subbasin, and of the direct runoff and baseflow, were used as input data for water quality simulation in streams, using the method in Brown and Barnwell [23]. For the baseflow component, organic-N and organic-P were excluded because these nutrients have greater adsorption characteristics to soil particles while moving downwards from the surface. Only soluble pollutants can infiltrate into aquifers from the land surface [29]. The QUAL2E model can only simulate carbonaceous biochemical oxygen demand (CBOD), and not bottle BOD₅ (five-day biochemical oxygen demand); therefore, bottle BOD₅ of the subbasin calculated by Equation (4) was converted to CBOD following Equation (5):

$$CBO{D}_u=\frac{bottleBO{D}_5-4.57(N{H}_3)(1-e^{-5\times knb})-{\mathsf{\alpha }}_{4}A(1-e^{-5\times \mathsf{\rho }})}{1-e^{-5\times kdb}}$$

(5)

where CBOD_u is the ultimate CBOD concentration (mg/L), hereinafter referred as CBOD, and bottleBOD₅ is the bottle BOD₅ concentration calculated by Section 2.1.1 (mg/L), NH₃ is the ammonia concentration (mg/L), A is the algal biomass (mg/L), $\rho$ is the respiration rate of algae (day⁻¹), knb is the nitrification rate coefficient in bottle BOD₅ at 20 °C (day⁻¹), and kdb is the deoxidation rate coefficient in bottle BOD₅ at 20 °C (day⁻¹).

Since initial dissolved oxygen (DO) and algae pollutant concentrations required in QUAL2E are not available in the EMC databases (Table 1), the initial concentration of DO was estimated using the method proposed by APHA (Washington, DC, USA) [36]. Changes in DO were calculated using the method in Brown and Barnwell [35]. The amount of algal biomass is closely related to chlorophyll a and can be derived by using its relationship with the concentration of chlorophyll a [37]. Chlorophyll a concentration was derived by using a simplified version of the exponential function proposed by Cluis et al. [37] (Equation (5)):

$$\begin{array}{l}Chl=0 , if Q_{DR,sub}<{10}^{-5}{m}^3/s, or (TP and TN < {10}^{-6})\\ Chl=\frac{0.5\times {10}^{2.7}}{Q_{DR,sub}} , if Q_{DR,sub}>{10}^{-5}{m}^3/s, and (TP and TN > {10}^{-6})\\ Chl=\frac{0.5\times {10}^{0.5}}{Q_{DR,sub}} , if Q_{DR,sub}>{10}^{-5}{m}^3/s, TP<{10}^{-6}and TN>{10}^{-6}\end{array}$$

(6)

where Chl is the chlorophyll a concentration in the direct runoff (μ/L), Q_DR,sub is direct runoff flow rate into the main channel (m³/s), TN is the total Kjeldahl nitrogen load (kmol), TP is the total phosphorus load (kmol).

The initial concentration of DO was estimated using Equation (7) proposed by APHA (Washington, DC, USA) [36]:

$$DO=\exp \left[-139.34410+\frac{1.575701\times {10}^5}{TempK}-\frac{6.642308\times {10}^7}{Temp{K}^2}+\frac{1.243800\times {10}^{10}}{Temp{K}^3}-\frac{8.621949\times {10}^{11}}{Temp{K}^4}\right]$$

(7)

where DO is the saturation concentration of dissolved oxygen (mg/L), and TempK is water temperature in Kelvin (273.15 + °C).

Water temperature was calculated using Equation (8) proposed by Stefan and Preud’ homme [38]:

$$TempK=273.15+(5.0+0.75T_{average})$$

(8)

where TempK is water temperature during the day (°C) and T_average is the average air temperature during the day (°C).

The initial concentrations of chlorophyll a and DO in the subbasin were used for input data of the instream water quality component of the watershed-scale L-THIA ACN-WQ model. The mechanism of instream water quality changes in the QUAL2E model is shown in Figure 2. More detailed description for each submodule of QUAL2E is explained below in details.

Weather data (e.g., solar radiation, duration of solar radiation and average air temperature) are required to simulate instream water quality changes using QUAL2E. A change in algae biomass in streams affects nutrient levels. The amount of algal biomass is simulated using solar radiation and duration as outlined by Brown and Barnwell [35]. The watershed-scale L-THIA ACN-WQ model was modified to use the weather data for instream water quality simulation.

In the watershed-scale L-THIA ACN-WQ model, the EMC data can be adjusted to consider local EMC characteristics. The subdivisions of TN and TP were added to the model to simulate the NO₃-N, NO₂-N, NH₃-N, organic-P and inorganic-P in streams (

Table 2. Description of event mean concentration (EMC), nitrogen, and phosphorus parameters used in instream water quality module.

Parameter Name	Description	Range	Default Value
Adj_EMCDR,N	Constant value for adjustment of nitrogen in surface	−0.9–0.9	1.0
Adj_EMCBF,N	Constant value for adjustment of nitrogen in aquifer	−0.9–0.9	1.0
Adj_EMCDR,P	Constant value for adjustment of phosphorus in surface	−0.9–0.9	1.0
Adj_EMCBF,P	Constant value for adjustment of phosphorus in aquifer	−0.9–0.9	1.0
TN_ratio1 1	Ratio of organic-N in total nitrogen	0.0–0.9	0.05
TN_ratio2 1	Ratio of NO3-N in total nitrogen	0.0–0.9	0.8
TN_ratio3 1	Ratio of NH3-N in total nitrogen	0.0–0.9	0.1
TP_ratio1 2	Ratio of organic-P in total phosphorus	0.0–0.9	0.5

Notes: ¹ Ratio of NO₂-N = 1 − (Sum of TN_ratio1, TN_ratio2 and TN_ratio3); ² Ratio of dissolved P = 1 − TP_ratio1.

).

Pollutant loads in the stream were simulated using the QUAL2E model and the parameters in QUAL2E were used in the watershed-scale L-THIA ACN-WQ model (

Table 3. Description of QUAL2E parameters used in instream water quality module.

Parameter Name	Description	Recommended Range in QUAL2E	Default Value
RS1	Local algal settling rate in the reach at 20 °C	0.15–1.82	0.3408
RS2	Benthic source rate for dissolved phosphorus in the reach at 20 °C	0.001–0.1	0.1
RS3	Benthic source rate for NH4-N in the reach at 20 °C	0.0–1.0	0.0
RS4	Rate coefficient for organic nitrogen settling in the reach at 20 °C	0.001–0.1	0.001
RS5	Organic phosphorus settling rate in the reach at 20 °C	0.001–0.1	0.08
RK1	Carbonaceous biological oxygen demand (CBOD) deoxygenation rate coefficient in the reach at 20 °C	0.02–3.4	0.3
RK2	Oxygen reaeration rate in accordance with Fickian diffusion in the reach at 20 °C	0.0–100.0	1.0
RK3	Rate of loss of CBOD due to settling in the reach at 20 °C	−0.36–0.36	−0.36
RK4	Benthic oxygen demand rate in the reach at 20 °C	0.0–100.0	0.0
BC1	Rate constant for biological oxidation of NH4 to NO2 in the reach at 20 °C	0.1–1	0.1
BC2	Rate constant for biological oxidation of NO2 to NO3 in the reach at 20 °C	0.2–2	0.2
BC3	Rate constant for hydrolysis of organic N to NH4 in the reach at 20 °C	0.2–0.4	0.03
BC4	Rate constant for mineralization of organic P to dissolved P in the reach at 20 °C	0.01–0.7	0.1
RTH	Algal respiration rate at 20 °C	0.05–5.0	0.05
TFAC	Fraction of photosynthetically active solar radiation	0.0–1.0	0
MMX	Maximum specific algal growth rate at 20 °C	1.0–3.0	1.0
IG	QUAL2E algae growth limiting option (1: multiplicative, 2: limiting nutrient, 3: harmonic mean)	1, 2, 3	1
A0	Ratio of chlorophyll-a to algal biomass	10.0–100.0	10
A1	Fraction of nitrogen algal biomass	0.07–0.09	0.071
A2	Fraction of phosphorus algal biomass	0.01–0.02	0.003
A3	Rate of oxygen production per unit of algal photosynthesis	1.4–2.3	1.4
A4	Rate of oxygen uptake per unit of algal respiration	1.6–2.3	1.6
A5	Rate of oxygen uptake per unit NH3-N oxidation	3.0–4.0	3.0
A6	Rate of oxygen uptake per unit NO2-N oxidation	1.0–1.14	1.0
Lam0	Non-algal portion of the light extinction coefficient	0–10	0
Lam1	Linear algal self-shading coefficient	0.006–0.065	0.006
Lam2	Nonlinear algal self-shading coefficient	0–1	0
KN	Michaelis-Menten nitrogen half-saturation constant	0.01–0.3	0.01
KP	Michaelis-Menten phosphorus half-saturation constant	0.001–0.05	0.001
KL	Light half-saturation coefficient	0.223–1.135	0.223
Knb	Nitrification rate coefficient in CBOD5	–	0.5
Kdb	Deoxidation rate coefficient in CBOD5	–	0.5
PN	Preference factor for ammonium nitrogen	0.0–1.0	0.0

). The default values of QUAL2E parameters were set in the model, as shown in Table 3, with reference to studies of Glavan et al. and Na [39,40].

2.2. Applications of the Watershed-Scale L-THIA ACN-WQ Model

To evaluate the watershed-scale L-THIA ACN-WQ model, estimated streamflow and pollutant loads were compared with observed data at two study watersheds (Dalcheon A and Pyungchang A in South Korea where the TMDL was carried out; Figure 3).

Dalcheon A is located in the Han River basin in South Korea. The average slope of this watershed (22.0%) is lower than that of Pyungchang A watershed (30.8%). Land cover types in Dalcheon A consist of forest (71.9%), agricultural (18.6%), urban (2.5%), bare land (0.9%) and others (6.1%).

Pyungchang A is also located in the Han-River basin. The land cover types consist of forest (80.5%), agricultural (12.8%), urban (1.7%), bare land (2.0%), and others (1.2%).

Information on the study watersheds including surface area, soil, climate, agriculture production type and water quality information are shown in

Table 4. Information of the study area.

Information of Study Area		Dalcheon A	Pyungchang A
Surface area (km2)		1200.33	1756.87
Average precipitation (2011–2014) (mm/year)		1315.6	1389.6
Average temperature (2011–2014) (degree)		11.7	12.8
Hydrologic soil group (%)	A	11.8	11.7
	B	16.0	15.9
	C	36.6	35.8
	D	35.6	36.6
Agriculture production type (%)	Upland crop	42.5	32.6
	Green house	14.7	20.7
	Ochard	21.6	22.6
	Paddy	21.2	24.1
Average water quality (2011–2014) (mg/L)	BOD	0.96	0.88
	TN	2.76	2.91
	TP	0.02	0.03

.

The watershed-scale L-THIA ACN-WQ model, which can estimate instream water quality, was applied to these study watersheds for the period 2010–2014. The hydrology model parameters and water quality-related parameters were calibrated using 8-day interval observed flow and water quality data (the only flow and water quality data available for these TMDL watersheds).

In addition to the adjusted EMC, various QUAL2E parameters such as algal respiration rate, algal settling rate, algal growth rate, nitrification rate of nitrite, and phosphorus settling rate were adjusted for a better fit with measured TN and TP data, when calibrating instream water quality of TN and TP.

3. Results and Discussion

3.1. Estimation of Pollutant Load

For assessing the model performance, parameters of the model were manually calibrated by comparing model simulations with observations for streamflow, TN, and TP. Pearson coefficient of correlation (R²) and Nash–Sutcliffe simulation efficiency index (NSE) [41] were used as the objective function, all of which are well known and generally used in evaluating the performance of models [39,42].

The calibrated CN values in Dalcheon A and Pyungchang A are 3% lower and 9% higher, respectively, than the default CN values (

Table 5. Calibration parameters for streamflow component of the L-THIA ACN-WQ model.

Watershed	Adj_CN 1	SLSUB 2	DRlag 3	αBF 4	aqfthr 5	Frconf 6	BFdelay 7	Mk1 8	Mk2 9	Mkx 10
Dalcheon A	0.03	1.5	4	0.7	20.0	0.10	1	0.05	0.95	0.2
Pyungchang A	0.09	1.0	8	0.5	30.0	0.05	5	0.25	0.75	0.5

Notes: ¹ Adjusted coefficient for CN; ² Adjustment for slope length; ³ Direct runoff lag coefficient; 4 Baseflow recession constant; 5 Threshold water level for baseflow contribution in the unconfined aquifer (mm); 6 Fraction of water flowing into confined aquifer; 7 Delay time in aquifer recharge following water infiltration at surface (day); 8 Weighting factor for influence of normal flow on storage time constant value; 9 Weighting factor for influence of low flow on storage time constant value; 10 Weighting factor for Muskingum method.

). The calibrated value for field slope length at Dalcheon A watershed is 1.5 times greater than that estimated from the relationship between average slope and slope length suggested by Wischmeier and Smith [43].

Threshold water levels (aqf_thr) in Dalcheon A and Pyungchang A were 20 mm and 30 mm, respectively.

Calibrated parameters of the watershed-scale L-THIA ACN-WQ model indicate that EMC adjustment coefficients affect instream water quality as well as pollutant loads from each subbasin (

Table 6. Calibration parameters for the water quality component of the L-THIA ACN-WQ model.

Parameters	Dalcheon A	Pyungchang A
Adj_EMCDR,N	0.80	0.70
Adj_EMCBF,N	0.90	0.90
Adj_EMCDR,P	0.60	0.40
Adj_EMCBF,P	0.40	0.35
TN_ratio1	0.03	0.03
TN_ratio2	0.95	0.95
TN_ratio3	0.01	0.01
TP_ratio1	0.40	0.40
RS1	1.00	1.00
RS2	0.001	0.001
RS3	0.001	0.01
RS4	0.10	0.001
RS5	0.10	0.1
RS6	2.50	2.50
RK1	0.50	0.50
RK2	50	50
RK3	0.36	0.36
RK4	2.00	2.00
RK5	2.00	2.00
RK6	1.71	1.71
BC1	0.10	0.55
BC2	0.20	2.00
BC3	0.02	0.4
BC4	0.01	0.01
TFAC	0.30	0.30
MMX	1	1
IG	3	3
A0	80	80
A1	0.09	0.09
A2	0.01	0.01
A3	1.60	1.60
A4	2	2
A5	3.50	3.50
A6	1.00	1.00
Lam0	1	1
Lam1	0.03	0.03
Lam2	0.054	0.054
KN	0.75	0.02
KP	0.020	0.025
KL	0.025	0.75
Knb	0.50	0.03
Kdb	0.03	0.045
PN	0.45	0.50

). This indicates that the EMC should be adjusted to account for differences in agricultural management (fertilizer/manure application rates, tillage, etc.) or urban impervious areas and population density at each watershed.

As shown in the calibration result, EMCs can be different for direct runoff and baseflow, and they should be separately calibrated.

In many other countries, including South Korea, it is difficult to predict the EMCs for baseflow by individual land cover type.

Therefore, if the L-THIA ACN-WQ model is calibrated with the EMC for baseflow in areas with various land cover types, the representative EMCs for baseflow by individual land cover type could be estimated. Furthermore, the representative EMCs for baseflow could be used for estimating a representative ratio of baseflow pollution loads to streamflow pollution loads that can contribute to estimating pollutant loads for baseflow in ungauged watersheds.

The parameters of the partitioning coefficients of organic N, NO₃-N, and NH₃-N in TN and organic-P in TP were defined based on the ratio of measured water quality data from outlets in the study watersheds (Table 6). Organic-N, NO₃-N, NH₃-N, and organic P (thus inorganic P) ratio values are similar in both watersheds, while NO₂-N concentrations are different. Thus, these partitioning values (organic-N, NO₃-N, NH₃-N, and organic P) from observed data were used as model inputs for both study watersheds, while these could be adjusted or calibrated based on local flow and water quality conditions. If no observed water quality data are available, these partitioning model parameters could be calibrated to reflect local conditions indirectly. Parameters related to water quality in streams (e.g., respiration rate of algae, growth rate of algae, settling rate for phosphorus, nitrification rate) were calibrated manually to better reflect stream characteristics.

3.2. Model Performance Compared to Observed Streamflow And Pollutant Loads

Estimated average streamflow using the watershed-scale L-THIA ACN-WQ model and observed streamflow at Dalcheon A watershed were 27.68 m³/s and 26.29 m³/s, respectively. The R² and NSE values were 0.79 and 0.78, respectively (

Table 7. Comparison between estimated and observed eight-day average streamflow in the period between 2011 and 2014.

Watershed	Average Streamflow (2011–2014) (m3/s)		R2	NSE
	Observation	Estimation
Dalcheon A	26.29	27.68	0.79	0.78
Pyungchang A	53.69	56.48	0.76	0.76

, Figure 4). Estimated average and observed streamflow at Pyungchang A watershed were 56.48 m³/s and 53.69 m³/s, respectively (Table 7, Figure 5).

Average estimated and observed TN pollutant loads in Dalcheon A were 5227.04 kg and 6077.44 kg, respectively, (

Table 8. Comparison between estimated and observed eight-day average pollutant loads in the period between 2011 and 2014.

Watershed	Pollutant Loads	Average Pollutant Load (2011–2014) (kg)		R2	NSE	Average Concentration (2011–2014) (mg/L)
		Observation	Estimation			Observation	Estimation
Dalcheon A	TN	6077.44	5227.04	0.81	0.79	2.76	2.70
	TP	242.13	195.75	0.79	0.78	0.02	0.03
Pyungchang A	TN	134,57.67	10,282.49	0.66	0.64	2.91	2.35
	TP	302.51	306.02	0.66	0.66	0.03	0.04

) and R² and NSE values were 0.74 and 0.72, respectively (Figure 6). Estimated and observed pollutant loads in Pyungchang A were 10,282.49 kg and 13,457.67 kg, respectively (Table 8). The R² and NSE values were 0.66 and 0.64, respectively (Figure 7).

Average estimated and observed TP pollutant loads in Dalcheon A were 195.75 kg and 242.13 kg, respectively. R² and NSE values were 0.70 and 0.64. Average pollutant loads estimated and observed in Pyungchang A were 306.02 kg and 302.51 kg. The R² and NSE values were both 0.66 (Table 8, Figure 8 and Figure 9).

The simulated and observed data for concentrations of TN and TP match reasonably well. For Dalcheon A, estimated and observed TN concentration values were 2.70 mg/L and 2.76 mg/L, respectively (Table 8). TP estimated and observed values were 0.03 mg/L and 0.02 mg/L, respectively. For Pyungchang A watershed, average simulated and observed values were 2.35 and 2.91 mg/L of TN, and 0.04 and 0.03 mg/L of TP, respectively.

These results show that the watershed-scale L-THIA ACN-WQ model can simulate not only streamflow but also pollutant loads at a watershed with less input data required than for complex watershed-scale models.

4. Conclusions

In this study, the watershed-scale L-THIA ACN model, developed by Ryu et al. [27], was enhanced by developing a field pollutant load estimation module and incorporating the instream QUAL2E model, in order to simulate pollutant loads at a watershed scale. Comparison of the simulated results of this L-THIA ACN-WQ model with observed values from two watersheds in South Korea reveal that the model can simulate not only streamflow but also pollutant loads reasonably well [44,45].

Eight parameters used for calculating pollutant loads in each land cover and 33 parameters for water quality changes in streams were added in the model. The model requires 51 parameters (ten for streamflow, eight for pollutant loads from the watershed, and 33 for instream water quality—default values for these 33 parameters, recommended by the US EPA, are available) in order to simulate streamflow and water quality.

Based on the calibrated parameters, the default values of the 33 instream water quality parameters could be used at ungagged watersheds in conjunction with ten (for streamflow) and eight (for watershed pollutant loads) model parameters prepared according to the localized conditions of the watershed in question. Further evaluations of the watershed-scale L-THIA ACN-WQ model under different weather, land use, and soil conditions are needed.

This study attempted to improve the shortcomings of the existing L-THIA model, which could simulate only pollutant loads for direct runoff. The model suggested in this study can estimate the total instream pollutant loads and water quality considering baseflow as well as direct runoff.

Furthermore, existing watershed models such as SWAT or HSPF involve complex hydrologic components and many parameters related to the watershed characteristics [46]. Consequently, users need to spend lots of time and effort to calibrate various model parameters.

On the other hand, the L-THIA ACN-WQ model can successfully estimate the total streamflow (direct runoff and baseflow) and instream pollutant loads with less parameters than those of SWAT or HSPF. In addition, the operating time of models for long-term simulations is significantly reduced when compared to the other catchment models.

In this regard, we believe that the L-THIA ACN-WQ model can be a useful tool for simulating streamflow and water quality for sustainable watershed management plans.

Nevertheless, there are still limitations in the developed model in that sensitivity and uncertainty analysis should be manually conducted. Moreover, the model has no function for simulating other water qualities (e.g., suspended soil) excluding the BOD, TN, and TP.

Thus, in the near future, an auto-calibration tool associated with various optimal-algorithms will be mounted on the model to conduct sensitivity and uncertainty analysis. Furthermore, new simple modules to simulate various water qualities such as SS and TOC will be added to the model.