Nonlinear model predictive control of a fixed-bed water-gas shift reactor: An experimental study

Compurers Printed

Vol. 18, No. 2, PP. 83-102,

them.

Engng,

in Great

Britain.

All

rights

1994

reserved

Copyright

0

0098-I 354/94 $6.00 + 0.00 1994 Pergamon Press Ltd

NONLINEAR MODEL PREDICTIVE CONTROL OF A FIXED-BED WATER-GAS SHIFT REACTOR: AN EXPERIMENTAL STUDY T.

G.

of Chemical

Department

(Received

16 November

f99.7; jinal

and T. F.

WRIGHT

Engineering,

The University

revision received

EDGAR

of Texas

2 June

at Austin,

Austin,

TX

78712,

1993; received for publication

U.S.A.

I6 June 1993)

Abstract-This paper describes new results on the experimental application of nonlinear model-predictive control (NMPC) to a fixed-bed water-gas shift (WGS) reactor. The development and experimental and how it impacts controller validation of an appropriate first-principles WGS reactor model, performance is discussed. The implementation of NMPC is computationally intense, requiring that a large nonlinear program (NLP) be solved at each sampling period. The significant computational burden

dictates that a relatively slow sampling rate be used. Infrequent sampling, however, diminishes disturbance rejection capabilities. To combat this problem, NMPC was implemented in a master-slave cascade

configuration where a low-level liner controller, having a significantly faster sampling rate, was employed. The control study was performed using a PC-based distributed control system (DCS) One of the three processors was dedicated to NMPC calculations. A complete and rigorous implementation strategy is described in the paper, and the performance of NMPC for set-point tracking of this nonlinear process is shown to be superior to adaptive or linear control. We also illustrate the ease with which NMPC accommodated feedforward control.

1.

algorithm

INTRODUCI-ION

may be difficult to achieve.

deficiency, The

severity

of the nonlinearities

cesses influences for successful of nonlinear ing point design.

the selection

control

is a standard

is constrained

to

nonlinear may

linear

perform

conditions

such

poorly.

of

continuous

of

batch

wide

also

system

for

range

start-up

back.

highly

this

of or

way

trajectory

prove

difficult

operating

One technique an adaptive and past

that attempts

associated feedback

operating

controller

can

technique

linear

process change

states

that

applications

In general,

thought

of

has

distinctly

two

variables

the

model-based

feedback states

the resulting

first

such

to

employ

advanced

been

via state or output

feed-

formations

required

a nonlinear

temperature

or

takes

that,

state

under

dynamic

system

control

of

successful

of

in

The

predictive

seek

often

computer models

A control

fail to broadly

capabilities for

real-time

strategy

power

control

strategy

objective

horizon

which

and

model-predictive

of an open-loop

a finite

this

a more

a repeated optimization over

was

While

level of robustness

of the increased

is nonlinear

which

strategy.

rigorous

and control.

advantage

(NMPC).

speed

control involves

performance

extending

from

the

current time into the future as done for linear model-

have been reported

control

use

was

the necessary trans-

for implementation

improvements

the

ex-

(1978)

ef al. (1983).

we

an

output

becomes

methodology,

control

transinto

or

Brockett

appealing,

Hence,

nonlinear

of computers

on very

concentration) many

error.

optimization

different

the parameters

model

involves system

this

is intuitively

permit

is the

have

a nonlinear

by Hunt

technique

Recent

an adap-

further

strategies

theory

linearization

of

feedback,

to

for which

control

geometric

system

applicable

uses

with time, whereas the

faster rate. While successful

as

operating

Presumably,

of adaptive

in the literature,

for the

law based upon current

vary slowly

(e.g.

at a much

to compensate linear controllers

control be

time scales.

the process model

with

conditions.

types of time-dependent different

controllers,

exist, or do not yield the desired

inadequacies

control

Global

further

control.

tive

years,

differential

eqivalent

the

tracking for

in nonlinear

Because of this

be a need

actly linear in a specified manner.

shut-down

and

recent

forming

chemical

in

to

used to effect linearization

operation

but

designed

in

processs

processes

plant

as nonisothermal

The

encountered

In

appears

for this research.

employing

when the process

region,

controllers

motivation

operat-

for control

when

a narrow

systems

reactors,

or

improvements

The linearization

is adequate

mild

pro-

algorithms

about a nominal

approach

This approach are

of control

of a process.

physical models

nonlinearities

in chemical

there

predictive

with an adaptive

et al., 83

control

1987;

(Cutler and Ramaker,

Eaton

and

Rawlings,

1980; Clarke

1991).

Rawlings

G. T. WRIGHT and T. F. EDGAR

84

and Muske lizing an

(1991)

infinite

equality

model

prediction

problem

programming plicit

and

using

Linear

NMPC

extends

with

primary

develop

objective

an advanced

fixed-bed

reactor

imentally.

Section

laboratory

scale

with computer

nonlinear

and

2 describes

are classified

with

maximum

to nonlinin an ex-

5, a first-principles for

egies

model

is developed.

A

parameters

are

Section

7, NMPC

pared

is evaluated

to more

traditional

suit-

control

strat-

and

the

reactor

to NMPC

issues,

and in Sec-

experimentally

and com-

control

FACILITY PROCEDURE

The

water-gas

shift

reaction

production

of

ammonia,

chemicals.

The

reaction

AND

better

implement

Catalyst

on a similar

used

a sulfur-tolerant

that is much (Fig.

1) consisted

processing

arises

hydrogen is

and

reversible

in the organic

and

mildly

exothermic: + H@(g) AH,, In typical feed

hydrogen, sulfur

Carbon

steam

A

is 40%

CO, ratio

condition

that

the

exceeds

four

depending

dioxide,

hydrocarbons dry-gas

reactor

insignificant

to dry-gas

satisfies

carbon

of

nominal

assuming

hydrocarbon, The

monoxide,

monoxide,

40%

CO,

and

may

vary,

the

steam

upon

effluent

Hr.

WGS

reaction

quired. tween

high

is run in either

reactor

conversion

Steam

quench

or in multiple of

carbon

streams

beds in a multi-reactor

configuration

bed behavior

varying

the reactor

dry-gas

ratio.

in the dry-gas

feed

stream

temperature

quately

rejected.

monoxide

are often

is primarily may

composition variations

and

be-

In either

influenced

by

or the steam

to

include that

is re-

located

configuration.

inlet temperature

Disturbances

was designed

sation

prior

obtained

to entering

from

fluctuations

flow, cannot

and

up-

be ade-

hydrogen

and

pressure-

controllers

(MFC).

building

header.

with

manual

the supply

gases as

A

3 m sections

water

The

tubing

A K-type

the assembly ture. The

heater

eral wool

was

from

the reservoir. in two-parallel The

parallel heater. 1.25 in.

packed with

with high

magnesium thermal

con-

was also placed

the cartridge

temperature

generator

I

in a 30.5 cm long,

thermocouple

input

20

a 750 W cartridge

insulator

to monitor

the power

of the steam

around

water, in

diaphragm

was vaporized

was housed

an electrical

ductivity.

stored

of 16 in., 3 16 stainless-steel.

was wrapped

conden-

Deionized

was

a

steam

quantity

to avoid

Chem-Tech flow

to

The

a known

the reactor. tap,

similar

(1991).

it enough

a building

deionized

lating

was

to vaporize

was used to regulate

The

2.3.

system

reservoir.

tube.

in series

gas

and the

from

from

used to direct

and to superheat

dia.

reactors

dioxide,

in conjunction

generating

assembly

adia-

reactor

to the reactor

used by Bell and Edgar

This

a single

carbon

valves

steam

to

composition

feed

the wet

Steam generator

oxide,

batic fixed-bed when

2.2.

tubing

ratio

gas

system.

via mass flow

were

but typically CO

generator,

air was obtained

valves

pump

goals. The

bypass

feed

20%

reactor

dry

desired.

polyethylene

of

(1991) catalyst

WGS

the

the fixed-bed

supplied

isolation

and

quantities

cata-

Previous

Edgar

the

parts:

the steam

cylinders

of water process

for

five

system,

were

generator the dry

shift

used.

and

to

a well-

2. I. Dry gas feed processing

generator

kcal/mol.

quantities

impurities.

composition

+ H,(g)

applications

carbon

small

CQl(g)

= -9.8

industrial

contains

-

was

the

chemistry

to characterize.

facility

gas processing

The

CO(g)

approx.

strategies,

by Bell

difficult of

system,

processing

effluent

of

cobalt-molybdenum

experimental

MFC (WCS)

reactor

more

the others.

shift catalysts

iron-oxide

Cl2-3-05)

