Traceability Technology for Sudden Water Pollution Accidents in Rivers

Abstract

The traceability technology for sudden water pollution accidents can be used for fast, accurate identification of a pollution source in the river. A correlation optimization model with the pollution source position and release time as its parameters is established based on hydrodynamic calculation and on the coupling relationship between forward concentration probability density and backward position probability density; and the solution of the model is realized by using a differential evolution algorithm (DEA). A coupled probability density method is to convert the traceability of a sudden water pollution accident into the optimization of two minimum values. This method is simple in principle and easy in solution, realizing the decoupling of parameter of the pollution source. The concept of gradient is introduced to the differential evolution algorithm, improving the efficiency of searching process. The proposed method of traceability was applied to the emergency demonstration project of the SNWDMRP. The results indicate that the model has good efficiency of traceability and high simulation precision and that traceability results have a certain guiding significance to the emergent regulation and control of sudden water pollution events in a river.

1. Introduction

Frequent occurrence of various kinds of water pollution events makes the research of pollutant traceability receive increasing attention, especially for long-distance water diversion projects including the SNWDMRP. The uncertainty of sudden water pollution events increases the difficulty in coping with pollution urgently. Finding out the source of water pollution is a prerequisite for realizing water quality prediction and water body pollution control. Therefore, the first task of urgently coping with sudden water pollution is to determine the source of pollution at the first time after such event, make a reasonable emergency handling plan on the basis of pollution source intensity and the occurrence place and time which have been determined, and meanwhile, provide preconditions for sudden water pollution prediction and early warning [1, 2]. The aforesaid traceability technology infers the position and time of occurrence as well as pollution source intensity through research on the transfer and conversion laws of pollutants in a river and on the basis of a monitored pollutant concentration process, realizing the reconstitution of the pollution event and playing an important role in the emergent regulation and control process of sudden water pollution events [3].

2. Basic traceability principle of sudden water pollution

Essentially, both pollutant concentration prediction and pollutant traceability fall into the process of pollutant traceability. To be more precise, the prediction of pollutant concentration distribution can be taken as the forward track process of pollutants, while the traceability of pollutants the process of backward track and traceability, as shown in Figure 1.

Figure 1.

Traceability process of sudden water pollution.

According to the information on a pollution event, the process of determining the concentration of the pollutant during the event by using the simulation technology in Chapter 2 is forward track of pollutant transfer and whereabouts, while the tack of the related information on the occurrence process and source of the event by using observed data including pollutant concentration distribution is the reversal of the forward track. Compared with the prediction of pollutant concentration distribution, the traceability of pollutant is obviously more complex, not only involving the process of track and backward reconstitution of the event but also having nonlinearity and ill-posedness as the reverse problem of prediction.

According to the connotation of sudden water pollution traceability, the definition of such traceability can be divided into five categories: The first category of traceability refers to reconstitution of unknown coefficients in a pollutant transport model according to some (observed) information of temporal and spatial distribution of a known pollutant, including longitudinal dispersion coefficient, lateral dispersion coefficient, and degradation coefficient, and the research is called a problem on parameter identification research; the second category refers to inference of the right pollution source (collection) item in the model by using observed values, including pollution source position and release intensity and time, and such traceability is also called an identification problem of the pollution source (collection) item (called a traceability problem for short); the third category is reverse inference of initial conditions on the basis of known information, and such traceability is also called a backward-time research problem; the fourth category means to infer the type or parameters of boundary conditions according to known information, namely, the reverse inference research of boundary conditions; and the fifth category is the combination of the four categories above. The pollutant traceability technology studied in this chapter is mainly about the second category traceability problem, that is, the identification problem of source (collection) items.

Attention has being paid gradually to the research of the traceability of sudden water pollution in recent years. There are two kinds of common research methods: deterministic method and probabilistic method [4, 5]. The deterministic method mainly includes a regularization method [6, 7], a trial and error method [8], a least square method [9, 10], and other optimization methods. These methods have a clear physical meaning. A single solution can be obtained by using them, and the error is often great if information is not accurate. However, the probabilistic method is a kind of random method based on Bayes inference and Markov chain Monte Carlo (MCMC) sampling, providing multiple reliable alternative results depending on the collection of distribution information on random variables [4, 11, 12]. On the basis of common traceability methods, this chapter introduces a new fast traceability method for sudden pollution in a channel, which combines the deterministic and probabilistic methods [1] and provides multiple reliable alternative results for urgently coping with sudden pollution without depending on distribution information on random variables.

3. Common traceability methods

3.1 The deterministic method

