Abstract

A two-dimensional water quality numerical model has been developed to analyse the influence of adsorption boundary on the variation of water pollutant degradation. The characteristic finite element method is adopted to solve the discrete pollutant concentration equations, so the time-monotonicity of concentration is ensured in the convection-dominated flow field with various adsorption conditions. The second-order Runge–Kutta method is adopted to track the upstream concentration source points along the streamline, so the precision of the calculated pollutant concentration is improved even under the streamline sharp-bending condition. The flow roughness and pollutant degradation coefficient are formulated by elements according to the actual situation, so the model can quantitatively analyse the effect of various ecological restoration schemes on purification capacity in water areas. The adaptability and robustness of the model were verified by a typical engineering case, where the calculated results coincide with the measured data.

INTRODUCTION

At present, the Chinese government has increased support for water pollution control and water ecological restoration to promote the construction of ecological civilization, with many water environment renovation projects having been implemented. The main purpose of water ecological restoration projects are either to promote the regional water environmental capacity, or to improve the self-purification capacity of a water area (Jiang et al. 2014), if only considering the landscape effect, without the quantitative analysis of water quality and the later intensive management of pollutant discharge, the project is likely to be unable to achieve the desired goal (Qin et al. 2002). Therefore, the quantitative analysis of water quality is very important for the determination of the layout of the water ecological restoration project and the corresponding water body operation scheme (Xu et al. 2017).

In a water ecological restoration project, in addition to taking measures to increase the water area and improve the replacement rate of water (Wang et al. 2017), various adsorption boundaries are generally used to purify the water quality, which has an important influence on the later maintenance of the ecological restoration project. For example, in order to improve the existing water purification capacity, an artificial wetland is constructed in the surrounding lake (Hao et al. 2017; Huang et al. 2017), by circulating water adsorption, to achieve the partial purification of lake water; for the lifting of the water landscape and improving water quality, the landscape planning department usually adopts floating plants and aquatic plants in the waters during the growing season (Guan et al. 2014; Lyu et al. 2017; Zhang et al. 2017), which can adsorb contaminants in water and enhance its purification capacity, however, the plants inevitably decompose due to death and decay (Tang et al. 2013; Qi et al. 2017).

With the development of science and technology some new measures appeared, such as artificial aeration, adsorption film, nano-zeolite bed, physical and chemical reagents, etc. (Chen & Lu 2007; Liu et al. 2011; Qiang et al. 2014; Wang et al. 2016; Wu et al. 2017). The common feature of these measures is that the pollutants are partially adsorbed or degraded when they pass through the above interface, thus, improving the self-purification capacity of the water. In addition, some adsorption boundary in the change of water self-purification ability will also change the water flow conditions at the same time, such as aquatic plants making the local roughness increase and the turbulent diffusion coefficient decrease, and it is easy to form a retention area or a backwater region (Li et al. 2017; Zheng et al. 2017) which has a negative impact on the overall efficiency of water replacement. To solve the above problems, a two-dimensional flow and water quality numerical model has been developed to quantitatively analyse the influence of an adsorption boundary on the flow field and pollutant degradation law in the water area, which establishes a foundation for future engineering application.

BASIC THEORY OF THE MATHEMATICAL MODEL

Governing equations

In the mathematical model, the two-dimensional water depth average flow equations are adopted:  
formula
(1)
 
formula
(2)
 
formula
(3)
Accordingly, the transport equation of pollutant concentration is:  
formula
(4)
where u and v are flow velocity, (m.s−1); h and are the water depth and bed elevation, (m); E is the turbulent exchange coefficient, (m2.s−1); C is the Chezy coefficient, , n is roughness; c is the concentration of pollutants in the water, (mg.L−1); is the diffusion coefficient for pollutants (m2.s−1); q is the source/sink discharge (s−1); is the concentration of pollutants in the source/sink discharge (mg.L−1). The pollutant degradation term is often described by the first-order kinetic model: , where k is the degradation coefficient of pollutant in water, (d−1).

