Fault Detection in Wastewater Treatment Systems Using Multiparametric Programming

In this work, a methodology for fault detection in wastewater treatment systems, based on parameter estimation, using multiparametric programming is presented. The main idea is to detect faults by estimating model parameters, and monitoring the changes in residuals of model parameters. In the proposed methodology, a nonlinear dynamic model of wastewater treatment was discretized to algebraic equations using Euler’s method. A parameter estimation problem was then formulated and transformed into a square system of parametric nonlinear algebraic equations by writing the optimality conditions. The parametric nonlinear algebraic equations were then solved symbolically to obtain the concentration of substrate in the inflow, Scin , inhibition coefficient, Ki, and specific growth rate, μo, as an explicit function of state variables (concentration of biomass, X; concentration of organic matter, Sc; concentration of dissolved oxygen, So; and volume, V). The estimated model parameter values were compared with values from the normal operation. If the residual of model parameters exceeds a certain threshold value, a fault is detected. The application demonstrates the viability of the approach, and highlights its ability to detect faults in wastewater treatment systems by providing quick and accurate parameter estimates using the evaluation of explicit parametric functions.


Introduction
Traditionally, wastewater treatment is a process of converting wastewater into bilge water that can be returned to the environment, and used for domestic and industrial applications.Nowadays, wastewater treatment plants (WWTPs) also focus on sustainability issues through recovery of energy and nutrients from wastewater [1][2][3].With the increasing number of WWTPs worldwide, and the increasingly stricter requirements for maintaining the quality of effluents, on-line process monitoring has become an important aspect for ensuring efficient operation and management of WWTPs [4].It involves a process of detecting faults and diagnosing their causes and locations.This is achieved by continuously monitoring the systems to detect any abnormal conditions, and then, evaluating and diagnosing the conditions with faults [5][6][7].
Fault detection methods for sensor faults in wastewater treatment (WWT) systems have normally used data-based methods, such as neural networks (NN) and principal component analysis (PCA).The neural network model was presented in Maier and Dandy [8] to model a wastewater treatment system.In Caccavale et al.'s work [9], faults in nitrogen sensors were detected by estimating the concentration of NO and NH using neural networks.Honggui et al. [10] showed how sensor faults can be diagnosed using the fuzzy neural network to estimate dissolved oxygen concentrations, pH, chemical oxygen demand (COD), and total nutrients.In Lee et al. [11], the kernel PCA was used to extract nonlinear relations in process variables, and it showed a better performance than linear PCA in process monitoring.Adaptive PCA was used in Baggiani and Marsili-Libelli [12] to compare the current plant operation with an exact performance based on a reference data set and sensor outputs.In Sanchez-Fernández et al. [13], a distributed PCA was applied to detect faults by minimizing the communication cost between the blocks in WWTPs.The classical principal component analysis was presented using Benchmark Simulation Model No. 1 (BSM1) in Garcia-Alvarez et al. [14], Chen et al. [15], and Carlsson and Zambrano [16].The combined use of PCA in data preprocessing and artificial neural networks has been presented in Gontarski et al. [17] to improve network performance.Besides that, fault detection in WWT has been discussed using an observer-based method in Fragkoulis et al. [18], where multiple actuators and sensors faults were detected.
In an aerobic WWT system, respiration rate is used as an indicator of biological activity for monitoring and control [19,20].The respiration rate is affected by the initial condition of biomass, substrate concentration in the inflow, and extrinsic growth behavior of the biomass on inhibitory substrates.In Wimberger and Verde [19], a fault detection and isolation was performed by evaluating the detectability and isolability of analytical-and signal-based methodologies using information from applying the sensitivity theory.However, respiration rate depends on intermittent aeration patterns, and the calculation can only be evaluated during air-off periods [21][22][23].In this work, we propose fault detection in the WWT system by detecting and monitoring the kinetic parameters of extrinsic growth behavior using multiparametric programming.The main idea is to detect faults by estimating model parameters and monitoring the residual of model parameters.In parameter estimation-based fault detection, faults can be associated with the specific parameters of the model.With this assumption, parameters of a system are estimated on-line repeatedly using well known parameter estimation methods.The presence of faults is indicated if there is a discrepancy between the values of estimated parameters and the 'true' parameters.An overview for fault detection using parameter estimation can be found in [24][25][26][27][28][29][30][31].
In our earlier work [32], the fault detection method based on parameter estimation by using multiparametric programming [33][34][35][36][37][38][39] was presented.In that work, nonlinear ordinary differential equations model was converted into algebraic equations using Euler's method.Then, a square system of parametric nonlinear algebraic equations was obtained by formulating Karush-Kuhn-Tucker (KKT) optimality conditions.The model parameters were then obtained as an explicit function of the measurements by symbolically solving the equations representing KKT conditions.The estimated model parameters were compared with the normal operation for fault detection.If the residual of model parameters exceeds a certain threshold value, a fault is detected.
In this work, the concentration of substrate in the inflow, inhibition coefficient, and specific growth rate were treated as model parameters and obtained as an explicit function of the measurements using multiparametric programming, and monitored for fault detection and diagnosis.
The rest of this paper is organized as follows: Section 2 presents the parameter estimation algorithm using multiparametric programming, and in Section 3, the wastewater treatment process reaction phase model is introduced.This section also includes detailed formulation for obtaining model parameters as an explicit function of measurements to detect faults.Section 4 evaluates the feasibility of parameter estimation using multiparametric programming in fault-free and faulty scenarios.Concluding remarks are presented in Section 5.

