Eigenvalue problems arising in the control of a distributed-parameter bioreactor

Control Eng. Prucfice, Vol. 4. No. 7, pp. 1015-1021, 1996 Copyright 0 1996 Elsevier Science Ltd Printed in Great Britain. All rights reserved 0967~0661/96 $15.00 + 0.00

PII:SO967-0661(96)00101-3

EIGENVALUE PROBLEMS ARISING IN THE CONTROL OF A DISTRIBUTED-PARAMETER BIOREACTOR M.T. Nihtila*, J.P. Babary*** and J.P. Kaipio** *Department of Computer Science and Applied Mathematics, University of Kuopio, Faculty of Natural and Environmental Sciences, P.O. Box 1627, FIN-70211 Kuopio, Finland ([email protected]) **Department of Applied Physics, University of Kuopio, Faculty of Natural and Environmental Sciences, P.O. Box 1627, FIN-70211 Kuopio, Finland ***Centre National de la Recherche Scientifique. Laboratoire d’Analyse et &Architecture des Systemes, 7. avenue du Colonel Roche, F-31077 Toulouse Cedex, France

(Received

October

1995;

infinalform

April 1996)

Abstract: A distributed-parameter model of a continuous-flow fixed-bed reactor is studied. The main emphasis lies on a structural property of the partial differential equation (PDE) system model. This property, which in lumped parameter systems is called the nonminimum phase property, has certain implications in the controller design. The controller design for the PDE model is constructed via semidiscretisation. The only space variable of the PDE model is discretised by using Galerkin’s finite element method (FEM). For some parameter values of the PDE model the linearising control of the semidiscretised model results in unstable behaviour in the sense that the zero dynamics of the model is unstable. This same unstable behaviour was also earlier observed in using the orthogonal collocation for the semidiscretisation. Connections between the location of the zeros of the original PDE model linearised around its steady-state solution and the stability/instability properties of the linearising control of the semidiscrete model are discussed in relation to the same issues in lumped-parameter differential system models.

Keywords: Distributed-parameter systems; Orthogonal collocation element method approximation; Linearizing control; Non-minimum dynamics; Linearized transfer function.

approximation; Finite phase property; Zero

compactness and high efficiency for biological (removal epuration of nitrogen and carbon their good integration in the compounds), environment (no noise, no smell), their low energy consumption and sludge production and the fact that they do not need secondary clarifiers.

1. INTRODUCTION The system studied here is a continuous-flow fixedbed bioreactor. It is aimed at removing harmful nitrogen (and carbon) compounds from the influent waste water. The goal is to decrease the overall concentrations at the outlet of the reactor under some prescribed allowable limits. Most conventional wastewater plants are based on the activated sludge process. The new type of biological wastewater treatment plants considered here are also called biofilter processes. Their main advantages are their

A submerged granular fixed-bed biofilter can be considered as a continuous-flow distributed-parameter bioreactor, in which the wastewater flows through the tubular reactor. The microorganisms fixed in the reactor absorb these harmful nutrients in such a way 1015

1016 that the substrate concentrations in the water that flows out.

M.T. Nihtila et al. 2. SYSTEM MODELS

have been decreased

Due to the clearly distributed nature of the bioreactor it is modelled by a set of two partial differential equations of parabolic type. The model is nonlinear due to the dependence of the specific growth rate of the microorganisms on the concentrations of nitrogen and biomass. The model has two independent variables: the time-to-go and the space variable, the distance from the inlet of the reactor. Some basic assumptions in modelling have been applied; for instance: the biomass and nitrogen concentrations are both described by one state variable, the gaseous phase is not taken into account because the gaseous hold-up may be ignored as compared to the liquid one, and the temperature and acidity (pH-value is close to 7) are constant. The modelling is based on mass balances and on the specific growth rate model of Contois for the biomass. The key technique in studying the behaviour of and designing a controller for the system is to approximate the infinite-dimensional PDE model by a nonlinear finite-dimensional ordinary differential equation (ODE) model. A method of orthogonal collocation (Tali-Maamar, et al., 1993) and a finite element method (FEM) of Galerkin type (Nihtilti, et 1994) have been applied to construct the al., approximations.

