Research on Water Quality Simulation and Water Environmental Capacity in Lushui River Based on WASP Model

: In recent years, the severe deterioration of water quality and eutrophication in the Yangtze River has brought much trouble to people’s lives. Because of this, numerous management departments have paid more and more attention to the treatment of the water environment. In order to respond to water environmental protection policy and provide management departments with a basis for reﬁning water quality, this paper takes the Zhuzhou section of Yangtze River-Lushui watershed as its research object. First, we used the Water Quality Analysis Simulation Program (WASP) model as a tool, and obtained the pollution load using the FLUX method formula. During the calibration process, the sensitivity analysis method, the orthogonal design method, and the trial and error method were used. Then, we veriﬁed the results by using water quality monitoring data published by Zhuzhou Ecological Environment Bureau. Following that, the water environmental capacity of the Lushui River in normal, wet and dry periods was calculated using the WASP model: the chemical oxygen demand (COD) was 14,072.94 tons/yr, 17,147.7 tons/yr and 10,998.18 tons/yr, respectively; ammonia nitrogen (AN) was 469.098 tons/yr, 571.59 tons/yr and 366.606 tons/yr, respectively; and total phosphorus (TP) was 93.8196 tons/yr, 114.318 tons/yr and 73.3212 tons/yr, respectively. The results show that the WASP model is applicable and reliable and can be used as an effective tool for water quality prediction and management in this area.


Introduction
During the past half-century, different models have been used to simulate water quality; some researchers show that models have become a primary management tool [1]. Therefore, in establishing water quality assessment scenarios, models have become an essential method for predicting water quality in a river [2]. The research progress and development trend of water quality models started relatively late, and the gap with foreign research mainly exists in the versatility and comprehensiveness of model development and the data required for modeling, but certain results have also been achieved.
In order to progressively protect the water quality environment, it is important to study the water environmental capacity, which indicates the amount of pollutants that can be contained in a body of water, according to a mathematical model that can obtain the flow volume, flow direction and water quality transport. Extensively, some models such as the MIKE  the Stream Water Quality Model (QUAL2K), as well as other models, can be used for the water simulation [3].
Ellina et al. applied the research of fuzzy implications via fuzzy linear regression in data analysis for a fuzzy model in 2020 [4]. Hu Kaiming et al. applied the WASP model to the simulation of sediment release items and found that endogenous phosphorus pollution had become a factor that could not be ignored in water quality changes [5]. Zhang Lingxia also used the fuzzy comprehensive evaluation model to evaluate the nutritional status of Boyang Lake in 2012, and used the WASP model to simulate and predict its water quality [6].
In addition, other research had been found in the United States and other countries to determine the impact of water assessment measures [3]. They have shown that river models play an essential role in the river basin because they can simulate the water quality in the river at different levels and provide results that could be useful for water assessment [4].
The purpose of this work was to use the WASP model to predict water quality in the Lushui River basin, using data from a monitoring point; to construct the Lushui River water quality model to simulate the concentration of COD, AN and TP; and to estimate the water environment capacity according to the water quality target required.
The research results aim to assist water quality management in addressing the uncertainty related to an ongoing effort by a citizens' group in monitoring water quality in the Lushui River basin, according to the cross-section water quality report and monitoring data provided by the Ecological Environmental Protection Technology and the Plan of Yangtze River Program in Zhuzhou City.
We selected the topic of researching water quality simulation and water environmental capacity in the Lushui River based on the WASP model to calibrate the water quality parameters of the Lushui River and then predict the water environmental capacity following the water environmental management requirements.