The deterministic method is developed mainly by reference to the traceability research of underground water. Perform research on pollutant concentration distribution in a simulation event by using a pollutant transfer and diffusion model, set up an optimization model which takes the error square sum between the simulation result and the result actually measured as an objective function, solve the objective function of the optimization model by using a deterministic algorithm, and then seek for, via iteration, the calculation result which best matches the values actually observed. Typically, the deterministic method includes the regularization method, trial and error method, correlation regression analysis, and the like, as well as some heuristic algorithms including genetic algorithm and simulated annealing algorithm. Regarding the tedious problems of the trial and error method, Han et al. [8] performed the inversion by using “local and basic expansion of a forward problem” and the optimal control solution by converting traceability into a minimum value problem. Jin and Chen [13] obtained a response relation between the concentration of pollutants in a channel via variation of an objective function established by using the Laplace transform. On such basis, they determine the pollution source intensity of an upstream cross section under environmental capacity control via reasoning. Chen et al. [14] transformed the calculation expression of instantaneous discharge concentration at a 1D single-point source, obtaining a linear regression model. And via regression analysis, they determined the release position, time, and intensity of the pollution source as well as the longitudinal dispersion coefficient of a channel. Min et al. [15] carried out, by using genetic algorithm, the research, respectively, on the identification problem of multiple parameters including the flow velocity, diffusion coefficient, and attenuation coefficient of a 1D river, as well as on the identification problem of right items of a 1D convection-diffusion equation. This part presents a simple introduction of the deterministic method by taking the correlation regression analysis method as an example.

When pollutants are released instantaneously, the analytic solution of the 1D convection-diffusion equation can be further determined via Fourier transform or dimensional analysis. Assuming that the occurrence position and time of pollution are $x_0$ and $t_0$, respectively, then we obtain

$$C=\frac{M}{\sqrt{4\mathit{\pi D}\left(t-t_0\right)}}exp\left(-\frac{{\left(\left(x-x_0\right)-u\left(t-t_0\right)\right)}^2}{4D\left(t-t_0\right)}\right)$$E1

where $M$ is the average intensity of the initial surface source of the released pollutant along a cross-section, $\mathrm{g}/{\mathrm{m}}^2$.

After taking logarithmic transformation of both ends of Equation (1), we obtain

$$lnC=\mathit{aX}+b$$E2

where $X={\left(\left(x-x_0\right)-u\left(t-t_0\right)\right)}^2$, $a=-\frac{1}{4D\left(t-t_0\right)}$, and $b=ln\frac{M}{\sqrt{4\mathit{\pi D}\left(t-t_0\right)}}$.

For a fixed pollution event, the occurrence position $x_0$ and the time $t_0$ are constants. Therefore, if the concentration of multiple cross sections can be obtained at the same time and the control time parameter $\left(t-t_0\right)$ is a constant, then it is easy to judge that $lnC$ and X (i.e., the observation position x) meet a linear relation; and a linearly dependent optimization model can be built by using a regression analysis method [14]:

$$R=\frac{{\sum }_{i=1}^n\left(ln{C}_i-\overline{lnC}\right)\left(X_i-\overline{X}\right)}{\sqrt{{\sum }_{i=1}^n{\left(ln{C}_i-\overline{lnC}\right)}^2}\sqrt{{\sum }_{i=1}^n{\left(X_i-\overline{X}\right)}^2}}$$E3

With the position of the monitored cross section as a variable, after selection of x₀ and t₀, the concentration measured C_i at different monitored cross sections x_i at the same time and corresponding to X_i, and at the monitored position corresponding to the same time are obtained via calculation. Theoretically, when x₀ and t₀ are selected as the true position and time of occurrence, R obtained via calculation should be 1.0.

From this, the problem of pollutant traceability is converted into determination of proper x₀ and t₀ via Formula (3), and the R meeting calculation is maximum, that is, 1.

The solution can be completed in five steps:

1. Firstly, calculate a proper $x^{\prime }=x_0+u\left(t-t_0\right)$, and it is selected to be 1.0 if meeting R. The calculation can be realized via derivation.

2. Calculate X_i corresponding to different observed cross sections at the same time according to known $x^{\prime }$, fit $lnC$ and $X$ in combination of Formula (2), calculate the corresponding slope a and intercept $b$, and therefrom, obtain $t_0$ via calculation.

3. Calculate the time parameter $t_0=t+\frac{1}{4\mathit{Da}}$ according to calculated slope $a$, then $x_0=x^{\prime }-u\left(t-t_0\right)$ in combination of known $x^{\prime }$ and finally the intensity of the initial surface source $M=\sqrt{4\mathit{\pi D}\left(t-t_0\right)}exp\left(b\right)$ in combination of known $b$.

4. The intensity of discharged pollutants calculated by using flow cross-section area $A$ at time $t_0$ at estimated pollution position $x_0$ is MA.