work

used have

was not to study

control

high-temperature,

The

than

but to use it as a model nonlinear

cobah-

catalysts

temperatures

understood, (United

oxide

detail

goal

lyst

Compressed

OPERATING

greater

and

the catalysts

as high-temperature operating

reaction,

regulated 2. THE EXPERIMENTAL

Iron

Since our principal

nitrogen

strategies.

in much

WGS

feed

is adopted,

6 is devoted

and implementation

a

3, 4 and

reactor

technique

estimated

of

data and

In Sections

solution

is validated.

for a

reactor

of acquiring the WGS

to

exper-

shift

model-based

model

development

strategy

the construction

for

model

was

strategy

water-gas

control.

use in real-time,

research

500°C.

oxide

are among

shift conversion.

studied

the

capable

advanced

to promote been

control

apply

fixed-bed

facilities

implementing able

to

copper-zinc catalysts

They

MPC

this

oxide,

the optim-

constraints

of

Iron

molybdenum

using quadratic

permit

directly

incorporates

models

manner.

The

tion

stabi-

controller

constraints

to be solved (QP).

systems

a nominally

predictive

horizon.

and inequality

ization ear

have developed

constrained

was varied

to the generator. was insulated

in

heater temperaby manipuThe with

exterior

2 in. min-

insulation.

Wet -gas feed processing

After

mixing

the dry gases with the steam,

gas mixture

was heated

to a temperature

165°C.

temperature

was maintained

This

the wet

of approx. using

PID

Nonlinear model predictive control

k

W8tW

Fig. 1. Schematic of the experimental facility for the water-gas shift reactor. control

to

mitigate

upstream

temperature

ations caused by the steam generator. 0.75 m of 0.25 in. tubing

fluctu-

Approximately

ran from the steam/dry-gas

mixing tee to the reactor inlet. A heat exchanger constructed ceramic

by wrapping

fiber-insulated,

this tubing 1.9 n/f

with 3.65 m of

nichrome

wire. Tem-

perature at the reactor inlet was controlied power

to the resistance

heater. This

was

by varying

assembly

com-

prised the feed preheater. 2.4.

The

reactor

The body and

was constructed

(0.035

in.).

(39 in.). taken

A

reactor.

was positioned

using a

total

3.175 cm

wall reactor

length

temperature.

The

0.089 cm

was

99.06 cm

points

were

along

the

heat

loss,

heater extended thermocouple and

extended

thermo-well tor. Power

were

two from

independent around

resistance

the reactor.

the inlet header

to the fourth

pair. The second continued to

the

last

heaters The first

thermocouple

from there pair.

The

housed

were in a

at the top of the reac-

heater was supplied

in a

to the steam generator heater. The

surrounding

the thermo-well

and extending

was packed

Pyrex glass beads. These facilitated

flow distribution

to

gas entry

served as catalyst

into

the active

support

bed.

They

when loading

with also

the reactor

bed. The active bed of the reactor was 36.6 cm long and extended

from

the seventh.

the first thermocouple

The

active

The void

The

cm/O.125

to particle diameter

fraction

exit section

of the reactor

normal

operation

rated plate located couple

ratio of approx.

Some

were

used

bed was 0.35.

was packed

and rested

the catalyst upon

at the reactor exit. Four

pairs extended

boundary

with iron-

in.), yielding a

in the catalyst

kaolin clay spheres. They supported during

pair through

bed was packed