Case Study
The Lushui River Basin is located in the Hunan and Jiangxi provinces, China. Its area is about 5675 km 2 , with 2278 km 2 located in Jiangxi Province and 3397 km 2 in Hunan Province. There are two sources of the Lushui River: south and north, with the south as the primary source. Lushui originates from the south of Qianla Mountain in Pingxiang City, Jiangxi Province, and enters Hunan Province from Jinyushi. It joins the northern source at Shuangjiangkou, and Xiangjiang River at Lukou Town (Lukou District). In Hunan, the water flows through Liuyang City and Zhuzhou City. In this study, Lushui River Basin in Zhuzhou City (Hunan Province) is selected as the research area (shown in Figure 1). The mainstream of Lushui River is 82 km long in Zhuzhou City, including 68 km in Liling City and 14 km in Lukou District. The Lushui River system is developed with many tributaries, which flow through twenty-three towns in Liling City, Lukou District, and Youxian County in Zhuzhou, with a total drainage area of 2788 km 2 . It also includes a total of six monitoring sections on the mainstream of Lushui River from upstream to downstream.
Lushui River Basin belongs to the humid climate of a subtropical monsoon, characterized by average annual precipitation of 1300-1600 mm, abundant rainfall, distinct seasons and various heat conditions. Northwest winds prevail in winter, yielding dry, cold weather; in summer, southerly winds blow. Consequently, the weather is hot, rainy, and prone to waterlogging and drought. In addition, the rain is irregular from September to October. From June to August, the Lushui River is in wet season; January, November, and December are dry seasons; other periods are normal water seasons [7]. Lushui River Basin belongs to the humid climate of a subtropical monsoon, characterized by average annual precipitation of 1300-1600 mm, abundant rainfall, distinct seasons and various heat conditions. Northwest winds prevail in winter, yielding dry, cold weather; in summer, southerly winds blow. Consequently, the weather is hot, rainy, and prone to waterlogging and drought. In addition, the rain is irregular from September to October. From June to August, the Lushui River is in wet season; January, November, and December are dry seasons; other periods are normal water seasons [7].
Lushui is the main river of Zhuzhou City and the mother river of Liling City. Lushui Basin has fertile land, abundant products, dense population and rapid economic development. Lushui is a local mother river established in Hunan Province, a comprehensive watershed management model in the demonstration area of regional cooperation along the border of Hunan and Jiangxi. Under the trend that the whole country pays increasing attention to the environment, it is a relatively complete small watershed with typical significance. In order to ensure the steady improvement of the water quality and the timely completion of the construction of the inter-provincial model river in Lushui, the WASP water quality model was used to simulate the water quality and water environmental capacity in Lushui Basin, in order to present the management department with a basis for refined water quality. Figure 1 shows the overview of the research area.

Determination of the COD, AN and TP
Chemical oxygen demand (COD) is an indirect measurement of the amount of organic matter in a sample. With this test, it is possible to measure virtually all organic compounds that can be digested by a digestion reagent. Chemical oxygen demand (COD) is important for determining the amount of waste in the water. Waste that is high in organic matter requires treatment to reduce the amount of organic waste before discharging into receiving waters. COD measures organic matter by using a chemical oxidant. It is critical that a strong enough oxidant is used to react with virtually all organic material in the sample [8].
Ammonia nitrogen is one of the major nitrogen forms in the nitrogen cycle, especially in natural waters [8]. Ammonia nitrogen consists of ammonia (NH 3 ) and ammonium (NH 4+ ) in natural waters. The common methods for quantitative analysis of ammonia nitrogen in natural water are the indophenol blue (IPB) spectrophotometric methods and ophthalaldehyde (OPA) fluorometric methods [9]. The IPB spectrophotometric methods are based on the Berthelot reaction [10].
Phosphorus (P) is an essential nutrient element that is used by all living organisms for growth and energy transport [11] and is often the limiting nutrient for primary pro- Lushui is the main river of Zhuzhou City and the mother river of Liling City. Lushui Basin has fertile land, abundant products, dense population and rapid economic development. Lushui is a local mother river established in Hunan Province, a comprehensive watershed management model in the demonstration area of regional cooperation along the border of Hunan and Jiangxi. Under the trend that the whole country pays increasing attention to the environment, it is a relatively complete small watershed with typical significance. In order to ensure the steady improvement of the water quality and the timely completion of the construction of the inter-provincial model river in Lushui, the WASP water quality model was used to simulate the water quality and water environmental capacity in Lushui Basin, in order to present the management department with a basis for refined water quality. Figure 1 shows the overview of the research area.

Determination of the COD, AN and TP
Chemical oxygen demand (COD) is an indirect measurement of the amount of organic matter in a sample. With this test, it is possible to measure virtually all organic compounds that can be digested by a digestion reagent. Chemical oxygen demand (COD) is important for determining the amount of waste in the water. Waste that is high in organic matter requires treatment to reduce the amount of organic waste before discharging into receiving waters. COD measures organic matter by using a chemical oxidant. It is critical that a strong enough oxidant is used to react with virtually all organic material in the sample [8].
Ammonia nitrogen is one of the major nitrogen forms in the nitrogen cycle, especially in natural waters [8]. Ammonia nitrogen consists of ammonia (NH 3 ) and ammonium (NH 4+ ) in natural waters. The common methods for quantitative analysis of ammonia nitrogen in natural water are the indophenol blue (IPB) spectrophotometric methods and ophthalaldehyde (OPA) fluorometric methods [9]. The IPB spectrophotometric methods are based on the Berthelot reaction [10]. Phosphorus (P) is an essential nutrient element that is used by all living organisms for growth and energy transport [11] and is often the limiting nutrient for primary production in terrestrial and aquatic ecosystems [4]. For the determination of total phosphorus (TP; unfiltered samples) and total dissolved phosphorus (TDP; filtered samples) it is necessary to convert all phosphorus-containing species into a detectable form. For natural waters, the most common method of detection is the molybdenum blue chemistry with spectrophotometric detection, which determines molybdate-reactive ortho-phosphate using either a batch or flow-based approach [11].