Problem Definition
Consider the following parameter estimation problem [40]: Problem 1.
subject to: dx j (t) dt = f j (x(t), u(t), θ, t), j ∈ J (2) where x(t) represents the J-dimensional vector of state variables in the given ordinary differential equations (ODEs) system, xj (t i ) represents the measurements of the state variables at the time points t i , u(t) represents the vector of control variables, and θ is the vector of parameters.An occurrence of fault can be attributed to changes in model parameters, θ, from the nominal values.A key difficulty with this approach for fault detection is that it requires an online solution of problem 1 at regular time intervals, which is computationally demanding and prone to failure of the numerical solver for Problem 1.These limitations can be overcome using multiparametric programming (MPP) to estimate the model parameters, θ, as an explicit function of measurements, xj (t i ), by treating θ as optimization variables and xj (t i ) as the parameters.The algorithm for parameter estimation using MPP to obtain a symbolic solution for model parameters is summarized as follows [32]: (i) The nonlinear ODEs model in Equation (2) was discretized using Euler's method to algebraic equations on the interval, t ∈ [0, t f ].The Euler's method provides where the step size is given by ∆t.(ii) Fault detection problem was formulated as a nonlinear programming (NLP) problem as follows: 2 (6) subject to: x j (0) = x 0 j , j ∈ J where h j is the set of nonlinear algebraic equations obtained by discretizing the ODEs given by Equation (5), and I = {0, 1} is considered in this work.(iii) For Problem 2, the Lagrangian function is given by where and λ j represents the Lagrange multipliers.The first-order Karush-Kuhn-Tucker (KKT) conditions are given by the equality constrains as follows (iv) The equality constraints corresponding to the KKT conditions given by Equations ( 12) and (13) were solved symbolically to obtain Lagrange multipliers and model parameters, θ(x), as an explicit function of measurements, x. (v) The solutions obtained in the previous step were examined and solutions with imaginary parts were ignored.(vi) The estimated model parameters, θ, were calculated using the measurements, x, by simple evaluation of θ(x).(vii) Faults were diagnosed by monitoring residual changes in model parameters.Any significant difference between estimated and observed model parameters may be attributed to occurrence of a fault.

Wastewater Treatment System
In this work, an aerobic sequencing batch reactor (SBR) model for fed-batch reactor operation mode was considered.The bioprocesses involved in the treatment used activated sludge and provided treatment for wastewater in five stages: Fill, react, settle, decant, and idle as shown in Figure 1 [41].During the fill stage, the wastewater was directed into the tank and mixed with the sludge from previous cycles.At the reacting stage, air was provided as function of the aeration process that consumes the waste as nutrition and produces carbon dioxide, nitrates, and nitrites.After sufficient time of reaction, the aeration process was stopped and the sludge was allowed to settle.At the decant stage, the treated wastewater was removed from the reactor and the sludge that remained was reused for the next cycle.The reactor then entered the idle stage which was used to prepare the SBR for the next cycle.
(iv) The equality constraints corresponding to the KKT conditions given by Equations ( 12) and (13) were solved symbolically to obtain Lagrange multipliers and model parameters, ( ) x , as an explicit function of measurements, x .
(v) The solutions obtained in the previous step were examined and solutions with imaginary parts were ignored.
(vi) The estimated model parameters, θ , were calculated using the measurements, x , by simple evaluation of ( ) (vii) Faults were diagnosed by monitoring residual changes in model parameters.Any significant difference between estimated and observed model parameters may be attributed to occurrence of a fault.