5. Considering the existence of fitting and observation errors, calculate $R$ via back substitution of obtained $x_0$ and $t_0$, then judge whether $R$ is 1. If $R$ is 1 (or the error is very small), then complete traceability calculation. Otherwise go back to Step (2), properly adjust the slope $a$ and intercept $b$, and then recalculate $x_0$ and $t_0$.

From the aforesaid solution process, we can see that the observed data of multiple cross sections are needed at the same time for the model to realize the traceability of the sudden pollution due to a point source.

3.2 The probabilistic method

Because a traceability problem has ill-posedness, when the deterministic method is used for solution, the error in observation or model calculation might result in great deviation in results and distortion in the result of traceability. For this, random devices are introduced to the study of traceability. Among such devices, the one mostly used is the probabilistic method based on Bayes inference and MCMC sampling. Bayes inference is a kind of method based on the theoretical foundation of probability theory, which can reflect the uncertainty of sudden water pollution events in a channel. It solves the posterior probability distribution of such parameters on the basis of making full use of a likelihood function and the priori information of parameters to be determined and then obtains the estimated values of all parameters of a pollution source by respective sampling. This method can offer a random distribution function of the traceability result of a water pollution event. Therefore, the methods based on Bayes inference are mainly used to estimate happening probability of a sudden water pollution event. By using them, we can obtain the posterior probability distribution of the traceability result instead of a single solution, quantitate the uncertainty of the result, and gain more traceability information on the event. In order to acquire the estimated value of the result, Bayes inference should be combined with related sampling methods, such as Markov chain Monte Carlo (MCMC), randomized quasi-Monte Carlo (Monte Carlo, MC), and other sampling methods [11, 12]. Among them, the MC method is an estimation method in which an initial value easily converges to a suboptimal solution no matter whether the initial value deviates from its true value. Therefore, the accuracy rate of the traceability result obtained via this method is not high. The defect in the MC method can be compensated by iterating the mode of Bayes inference’s combination with the MC or MCMC method into the obtained distribution function of the traceability result. The MCMC method is a sufficiently long Markov chain acquired via random walk, only by which can the proximity of the sampling result to the posterior distribution of the traceability result be ensured, that is, the limit distribution of Markov chain is used to express the posterior probability density function of the traceability result. For this, the MCMC method expands the application of Bayes inference in the traceability research of environmental pollution events.

Taking Bayesian MCMC method, for example, this part provides a simple introduction of the traceability method in probability statistics. The method takes all variables in a traceability research model as random ones and considers that the solution to the traceability problem of a sudden water pollution event is a probability distribution. Firstly, it converts the priori information of the solution into a priori probability distribution by using Bayes method, then combines with observed data, and finally obtains the posterior probability distribution of the solution to be determined from the random sampling process of Markov chain by using the likelihood function between observed values and the values calculated via simulation. The probabilistic method includes the following main processes:

1. Model a basic traceability problem, and quantitate the priori information quantity of related traceability solutions (the intensity, occurrence position, and time of a pollution source) by adopting a reasonable probability distribution function.

2. Select and establish a reasonable likelihood function based on a channel water quality coupling simulation model and in combination of the information related to the site of a sudden water pollution event and to the hydrologic data of water space of the event.

3. Obtain a posterior probability distribution of traceability solutions based on the prior probability distribution and the likelihood function.

4. Perform sampling of the posterior probability distribution and acquire the estimated values of traceability solutions to the event.

Firstly, express the probability distribution table of the priori information on the known parameter vector $\theta$ before collecting observed data into $p\left(\theta \right)$. After obtaining observed data, the posterior distribution of the known parameter acquired via Bayes inference is $p\left(\theta \left|d\right)\right)$, meeting

$$p\left(\theta \left|d\right)\right)=p\left(\theta \right)p\left(d\left|\theta \right)\right)/p\left(d\right)$$E4

where $\theta$ stands for three parameters of the pollution source; d is the measured value of pollutant concentration.

$p\left(d\right)$ is the probability distribution of the measured value. $p\left(d\right)$ obviously has no relation with parameter $\theta$. In the case with a known concentration actually measured, $p\left(d\right)$ can be understood to be 1. From this, we obtain

$$p\left(\theta \left|d\right)\right)=p\left(\theta \right)p\left(d\left|\theta \right)\right)$$E5

The prior distribution $p\left(\theta \right)$ of parameter $\theta$ in Formula (5) can be deemed a uniform distribution of a respective parameter within a prior range of values.

$$p\left(\theta \right)=\prod U\left({\theta }_i\right)$$E6