Model Description and Application
WASP model has various modules such as organic toxicants, simple toxicants, eutrophication, mercury, and non-ionic organic toxicants. For this research, we used the Water 2021, 13, 2819 4 of 20 EUTRO module of WASP; it is more applicable than the traditional one. In fact, according to its ability, it is able to combine several interacting systems comprising of biochemical oxygen demand (BOD), organic nitrogen, DO, ammonia, nitrates, bacteria, phytoplankton, phosphates, pH and solids [6]. The EUTRO module simulates the water quality in terms of nutrients, bacterial decay, DO, and reactive pollutants in the water.
Based on the input data, the specific chemical kinetics equations and general mass balance are used for the simulation [5]. Equation (1) shows the mass balance around an infinitesimally small fluid volume: where C is the concentration (mg/L); t is the time in days; u x , u y and u z are the longitudinal, lateral, and vertical adjective velocities (m/day), respectively; E x , E y , and E z are the longitudinal, lateral, and vertical diffusion coefficients (m 2 /day), respectively; S L is the direct and diffuse loading rate (g/m 3 day); S B is the pollution load of the boundary rate (including upstream, downstream, benthic and atmospheric) (g/m 3 day); and S K is the total kinetic transformation rate (g/m 3 day). The basic concept of writing a mass balance equation for a body of water is to account for all the material entering and leaving the water body via direct addition of material (runoff and loads), via physical, chemical and biologic transformations, and via dispersive and advective transport mechanisms [4].

Model Input
The WASP model is empirical, so its establishment needs input data such as hydrology, water quality, topography, etc. This step also poses as one of the difficulties in the establishment of the model. This part gives a detailed introduction to the collection, processing, and input processes of the primary data, such as basic model information, segmentation, control unit body data, initial and boundary conditions, and pollution load.

Basic Information
The basic model information (data set) mainly included the start and end time of the simulation, the sub-module, the hydrodynamics module, etc. Euler's equations (Euler) are used to solve this research. Moreover, to determine the accuracy of the simulation, the maximum time step (max time step) in this study was set to 1. As this research mainly studied the 1D water quality model of the mainstream of the Lushui River, the 1-D Network Kinematic wave module was selected. The meteorological data included the environmental parameters, such as solar radiation, rainfall, wind speed, segment, air temperature, etc. Furthermore, the initial concentrations of the total phosphorus (TP), chemical oxygen demand (COD) and ammonia nitrogen (AN) in the surface water were taken as their value according to the concentration values monitored on 1 January 2019. The fraction of pollutants dissolved (fraction dissolved) used the model default value of 1.

Segmentation and Division of Control Units in the Study Area
According to the above basic principles and requirements, the research steps for the Lushui Basin control unit were as follows: 1.
Collection of primary information data of natural geographical features in the Lushui area, including digital elevation (DEM), topographic maps, tributaries, small watershed boundaries, county and city administrative boundaries, environmental function areas, urban built-up areas, companies, drinking water sources, etc.; 2.
Processing of geographic information data based on ArcGIS, and the above primary geographic information data were overlaid and analyzed. Afterward, a preliminary sketch of the control unit was obtained; 3.
Division of the control units: after understanding and processing the collected information, the preliminary control unit sketch was fine-tuned according to the control unit's principles of watershed division and clean boundary. Extraction and confirmation were performed upon the boundaries of different water catchment units and control units on different environmental functions in the same administrative region, to achieve a one-to-one correspondence between pollution source water quality and strong operability; 4.
Naming of the control units: in the absence of nationally unified naming standards and norms, naming was undertaken as a process of trying to understand and conform to the actual situation of each control unit. Therefore, the Lushui River was discretized into six control unit bodies. Each control unit body was given detailed geometric attributes. Some data were obtained from monitoring points during segmentation and generalization. In addition, other geometry parameters were obtained by the formula on WASP 8 User's Guide. The measured river geometries and water velocities data were used to obtain the value of the exponent and coefficient of velocities and the depth of each segment, as shown in Equations (2), (3), (4), (5) and (6) [12]: where  Figure 2 and Table 1.

Pollution Load
The pollutant loads could be calculated through the water quality, and the collected flow data could be calculated, since the monitoring point's data represented the information from all watershed sources upstream of the monitoring point. Various approaches, methods, and formulas have been developed to calculate pollutant loads using monitoring data. For this research, based on the FLUX method, the pollutant load was calculated by multiplying the concentration with the flow of each indicator in the river. FLUX method was developed by William Walker (1996) for the United States Army Corps of