The bioreactor of the denitrification process is modelled by using partial differential equations, due to the clearly distributed nature of the system. Different variants of the basic model structure, which include compartments for the biomass, which is fixed, and for the substrates, which flow with the flowing liquid medium, have been developed. A fourdimensional model is studied in (Babary, et al., 1996). Two-dimensional models with and without diffusion were considered in (Dochain, et al., 1992) and (Nihtila, et al., 1992). Here the system is described by a two-dimensional model with one state variable for the biomass and the other for the substrate. The overall equations:

model

is

given

by

the

following

(1)

t

ElO,T[ ;

zEGn]O,l[

Pow2) = PLmU21&J4,U, +%I Direct application of the linearizing control technique in the controller design of the approximated (semidiscretised) systems results, however, in unstable behaviour of the internal modes of the controlled system for some values of the process parameters (c.f. (Tali-Maamar, 1994)). It is a known fact that in nonlinear systems which are of non-minimum phase character, i.e. their zero dynamics is unstable, the application of linearizing control techniques leads to internally unstable closedloop systems (Isidori, 1989). The origin of this study, consequently, lies in the question: Is the nonminimum phase property of the ODE systems obtained via semidiscretisations also the property of the original infinite-dimensional PDE system, or is it only due to the approximation methods applied? The organisation of this paper is as follows. The system models - the original PDE system and its approximation obtained via FEM-semidiscretisation are presented in Section 2. The linearising control idea and the algorithm for the semidiscretised model with some simulation results are given in Section 3. The transfer function of the PDE system model linearised around the steady-state solution function is presented and discussed in Section 4. Some concluding remarks complete the paper.

(3)

in which u, and u2 denote the concentrations

of the

biomass and nitrogen, respectively. The specific growth rate model, Eq.(3), is the model of Contois (1959). The boundary conditions are given by

(4)

f&Q)=o. Functions

c(t)

(5) and y(t) = ~Jl,t)

are the control

input, defined as the flow rate per the cross-sectional area of the bioreactor, and the output, correspondingly. The initial condition functions which are obtained for constant values of the control and the disturbance variable s,(t) are chosen as the steadystate solutions of the problem.

u,(z,O)= qz>

(6)

u,(z,O) = qz)

(7)

For the chosen specific growth rate model of Contois (3) these are given by the equations

Eigenvalue Problems in a Distributed-Parameter Bioreactor