In the pollutant concentration equation, the effects of the adsorption boundary on the distribution of pollutant concentration are mainly described in the pollutant degradation term, where the degradation coefficient k value is dependent on the water temperature, velocity, adsorption boundary characteristics and the initial concentration of pollutants (Guo et al. 2008). For this reason, the finite element method is used to solve the discrete equations, so that the parameters of the flow and concentration can be specified separately for each numerical element, and by modifying the code, we can also use a specific pollutant degradation function.

Under the action of the flow field and the adsorption boundary, the pollutant concentration should be attenuated along the streamline and, for any point in the water, the pollutant concentration increases with time and then tends to the equilibrium value. Thus, a characteristic finite element method is adopted to guarantee the monotonicity of concentration under the strong convection flow condition and the adsorption boundary.

The flow boundary conditions are usually the flow velocity of the nodes at the inflow boundary and the water level of the nodes at the outflow boundary. The pollutant concentration boundary conditions include two situations: one is to set the pollutant concentration of the nodes on the inflow boundary, which can simulate the transport and degradation of upstream pollutants along the river flow. The other is to set the concentration of the nodes on the sidewall boundary, which can simulate the discharge process of pollutant sources along the river. The water flow can be either constant or unsteady, but in any case, the movement of pollutant is always an unsteady process.

Characteristic finite element method

For the convection-dominated transport problem, the variation of pollutant concentration at any point in the flow field is mainly relative to upstream pollutant transportation besides local concentration diffusion. However, the traditional finite element method is equivalent to the central difference, which needs to increase the additional numerical viscosity to suppress the oscillations that easily cause local distortion. Therefore, a characteristic finite element method is adopted in this paper. It only needs to source the characteristic point along the streamline and to simplify the transport equation. The specific methods are as follows.

In the diffusion equation of pollutant transport, the full derivative operator of the time derivative term and the convection term can be replaced by the following implicit schemes:  
formula
(5)
where is the characteristic concentration, the value of the characteristic point K at the time, and the position coordinates of the point can be expressed as:  
formula
(6)
The transport equation of pollutant concentration is discrete as follows:  
formula
(7)
By substituting Equation (5) into the above equations, the finite discrete element equation in each element is obtained:  
formula
(8)
 
formula
(9)
 
formula
(10)

In the equations above, the shape functions in the element are expressed as ; the concentrations of the three nodes are ; the characteristic concentrations of the nodes are ; and the source sink node concentrations are . Through numerical integration, the element coefficient matrix and the right end vector are obtained for each element and assembled into the global matrix to calculate the node concentration in the whole domain.

The determination of characteristic points

For each node in the element, the position coordinates of the characteristic points K should be obtained through the integrated solution in Equation (6). The specific method is shown in Figure 1, starting from the node position, tracked upstream along the streamline and finding the characteristic point. Considering the error caused by streamline bending, the second-order Runge–Kutta method is used:  
formula
(11)
where and are the position coordinates of the concentration characteristic point K; and are the position coordinates of the node P; and are flow velocities for the node; and are velocities of the halfway characteristic point M. Using the upstream tracking along the streamline to find the characteristic points, the accuracy of the characteristic concentration is improved by one order.
Figure 1

Sketch of the characteristic point position.

Figure 1

Sketch of the characteristic point position.

After obtaining the coordinates of the characteristic points (also including its halfway characteristic points), we need to find the element in which the points are located, better known as the characteristic element. And so, the characteristic concentration value can be interpolated according to the shape function and the three node concentrations at the time . In the program design, the numbers of characteristic elements for all nodes are stored in the first step. At the next step, the program initially starts the characteristic element check from the pre-stored element. If the characteristic point is not in this element, the program performs a bi-directional search from the stored number of characteristic elements to accelerate the computing speed. For the constant problem, the search for the characteristic element number is only executed once.

MODEL VERIFICATION AND ANALYSIS