Pollution Load
The pollutant loads could be calculated through the water quality, and the collected flow data could be calculated, since the monitoring point's data represented the information from all watershed sources upstream of the monitoring point. Various approaches, methods, and formulas have been developed to calculate pollutant loads using monitoring data. For this research, based on the FLUX method, the pollutant load was calculated by multiplying the concentration with the flow of each indicator in the river. FLUX method was developed by William Walker (1996) for the United States Army Corps of Engineers [13]. It can estimate loads of nutrients, sediment, and other water quality constituents using flow and concentration data. This technique was developed to model eutrophication in reservoirs and is now also commonly used for lake modeling [13]. FLUX method can estimate the monitoring period and loads using daily flow values and periodic water quality concentration values based on the relationship between concentration and flow, either as a whole, or by dividing the data into groups by some stratification. The water quality data was obtained from the Zhuzhou City Ecological Environment Bureau. The pollution loads in Tables 2 and 3 are the total pollution load value for all sectors such as point and non-point source pollution load, namely, agriculture, industry, domestic sector, and livestock, etc. To obtain the pollutant load value of every sector, it is essential to know the Water 2021, 13, 2819 7 of 20 percentage of pollutants from each sector. The above-mentioned scenarios can be rewritten as Equation (7) [14,15]: where W t is the total pollution loads; W p is the pollution load of the point source; W np is the pollution load of the non-point source; t is the time; C p (t) is the pollutant concentration of the point source at time t; and C np (t) is the pollutant concentration of the non-point source at time t. Due to the lack of continuous data obtainable from the river water quality monitoring point, the integral Equation (7) should be changed into the discrete Equation (8): We can obtain the pollutant concentration data and total pollution load with the monthly flow data, using Equations (9) and (10): where C i is the concentration from monitoring data in the ith month (mg/L); Q i is the average flow in the ith month (m 3 /s); and ∆t is the period of the ith month. W i is the average total pollution loads in the ith month (kg/d): α is the conversion factor and its calculation formula is shown in Equation (11): The results of annual average total pollution load of COD, AN and TP in Lushui River in 2019-2020 are shown in Tables 2 and 3.

Discrete Exchange
The discrete exchange is a table integrating WASP, consisting of times, dates, and dispersion coefficient values (m 2 /sec). During data input, these points are interpolated by the dispersion coefficients based on the dispersion function interpolation option selected. The river longitudinal dispersion coefficient represents the area (m 2 /s) of the pollutants in the river. The pollutants are longitudinally dispersed along the river direction per second. The integral formula and empirical estimation formula could determine it. However, the integral formula is suitable for river hydrological data and river section data, therefore, for this research, we used the empirical formula for estimation. To date, the well-known empirical formulas mainly include formulas established by Fischer [16], , Iwasa [18], and Seo [17]. According to the completeness of the collected data, Mcquivey-Keefer appears as the most straightforward and most feasible to calculate, and its calculation is noted in the formula (Equation (12)): where, E x is the longitudinal dispersion coefficient in m 2 /s, Q is the flow of the river in m 3 /s, B is the average water surface width in m, J is the hydraulic slope %, g is the acceleration of gravity m/s 2 , H is the average water depth in m, and u is the average flow velocity of the segment m/s. According to the Mcquivey-Keefer formula, the longitudinal dispersion coefficient corresponding to each river section is calculated, as shown in Table 4.

Initial and Boundary Conditions
The initial boundary conditions determine the initial value of water bodies before running the model [19]. Taking into account the influence of tributaries, we used Segment 1, Segment 4, and Segment 5 as the boundary conditions of the study area during segmentation and generalization. In establishing the quality model, it was necessary to input the flow data of the boundary conditions in the flow option (Flows) in advance; otherwise, the boundary concentration could not be inputted into the boundary condition option (Boundaries). The boundary flow and boundary concentration data were derived from the monitoring data of hydrological stations and provincial monitoring sections. Table 5 shows the initial boundary concentration.

Sensitivity Analysis and Calibration of Parameters
Parameter sensitivity analysis is an integral part of the water quality model, establishment and process. Its core purpose is to determine the sensitivity of each parameter of the model and provide a reference for parameter calibration. Parameter calibration was obtained based on sensitivity analysis. The first step was to assume a set of parameters, and then substitute them into the water quality model to output the simulated value. The next step was to compare the simulated value with the monitored value. If the difference between the simulated value and the monitored value was not too significant, then the assumed parameters could be used as the water quality model parameters. Conversely, if the simulated value was significantly different from the monitored value, the parameters could be re-adjusted until the error between the simulated value and the monitored value met a sure accuracy. During this research, the 2019 data year from January to December was used for the calibration, and the year 2020 for validating the parameters.

Sensitivity Analysis of Parameters
To calibrate the parameters for this research, we selected the local sensitivity analysis method, and the sensitivity analysis method based on the orthogonal design method.