Qz) = t qcosh[q(l-41 +Pwdl

(8)

-zwpwz)l

[qcosbq+ (a ‘PI dnhgl

&(~,,,-k,Y@&,) The coefficients

(9)

p and 4 are defined by (10)

a =k,(~r,,,-k~)I(&~

(11)

4=/m

(12)

A firm mathematical basis for the problem is given in the form of the following existence theorem. Theorem

1.

Suppose that c and s,(t)

C*([O,-[ ) functions

Finite element methods are mathematically well proved and justified spatial discretization methods for PDE systems (Thomee and Wahlbin, 1975; Thomee, 1984). Their error estimates and convergence proofs as a function of the number of discretisation points are based on the abstract theory of the solutions in a weak sense in suitable Sobolev spaces. Here, the solutions of the initial/boundary-value problem are approximated by finite series

; p=P/2

P=FjD

(15)

in which the base functions $i are triangular, positive, piecewise linear functions. Application of the series (15) in the weak form of the PDE system results in the ODE system

are positive Mi = -DKx-c(t)[Rn+w(t,x)]

with the property

ut E

solution

u = (ut, ZQ such that

C1(14T[,C(~))nC<[o,T[,C(G))

nCW,TLC2(&). The proof is based on the abstract formulation problem in the form

where

pT]T

+A(t)v = P&v)

; v(0)

= v,,

v = (v,,v,)

= (ut, $ -so)

where the 2N-vector x and the coefficient matrices M, K, and R, and the vector-valued functions w and F can be found in (Nihtill, et al., 1994). The orthogonal collocation method, studied for example by Villadsen and Michelsen (1978) offers another alternative for spatial discretisation. Then the base functions & are obtained via Lagrange interpolation polynomials. Details of the application can be found in (Tali-Maamar, et al., 1993; Tali-Maamar, 1994).

u2 E C’(lO,T[,C(~))nC<[o,T[,C(G))

dv/dt

= [d

(13)

; s,(O) =i,

then the problem of the equations (l-3) with the boundary and initial values of (4-9) has a unique classical

+F(x) (16)

x(0) c(0) = c

1017

of the The system (16) can be put into the basic control theoretic form in the case of the constant disturbance (14)

is in a suitably

defined Sobolev space. The operator A(t) , which generates an analytic semigroup for any fixed t, and E is a uniformly Lipschitz continuous function with respect to v. It is obtained from the Contois model (3). The proof is outlined in (Nihtill, et al., 1994) and given in detail in (Goldstein and Tervo, 1996).

s,(t),

i

i.e. (17)

=f(x) +g(x)40

(18)

r=W) where the output mapping

is given by

h(n) = [0 ... 0 11.x.

(19)

2.1 Semidiscretisation The input-output behaviour of the PDE system was simulated in the open-loop and closed-loop cases. Galerkin’s finite element method was applied to convert the PDE system into a set of nonlinear ODES. Details of the semidiscretisation are presented in (Nihtila, et al., 1994).

3. LINEARISING CONTROL OF THE SEMIDISCRETE PROBLEM The relative degree of the model (17-18) is obtained by differentiating the output function h(n) with respect to the time. Then by taking into account (17) the control c(f) appears in the diffentiated giving the relative degree r = 1.

output

1018

M.T. Nihtila et al.

The asymptotic output tracking of the reference input ya is obtained by applying the linearising control law (Isidori, 1989; p.191, Eq(5.2)) c(t) =

+-$-q(x) 8

+ 3,

- G-J[W) - YRI1

Then the closed-loop transfer function

G&)

(20) ’

system (from y, to y) has the

4. ON THE TRANSFER FUNCTION THE LINEARISED PDE

A possible way of studying the qualitative behaviour of the linearising control of the semidiscretised system is to study the zeros of the (transcendental) transfer function of the linearised system. Now, the linearisation of the PDE system means linearisation around the steady-state

(21)

= l/CA +co).

Typical parameter values, in the case where N = 8, which gave the stable linearising control for the semidiscrete (ODE) system, are given by R, =

0.4

kd = 0.05

p, = 0.35

D = 0.05

c

say h(t),

the small changes,

small changes, say 6y(t), given by the equation au(A)

disturbance

given

of the input on the

of the output is formally

(22)

A(r(A))

’

After some tedious calculations, including the solving of the non-constant coefficient linear ODES, the numerator and the denominator functions are obtained in the following form

A(r) = [(p2 +r2)sinhr + 2prcoshr](r2

-q2)

(23)

= 0.1 . B(r) = b,rcoshr+b2sinhr+r3

A simulation

ii = (&, &)

= j+(U)

WA)

0.4

function

by Eqs (8-9). The transfer function (as a function of the Laplace variable h) of the system linearised for

In order that the linearising control law also results in an internally stable system, the zero dynamics of the model (17-1 S), with respect to the (constant) reference input y,, must be stable.

k, =

OF

-b,r

(24)

run is depicted in Fig. 1 for c,, = 1. The so changed

from 5.0 to 4.5 at t = 1.66,

and the reference ya from 0.77, that corresponds to the value 0.1 of the constant control, to 0.87 at t = 10.

b, = 2pqlsW b2 = 2p2q/sinh4

(25)

b3 = q2 +2p2 +2pqcoshq@hq K = (alD)s,W[qcosh4+(a

+P)sWI

in which the inner function r is given by the equation Disturbance

4’

i

0

/

5

r(h)2 = (h2+dll

10

I 20

15

Output and referent y

+4)/[0(1

+a,)]

(26)

d, =a, +k,a,+c7./(4D) 4 = k,kda, +a,F2/(4D)

(27)

a, = (I -kdlcL,)kd az=(I-kd/cL,)(cL,-k,>lkz. It has to be noted that when considering

the relation

between the location of the zeros of B(r( ii)) on the complex h-plane and the stability of the corresponding transfer function l/B(r( h)), the original overall transfer function can also be given by Fig.‘l.

Simula&n the disturbance

of step Responses fb”, the steps2& s, and in the reference

ya

by(h)

= Kf(r(“))

WA)

W(V)

(28)

1019

Eigenvalue Problems in a Distributed-Parameter Bioreactor

A(r) =[p2f(r)+rsinhr B(r)

(29) (30)

=bphr+b~(r)+r~-b3 = sinhr/r

f(r)

+2pcoshr](r~-q2)

By considering functions

f(r)

(31)

. the power and

series expansions

of the

It2 +(dl -q2D)12

oo&r

(32)

@%+~+~+;+... .

.

.

Mshr=l+~+~+~+... . . .

a =-a,

(36)

which has the zeros )c,=o;

functions A (r( I)) and B( r( J.)) as functions of h have an essential singularity at the singularity of the r(A),

=0

-(&-l&l

A, = (&/&

- 11% .

(33)

it is seen that B(r) in fact depends only on fl but not on r itself. The Riemannian cutting of the complex plane, due to the square root in the expression of r in calculating the inverse transforms may then be avoided. These issues need, however, a more detailed and careful analysis, especially to be able, if it is at all possible, to connect the location of the zeros and the stability/instability of the linearising control. Furthermore, the numerator and denominator

function

(35)

=0

The former equation, due to the fact that the coefficients d, and d, are strictly positive, gives the result that the both roots have negative real parts. The latter reduces to the equation h[l

r

+d, -q2Da,

i.e. at the point

However, based on these singular zeros, for k, < k,, nothing can be concluded about the stability of the FEM-approximated semidiscrete system. The abovementioned common pole-zero group seems to be the result of the uncontrollability of the output of the original PDE system. It is clear that the output always stays in the following interval $(I,0

E

W,s,l

for a constant disturbance, uncontrollability property.

showing

a kind of output

A subset of the ordinary zeros in the r-plane and the corresponding h-plane are given in Table 1. Furthermore, it can be concluded that for the sequence of r-zeros, which are pairwise ri and ii,

(
This essential singularity cannot be removed by multiplying the numerator and denominator functions by some suitable power

iF+,

Iri I = ma

The corresponding equations as is evident by considering the substitution of the power series (32-33) to the equations (29-30).

2

I2

+d,l

(37) h-zeros

are obtained

from

+d,

(38)

ri = The parameter values which also gave stable linearising control for the semi-discrete system, and which were applied in calculating the system zeros, are given by (cf. the list below equation (21)) k, =

0.4

k, =

0.4

kd =

0.05

p, =

0.35

D =

0.005

F

0.2

s, =

5.0 .

=

A specific feature of the transfer function

the

D(A + a,)

7; =

A2 +d,3c +d,

(39)

D(A + a,) Due to the infinite limit (37) at least two of the corresponding four h-zeros (pairwise complex conjugate values) have the limits lim li = *lim ii = -a,, i++CS r-.+m

is that the

functions A(r( A)) and B(r( A)) have three common zeros in the r-plane: r,, = 0, r,,2 = 2 q. These correspond to the four zeros in the k-plane obtained from the equations (c.f. (26)) (34)

as is evident from the root locus principle. Actually, this is also seen in Table 1, because in the last row there are Z?e(l& = -0.04289;

Zm(Q

=O (q

=0.04286).

1020

M.T. Nihtila et al. Table 1. Some zeros of the linearised

transfer function

- r-zeros: rk = r, +

-21 -21 -21 -21 -22 -22 -22 -22 -22 -23 -23 -23 -23 -23 -23 -23 -24

ar2

;

Fk = r, - ir, ;

6 12 18 24

2590432 5060987 7453690 9867039 2357998 z7' 4939537 43 7602692 50 0331573 56 3110428 62 5925789 68 8767166 ii 1626663 81 4498291 87 7377930 94 0262299 100 3149414

5. SOME SPECULATIONS, AND CONCLUDING

0000000

-.0430097 -.0430066 -.0429933 -.0429778 -.0429626 -.0429488 -.0429368 -.0429268 -.0429177 -.0429107 -.0429049 -.0428993 -.0428943 -.0428910 -.0428866

OPEN QUESTIONS REMARKS

The infinite-dimensional nature of distributedparameter systems is the main obstacle in the efforts to transfer concepts and ideas that are directly usable in finite-dimensional, even if nonlinear, ODE systems. The concept of linearising control may have its counterpart in PDE systems. Especially, the concept of zero dynamics can be formulated (to be done later on) as the inverse problem of the original control problem of Eqs(l-7). On the other hand, the boundedness of the solution of the inverse problem as the time goes to infinity is not necessarily directly related to the location of the zeros of the original PDE system linearised around its steady-state solution, as it is in the case of finite-dimensional ODE systems. The (transcendental) transfer function of the linearised PDE system has four common pole-zero pairs. This is seemingly due to the evident uncontrollability of the output, which corresponds to the analogous case in linear ODE systems. The practical treatment, e.g. the controller design for PDE systems having pointwise boundary control and output, can be based on some semi-discretisation, i.e. spatial, procedure. Then the system model is approximated by an ODE system. Several theoretical questions, however, remain unanswered in this paper: What is the relation of the instability of the linearising control in the ODE approximations obtained via the orthogonal collocation and FEM for some parameter values and the location of zeros of

set.

A. = A1 + il, ; ii%_, = A, -il., 2k-1 = 1, +ih, ; Xx = A, -il, %k

0000000

4476109 5126419 6811600 9022617 1372890 3663654 5817642 7816715 9666996 1382580 2978764 4469700 5867672 7183132 8424950 9600658 0716686

- h -zeros:

from the infinite

0006409 0003902 0002611 0001834 0001334 0000996 0000757 0000586 0000464 .0000367 0000295 :0000250 0000208 :0000174 .0000153

-:1389555 -.6886393 -1.6153923 -2.9283857 -4.6341267 -6.7363396 -9.2368336 -12.1363010 -15.4348450 -19.1322803 -23.2282982 -27.7225494 -32.6146812 -37.9043808 -43.5913429 -49.6753235

0000745 -113478414 -2.7121081 -4.1060500 -5.5316396 -6.9864960 -8.4669294 -9.9694204 -11.4910402 -13.0294533 -14.5827875 -16.1495399 -17.7284737 -19.3185501 -20.9188995 -22.5287647 -24.1474953

the linearised transfer function of the original PDE system? Is it possible to formulate a general theory of linearising control for (a certain type of) PDE systems?

REFERENCES Babary, J.P., Bourrel, S.V., Nihtila, M.T., and Dochain, D. (1996). Sur la representation d’etat des systemes a parambtres rtpartis, in the book: Les Systhes de Rkgulation, A.Rachid, ed., Masson, France, 23 pp. (in print) Contois, D. (1959). Kinetics of bacterial growth relationship between population density and specific growth rate of continuous cultures, Journal

of General

Microbiology

21, 40.

Dochain, D., Babary, J.-P., and Tali-Maamar, N. (1992). Modelling and adaptive control of nonlinear distributed parameter bioreactors via orthogonal collocation. A utomatica, 28(9), 873-883. Goldstein, J.A. and Tervo, J. (1996). Existence results of solutions for a system of nonlinear partial differential equations related to bioreactor. Dynamic Systems and Applications, 10 pp. (in print) Isidori, A. (1989). Nonlinear Control Systems, 2nd Ed., Springer, Berlin. Nihtila, M., Tali-Maamar, N., and Babary, J.P. (1992). Controller design for a nonlinear distributed parameter system via a collocation approximation. Proceedings IFA C Symposium on Nonlinear

Control

Bordeaux, France, ed., pp. 447-452.

System

Design,

June 24-26,

NOLCOS92,

1992, M.Fliess,

Eigenvalue Problems in a Distributed-Parameter Nihtilti, M., Tervo, J., and Kaipio, J. (1994). Simulation of a nonlinear distributed parameter bioreactor by FEM approach. Reports, A5, Department of Computer Science and Applied Mathematics, University of Kuopio, 19 pp. Tali-Maamar, N. (1994). Modblisation, Analyse et Commande d’un Pro&de Biotechnologique a Gradient Spatial de Concentration. These de Docteur de l’Universin5 Paul Sabatier de Toulouse, Rapport LAAS no: 94079, 106 pp. Tali-Maamar, N., Damak, T., Babary, J.P., and Nihtill, M.T. (1993). Application of a collocation

method for simulation bioreactors. Mathematics lation,

1021

Bioreactor

of distributed and Computers

parameter in Simu-

31(1 l), 303-319.

Thomee, V. and Wahlbin, L. (1975). On Galerkin methods in semilinear parabolic problems. SIA M .I. Numer. Anal. 12, 378-389. Thomee, V. (1984). Galerkin Finite Element Methods for Parabolic Problems. Lecture Notes in Mathematics. 1054, Springer, Berlin. Villadsen, J. and Michelsen, M.L. (1978). Solution of Differential Approximation.

Equation

Models

Prentice Hall,

by Polynomial

London.