The Luan-he river regulation project of Chengde city in the Hebei Province of China was taken as an example. The layout of the regulation project is shown in Figure 2. The main stream is about 4.7 km long, and averages 300 m wide. The tributary is about 1.6 km long and averages 200 m wide. The research scheme includes two cases. The first case is to analyse the variation of the flow field and pollutant distribution after the traditional river course has been carried out. The second case analyses the effects of the ecological restoration projects in specific regions on the purification capacity of the water area. In the calculation, using the river organic pollutant chemical oxygen demand (COD) as research index, according to the annual environment bulletin of water quality in Chengde city in 2011–2012, the water quality at the monitoring sections of the upstream Shang-ban-cheng bridge at inflow boundary I and the downstream Wu-long-ji bridge at the outflow boundary are all Grade IV–V, so the concentration of COD at inflow boundary I of the main stream is 40 mg.L−1, which corresponds to Grade V water quality. The water quality of the tributary is clearly better than the main stream because of the completed ecological restoration, so the pollutant concentration at inflow boundary II of the tributary is 2.0 mg.L−1. The flow of the main stream at inflow boundary I is 500 m3.s−1, the tributary flow at inflow boundary II is 100 m3.s−1, and the water level at the downstream outflow boundary is controlled at 275 m, so that both flows from the main channel and tributary will be mixed at the junction area and diluted and attenuated at the lower reaches of the river channel. The computational domain is divided into four regions where the different parameters are provided for water flow and water quality calculations according to the design schemes. A series of formulas for the turbulent exchange coefficient were evaluated by Najafzadeh & Tafarojnoruz (2016), wherein the Kashefipour & Falconer (2002) formula is adopted to calculate the turbulent exchange coefficient, and the diffusion coefficient of pollutant is 0.9 times the turbulent exchange coefficient. River roughness and pollutant degradation coefficient are selected in different parts of the river (Guo et al. 2008; Zhang 2012), with the detailed parameter settings shown in Table 1. The calculation conditions above remain constant.

Table 1

Calculation parameters of river regulation project

CaseRiver functional areaTurbulent exchange coefficient (m2.s1)RoughnessPollutant degradation coefficient (d−1)Case description
(1)Beach land 0.9 0.035 0.18 No ecological restoration measures, the degradation coefficient of pollutant remains unchanged 
(2)River band 0.6 0.04 0.18 
(3)Main channel 1.2 0.022 0.18 
(4)Beach land 0.9 0.035 0.18 
(1)Beach land 0.9 0.035 0.18 The degradation coefficient of pollutant was increased due to the planting of aquatic macrophytes in region (4) 
(2)River band 0.6 0.04 0.18 
(3)Main channel 1.2 0.022 0.18 
(4)Ecological restoration 0.6 0.08 0.35 
CaseRiver functional areaTurbulent exchange coefficient (m2.s1)RoughnessPollutant degradation coefficient (d−1)Case description
(1)Beach land 0.9 0.035 0.18 No ecological restoration measures, the degradation coefficient of pollutant remains unchanged 
(2)River band 0.6 0.04 0.18 
(3)Main channel 1.2 0.022 0.18 
(4)Beach land 0.9 0.035 0.18 
(1)Beach land 0.9 0.035 0.18 The degradation coefficient of pollutant was increased due to the planting of aquatic macrophytes in region (4) 
(2)River band 0.6 0.04 0.18 
(3)Main channel 1.2 0.022 0.18 
(4)Ecological restoration 0.6 0.08 0.35 
Figure 2

Partitions of the computing grid in the river regulation area.

Figure 2

Partitions of the computing grid in the river regulation area.

Figures 3 and 4 show the river flow field and the distribution of pollutant concentration after the conventional river regulation scheme. In the downstream of the junction area a high-velocity zone appears with a velocity magnitude of approximately 0.14 m.s−1 and is mainly located in the middle of the river. Due to the clear water discharging in the tributary, there appear obvious concentration demarcations in the junction area. The main pollutants in the river are on the right side, and along the left bank the water quality is better. This is until the downstream outlet section where the peak value of pollutant concentration is up to 30–36 mg.L−1, which still belongs to the Grade IV–V water.

Figure 3

Distribution of flow field after river regulation (m.s−1).

Figure 3