Therefore, to determine the posterior distribution $p\left(\theta \left|d\right)\right)$ of pollution source parameters, we firstly need to determine distribution $p\left(d\left|\theta \right)\right)$, which can be defined as the likelihood function of a measured and a predicted value. Let $d_i$, $C_i\left(x,t\left|\theta \right)\right)$, and $p\left(d_i\left|\theta \right)\right)$ to be the measured value, predicted value, and likelihood function of no. $i$ measuring point, respectively, and ${\epsilon }_i=d_i-C_i\left(x,t\left|\theta \right)\right)$ to be an error of measurement, $i=1,2,\dots ,N$. Assuming that obedience mean of error ${\epsilon }_i$ is 0, standard deviation is the normal distribution of ${\epsilon }_i$, and all measuring points are independent from each other, then

$$p\left(d\left|\theta \right)\right)=\prod _{i=1}^{N}p\left(d_i\left|\theta \right)\right)=\frac{1}{{\left(\sqrt{2\pi }\sigma \right)}^N}exp\left(-\sum _{i=1}^N\frac{{\left(d_i-C_i\left(x,t\left|\theta \right)\right)\right)}^2}{2{\sigma }^2}\right)$$E7

Thereby, we obtain the posterior distribution of pollution source parameters $p\left(\theta \left|d\right)\right)$:

$$p\left(\theta \left|d\right)\right)=p\left(\theta \right)p\left(d\left|\theta \right)\right)=\prod U\left({\theta }_i\right)\cdot \frac{1}{{\left(\sqrt{2\pi }\sigma \right)}^N}exp\left(-\sum _{i=1}^N\frac{{\left(d_i-C_i\left(x,t\left|\theta \right)\right)\right)}^2}{2{\sigma }^2}\right)$$E8

However, complicated $C_i\left(x,t\left|\theta \right)\right)$ or great model parameter space and number of dimensions make $p\left(\theta \left|d\right)\right)$ very abstract and difficult to be expressed visually. Therefore, direct use of Bayes method almost cannot solve an actual problem directly, which, however, can be solved by means of the MCMC method. According to different transition probability matrixes making up the Markov chain, the MCMC method mainly has the following sampling algorithms: Gibbs sampling algorithm, Metropolis-Hastings algorithm, and self-adaptive Metropolis algorithm. The self-adaptive Metropolis algorithm can converge to a target distribution for any prior distribution of $\theta$, so this algorithm is selected for sampling.

According to the aforesaid derivation process, the main steps for traceability of a sudden pollution event based on Bayes inference and the MCMC method are as follows:

1. Suppose i = 0, and then initialize different variables.

2. Produce and accept random variables, and make up a Markov chain.

   1. The priori parameter $\theta$ producing a uniform distribution equals $\theta \left(m,x_0,t_0\right)$.

   2. From the produced pollution source parameters, calculate the concentration value at the observing point by using the methods in Chapter 3.

   3. Calculate the likelihood function $p\left(\theta \left|d\right.\right)$ by using Formula (8).

   4. Calculate the acceptance probability of the Markov chain by using the following formula:

      $$\alpha =min\left\{1,\frac{p\left({\theta }^{\ast }\left|d\right)\right)}{p({\theta }_i\left|d)\right)}\right\}$$

   5. A random number $R$ is produced which uniformly distributes within 0–1. If $R<\alpha$, then accept the parameters for this test and make ${\theta }_{i+1}={\theta }^{\ast }$. Otherwise let original parameter ${\theta }_{i+1}={\theta }_i$.

3. Repeat Steps (1)–(5) until reaching the number of predefined iteration or acquiring the posterior sample number preset for pollution source parameters, and then perform statistics of the posterior distribution law of all parameters and complete the traceability calculation.

From the aforesaid process, we can see that the probability statistics method depends on the priori range of values of parameters and the information of calculation error distribution and that the convergence rate of producing parameter values randomly by using uniform distribution is very slow if the priori range of values is large. Moreover, the probabilistic method reaches solution based on Markov chain Monte Carlo random sampling, without directivity for making no optimal judgment of existing results, also resulting in a lower efficiency of traceability.

4. Coupled probability density method

4.1 Basic principles and method

According to the introduction and analysis of the common methods mentioned above, the deterministic method has a clear physical meaning but a single solution, resulting in a great error in the case with inaccurate information. Although the probability statistics method can offer multiple reliable alternative results, it depends on the acquisition of distribution information on random variables. For the point of view of research trend, a new generation of traceability method which integrates the deterministic method and the probability method will become the future trend of pollutant traceability research. The coupled probability density method proposed in this chapter combines the deterministic method with the probability method and sets up a sudden water pollution traceability model based on a coupled probability density function (C-PDF). Based on hydrodynamic calculation and considering the observation error of a system, an optimization model is built by taking pollution source position and release time as parameters and through correlation analysis of the forward concentration distribution probability density and backward position probability density of pollutants, and then the solution of such model is realized by using DEA. On such basis, a minimum value optimization model is established based on the forward concentration distribution probability density function of pollutants to determine the intensity of the pollution source.

