Dynamic estimation of temperature and concentration profiles in a packed bed reactor

Ckemral Engineering Science. Vol. 44, No. 10, pp. 2087-2106 Printed m Great Rrilain

1989.

OOiS2509/89 S3oO+O.W 0 1989 Pergamon Press plc

DYNAMIC ESTIMATION CONCENTRATION PROFILES

Department

LARRY C. WINDES,’ of Chemical Engineering, (Received

OF TEMPERATURE AND IN A PACKED BED REACTOR

A. CINAR* and W. HARMON University of Wisconsin, Madison,

15 December

1987; accepted 22 December

RAY* WI 53706, U.S.A.

1988)

Abstract-A two-dimensional nonlinear distributed-parameter optimal state estimation algorithm has been develoljed and applied to an exothermic packed bed catalytic reactor. Temperature and concentration profiles were estimated from limited temperature and/or concentration measurements. Orthogonal collocation was used to solve the differential sensitivity and filtering equations. Both simulations and on-line experiments with a pilot scale packed bed reactor for the partial oxidation of methanol have been used to study the effects of measurement choice, model parameter errors, and unmeasured process disturbances. Excellent tracking of the reactor profiles was achieved with only one or two properly placed temperature

systems has been slow. The theory of DPS estimation

INTRODUCTION

In order to exert the necessary control action for safe and profitable operation of an exothermic packed bed reactor, it is useful to have real-time estimates of temperature and composition profiles. However, this knowledge of the reactor state variables is limited by the type and number of measurements and the time delay in obtaining them. For example, a gas chromatograph analysis of a single concentration measurement often takes longer than the time transients in the reactor and cannot be used directly to control reactor dynamics. Many industrial reactors have limited measurements consisting of only a few thermocouples. Thus algorithms for providing reliable estimates of reactor temperature profiles, conversion, and selectivity for safety and process control purposes would be welcome tools in the process industry. State estimation, or filtering, involves generating the best current estimate of the state variables given a dynamic system with a limited number of noisy measurements over some past time interval. Errors occur in these measurements due to sampling or instrumentation problems, while errors arise in the process representation because of structural deficiencies or oversimplification of the mathematical model, uncertainty in the parameters for the model, and unknown internal and external disturbances in the operation of the physical system. Errors in the initial conditions are considered due to uncertainties in the state variables at the time when filtering begins. Since the tubular packed bed reactor is a highly nonlinear distributed parameter system, previous studies in these areas are of interest. The theory of state estimation for distributed parameter systems (DPS) was developed extensively during the 197Os, but the progress in application of these ideas to real +Present address: Eastman Kodak Chemicals Division Research Laboratones, Kingsport, TN 37662, U.S.A. $Present address: Department of Chemical Engineering, Illinois Institute of Technology, Chicago, IL 60616, U.S.A. $Author to whom correspondence should be addressed. 2087

has been reviewed by Tzafestas (I 978) and Sawarogi et

al. (1978), and applications were reviewed by Ray (1978). Several results and tables of the state estimation equations for linear DPS systems are presented by Bencala and Seinfeld (1979). The approaches that have been used for nonlinear estimation of distributed parameter systems have been reviewed by Tzafestas (1978). In this work we will develop a state estimation algorithm based on the least squares approach. Previous work in this direction with nonlinear DPS includes an approximate nonlinear filter by linearization (Seinfeld, 1969; Seinfeld and Hwang, 1970), a finite difference approximation (Seinfeld et a!., 1971): and estimation of a system parameter simultaneously with state estimation (Hwang et al., 1972). Following these basic filters, further work has included (1) discrete measurement points, discrete measurement times (Ajinkya et al., 1975), (2) estimation with time delays (Yu et al., 1974), (3) coupled ODES and PDEs (Yu et al., 1974; Ajinkya et al., 1975), (4) stiff systems having a small parameter (Soliman and Ray, 1979a). A comprehensive development of nonlinear filtering theory is given by Lausterer (1977) who built upon the earlier work of Yu (1973). A more detailed literature review of distributed parameter state estimation as it relates to this work is given by Windes (1986). There has been no distributed parameter state estimation of the dynamic temperature and composition profiles in packed bed reactors other than a few test simulations with a highly simplified reactor model (Soliman and Ray, 1979b). The early work in estimation for chemical reactors either assumed a linear, lumped model (Goldman and Sargent, 197 1) or used a finite difference technique in early lumping (Pell and Aris, 1970; Balchen et al., 1971). Some investigators have estimated the catalyst activity in the reactor

2088

LARRY C. WINDES ez ul.

while assuming the temperatures and concentrations were at quasi-steady state (Joffe and Sargent, 1972; Ajinkya et al., 1974; Kuruoglu et al., 1981). Kuruoglu and co-workers (1981, 1982) used DPS theory to develop the quasi-steady-state filter for the catalyst activity in order to generate temperature profiles. For estimation of temperature and concentration profiles in tubular reactors, the vast majority of investigators have followed these steps:

(1) application

of orthogonal collocation to the reactor model to generate one ODE per state variable at each collocation point (early lumping), (2) linearization of the nonlinear ODES, of the “standard” lumped par(3) implementation ameter Kalman filter. The first to use this procedure was Vakil et al. (1973), and many other examples have been reported since that time (e.g. Jutan et al., 1977; Jorgensen and Clement, 1977; Sorensen, 1977; Kumar and Seinfeld, 1978; Silva et al., 1979; Sorenson et al., 1980; Clement et al., 1980; MacGregor and Wong, 1980; Wallman and Foss, 1981). One case of noniinear estimation was the use of an extended Kalman filter in simulated tests (Rutzler, 1980). The alternative used in this paper is to apply distributed parameter estimation theory first and to solve the full distributed estimation equations by orthogonal collocation (late lumping). This has an advantage of retaining more of the model structure further in the estimation solution process. Early Iumping and late lumping for linear distributed parameters systems have been directly compared (Cooper et ol., 1986). The late-lumped fitters were found to

converge and track more efficiently. The full covariance matrices are retained in the final form of the late-lumped filter as opposed to diagonal covariances in the lumped Kalman filter. The distributed parameter filter incorporates a spatial functionality for the error weighting matrices which accounts for interactions between neighboring spatial locations. By contrast, the error weighting matrices from the early lumping method do not account for the spatial location. Finally, we use the full nonlinear process model in the equations for the state estimator. ln the present work, we have developed the theoretical estimation equations for two-dimensional distributed-parameter systems having a separation of time scales. These general equations arc then applied to the estimation of axial and radial tctnperature and concentration profiles in a packed bed reactor. The particular reactor studied is the pilot scale, jacket cooled, tubular reactor in our laboratory shown in Fig. 1. The reactor is designed to carry out the partial oxidation of methanol to formaldehyde and suffers from a secondary reaction of formaldehyde oxidation to carbon monoxide and water. To achieve the optimal yield of formaldehyde, one must have close control of the temperature and concentration profiles in the reactor; thus, on-line state estimation of these profiles is important. The work reported here includes the development of new nonlinear estimation results and the first application of rigorous state estimation theory to an experimental packed bed reactor. DEVELOPMENT OF A NONLINEAR DISTRIBUTED PARAMETER STATE ESTIMATOR The theoretical equations estimator have been derived

“E .1.

Fig. I. Schematic of the packed bed reactor pilot plant.

for a nonlinear in general terms

state for a

Estimation

of temperature

two-dimensional nonlinear distributed parameter system containing “slow” and “fast” states. This derivation drew upon the general discussion of Ray (1981) and is based upon earlier work dealing with nonlinear filters (Lausterer, 1977), linear two-dimensional systems (Lausterer et al., 1978; Lausterer and Ray, 1979), and nonlinear systems with a small parameter multiplying the time derivatives of some state variables (Soliman and Ray, 1979a, b).

General results The model equations take the form of a general nonlinear parabolic PDE: Kx,(v, t)= f (x, x,,, x,,, v, 0 + W, f)

(1)

with boundary conditions: O=bi(x,x,,

i=O,i=l

t)+&(t)

(2)

where v is the independent spatial variable vector in the N-dimensional domain Q; here x represents n state variables which are functions of space and time, x, and x,, are the first and second spatial derivatives of x, and the n x n diagonal matrix K incorporates the small parameter E in the form: “I xn* 0