•
Local sensitivity analysis: Local sensitivity analysis is the most used method in parameter sensitivity analysis. Its technical route is practical and straightforward. It can effectively investigate the influence of parameters on model output [20]. The basic steps of the local sensitivity analysis of the parameters in the Lushui River are as follows: 1.
Determine the analysis method: At present, standard local sensitivity analysis methods at home and abroad include the Lenhart sensitivity analysis method [21], Morris classification screening method [22], improved Morris classification screening method [3], and other methods. In this research, the Lenhart sensitivity analysis method was selected for local sensitivity analysis due to its simple principle, strong practicability, and wide application. The formula to calculate parameter sensitivity (S S ) [23] is as follows: where S S means sensitivity index, X is the parameter value, Y is the simulation value, ∆Y is the relative change of the pollution caused by the change of parameter, and ∆X is the relative change of parameter value. Table 6 shows the grading standard and level of the sensitivity. 2.
Determine the initial value of the parameter: When performing local sensitivity analysis, we determined the initial value of each parameter for analysis. For this research, the results were the parameter's initial value drawn from previous research on the WASP model.

3.
Determine the range of change (−50%, +50%): According to the model range among parameters and the research purpose, the variation, such as COD, AN, and TP, was uniformly determined to be ±50% of the parameter's initial value of study, to better compare each parameter's sensitivity. When calculating a specific parameter, we maintained other parameters as unchanged. Then, we ran the model to obtain the simulated output concentration when the parameter changed by ±50%. Finally, we calculated the sensitivity of the parameter to the simulated output concentration according to Equation (13).

•
Sensitivity analysis based on orthogonal design method: The orthogonal design method uses multiple factors and levels through the test to determine the optimal parameter combination, and fewer experiments are required to finish the test. The best advantage of this method is the more significant number of factors and levels. The local sensitivity analysis method only considers the impact of a single parameter on the model's output, ignoring the interaction between model parameters. Various parameters in the water quality model often affect each other, and there is a general phenomenon of different parameters with the same effect [20,24]. In order to better examine the influence of each parameter on the model output and the influence of the interaction between each parameter on the model output, and to determine the theoretical optimal parameter combination of COD, AN, and TP, this research used the orthogonal design method to further the parameter sensitivity analysis. There are two main methods for analyzing the results of orthogonal design: range analysis and variance analysis. In this research, we used the range analysis method to analyze the orthogonal design results, and the calculation formula of the range value is shown in Equation (14) [25]: where, R j is the range value of the j-th parameter, and K ij is the average value of the simulation results of the parameter j at the i level.
In this study, the specific steps of the sensitivity analysis of the parameters based on the orthogonal design method were as follows: 1.
Determination of the parameter level: according to the value range of each parameter, five parameter levels were set uniformly, namely 80%, 90%, 100%, 110% and 120% of the initial value of each parameter; 2.
Determination of the orthogonal design scheme: this orthogonal design was completed using the same level orthogonal table. The orthogonal table was the source of the orthogonal design method, which forms as follows: where L means the symbol of the orthogonal table, n is the number of trials arranged, m is the number of factors arranged, and k is each factor level number. According to Equation (15), there are many forms of standardized orthogonal tables, such as L 4 (2 3 ), L 8 (2 7 ), L 9 (3 4 ), L 16 (4 5 ), L 25 (5 6 ), etc. Additionally, according to the local sensitivity result and the selection principle of the orthogonal table, six factors and five levels were used to complete the orthogonal table, and arsenic was five factors and five levels. Therefore, the selection of L 25 (5 6 ) orthogonal table met this requirement. On this basis, the various levels of each parameter were numbered 1, 2, 3, 4, and 5 respectively.

3.
Performance of range analysis: this orthogonal design used the L 25 (5 6 ) orthogonal table for design, so each indicator required a total of 25 tests; that is, COD, AN, and TP all need to run the model for 25 simulations. During the simulation, the following steps occurred: collection of the output concentration values of COD, AN, and TP, inputting of the corresponding positions into the orthogonal table, calculation of the range value of each parameter according to the range analysis method, and finally, obtaining the range analysis results of each parameter. The orthogonal design method is a qualitative analysis, which mainly determines the model's critical parameters by comparing the relative magnitude of each parameter range, so the smaller range value does not affect the results of the analysis model. Therefore, the orthogonal design method was consistent with the local sensitivity analysis results. Results were determined with both analyses where the key parameters affecting the simulated output concentration were nitrification rate constant at 20 • C (K 12 ), nitrification temperature coefficient (Θ 12 ), dissolved organic phosphorus mineralization rate constant at 20 • C (K 83 ), dissolved organic phosphorus mineralization temperature coefficient (Θ 83 ), COD decay rate constant at 20 • C (K d ), and COD decay rate temperature correction coefficient (Θ d ).