Wastewater Treatment System
In this work, an aerobic sequencing batch reactor (SBR) model for fed-batch reactor operation mode was considered.The bioprocesses involved in the treatment used activated sludge and provided treatment for wastewater in five stages: Fill, react, settle, decant, and idle as shown in Figure 1 [41].During the fill stage, the wastewater was directed into the tank and mixed with the sludge from previous cycles.At the reacting stage, air was provided as function of the aeration process that consumes the waste as nutrition and produces carbon dioxide, nitrates, and nitrites.After sufficient time of reaction, the aeration process was stopped and the sludge was allowed to settle.At the decant stage, the treated wastewater was removed from the reactor and the sludge that remained was reused for the next cycle.The reactor then entered the idle stage which was used to prepare the SBR for the next cycle.The aerobic system involves aerobic growth and endogenous respiration reactions given by: Figure 1.The sequencing batch reactor stages (adapted from Reference [41]).
The aerobic system involves aerobic growth and endogenous respiration reactions given by: Endogenous respiration : where S c represents the concentrations of organic matter, S o is the concentration of dissolved oxygen, and X represents the concentration of biomass.The mathematical model of the process is given by the following equations [42]: ) where µ is specific growth rate, q in is inlet flow rate, V is volume, k 1 and k 2 are yield coefficients, b is endogenous respiration kinetic constant, S c in is inlet organic matter, S o in is dissolved oxygen concentrations, k La is transfer coefficient, and S o s is oxygen saturation concentration.The specific growth rate, µ, is represented by the Haldane model and is given by Equation (20).The parameter values for the wastewater treatment process reaction are shown in Table 1 [19,42].

Fault Detection Problem for the Wastewater Treatment Process
In this work, a method to estimate and detect faults in wastewater treatment is presented.The concentration of substrate in the inflow, S c in , inhibition coefficient, K i , and specific growth rate, µ o , were treated as model parameters and obtained as an explicit function of measurements using multiparametric programming and monitored for fault detection and diagnosis.By monitoring the estimated model parameters, process faults can be detected and diagnosed.We took this approach because the respiration rate indicates the biological activity, which is used for monitoring and control, and the respiration rate is affected by S c in , K i , and µ o .
Thus, the objective of this fault detection problem is to estimate the model parameters, S c in , K i , and µ o by minimizing the error of parameter estimate, ε FD , between the measurement of state variables and model predicted value of state variables as shown in problem 3.
The formulation and solution of the parameter estimation problem using MPP for WWTP are presented as follow: (i) The nonlinear ODE model in Equations ( 16)-( 19) is discretized using explicit Euler's method and reformulated as the following algebraic equations: (ii) The discrete-time fault detection problem is formulated as the following NLP: Problem 4.
subject to: (iii) Equations ( 27)-( 30) are substituted into Equation ( 26) to obtain: The derivative of g with respect to model parameters θ is then obtained and equated to zero.The resulting equality constraints are then solved analytically.Hence, the gradient of g with respect to model parameters, S c in , K i , and µ o , is given by (iv) The equality constraints in Equations ( 36)-( 38) are solved analytically in Mathematica, and the symbolic solution for model parameters is given by (v) Equations ( 39)-( 41) represent the symbolic solution for model parameters, S c in , K i , and µ o , obtained as explicit functions of the state variables, X, S c , S o , and V.In this case study, single fault was assumed to occur at any time.Hence, as an example in Equation ( 39), the solution of S c in was obtained in terms of model parameters, K i and µ o , and state variables, X, S c , and V.
In Equation ( 40), the solution of K i was obtained in terms of model parameters, S c in and µ o , and state variables, X, S c , S o , and V.The solution of µ o was obtained in terms of model parameters S c in and K i , and state variables, X, S c , S o , and V, as shown in Equation (41).Simple function evaluation was performed to evaluate the model parameters without the need for solving the online optimization problem.Then, the fault detection was performed by monitoring the residuals of model parameters.Any substantial discrepancy between estimated and observed model parameters indicates changes in the process and may be attributed to a fault.