oxide catalyst pellets (0.3175 reactor diameter

thermocouples

inlet section

heater,

to the first pair of thermocouples

9.4.

of

located

to the cartridge

manner analogous volume

cartridge

assembly

wall

first pair

tightly

using a 300 W

of a

11.43 cm (4.5 in.) from the top of the reactor

were wrapped

heated

consisted

and the last at 72.39 cm (28.5 in.). In order to minimize

Insulation.

feed gases and the reactor

and a corresponding

pair of measurements temperature

321

of

measurements

at 11 equally-spaced

Each

vertically

(1.25 in.),

thickness

pair of temperature

centerline-bed located

with

The

axially

The

prior

of the reactor

stainless-steel

reactor was then insulated with several layers of thick Cerablanket

to

conditions

into

the reactor

examine discussed

bed

a perfothermo-

exit section.

the validity later. The

section was 51.05 cm (20.1 in.) in length.

with

of model total

exit

G. T. WRIGHT and T. F. EDGAR

86 2.5.

Efluent

Upon

exiting

through water

gus processing the reactor,

a series vapor

effluent

of

content

dry-gas

the

of the gas and to determine

the

composition. effluent

This

coupled

information

for

each

species

was

Infrared

A

bus

DCS

with

16 MHz

functionally

were

necessary provide

for

The

available

platforms.

software

DMACS

ing easy

including

execution,

The

uler which

time

second

station

stations

execution

were

device

capable the

are, in essence,

of accomplishing

host

stations

computer

an RS-422/485 were 2.7.

Reactor

those

intelligent

for

suggested

used

from

shut-down. performed

by the catalyst

Total

inlet

reaction

catalyst

bed When

the reducing

process

typically

gas flowrates

similarly An

for

control.

performed

in depth

experiments ranging

from

from

10 to by

temperature

This

rejecting

is

were performed

inlet

experiments.

facility

and control

ranged

experiments

active-bed

method

reaction

collection

temperatures

value using PID

and replicating

proved

at

a

to be an

upstream

disturbances

Closed-loop

experiments

using

discussion

a cascade

control

of this is given

later

in the text.

3. DYNAMIC

of

MODELING

which The of

Optomux over

brain boards

catalysts

modeling

manufacturer,

United

fixed-bed

must be given

ranging

from mass

siderations

WATER

GAS

level of model

the

availability

mation

for

among

the

such

con-

of phenom-

intra-particle partial

simple

and

These

con-

differential

reaction

schemes.

constrained

physical

property

products,

accurate

of catalyst

by

inforrate

characteristics,

things. model

is

to

in the modeling

reactor

complex,

and

and knowledge

other

to

is ultimately

reliable

reactants

expressions

for

detail

of

reactors,

transport.

complex

even

The

how

energy

to

models,

catalytic

to a multitude

fluid and

lead

equation

How

were

OF THE FIXED-BED REACTOR

SHIFT

bed the

entire

The

465°C.

and closed-loop

at inlet

the

interphase

operation.

treat

that

stabilized,

data

Open-

Steady-state

constraint to

through-

as

approx.

of the WGS

by start-up,

were

ena

multiplexers,

series

was con-

15-20 h.

and

effective

gas

to limit the rate

certain

The

a

was

to reaction

downward

ultimately

operation

16 SLM.

process

phases

exceeded

was complete.

Optomux a mounting

operation

procedures

procedure

sideration

the host computer

multidrop

the temperatures

wet

of

flowrate

in the bed to 65°C.

to make

never

Steam

was exercised

adjusted

proceeded

temperatures

scheme.

tasks independent

The

with

was

pro-

computer.

serial link. The Optomux

configured

The

system.

communicated

slowly

and anlog

is a controller

most

increase

user

a sched-

Each

to a host

care

of temperature

When

board

Special

a 1:2 dry-gas

was used

digital and

The

reach-

by intro-

consisted

temperature

to

continuous

independent

board

CO.

SLM. bed

Upon

with

to

mixture.

was initiated mixture

Hz,

and reducing

were

interval

using

employed.

brain

at 68

(300°C).

Sev-

intervals.

rack. The Optomux boards

fixed

of COI,

gas

6-8 SLM.

dry-gas

state. was

ambient

steam

gas mixture

The

out the heating

desired

allow-

code.

and

and control,

operates brain

ratio.

I mixture

ditions

to

iron a gas

its oxidized

it from

activation

UC1

procedure

by heating

a wet process

steam

the

employing

from

was typically

maintaining

on several

techniques

to spawn

The

executing

execution,

of a brain

as a slave

flowrate

270 to 300°C.

reaction

architecture for

to

a commer-

user-written

of these

was designed

consisted

Steam

Normal

Ethernet

together.

and was implemented

data acquisition

Optomux

coaxial

the catalyst

a 2: 3 nitrogen

categorized

adaptor-

computers

available

provided

as user-specified

For

the

DMACS,

for

event-based

predominantly grams

were

specific

execution.

Ethernet

package

access

mechanisms

tasks

intensive

has an open

database

using

took

optimizations

for the WGS

FIX

for

for highly

the computers

Intellution’s

was

33 MHz

Thin

the DCS

im-

The

of

communication.

was

150°C

of

entailed

in the activation

the catalyst

temperature

several

storage,

as the

in each

compat-

machine

interface.

such

32-bit and one

among and

calculations.

installed

used to drive

facility

step

used to raise the

architecture

16 MHz

collection

was used to wire

software

prepare

to

traditionally

used as platforms

of NMPC were

eral

tasks

computer

data

Intel

This

host

computations

boards

first

Activation

basically

to reduce

The

ducing

control

IBM-PC

the

computers.

numerical

deter-

computers

research.

and as an operator

machines

cially

33 MHz

80386/80387

this

primary

trending

cable

two

divides

for

medium

(UCI).

ing this temperature

balance.

distributed

compatible

in

catalyst

Cl2

maintained

by a single

autonomous used

flowrates

to completely

by material

architecture,

16-bit

used

plemented

mass

Inc.

oxide

2: 2:

IBM-PC

Intel

was

inlet

sufficient

employing

80386/80387 ible

were

equipment

global

system

analyzers

CO and CO2 compositions.

the exit composition Digital

2.6.

designed

gases passed

to eliminate

used to determine

mine

the product

units

Catalysts,

model

would

but also unsuitable

as optimization

and

be

used

effort.

is A

prove

to

be

for on-line nonlinear

another

rigorous not

key fixedonly

applications

control.

What

Nonhnear

model predictivecontrol

follows is a discussion of the assumptions used to develop a sufficiently simple, yet accurate, model for real-time implementation. The WGS reactor constructed for this study was designed to operate under nearly adiabatic conditions. This was achieved by adding guard heaters to the reactor to minimize radial heat loss. The reactor was also insulated with several inches of fiberglas insulation. While these efforts did not eliminate heat loss entirely, an adiabatic assumption was employed for model simplification. Although fixed-bed reactors are heterogeneous systerns with both fluid and solid phases, it is often reasonable to assume that the mass within a volume element can be characterized by a single bulk temperature, pressure and composition. The pseudohomogeneous assumption is a valid approximation provided composition and temperature gradients between the fluid and solid phases are small. This situation prevails when reaction resistance is large relative to mass and heat transfer resistance. Windes e? nf. (I 989) compared one- and two-phase models for the oxidation of methanol. They concluded that qualitatively these models compare favorably under most circumstances. They further concluded that even if the pseudo-homogeneous assumption were not strictly applicable, the one- and two-phase models compare well quantitatively with some parameter adjustment. Bell and Edgar (1991) employed this assumption in a WGS reactor model for a system similar to the one constructed for this research. Their results confirmed that assuming homogeneity is a practical simplification yielding good results for the WGS reaction. The work of Ampaya and Rinker (1977) and Bonvin (1980) further supported this conclusion. The spatial dimensionality of a fixed-bed reactor model may profoundly affect the model’s capacity for accurate prediction. For small diameter reactors running under adiabatic conditions radial gradients can often be ignored, but for nonadiabatic exothermic reactions where radial gradients can be large, failure to model the radial dimension may render the model useless. As stated earlier, the reactor in this research was designed to minimize heat loss and thus large radial gradients as done by Bell and Edgar (1991). A 1-D model was therefore developed, which also minimized the number of states upon discretization of the distributed parameter system. The simplifications mentioned thus far had the greatest impact upon the size of the model and were implemented primarily for this reason. Other assumptions described below were based solely upon physical considerations. Because the reactor was operated at low pressures, the process gases were assumed to

87

behave ideally. In addition, pressure drop across the reactor was small, obviating the need for equations describing this effect. The residence time for the process gas was on the order of 1 s. This was small compared to the time-scale for changes in catalystbed temperatures. Therefore, the quasi-steady state assumption that concentration changes are instantaneous relative to temperature changes was adopted. Dispersion is negligible when the ratio of reaction length to particle diameter exceeds 100 (Carberry and Wendel, 1963, Rase, 1977). It was included here for numerical conditioning as suggested by Windes (1986). The final model consisted of two partial differential equations with axial and temporal dimensions. Danckwerts (1953) boundary conditions were used at the reactor inlet for the catalyst bed balances, and zero gradient boundary conditions were applied at the exit. The dimensionless model is presented below. Model symbols are defined in the Nomenclature. 3.1. Carbon 0=

monoxide

balance

-$+i[+]-Da(-i,,) pe,

3.2. Catalyst

(1)

bed energy balance

z= hI df2

af a?= Le x=z+&

1

- %w@

-

fw) + fit-fco).

(2)

The WGS reactor model was nondimensionalized using the following definitions: f=$, rsf

i=Z

L’

‘=I. rcf

The reference time t,, was chosen to be the “residence time” based upon the initial gas velocity L/v,. The reference temperature was taken to be the reactor inlet temperature TO . These variable definitions led to the dimensionless groups given in Table 1. The dimensionless boundary conditions were: af

= Pe,(? ai

i =o:

3Y, ca2 i=

1:

= Pe,(Yco K?Yco aidi

-

fO),

(3)

- YO),

(4)

-O_

(5)

G. T. WRIGHT and T. F. EDGAR

88

model,

resulting

states. DAEs integration

since

conditions initial

conditions

boundary

condition

employed

the reactor exit for the above equations approximation However, cannot support,

from

becomes

the

similarly

nondimensionalized

and

profiles, Rinker

Newsome

expressions

special

that

Lee

(-rco)

in

two

(1962), Bonvin

Many

rate

have been proposed,

but

for

second

was given

reaction

equilibrium

Inc.

was employed.

- JJc,,-YH*I&l(~)l

has the advantage

of directly incor-

the effects of steam. This is especially

useful

if the flow of steam to the reactor is adjustable. rate constant temperature

to

order rate expressions

Catalysts,

= kt&J.Yrilo

This expression

were

Moe

(1980),

consideration

account

by United

so

and Bell and Edgar (1991)

for the reaction

The following

porating

resulted

an array of shift catalysts.

expressions

systems

conditions

including

(1977),

(1980),

in this research

provided

the reactive

bed to the inert

one for each equation.

of investigators

have studied

As

the effluent concen-

initial

and

effects.

Thus,

behavior.

k was assumed

The

to have an Arrhenius

dependence.

and

set of in

initial

consistent

a way

dynamic

that

start-up.

using an explicit DAE,

may be initialized

but

similarly:

k = f(x, Y),

(6)

0 = g(x, Y).

(7)

states

of

the

DAE,

x were

given

values xt,, determined by some initial (perhaps arbiso that trary) profile. y0 was then determined equation

(7) was satisfied.

the differential

even when it

and for adiabatic The

A number (1980),

ceases. fixed,

temperature.

dimensionless Ampaya

verified.

the active catalyst

reaction

tration

observed

assumption

be experimentally

gases move

does

to experimentally

this is a common

at

is a numerical

determined

is outlined

differential

14 algebraic

a consistent

steady-state

DAEs

and

difficult to initialize for

In this research,

were

both

The technique implicit

The zero gradient

finding

is nontrivial.

permitted

The

in 12 differential

are notoriously

were

determined

course, ant,

Having

and algebraic to

permitted

the same

satisfy

dynamic

technique

state start-up.

How

determined

states, equation