4.
The analysis of theoretical and optimal parameters: in parameter sensitivity analysis, each initial parameter value is not the fixed value of the parameter. Consequently, the primary purpose of the combined analysis is to find the parameter combination that is closest to the monitored value as the initial combination of parameter rate timing, by improving the efficiency of parameter calibration. For example, the annual average values of the monitoring concentrations of COD, AN, and TP in Segment 5 in 2019 were 13.33 mg/L, 0.18 mg/L, and 0.0525 mg/L, respectively. The concentrations were compared with the output simulated value, and the closer monitoring value was the corresponding parameter combination. At that point it was considered to be the theoretical optimal parameter combination. After comparison, the theoretical optimal combination of COD, AN, and TP was determined, as shown in Table 7.

Calibration of Parameters
After the local sensitivity method and the orthogonal design method, there are four main methods currently available for parameter calibration: theoretical, experimental, empirical, and trial and error. For this research, we used the trial and error method to calibrate the model parameters. When calibrating a specific parameter, it is essential to keep the other parameters of the model unchanged. Therefore, we continuously selected values for it within its value range, and observed the error between the simulated output value and the monitored value. The trial and error method aims to reach a good result by trying out various means until errors are satisfactory and reduced to negligible [26,27]. It has frequently been used in water quality models and has accomplished sound efficacy in recent years [28]. Generally, after the simulation, if the error is less than 15%, then the calibration of the parameter is considered completed. Otherwise, the calibration will continue. The steps above are repeated until all parameter calibrations are completed, with the calibrated parameter combination fine-tuned, and the parameter calibration process completed.

•
Model verification: After completing the parameter calibration, model verification requires rationality of the model parameter selection and applicability of the model test. If the verification result is satisfactory, the model has a predictive function. If the verification result is unsatisfactory, the model structure needs to be reconsidered with a re-rating of the set parameter values until the result is ideal. For our model verification, the monitoring data and simulation output data for all segments was selected from January to December 2020, and then the calculation formula of relative error [27] was used, as shown in Equation (16): In the formula, f is the relative error between the measured value and the estimated value; C M is the measured value of each water quality index; and C S is the simulated value of each water quality index. We used the basis of general performance ratings to evaluate the model [29], as shown in Table 8. The model was verified by comparing the simulated data and measured data from the calibration on WASP. We used the statistical estimators to verify the accuracy of the simulated results suggested by the previous studies [24][25][26][27][28]. The coefficient of determination (R 2 ) determines the goodness of fit between simulated and observed data. The range value of R 2 was from 0 to 1. If the value of R 2 was close to 1, the model simulation fit well with the monitored data. Table 8 shows the general performance ratings of the coefficient of determination (R 2 ).

Water Environment Capacity
From the standards values of surface water environmental quality (GB3838-2002), we analyzed the situation to determine the division of water environment function areas and water quality targets of the Lushui River. As a result, three total capacity control indicators (COD, AN, and TP) were determined, as shown in Table 9. At present, there are three methods to calculate water environment capacity: analytical formula method, system optimal analysis method, and trial and error method, each having advantages and disadvantages [29]. First, the water environment capacity was calculated using the calibrated parameters on the WASP model. Then, through the trial and error method [14], the pollution loads value of the COD, AN, and TP were adjusted by trial and error until the water quality prediction results met the water quality objectives [29]. The steps of the simulation were as follows: 1.
Water quality aims were to determine the base level of the water environmental management conditions of the Lushui River. Due to the water environmental assessment conditions of Zhuzhou City, the goal of present environmental work in Zhuzhou to upgrade the water quality of the Lushui River to Grade II standard. Water quality standards [30,31] are shown in Table 9.

2.
To obtain various pollutant environmental capacities and simulate the pollution loads of Lushui River, we determined the water quality requirements for water quality control in the Lushui River basin as the aims of the simulation, and adjusted the input pollution loads of each indicator, such as COD, TP and AN.

3.
The water quality of the Lushui River was set to Grade II to obtain the water environmental capacity of the Lushui River, and we also considered the pollutants discharged into the Lushui River from both shores.

Parameter Sensitivity Analysis
The primary parameters were selected to calculate the sensitivity to each water quality index. The parameter change was ±50% [8], which did not include the default value of the model. The default value of the model and the formula option parameters are shown in Table 10:  According to the model parameter sensitivity analysis: 1.
For AN: only Θ 12 was high sensitivity, K 12 was medium sensitivity, and the sensitivity of the remaining parameters were small to negligible; 2.
For COD: Θ d was high sensitivity, K d was medium sensitivity, and the sensitivity of the remaining parameters were small to negligible; 3.
For TP: Θ 12 was very high sensitivity, Θ 83 was high sensitivity, K 83 was medium sensitivity, and the sensitivity of the other parameters was small to negligible.
In addition, there were situations where multiple water quality indicators were affected by the same parameter such as Θ 12 , which is a sensitive parameter for ammonia nitrogen and total phosphorus. This caused the four reaction systems of the EUTRO module [8] (that is, the N cycle, P cycle, and DO balance, and algae growth system) to interact with each other.