K=

0 cl n2xn2

1

(3)

Thus, x is partitioned into nl “slow” states x1 and n2 “fast” states x2; the other variables and functions are partitioned similarly. The discrete measurement operator with m measurements made at location vo and time tj is: Y(cj)=h(x(vY,

rjh

x(vi3

rj)t

. . 3 x(viT

tj)+V(tj)*

t4)

Because of the nonlinear character and unavailability of rigorous mathematical results, a “formal” approach is adopted where the state estimation problem is recast as a deterministic optimal control problem and is termed optimal estimation. For optimal state estimation, it is desirable to minimize the errors over time in a weighted least squares sense: *=

lr 5 0 [SJ x

SI

[Kx,(u’, t)-f(x,

+Ff i

X

cl

[Kx,(v, c)-f(x,

(6)

comes directly from the resulting filter equation for a linear system. However, for a nonlinear system, P is not so easily computed, and we must make use of differential sensitivities as filter parameters. The filter algorithm requires solving two sets of equations-one for the differential sensitivities and another for the state estimates. The equations for the state estimates can be viewed as model equations coupled with a proportional controller with space and time varying gains(G) forcing the estimates toward the measurements. For a two-dimensional distributed model where “r” and “z” represent the spatial coordinates, the state estimation equations take the form: KW,

Z, t) =f(%

fi,, %,, Et,, g,,, r, Z, t)

+ G(r, Z, t) (y(t) -h(%(r*,

z*, t))). (7)

The matrix of gains G has m columns (corresponding to m measurements) which are determined by the differential sensitivity calculations. The gain is the product of the differential sensitivity evaluated at the measurement location, the partial derivative of the measurement operator, and the measurement weight matrix Q, i.e.: G(r, z, t)=

P(r, r*, z, z*, t)

The differential

sensitivity

ah(qr*,

z*, t)) =

iX(r*, 2*, t) > is obtained

Q(t).

(8)

by solving.

KP,(r, s, z, w, t)K=f&-, z)PK+f* r(r, z)P,K I,(r, ~)p,~K+f~(r, z)PK

KPf;(s,

w) +

KP,f$s, W)

+ KP,f@, 4 + KP,,f;Js, 4 meaS - C W, rZ, z, zZ, t)V,P(rZ, s, zZ, w, t) In

hb(v:>I)>f)l

1 dt

R(v, Y, t)R’(v, v’, t)dv=

R

P(v, v’, t) = E[%(V, t), %T(v’, t)]

+ KP,,f;_(s, w)+ KPf;(s, w)

v,o, 0

(5)

+ R(r, s, z, W, t1

(9)

where V, is related to the measurement weighting matrix and the measurement operator, i.e.

where:

s

By minimizing the least squares objective (5), we wish to obtain the optimal estimate f of the state vector x. If we let the error in the estimate be %=x -2, then the covariance of the estimate error

+

+ bX(x> xv, OR,(t)b,b, x,, ~1 +b:R,(t)b,

Q-’ _ covariance of measurement noise R -covariance of process noise R, I, R; ’ - covariance-of boundary noise.

+ fn2(7, z)P, K + fnsz(I, z) P,, K

v’, t)] dv dv’

k

[Yktf) -

2089

profiles

The symmetric positive definite weighting matrices 0, R, R,, and Ri have no strict statistical significance for this nonlinear case, but by analogy with linear problems having Gaussian white noise:

+f9

v, t)-JTR’(v, v’, t)

h,(t) - h(x(o:, t), WCW,

and concentration

16(v-v’).

(10)

2090

LARRY

c.

WINDES

The differential sensitivity equations are nonlinear due to the measurement terms. Note that large process error weight R increases P while large measurement weight (corresponding to small measurement error) decreases P. These general state estimation equations can be reduced to specific equations for fast states and slow states. This effectively involves reformulating the equations as s-+0 in K defined in eq. (3) in a similar way as Soliman and Ray (1979). The details of this general result may be found in Windes (1986). Packed

bed reactor

state estimator

The general state estimation equation was to yield two-dimensional state estimates packed bed reactor in our laboratory. The system, shown in Fig. 1, is designed to carry

applied for the reactor out the partial oxidation of methanol to formaldehyde. There is an undesirable second reaction step in which formaldehyde further oxidizes to CO and H,O. Careful control of the reactor temperature and concentration profiles provides valuable information for on-line optimization and control. The reactor system and model are described in detail in Windes (1986), Windes PCal. (1988), and Schwedock et al. (1988): thus we will provide only a brief summary here. The reactor modelling equations used by the state estimator are for a two-dimensional heterogeneous packed bed reactor model and have the form given below. Material co,

balances

on gaseous

species

yj (methanol,

etc.).

(&2+& 7nz(/I

I

g!Lj mr c

22 ar > -_(I --E~P~R;P+ 7)

(11)

where the reaction rate expressions, effectiveness factor calculations, catalyst decay, etc. which come into Rj are described in great detail in Windes et al. (1988) and Schwedock et al. (1988).

et

d.

tions at the wall are:

aT,

-F=Bi,,(Tw-T7:);~=0

while at: t-=0;

balance for the fluid

ai-,

ZT,

‘?yj

ar

?r

&

-~_---__-_~o.

(1) Partition

the differential sensitivity based on fast and slow variables. This scalar FL1 for solid temperature and PZ2 for the fluid variables. There are terms Prl and Pi’. Hence,

where (P”)r

= P12.

(17)

(2) Introduce

the pseudo-steady-state assumption which allows the dynamics of the fast states to be very rapid compared to the slow states. Thus, EPI’AO, &P”+O, and .sP” is redefined as Prz. In this way we track only the slow states and assume the fast states are always at their quasisteady state. (3) Insert the packed bed chemical reactor model into the general state estimation equations with

]=

----_;;__-[methanol]

[

Energy

balance

for

the solid

phase.

(13) The exothermic heat of reaction is removed by a circulating oil bath which maintains the tube wall at a constant temperature (r,,,), and the boundary condi-

equations involves a a matrix also cross

(16)

phase.

(12)

(15)

The heat transfer and kinetic parameters have been identified experimentally and through published correlations and data (cf. Windes et al., 1988; Schwedock et al., 1988). Note that all gas phase concentrations and temperature are “fast” variables in quasi-steadystate, while the solid phase temperature is a “slow” dynamic variable. To obtain the final equations for the state estimator, the following procedure was used:

+-Energy

(14)

Y

cc01

1

where the slow state is the catalyst pellet (solid phase) temperature, and the fast states are the gas temperature plus the composition of each independent species (methanol and carbon monoxide). The small parameter E is the ratio of heat capacities of gas and solid within the reactor. structure. Here, we con(4) Choose a measurement sider NTtemperature measurements T; at positions [r;, ~$1 and NY composition measurements ML(CH,OH) and C;;(CO) at positions [r;;, z;;].

Estimation

of

temperature and concentration profiles

(5) Derive the measurement

operator and determine V,. Here the temperature measurement is assumed to be a weighted average of the solid and fluid temperature: T;=ef$+(l -e)ff+ qr-. The measurement error is assumed independent of location so that:

where

(23) Carbon monoxide

e”QT e(l-e)Q,

V,=

e(l-e)QT

0

0

(~-c)~Q~

0

0

0

0

0

0

QM 0 0 Qc1

at

x [Tk-ef&.,

+F

k

{p^"( 23

+ P^$$(r,

z;. t)-(1 II rk, r, z:,

ri, z,

v,z;, z)Qr -e)Ff(r;,

z)QMCW-&L.,

zi)Q,.Ci

z;, t)]} 4,

- k(r2, z;, t)]}

01 (24)

where:

Solid temperature

afAr32, t) =f^l

concentration