start-up.

(6).

This,

Equally

was employed

best to choose

both

the derivatives of

import-

for steady-

x0 was addressed

as follows. For the WGS

reactor model,

the differential

corresponded

to spatially

temperatures,

which were measurable.

tributed,

radially

distributed

averaged

prised the algebraic

bed

the reactor

mate

model

were used

temperatures

nonlinear

collected

to determine

at the spatial

nodes

when implementing

equation

sets that

the above procedures

a subroutine

from

libraries

et al.,

Caracotsios’

the MINPACK DASAC

was used to integrate the model in time and to

evaluate to

1980).

arise

were solved

using HYBRD, (1986)

via

applications.

algebraic

(More

axially approxi-

These were used to initialize the

for experimental

The

The seven equally-

measurements

linear interpolation. model

discom-

states, but only the composition

temperature

along

Spatially

compositions

at the reactor exit was measurable. spaced

states

centerline-bed

the

state

and

sequence

DASAC

output of

is based upon

gration

algorithm,

sensitivities

manipulated

with

respect

variable

moves.

the predictor-corrector

DASSL,

developed

by

intePetzold

(1983). 4.

The

model

section finite

were

TECHNIQUES

equations solved

element

piecewise-simple suitably

SOLUTION

chosen

developed

numerically

technique.

A

polynomials partition

model

solution.

equation

scheme.

were

to spatially

respect

Twelve

of

to some con-

to the true

ultimately (DAE)

set,

led to which

predictor-corrector

piecewise linear elements

discretize

FSTIMATION VERIFICATION

AND

MODEL

a Galerkin

combination

approximation

using an implicit,

integration used

with

This approximation

a differential-algebraic was integrated

using

linear

of the spatial domain

stitute a finite-dimensional

5. PARAMETER

in the previous

the WGS

reactor

For the WGS to be fitted

reactor model,

consisted

l

E,, activation energy. A, pre-exponential factor.

l

CP‘Y

l

These poorly

heat capacity parameters known

the parameter

of the following

set of

parameters:

of solid medium.

were

and could

chosen not

because accurrately

they

were

be deter-

mined a ~riori. The first two of these parameters

are

model predictive

Nonlinear

however, was how best to choose T, for optimal conditioning of the estimation problem. The parameter estimation problem was solved for several values of T,,, . The value that was ultimately employed yielded the smallest 2 - d intervals for the parameter estimates as determined by GREG (Caracotsios, 1986), a parameter estimation package developed by Caracotsios. A weighted least squares cost function was used as the measure of plant-model mismatch in this research:

clearly kinetic rate parameters and the last was used primarily for fitting reactor dynamics since it appears only in the Lewis number which is the coefficient for the time derivative of dimensionless temperature. The need to estimate heat transfer parameters was removed by virtue of the adiabatic assumption which eliminated the reactor wall energy balance. This assumption proved to be quite reasonable since heat losses to the surrounding were small relative to the heat generated by reaction. Seven temperature measurements located axially along the centerline of the active bed were used for parameter estimation. The measured effluent CO composition was also employed initially, but subsequent studies proved the parameters to be insensitive to this value when used with the multiple temperature measurements. The model states were most sensitive to the rate parameters on the interior of the active bed for high-temperature operation and at the exit for lowtemperature operation. The state sensitivities with respect to the activation energy and the pre-exponential constants varied roughly in proportion to one another. Linearly dependent sensitivities lead to virtually dependent first-order necessary conditions for the parameter estimation problem and an ill-posed estimation problem. For this reason the “centering” technique described by Bates and Watts (1988), which improves the conditioning of the estimation problem, was employed. The Arrhenius expression:

The cost function is equivalent to the maximum likelihood estimator when the measurement errors are uncorrelated and normally distributed and their variances are constant. Since these assumptions do not apply to our data, we are content to interpret the results simply as weighted least squares estimates. Because each temperature measurement was assumed to be similarly accurate, the weights u, were each given a value of unity. For measurements of varying qualities, however, the weights can be- adjusted to reflect measurement confidence. The weights can also be used as scaling factors for measurements of different magnitudes. This was unnecessary here since the equations and the data were scaled via nondimensionalization. The parameters were estimated using eight steadystate data sets and three dynamic data sets. The dynamic experiments consisted of perturbations to the active bed inlet temperature usually in the form of first-order exponentially filtered steps. Flowrates ranged from approx. 9.6 to 12.0 SLM, and the inlet conditions varied as indicated in Table 2 for the steady-state experiments. Experiments rb1028a,b, and c were performed to verify reproducibility of results from several months earlier. Parameter estimates obtained from steady-state data for the activation energy and the pre-exponential constant are given in Table 3. Also listed are the 2 - cr intervals associated with each parameter estimate. The 2 - D interval is a simple measure of the quality of the parameter estimate-small values relative to

was rewritten as

, where A’=Aexp(-&).

The mean temperature T, was chosen to lie within the range of observed bed temperatures. The primary effect of centering was to reduce the collinear dependence between the sensitivities. What was not clear, Table

Experiment rbO7231a rbO723lb rbO7241a rbO724lb rbO724lc rb10281a rbl028lb rb10281c

ID

2.

Operating

conditions

Inlet tcmperaturc “C 285.00 280.00 28 I .24 285.00 290.00 281.24 285. I7 290.00

89

control

for

steady-state Inlet

mol

experiments fraction

co

H>O

CO1

HZ

Total flowrate SLM

0.154 0.154 0.167 0.167 0.167 0.165 0.165 0.165

0.534 0.534 0.588 0.588 0.588 0.585 0.585 0.585

0.156 0.156 0.167 0.167 0.167 0.166 0.166 0.166

0.156 0.156 0.084 0.084 0.084 0.084 0.084 0.084

9.6 I 9.6 I 11.95 11.95 11.95 I i.97 11.97 Il.97

G.

90 Table

3. Optimal

parameter

estimates

from

steady-state

Estimated value

Parameter

E. T, = 297°C SE of residuals

WRIGHT and T. F. EDGAR

experiments

2-u Interval

1.055 x 10-s 2.783 x IO+’

A’

T.

4.258 x lO-7 1.189 x lo+’

7.199X 10-3

the estimates are preferred. These intervals are strictly valid only if the parameter estimates are independent and normally distributed. Figure 2 illustrates the good experimental data/model agreement for three of the eight steady-state experiments. The sample given represents low; medium- and high-temperature operation. As stated earlier the heat capacity of the solid phase was estimated to obtain a good dynamic fit. The rate parameters were also re-estimated using the steady-state parameter estimates as initial guesses. Table 4 lists the parameter estimates obtained when the dynamic experimental data were employed. Because the rate parameters varied only slightly from the values obtained using steady-state data and since all parameters were well determined, we may conclude that the good dynamic fit illustrated in Fig. 3 for the seven equally spaced axial centerline bed temperature measurements was primarily achieved via the solid heat capacity estimate. Note that the dynamic data led to smaller 2 - o intervals for the kinetic parameters. Figure 3 shows experiment rb7241d. The remaining dynamic experiments behaved similarly. We may also conclude that the model was valid over the nonlinear operating space given by the conditions in Table 2.

6.

CONTROLLER

DEVELOPMENT

There are two ways of performing model-predictive control calculations. The first method is sequential and employs separate algorithms to solve the differential equations, and carry out the optimization. First, a manipulated variable profile is guessed, and the differential equations are solved numerically to obtain an open-loop variable profile. Based upon the numerical solution, the objective function is evaluated. The gradient of the objective function with respect to the manipulated variable is determined either by finite differencing or by solving sensitivity equations. Finally, the control profile is updated using some optimization algorithm, and the process repeated until the optimal profiles are obtained. This constitutes a sequential solution and optimization strategy, and recent versions of this strategy have been reported by: Asselmeyer (1985), Morshedi (1986), Jang et al. (1987), Kiparissides and Georgiou (1987) and Peterson er al. (1989). The availability of accurate and efficient integration and optimization packages permits implementation of this method with little programming effort. However, constraint handling is poorer than in an alternative method which uses a simultaneous solution and optimization strategy. When the second or simultaneous approach is adopted, the model differential equations are discretized, and along with the algebraic model equations are included as constraints in a nonlinear programming (NLP) problem. The optimization of the objective function is performed such that

460 model model model

440

420

E 3

z

P p1

380

-

rb7241a rb724lb rb7241c rb724la rb7241b rb7241c

----. ----l

+ D

t

360

340

t

320

I

I

0

0.2

Fig. 2. Steady-statemodel predictionsand

0.4 Normalized

Axial

0.6 Length

experimental observations rb724 lc.

0.8

.

1

for experimentsrb724I a, rb724I b,

Nonlinear model predictive control Table

4. Optimal

parameter

estimates

Parameter A’

dynamic

Cp* 7-m=297”C SE of residuals

3.210 9.953 6.123

7.282

the discretized

model

differential

isfied and other constraints variables

are

this method

have

Asbjornsen

(1977),

Biegler

(1987),

(1988,

1989,

Renfro

Having

are sat-

reported

by Hertzberg

and

In this work,

a sequential

for upon for

experimental the the

dimensionality

two

equations taneous directly

arising

from

optimization

simultaneous

and

of the sampling

NLP

number

sufficiently

large, the computational horizon

associated

with

sequential

the

strategy.

the computation

is extended than

intense

integrations

varies

an integer order than

required

when

is smaller

for

for the simultaneous

path strategies

is

the

of

objective

infeasible.

Another

solving

the NLP

ant for real-time the sampling

path

advantage

output

of feasible

the luxury

suboptimally.

Experimental Model

Data Data

path

some fixed time, optimization

----.

380

t 360

340

320

300

280 0

Fig. 3. Dynamic

model

50

100

150

Time

predictions and experimental experiment number

200 250 (minutes)

of

usually

If the time limit is approached

in the course of solving the NLP,

2

H d

tech-

import-

400

::

fail,

since the solution

420

-z

can be

of intentionally

460

440

the

is usually

This becomes

implementation

interval.

at

strategies

since the model

must not exceed

if the

need only comfunction

the controller

niques is that one enjoys

offers

is that

to the value at failure.

infeasible

there is no direct recourse

the

to op-

When

While

to be computation-

fail, the controller

of the optimization

either

strategy.

first of these

the

sol-

but the

path technique

values

each NLP

sequential

appear

should

the prediction

approach

path solution

optimizer

If this value improves,

and

equations),

pare

the

while seeking

strategy can employ

The

but

simul-

a feasible path technique

the feasible

advantages.

implemented.

that

in

infeasible

beginning

the relative increase in

time required

horizon

arose

burden of a large

more

Moreover,

strategy

strategy (PH),

or infeasible

several

constraint

If the model

prediction

becomes

which of

a feasible

strat-

on the contrary,

optimization

of the model

optimization

ally more efficient,

in the simul-

solution horizon

time.

sol-

this choice

predominantly

discretization

with the prediction

multiple

the

and

(at least in terms

sequential

path

be satisfied

at each iteration,

at each iteration

et al.

Rawlings

The

and sol-

the constraints

path strategies,

Patwardhan

optimization

The

Feasible

ution strategy is necessarily

based

of

strategies.

and satisfy

and

We justified

application

be solved)

equations

Cuthrell and

Infeasible

(the model

satisfy the constraints

( 1990). ution strategy was employed.

strategy.

that the constraints

the optimum.

Eaton

optimization

for a feasible path approach

and

(1984),

et al. (1987),

1991)

path

find the optimum

arising

is indepen-

horizon.

we opted

an infeasible

constraints

and solution

egies do not require

taneously.

employing

of NLP

selected a sequential

ution strategy, vs

results

number

optimization

dent of the prediction

x IO-’ x IOf’ x IO-’

equations

The

from sequential

x 1O-3

Key

Biegler

1990,

timization.

on the states and manip-

met.

been

experiments

2-a Interval

1.070 x 10-S 2.734 x 1O+4 4.480 x 10-l

6

ulated

from

Estimated value

91

300

observations for experiment 3 for conditions).