Model Calibration
Following the calibration process using the data from 2019 and the results of the sensitivity analysis of parameters using local sensitivity analysis, orthogonal design, and the trial and error method, the parameters K 12 , Θ 12 , K 83 , Θ 83 , K d and Θ d were calibrated, and the remaining parameters were used as the default values. The Table 11 shows the model parameter calibration results.

Simulation Results
Based on the WASP model parameters calibrated, monitoring data of AN, COD and TP from the year 2020 were used for the validation. The simulation results of each segment of the Lushui River in wet, dry and normal water seasons are shown in Table 12, and Figure 3 shows the comparison of the simulation results and measured data of COD, AN and TP in Segment 5.

Model Calibration
Following the calibration process using the data from 2019 and the results of the sensitivity analysis of parameters using local sensitivity analysis, orthogonal design, and the trial and error method, the parameters K12, Θ12, K83, Θ83, Kd and Θd were calibrated, and the remaining parameters were used as the default values. The Table 11 shows the model parameter calibration results.

Simulation Results
Based on the WASP model parameters calibrated, monitoring data of AN, COD and TP from the year 2020 were used for the validation. The simulation results of each segment of the Lushui River in wet, dry and normal water seasons are shown in Table 12, and Figure 3 shows the comparison of the simulation results and measured data of COD, AN and TP in Segment 5.   The comparison of the simulation results and measured dat Segment 5 (2020) are shown in Figure 3.

Water Environment Capacity
The calibrated WASP model was used to calculate the pollu tion of each river section as the upstream inflow water concentrati section. The results for the water environmental capacity of amm total phosphorus in wet, dry, and normal water seasons are show

Generalization and Polltion Load
The results from the annual average total pollution load Lushui River from 2019-2020 showed that loads in Segment 3 and than other segments (shown in Tables 2 and 3). Additionally, it w The comparison of the simulation results and measured data of COD, AN and TP in Segment 5 (2020) are shown in Figure 3.

Water Environment Capacity
The calibrated WASP model was used to calculate the pollutant effluent concentration of each river section as the upstream inflow water concentration of the following river section. The results for the water environmental capacity of ammonia nitrogen, COD, and total phosphorus in wet, dry, and normal water seasons are shown in Table 13.

Generalization and Polltion Load
The results from the annual average total pollution load of COD, AN and TP in Lushui River from 2019-2020 showed that loads in Segment 3 and Segment 6 were smaller than other segments (shown in Tables 2 and 3). Additionally, it was shown that the load value in Segment 4 was higher than other segments because there was an intersection of Lushui River with the primary river (Yangtze River). The loads reduced in the middle segments of the river where the area is characterized by urban land use and small-scale urban agriculture. The greatest improvement was observed at the most polluted segment of the Lushui River, which is situated precisely in Segment 1. It could therefore be contributing to the high change in the pollutant load.

Parameter Sensitivity Analysis
The sensitivity analysis conducted on WASP (Table 10) was performed by varying the parameters by (-50%, +50%). The sensitivity analysis result revealed that six parameters (K 12 , Θ 12 , K 83 , Θ 83 , K d and Θ d )were identified as medium to highly influential. The sensitivity index of all these parameters were above 0.05, indicating that they are the key parameters affecting the simulated output concentration, according to the model parameter sensitivity analysis (Tables 6 and 10). The high sensitivity of these parameters was also reported by Guo et al. [25].

Simulation Results and Model Verification
We used WASP to simulate and analyze the water quality of the Lushui River in 2020 during the normal, dry, and wet water seasons (Table 12, Figure 3). The concentration change trend of each water quality index in the three water seasons was the same. The concentration of COD was normal water seasons > dry water seasons > wet water seasons, and the concentration of AN was normal water seasons > wet water seasons > dry water seasons. The concentration of TP did not change much with the three water seasons. At the same time, from the perspective of seasonal characteristics, the concentrations of COD, TP and AN in different water periods differed mainly due to the large precipitation during the high water period, which led to the release of endogenous pollution in the river sediments. In this regard, there have been related research reports abroad. For example, Hu et al. used the WASP model for the simulation of sediment release items and found that endogenous phosphorus pollution had become a factor that could not be ignored in water quality changes [22]; whereas ammonia nitrogen and COD had different water periods. The main reason was the change of river flow. Rainfall is abundant during the high water period and the river flow is relatively large, which had a certain dilution effect on the pollutants in the river.