Distribution of flow field after river regulation (m.s−1).

Figure 4

Distribution of pollutant concentration after river regulation (mg.L−1).

Figure 4

Distribution of pollutant concentration after river regulation (mg.L−1).

Figures 5 and 6 are changes of river flow field and pollutant concentration distribution after the implementation of ecological restoration in a specific region. In the ecological restoration area on the right bank, the aquatic plants are densely planted to increase the local roughness, which causes the high-velocity area on the left bank side, that subsequently increases the pollutant mixture of mainstream and tributary. At the same time, the degradation coefficient increases and the pollutant concentration decreases in the ecological restoration area. The peak value of pollutant concentration decreases to 25–30 mg.L−1 at the downstream outlet, which belongs to Grade III–IV water.

Figure 5

Distribution of flow field after ecological restoration in the specified region (m.s−1).

Figure 5

Distribution of flow field after ecological restoration in the specified region (m.s−1).

Figure 6

Distribution of pollutant concentration after ecological restoration in the specified region (mg·L−1).

Figure 6

Distribution of pollutant concentration after ecological restoration in the specified region (mg·L−1).

Following comprehensive comparison, ecological restoration in the second case was performed by design in 2012. According to the environmental bulletin of Chengde city in 2015–2016, the water quality monitoring section of the downstream Wu-long-ji bridge was of Grades IV and III respectively, wherein the COD decreased to 20–30 mg.L−1. This tendency shows consistency with the conclusion of the model prediction.

CONCLUSIONS

Adsorption boundaries exist widely in nature and bring great differences in river pollutant degradation. Now, in many river regulation projects of China, some ecological or artificial restoration measures with adsorption or absorption characteristics have been used to enhance water purification capacity. The quantitative study of these adsorption boundary effects is very important for ecological environment preservation in the future. Therefore, a two-dimensional flow water quality mathematical model is developed to analyse the influence of an adsorption boundary on the concentration degradation of water pollutants. By using the characteristic finite element method, the model can ensure the time-monotonicity of the calculated concentration, even in the convection-dominated flow field with various adsorption conditions. The second-order Runge–Kutta method is used to track the location of the characteristic points along the streamline, and the precision of the calculated concentration is improved under the streamline bending conditions. By specifying the flow roughness and pollutant degradation coefficient in each element, the actual action of various sub-regional water quality improvement measures can be considered. The ecological restoration case, an actual river regulation project, is analysed and the calculated results show that the COD concentration at the downstream outlet decreases from the original 30–36 mg.L−1 to 20–30 mg.L−1 after implementing ecological restoration, which shows consistency with the measured data. In practice, the model can be used to quantitatively evaluate various ecological restoration measures in the water area, and provide support during the decision-making process.

ACKNOWLEDGEMENTS

The work described in this paper was supported by a grant from the National Key Research and Development Program of China (No. 2016YFC0401708) and the National Natural Science Foundations of China (No. 51709278 and No. 51679261).

REFERENCES