350

rb7241d

400

450

(see Table

2,

can be

G. T. WRIGHTand T. F. EDGAR

92

halted and the solution from the last complete ation can be implemented. 6.1. Sequenti& solution and optimization

iter-

strategy

The nonlinear model predictive controller implemented in this research was formulated as the following NLP: mip @[x(ri), u(zi)]

i = 1,2,.

_ , PH,

subject to satisfying: 1. Model

differential and algebraic equations: E(t)*

= f]x(O,

uU)l>

(8)

u(t) E 4t”, E(t) E W” X“, where x(2) E W”, rank[E(t)]
u(r) =

u1,

r,< t < t,,

u2r

1, < t < t2,

.

1: ucn I

(9)

kn - I G t -= q-h,

3. Initial conditions: x(r) =x0. 4. Simple bounds on differential and algebraic state variables over the prediction horizon: x, c x(ti) < X”,

i=l,...,PH.

5. Simple bounds on manipulated u,
j = 1,. . , CH.

bounds on manipulated

Iu~-u,+~I
variables:

j=l,...,

variables: CH-1.

The resulting NLP was solved using the generalized reduced gradient (GRG2) (Lasdon and Waren, 1986), and u, was implemented. This process was repeated at every sampling instant. The effect of modeling error and unmeasured disturbances was treated as an additive, unmeasured disturbance, and was estimated at the kth sampling instant in a manner similar to DMC: xr,,rcdic,cd.& = Xmx,d.k dk

=

xmeasud.k

Real-time solution of the NLP arising from the controller formulation required special care. As stated earlier, an NLP was solved at each sampling interval. The solution times ranged from little computation time near steady-state to very large computation time for the initial stages of set-point tracking. While it was not possible to determine an absolute upper limit on the computation time required to solve a given NLP, a bound of 5 min (corresponding to the sampling time) permitted completion of the solution process for most of the NLPs which arose. If, based upon the dynamic characteristics of the system, this sampling rate had been too slow, this control strategy would have been abandoned. For the occasional NLP which, if allowed to go to completion, would have required more than one sampling interval to converge to the optimum, the following approach was adopted. The optimization code, GRG2, was modified to call the system clock at the beginning and end of each major optimization iteration. The time required for solution of the NLP. Before each major iteration was initiated, the elapsed time was subtracted from the sampling time to determine the maximum time available for continued computations. The times for previous major iterations were used to extrapolate an expected iteration time for the current iteration. This value was compared with the remaining available time for computation. Based upon the comparison, optimization was either terminated or continued. When terminated, the results from the previous major iteration was taken to be the solution to the NLP, provided the objective function had improved relative to its initial value. 6.1.1. Obtaining gradient information. The computational burden associated with integration of the system equations coupled with a need for gradient information constitutes the greatest impediment to the on-line implementation of sequential optimization and solution strategies. These problems are exacerbated when finite differencing is employed to calculate the gradient information. Let us define v to be the partitioned vector: v = [UT - UT] . . J&IT. Then the gradient respect to v is:

+a,_,,

-

v,o xprcdid,k

This constituted the feedback portion of the algorithm which distinguishes it from the open-loop, optimal manipulated variable profile calculation methods of Biegler (1984) and Renfro er al. (1987). When a perfect process model is used d becomes the additive disturbance in the process output.

of the objective

(10) function

with

= [V,, @I . . . IV,,@],