In case of a sudden event, the majority of pollutants enters a channel in the form of a point source; move and convert with water flowing in the channel. Based on traceability research at home and abroad and the transfer and conversion law of pollutants in a channel, this paper sets up a 1D traceability model for sudden river water pollution depending on a coupled probability density function (C-PDF) and infers the pollution source in the light of observed concentration of pollutants, thus realizing the reconstitution of pollution event. Based on 1D hydrodynamic calculation and considering the observation error of a system, an optimization model is built by taking pollution source position and release time as parameters and through correlation analysis of the forward concentration distribution probability density and backward position probability density of pollutants, and then the solution of such model is realized by using DEA. Meanwhile, a minimum value optimization model is established based on the forward concentration distribution probability density function of pollutants to determine the intensity of the pollution source.

The magnitude of pollutant concentrations at different cross sections in a channel can be expressed as the statistic at such cross sections where tiny substance particles appear at a respective time. Statistically, the high or low probability that substance particles appear at some position can be expressed by a probability function, and the concentration distribution of pollutants in a channel can be also described by some probability density function. On the contrary, without knowing the pollution source, the pollutants observed at some cross section in a river might come from any upstream place, whose probability also can be described by a probability density function. Considering that pollutant traceability is the reverse problem of pollutant concentration distribution, a probability density function used to describe such distribution is called a forward concentration probability density function. However, a probability density function used to describe the probability that the pollution source is at different positions is called a backward position probability density function. Depending on the definition of the probability density function, a forward probability density function can be obtained after normalization of concentration.

According to $\underset{x}{\int }C\left(x,t\right)\mathit{dx}=m_0$, the normalization of $C\left(x,t\right)$ can be done, as shown below:

$$c\left(x,t\right)=C\left(x,t\right)/m_0$$E9

where $c\left(x,t\right)$ is the forward concentration probability density value corresponding to $C\left(x,t\right)$, with a dimension of m⁻¹, indicating the probability that pollutants appear at cross section x at time t.

According to Neupauer and Wilson’s derivations, forward concentration transport and backward position traceability are adjoint processes to each other [16, 17]. For a 1D channel, $P\left(x_s,t^{\prime }\right)$ is used to express the probability, judged from an observed cross-section x_d, that a pollution source is at x_s at time $t^{\prime }$ (i.e., the probability that pollutants are transported from cross-section $x_s$ to cross-section $t_d-t^{\prime }$ within time $x_d$, then $P\left(x_s,t^{\prime }\right)$ meets the adjoint state equation of a convective diffusion Eq. (10)) and normalization conditions ($P\left(x_s,t^{\prime }\right)$ also has a dimension of m⁻¹), as shown in the following formula:

$$-\frac{\partial P\left(x_s,t^{\prime }\right)}{\partial t}+\frac{\partial \left(\mathit{uP}\left(x_s,t^{\prime }\right)\right)}{\partial x}+D\frac{{\partial }^{2}P\left(x_s,t^{\prime }\right)}{{\partial }^{2}x}=0$$E10

$$P\left(x_d,t_d\right)=1$$E11

where $t^{\prime }$ is the time point for inverse calculation and $t_d$ is the time point for observing pollutant concentration. Formula (11) indicates that the pollution source only can be at the observed cross section when pollutants are not transported ($t^{\prime }=t_d$) but appear at such cross section.

Like concentration transport process, Formula (11) can be considered a position probability transport process. Similarly, when u is constant, an analytic solution can be obtained, as shown below:

$$P\left(x_s,t^{\prime }\right)=\frac{1}{\sqrt{4\mathit{\pi D}\left(t_d-t^{\prime }\right)}}exp\left(-\frac{{\left(x_d‐x_s-u\left(t_d-t^{\prime }\right)\right)}^2}{4D\left(t_d-t^{\prime }\right)}\right)$$E12

We can see that $P\left(x_s,t^{\prime }\right)$ and $C\left(x,t\right)$ have the same form. Actually, the relation between Formulas (1) and (12) decides that notwithstanding flow field, the relation between $P\left(x_s,t^{\prime }\right)$ and $C\left(x,t\right)$ can be determined by Figure 2, on which the arrows show the direction of transport (advection item).

Figure 2.

Processes of forward concentration transport and backward position probability transport. a and b are the concentration probability distributions along the river at release time t₀ and after Δt, respectively; c and d are the position probability density distributions at monitoring time t_d and of backward traceability for Δt, respectively.

Pollutant concentration distribution and traceability are inverse problems to each other, but both of the problems follow the same physical law basically. Figure 2 shows that when $t-t_0=t_d-t^{\prime }$, the probability $c\left(x_d,t\right)$ that pollutants move from source $x_0$ to cross-section $x_d$ after time $t-t_0$ equals to the probability, judged by an observer at cross-section $x_d$, that pollutants move from cross-section $x_0$ to cross-section $x_d$ after time $t_d-t^{\prime }$. Namely, if $t-t_0=t_d-t^{\prime }$, the following Formula (13) holds