aC(r, z, f) p-=*=f$+~{(l--e)PI:(r;, St k

This procedure leads to the following packed bed reactor state estimator equations:

__-

2091

+ “c’ { [eP’ yr, r;, 2, z$ t) k

+(I -e)Pi2(r,

G. 2, z$ KiQTCTk-efss(r$

-(I

A, 01)

-e)ff(&

A, t)

The state estimator boundary conditions are the same as those of the model. The differential sensitivities given by eq. (16) may be calculated from rather complex equations (cf. Windes, 1986). To illustrate the form, the equation for PI’ is given in the Appendix. Equations for the remaining sensitivities PI’, Pzl, Pz2 are tedious and omitted to conserve space. Complete results and the theoretical derivation may be found in Windes (1986).

t “c’{ p:yr, r;, z, z;r’,t) Q,,JM; - A&r;, z;, t)] k

IMPLEMENTATION OF THE DISTRIBUTED PARAMETER STATE ESTIMATOR In carrying out the real-time implementation of state estimation algorithm described in the previous section, there are two tasks:

FIuid temperature 8 TJ (r, z, t) at

=O=j-:

-e)Fff(r,

+?{(I

r;. z, z;)QT

k

x[Tk-ef&,, + y{pti(r,

2;. t)-(1-e)Tf(r;,

(1) to solve the differential sensitivity equations for P in order to calculate estimator gains G, (2) to solve, in real time, the state estimation equations for 8.

z;.t)]}

r:, z, z’JQM[MP -A&$,

2;. t)]

k +

sf:(r,

r;, z, zL)QJC;i - d(ri, zt, t)]}

(20)

where: ii=st,ci,-i,,Z+&Z

+&(?&r%). Methanol

(21)

concentration

aA(r, z, t)

NT

=O=.r^:+x

at

{(l-e)P::(r;,r,~;,z)Q~ k t)-(1

x[T;-eeTs(r;,z;,

-e)ff(r;,z;,t)])

NY +C{G(r,

4,

z,

zi)Q,,,[ML

-

J@ri,

z;, t)]

k +

P::(r, ri, z, z;‘)Q,[C;

-d

(r;, zl,

t)]}

(22)

We will first describe the numerical methods used to accomplish those tasks and then present some results for the packed bed reactor. Solution method For the packed bed reactor, the state estimates and differential sensitivities require the solution of the 14 partial differential equations presented in the last section. Applying orthogonal collocation to the spatial variables reduces the problem to the solutions of ordinary differential equations in time for the solid phase temperature T, and its differential sensitivity PI1 plus algebraic equations for the remaining variables. The solution of this large number of equations would be too timc-consuming for real-time, so the differential sensitivities are solved in advance off-line. The time-varying values of the differential sensitivities can then be represented by fitting the results to simple

2092

LARRY C.

WTNDES et al.

interpolation equations so that values can be quickly generated at any time during on-line solution of the filter equations. The procedure for solving the on-line filter equations is similar to the solution method for the model equations (cf. Windes et al., 19&7), except that values of measurements and estimator gains must be entered at each time step. In principle, the filter equations are of the form