where V”,@ = y %,x,(0+@,, i--l

j = 1,. . . , CH.

are row vectors and xi is the solution to the DAT at time li. The sensitivities of the differential and algebraic states with respect to each vector ui in the

93

Nonlinear model predictive control sequence

of piecewise constant

x,,. In order to determine function

with respect

to v, the sensitivity

for the states of the DAE Let W, represent vector

inputs are denoted

,__......................,

equations

must be determined.

the sensitivity

with respect

Feed

by

the gradient of the objective

matrix

of the state

to u,: w, = x*,.

Then

the dynamic

evolution

of W, is determined

R

by:

e E(t)W,

i = 1,

= J(x, v)W, + B,(x, v),

a

. . . , CH,

:

(11)

0

where

zj-, < t -c t,,

f” (x. q 1,

Bj(Xv VII =

3

and J(x, v) is the Jacobian equations

are subject

to the initial conditions:

CONTROL

a consequence with

of

a slow

required to accommodate tation

times

Because

evolving

this factor

rejection plemented

in

configuration

As described was

power

varying

the power

to

impact

NMPC

directly

determined

active bed inlet temperature lead

to

active

the bed

boundary This

desired inlet

input to the heater.

The

WGS

was influenced inlet

ranging from 3 min.

generated

Power

value

that would

behavior.

time

interest,

constituted

set-point

the

that the

energy

reactor

for

the inlet

balance.

reactor

reactor

in

9 to

was

be feed-gas

quite

13 SLM.

linear

Time

for

response,

variation

static

of 15%.

gain

was

2.8

with

using

reactor

for

the

was

of

con-

the reactor

inlet

by a model-based

PID

ITAE

controller tuning

5 is typical

constant 6 min

for

and

was

rule

for

of the closed-

had the appearance

time

control

the PID

transfer function

approx.

flowrates

cascade that

track

The an

Figure

which

We have already concluded the relationship

and

flowrate.

defined

flowrates

easily measured

varied

a maximum was

model that

subsequent

The relationship

a function

Since the reactor inlet section

diffused

master-slave,

master).

closed-loop

defined.

10.3 to 14.3 min and the dead time was approx. The

little

effectively

tracking.

behavior

constants

bed

impacting

target value computed (the

perature

foop

inlet temperature

predominantly

behavior

from

inlet temperature

the

therefore,

interest

in the active

of the heat

was used to close the loop.

order plus dead-time The

a portion

of a first-

step response.

the flowrates

the dead-time

of was

3 min.

the NMPC illustrated

loop

(e.g.

CO),

it was imperative slave)

controller

of

of a reactive gas mixture

varied

of

(the

independent

feed

reaction

control

configuration,

tuned,

when

reaction and cool-

the reactor

contained

delay

light

no

heating

Because the static gain, time constant

PID

output.

presumably

Recall

for the reactor

its

beads,

to the inlet, marginally

inlet behavior.

In

However,

from

upstream and

glass

inlet was virtually

was composed

the gas mixture

composition,

The open-loop The

stream

temperature

upon

Pyrex

feed composition.

specifically

Fig. 4. 7.1.

ing of the reactor

troller

a target

for the WGS

While

inert

in this region. Therefore,

all bed temperatures

was used to construct

strategy

was

employed.

as the NMPC

temperature

condition

relationship

control

bed

im-

control

bed inlet temperature.

was never determined Instead,

the primary

the inlet feed

we focused

on the active

cascade rate was

level affected

compositions,

was

linear controller

sampling

previously,

reactor

process.

NMPC

master-slave

faster

long compu-

the disturbance

NMPC,

where a low-level

burden rate was

solution

diminished

of

a

a significantly

and

also

capabilities

the

with

occurred

sampling

the relatively from

*

packed

STUDIES

the computational

NMPC,

----------

Fig. 4. Cascade control configuration for implementing model-based control strategies.

= 0.

7. EXPERIMENTAL

As

l-r

of f(x, v). The sensitivity

W,(r,)

associated

-------a

otherwise,

0

of power

and subject

PID

first-order

plus

described

above.

between bed

dead-time

well as

was poorly

which

Fortunately,

inlet was eliminated

controller

was

bed behavior

to the inlet heater to disturbances

and

inlet bed tem-

behavior

governing

or modeled.

the reactor the

that given flowrate

were not

the need to by assuming

consistently

generated

a

closed-loop

response

as

G. T. WRIGHT and T. F. EDGAR

288

287

286

285

284

1

48

46

44

’

1

0

20

40

60 Time

80

100

(minutes)

Fig. 5. Close-loop response of the reactor inlet temperature for a 5°C set-point increase and a Aowrate of 13 SLM.

7.2.

Closed-loop

NMPC

would

experiments

have made

the problem

feasible for real-time For CH

all

NMPC

was

variable

unity, move

prediction

experiments permitting

over

horizon

120 min. Aggressive

the

the control

only

entire

one time

PH was 24 sampling control

is generally

large CH and small PH. Our objective, not

to demonstrate

aggressive

consistent

control

Therefore,

the tuning

propriately.

over

control,

a broad

Furthermore,

horizon. intervals achieved however, but

operation

parameters a larger

horizon

manipulated

NMPC

for

ution

was

region. ap-

horizon

it would number

Although

or

smooth

were selected control

The

because increased

have

been

on the WGS

every

attempt

was

an

made

to

make

efficient,

the

the sol-

from

3 to 4.5 min.

In

light of this, the sampling

interval

T, for control

of

the sixth

was chosen

bed temperature

This value complied sented

required

by

variables.

computationally

NLP

in-

system

accompanied

of optimization

algorithm of each

computationally

application

by

Wittenmark

Seborg (1990),

with established et

al.

who

to be 5 min.

guidelines

pre-

(1989)

and

Astrom

suggest

that

the sampling

and

95

Nonlinear model predictivecontrol rate be less than a tenth of the dominant time constant or that the ratio of the sampling rate to the time constant lie between 0.1 and OS. The dominant time constant for the system was approx. 55 min and the dead time, 30 min. A nominal dry-gas inlet composition of 40% CO, 40% CO2 and 20% H2 was used for all experiments unless otherwise stated. The volumetric dry-gas to steam ratio was 0.625, and the total gas flow was 13 SLM.

The first experiment was intended simply to demonstrate that NMPC handles set-point tracking smoothly and efficiently. Figure 6 depicts the closedloop response of bed temperature 6 to a 16.5% step-change in its set-point. For clarity, we reiterate that the manipulated variable for the NMPC loop was the inlet temperature set-point, and that the manipulated variable for the PID loop was power to the inlet heater. NMPC was permitted to change the

325

I I ~______________________________._____~

320 G .m

t

set

-Point

-----

Outpuf

-

____________I 4

1

300

240 288 286 284 282

set-POi

nt

----

output

280 278 276

I

3

44.0

2

42.0

1

H

c.

I

0

50

150

100 Time

Fig. 6. NMPC

200

250

300

(minutes)

experimentnumber I: WGS reactor responseto a 16.5”C step increasein the set-point for bed temperaturenumber 6.

G. T. WRIGHT

96

and

inlet temperature set-point by no more than +2.5”C per control interval. Absolute limits of 270 and 300°C were also enforced. Although small oscillations of f 1°C persisted, it is apparent that NMPC achieved the desired set-point. These minor oscillations were the direct result of small inlet temperature oscillations about the inlet temperature target value. Since NMPC is modelbased, dead-time compensation is inherent, provided 360

I

350

the model accounts for it. For the WGS reactor, a temperature variation at the inlet initiates a thermal wave, which amplifies as it propagates through the active bed. This phenomenon effectively creates a lag that would be modeled as pure delay in a transfer function representation of the system. Notice that NMPC required only three sampling intervals to determine the inlet temperature that would drive bed temperature 6 to set-point. Furthermore, once this

1

Set-Point ourput

----’

Set-Point

----

T. F. EDGAR

1

1

-

340

330

320

300

290

output

-

286

284

42.0

36.0 0

100

200 Time

Fig. 7. NMPC temperature

experiment number 2: WGS reactor number 6, which spans the operating

300 (minutes)

400

500

response to a sequence of set-point increases for bed space for nominal feed composition and flowrate.

Nonlinear model predictivecontrol value had been determined, manipulated variable changes virtually ceased, despite the initial absence of response of temperature 6. This example and others that follow powerfully illustrate the inherent deadtime compensation of NMPC. The second experiment was intended to illustrate the ease with which NMPC progresses from a state of virtually no reaction to a state of almost complete reaction when applied to the WGS reactor. This example highlights the effective use of NMPC for plant start-up. Figure 7 shows a sequence of set-point changes. The first was an 18S”C stepincrease from 306.5 to 325°C and the second a 25°C step-increase to 350°C. A velocity constraint permitted NMPC to manipulate the inlet temperature set-point by no more than + l.O”C per control interval. As with the previous example delay time was easily accommodated as evidenced by the absence of excessive manipulated variable movement. More significant, however, was the successful handling of the static gain variations. Figure 7 clearly illustrates that for an 18.5”C change in the output, an inlet temperature change of approx. 5°C was required. For the subsequent 25°C output increase, which occurred at higher CO conversion, an inlet temperature change of approx. 2.4”C was required, a tripling of the static gain. NMPC inherently recognized these gain varitions and responded accordingly. For the same operating conditions, Fig. 8 illustrates the poor simulated response achieved using a PID controller, tuned with ITAE rules for set-point tracking.

97

We note at this point that when Ziegler-Nichols tuning rules are adopted for PID tuning, positive static gain variations should not exceed approx. lOO%, and even this value is borderline. While more advanced PID tuning strategies have been developed more recently, this rule of thumb still loosely applies. The third and final NMPC experiment dealt with the disturbance rejection capabilities of NMPC. First, steady-state was achieved with an output setpoint of 310°C. At time equal to 30 min, the dry-gas flowrate was decreased by 10% to 4.5 SLM. The dry-gas composition and steam flowrate were not altered. Since flowrate and composition measurements appeared as parameters in the NMPC model, an inherent feedforward action caused an immediate drop in the inlet temperature set-point (Fig. 9, arrows mark the flowrate decrease). This occurred before any significant response in bed temperature 6, the feedback variable. In fact, the output only began to respond approx. 10 min later. This disturbance rejection example highlights a flaw of the NMPC control implementation. When constructing the controller, it was assumed that the inlet temperature loop had a perfectly consistent first-order response for set-point tracking. Figure 9 clearly demonstrates that this assumption was violated in the presence of a disturbance. In the next section we discuss the ramifications of this assumption by examining the model states with and without inlet temperature feedback.

370 I 360

-

350

-

z

340

-

2 :: 1

330

-

iz

300 0

I 100

1 200 Time

300 (minutes)

400

500

600

Fig. 8. Simulated WGS reactor response to a sequence of set-point increases for bed temperature number 6, which spans the operating space for nominal feed composition and ilowrate.

G. T. WRIGHTZUI~

98

T. F.EDGAR

320 318

set-Polnc output

----' _

Set-Point output

----. -

316 314

312 310 308 306

286.0

265.0

284.0

282.0

281.0

260.0

40.0 39.5 39.0 36.5 36.0 37.5 37.0 36.5 36.0 35.5 35.0 1

1

1

0

so

150

100

Time

200

lminucesl

Fig. 9. NMFC experiment number 3: WGS reactor response to a 10% step decrease in the nominal dry-gas flowrate.

7.3. Comparison

of

plant

and

model

outputs

for

NMZ’C Experiment three clearly demonstrated a violation of the assumption that the closed-loop behavior of the inlet temperature loop was first-order. The unexpected response of bed temperature 1 was a consequence of the 10% step-decrease in the dry-gas flowrate. Figure 10 compares the output response

(bed temperature 6) actually experienced in the plant to the model response. Notice that the model temperature increased slightly due to the decreased drygas flowrate (arrows mark the flowrate decrease). When compared to the actual temperature increase in the plant, however, the model temperature increase was small. The temperature increase in the plant output was the cumulative effect of increased residence time,

Nonlinear model predictive control

99

I

330

Plant Model Model

__-_. --.--'

vf Inlet Temp. Feedback w/o Inlet Temp. Feedback

325

310

0

50

150

100 Time (minutes)

200

Fig. 10. Comparison

of plant and model states with and without inlet temperature feedback when a 10% step decrease in the nominal dry-gas flowrate is applied to the reactor system.

which

permitted

temperature removal.

and increased

be a slower

Incorportion

have required

inlet

Such

a model

be otherwise

would

accounted where

inlet section

only

size while adding

simulation

increase

for. Figure

10 shows a closed-

the actual

reactor

almost

derivative by the plant.

was

A small,

varying bias ranging from 5 to 10°C persisted, mechanism

of NMPC 380

370

is designed 1

with that

this

phenomenon.

better

We

plant-model

the model

Inlet

feedback

conditions

tively achieve a feedforward

control

would

-

u

350

-

E J E

340

-

P 8

330

-

Figures ture

pre-

was

slowly

no

output

but the to effec-

11 and 12 illustrate

feedback

set-point

strategy (Wright,

would

tracking.

be

In fact,

substantial

that an inlet tempera-

much

distinction

with or without

less

significant

for experiment

feedback.

between

the

temperature 1

operation.

In

experiment

one,

__----_ --._ ,,~~:___---~-~'---r..__ ../-9' ,.;? :,s :*

w/ Inlet Temp. Feedback w/o Inlet Temp. Feedback

_____ ----_' -

Fig.

11. Comparison

i 100

model

As with experiment

l1 0

for

two there

three, the bias varied slowly here, increasing with high

Plant Model Model 300

to

effec-

1992).

-

360

is

is used as an

temperature

boundary

conclude,

agreement

if the actual inlet temperature

to the model.

reset

inlet tem-

the model-output

time

input

that can

In this case,

cisely that experienced

achieved

the overall

information

deal

therefore,

be mod-

perature is used as an input to the model.

feedback

tively

rate of heat

of the second effect would

that the reactor

reactor model loop

reaction,

simply

Only the first of these effects was considered

in the model. eled.

more

caused

200 Time

300 (minutes)

400

of plant and model states with and without inlet temperature sequence of set-point changes was applied to bed temperature 6.

500 feedback

when

a

use

of

G.

325

T.

WRIGHT and T. F. EDGAR

-

-

300

1 0

I 50

Model

w/

Model

w/o

100 Time

Inlet Temp. Feedback Inlet Temp. Feedback

150 (minutes)

1 200

__--. -----'

1 250

300

Fig. 12. Comparison of plant and model states with and without inlet temperaturefeedback when a 16S”C step-pointincreasewas appliedto bed temperature6.

inlet temperature as a model input had a marginally greater effect than for experiment two. Even so, the output behavior in the plant could not be captured for the time interval ranging from 100 to 175 min. This example illustrates, however, that NMPC is robust to plant-model mismatch. It is clear, that in the absence of repeated disturbances, both techniques lead to the same model output, but inlet temperature feedback may significantly affect transient behavior. 7.4. Comparison

with closed-loop

GPC

Adaptive GPC was implemented using the control structure outlined for NMPC. Unlike the NMPC experiments, however bed temperature 4 was taken to be the controlled variable. Because the computation time required for adaptive GPC was small, a sampling time of 2 min, based solely upon the openloop dynamics of bed temperature 4, was adopted. For an inlet temperture of 280°C and nominal values of composition and flowrate, the dead-time was approx. 14min, the time constant 25 mm, and the gain 1.5. The gain increased by approx. 120% from this low reaction state to a state of complete reaction. The control experiment described below used a prediction horizon of 20 sampling intervals, and a control horizon of unity. A move supression factor of 10 was also employed, and GPC was permitted to change the set-point by no more than _t 1°C per control interval. The recursive least squares estimation algorithm of Chen and Norton (1987) was employed for parameter estimation. The model took the form: A (4 -‘ly(t)

= q -‘B(q -‘)24(t) + ci,

(12)

where A (q -‘) and B(q -‘) are polynomials in the backward shift operator of orders 1 and 3, respectively. Figure 13 illustrates a sequence of three 5°C stepincreases in the target value for bed temperature 4. While the close-loop response for the first and second increments were satisfactory, it is clear that the response became progressively worse with increasing operating temperature. The oscillatory behavior was obtained despite parameter adaptation. In addition, this control problem was less challenging than the problem to which NMPC was applied. We concluded, therefore, that traditional adaptive control is not well suited for WGS reactor start-up, since the linear model, even with parameter adaptation, does not adequately reflect the rapidly changing nonlinear dynamics of the system. 8.

CONCLUSlONS

The primary goal of this research was to develop an advanced nonlinear control strategy for fixed-bed catalytic reactors. The control method was applied experimentally using a laboratory-scale water-gas shift WGS reactor. The following conclusions may be drawn from the results of this work. An adiabatic, pseudo-homogeneous WGS reactor model represented the physical system well over the operating space of interest. The physically reasonable, simplifying assusmptions that were adopted proved useful in developing a low-order model, suitable for implementation in an NMPC framework. The Galkerin technique on finite elements with piecewise linear polynomial approximations proved not to

Nonlinear model predictive control

306

302

296

Set-Point

output

----’

-

286

284

46.0

3

50

100

150

200 Time

250 (minufcr)

300

350

400

450

Fig. 13. Adaptive GF’C experiment: WGS reactor response to a sequence of set-point increases for bed temperature number 4 for nominal composition and flowrate.

be susceptible

to oscillatory

behavior

over the spatial

domain when 12 nodes were employed for discretization. Estimation of “dynamic” and steady-state” parameters were efkctively decoupled for the pseudohomogeneous WGS reactor model. Furthermore, information-rich dynamic data yielded good parameter estimates with less experimental effort. The control experiments demonstrated that absolute plant-model agreement was not imperative for

good control using NMPC. However, temporal firstderivative information, consistent with plant behavior, was crucial to good performance. NMPC was better suited for feedforward dynamic compensation than linear techniques since the nonlinear model has an inherent characterization of the feedforward mechanism. Feedforward control significantly enhances NMPC performance. Finally, NMPC was effectively used to start-up the WGS

G. T. WRIGHT and T. F. EDGAR

102

NMPC

system. ditional operating control

was superior

control

techniques

regions appears

were to

estimates

varied

successful

parameter

as

in this regard since

broad

traversed.

be unsuitable, rapidly

as

adjustment

Adaptive since

the

linear

parameter

state,

extremely

to tra-

nonlinear

making difficult.

REFERENCES

Ampaya J. P. and R. G. Rinker, Autothermal reactors with internal heat exchange. J. Chem. Engng Sci. 32, 1327 (1977).

Asselmeyer B., Optimal control for nonlinear systems calculated with small computers. J. Opt. Theory Applic. 45, 533 (1985).

Astrom

K. J. and B. Wittenmark, Computer Controlled 2nd Edn. Prentice-Hall. New

Svstems: Theorv and Des&n. Jkrsey (1990). _ I

Bates D. M. and D. G. Watts, Nonlinear Regression Analysis and Its A&ications. Wiley, New York (1988). Bell N. H.-&d T. F. Edgar, Modeling. of a fixed-bed water-gas shift reactor 1: Steady-state model verification. J. Process Confrol 1, 22-31 (1991). Biegler L., Solution of dynamic optimization problems by successive quadratic programming and orthogonal collocation. Computers them. Engng 8, 243 (1984). Bonvin D., Dynamic modeling and control structures for a tubular autothermal reactor at an unstable state. Ph.D. Thesis, University of California at Santa Barbara (1980). Brockett R. W., Feedback invariants for nonlinear systems. In Proc. 7th ZFAC World Congr., Helsinki (1978). Caracotsios M., Model parametric sensitivity analysis and nonlinear parameter estimation-theory and applications. Ph.D. Thesis, The University of Wisconsin at Madison (1986). Carberry J. C. and M. M. Wendel, A computer model of the fixed bed catalytic reactor: the adiabatic and quasi-adiabatic cases. AZChE Jf 9, 129-133 (1963). Chen M. J. and J. P. Norton, Estimation technique for tracking rapid parameter changes. ht. J. Comrol 45, 1387-1398 (1987). Clarke, D. W., C. Mohtadi and P. S. Tuffs, Generalized predictive control-part i the basic algorithm. Automatic 23,

137-148

(1987).

Cuthrell J. E. and L. T. Biegler, On the optimization of different/algebraic process systems. AZChE JI 33, 1257 (1987).

Cutler C. R. and B. L. Ramaker, Dynamic matrix controla computer control algorithm. In Joinr Automatic Control Conf. Proc. Paper WP-5B San Francisco, CA (1980). Danckwerts P. V., Continuous flow systems-distribution of residence times. Chem. Engng Sci. 2, 1-13 (1953). Eaton J. W. and J. B. Rawlings, Feedback control of chemical processes using on-line optimization techniques. Computers them. Engng 14, 169-479 (1990). Eaton J. W. and J. B. Rawlings, Model predictive control of chemical processes. Chem. Engng Sci. (1991). Hertzberg T. and 0. A. Asbjomsen, Parameter estimation in nonlinear differential equations. In Computer Appiications in the Analysis of Data and Plants. Science Press, Princeton (1977). Hunt L. R., R. Su and G. Meyer, Global transformations of nonlinear systems. IEEE Trans. Automat. Contr. AC28,

24 (1983).

Jang S., B. Joseph and H. Mukai, On-line optimization of constrained multivariable processes. AZChE Jf 33, 26 (1987). Kiparissides C. and A. Georgiou, Finite-element solution of nonlinear optimal control problems with a quadratic performance index. Computers them. Engng 11,77 ( 1987). Lasdon L. S. and A. D. Waren, GRG User’s Guide. Technical Report, University of Texas at Austin (1986). Lee A. L., Effect of raw product gas composition on shift catalysts. Technical Report Project 61015, DOE Contact DE-AC2I-78ET11330, Institute of Gas Technology. Chicago (1980). Moe J. M., Design of water-gas shift reactors. Chem. Engng. Proc. 58, 33 (1962). Mor& J. J., B. S. Garbow and K. W. Hillstrom, User Guide for Minpack- 1. Technical Report ANL-80-74, Argonne National Laboratory (1980). Morshedi A. M., Universal dynamic matrix control. In Chemical

Control

Conf.

IZI, Session

VI,

Paper

No

2,

Asilomar, CA (1986). Newsome D. S., The water-gas shift reaction. Cufal. Rev.Sci. Engng 21, 275 (1980). Patwardhan A. A. and T. F. Edgar, Nonlinear modelpredictive control of a packed distillation column. In Am. Control Conf., Boston, MA (1991). Patwardhan A. A., J. B. Rawlings and T. F. Edgar, Model predictive control of nonlinear processes in the presence of constraints. In AZChE JI Annual Meeting, Washington, DC (1988). Patwardhan, A. A., J. B. Rawlings and T. F. Edgar, Model predictive control of nonlinear processes in the presence of constraints. In ZFACJZEEE Symp. on Nonlinear Control Systems Design, pp. 45&460 (1989). Patwardhan A. A., J. B. Rawlings and T. F. Edgar, Nonlinear model predictive control. Chem. Engng Commun. 87,

123 (1990).

Peterson T. E., E. Hernandez, Y. Arkun and F. J. Schork, Nonlinear predictive control of a semi-batch polymerization reactor by an extended dmc, Paper No. tp4. In Am. Confrol. Conf., Pittsburg, PA (1989). Petzold L. R., A description of DASSL: a differentialalgebraic system solver. In Scientific Computing (R. S. Stepleman, Ed.), pp. 65-68. North-Holland, Amsterdam. Rase H. F., Chemical Reactor Design for Process Plants, Vol. 1. Wiley, New York (1977). Rawlings J. B. and K. R. Muske, The stability of constrained receding horizon control. Submitted to IEEE Trans. Auto. Contr., Sept. (1991). Renfro J. G., A. M. Morshedi and 0. A. Asbjornsen. Simultaneous optimization and solution of systems described by differential/algebraic equations. Computers them. Engng 11, SO3 (1987). Seborg D. E., T. F. Edgar and D. A. Mellichamp, Process Dynamics and Control. Wiley, New York (1989). United Catalysts

C- I2 Operating

Instructions.

Windes, L. C., Modeling and control of a packed bed reactor. Ph.D. Thesis, The University of Wisconsin at Madison (1986). Windes, L. C., M. J. Schwedock and W. H. Ray. Steady state and dynamic modeling of a packed bed reactor for the partial oxidation of methanol to formaldehyde-I, model development. Chem. Engng Commun. 78, 143. (1989).

Wright G. T., Modeling and nonlinear predictive control of a fixed-bed water-gas shift reactor. Ph.D. Thesis, The University of Texas at Austin (1992).