$$P\left(x_0,t^{\prime }\right)=c\left(x_d,t\right)$$E13

From Figure 2 and Formula (13), we can see that forward concentration transport process is highly coupled with backward position probability transport process. The two processes are exactly the same except for backward-time calculation direction. Because $P\left(x_s,t^{\prime }\right)$ has no relation with the intensity of the pollution source, $P\left(x_s,t^{\prime }\right)$ can be directly calculated by using Formula (10) or (12) with $x_s=x_0$ and $t^{\prime }=t_d+t_0-t$ given. Therefore, a linearly dependent model can be established by substituting $P\left(x_s,t^{\prime }\right)$ for $C\left(x_d,t\right)$ based on such coupling relation in Formula (13) to realize traceability calculation.

Assuming that observed concentration series is $C_i$ and that the corresponding calculated position probability density series is $P_i$, $i=1,2,\dots n$, according to the aforesaid derivation, we can obtain the expression of the correlation coefficient r of the two series is as follows:

$$r=\frac{{\sum }_{i=1}^n\left(C_i-\overline{C}\right)\left(P_i-\overline{P}\right)}{\sqrt{{\sum }_{i=1}^n{\left(C_i-\overline{C}\right)}^2}\sqrt{{\sum }_{i=1}^n{\left(P_i-\overline{P}\right)}^2}}$$E14

where $\overline{C}$ and $\overline{P}$ are the arithmetic averages of C_i and P_i, respectively.

According to Formulas (13) and (14), we can set up the following function:

$$min\left(\mathit{abs}\left(1-r\right)\right)$$E15

The constraint condition is the value range of x₀ and t₀, which is given by priori information and generally is obtained via field survey or existing information. See Formulas (16) and (17):

$$x_{0min}\le x_0\le x_{0max}$$E16

$$t_{0min}\le t_0\le t_{0max}$$E17

Through solution of the aforesaid optimization model, we can obtain the release position x₀ and time t₀ of the pollution source.

The optimization model realizes the decoupling of pollution source intensity and time parameter, so determine position and time parameters firstly and then the intensity of the nonpoint pollution source m₀. Considering that the occurrence position and time of the pollution source have been determined and that the range of pollution source intensity m₀ can be roughly determined by the previous model, the forward probability density function can be used to build an optimization model to calculate pollution source intensity. See Formulas (18) and (19).

$$\text{Target function}min\left(\sum {\omega }_i{\left(m_0\ast c_i-C_i\right)}^2\right)$$E18

$$\text{Constraint condition}{m}_{0min}\le m_0\le m_{0max}$$E19

where coefficient ${\omega }_i=\frac{1}{{\left(C_i+1.0\right)}^2}$ is used to remove small concentration error loss arising from a large concentration difference.

After verification and analysis, the optimization models (18) and (19) can not only solve the intensity of pollution source m₀ rapidly but also have an advantage in observed error filtration.

4.2 Model solution

The model uses a smart differential evolution algorithm (DEA) method for solution, which is similar to the genetic algorithm [1, 18]. Solving a traceability optimization model by using the DEA includes population initialization, mutation, recombination, and selection. Take the first optimization model made up by Formulas (14) and (16), for example. The DEA defines $x_0^{\prime }$ and $t_0^{\prime }$ as property elements and produces basic evolutionary individuals, with $F=\mathit{abs}\left(1-r\right)$ as a fitness target function. The detailed steps of solution are as follows:

1. Population initialization: let the scale of population equal to NP, and then produce the first individual $X_i^{\mathrm{G}\mathit{en}}\left(x_0,t_0\right)$ from Formula (20) ($i=1,2,\dots ,\mathit{Np},\mathit{Gen}$ are the individual and evolution generation numbers, respectively):

   $$X_i^{\mathit{Gen}}=X_{\min }+\mathit{rand}\left(0,1\right)\left(X_{\max }-X_{\min }\right)$$E20

   where $\mathit{rand}\left(0,1\right)$ stands for a random number within [0,1]; $X_{\min }$ and $X_{\max }$ are the minimum and maximum values of ${X_i}_{\mathit{Gen}}^{}$, respectively. After initialization, calculate their fitness values $F\left(X_i^{\mathit{Gen}}\right)=\mathit{abs}\left(1-r\left(X_i^{\mathit{Gen}}\right)\right)$, respectively.