Fault-Free Scenario
In the fault-free scenario, the simulated measured values and model predicted values for concentrations were obtained using the model parameters listed in Table 1 with initial values, X(0) = 4734 mg/L, S c (0) = 0 mg/L, S o (0) = 6 mg/L, and V(0) = 3 L; the state profiles thus obtained are shown in Figure 2. In this system, noise was added as random data to evaluate the effectiveness of the proposed method.The estimated model parameters, S c in , K i , and µ o , were calculated using Equations ( 39)-( 41) with step size, ∆t = 0.001 h and are shown in Figure 3.The result shows that the estimated model parameter is close to true model parameters.The diagnosis of fault was carried out by monitoring value of the residual value of model parameters and is shown in Figure 4.For each parameter no fault was detected as the residual was less than the threshold.The threshold was chosen as 10% from the nominal value.These results indicate that the technique proposed in this work can accurately estimate the model parameters, and this was achieved by carrying out simple function evaluations that alleviate the computational burden required for online implementation.4. For each parameter no fault was detected as the residual was less than the threshold.The threshold was chosen as 10% from the nominal value.These results indicate that the technique proposed in this work can accurately estimate the model parameters, and this was achieved by carrying out simple function evaluations that alleviate the computational burden required for online implementation.

Faulty Scenario
Investigation for the faulty scenarios was implemented where three faulty scenarios have been considered, where the percentage of change kinetic parameters is given in Table 2.In this case study, single fault was assumed to occur at any time.The step size was given as 0.001 t ∆ = h.For faulty scenario 1, the estimated concentration of substrate in the inflow is shown in Figure 5.We can see

Faulty Scenario
Investigation for the faulty scenarios was implemented where three faulty scenarios have been considered, where the percentage of change kinetic parameters is given in Table 2.In this case study, single fault was assumed to occur at any time.The step size was given as ∆t = 0.001 h.For faulty scenario 1, the estimated concentration of substrate in the inflow is shown in Figure 5.We can see that the estimated parameter S c in decreased from 168 mg/L to 118 mg/L after 3 h.The residual of S c in was monitored for fault detection, and is shown in Figure 6.This figure shows that the fault was declared at 3 h as the residual was more than a threshold value of 10%.For faulty scenario 2, the estimated half saturation coefficient, K i , is shown in Figure 7.We can see that the estimated parameter K i increased from 3.753 mg/L to 4.878 mg/L after 3 h and the residual of K i was monitored for fault detection.Figure 8 shows that after 3 h, the residual of K i was more than 10% of threshold and therefore, the fault was declared at 3 h.Figure 9 shows the estimated parameter for µ o where the estimated model parameter was decreased from 0.1916 1/h to 0.1341 1/h.The fault detection was monitored using residual, and the result is shown in Figure 10.The residual of µ o increased to 30% and indicates that faulty scenario 3 occurred at 3 h.These results indicate that there were faults at specified scenarios, and provided quick and accurate fault detection using explicit parametric functions.was monitored for fault detection, and is shown in Figure 6.This figure shows that the fault was declared at 3 h as the residual was more than a threshold value of 10%.For faulty scenario 2, the estimated half saturation coefficient, i K , is shown in Figure 7.We can see that the estimated parameter i K increased from 3.753 mg/L to 4.878 mg/L after 3 h and the residual of i K was monitored for fault detection.Figure 8 shows that after 3 h, the residual of i K was more than 10% of threshold and therefore, the fault was declared at 3 h.Figure 9 shows the estimated parameter for o µ where the estimated model parameter was decreased from 0.1916 1/h to 0.1341 1/h.The fault detection was monitored using residual, and the result is shown in Figure 10.The residual of o µ increased to 30% and indicates that faulty scenario 3 occurred at 3 h.These results indicate that there were faults at specified scenarios, and provided quick and accurate fault detection using explicit parametric functions.

Concluding Remarks
We have proposed a fault detection methodology for the wastewater treatment system using multiparametric programming where the model parameters were efficiently calculated by performing simple function evaluations without solving the online optimization problem.In this work, the related kinetic parameters for the faulty process that were investigated were in c S , i K , and o µ , which affected the respiration rate; these kinetic parameters were obtained as an explicit function of measurements.The estimation of kinetic model parameters in faulty and fault-free scenarios has shown good performance in the accuracy of parameter estimation-based fault detection.This demonstrates the advantages of multiparametric programming-based parameter estimation for detecting faults in wastewater treatment plants quickly and accurately, and reducing the online computational burden.Future work will focus on investigating the case when more than one fault simultaneously occurs.
Author Contributions: E. C. M. and V. D. contributed to the design and implementation of the research, to the analysis of the results and to the writing of the manuscript.
Funding: E.C.M. would like to thank the Ministry of Education (MoE) Malaysia, and University Malaysia Perlis (UniMAP) for the financial support.