Di&-ential sensitivity and estimation gain calculations The differential sensitivity calculations were used to provide estimator gains for each measurement. Because there are 4 states, there are 4 time and spatially varying gains for each measurement. For exampIe, the [solid temperatureltemperature measurement k] estimator gain at the collocation point [r,, z,] of the reactor model is: G r,lr;(r,,

2.))= (eP’ I(?,, &, zyr.&) +(l -e)PY(r,,

Thus the rate of change of the model states are adjusted by a measurement feed-back term which tends to force the measurement estimates toward the actual measurements. It was found in our work that 2 radial and 8 axial collocation points were adequate for solving the state estimator in real-time. To calculate the filter gain, one must have values of the differential sensitivities. Because the differential sensitivity equations are functions of the process state variables, solving them off-line requires the use of nominal values of the state variables in their solution. For the packed bed reactor studied here there was only a very small error introduced by this approximation. In facl, one set of nominal differential sensitivities were adequate for state estimation under a wide range of operating conditions. Although large shifts in the position of the steady state profiles caused performance to deteriorate somewhat, it was only in cases of significant changes in catalyst activity that the differential sensitivities had to be recalculated. The differential sensitivity equations are solved in two steps. The six sensitivity equations p*’ are solved first because they do not depend upon PI1 or PI’. In the second step, the equations for P1’ were integrated forward in time together with the three PI’ equations at pseudo-steady stare. Because the sensitivity equations for the packed bed reactor contain 4 spatial variables, 4 dimensional collocation methods were required in their solution. Thus values of P were calculated at the collocation points r., sg, zy, wti The number of axial and radial collocation points used is a compromise between accuracy and computational load. After reduction of the number of equations by expressing the boundary collocation values in terms of interior collocation points, the solution for a single differential sensitivity equation (e.g. for PI’ in the Appendix) with NR radial and NZ axial collocation points requires (NR - 1)’ x (NZ -2)’ equations. Tn our work, 2 radial and 9 axial collocation points were found to provide adequate accuracy for the differential sensitivities in all situations. For some cases, even one radial and 5 axial collocation points proved sufficient. In the interest of space, the collocation equations will not be derived here; however, a full detailed description may be found in Windes (1986).

rL zyi 4)Qp

(27)

These gains are stored at the collocation points, and the collocation interpolation weights allow the gains to be expressed at the collocation points of the estimator equations. When solving the differential sensitivity equations, several parameters must be selected by the system designer. These are:

(1) Number and location of measurements. (2) Initial values of P, ,(r, s, z, w, 0). Generally, only

steady state gains are used; therefore, initial P, 1 values are not critical. (3) Measurement error covariance (Q-l). It is assumed that all measurements are independent, and the variances are not a function of location. (4) Magnitude and spatial correlation of the process error covariance (R). The shape of the R matrices was selected to be: R(r, S,I, w)=p(l--(z-~w)~)“~(l

-(r-~)~)~~(rs)~~‘. (28)

This functional form was chosen to be symmetric and decreasing with distance between (T, z) and (s, w). The exponent “prl” determines the radial dependence of the process error. The values of these parameters used in the simulations are given in Table 1. To illustrate the spatial dependence of the differential sensitivities and estimator gains, four different realistic measurement situations were considered: (A) 1 bed ment. (B) 1 bed urement. (C) 1 bed (D) 2 bed

temperature,

1 methanol + CO measure-

temperature,

1 carbon monoxide

meas-

temperature measurement. temperature measurements.

For measurement combinations A, B, C, D above, the user-selected parameters are given in Table 1. The measurement errors assumed correspond to standard deviations of 5°C in temperature, 0.02 in methanol conversion, and 0.003 in methanol conversion to carbon monoxide. The model parameters were obtained through off-line parameter fitting of data taken when the catalyst had not yet significantly deactivated (cf. Windes et al., 1987; Schwedock et al., 1987). The base steady state for the differential sensitivity calculations was chosen as follows:

Estimation Table

of temperature

1. Parameters

Case:

A: IT;

B: 1T; 1c

0.005 0.0 0.002 1.0 0.5 0.0

feed and coolant temperature: feed mole fraction methanol: flow rate: 1.4 g/s

The

values

of differential

C:IT

of

sensitivities

D:2T

0.0 1.E4

1.E5 0.5 0.0025

0.5 0.0025 0.001 0.001 0.001 0.0 1 0.0 0.0 0.005 0.0 0.002 1.0 0.5 0.0

Ml01

0.00 1 0.001

0.0 1 0.0 0.001 0.0 2.5E-5 1.0 0.5 0.0

R given

2093

profiles

0.25; 0.7

0.25 0.8 0.0 1.0 1.E4

0.0

0.0 0.0

Note:

for calculation

1M and C 0.25 0.6 0.0 I.0 1.E4 2.5E3 l.E4 0.5 0.0025 0.00 1 0.001 0.00 I 0.01

collocation

and concentration

represent

the

mean

values

for all

points.

250°C 0.05

The differential sensitivities were solved using 1 radial and 5 axial collocation points. The number of axial collocation points was increased to 6, 7, and 8 points, and 2 radial collocation points were also tested. Increasing the number of collocation points in the differential sensitivity computation had only minor effects on the estimator gains. However, if a particular radial dependence of the process or measurement errors is assumed, then two radial collocation points is preferable. If there are three or more axial temperature measurement locations, then 7 to 9 axial collocation points are recommended. In case B, several combinations of weighting between CO measurement error and process error were tested. The case shown here assumes reasonably accurate measurements, small CO process error in R::, and larger R$* (correlation between temperature and CO process errors). This means the CO measurement has considerable effect on the temperature estimates, which in turn drive the CO estimate toward the measured values. The corrected temperature estimates have the added benefit of iniluencing the methanol estimates toward the true state if the production of methanol and CO is correlated. Since the differential sensitivity is analogous to the estimate error covariance for linear equations, the magnitude of the differential sensitivity is interpreted as an indicator of the quality of the estimates. Thus, for a given nominal steady state, alternative measurement locations can be compared for a given set of error parameters, and a lower overall diflerential sensitivity profile indicates a better measurement set. In addition,

regions where the differential sensitivity and thus the estimator gain tends to be large are found to be good choices for measurement locations. The radial location for temperature measurements was chosen at the tube centerline. Reasons for centerline sensor pIacement are: (a) the temperature-related differential sensitivities are the largest, (b) greater measurement accuracy is expected because the radial temperature derivative is equal to zero, (c) the temperatures are the highest, so centerline sensors focus on the hottest and most parametarically sensitive regions of the reactor. In the linear case, as trace P” becomes small, the estimate variances become small. The best measurement locations can be located through minimization of the trace of the covariance (e.g. Colantuoni and Padmanabhan, 1977). For the simplest case of 1 temperature measurement, the choice of axial location was made by observing the diagonal terms of the discretized P, 1 matrix. These terms are indicative of the variances of the solid temperature estimates as a function of position, From Figs 2 and 3, the “best” measurement location is shown to lie between z = 0.2 and z =0.3 because the resulting diagonal differential sensitivity profiles are the smallest, and the estimator gains are the largest at these points. Note that there is a sharp deterioration of performance if the thermocouple is placed in the downstream section of the reactor. For case A, the optimal location for a composition sensor usually involves some trade-off between the best position for methanol measurement and the best

2094

LARRY

I 0.0

C. WINDES

et al.

I

I

I

I

I

I

0.2

0.4

0.6

DIMENSIONLESS

AXIAL

I

1

I

0.8

1.0

DISTANCE

Fig. 2. Effect of axial measurement location upon differential sensitivity !‘I’. Curves denote values elements of P ’ 1and are computed with a single temperature measurement at: - ~~ z’ = 0.2;. ,f=o.3; - -. f=@,4; ~ - - z’zO.5; ~ Z’EO.7

of

diagonal

I

l-

:

I

0.6

DIMENSIONLESS

AXIAL

--A 1

1:o

DISTANCE

Fig, 3. Effect of axial measurement location upon estimator gain for the solid temperature. Curves are for a ~~ :‘=0.7. single temperature measurement at: ~ _~_--‘=0,2; ..,, ~‘~~.3;~.-.;‘~~.4;--~~‘~0.5;

position for carbon monoxide measurement. The best position for measuring methanol seems to be near the infection point of the axial methanol conversion curve, and the combined measurement was tocated at z”=O.6. Due to the experimental constraint that a long delay exists in obtaining methanol measurements from a gas chromatograph, case A cannot be implemented in real-time and is possible only in simulation. The carbon monoxide measurement is best made near the reactor exit, and the location z”=O.S was found to give good results. Rapid measurement of

CO is possible using an infrared detector, but case B was tested only by simulation. For multiple temperature measurements, a reasonably uniform differential sensitivity and overall gain profile is desirable while maintaining low interactions between measurements in the estimation equations. Figure 4 shows the differential sensitivity and estimator gains for thermocouple positions z’[ = 0.25 and z\=O.7. The estimator gain for each measurement position selected is the largest in the vicinity of its location.

Estimation

of temperature

and concentration

profiles

2095

T,:T I

0.0

0.2

I

-

0.4

Dimension

I

1

0.6

0.8

less Axial

Distance

Fig. 4. Differential sensitivities and estimator gains for two temperature measurements at r’= 0, z’= 0.25 and pxs. - _ PlZ. . 0.7. (A) Diagonal elements of PI’, (B) diagonal elements of PI’: ~Pi’, (C) estimator gains for the solid temperature equation: (1) for z’=O.2;, 12) for z’zd.7.

and 60cm. The estimator gains were calculated using one radial and 9 axial collocation points. These estimator gains will be used to obtain the experimental results for on-line estimation shown in later sections and also for control and optimization experiments to be reported in future papers. STATE ESTIMATOR

DIMENSIONLESS AXIAL DISTANCE

Fig. 5. Estimator gains with three centerhne temperature sensors. Solid temperature equation. Dots indicate sensor locations.

The performance of the state estimator applied to ’ our packed bed reactor system was tested in simulation and through on-line implementation to our pilot scale reactor. A schematic of the procedure is shown in Fig. 6. The following aspects have been investigated: (1) Influence of model accuracy-model

(2) In preparation for on-line experiments in state estimation and control, the estimator gains were calculated at a nominal steady state from a model based on off-line parameter fitting of open-loop experiments. The catalyst was substantially deactivated near the entrance of the bed. Only temperature measurements were used in on-line estimation due to the time delays for the composition measurement The estimator gains for solid temperature are shown in Fig. 5 for the three temperature measurement positions at 25, 40

PERFORMANCE

(3) (4) (5) (6)

parameter errors. Accuracy of collocation approximation-number of differential sensitivity collocation points. Number and location of sensors. Influence of user-selected filter parameters. Initial estimate of state variables. Effect of unmeasured disturbances.

The effectiveness of the state estimator is based on comparing the estimates to the actual or simulated data. The temperature at the measurement location, the maximum temperature, the centerline temperature

2096

LARRY C. WINDES rl ~1.

Simulation testing Using the estimator gains cakulated as described above, the state estimator was tested by simulating the reactor, generating “measurements” by adding random error to the simulation results, and using the artificial measurements as input to the state estimator. Two types of situations are simulated:

PROCESS INPUTS

STATE INPUl

TWX

YIELD

-

-

(1) Poor initial condition during start-up with various levels of error in the parameters’ (a) no parameter error, (b) erroneous heat transfer parameters, (c) erroneous kinetic paramctcrs (preexponential factor). Unmeasured disturbances in operating conditions: (a) wall temperature disturbances, (b) inlet temperature disturbances.

IMERPOL4TION

ZATT,

Fig. 6. Schematic for state estimation.

(2) profile, and the compositions at the measurement location and the exit are examined. In the case of simulation, the true states are completely known, and the comparison is straightforward. In the case of reactor data, numerous thermocouple readings are recorded, and the temperatures are compared with the results of the state estimation based on a few selected measurement locations. The results are equivalent for on-line and off-line state estimates, but on-line implementation is constrained by computing time and requires software to gather and process the data in real-time. 325

For test case la, an important criterion for the performance of the estimator is the speed of convergence to the correct temperature profile. The initial temperature profile estimate is 50°C too high. The temperature estimates must decrease rapidly toward the actual temperatures, and then they must increase as the reactor temperature rises during the start-up simulation.

values

Temperature

are given

profiles

in Fig. 7 and

and

composition

illustrate

the

,=,.5 min

A T

250

/ 0.0

I

14

I

I

28

AXIAL

/

I

I

42

I

56

1

70

DISTANCE (cm)

TIME (min)

Fig. 7. Simulation--convergence of temperature estimates given poor initial guess. Two temperature measurements at z’=O.25 and 0.7; initial condition: estimator, 300°C; actual, 250°C; --.~ ~~~~actual, estimates. (A) Centerline temperature profiles. time: 0.5 min, 1.0 min, 1 5 min, (R) exit methanol conversion; exit carbon monoxide.

rapid

Estimation of temperature convergence downstream

of the estimates. of the hotspot

composition)

greatly

and concentration

profiles

2097

