Identification and control of anaerobic digesters using adaptive, on-line trained neural networks

Computers chem. Engng Vol. 21, No. 1, pp. 113-143, 1997

Pergamon 0098-1354(95)002,44-8

Copyright © 1996 Elsevier Sciencc Ltd Printed in Great Britain. All rights reserved 0098-1354/96 $15.IX)+ 0.(XI

IDENTIFICATION A N D CONTROL OF A N A E R O B I C DIGESTERS USING ADAPTIVE, ON-LINE T R A I N E D N E U R A L NETWORKS C. EMMANOUILIDESand L. PE'rROU Aristotle University of Thessaloniki, GR-540 06 Thessaloniki, Greece (Received 19 January 1994; final revision receioed 9 June 1995)

Abstract--This paper introduces an anaerobic digestion identification and control scheme, based on adaptive, on-line trained neural networks. Anaerobic digestion is a complex, nonlinear biochemical process, widely used for the treatment of organic sludge in municipal wastewater treatment plants. Conventional control schemes usually fail to overcome the typical difficulties encountered in systems with complex nonlinear dynamics and difficult-to-measure or time varying parameters. It is shown by simulation results that, under a predictive control approach, adaptive on-line trained neural networks are successful in tackling such problems, in the case of anaerobic digestion. The proposed control scheme features desired tracking, regulation and robustness properties in various anaerobic digestion control tasks, including set points or process inputs variations, even in the presence of measurement noise or in cases of process parameter changes. In addition, the performance of three training algorithms, the back-propagation and two different random optimisation techniques, is examined over the neural controller training task. In all cases the random optimisation techniques converge much faster than the back-propagation algorithm. Copyright © 1996 Elsevier Science Ltd

I. INTRODUCTION

Anaerobic digestion is a complex, nonlinear, biochemical process, widely used in municipal wastewater treatment plants for the treatment of organic sludge. The process can be considered to take place in three stages. In the first one, the hydrolysis stage, the complex organics are decomposed to simpler ones. In the acidification stage, the acidogenic bacteria activity causes the conversion of the simpler organics to volatile acids (acetic, propionic, butiric), carbon didoxide and hydrogen. The final, most important and sensitive stage is the methanogenesis stage, where methanogenic bacteria consume the volatile acids, producing methane and carbon dioxide. The interaction of the carbon dioxide and the hydrogen produced in the second stage results in the production of some extra methane. The growth of the methanogenic bacteria is slower than the corresponding one of the acidogenic ones. In addition, methanogenic bacteria are very sensitive to sudden temperature changes and the presence of oxygen or toxic materials. However, in order to establish stable operation in an anaerobic digester, equal rates of volatile acids production and methanogenic bacteria growth must be achieved. In practice, both the observation and control of the anaerobic digestion have proved to be quite difficult tasks. The main reasons for that is the nonlinear nature of the process and the lack of reliable on-line measurement sensors, capable of monitoring some

of its basic variables and the variation of some parameters of the process through time. Such difficulties are usually met in many bioprocess systems. The use of neural networks appears to be of great value within the area of process identification and control. The ability to achieve accurate nonlinear mappings from input-output pairs of data, without knowing the functional relationship that relates them, is perhaps the most important of a series of features that neural nets share. In addition (Hunt et al., 1992), their learning and adaptation capacity, data fusion ability, applicability to multivariable systems and the convenience they present to parallel and hardware implementation make a strong framework, useful for identification and control purposes. Recent research have shown the applicability of neural networks to various process identification and control tasks. Bhat and McAvoy (1990) used the back-propagation algorithm to model the dynamic response of pH in a continuous stirred tank reactor (CSTR). Ydstie (1990) applied an indirect predictive neural controller to a reactor temperature control problem and developed a direct adaptive controller based on the identification of the system inverse. The effective use of neural networks to process modeling was also shown by Willis et al. (1991), where a neural network model was used to provide estimates of difficult-to-measure controlled variables, by inference from other, easily measured ones. They also proposed a neural network based 113

114

C. EMMANOUILIDESand L. PETROU

inferential control scheme to regulate the product composition in a pilot plant distillation column, Furthermore, they applied a neural predictive controller to a CSTR, in order to enable the reactor to operate at a near maximum feed conversion conditions, while preventing thermal runaway. Instead of the popular back propagation algorithm they used the chemotaxis algorithm, a kind of random optimization technique. Chemotaxis appeared to be superior to the back-propagation in terms of its simplicity, ability to avoid local minima and speed of convergence. The chemotaxis algorithm was also applied by Di Massimo et al. (1992) to provide on-line biomass and penicillin estimations in an industrial fermentation, while Pollard et al. (1992) investigated the capability of neural networks to provide accurate estimations of tray temperature in a distillation column. A common difficulty encountered when employing artificial neural networks to control biochemical processes is the selection of the initial training data. The training data set should span the whole range of operation in order to encompass information about all the dynamics of the process, i.e. the data set should be persistently exciting. Obtaining such a data set in the case of anaerobic sludge digestion from real measurement data is not realistic. Moreover, even if the adequacy of the initial training data set could be guaranteed, the fact that the anaerobic digesion is a process with time-varying parameters indicates the necessity of developing a process model capable of monitoring the parameter variations through time. Yang and Linkens (1994) have shown that by applying on-line adaptive neural networks to control processes with slow dynamics, the selection of the initial training data set becomes trivial and parameter variations through time can be successfully monitored. The aim of this paper is to show that effective anaerobic sludge digestion and control can be achieved by employing on-line trained adaptive neural networks under a predictive control strategy. This approach has some clear benefits compared to conventional control or non-adaptive neural network based schemes. First of all, it successfully overcomes difficulties arising from the nonlinear nature of the process, the lack of reliable on-line sensors and the fact that process parameters may change with time. On the other hand, the on-line estimation of physical process parameters, which is essential in a traditional parameter adaptive control approach (Polihronakis et al., 1993), is avoided, since a neural model can directly provide estimates of the process states. In addition, the presence of a persistently exciting initial training data set becomes

unnecessary, as the on-line adaptation of the neural networks adds flexibility to the overall control scheme. The predictive control scheme presented here, based on adaptive on-line trained neural networks, features desired tracking regulation and robustness properties, in various control tasks, as is shown by the extensive simulations carried out, covering a wide range of the process operation, even under the presence of measurement noise or parameter changes. All the simulations have been carried out using the model proposed by Buhr and Andrews (1977), which has been verified in the past by experimental data (Carr and O'Donnel, 1977; Polihronakis et al., 1993). The neural controller has been trained using three algorithms, the back propagation and two different optimisation techniques, the chemotaxis and a random search algorithm suggested by Solis and Wets (1981). The random optimisation techniques appeared to be superior to the back propagation both in terms of speed of convergence and computational simplicity.

2. NEURAL NETWORKS

A variety of neural network architectures and training algorithms are referred to in the literature (Hertz et al. 1991). A neural network consists of a structure of computing elements, called nodes or neurons. Each one of these nodes is a simplified model of the human brain cells. Neurons are interconnected with unidirectional synapses and each synapse has its weighting factor (synaptic weight or simply weight). Neural networks may have layered structure, while different interconnection schemes, such as intra-layer, inter-layer and recurrent connections can exist. The first neuron layer is called the input layer, while the last one is the output layer. Any other layer is called hidden layer. The behaviour of the neural network depends on its structure, and on the basis and activation function employed. Given a set of inputs and synaptic weights, the basis function determines the net value of the neuron. Usually a linear-basis function (LBF) or a radialbasis function (RBF) is employed. The activation function is applied to the net value to form the neuron output. The activation functions most commonly used are hard limiter, threshold, sigmoid and Gaussian functions. In this work only layered feedforward networks with linear basis function and sigmoid activation function are employed. Specifically, the activation function is: 1

g(h) - 1 + exp( - 2flh)'

(1)

Anaerobic digestors where h is the net value and fl a parameter that determines the slope of the activation. The determination of the appropriate weights that correspond to an accurate input-output mapping requires the training of the network. This training may involve the presence of a set of input-output pairs or patterns (supervised learning). The training set and the network architecture must be suitably chosen to ensure the desired generalisation. A variety of training algorithms can be found in the literature. Three of them are employed here and a brief discussion of them follows. This discussion is restricted to the main features of the employed variant, for each training algorithm.

Back propagation

Based on the gradient descent optimization, back propagation is probably the most popular training algorithm for feed forward networks. The basic back propagation algorithm has several drawbacks. The most important ones are slow convergence, the possibility to become stuck in local minima and the computational complexity. Many variations of the basic algorithm that improve its performance have been suggested. The use of a momentum term generally speeds up the convergence and smooths the trajectory of the weights during the update procedure. During training, both learning rate and momentum can be modified, in order to improve convergence. Chauvin (1989) suggested the use of an energy term in the cost function. It is argued that the presence of such a term can help convergence and avoidance of local minima. Chauvin indicated that, for energy functions of order higher than two, the behaviour of the algorithm does not significantly differ from an algorithm that employs a second order energy function. The algorithm suppresses the hidden nodes that have constant activation whatever pattern is inserted in the network and thus ~'optimizes" the use of hidden units. Applying gradient descent to this cost function, results in a delta rule that differs from the basic delta rule of the back propagation, in the presence of an extra energy delta for each hidden unit, whose energy is to be minimised.

decreases the cost function, then it becomes accepted and is repeated until it no longer causes an improvement to the cost function. This algorithm has clearly several attractive features, like programming and computational simplicity and flexibility to avoid local minima in some cases. However, if the training is based upon large number of patterns, then the evaluation of the cost function over the entire training set significantly slows down the procedure. It is also worth mentioning that, once a weight update becomes accepted, it is repeated until it no longer decreases the cost function, even when the improvement is too small. The algorithm can be modified so that a weight update would not be allowed to go on as long as it decreases the cost function, but only as long as the sequence of positive steps does not reach a maximum number of allowable steps. Thus the possibility of a repeated tiny positive step will be eliminated, with the cost of losing a highly positive step when it still can improve the cost function. This modification adds some flexibility to the basic chemotaxis algorithm and reduces the possibility of long training time, while, on the other hand, slightly slowing down the speed of the algorithm in some cases.

Random search

In contrast to chemotaxis, the random search optimization algorithm performs a random search at every iteration. Thus, the weight update procedure becomes more flexible. The algorithm proposed by Solis and Wets (1981) has computational and programming simplicity similar to those of chemotaxis. This similarity applies also to the remark that the larger the training set is, the more computational effort is needed to calculate the cost function at every step. The procedure of the algorithm is: 1. Choose W ( 0 ) • R". Set k = 0 , s = f = 0 and b(0) : 0 • R". Fix p( - 1), Plb, ex, ct, Sex, f,.,. 2. Set - p(k- l).ex p(k)=

I

p(k-1).ct

_ p ( k - 1)

Chemotaxis

The chemotaxis algorithm (Willis et al., 1991) empoys a random weight update procedure, instead of the deterministic one involved in back propagation. If a randomly generated weight update vector

115

ifs>~s~, iff>~L,

otherwise

(2)

if p(k)<~l~b then stop; otherwise generate a Gaussian distributed random vector AW(k) with mean value b(k) and standard deviation

[p(k)l.

116

C. EMMANOUILIDESand L. PETROU

training pattern set for the initial training of the anaerobic digestion model, which will be explicitly discussed later. The equations describing the dynamics of the two main variables of the process, the substrate concentration and the methane rate production, are:

3. Set

w(k + 1) = W(k) + AW(k)

if E{W(k) + AW(k)}

>e{w(k)}, sets=s+l,f=0and

~3T= - U ' S T - K" Y+ U'So y2 ? = (u - k~). Y - K . - -

b(k + 1) = 0.4. AW(k) +0.2.b(k)

ST

W(k) - zXW(k) if E{W(k) - AW(k)} +

< E{W(k)}, set s=s+l,[=O and b(k + 1) = b(k) -0.4.AW(k) W(k)

otherwise, s e t s = O , f = f + 1, b(k + 1) = b(k)/2

(4)

(3)

and return to step 2 with k = k + 1. where k is iteration number; thb is the lower bound for p(k). If p(k) reaches this value, after a long series of unsuccessful choices of AW(k), then it is useless to continue searching for better W(k), thus current W(k) can be considerd to be optimum, s and f are, respectively, the number of successive suecesses and failures in finding better W(k). If s, f reach respectively Sex and f¢, then p(k) is expanded (by a factor ex) or contracted (by a factor ct). W(k) is the network weights vector at the k-th iteration. E{. } is the cost function. 3. ANAEROBIC DIGESTION DYNAMICS

Among the mathematical models suggested for the anaerobic digestion process the one developed by Buhr and Andrews (1977) was selected for three main reasons: (a) the model has been verified by experimental data (Carr and O'Donnell, 1977; Polihronakis et al., 1993), (b) it is descriptive, since it provides equations for all the main variables of the process, (c) it is rather simple, as it considers only methanogenesis as the rate limiting phase. In this work, all the simulation tests for the anaerobic digestion process have been carried out using the Buhr and Andrews model. However, this model is too complicated for control purposes. Polihronakis et al. (1993) have derived a simplified model by considering a subsystem of the model originally proposed by Buhr and Andrews. The basic statespace equations describing the bacterial growth were modified to exclude the organism concentration, which, in practice, is an unmeasurable variable. This model is adopted here in order to obtain a

L\S

v + , .Xo

/

1. u,

K-V,,,s_l

(5)

where ST is total substrate (volatile acids) concentration, (mol/l, as equivalent acetic); U is influent dilution rate, (l/day); So is influent substrate concentration, (mol/l, as equivalent acetic); Y is methane production rate, (l/day); p is specific growth rate, (1/ day); kd is decay coefficient, (l/day); Yx/s is yield of organisms per mole of substrate consumed (Buhr and Andrews, 1977); X0 is influent organism concentration, (mol/I, as equivalent CsH7NO2); ST, ? are the derivatives of ST, Y, respectively; K is the yield coefficient, K = 1/(D. V. Yc.41x" Yx/s) (Polihronakis et al., 1993); D is gas volume conversion factor (1/ tool); Yell,Ix is yield of methane per mole of organisms formed (Buhr and Andrews, 1977); V is digester liquid volume (1). The control scheme proposed in this work, considers the influent dilution rate (U) as the manipulated input, expressed in terms of sludge volume to digester volume per time unit, or equivalently as (l/day). The controlled variables are the total substrate concentration (ST) and the methane production rate (Y). The influent substrate concentration (S,) and the influent organism concentration (X0) are considered as external disturbances. The influent substrate concentration (So) is measured indirectly and ()to) is usually considered equal to zero. However, in cases of digested sludge recirculation, X, should be taken into account, thus a nonzero value for X0 is considered in this work. By employing first order Euler approximation, equations (4) and (5) are discretised as follows: S.r(n + 1) = Sv(n) - H. {U(n). [S0(n) - ST(n)] -- K ( n ) .

Y(n)} + e,(n)

(6)

Y(n + 1) = Y(n) y2(n)'~

+/4. f~[l, - k,]. Y(,,) - K(,,). S+--~S~)J +

H ' ~ [ S°(n)

[ [ST(n)

-

x U(n) + e2(n),

2] Y •

I~(n)'Xo(n)) K(n). Yx/s J

(n)+---

(7)

117

Anaerobic digestors where (n) is the discrete time index (n = 0, 1, 2 . . . . ), H the sampling interval, and el(n) and e2(n) represent the cumulative effect of noise, modeling and discretisation errors. Considering that/~(n) can be determined by the current values of pH and temperature (0), equations (6) and (7) indicate that for a given set of measurements for ST(n), S0(n), Y(n), pH(n), O(n) and U(n) (which is the control input), a prediction for the future values of ST(n + 1) and Y(n+ 1) can be derived. Instead of using a parameter estimation algorithm for the parameters/~ and K, an attractive alternative is to train a neural network to provide accurate predictions of Y(n + 1) and Sv(n + 1), when ST(/'/), So(n), Y(n), pH(n), O(n) and U(n) are applied at its input layer, Then, on the basis of the model-network predictions, a predictive control strategy for the anaerobic digester can be implemented. 4. ANAEROBIC DIGESTION CONTROLLER

4.1. General description In many cases of nonlinear dynamical systems control, the availability of a neural model of the system can ease the burden of designing an effective control scheme. Such a model can be directly employed in a model based control strategy. Among the control strategies suggested in the literature, the predictive control scheme seems to offer stable solutions for nonlinear systems. In a typical neural network predictive control strategy, a neural model is used to provide predictions of the future plant response over a specified predictive horizon. Then an optimisation routine utilises these predictions in order to minimise a cost function, usually according to the generalised minimum variance control law. It

is crucial in this procedure to select the right length of the predictive horizon as well as to balance the importance of minimising the control action variance together with the variance of the controlled variables from their set points. Another important point of consideration, when predictive control is employed, arises from the fact that the optimisation criterion is not the minimisation of the difference between the system output and the set point, but the difference between the neural model output and the set point. When an accurate neural model is available the neural controller derived is very close to the inverse of the system. However, when there is a considerable mismatch between the real system and the neural model the overall control performance is poor. The neural model accuracy can be kept at acceptable level when adaptive on-line trained neural model networks are employed. In particular, Yang and Linkens (1994) employed a moving-time-zone, fixed length training data set in order to achieve on-line adaptation of the neural model. Thus, in systems with slow dynamics, the adaptive neural model can successfully track the variations of the controlled system parameters. The same on-line training procedure for the neural model is followed here in the anaerobic digestion control scheme shown in Fig. 1. Given a set of measurements X of the process input and output variables, the predictive accuracy of the neural model is checked. If necessary, the online model training algorithm is activated to improve the neural model predictive performance. Then, the controller network, which receives at its input the measured values X and the set point values Ysp, STsp of Y, Sv respectively, offers at its output a suggested value for the manipulated input U, i.e. the influent dilution rate. Based on these values of X, U, the

'l

Plant + ..............

X

/ /"

/

sT,p

_T

Con~oHex u

II

Model N~,ml

Network /' !"

Network ,

Fig. 1. Block diagram of the on-line adaptive neural network predictive controller.

ll8

C. EMMANOUILIDESand L. PETROU

neural model offers an estimate of the future values 1;', ~T of ]I, ST, respectively. If these values are not close enough to the set point values, the controller is adequately trained, until it provides an acceptable value for U. The objective of the proposed controller is to ensure that, regardless of input, load or environmental disturbances, the controlled variables remain close to their prespecified set point values, even under the presence of measurement noise or anaerobic digestion parameter variations, while, at the same time, the stable operation of the digester is secured. The controller presented here can regulate either one of the process outputs (i.e. the methane production rate, Y and the substrate concentration, ST) separately, or both of them, in a unified control approach (combined control strategy).

4.2. Neural modeling of the anaerobic digestion The anaerobic digestion predictive control algorithm presented here requires that an initially trained fairly accurate neural model of the process exists. This model is dervied using available a priori information for the anaerobic digester. The source of this information may be a data set of real life measurements from a well constructed series of experiments or a pattern set directly derived from an existing mathematical model of the process. Even though the former set may not contain all the necessary information of the process dynamics or may be contaminated with noise and the latter set may be inaccurate due to modeling or discretisation errors, the involvement of on-line neural model adaptation in the control procedure, can improve the predictive accuracy of the model. In other words, it is unnecessary that the a priori information should consist of a persistent exciting data set. In order to proceed with the development of a neural model for the anaerobic digestion process, some practical difficulties concerning the procedure of measurement acquisition, should be stated. Since the set of measurements required, involves time consuming laboratory analysis, a sampling interval of less than a few hours is not realistic. Thus, a sampling interval H = 6 h is chosen. For such a sampling interval, the predictions of the methane production rate and substrate concentration, obtained from equations (6) and (7), are not accurate enough, since the assumption that the higher order derivatives can be ignored is invalid. Furthermore, the large number of parameters involved in the process dynamics, makes the neural network initial training task rather tedious if a neural predictor with a wide predictive horizon should be developed. Thus it can be reasonably

questioned, whether the development of a neural model for the process with wide predictive horizon, could be cost efficient. Instead, it is considered more preferable to utilize a neural model which acts as a one step ahead predictor and employ on-line training to further improve it's predictive accuracy. The initial training of the neural model for the anaerobic digestion process was based on a training set of 48,000 patterns of input--output pairs. Each pattern consists of the current values of the input variables, i.e. the temperature {0(n)}, total substrate concentration {ST(n)}, methane production rate {r(n)}, influent substrate concentration {S0(n)} and pH {pH(n)}, together with the value of the control variable, i,e. the influent dilution rate {U(n)}, as well as the values of the output variables, i.e. the future values ST(n + 1) and Y(n + 1), obtained from the simplified model (equations (6) and (7). All the values involved in neural network calculations throughout this paper have been normalized according to: X

xN = rain + (max - min)

-

-

-

Xmi n

-

(8)

Xma x -- Xmi n

where XN, X are the normalized and non-normalized values, respectively, min and max are selected to be 0.1 and 0.9, respectively, for the output variables and - 1 , 1, respectively, for the input variables, while Xm~, and Xmax are the extreme values corresponding to each one of the process variables. The extreme values for all process variables are shown in Table 1. In the pattern generation procedure the sampling interval was taken equal to 0.25 days, X,=0.001 (mol/I) for a 101 digester volume. The normalized measured values are inputs to the network and the normalized future values of Y and ST correspond to the output patterns. The activation function employed is the one given by (1) with fl = 0.3. The basic back propagation algorithm was modified to allow the adjustment of the learning rate and momentum in order to improve convergence. A quadratic cost function with the addition of an energy term was selected:

E=/~:,.E:r+g:,.E:.,

(9)

where E is the cost function, Ecr is the standard quadratic cost function of the back propagation algorithm, i.e. the summed squared error over the Tablc 1. Minimumand maximumvaluesof thc processvariables ST Y (mol/I) (I/day) pH Maximum 0.035 Minimum 0.0005

30.0 0.5

6.9 6.0

U St, (I/day) (tool/I) 0 (*C) 0.2 0.01

1.5 0.05

39.0 35.0

Anaerobic digestors whole set of patterns, E~. is the sum of the energy "spent" by the neurons over the whole set of patterns and/~er and/~e, weighting factors. The weighting factor/~¢r has been set to 1, while/~¢, has been initialized at 0.05 and is gradually reduced during the training procedure in order to improve convergence. The energy function, (e,), of the neuron activation was: 0 2 %(02) = l + 02,

(10)

where O is the neuron output. After several training experiments a network with two hidden layers and six neurons per hidden layer appeared to be a good compromise betwen the predictive accuracy achieved and the network size. The patterns have been randomly selected, so that the set of the input variables follows a multivariate uniform distribution with extreme values, the ones shown in Table 1. Both the stable operation region of the digester and a large region of input variables combinations that lead to system failure are included. Note that following this procedure, the input set of each pattern, is not necessarily compatible, i.e. it does not correspond to a realistic stage of an anaerobic digester. This is due to the fact that there is a strong interdependence between the input variables. However, the possibility of selecting a pattern, that corresponds to a non-realistic state, can be reduced. This can be achieved, by rejecting combinations of input variables leading to either Y(n+ 1) or ST(n+ 1) values outside the region of Table 1 or causing an extremely large change in ST and Y. The latter condition can more precisely be expressed as follows: -18.0
(11) (12)

This constraint has been selected by simulation tests corresponding to extreme load or environmental changes. In particular, even in cases of extreme change which may appear in a real anaerobic digester, the resulting changes in Y and Sx fall inside the area described by equations (11) and (12). Thus, such a constraint should be considered more as a result of the interdependence between the process variables than as an arbitrary assumption, Moreover, the presence of "outliers" in the pattern set can be a serious obstacle to the training procedure. In order to suspend the inhibitory effect of these "outliers" a robust error suppressor may be used (Kosko, 1992). In the present training, after 6,000,000 iterations, an error suppressor of the form:

119 a°e

e' =

(13)

(b.e) 2 + 1

was employed, where e' is the suppressed and e the non-suppressed error, between the normalised values of the desired output and the corresponding network output. The choice of the values of the parameters a, b that appear in equation (13) depends on the particular training problem and is made upon the basis of physical limitations on the real process. The complete training procedure required approximately 10,000,000 iterations. The neural model has been tested by several simulation tests, covering a wide operating region of the digester. The initial values of the variables used in the simulations are (ST=0.004495 mol/l as acetic), (Y--- 7.434232 l/day), (U = 0.066 I/day), (S0=0.43mol/i as acetic), (0=38.0°C) and (pH = 6.679287). The predictive performance of the network has been examined by carrying out six different simulation tests, covering various environmental and load conditions. In all cases the Buhr and Andrews (1977) model has been considered as the "real plant". The step changes imposed in the input variables are briefly presented and commented on in Table 2. Case A is shown in Fig. 2. The network predictions of both ST and Y are quite accurate. Figure 3 shows cases B and C. In these cases of extreme shocks, larger prediction errors appear. However, the neural model predictive performance is improved when the on-line model training is applied, as will be shown in the control simulations presented in the following sections. Cases B and D are shown in Fig. 4. Finally, for cases E and F similar predictive performance of the network is obtained, as shown in Fig. 5. A moving-time-zone fixed-length training data set is employed for on-line neural model training (Yang and Linkens, 1994). The data set at the n-th time instant has the form: Din(n) = {[Xm(n -- l m ) , y m ( n -- l m ) ] ,

[Xm(n - l m + 1), ym(n -- lm +

1)l . . . .

[Xm(n), ym(n)]},

(14)

where Xm(n) = {Y(n - 1), ST(n -- 1), pH(n - 1), U(n-1),So(n-1), 0(n-l)}

(15)

ym(n) -----{Y(n), ST(n)},

(16)

and

120

C. EMMANOUILIDES a n d L. PETROU Table 2. Verification tests for the neural model

Case

Changes in the input variables So = 0.6 mol/l U = 0.0924 I/day 0 = 36.0°C U = 0 . 0 6 6 I/day S~= 0.43 mol/l 0=38.0°C

A series of step changes in both So and U combined with sudden temperature changes.

2

An extreme case of organic overload; the digester is driven very close to failure, but the system remains stable.

3, 4

t = 2 days: t = 6 days: t = 10 days: t = 14 days: t = 18 days: t = 2 2 days:

B

t = 2 days:

C

t = 2 days: S0= 1.35 mol/I

Similar case to B; the minimum organic shock that causes system failure.

3

D

t = 2 days: So = 1.3 mol/I t = 3 days: 0 = 36.0°C

Same case as B with the addition of a sudden drop in temperature that leads to system failure.

4

E

t = 2 days: U = 0 . 1 8 I/day

A hydraulic overload; system remains stable.

5

F

t = 2 days: U = 0 . 1 8 7 I/day

Similar case as E; the minimum hydraulic shock that causes system failure.

5

So = 1.3 mol/I

1 tm Em = ~ ~ {[Y(n - k) - ~'(n - k)] z + [ST(n -- k ) - $T(n -- k)]2},

(17)

where ~(n), ~T(n) are the neural model predictions when the input vector Xm(n) is fed to the input layer and index m refers to the model neural network. The back propagation algorithm with adaptive learning rate and momentum is employed in the minimisation of the above cost function. To reduce the effect of noisy measurements and improve the robustness of the neural model Yang and Linkens (1994) suggested that the weight adaptation procedure should be carried out according to: 1) =

Fig.

A

while lm is the fixed data length of the training data set Din. The cost function to be minimised is:

W m ( n 4"

Comments

am" W i n ( n ) -k (1 - a m ) "

W ' ( n + 1), (18)

where Win(n) is the neural model network weight vector at the n-th time instant, W ' ( n ) the neural model weight vector after the minimisation of equation (17) is completed and am is the model adaptation rate. The effect of the parameters/.1 and am is important for the on-line training of the neural model. A suitable value of lm should be selected, in order to ensure that the information included in the training data set Dm is rich enough for the neural model training. However, if an excessively large value Oflm is selected, then the neural model may be overtrained. In practice, the time zone should be selected on the basis that Dm contains all the information that is necessary for the description of the

process response. In particular, Yang and Linkens (1994) suggested that lm should be proportional to the primary time constant of the system. However, there is no existing rule for the precise determination of Ira. For the anaerobic digestion controller presented in this paper, a value lm= 40 that corresponds to a time zone of 10 days is found a suitable selection, since this time period encompasses the typical behaviour cycle of the process dynamics. On the other hand, the selection of am should be made with criterion the predictive performance of the neural model. A small am is more suitable for a system with relatively rapid dynamics, while a larger value can increase the model robustness and reduce the effect of noise. A value of 0.2 was found suitable for the anaerobic digestion neural model training and is kept fixed during all the control simulations carried out throughout this paper. Finally, a maximum allowable number of 2000 iterations is imposed for the on-line neural model training to ensure that the time consumed should not become too large for real-time calculation. 4.3. C o n t r o l l e r d e s c r i p t i o n

The controller neural network shown in Fig. 1 represents the inverse of the system, in the sense that when the input layer is fed at the n-th time instant with the measured values X(n), i.e. the methane production rate {Y(n)} the substrate concentration {ST(n)}, the pH {pH(n)}, the temperature {0(n)}, the total input substrate concentration {Son)}, as well as the set points Ysp(n), STep(n) for Y and ST,

Anaerobic digestors

121

"~ 0.010 -

Plant ..........Prediction -

0.008 O

0.006 0

0.004

0.002 0q

T i m e (Days) 0.000

t

I

I

t

0

I

I

I

I

I

5

I

t

I

I

I

I

I

f

I

'

I

,

,

20

15

10

I

'

l

,

I

25

30

(a) 0~

16 -

Plant

-

.......... Prediction

0

....... /F

0

4

I

0

I

5

I

I

I

F

I

r

~

I0

I

r

]

15

I

I

F

I

I

20

I

I

Time

(Days)

I

I

I

I

25

I

I

I

30

(b)

Fig. 2. Predicted response of the process with a series of step changes in So and U combined with sudden temperature changes (case A). (a) Substrate concentration response. (b) Methane production rate response.

respectively, it provides at its output a suggested value for the manipulated input, i.e. the influent dilution rate. This value, if it becomes accepted, should drive the system outputs to the set point levels. The neural model receives this value together with X(n) and offers the predicted values ~'(n + 1) and ~v(n + 1). If these values are not close enough to the prespecified set points, the learning algorithm adjusts the controller network weights in order to minimize a prespecified cost function. A maximum allowable number of 500 iterations is also imposed for the on-line controller training in order to compromise with real-time calculation limitations.

It should be noted that, since the control scheme acts as a one-step-ahead predictive controller, the above minimisation procedure can result in excessive control effort and large overshoots. However, this can be easily avoided, as will be shown later, either by using the generalized minimum variance as the cost function to be minimized, or by employing a simple low pass filter to smooth the control action suggested by the controller network. In the former case, a tuning factor has to be suitably selected to determine the importance of the control effort minimisation, while in the latter case it is a matter of appropriate filter coefficients selection. For the case

122

C. EMMANOUILIDESand L. PETROU 0.035 Plant .......... Prediction

"-" 0.028 e. O .,=~

,~ 0.021

,,,o

~ 0.014

0.007 Time 0.000

I

,

i

i

I

'

'

'

I

'

5

0

'

r

r

I

I

10

t

(Days) ,

,

,

15

20

(a) 30

Plant ...........Prediction -

......................................................

-

I---,

,

¢¢

D

20

~ 10-

Time ~

0

I

0

r

I

I

i

I

i

5

I

I

10

I

I

I

I

I

15

i

(Days) I

I

I

20

(b)

Fig. 3. Predicted response of the process with organic overload (cases B and C). (a) Substrate concentration response. (b) Methane production rate response. of training the controller by the back propagation algorithm, the use of such filter is further justified, as will be explained later in this section. Provided that the neural model is accurate enough, it is unnecessary to have an initially trained controller neural network, since the controller training is goal directed, i.e. it will be completed when the controller output has acquired a value Uc(n) that, according to model predictions, will drive the plant response to the set point level. This observation simply means that the availability of a prop e d y initially trained controller network is not a prerequisite for the control procedure. However, when a priori information is available, it is worth developing such a model, as the on-line controller training will become much easier. In particular, a

data set obtained in exactly the same way with the set employed for the neural model initial training can be used for the initial controller training. The difference here is that the future values of Y(n + 1), ST(n + 1) are considered as the set point values of Y, ST respectively, while U(n) is considered as the desired controller network output that would drive the system to the set point values. Thus, each pattern consists of the current values of the input variables, i.e. the temperature {0(n)}, total substrate concentration {ST(n)}, methane production rate (Y(n)}, influent substrate concentration {So(n)} and pH {pH(n)}, together with the future values Sx(n + 1) and Y(n + 1), obtained from the simplified model (equations 6) and (7), as well as the value of the control variable, i.e. the influent dilution rate

123

Anaerobic digestors *~ 0.035

,.•

Plant

0.030 0.025

0.020 0.015

3 v 0.010

0.005 0.000

1

Time I

I

I

0

I

I

I

I

I

5

I

I

I

I

I

10

I

I

(Days) t

I

I

15

20

(a) 30 13

:""""""("',

t~

......................................................... ,ol

20

.2 0

0

~ lO

_

~

-

Plant

.......... P r e d i c t i o n

q) ~

Time (Days) 0

I

0

i

I

I

I

I

I

5

I

I

10

I

I

I

I

I

15

l

I

i

i

20

(b) Fig. 4. Predicted response of the process with organic overload followed by a sudden drop in temperature (cases B and D). (a) Substrate concentration response. (b) Methane production rate response. {U(n)}, which is now considered as the network output. Based on such a data set, consisting of 48,000 patterns, a controller network of two hidden layers, with seven nodes in the first hidden layer and three in the second and employing the activation function of equation (1) with f l = 0 . 5 , has been trained using the same variant of the backpropagation algorithm as the one empoyed for the initial neural model training. However, since the controller network has only one output, the training procedure was shorter and completed within 5,000,000 iterations, including approximately 1,000,000 iterations with the robust suppressor of equation (13) employed. Apart from the on-line controller training described above, the neural

model on-line training is activated, when the overall predictive performance of the model network over the moving-time zone fixed-data-length model training set Dm becomes poor. In practice, the permitted change in the control input strongly depends on the preceding stage of the anaerobic digester, i.e. the sludge thickener. In the simulations performed, this dependence has been taken into account as a constraint. Thus, per each sampling interval, changes in the control input higher than 15% of its current value are not allowed. Without such a constraint, improved tracking and regulation could be achieved. However, this would be of no real value, since it could not be implemented in a real digester.

124

C. EMMANOUILIDESand L. PETROU 0.05 Plant .......... P r e d i e t i o n

- -

"-" 0.04 o

~ o.oa

o 0.02 Co

0.01 Time

0.00

i

i

i

0

I

i

i

i

i

5

I

i

I

i

I

(Days) I

I

r

15

10

20

(a) "O

3O

2O, F

~ 10-

Plant ..........Prediction -

Time ~

0

I

0

1

I

I

i

I

i

5

[

I

i

i

i

i

10

I

15

I

(Days) I

i

20

(b) Fig. 5. Predicted response of the process with hydraulic overload (cases E and F). (a) Substrate concentration response. (b) Methane production rate response.

Several different control tasks, including set point changes and disturbances have been carried out, even under the presence of measurement noise or process parameters changes. Cases that lead the open loop system to failure have also been considered. Each one of these tasks was carried out using three different training algorithms for the online controller training, the back propagation, the chemotaxis and the random search optimization algorithms described in Section 2. The back propagation algorithm has been applied without energy term in the cost function, but with adaptive learning rate and momentum. In order to reduce the control effort the generalized minimum variance law has

been used, which in the present case suggests that the following cost function is to be minimized: ny

Ec(n) = 7 [Ysp(n + 1) - ~'(n + 1)] 2 ns

+ ~- [ST,(n + 1) -- ~T (n + 1)]2 q + --" [Uc(n) - U ( n - 1)] 2, 2

(19)

where n~ and n~ are weighting coefficients that allow the output variables to be separately or jointly controlled (combined control strategy), q is a scalar factor, U~(n) is suggested by the controller value for

Anaerobic digestors the influent dilution rate and index c refers to the controller neural network. When no filtering is applied to the value that the neural controller finally suggests as the appropriate control action, the adopted influent dilution rate value at the n-th time instant is U ( n ) = Uc(n). In order to employ the back propagation algorithm to the controller training, the gradient of the cost function with respect to the control input must be available. In the case of the squared error function, this gradient is calculated by propagating the error through the neural model. However the neural network controller can not be trained by the back propagation algorithm when the cost function of equation (19) is employed, with q#:0, since the gradient of this cost function with respect to the control input is not available. The choice of the squared error cost function (q = 0) allows the calculation of this gradient. In this case, a low pass filter smooths the control action suggested by the controller. A suitable choice was found to be: U(n) = f ( 0 ) . Uc(n) +f(1)" U(n - 1),

(20)

where f(0) and f(1) are equal to 0.2 and 0.8, respectively, for the methane production rate control case, and 0.25 and 0.75, respectively, for the substrate concentration control scheme. A maximum allowable value of 16 (I/day) is imposed for the methane production rate set point. This choice is justified by the fact that the digester operation at higher rates of methane production can easily turn out to be unstable, even with small input disturbances. When Y approaches 16 (I/day) the control is gradually directed towards the methane production rate, by modifying the weighting factors ny and ns, on the basis of the prediction ~'(n + 1): 1

ny(n + 1) = 1 + exp[ - 8. (lZ(n + 1) - 16)] (21) n~(n + 1)= 1 -ny(n + 1).

(22)

The values for ny and n~ given by equations (21) and (22) become accepted only if ny is greater than its prespecified value. Since the substrate concentration and the methane rate production are not independent variables, the selection of the set points Ysp and STso cannot be arbitrary. Thus, when combined control strategy is applied, instead of choosing arbitrarily the set points, Polihronakis (1992) suggested an attractive alternative. According to this, an initial selection of ST~o is made. Considering the steady state of the anaerobic digestion, an indicative relation between Y and ST is derived by equations (4) and (5).

Uss

y,, = ~ . [S,~ - ST~,],

(23)

125

where index ss denotes steady state values and for a 101 digester: K~-~ 11.748/(273 + 0),

(24)

according to the Buhr and Andrews (1977) model. Applying (23), the set point for Y is derived according to the equation:

us~

Ysp(n) = K~--~" (So(n) - Svsp(n)),

(25)

where U~ is the desired value of U for the steady state operation of the digester, which is taken as Us~=0.066 (1/day), while K~s is given by equation (24). When the automatic set point selection for Y is employed, for values of Ysp lower than 16 (I/day) the value given by equation (25) becomes accepted, otherwise y,p is set to 16. 4.4. Simulation tests All simulation tests have been performed considering the cost function of equation (19) either with q g:0 (generalized minimum variance law), or with q = 0 (simple quadratic criterion). In the former case the tuning factor q is regulated to a value equal to 1, when the control main target is towards the methane production rate, or to 0.25, when the control target is driven towards the substrate concentration. These values for q have been selected by a trial and error procedure, in order to ensure stability, without vigorous changes in the control action. Thus upsets to the thickeners can be reduced. They also seem to offer an improved plant response, since higher values of q result in smaller rise time, while lower q values invoke an overdamped response. The generalized minimum variance law involves higher training time in exchange of improved control effort for the same mean error values. In all cases the training algorithms have successfully converged, providing the influent dilution rate value. However, both the chemotaxis and the random search algorithms have been found to converge much faster than the back propagation. The simulations have been carried out using the Buhr and Andrews (1977) model as the "real" plant. The anaerobic digester is considered to start up with the initial conditions mentioned in the modeling section, while the initial set point values are (Sv~p= 0.004495 tool/l, as acetic) and (y,p = 7.434232 l/day). A list of the control task simulations performed, is given below. All tests have been carried out for all three training algorithms. However, no significant differences appear between the plots obtained using the three different training algorithms, with the exemption of the slightly reduced control effort spent, when the generalised minimum variance control law

126

C. EMMANOUILIDES and L. PETROU

is employed. This improved control action is achieved at the expense of an increase at the required CPU controller training time. Thus, only the random search training technique results, with the generalised minimum variance control law employed, are shown in the following. A. A square wave point change followed by step input changes: t = 2 days: 20% increase in I(v, S~p. t = 40 days: Y~p,S~p are set back to initial levels. t = 8 3 days: S0=0.6 moi/l. t = 139 days: So = 0.43 moi/l. t = 195 days: 0 = 35.5°C. t = 219 days: 0 = 38.0°C.

Different pairs of weighting coefficients ny and n~ are chosen in the subcases considered: (ny = 0, n~= 1), (ny= 0.02, n~=0.98), (ny=0.1, n~= 0.9), ( n y = l , n~=0), (ny=0.65, n~=0.35). The substrate concentration and methane production rate response and prediction plots corresponding to the simulation tests where the main control target is ST (ns>ny) are shown in Fig. 6(a-d). It can be seen that the substrate concentration (ST) follows very closely the set point in the case that ny = 0, whereas in the other two subcases, where ny becomes significant, the steady state error of (St) is increased, but the (Y) error is now reduced. In all cases the model predictions are very accurate, owing to the on-

0.006 Set Point ..........ny=O, no=1 ny=O.02, nt=0.98 ny=O.l, n°=0.9

- -

0

0.005

!

o 0.004 fl

Time

(Days)

m 0.003 0

25

50

75

100

125

150

175

200

225

250

(a) 0.0060

"~ ~ "-4

- Plant .......... P r e d i c t i o n (ny=0,

n,=l)

:~ 0.0050

0.0040 .,a

Time

(Days)

200

225

o~ 0.0030 0

25

50

75

100

125

150

175

250

(b) Fig. 6. A square wave set point change followed by step changes in S0 and 0. Substrate concentration control. (a) Substrate concentration response. (b) Substrate concentration predicted response. (c) Methane production rate response. (d) Methane production rate predicted response. (e) Control action. (f) pH response. (g) Variation of the specific growth rate.

Anaerobic digestors

127

13-Set P o i n t ..........nr=O, n,=l

~=~ v

nr=O.02, n,=0.98

ny=0.1, n,=0.9 /-'"'"-'"'i

@ o

9

o

t

,','

r

~ "k..

~,,

7

Time

~

(Days)

5 25

50

75

100

125

150

175

200

225

250

(e) 13 Plant .......... P r e d i c t i o n (ny=O, n , = l )

-

.t=J

-~ 11

@ .#J

0

9

0 7"

Time (Days) 5

Illill~lllllllilliltilllllltllllllliillililllIIl~ 25 50 75 100 125 150

175

200

225

250

(d) O. lO t~

^

.......... ,n,,,=- uO: O 2 . 'nn _=z-__0 . 9 8

eA

"~ 0.08 @

o.o8 e~

Time (Days) 0.04

, , , i i ,i

25

J i i~i

50

,II,i

75

,Ill

i i l l I I , i l l

lO0

1.25

Fig. 6

(continued).

(e)

l I I i l

150

175

iII

Jill,

200

I i i II

225

250

128

C. EMMANOUILIDESand L. PETROU 6.72 ..........

ny=0, no=l ny=0.02, n . = 0 . 9 8 ny=0.1, no=0.9

6.68

6.64

Time 6.60

fliilrl 0

J II 25

50

Illl

II

75

II

II

II

100

IllllJl

IIIIII

125

150

175

(Days)

I I JIIII

I I'i

200

225

I II

250

(f) ~ , 0.13 ny=0, n , = l .......... ny=0.02, n,=0.98 . . . . . ny=0.1, n,=0.9 0.11

-~ 0.09 0

¢)

t t

~, 0.07 ~U

Time 0.05

IIrlll 0

25

irll]r]liillllii[l~tllilil[lilirllltr

50

75

t00

(Days) till

125

150

175

200

225

250

(g) Fig. 6 (continued).

line model training. The mean value of the absolute prediction error is kept below 0.3% for both Y and ST. The influent dilution rate (U), pH and specific growth rate (/~), corresponding to each one of the ST targeting control subcases, are illustrated in Fig. 6(e-g). Similar results were obtained for the next two subcases, where the main control target is Y. These results are illustrated in Fig. 7. B. An extreme case of step input changes that would drive the open-loop system to failure: t = 3 days: S, = 1.3 mol/l. t = 7 days: 0=36.6°C. t = 7 1 days: 0=38.0°C. t-- 103 days: S0 -- 0.43 mol/1.

Sx~pis set to 0.004495 moi/l, while Y~pis automatically set according to equation (25). The following subcases are considered: (ny=0, n~= 1), (ny=0.02, n~= 0.98). The simulation plots of case B are shown in Fig. 8. As long as (Y) remains below 16 l/day, (ST) remains very close to the set point. When Y approaches 16 i/day the control is gradually directed towards Y, preventing any further increase in (U). Thus the digester is driven to a lower (St) steady state value and the plant remains stable. When (S0) is set back to the initial level the digester settles back to the normal operation. The modifications introduced by equations (21) and (22) ensure that no abrupt changes in U occur, when (Y) is around 16 (l/day). Note the

Anaerobic digestors

where k=O, 1,2 . . . . . and A0(k), AS0(k) are Gaussian distributed random variables with zero mean and standard deviation 0.04°C for the former and 1.5% of the current So value for the latter. Two subcases are considered:

initially poor predictive performance when So is set back to 0.43 i/day. However, the improvement in the predictive accuracy, due to the online model training involved, becomes obvious very soon. The mean value of the absolute prediction error is kept below 0.65, and 2.3% for Y and ST respectively. C. A more realistic case; a wide operating region of the digester is covered. A series of random disturbances A0(k), AS0(k) are imposed on 0 and So, respectively, every 6 h, according to:

O(n+l)=O(n)+AO(n)

(ny = 0, ns = 1), (ny = 1, n~ = 0). The set points changes are: t = 2 days: 20% increase in Y~p, S~p. t = 40 days: S~p is set back to the initial level. t = 120 days: Y~p is set back to the initial level.

(26)

So(n+ 1)=So(n)+ASo(n),

129

Case C is illustrated by Figs 9-11. The random disturbances imposed on (So) and (0) are shown

(27)

c¢11- -

Set Point ..........ny=l, n,=0 ..... ny=0.65, n°=0.35

v

90

L

t

..........

o 7-

e, o~

Time

~D

(Days)

5-

0

25

50

75

100

125

150

175

200

225

250

(a) II -

Plant ..........Prediction -

(ny= 1, n,=O)

o O

k

o h gz. O~

,= ¢0

25

50

75

100

125

150

175

Time

(Days)

200

225

250

(b) Fig. 7. A square wave set point change followed by step changes in So and 0. Methane production rate control. (a) Methane production rate response. (b) Methane production rate predicted response. (c) Substrate concentration response. (d) Substrate concentration predicted response. (e) Control action. (f) pH response. (g) Variation of the specific growth rate.

C. EMMANOUILIDESand L. PETROU

130

J.q

0.0065

-

o

Set Point .......... n y = l , n.=O ..... n~=0.65, n,=0.35

o

0.0055

.,J l

-'-.

0.0045 ,,'"

U/

,Y ,,x

0.0035 CD

I

Time

(Days)

5~ 0 . 0 0 2 5 25

50

75

t00

125

150

175

200

225

250

(c) 0.0065 -

-

..........

^ o

Plant Prediction (ny= 1, n.=O)

0.0055 t,, 0.0045

0

2

o.oo35

~q

u~ 0 . 0 0 2 5

I I i I I i I i

25

I i

50

i I i I i

75

i

I ] t I i

i i I i I I I i i

I00

125

150

~ J I I I I i

175

I I i

200

t I I i

F i I

225

250

(d) ~-. 0 . 0 9 ny=l, n,=O .......... n y = 0 . 6 5 , n o = 0 . 3 5 - ~ 0.08

1

0.07 o

2

..

_

L

0.06

,,,o

e, 0.05 Time O.O4

(Days)

I I I I I r I ~ J I I I I I t t I I I t I 3 I r t i I t I I I I r { I r t I [ I I I r {

25

50

75

IO0

Fig. 7

I25

(continued).

t50

175

200

225

r I r i

250

Anaerobic digestors

131

6.70 n~,= 1, ..........

nm=O

ny=0.65, n,=0.35

6.68

=: :,, 6.66

6.64 Time 6.62

I I I I I I'1

25

(Days)

~ I I I I I I I I I ( I I I I I I I ~ I I I i I I I I I I I I I I I I I I I I I I I

50

75

tO0

125

150

175

200

225

250

(f) 0.13 .q

~0

ny=l,

I)

n.=O

.......... ny=0.65, n,=0.35

I

O. l l

-~ 0.09 0L , L~ ¢.)

;.7, 0.07

.F,,t

0.. O"2

Time 0.05

(Days)

lllilillillll;ll;l~llllllllllIlllllllilllllllllll

25

50

75

I00

125

150

175

200

225

250

(g) Fig. 7 (continued).

in Fig. 9. In Figs 10(a) and (b), and l l ( a ) and (b) the responses for the substrate concentration and the methane production rate, as well as the model predictions are shown, respectively. The controlled variables remain very close to the set points and the mean value of the absolute steady state error remains well below 1.7% in both subcases. The predictive performance of the neural model is also excellent, since the mean value of the absolute prediction error is kept below 0.35% for both Y and ST. In Figs 10(c) and 11(c) the corresponding control actions are given. D. Process parameters change through time together with the series of events of case A. The following parameter changes are introduced:

t = 0 days:

the specific growth rate (H) specified in the Buhr and Andrews (1977) model, is scaled by 1.1, while the yield coefficient (K) is scaled by 0.9. t = 70 days: the above scaling of g reverts to unity. t = 120 days: scaling of K is also set back to unity. t = 170 days: K is scaled by 1.1 and H by 0.9. The proposed scheme appear to be effective in controlling the anaerobic digester operation when significant process parameter changes through the time occur, as is illustrated by Figs 12 and 13. The process parameter changes are

132

C . EMMANOUILIDES

well monitored, due to the on-line neural model adaptation. Within a period of less than 10 days after the parameter changes, the model predictions become very accurate and the controller is capable of providing the appropriate control action. It should be noted that the parameter changes introduced are very abrupt and most unlikely to be met in a real digester, where slow parameter variations through time may appear. Since the proposed control scheme can track so significant variations of the process parameters in a relatively short time period, it can be reasonably argued that more realistic cases of slow parameter variations can also be tackled.

and L. PETROU

E. Random measurement noise is introduced in Y and ST measured values with noise-to-signal ratio of 5 and 10%, respectively. Two subcases are considered: (ny = 0, n s = 1), (ny = 1, n~ = 0). The control performance under random measurement noise is given in Figs 14 and 15. The mean value of the absolute error is approximately 3.5% for ST and 2.6% for Y. The corresponding mean values of the absolute prediction error is approximately 3.3% for ST and 2.5% for Y. In other words, even under the presence of

"~ 0.015 - Set Point .......... n y = O , n,= I ..... n y = O . 0 2 , n,=0.98

"-" 0.012 O

0.009

¢...)

0.006

] 0.003

o~

i' u~ 0.000

p J i

0

I i

i

r I ,

20

,

,

40

I

f

r i

60

[

~ r i

80

T i m e (Days)

I J i

i

100

I

i

f

i

120

I i

i

i

140

160

(a) -~ o.o15 -

-

Plant

.......... P r e d i c t i o n

"-~ 0.012

(nT=O,

n.=l)

o 0.009

o~ 0.006 ¢..) .........

0.003

Time (Days) 0.000

,

0

i

i

I

20

'

i

i

I

40

i

i

i

i

60

I

i , l i

i 80

t

I

100

i

i

I

I

120

i

I

i

I

t40

I

I

I

160

(b) Fig. 8. A n e x t r e m e c a s e o f input c h a n g e s that w o u l d d r i v e the o p e n - l o o p s y s t e m to failure. (a) S u b s t r a t e c o n c e n t r a t i o n r e s p o n s e . (b) S u b s t r a t e c o n c e n t r a t i o n p r e d i c t e d r e s p o n s e . (c) M e t h a n e p r o d u c t i o n rate r e s p o n s e . (d) M e t h a n e p r o d u c t i o n rate p r e d i c t e d r e s p o n s e . (e) C o n t r o l a c t i o n . (f) p H r e s p o n s e . (g) V a r i a t i o n o f the specific g r o w t h rate.

Anaerobic digestors

133

20 ny=O, n,=l ..........ny=O.02, n,=0.98

/ J

-~ 16

o

o

Time (Days) ~

4

I

[

I

I

I

I

20

[

i

I

I

40

I

'

'

I

i

60

i

i

80

I

i

i

i

100

I

~

I

I

120

I

'

l

140

160

(o) ~.~ 20

Plant ..........Prediction

- -

(ny=0, n,=l)

Y

O

~

4-

~

0

Time (Days) i,

0

,

I

i,,

20

i

i

i

40

i

i

i

i

i

i

i

,

80

60

I

i

i

s

100

I

,

i

s

120

I

,

,

s

140

160

(d) 0.16 ny--0, ..........

ny=O.03,

n,= 1

n,=0.98

0.12

e~ o 0.08

,/

0.04 ~9

Time (Days) 0.00

'

0

20

40

60

80

(e) Fig. 8

(continued).

I

100

J

'

q

I

120

i

i

,

I

140

'

'

I

160

134

C. EMMANOUILIDESand L. PETROU 6.8 .......... 1"17=0' ny=0.02,rlt=l n,=0.98

l

6.7

~6.6

6.5

' r ' , ' ' ' , ' ' ' , ' ' ' l ' '

6.4

0

20

40

'''

I

60

80

100

I YSl

120

140

160

(f) 0.20

0.15

2

1

lay=O, 1"1.,=1 .......... ny=O.02, n,=0.98

0.10 O

0.05 T i m e (Days} 0.00

I

0

I

t

I

20

I

I

I

I

40

l

I

I

I

I

60

i

I

I

i

80

I

I

I

100

I

I

l

120

I

I

I

I

140

I

I

I

160

(g) Fig. 8 (continued).

significant measurement noise the control performance is acceptable. All the above simulations involved on-line model and controller training. The mean value and the standard deviation of the CPU computational training time (HP Apollo 735) consumed for the model training were below 0.3 and 0.45 s, respectively for all the simulations with the exemption of the noisy measurement cases that involved 1.22 and 0.01 s, respectively, i.e. the maximum allowable training time. The on-line controller training have been carded out either employing the generalised minimum variance law or the simple quadratic criterion together with the linear filter of equation (20). In the former case only chemotaxis and random search

algorithms were employed, while in the latter the back-propagation algorithm was also involved. In all cases the random optimisation techniques converge much faster than the back-propagation. An indicative chart is shown in Fig. 16, where the CPU time required for some of the control tasks of the cases A and C is illustrated. It is seen that the use of the simple quadratic criterion involves very fast training. This happens because the initial control network has been already trained to provide the correct control action and no significant on-line training is required. On the other hand, the superiority of the random optimisation techniques, in terms of CPU training time, is clear. The algorithms employed exhibit similar CPU training time performance in all the simulations tests performed.

Anaerobic digestors

135

~" 0.9

~ 0.7 0

.~ 0.5

eL

Time 0.3

IlIllllllllllllllllllllll

25

50

(Days)

Jlllllllllllllllllllllll

75

100

125

150

175

200

225

250

(at 39.0 38.8 (.J o v

38.6 38.4

t~ O9

38.2

~9

38.0 37.8 Time (Days) 37.6

II'l~lllllllll'llllJlllll~llllllllJlll,,IJllll,

0

25

50

75

I00

125

150

175

200

225

250

(b) Fig. 9. Random input disturbances. (a) Input substrate concentration. (b) Temperature.

5. D I S C U S S I O N A N D C O N C L U S I O N S

This work demonstrates the use of adaptive, online trained neural networks for identification and control of a complex, nonlinear bioprocess, i.e. of an anaerobic digester. The process has been simulated with a valid mathematical model. It is shown that under a predictive control approach the proposed neural network based controller has some clear advantages over traditional control strategies, mainly due to the on-line adaptation of the neural networks, which is based on the latest available measurement information. First it exhibits desired tracking, regulation and robustness properties in cases such as set points or process input variations. Then it successfully overcomes problems arising

from the nonlinear nature of the process, the presence of measurement noise or from process parameter changes through time. In addition the presence of a well-constructed training data set that ensures persistent excitation is unnecesary, since the on-line nature of the networks training improve modeling accuracy and control performance, by adding significant flexibility to the control scheme. All the above are accomplished with such a short computational time involved, that poses no problem to real time implementation, since anaerobic digestion is a process with slow dynamics. The initial model and controller training have been carried out using the back-propagation algorithm, but options like adaptive learning rate, momentum, introduction of an energy term in the

136

C. EMMANOUILIDESandL. PETROU

0.0065

~

.......... S e t P o i n t _

0.0055

O.0045 o

2

0.0035

Time 0.0025

I

llllll~lrl~llll~llliJlJlllIF~lll 25 50 75 100

125

I I

I

i

175

150

I

I I I

(Days) I

I

I

r I

I I

200

225

250

Time

(Days)

200

225

(a)

g 9

6

25

0

50

75

100

125

150

175

250

(b) 0.10

"

0.09

a~ ~: 0.08 o 0.07 e~

0.06 Time 0.05

I I

I I

i

25

I

I l

I

]

50

I

I J I

I

75

I

I I

r

i

r i

100

I

r ]

I I

125

I I i

I

150

I

I I

I

175

i

I

i

I I

200

(Days)

I I i

[

225

i

I I

J

250

(c) Fig. lO. Substrate concentration control under the presence of random disturbances in S, and O. (a) Substrate concentration response. (b) Methane production rate response. (c) Control action.

Anaerobic digestors

137

0.006 @

- -

fly=l,

11,=0

0.005 o 0.004

@ 0.003

0.002

Time (Days) 50 0 . 0 0 1

I!~l,l~lllrlllJ~,,llllll,,llllllll,,lll,llllll,I 25

50

75

100

125

150

175

200

.......... S e t

- -

225

250

Point

ny=l,

n,=0

O ¢9

=

8

0 la, 7t~

Time It

I tilt 25

l l i l l i l [ l i t i J l l l i l l l i t l l l i l l l l t 50 75 I00 125 150 175

(Days)

ilililillll 200 225

250

(b) ,-., 0 . 0 9 •~1

-

-

ny=l

ni=O

o.o7 o

2 0.05

2 Time (Days) 0.03

i Fi i i i ~ ~~ i i i j i I 0 25 50 75

t f i i i i ~ i i i r i i ~ t i i [ ~ i i i i i i i i [ I ~i i 100 125 150 175 200 225 250

(c) Fig. l l . M e t h a n e p r o d u c t i o n r a t e c o n t r o l u n d e r t h e p r e s e n c e o f r a n d o m d i s t u r b a n c e s in S0 a n d 0. ( a ) S u b s t r a t e c o n c e n t r a t i o n r e s p o n s e . ( b ) M e t h a n e p r o d u c t i o n r a t e r e s p o n s e . (c) C o n t r o l a c t i o n .

138

C . EMMANOUILIDES a n d L . PETROU

"~ 0 . 0 0 6 .......

Set Point

- -

n,=O,

n°:l

e" O :~ 0 . 0 0 5


@ 0.004 .,J oO

Time 0.003

I I I I I I i

(Days)

I I I I I t I ~ I I I I I ~ I I I I J I I I I I I I } I I t t I I t I I I i

25

50

75

100

125

150

175

200

t r T a

225

250

(al 0.006 Plant ....... Prediction (n~,=O, n , = l )

-

o

-

0.005

0 0

0.004 L

Time 0.003

u r r r I J t ~ I

25

I I I r i I I t u I I r t I I [ I t I t [ I I t I J I t t

50

75

LO0

125

150

175

(Days)

I [ t I I I [ t t I I

200

225

250

(b) ..... 0.12 ¢¢ ny=O, n , : l 0.10

e. o 0.08

\

0.06

Time 0.04

Illrlll

0

25

50

75

fllrll

100

125

(Days)

~lll~llll~llllllllrll

150

175

200

225

250

(c) Fig. 12. Substrate concentration control under process parameters variation through time. (a) Substrate concentration response. (b) Substrate concentration predicted response. (c) Control action.

Anaerobic digestors

139

111_=

....... Set Point - -

10

8C

ny=l,

n,=O

9o 0

Q

~

~,___._._J~__ ....

7

6

Time

,,,o

5 0

r t t t It

t t t It

25

50

t t ~ll

t t ~]l

75

l l J t l l l l t

125

100

(a)

111 ~

~

t i t ] t i

150

tint

175

(Days)

lit['tll]

200

225

250

Plant ....... Prediction (ny= 1, n . = O )

- -

_

9

o o

8

l

o

e, ~

6

Time ~

5

t t t t [ t t [ t I t I I t t t t't

0

25

50

75

t I t

(Days)

I i i I I I I [ I ; I I J i I I I I I I I I I ] I I I

lO0

125

(b)

150

175

200

225

250

~ . 0.09 ny=l,

n°=0

~9

0.07

f__l

o

¢~ 0.05

e-

Time 0.03

[I

0

25

50

75

LO0

Ill,Ill 1.25

150

~lll] 175

I]ll

~lll 200

(Days) llll~ 225

250

(c) Fig. 13. Methane production rate control under process parameters variation through time. (a) Methane production rate response. (b) Methane production rate predicted response. (c) Control action.

140

C. EMMANOUILIDES a n d L. PETROU

*~ 0.008 ....... S e t P o i n t - ny=O, n.:l

@

0.006 0) C)

0.004

T i m e (Days) 0.002

I i

i

t [ i

r i I I t i i

25

50

i [ I I i [ ] I i i I [ i I r I [ I I [ I [ t I 1 I ] t I i I

75

100

125

150

175

200

[ I i

i I

225

250

:>-, "~ 13 .,.J

- -

ny=O,

n,=l

11 o 9 @ 7 T i m e (Days) Q)

5

r J l l l l l l l l l l l l r l f l l l r l l l l r l l l [ l l l l [ l l l l l l l l l l l l l ~

25

50

75

ZOO

125

150

175

200

225

nr=O,

n,=l

250

(b) . ~ 0.12 t - -

"-~ 0. I 0

o

-

0.08

0.06

T i m e (Days) 0.04

Jl~tllr

0

II

25

[ l l l l i l l r l i i l l l l l t l l l l i i l l l r

50

75

100

125

lllll

150

175

200

Jlltllt

225

250

Fig. 14. Substrate concentration control under the presence of random measurement noise. (a) Substrate concentration response. (b) Methane production rate response. (c) Control action.

Anaerobic digestors

141

"~ 0 . 0 0 7 - -

Ily=I,

n,:O

0.006 0 L

0.005 @

-~ 0 . 0 0 4

Time

(Days)

0.003 25

50

75

100

125

150

175

200

225

250

(a) 13 .......

Set

- -

ny=l,

Point n,=O

3 0

9. o e~ 7-

5 0

;llltllltlll;tll*lllll;Elll;lll;lll;ll;I;l~llllll 25 50 75 100 125

150

175

Time

(Days)

200

225

250

(b) 0.12 - -

fly=

ns=O

1,

0.10 e~ e" o 0.08

r~

0.06 e"

Time 0.04

t t I I I I ; I I I I t It

25

50

I I I FI[I

75

100

t~

Jl

Ill

125

I ] I I I I I i I FI]t

150

175

200

ill

(Days) II

225

lii

250

(o) Fig. 15, Methane production rate control under the presence of random measurement noise. (a) Substrate concentration response. (b) Methane production rate response. (c) Control action.

C. EMMANOUILIDESand L. PETROU

142

CPU Training Time (see) 0.1

0.08

[ ] Chemomxis & GMV

i

0.04

0.02

• Random Search & GMV

I

0.06

dI ny=0,ns=l ny=l, ns=0 RandomInputDisturbances(CaseC)

• Random Search & Filter Illl Chemotaxis & Filter [ ] Backpropaption & Filter

m:i ny=0,ns=l ny=l, ns=0 SquareWaveSet PointChanges Followedby Sto Changesin So, O (ca A)

Fig. 16. Backpropagation, chemotaxis and random search performance in terms of speed of convergence.

cost function a n d the use of an e r r o r suppressor f o r m u l a h a v e b e e n employed. T h e on-line m o d e l t r a i n i n g is b a s e d o n a m o v i n g - t i m e - z o n e , fixed-data length m e a s u r e m e n t data set, w h e r e a b a c k p r o p a g a tion a l g o r i t h m v a r i a n t has also b e e n employed. O n the o t h e r h a n d , t h r e e different algorithms, i.e. back p r o p a g a t i o n , c h e m o t a x i s a n d a r a n d o m search technique were e m p l o y e d in the neural controller online training task. F u r t h e r m o r e , in the last two cases, the generalised m i n i m u m variance law has b e e n c o n s i d e r e d , in o r d e r to reduce the control effort. It is s h o w n t h a t , for the same tracking and regulation p e r f o r m a n c e , r e d u c e d control effort is a c h i e v e d w h e n the generalized m i n i m u m variance law is e m p l o y e d , at the expense of increased C P U t r a i n i n g time. In all cases b o t h the c h e m o t a x i s and the r a n d o m search algorithms converge much faster t h a n the back p r o p a g a t i o n algorithm.

NOMENCLATURE

D = gas volume conversion factor (I/mo 0

Dm= Moving-time-zone fixed-length model training data set E{.} = Cost function (employed in the random search algorithm) Ec = Cost function employed for the on-line controller network training Em = Cost function employed for the on-line model network training E,r = Error term of cost function E E¢. = Energy term of cost function E e~.(-) = Energy function of neuron activation f(.) = Coefficients of the linear filter employed at the neural controller output g(.) = Neuron activation function h = Neuron net input H = Sampling Interval K = Yield coefficient, K = I/(D. V. YCH~IX"Yx/s) (Polihronakis et al. 1993) ka = decay coefficient (l/day) lm= Fixed data length of Dm ny, n~= Weighting coefficients that determine the weighted influences of Y and ST error in the evaluation of Ec O = Neuron output q = Generalised minimum variance control law tuning factor S, = lnfluent substrate concentration, (mol/I, as equivalent acetic) Sr = Total substrate (volatile acids) concentration, (mol/I, as equivalent acetic) S'r. = Set point value of Sr

Anaerobic digestors = Derivative of S.r ~v = Estimate of the future value of St, offered by the neural model t = Time (days) U = lnfluent dilution rate, (I/day) U~= Suggested by the neural controller value of U V = Digester liquid volume (I) Wm= Model neural network weight matrix Xm(n) = Input pattern vector of Din, consisting of the measured values of Y, S v, pH, S., 0 and the manipulated control input value U at the (n-I)-th time instant X0 = Influent organism concentration, (tool/I, as equivalent CsHTNO2) Y= Methane production rate, (I/day) ~"= Derivative of Y Yx/s = Yield of organisms per mole of substrate consumed (Buhr and Andrews, 1977) Ycn4/x = Yield of methane per mole of organisms formed, (Buhr and Andrews, 1977) y,,(n) = Output pattern vector of D,,, consisting of the measured values of Y, S~r at the n-th time instant '~'= Estimate of the future value of Y, offered by the neural model y , = Set point value of Y am = Model adaptation rate t5 = Neuron activation slope parameter /~ = specific growth rate, (l/day) kt,,/~e. = Weighting factors of E~r and E~.,, respectively.

REFERENCES

Bhat N. and McAvoy T. J., Use of neural nets for dynamic modeling and control of chemical process systems. Computers chem. Engng 14, 573-583 (199(I). Buhr H. O. and J. F. Andrews, The thermophilic anaerobic digestion process. War. Res. 11, 12%143 (1977).

143

Carr A. D. and R. O'Donnel, The dynamic behaviour of an anaerobic digester. Prog. Wat. Technol. 9, 727-738 (1977). Chauvin Y., A back-propagation algorithm with optimal use of hidden units. Advances in Neural Information Processing Systems, Vol. 1 (Edited by D. S. Tourtezky), California, Morgan Kaufmann (1989). Di Massimo C., G. A. Montague, M. J. Willis, M. T. Tham and A. J. Morris, Towards improved penicillin fermentation via artificial neural networks. Computers chem. Engng. 16, 283--291 (1992). Hertz J., A. Krogh and R. G. Palmer, Introduction to the Theory of Neural Computation. Addison-Wesley, New York (1991). Hunt K. J., D. Sbarbaro, R. Zbikowski and P. J. Gawthorp, Neural networks for control systems--a survey. Automatica 28, 1083--1112 (1992). Kosko B., Neural Networks and Fuzzy Systems, Prentice Hall, Englewood Cliffs (1992). Pollard J. F., M. R. Broussard, D. B. Garrison and K. Y. San, Process identification using neural networks. Computers chem. Engng. 16, 253-270 (1992). Polihronakis M., Control methods of the anaerobic digestion process. Ph.D. Thesis, Aristotle University of Thessaloniki, Greece (in Greek) (1992). Polihronakis M., L. Petrou and A. Deligiannis, Parameter adaptive control techniques for anaerobic digesters-real life experiments. Computers chem. Engng. 17, 1167-1179 (1993). Soils F. J. and R. J-B. Wets, Minimization by random search techniques. Math. Op. Res. 6, 1%30 (1981). Willis M. J., C. Di Massimo, G. A. Montague, M. T. Tham and A. J. Morris, Artificial neural networks in process engineering, lEE Proc.-D, 138, 256-266 (1991). Yang Y. Y. and D. A. Linkens, Adaptive neural-networkbased approach for the control of continuously stirred tank reactor, lEE Proc.-Control Theory Appl. 141,341349 (1994). Ydstie B. E., Forecasting and control using adaptive connectionist networks. Computers chem. Engng. 14, 583599 (1990).