Model Verification and Model Accuracy Evaluation
Segment 5 was selected as the model verification segment. The simulation results of this segment from January to December 2020 were used, and then compared with the measured value. The relative error (f) results are shown in Table 14, and the WRDB on WASP was used to perform linear regression analysis on the simulated and measured values to discuss the degree of agreement. For details, see Figure 4.    According to the error analysis results (Table 14), the average relative error between the simulated value of ammonia nitrogen and the measured value was 7.418%. There were 11 groups (12 groups in total) with a relative error of less than 15%, accounting for 91.66%. The average relative error of the simulated COD value and the measured value was 6.086%, of which 11 groups (12 groups in total) had a relative error of less than 15%, accounting for 91.66%. The average relative error of the simulated total phosphorus value and the measured value was 8.637%. The relative error was that there were 10 groups (12 groups in total) less than 15%, accounting for 83.33%. Table 14 shows that the WASP model can be used as an effective tool for water quality prediction and management in the area. The linear regression equation and the accuracy evaluation statistics of the simulation results of COD, AN and TP in Segment 5 (2020) are shown in Figure 4 and Table 15. The coefficient of determination (R 2 ) sets the fit between monitored and simulated data, and its range value is from 0 to 1. If R 2 was close to 1, the model simulation fit satisfactorily with measured data. It can be seen from Figure 4 and Table 15 that the determination coefficients R 2 were all above 0.85 and the Index of Agreement (r) were all above 0.90, indicating that the simulated values were in good agreement with the monitored values and had a reasonable correlation. Similar findings relating to the overall effectiveness of water quality simulation studies in other domestic basins were reported by Guo et al. [25].

Water Environmental Capacity
The water environment capacity of AN, TP and COD in most river sections during dry seasons was smaller than that in normal and wet seasons. This was mainly due to the abundant precipitation in wet seasons, the larger river flow, and the larger water dilution environment capacity. In addition, from Table 6, it can be seen that the environmental capacity of each segment based on the water environment function was quite different, which was mainly affected by the water environment function zoning. When the water environment function was landscape recreation water area and industrial water area, the water environment capacity was relatively large; when the water environment function was a drinking water source protection area, the water environment capacity was relatively small. On the whole, the environmental capacity of the control unit based on human health benchmarks was much smaller than that of the control unit based on the water function zone standard. It can be seen that rivers frequently alternate functional zoning, protection operations are difficult, and management is inconvenient. Although the water quality can meet the standard, it does not guarantee human health very well.
The reductions in pollution load required to meet the water quality objectives were calculated by comparing the pollution load emissions from the water environmental capacity of the Lushui River [11]. Negative values indicated that the environmental capacity remained surplus to the pollution load and has no need to be reduced. Positive values of reduction indicated that the pollution load exceeded the environmental capacity and needs to be reduced. From the perspective of seasonal characteristics, according to the results from the environmental capacity and the pollution load emissions, the reduction in COD, TP and AN in different water periods was around 10 to 15%, meaning that there is a need for reduction. Therefore, in order to achieve the water quality target in the Lushui River, further research is needed to strengthen pollution source prevention and supervision.

Conclusions
In order to respond to the national protection policy for the Lushui River, a onedimensional mathematical model of the water environment was established to study the water environment capacity of ammonia nitrogen, COD and total phosphorus of the Lushui River. In this research, the WASP model was calibrated and validated using data for the years 2019-2020. The calibration results showed that, for ammonia nitrogen, only Θ 12 had high sensitivity, K 12 had medium sensitivity, and the sensitivity of the remaining parameters were small to negligible. Additionally, for the parameter sensitivity of COD, Θ d had high sensitivity, K d had medium sensitivity, and the sensitivity of the remaining parameters were small to negligible. Furthermore, for total phosphorus, Θ 12 had very high sensitivity, Θ 83 had high sensitivity, K 83 had medium sensitivity, and the sensitivity of the other parameters were small to negligible.
After the calibration process, the simulation results in the year 2020 were verified by all of the determination coefficients R 2 being above 0.85, and all of the correlation coefficients r being above 0.90. Therefore, this information indicates that the WASP model is applicable and reliable.
The water environmental capacities were obtained by calibrated parameters on the WASP model using the trial and error method. During the normal, wet and dry water seasons, the COD values were 14,072.94 tons/yr, 17,147.7 tons/yr and 10,998.18 tons/yr, respectively. However, during the normal, wet and dry water seasons, the results for ammonia nitrogen were 469.098 tons/yr, 571.59 tons/yr and 366.606 tons/yr, respectively. Moreover, total phosphorus was 93.8196 tons/yr, 114.318 tons/yr and 73.3212 tons/yr, respectively. Therefore, the results from the WASP model are applicable and reliable and can be used as an effective tool for water quality prediction and assessment in the Lushui River. It can also provide decision-making references for water environmental protection and management in the Lushui River.