REFERENCES
Chen
F. M.
&
Lu
Q.
2007
Preparation of water purifying agent from analcite rocks and study on its fluorine removal
.
Chinese Journal of Environmental Engineering
1
(
9
),
30
34
.
Guan
Y.
,
Feng
Z.
,
Zhang
C.
,
Yang
Z.
&
Qin
Z.
2014
Applications of aquatic plants in river ecological restoration projects
.
Environmental Science and Management
39
(
5
),
149
151
.
Guo
R.
,
Li
Y. B.
&
Fu
G.
2008
Controlling factors of degradation coefficient on organic pollutant in river
.
Journal of Meteorology and Environment
24
(
1
),
56
59
.
Hao
M. X.
,
Huo
L. L.
&
Wu
S. S.
2017
Research progress on water purification of plants in constructed wetland
.
Environmental Engineering
35
(
8
),
5
10
.
Huang
G. D.
,
Du
J. Q.
,
Zhang
Y.
,
Shi
B.
,
Xing
B. Y.
&
Feng
T.
2017
Application of constructed-wetland based water treatment and optimization technologies in Taihu lake basin
.
Wetland Science and Management
13
(
3
),
10
15
.
Jiang
X.
,
Wang
S. H.
,
Yang
X. F.
,
Wang
W. W.
&
Li
J. L.
2014
Change in water quality and ecosystem of Lihu lake before and after comprehensive water environmental improvement measures
.
Research of Environmental Sciences
27
(
6
),
595
601
.
Kashefipour
S. M.
&
Falconer
R. A.
2002
Longitudinal dispersion coefficients in natural channels
.
Water Research
36
(
6
),
1596
1608
.
Li
Z. J.
,
Wang
D.
,
Zhu
M. M.
,
Nie
R. H.
&
Yang
K. J.
2017
Flow characteristics of emergent vegetation patch in channel
.
Advanced Engineering Sciences
49
(
2
),
86
91
.
Liu
X. W.
,
Xie
D. P.
,
Li
K. M.
,
Jin
Z.
,
Jiang
D.
&
Li
M. Y.
2011
Research on the impact mechanism of different aeration level on biogeochemical cycling of nitrogen in sediments
.
Ecology and Environmental Sciences
20
(
11
),
1713
1719
.
Lyu
J. Z.
,
Zhang
S. T.
,
Li
C. B.
,
Nong
Y. J.
&
Ren
Y. Z.
2017
Effect of aquatic plant on water purification
.
Environmental Science and Technology
40
(
S1
),
191
195
.
Qiang
L.
,
Shi
Y. L.
,
Xu
J.
&
Liao
W. J.
2014
Research and application progress of nanomaterials in water treatment
.
Environmental Science and Technology
37
(
120
),
252
256
.
Qin
B. Q.
,
Wu
Q. N.
,
Gao
J. F.
,
Fan
C. X.
,
Xu
G.
,
Chen
W. M.
,
Mao
R.
&
Chen
Y. W.
2002
Water environmental issues in Taihu lake of China: problems, causes and management
.
Journal of Natural Resources
17
(
2
),
221
228
.
Tang
J. Y.
,
Cao
P. P.
,
Xu
C.
&
Liu
M. S.
2013
Effects of aquatic plants during their decay and decomposition on water quality
.
Chinese Journal of Applied Ecology
24
(
1
),
83
89
.
Wang
X. M.
,
Zhai
S. H.
,
Zhang
H. J.
,
Hu
W. P.
,
Li
Q. Q.
&
Han
T.
2017
Research on appropriate hydraulic retention time on basis of water quality improvement of Lake Taihu
.
Journal of Lake Science
29
(
1
),
9
21
.
Wu
P.
,
Lu
S. J.
,
Xu
L. Z.
,
Liang
Q. Q.
&
Shen
Y. L.
2017
Efficiency and mechanism of nitrogen and phosphorus removal in modified zeolite wetland
.
Environmental Science
38
(
2
),
580
588
.
Xu
B.
,
Yang
Y. S.
,
Wang
Y.
,
Zhang
G.
,
Fan
W.
,
Lu
Y.
&
Gao
C. P.
2017
Simulation and application of water quality model in a polluted river under eco-restoration engineering scheme
.
Chinese Journal of Applied Ecology
28
(
8
),
2714
2722
.
Zhang
B. W.
2012
Calculation and value selection method of the roughness of natural river channel
.
South-to-North Water Diversion and Water Science and Technology
10
(
1
),
25
28
.
Zhang
H.
,
Wang
P. H.
,
Zhang
W. J.
,
Zuo
J.
,
Zou
X. X.
,
Ji
J.
,
An
S. Q.
&
Zhao
D. H.
2017
Improving nitrogen removal by phenophase arrangement of aquatic macrophyte species in the floating-bed constructed wetland
.
Journal of Lake Science
29
(
3
),
575
584
.
Zheng
S.
,
Wu
Y. H.
,
Baiyin
B.
&
Mu
X. P.
2017
Experimental research on Manning's roughness coefficient of an open channel with aquatic vegetation
.
Journal of Hydraulic Engineering
48
(
7
),
874
881
.