2. Mutation: select three individuals $X_{r1}^{\mathit{Gen}}$, $X_{r2}^{\mathit{Gen}}$, and $X_{r3}^{\mathit{Gen}}$ from population $\mathit{Gen}\left(1<\mathit{Gen}<max\mathit{Gen}\right)$ via uniform sampling, and then produce a mutation individual $V_i^{\mathit{Gen}}$ from Formula (21):

   $$V_i^{\mathit{Gen}}=X_{r1}^{\mathit{Gen}}+\mathit{CF}\left(X_{r3}^{\mathit{Gen}}-X_{r2}^{\mathit{Gen}}\right)$$E21

   where $\mathit{CF}$ is a scaling factor, generally $\mathit{CF}\in \left[0.5,1\right]$, $r1,r2,r3=1,2,\dots ,\mathit{NP}$, and $r1,r2,r3$ are all not i. If $V_i^{\mathit{Gen}}\notin \left[X_{\min },X_{\max }\right]$, then regenerate a mutation individual $V_i^{\mathit{Gen}}$ from Formula (20).

3. Recombination: before and after mutation, each property element of an individual undergoes recombination and produces a new individual $U_i^{\mathit{Gen}}$. The rule of recombination is as follows:

   $$U_i^{\mathit{Gen}}\left(x_j\right)=\left\{\begin{array}{l}{V}_i^{\mathit{Gen}}\left(x_j\right),\mathit{rand}\left(0,1\right)<\mathit{CR}\mathit{\text{or randn}}=j\\ X_i^{\mathit{Gen}}\left(x_j\right),\mathit{\text{others}}\end{array}\right.,j=1,2$$E22

   where CR is a recombination constant, generally CR ∈ [0.8,1], $U_i^{\mathit{Gen}}\left(x_j\right)$ refers to the property value of $x_j$ in $U_i^{\mathit{Gen}}$, and $\mathit{\text{randn}}$ is a random value of 1 or 2.

4. Selection: superior individuals replace an original individual and go into the next generation. Select superior individuals according to Formula (23):

   $$X_i^{\mathit{Gen}+1}=\left\{\begin{array}{l}{U}_i^{\mathit{Gen}},F\left(U_i^{\mathit{Gen}}\right)<F\left(X_i^{\mathit{Gen}}\right)\\ X_i^{\mathit{Gen}},F\left(U_i^{\mathit{Gen}}\right)>=F\left(X_i^{\mathit{Gen}}\right)\end{array}\right.$$E23

According to the given steps, perform cycle evolution and iterative computation until adaptability meets requirements (the value of target function meets limited conditions, for example, it can be controlled by $r\ge 0.95$) or ends when it evolves to a maximum algebra $max\mathit{Gen}$. After that, select the individual having the minimum value of fitness function (the parameter values that individual elements stand for are x₀ and t₀ to be determined), and the calculation is completed.

4.3 Improvement to coupled probability density method

The model uses the DEA for solution, the result of calculation has a great randomness, and some optimized results might appear repeatedly. Therefore, some optimized rules should be added to control the value of optimized fitness function to make it change unidirectionally. For this, this chapter puts forward the application of gradient concept in a differential evolution process to provide an optimization direction and improving optimization efficiency.

Firstly, add two new attributes signifying gradient into an individual, namely, gradient feature factor $\mathit{gc}\left(X\right)$ and gradient direction factor $\mathit{gp}\left(X\right)$. A population individual newly defined is $S_i^{\mathit{Gen}}=S_i^{\mathit{Gen}}\left(X_i^{\mathit{Gen}},\mathit{gc}\left(X_i^{\mathit{Gen}}\right),\mathit{gp}\left(X_i^{\mathit{Gen}}\right)\right)$, where the gradient feature factor $\mathit{gc}\left(X_i\right)$ and gradient direction factor $\mathit{gp}\left(X_i\right)$ are determined from Formulas (24) and (25).

$$\mathit{gc}\left(X_i^{\mathit{Gen}}\right)=\left\{\begin{array}{ll}1,& F\left(X_i^{\mathit{Gen}}\right)\le F\left(X_i^{\mathit{Gen}-1}\right)\\ 0,& F\left(X_i^{\mathit{Gen}}\right)>F\left(X_i^{\mathit{Gen}-1}\right)\end{array}\right.$$E24

$$\mathit{gp}\left(X_i^{\mathit{Gen}}\right)=\left\{\begin{array}{l}\left[\mathit{sgn}\left(X_i^{\mathit{Gen}}\left(x_1\right)-X_i^{\mathit{Gen}-1}\left(x_1\right)\right),\mathit{sgn}\left(X_i^{\mathit{Gen}}\left(x_2\right)-X_i^{\mathit{Gen}-1}\left(x_2\right)\right)\right], \mathit{gc}\left(X_i^{\mathit{Gen}}\right)=1\\ 0, \mathit{gc}\left(X_i^{\mathit{Gen}}\right)=0\end{array}\right.$$E25

where sgn is a symbolic function.

After setting gradient attributes, the trend information that facilitates population evolution can be kept during differential evolution selection. For example, supposing an individual $X_i^{\mathit{Gen}}$ is superior to its previous generation of individuals $X_i^{\mathit{Gen}-1}$ and $\mathit{gp}\left(X_i^{\mathit{Gen}}\right)=\left[1,-1\right]$, then this means that the transient superior evolution trend of this individual is “increasing position parameter $x_0$ and reducing time parameter $t_0$.” The information on such trend can be used to guide the evolution of the next generation of individuals. When making improvement to DEA optimization in actual applications, we can utilize and obtain gradient information as an additional rule to control evolution direction and guide the evolution of each generation of individuals.

5. Application examples

As shown in Figure 3, a emergent treatment demonstration project of sudden water pollution events was carried out from the check gate in the Fangshui River at Jingshi Section (Pile 100 + 294.750) to the check gate in the Puyang River (Pile No. 113 + 492.750) of the SNWDMRP on March 22, 2014. Cane sugar was used as tracer, and the concentration of cane sugar solution, as a water quality detection index. The experimental channel has a trapezium cross section, whose parameters are bottom width = 18.5 m, side slope = 2.5, average water depth = 4.0 m, channel bottom slope = 0.00005, manning roughness = 0.015 (design value), the experimental estimated value of dispersion coefficient = 3.43 m²/s, and constant flow ≈ 8.0m³/s. About 800 kg of cane sugar was put instantaneously into the channel at 9:00 AM on March 22, and four monitoring cross sections were set downstream. In order to ensure successful instantaneous putting of 800 kg cane sugar, the cane sugar was dissolved in heated clear water of 1.0 m³ (about 80.0°C) timely and in advance. When the experiment started, a floating bridge already built was used to directly pour dissolved cane sugar solution into the middle of the channel. The concentration monitoring process at 1508 m downstream of the cross section where the solution was poured is shown in Figure 4, and the start time of monitoring was 9:30.

Figure 3.

Schematic diagram of experiment channel selected for demonstration project.

Figure 4.

Actual concentration monitoring process at monitored cross section.

Supposing the pile number of the pouring cross section is 0 m, then the pile number of the monitored cross section is 1508 m. The width and water depth at the experimental section of the channel were very small with regard to the length. After checking calculation, the pollutants 1508 m from the monitored cross section were uniformly mixed roughly in the lateral and vertical directions. The time for urgently coping with sudden pollution is short. Without considering degradation, the transport of pollutants can be described by using a convective diffusion equation at a 1D constant flow. A traceability model is applied for traceability calculation, which was established by using the probability statistics method and the proposed coupled probability density method among conventional traceability methods. The results are shown in Tables 1 and 2.

	M (kg)	x0 (m)	t0
Actual value	800.00	0.00	9:00
The value calculated by the probabilistic method	781	−570	9:14 AM
The value calculated by the coupled density method	731.13	−28.42	8:48 AM

Table 1.

Table of results of traceability by all methods.

	M (kg)	x0 (m)	t0
The error calculated by the probabilistic method	2.38%	570	14 min
The error calculated by the coupled probability density method	8.59%	28.42	12 min

Table 2.

Analysis of errors of traceability by all methods.

From Tables 1 and 2, we can see that the error in the position of a pollution source calculated by using the coupled probability density method is not up to 30 m, time error about 12 min, and the intensity of pollutants not larger than 10%. However, considering the calculation done by the probability statistics method, although the precision of its pollutant intensity is higher, the position and time are not calculated more accurately, especially the error of the position calculated is above 500 m. Furthermore, to solve the same problem, it takes much more time for calculation by the probability statistics method than by the probability density method. As a whole, when a sudden water pollution accident happens during water diversion, the adoption of the traceability technology based on the coupled probability density method is recommended more.

6. Summary

Due to uncertainty of sudden water pollution events, after occurrence of an event, it is very difficult to determine the location and release time by experience. Therefore, the first task for emergent treatment is to determine the pollution source. The main content of fast traceability of sudden pollution is to determine the location, release time, and the intensity of the pollution source. Regarding fast determination of the conditions of a pollution source in case of sudden water pollution under SNWDMRP, this chapter presents detailed introduction and analysis of the fast channel traceability technology including conventional traceability, provides practical effect of the traceability technology via analysis of example applications, and recommends the application of the coupled probability density method to traceability calculation according to traceability results, in order to provide support for fast traceability and emergent regulation and control of sudden water pollution under long-distance water transfer projects.

However, the traceability model can have good application results for single-point source of sudden water pollution accident in a river. The traceability model for multipoint sources of sudden water pollution accidents needs deeper research and can refer to the ongoing research by Wang [19].