Concluding Remarks
We have proposed a fault detection methodology for the wastewater treatment system using multiparametric programming where the model parameters were efficiently calculated by performing simple function evaluations without solving the online optimization problem.In this work, the related kinetic parameters for the faulty process that were investigated were S c in , K i , and µ o , which affected the respiration rate; these kinetic parameters were obtained as an explicit function of measurements.The estimation of kinetic model parameters in faulty and fault-free scenarios has shown good performance in the accuracy of parameter estimation-based fault detection.This demonstrates the advantages of multiparametric programming-based parameter estimation for detecting faults in wastewater treatment plants quickly and accurately, and reducing the online computational burden.Future work will focus on investigating the case when more than one fault simultaneously occurs.
Processes 2018, 6, x FOR PEER REVIEW 9 of 16 effectiveness of the proposed method.The estimated model parameters, in c S , i K , and o µ , were calculated using Equations (39)-(41) with step size, 0.001 t ∆ = h and are shown in Figure 3.The result shows that the estimated model parameter is close to true model parameters.The diagnosis of fault was carried out by monitoring value of the residual value of model parameters and is shown in

Figure 2 .
Figure 2. State variables profiles for the fault-free scenario.Figure 2. State variables profiles for the fault-free scenario.

Figure 2 .
Figure 2. State variables profiles for the fault-free scenario.Figure 2. State variables profiles for the fault-free scenario.

Figure 3 .
Figure 3.Estimated model parameters value for the fault-free scenario.(a) Concentration of substrate in the inflow, in c S ; (b) inhibition coefficient, i K ; and (c) specific growth rate, o µ .

Figure 3 .
Figure 3.Estimated model parameters value for the fault-free scenario.(a) Concentration of substrate in the inflow, S c in ; (b) inhibition coefficient, K i ; and (c) specific growth rate, µ o .

Figure 4 .
Figure 4. Residual of estimated model parameters for the fault-free scenario.(a) Concentration of substrate in the inflow, in c S ; (b) inhibition coefficient, i K ; and (c) specific growth rate, o µ .

Figure 4 .
Figure 4. Residual of estimated model parameters for the fault-free scenario.(a) Concentration of substrate in the inflow, S c in ; (b) inhibition coefficient, K i ; and (c) specific growth rate, µ o .

Fault Kinetic Model Parameter Fault 1 Fault 2 Figure 5 .
Figure 5.Estimated model parameter in c S for faulty scenario 1.

Figure 5 . 16 Figure 6 .
Figure 5.Estimated model parameter S c in for faulty scenario 1. Processes 2018, 6, x FOR PEER REVIEW 13 of 16

Figure 6 .
Figure 6.Residual of the estimated model parameter S c in for faulty scenario 1.

Figure 6 .
Figure 6.Residual of the estimated model parameter in c S for faulty scenario 1.

Figure 7 .
Figure 7.Estimated model parameter i K for faulty scenario 2.

Figure 8 .
Figure 8. Residual of the estimated model parameter i K for faulty scenario 2.

Figure 7 .
Figure 7.Estimated model parameter K i for faulty scenario 2.

Figure 6 .
Figure 6.Residual of the estimated model parameter in c S for faulty scenario 1.

Figure 7 .
Figure 7.Estimated model parameter i K for faulty scenario 2.

Figure 8 .
Figure 8. Residual of the estimated model parameter i K for faulty scenario 2.

Figure 8 . 16 Figure 9 .
Figure 8. Residual of the estimated model parameter K i for faulty scenario 2. Processes 2018, 6, x FOR PEER REVIEW 14 of 16

Figure 10 .
Figure 10.Residual of the estimated model parameter o µ for faulty scenario 3.

Figure 10 .
Figure 10.Residual of the estimated model parameter µ o for faulty scenario 3.

Table 1 .
Model parameters for the wastewater treatment process.

Table 2 .
Faulty scenario for the wastewater treatment system.
in c S

Table 2 .
Faulty scenario for the wastewater treatment system.