A second measurement (either

improves

temperature

the estimates

of the reactor. test case la, one expects that would eventually converge to the actual

or

in that

section For

the estimates states because

the model

is perfect; however, even in test case 1b (heat too low), good conversion, selectivity, and temperature estimates are obtained when the model has a biased error (cf. Fig. 8). For case lc (when there are kinetic parameter errors and only temperature transfer

measurements

are availahle)

the estimated

tempera-

ture profile converges towards the temperature urements, but the correct reactor temperature gives

incorrect

composition compensate

conversion

and

selectivitics.

measprofile Thus,

measurements are needed to accurately for poor kinetic data. Figure 9 shows an

example of estimates in which the mbdel assumes kinetic constants (preexponential factor) 10% too large. The second test case situation is given by a series of unknown changes in the wall temperature, feed temperature, or flow rate. The estimator assumes the reactor is operating at T,,, = T, = 25o”C, and the initial conditions

are for the corresponding

case 2a, the measurements are from reactor operation

provided

steady

state.

0.0

14

42 AXIAL

DISTANCE

5b

70

(cm)

Fig. 8. Simulation-estimation with errors in initial guess and model parameter. Radial fluid heat transfer in estimator 15% too low. Two temperature measurements at z’=O.25 and 0.7; initial condition: estimator, 300°C; actual, 25D”C; actual states, - - - estimates of centerline temperature profiles.

In

to the estimator

with fluctuating wall temin Fig. IOF. The estimator must perform even when the wall temperature provided to the estimator is constant. A single temperature measurement at Z’ = 0.25 gives good estimates of the changing reactor temperature (cf. Fig. IOA-E). Since the temperature measurements are made at the centerline perature

250

as shown

and the disturbance

is propagated

from the wall, there

0.0

14

28 AXIAL

42 OISTANCE

56

70

(cm)

is a delay in the temperature response when compared with the upsets in conversion/CO production. Therefore, a CO measurement greatly helps to speed the response during the rapid dynamics of this test case. Step changes in the inlet temperature comprised the second test of estimation with unmeasured disturbances in the operating conditions, 2b. This case is more difficult than having wall temperature disturbances because the hotspot shifts axially. In order to correctly track this shift in the temperature profile, more than one temperature

measurement

cause the temperature

must

be available

be-

increases in one portion of the in other portions. The reactor

reactor and decreases model with uniform catalyst activity was used in the estimation simulation with three temperatures measurement locations at 25, 40 and 60cm. The inlet temperature disturbance to the reactor and the tracking of maximum reactor temperature is shown in Fig. 11. The estimates of the changing temperature profiles (Fig. 12) are excellent except for the first 1&I 5 cm of the reactor. Also, as shown in Fig. 13, the exit composition from the reactor is correctly tracked. This test conclusively illustrates that the estimator can trace

shifting

temperature

three thermocouples. It is obvious that

lations,

profiles

with dnly two

the differential sensitivity even when carried out off-line, require

or

calcuexten-

TIME

(min)

Fig. 9. Simulation+stimation with errors in initial guess and model parameter. Kinetic constants in estimator 10% too high Initial condition. estimator, 300°C: actual, 250°C. Step change in methanol feed mole fraction from 0.05 to 0.04 actual states, ~ one temperature at time = 2 min; ~ measurement at r’=0.25 and one CO measurement at z” =0.8, two temperature measurements at z’=O.25 and 0.7. (A) Centerhne temperature profiles, (B) exit methanol conversion, (C) exit carbon monoxide mole fraction. sive computational effort to generate the optimal estimator gains. Thus a study was carried out to determine the performance of an observer with estimator gains generated in an ad hoc fashion. The two ad hoc procedures were:

2098

LARRY

C

WINDES

et al

Fig. 12. Simulation-temperature measured disturbance

Fig. 10. Simulationestimation with unmeasured disturbances in wall temperature. One temperature measurement at z’ =0.25; ~ actual state; 0 L L measurements;,estimates. (At(D) Centerline temperature profiles, time: (A) 0.5 min, (8) 1 min, (C) 2 mln, (D) 3 min. (E) Temperature at r=O, z=O.25. (F) Unmeasured disturbance in wall temperature.

=o+iF 3.80-O

-~

1s

io

$5 _ -

30 B

actual -

5

10

T

15

I ME

c8timates

20

30

(min)

Fig. 11. Simulation-estimation with unmeasured in Inlet temperature. (A) Feed temperature ance pattern. Estimated T,=25O’C. (8) Estimates mum temperature.

ances

25

E

Fig. awes

SO330+ 0

t

disturb-

disturbof maxi-

(1) use of constant estimator gains with the same mean as the optimal ones, (2) use of estimator temperature gains shaped like the reactor temperature profile.

0.98

METHANOL

profile estimates

with un-

in inlet temperature.

,080

A

CARBON

MONOXIDE

13. Simulation--estimation with unmeasured in inlet temperature. (A) Methanol conversion (B) fraction CO production.

B

disturb at exit

Figure 14 shows the resulting estimated reactor tern perature profiles for test case la with no measurement errors added. These results show the optimal gains to be significantly better than the ad hoc procedures even though the rough magnitude of the gain was known. For more numerous measurements of different types (i.e. temperature and composition), guessing the appropriate combination of gains becomes much more difficult than for the single temperature measurement, and detail calculations are even more helpful. Although “tuning” the estimator to achieve adequate results is possible for these simulation cases, the actual states are not available for a real process, and the true performance of the estimator is difficult to assess.

Estimation 325

of temperature

,

and concentration l

0

14

0

i4

42

28 AXIAL

DISTANCE

56

i”

56

70

(cm)

Fig. 14. Simulation-estimation using ad hoc estimator gains-- temperature convergence given poor initial guess. (Compare with Fig. 7.) Initial condition: estimator, 300°C; actual, 250°C; ~~- actual, - - - estimates. (A) Constant gains of 10, (B) gains with axial shape matching the steady state temperature profile.

Thus, the optimal estimator gains seem worth the computational effort. Based on the results of simulation testing, quite good estimates are possible for most estimation situations using optimal estimator gains and a single thermocouple located somewhat upstream from the reactor hotspot. However, an additional measurement of either downstream temperature or carbon monoxide near the exit enhances the quality of estimates considerably. A second measurement of carbon monoxide concentration is preferable to a second temperature measurement, especially if there are rapid fluctuations in yield. For situations when the hotspot shifts, such as during flow rate or inlet temperature changes, multiple

temperature

sensors

are

necessary.

To

ac-

comodate these disturbances, modelling errors, and catalyst deactivation which gradually shifts the temperature profile, three temperature measurements are recommended. In comparing different choices of estimator gains, “optimal” estimates (from the differential sensitivity calculations) give better results than ad hoc observers with arbitrarily selected estimator gains. Real-tima experimental testing State estimation experiments

with

a piIot

scale

packed bed reactor were used to demonstrate the online capabilities of the state estimation algorithm.

2099

profiles

Since no signals are fed back to the reactor, the same results could be obtained by collecting the data and running the program off-line. However, the program must function on-line and produce state estimates in real-time for use in parameter estimation and control. Thus in our work, on-line state estimation was carried out during reactor start-up, with disturbances in the operating conditions, and with poor initial conditions and model parameters. The overall strategy for experimentation follows the rationale given in the previous section for simulation testing. In these experiments, however, only reactor temperatures are used as measurements. Temperature measurements were made at six second intervals, and the state estimator had no difficulty in carrying out the real-time computations necessary to produce state estimates at each time step. A PDPl l/55 minicomputer was used for implementing the state estimator as well as data logging. Not all of the temperature measurements taken in the experiment are utilized in the state estimation algorithm. Thirty reactor temperature measurements are taken, but only l-3 of these are utilized in the state estimation. The remainder are used to indicate the actual temperature profile in the reactor which is compared with profile generated by the state estimator. Thus, the advantage of this highly instrumented reactor is that the performance of the distributed parameter state estimation algorithm is evaluated not only at the measurement locations for the estimator, but also at other “observed” locations. An overview of the five state estimation experiments (El-E5) is shown in Table 2. These experiments demonstrate performance for:

(1) Factors

(2)

in design of the estimator: (a) changing numbers of temperature measurements (El, E2), (b) poor heat transrer parameters (E3), (c) poor initial conditions (El, E2, E3). Factors in reactor operation: (a) start-up conditions (El, E2, E3), (b) known change in feed mole fraction methanol (El, E2, E3), (c) known change in total flowrate (E4), (d) unknown disturbance in wall temperature (E5).

Experiments El, E2 and E3 follow the same experimental procedure: start-up after 1 min with 5% methanol feed, change from 5% to 3% methanol feed after 7 min. Experiment E4 involves beginning at steady operation and performing equal and opposite step changes in the reactor flowrate during a 20 min time period. Ramp changes were made in wall temperature during experiment

E5. The wall temperature

measure-

ment was not given to the state estimator, and wall temperatu-e changes represented unknown disturbances. All dynamics of the state estimates for experiment E5 are measurement

driven

by the

one

temperature

LARRY C.

2100 Table

2. State estimation

rt

ai.

experiments Mcas. Lot. (cm)

Comments

0.7 0.7 0.7

25 25 25

initial est. temp. = 300°C reactor start-up feed mole fraction change

0.00 0.05 0.03

0.7 0.7 0.7

25, 40, 60 25, 40, 60 25, 40, 60

initial est. temp. = 300°C reactor start-up repeat El

260 260 260

0.00 0.05 0.03

0.7 0.7 0.7

25, 40, 60 25, 40, 60 25, 40, 60

initial est. temp. = 300°C radial heat transfer 25% low repeat El

O&3.3 3.3-l 3.7 13.7-20.0

260 260 260

0.05 0.05 0.05

1.4 0.5 1.4

25, 40, 60 25, 40, 60 25, 40. 60

initial cond. at S.S. step decrease in flow step increase in flow

O-25

260’

0.05

0.7

25

Expi

Time (min)

T,.=T, (“C)

El El El

&1 l-7 7-l 2

260 260 260

0.00 0.05 0.03

E2 E2 E2

@l 1-7 7-12

260 260 260

E3 E3 E3

61 1-7 7-12

E4 E4 E4 E5 +The actual

WINDES

wall/feed

FlOW (g/s)

Cl-I&I

temperature

Trn,,,

unmeasured

temp. disturbance

makes ramp changes between 252 and 268°C.

located. Profiles 2-5 show snapshots of the dynamic response of the estimated profiles compared with the experimental profile4 after reactor start-up at 1 min. The estimates follow the reactor profiles adequately and are especially good in the upstream half of the reactor. The estimates in the downstream half of the reactor are high during the lirst 3 min but improve as the reactor approaches steady state.

0-l -,min 0 - 2 - l.Srnl” d - 3 - Jrnl” q -•-4mtn x7-S-smlmln

‘--‘--‘._*_____~

--_..__*_-o__.,

“i5Oj

,

.

,

I

.

~_..~_~__&-..---b._-_-__

,

,

.

,

,

0

A X $I ZJx

I A L

0

I 5”;

A

N C

\

E

(cm)

,

70

1

Fig. 15. Stale estimation experiment El-reactor start-up. Comparison between measured and estimated temperatures. One temperature measurement at z=25 cm. Initial condotions: reactor, 260°C; estimates, 300°C. (A) Temperature profiles, (B) temperature at the sensor location.

Reactor start-up with one temperature measurement The results of experiment El are presented in Fig. 15. The state estimator used one temperature sensor located at the reactor centerline and 25 cm from the reactor entrance. An initial guess of a uniform reactor temperature of 300°C in the state estimator converged rapidly toward the true ractor conditions of 260°C during the first minute. Profile 1 of Fig. 15 indicates that convergence is most rapid in the upstream half of the reactor where the measurement is

Reactor start-up with three temperature measurements Experiment E2 repeats experiment El with three temperature sensors at 25, 40 and 60 cm. The performance of the estimator, illustrated in Fig. 16, shows significant improvement over the case of one temperature measurement. The estimator converges to the data within 1 min, and the estimates track the entire centerline temperature profile accurately through the entire reactor start-up. The time response of the estimates is shown at each measurement location, and the dynamics agree well with the data. The composition profiles were also estimated (cf. Fig. 17) based on the estimated temperature profiles and the reactor model even though no composition measurements are made during these short time spans. Estimation with model parameter errors Experiment E3 repeats experiment E2 with incorrect values of the radial fluid heat transfer coefficient (Pe,,J and catalyst heat capacity (Lewis number). The radial heat transfer was chosen 25% lower than the best available value obtained from off-line parameter fitting. The Lewis number was set 20% low compared with the value from off-line data analysis. This error caused a dynamic response faster than the real reactor. The estimated temperature profiles shown in Fig. 18 give reasonably good results even for these rather severe parameter errors. However, the hotspot region of the reactor increased too fast and is about 5’ too hot at steady state. The response at the sensor locations (Fig. 19) indicates good estimates of

Estimation

of temperature

and concentration Estimation

-lOWtkll

AXI

T I M E

(min)

A L

A N C E

D I S”;

T

IYE

(min)

12

(cm)

‘““m’ T E I H

(mm)

Fig. 16. State estimation experiment E2-reactor start-up. Comparison between measured and estimated temperatures.

Three temperature measurements at z = 25,40, 60 cm. Initial conditions: reactor. 260°C; estimates, 300°C. (A) Temperature profiles, (B), (C), (D) temperatures at the sensor locations.

the reactor dynamics, but the steady state results show more steady state offset than for experiment E2. The most important comparison to be made is between the estimator and the model alone in absence of state estimation. The estimator gives greatly improved performance over the model when there are significant modelling errors. In the case of errors in the Peclet number as in this experiment, the steady state error in the temperature estimates was reduced by 15-20” when compared with the model above, and the improvement in the dynamic response was dramatic.

E-D START-UP

I

EL.wsm llME AFIER AT TIME = 1 min

SVRT-UP

0.5min

70

AXIAL DI&&CE

44/10-H

changes

6min

0

CES

flow

0.06 ,

TIMEAFTER

1 -

Fig.

during

Experiment E4 tests the estimator for the more difficult dynamic responses to step changes in the total flow rate. The estimates converge quickly to their steady state values before the first step change. The dynamic responses are well-estimated, but some steady state otiset exists. The estimates tend to be too high near the reactor entrance and too low near the reactor exit. This characteristic of the model is not completely compensated for by the estimator. The temperature profile estimates compared with the reactor data for a decrease in flowrate are shown in Fig. 20. The hotspot of the reactor moves upstream by -20 cm, but the temperature changes very little at 42 cm. This temperature invariant location is seen in both the data and the temperature estimates. The estimator does an excellent job of determining the shifting temperature profile. Figure 21 shows the temperature profile estimates for an increase in the flow rate. Since the final flow rate is much larger, the dynamics are much faster than for the flow decrease. The reactor hotspot shifts -20 cm downstream, the upstream temperatures decrease rapidly, and a rapid increase in the downstream temperature lags somewhat behind. The region at 40 45 cm shows only small changes in temperature. All of these dynamic features are nicely represented by the estimated temperature profiles. The composition profiles for flow rate changes were also estimated, but no composition measurements were made during these rapid transients. An important aspect of the estimation of the temperatures is the correct estimation of the maximum temperature in the reactor. The maximum temperature estimate must show realistic dynamic behavior even if the model is in error. One of the worst deficiencies of the constant parameter reactor model used for estimation and control is inaccuracies in representation of substantial flow rate changes. The temperature maximum as determined from the data, the estimates, and the reactor model alone are shown in Fig. 22. The estimates show all of the essential

AT TIME = 1 m’m

4 -

2lOl

profiles

17. State estimation

experiment

(cm) E24stimated

70 ’

AXIAL

composition

DI&NCE profiles

(cm)

during

reactor

start-up.

2102

LARRY

C.

et al.

7b

0

AXIAL

DIS’TANCE

(cm)

Fig. 18. State estimation experiment E3-reactor starl-up. Comparison between measured and estimated temperature profdes. Three temperature measurements at z=ZS, 40: 6Ocm. Radial heat transfer parameter erroneously low. Initial conditions: reactor, 26O’C; estimates, 300°C.

B

WINDES

2

.

6

T I M E

e

A X I A L

Cl

I S-7 A N C E (cm)

Fig. 21. Stateestimation experiment E4k- --flow change from 0.5 to 1.4 g/s. Comparison between measured and estimated temperature profiles. Three temperature measurements at z = 25. 40, 60 cm.

I@

(min:

Fig. 19. State estimation experiment E3--reactor start-up. Comparison between measured and estimated temperatures al the sensor tocations. Three temperature measurements. Radial heat transfer parameter erroneously low. Initial conditions. reactor, 260°C; estimates, 300°C.

Fig. 22. State estimation experiment E4 flow rate disturbances. Maximum reactor temperature estimates with three temperature measurements.

Estimation

during

wall temperature

unmeasured

dis-

turbances

Experiment E5 demonstrates state estimation in the case of the unmeasured wall temperature disturbances shown in Fig. 23B. The estimates at the sensor location track the measurements nicely while the model alone in this case exhibits a constant biased output. The impressive feature of these results is the speed at which the estimator responds to the rapid changes in the reactor while driven entirely by the measurement. Not only is the temperature at the measurement location estimated correctly, but also the entire temperature and composition profiles are adjusted based on A X

I A L

D I S’;

A N

C E

(cm)

Fig. 20. State estimation experiment E4a-Bow change from 1.4 to 0.5 g/s. Comparison between measured and estimated temt)erature

profiles.

Three temperature z = 25, 40, 60 cm.

measurements

at

the

gains.

Figure

24 gives

the tempera-

at six different times during temperature profile estimates

the experare good,

especially at the reactor hotspot. The exit composition estimates in Fig. 25 show an oscillating behavior directly corresponding to the changes in the true reactor

features of these two flow changes, while the model exhibits incorrect dynamic responses and is especially poor for the flow increase. The steady state agreement with the data is also much better when the state estimator is functioning.

estimator

ture profile iment. The

temperature.

dynamic

Summary und conclusions

Two-dimensional, nonlinear, distributed-parameter state estimation has been applied in real-time packed

bed chemical

reactor

system

to a

for the first time.

Estimation

ol temperature

pnd concentration

profiles

TIME (min)

2103

TIME (min)

Fig. 25. State estimation experiment ES---exit composition estimates for unmeasured wall temperature disturbances. One centerline temperature measurement at z = 25 cm. (A) Fractional methanol conversion at exit, (B) carbon monoxide produced as fraction of methanol in feed.

T I M E

(min)

Fig. 23. State estimation experiment ES-unmeasured wall temperature drsturbance. One temperature measurement at z = 25 cm. (A) Comparison between measured and estimated temperatures at the sensor locations, (3) wall temperature

system because of the close coupling between solid and fluid temperatures and the fact that no separate measurements are made of solid and gas temperatures. The differential sensitivities were solved off-line using four-dimensional orthogonal collocation which produced estimator gains for the state estimate equations Some simplifications were found to be appropriate for the differential sensitivity equations:

disturbance.

(1) fewer discretization points than for the model, (2) a nominal steady state, (3) time-invariant estimator gains during continuous reactor operation. The choice of the nominal steady state was not critical so long as it was reasonable, but the measurement locations were important. Good measurement locations can be selected by examining the differential sensitivity profiles of candidate measurement sets. The recommended temperature measurement locations are:

Fig. 24. State estimation experiment E5--unmeasured wall temperature disturbance. Comparison between measured and estimated temperature profiles. One temperature meas-

urement at z = 25 cm.

These algorithms form the cornerstone of a control strategy for the packed bed reactor to be described in later papers. The derivation of the state estimation algorithm yielded two sets of partial differential equations: (1) state estimates and (2) differential sensitivities. A key element is the use of the full twodimensional nonlinear model in the state estimate equations. This procedure utilizes all of the power of the process model in order to estimate the temperature and composition profiles from the available measurements. These results were based on a two-dimensional heterogeneous model, but a pseudo-homogeneous model would also be adequate for this particular

(1) slightly upstream of any anticipated reactor hotspots, (2) as far downstream as any expected hotspot, (3) at the hotspot of the nominal operating point. Based on the results of simulation testing, remarkably good estimates are possible for most estimation situations using a single thermocouple located somewhat upstream from the reactor hotspot. However, an additional measurement of either downstream temperature or carbon monoxide neat the exit enhances the quality of these estimates. Measurement of carbon monoxide is preferable&o a second temperature measurement, especially if there are rapid fluctuations in yield. For situations when the hotspot shifts, such as during flow rate or inlet temperature changes, multiple temperature sensors are necessary. To accomodate such disturbances, substantial modelling errors, or catalyst deactivation (which gradually shifts the temperature profile), three temperature measurements are recommended.

LARRY C. WINDER ef cd

2104

The state estimation algorithm gave good results in tracking of the reactor behavior during simulation and on-line experiments. The entire two-dimensional temperature and concentration profiles were estimated using substantially less than real-time for online computations. The estimation of the maximum reactor temperature was reliable, and dynamics of reactor operating changes were accurately represented. The estimator rapidly converged when given a poor initial guess and withstood the effects of erroneous model parameters or model mismatch with the experimental system. In comparing different choices of estimator gains, the “optimal” estimator gave much better results than ad hoc observes with arbitrarily sclectcd estimator gains. Additional desirable features which could be easily included in the state estimator would be the incorporation of delayed composition measurements [cf. Ray (1985)] and the simultaneous estimation of uncertain parameters such as the radial heat transfer Peclet number along with the state variables. NOTATION

a b BK c Da E e f G H h K M Mr Mz m NT NY n JQJr, Per, I%,, P%ZZ Pe,,, Pesr P a R R,, R, R,

Ri

length to radius ratio of reactor nonlinear boundary operator dimensionless reactor temperature rise dimensionless mole fraction carbon monoxide dimensionless Damkdhler number, ratio of reaction rate to flow expectation operator measurement weighting for solid temperature nonlinear function estimator gain linear measurement matrix nonlinear measurement operator diagonal matrix containing Is and cs (a small parameter) dimensionless mole fraction methanol number or radial measurements number of axial measurements number of measurements number of temperature measurements number of composition measurements number of states radial and axial Peclet numbers for fluid heat transfer radial and axial mass Peclet numbers radial and axial Peclet numbers for solid heat transfer differential sensitivity matrix weight for measurement error weight matrix for process model error weight matrix for boundary error effective rate of reaction of key species which releases heat equivalent to all reactions occurring: R, + (AH,/AH,)R, dimensionless rate of reaction of species j; j= M, C; based on reactor volume

Rj.k RJ.Ts i-

s ‘-h

T TS t “, V

u’, VII w X

Y 2

derivative of rate of reaction of speciesj w.r.t. species k derivative of rate of reaction of species j w.r.t.

T,

radial coordinate dummy radial independent variable Stanton number, dimensionless interphase heat transfer temperature solid temperature time measurement operator multidimensional independent spatial variable vector dummy independent variables dummy axial independent variable state variable vector measurement vector axial coordinate

Greek letters Dirac delta function 6 a small parameter, the ratio of fluid to solid E thermal capacitance measurement noise process or model noise :! Y objective function 52 domain of functional space Subscripts derivative

with respect

F

fluid

A4 I s

derivative derivative derivative

T t

temperature derivative w.r.t. time

to carbon

monoxide

with respect to methanol w.r.t. radial coordinate w.r.t. s

derivative w.r.t. u derivative w.r.t. w derivative w.r.t. state variable x derivative w.r.t. axial coordinate z a, 8, y, 6 collocation indices V

W

Superscripts temperature r, composition

measurement measurement

estimate 0, * T

1 2

error general

measurement

location

transpose of a matrix slow states fast states REFERENCES

Ajinkya, M. B., Ray, W. H. and Froment, G. F., 1974, On-line estimation of catalyst activity profiles in packed-bed reactors having catalyst decay. I & E C Proc. Des. Deu. 13,

107-112. Ajinkya, M. B., Ray, W. H., Yu, T. K. and Seinfeld, J. H., 1975, The application of an approximate non-linear filter to

systems governed by coupled ordinary Int. J. Systems SC;. 6,

cntial equations.

and partial 313-332.

differ-

Estimation

of temperature

Balchen, J. G., Fjeld, M. and Olsen, T. O., 1971, Multivariable control with approximate state estimation of a chemical tubular reactor. Proc. IFAC Symp. on Multivariable Systems, Dusseldorf, Paper 51.1. Bencala, K. E. and Seinfeid, J. H., 1979, Distributed parameter filtering: boundary noise and discrete observations. Inc. i. Systems Sci. 10, 493-512. Clement, K., Jorgensen, S. B. and Sorensen, J. P., 1980, Fixed bed reactor Kalman filtering-experimental investigation of discrete time case without stochastic disturbances. Chem. Engng Sci. 35, 1231-1236. Colantuoni, G. and Padmanabhan, L., 1977, Optimal sensor locations for tubular flow reactors. Chem. Engng Sci. 32, 1035. Cooper, D. J., Ramirez, W. F. and Clough, D. E., 1986, Comparison of linear distributed-parameter filters to lumped approximants. A.1.Ch.E. J. 32, 186. Goldman, S. F. and Sargent, R. W. M., 1971, Applications of linear estimation theory to chemical processes: a feasibility study. Chem. Engng Sci. 26, 1535-1553. Hwang, M., Seinfeld, J. H. and Gavalas, G. R., 1972, Optimal least-square filtering and interpolation in distributed parameter systems. J. Math. Analysis Appl. 39, 49-74. Joffe, B. L. and Sargent, R. W. H., 1972, The design of an online control scheme for a tubular.catalytic reactor. ~rans. Inst. Chem. Engrs 50, 27@282. Jorgensen, S. B. and Clement, K., 1977, Experimental investigation of reaction zone control in a fixed bed reactor. From Chemical Engineering with Per Soltofi, pp. 69-85 (Edited by K. Ostergaard and A. Fredenslund). Teknisk Forlag, Copenhagen. Jutan, A., Wright, J. D. and MacGregor, J. F., 1977, Multivariable computer control of a butane hydrogenolysis reactor: Part Ill. On-line linear quadratic stochastic control studies. A.1.Ch.E. J. 23, 751-758. Kumar, S. and Seinfeld, J. H., 1978, Optimal location of measurements in tubular reactors. Chem. Engng Sci. 33, 1507-1516. Kuruoglu, N., Ramirez, W. F. and Clough, D. E. 1981, Distributed parameter estimation and identification for systems with fast and slow dynamics-a tubular, fixed-bed catalytic reactor to form styrene. Chem. Engng Sci. 36, 1357-1364. Kuruoglu, N., Ramire? W. F. and Clough, D. E., 1982, Use of a distributed mini/micro computer system for tubular reactor profile and catalyst activity identification. AIChE Annual Meeting, Paper 1 lc, Los Angeles. Lausterer, G. K., 1977, State estimation and optimal control of distributed parameter systems, PhD Thesis. SUNYBuffalo. Lausterer, G. K., Ray, W. H. and Martens, H. R., 1978, Real time distributed state estimation applied to a two dimensional heated ingot. Automatica 14, 335-344. Lausterer, G. K. and Ray, W. H.. 1979, Distributed parameter state estimation and optimal feedback-control-an experimental study in two space dimensions. IEEE Trans. Auto. Control 24, 179-190. MacGregor, J. F. and Wong, A. K. L., 1980, Multivariate model identification and stochastic control of a chemical reactor. Technometrics 22. Pell, T. M. and Aris, R., 1970, Problems in chemical reactor analysis with stochastic features--control of linear distributcd systems on discrete and corrupted observations. I & E C Fund. 9, 15-20. Ray, W. H., 1978, Some recent applications of distributed parameter systems theory-a survey. Automaticn 14, 281-287. Ray, W. H., 1981. Advanced Process Control. McGraw-Hill, New York. Ray, W. H., 1985, Polymerization reactor control. Proc. 1985 ACC. Rutzler, W., 1980, Nonlinear estimation methods for tubular reactors. PhD Thesis, University of Minnesota. Sawarogi, T., Soeda, T. and Omatu, S., 1978, Modelling,

and concentration

2105

profiles

Estimation, and Their Applications fbr Disrriblrted Parameter Systems. Springer-Verlag, Berlin. Schwedock, M. J., 1983, Modelling and identification of a catalytic packed bed reactor. PhD Thesis, University of Wisconsin. Schwedock, M. J., Windes, L. C. and Ray, W. H., 1988, Steady state and dynamic modelling of a packed bed reactor for the partial oxidation of methanol to formaldehyde. II. Experimental results compared with model predictions. Chem. Engng Commun. in press. &infield, J. H., 1969, Non-linear estimation for partial differential equations. Chem. Engng Sci. 24, 75-83. Seinfeld, J. H., Gavalas, G. R. and Hwang, M., 1971, Nonlinear filtering in distributed parameter systems. J. Dynamic Syst. Meas. Control 93, 157-163. Seinfeld, J. H. and Hwang, M., 1970, Some remarks on nonlinear filtering in distributed parameter systems. Chem. Engng Sci. 25, 741-743. Silva, J. M., Wallman, P. H. and Foss, A. S., 1979, Multi-bed catalytic reactor control systems: configuration development and experimental testing. I & E C Fund. 18, 383-391. Soliman, M. A. and Ray, W. H., 1979a, Nonlinear state estimation of packed-bed tubular reactors. A.I.Ch.E. J. 25, 718-720. Soliman, M. A. and Ray, W. H., 1979b, Non-linear filtering for distributed parameter systems having a small parameter. Int. 1. Control 30, 757-772. Sorenson, J. P., 1977, Experimental investigation of the Optimal control of a fixed-bed reactor. Chem. Engng Sci. 32, 763-774. Sorensen, J. P., Jorgensen, S. B. and Clement, K., 1980, Fixed bed reactor Kalman filtering and optimal control---l. Computational comparison ofdiscrete vs. contmuous time formulation. Chem. Engng Sci. 35, 1223-1230. Tzafestas, S. G., 1978, Distributed parameter state estimation, Chap. 3, pp. 135-207, Distribyted Parameter Systems (Edited by W. H. Ray and D. G. Laintiotis). Dekker. Vakit, H. B., Michelsen, M. L. and Foss, A. S., 1973, Fixed bed reactor control with slate estimation. Z & E C Fund. 12, 328-335. Wallman, P. and Foss, A. S., 1981, Experiences with dynamic estimators for fixed bed reactors. I & E C Fund. 20, 234239. Watanabe, K., Yoshimura, T. and Soeda, T., 1980, Optimal non-linear estimation for distributed pirameter systems via the partition theorem. Int. J. Systems Sci. 11, 11131130. Windes, L. C., 1986, Modelling and control of a packed bed reactor. PhD. Thesis, University of Wisconsin-Madison. Windes, L. C., Schwedock, M. J. and Ray, W. H., 1988, Steady state and dynamic modelling of a packed bed reactor for the partial oxidation of methanol to formaldehyde. 1. Model development. Chem. Engng Commun. in press. Yu, T. K., 1973, Optimal filtering for systems governed by coupled ordinary and partial differential equations. PhD Thesis, California Institute of Technology. Yu, T. K., Seinfeld. J. H. and Ray, W. H., 1974, Filtering in nonlinear time delay systems. IEEE Trans. Auto. Control 4, 324-333. APPENDIX The equation $(P’Y

for differential

r, s, 2. w, 0)=2B,{R,,.&r,

-St,{4P”

-2P:+,

+ 2B,{R,,(r,

sensitivity

P’ ’ is:

z)+R,,.,(s,

s, z, w, t)-2P;qs,

w)}P” r, w, z, t)}

z)P:*(s, r, w, 2, tJ

+ R,,&,

z)P:‘(s,

+R,..&,

w)P:‘(r,

r, w, z, t) s. z, w, r)+R,&-,

w)P:‘(r.

s, z, w, t)}

LARRY C. WINDES

2106

+${IJ::+P;:: $!‘+L:l+P:;L+P;s ,ri

++

3

NT

-c

z, 4,

OQTP"(r;, s, z;, w, t)

k

+e(l -e)P”(r, +e(l --e)P:‘(r,

r:,r, zf,t)QTP;*(s, r;. w. zi, t) r:,z, z:,t)QrP1’(r;I, s, z;,W, t)

--e)*P:2(i-,i-1,z, z:,r)QTP:2[S,r;, w, z;, t)}

+(1

NY NY

-c

NT

1 {e2P”b-, r’,, z

11

et ul.

2 {P:‘(r, r:‘, z, z:‘, i

+P:‘(r,

r:-‘, z, z;‘,

+R'l(r s 3,.

Other

r)Q,P:‘(s.

r;. w. z;‘, t)

x

z

differential

w,

t)Q&'(s. r;. w, z;, r)}

t).

sensitivities

have a similar

structure.