Inverse Simulation for Nonlinear Systems Analysis

Co pyright

©

IFAC Simulation of

Co ntro l Syste m s . Vie nna. A.u stria . 19Rf)

INVERSE SIMULATION FOR NONLINEAR SYSTEM ANALYSIS D. Kraft I lISlill/ l f

fiJl" F lighl

.'l\"sII'IIIS D\" I/ ({lIIirs, GI'I"IIWII Aems/}(IU' R l'sl'a r r h E sI({/Jlisillll I'll I .

0 111'1#;1/"111101"11, D·8031

W esslillg.

FRG

Abstract. A heuristic analysis tool for nonlinear dynamical systems is described, which is based on the solution of a sequence of optimal control problems. These are solved by a direct shooting method, which uses information of system simulation to update the control parameters into a desi red di rection. It is possible to consider vector optimal control problems. Key wo rds • Opt ima 1 cont ro 1, systems, system analysis.

non 1 i nea r

programmi ng,

vector

opt imi za t i on,

non 1i nea r

I NTRODUCTI ON The analysis and design procedure outlined above will be exemp 1 if i ed on the problem of cont ro 11 i ng the cryogenic windtunnel of the DFVLR at Cologne (Kraft, 19 86).

A necessary prerequisite for the design of control systems for nonlinear dynamical processes is a careful analysis of the time-response of the process due to the control inputs. Choice and dimensioning of the controls as well as the corresponding investment costs are based on these analysis results. Such investigations also give informations on the reachable operation region of the process and on deviations of certain output values from their required values resulting from systems disturbances.

ANALYSIS PROBLEM In analyzi ng dynamic systems we have the foIl owi ng primary problem : determine the set of controls within a given admissible control domain which carry the system states from certain initial values at the beginning of the process to prescribed final values at the end of the process while keeping them within a feasible region.

A widely used desi gn methodology consists in system linearization and eigenvalue/eigenvector assignment by pole-shifting or optimization techniques. The results obtained for the linear model are then tested on the nonlinear model by system simulation. To this end special simulation software too 1s have been deve loped in t he past (ACSL, CSMP, ECSSL, etc . ). Common to all of them is that the user has to specifiy certain input functions and initial values of the system state variables; and within the simulation the initial value problem with differential equations is solved yielding the time history of the state variables.

The secondary problem consists in minimizing a cost functional which depends on the states and/or the controls within the entire time interval of the process. While optimizing the process this secondary task guarantees local uniqueness of the primary problem. The third problem is to determine the sensitivity of the cost with respect to changes in the admissible control domain and to find a best compromi se between the cost associated with the process and that in realizing the control.

This procedure is somewhat unsystematic as the analyzi ng engi neer is not interested primarily in a state trajectory resulting from some predetermined input. He woul d rather 1i ke to specify the output and determi ne the (eventually constrai ned) input tha t is necessary to generate the des i red output. This means that a solution for the inverse simulation problem is searched. An overview of this idea is given in Fig. 1.

Examples of criteria in optimum engineering design are for control and state region specification and payoff, respectively: values and rates of mass flows in valves, of rotational speed of fans and servomotors, of deflection of aerodynamic devices, accuracy, safety and velocity of the process,

The adequate mathematical tool for this is the solution of constrai~ed boundary value problems, together with the minimization of some cost functional to give a (locally) unique solution. This optimal control problem will be solved sequentially for systematically varied constraints to gain insight into the system sensitivity, as is shown in Fi g. 2.

consunption of material, cost and time of the process. The payoff may consist of a vector of cost functionals the elements of which may behave in a conflicting manner. The statements above can be forma 1 i zed as follows : determine piecewise continuous control functions m u(t) £ R , t £ [ to,tf l and parameters p E RP to

Generally the boundary values for the problem are given by the operation phases of the process while the constraints are defined by the operation region. It should be clear that state constraints as well as control constraints are permitted.

solve the constrained boundary value problem

155

156

D. Kraft

x-

o £ Rn, o £ Rq ,

f(x ,u,p)

b(x(t ) ,x(t f )) o c(x,u,p) cl .; u(t)

£

.; C

u

£

t

£

[to,tfl,

r R , t

£

[to,tfl,

U C Rm.

bi (y)

M(x(t f ) ,p).

Then by systematically vari ing cl (t) and cu(t) determine the sensitivity of the cost vector Wlth respect to these variations. Note that for every variation a complete optimal control problem is to be solved. It should be clear that an efficient and robust numerical solver of this problem has to be the core of this analysis methodology. ITERATIVE NUMERICAL SOLUTION A direct shooting method is applied to the approximate solution of the above prob1 em (Kraft, 1985). For thi s reason the i nterva 1 of the independent variable is discretized by the knot sequences i i i "' i := {to
t j , the knots themse1 ves, and the parameters p of the above problem: 1 m {Yi } = {u1 (to),u1(t1),· · ·,um(tf_1),um(tf)'

(u~,t~) approximate the control func-

tions in the entire interval [ to,tf l by, for example, interpolat i ng cubic sp1ine functions 1

S"' i := [ to,tf l .. R (i) (ii)

£

0,

2 C [ to,tf l , 0, •.• ,f-1,

S", i

l, ... ,m,

(i i i) S", ij

O, ••• ,f-l.

By this procedure the pair (u,p) in the optimal control prob1 em is replaced by a finite control parameter set, as is shown in Fig. 3 for m=l. Next the initial value problem is solved with appropriately chosen initial conditions (where Xl (to) augments the set {Yi } if the initial state is not prescribed for some index 1). Thus x(t) is known once y is given, and the remaining optimal control problem reduces to a non1inear vector programmi ng problem £

I . = {i : i=1, ••• ,l}

£

Y : = { y:

i =1, ••• ,q, i =1, ••• , r,

j=1, ••• ,s},

where the constraints are imposed on a finite set of communication knots "' s := {t o < t1 < ••• < t s _1 < t s } only, as is shown in Fig. 4 for infeasible xk and feasible x*. This problem is solved by a constrained variant of the method of Krei Bel mei er and Stei nhauser (1979), which finds the efficient or Pareto-optima1 points (Jahn, 1985) with the help of the fo11 owi ng interactive a1 gorithm : a)

choose a suitable starting approximation o yO, together with an upper bound u on J such that

with the notations a
id c)

The pairs

Y

d1 .; d(Y) .;

and to simultaneously minimize the cost vector J (u ,p)

subject to

id

k k k k-1 • choose u such that J(y )< u < u

Remarks: Step b) is solved by the minimax algorithm of Murray and Overton (19 80). In step c) the deci si on-maker interact i ve 1y reduces the sequence k of upper bounds {u } in such a way that the best possible compromise betv.een the conflicting cost elements is eventually found. The main analysis step now consists in the systematic variation of the upper and lower bounds specifying the set Y and in the determination of the sensitivity matrix aJi/ ac j . From this information the decision-maker chooses specifications for the performance of control actuators (Kraft, 1986). ANALYSIS EXAMPLE As an example for the application of the analysis procedure we choose the optimal control of a cryogenic wind tunnel (Kraft , 1985a,1986). Figure 5 shows a cross-section of the cryogenic wi nd tunnel of the DFVLR at Cologne (KKK = KryoKana1 Ko1 n) together with the four control elements : fan, liquid and gaseous nitrogen injection and gaseous nitrogen ejection. These control the following state variables of the tunnel gas: velocity w, temperature T, and pressure p, where the values are those in the tunnel test section. The mathematical model describing the connection of the control inputs with the state variables which are to be controlled within their constraints is established through thermo-f1uidmechanica1 balance principles. The balance of mass flows yields the continuity equat i on

Inve rse Simulation fo r Nonlinea r System Analysis

m where m is the total ma ss flow and i ndi ces 1 to 3 indicate the mass flows of liquid nitrogen injection, gaseous nitrogen injection, and gaseous nitrogen ejection, respectively. The caloric part of the energy equation gives an expression for the change of temperature

mT

where

m2 (t-T2)[T-rT U]

+

+

m3(t- T 3)[ r - 1] T 3 kpw A/2

-

a wAw(T- Tw)

-

aMAM(T-TM)'

mi(t-Ti)

indicates

We will consider minimLl11 time problems other cost effectiveness factors like mass and consumed el ectric energy are importance w.r.t. optimization because nearly compl etely prescribed by the conditions (boundary values). This is Kraft (1986). the

influence

of

the

transportation time lag in the controls with respect to their action on the gas in the test section. A is the test cross section, p the density of nitrogen, r the evaporation enthalpy, c the speclflc heat at constant volume, n the degr~e of efficiency, and a the heat transfer coefficient. Indices Sand U indicate saturation and ambient condi t ions, respect i ve ly, and i ndi ces Wand M represent the energy flow from the wall and the metal parts within the tunnel. The kinetic part of the energy equation relates the change of velocity of the tunnel gas to the external work introduced by the fan and the head loss of the gas 2 2 plw g(T)n - kpw /2 with fan rotational characteri st i c.

velocity nand

g(T)

as

fan

The gas pressure can be derived from the introduc ed quantities and with R as specific gas constant via the state equation for a perfect gas

p

pRT,

=

or in di fferent i a 1 form Pip

=

m/m +

T( l/T+PT/ P )

+ 'ilpw/ p + uPa / p ,

a

is the time derivawhere p~ denotes ap / ax , and tive or the angle-of-atta ck of the model to be tested in the tunnel. Kraft (19 86) gives further information on the development of the above model equations. The nonlinear differential eq uations are strongly coupled, and the control and state variables entering .them have to satisfy accuracy and safety constralnts (termed as operating conditions) in the entire time interval as can be seen from Tables 1 and 2. It is therefore evident that careful control analysis must precede operation of the tunnel and tha~ the tunnel gas has to be controlled by an automatlc control system in the entire operating regl0n

o~

w

~

The control design invol ves a trade-off between t\\Q conflicting demands, namely the above mentio ned satisfaction of the operating conditions and the requi rement for high cost effectiveness and large tunnel productivity. For tighter specifications the ope r ation time required to achieve them is longer and therefore the tunnel productivity in t~rms of possible test runs will be lower, and Vl ce versa. This control trade-off will be illustrated for cryogenic temperatures (T = 100 K) and for the operating phase "acceleration" of the tunnel \\Qrking gas from idle to ma x imlJ1l test velocity at thi s temperature (0 m/ s ~ w < 70 m/ s).

m1 (t-T 1)[ T-rTs+r/cv ] +

157

100 m/ s

100 < T < 300 K and in a11 operating phases (cool down, stand-by, accelerate, test, decelerate, warm-up).

on ly, as injected of mi nor they are ope rat i ng shown in

Th~

.a nalysis is organized in the following steps. Inltlally the state tolerances are varied at nominal control constraints and the influence of this variation on the operation time is checked. In the second step the control tolerances wi11 be varied at nominal state constraints. Again operation time is the criterion of interest. ~nally the influen ce of valve location is investigated for nominal operating conditions. The results of the first step are slJ1lmarized in Fig. 6, where the operation time t f is given for increaSing tolerances ~ T and 6 p of temperature and pressure, respectively. In the lower graph (x) H and ~ p are linearly related with the indicated values. In both the other graphs the state variables H (ll ) and li p ( D) are not a11 owed to deviate from thei r nomi nal values (T = 100 K, P = 1000 hPa). In these cases as well as in the limiting case li T = 0 and 6p = 0 the corresponding differential equations change to nonl inear algebraic equations from which the corresponding control functions can be determined. Thus a considerable reduction of computer time can be achieved. A quantitative derivation of these results is given in Kraft (198 6). If the constraints on the state variables are tightened by decreas i ng the to I erances more dynami c freedom is taken from the process and thus in all three cases the operation time increases up to its limit i ng value where 6 T and 6 p vanish simultaneously. In step two of the control analysis we concentrate on t he v~ ri at i on of ~ he contro I to I erances of 1fl2 because lt has a strl klng influence on the operatlon duratl0n. The remaining control tolerances are k ~pt at their nominal values. As in step one, restrlctlng the control takes freedom from the process and prolongs its duration. This is indicat~d in Fi g: 7 for stat i ona ry (~b = 0) and nonstatl0nary (xh 0) boundary conaitions. For the former condltions the time gained when injecting gaseous nitrogen with ambient temperature at nominal values (m 2 = 4 kg/s) as opposed to noni nject-

*

ing is 250 percent. Fina11y, for this analysis step, the control functions 1fl2 are constrained to the values shown in Fig. 8. The bang-bang structure of the graphs can be expected from the bil inearity of the system in these components At the left and right boundaries the control r~tes are also constrained.

158

D. Kraft

The last step in the control analysis demonstrates an alternative possibl il ity of this analysis approach. Instead of vari ing the boundaries of the controls and states certain design parameters that defi ne the process are changed and the i nfl uence of this on the cost is observed and evaluated for the design specifications. This possibility is exemplified by the variation of the valve location for liquid and gaseous nitrogen. The in fl uence of the di stance 1 of the liquid nitrogen valve from the test secti6n on the process duration is shown in Fig. 9 for nominal and optimal location 12 of the gaseous nitrogen valve. For extended physical interpretations of all observed phenomena to Kraft (1986) is referred.

g;'1 2:1 input

find oulpul SIMUL ATION

Slocc ced

Protlem

U

I~ r-

MC~lL L

i " I x. ul

I

OPTIMIZATION In verse Problem

ut ~

~ I

gi ven output

find Input

CONCL US IONS The analysis results of the last section often show surprising results in the behaviour of nonlinear systems. These effects can rarely be obtained by trial and error simulation. Instead of this a systematic analysis tool which combines simulation and optimization and which has been described in this contribution is needed in practical system analysis and design.

Fig. 1.

The inverse simulation problem.

REFERENCES Jahn,

J. (1985). Some characterizations of the optimal solution of a vector optimization problem. OR Spektrum, 7, 7-17.

Kraft, D. (1986). Optimalsteuerungen ein systematisches Hilfsmittel zur rechnergestUtzten Erforschung der dynamischen Moglichkeiten ei nes Ti eftemperaturwi ndkana Is. To appear as DFVLR FB 19B6. Kraft, D. (1985). On converting optimal control problems into nonlinear programmlng problems. In K. Schittkowski (Ed.), ComputatioSpringer, nal Mathematical Programming. Berl i n. Kraft, D. (1985a). Optimal control of a high performance wind tunnel. In G. di Pillo (Ed.), Control Applications of Nonlinear ProgrammOptimization. Pergamon Press, ing and Oxford.

Fig. 3.

Cross section of the cryogenic wind tunne 1.

Krei(3elmeier, G. and Steinhauser, R. (1979). Systematische Auslegung von Reglern durch Optimierung eines vektoriellen GUtekriteriums. Regelungstechnik, 27, 76-79. Murray, 101., and Overton, M.L. (1980). A projected Lagrangian algorithn for nonlinear minimax optimization. SIAM J. Sei. Stat. Comp., 1, 345-370.

U;

t;

Fig. 4.

Control parameterization.

T;

tf

Fig. 5.

Feasible and infeasible solution.

159

Inverse Simulation for Nonlinear System Analysis

/

1/ / 6u,

6x j

-+-

of a;= Ir= SEN SiiIViiY . ~' arc;

- "-

"'\

I,

\, / //

10

'-

Fi g. 2.

6x,6 u

System sensitivity analysis.

LilT i t ' ime

~

1\ ~ "'~ \ """" '\

In

~

~

--s-.

~

1\

~ t---.

~ h-

t.p ~nd

t.T

r---I><1.00

f

xb

a

1.50

2 . 00 (hPa) ,

t.T

(K)

E0 -7

Influence of state constraints on operation time.

!\

1\

~ \ '" 2

~~~

r---- ~ ~ ~

.0

.1

,, l?

"

= c pti

= r ami al

~ ""

.~ r--

o~

l---

~ I--"

~

.0.1.2.3.4.5.5.7.8.91.flE3 II (m) -7

Fig. 8.

Influence of valve location on ope rat ion time. TABLE 2 Actuator

(in

fr e

.2

.3

.4

- ea~

.5

l----E>

~

ibl )

.5

.7

.8

Fig. 7.

mi n

State Constraints lI w( %)

lI T(K)

li p (%)

Cooling

L

F

±2.0

Stand-by

L

±5.0

±2 .0

Accelerate

F

±2 .U

±2.U

± 1. 5 (± U. 5)

(:t 0.05)

Testing

±1.0 (± 0.5)

±1.U

Decelerate

F

± 2.0

±2 .U

Via rm-up

L

F

±2.U

(F

=

Optimized reference,

L

Control Rate Constraints

ma x

min

max

-

5 Rpm/s

+ 5 Rpm/s

12 kg / s

-

U.375 kg / s 2

+ U.375 kg / s

0 kg / s

4 kg / s

U.4UU kg / s

0 kg/s

12 kg/s

-

Fa n

0 Rpm

Nit rogen Injection

0 kg / s

Bl ow-i n Ejection

500 Rpm

1.0El

Influence of constraints on gaseous nitrogen injection on operation time.

Control and Control Rate Constraints Control Constraints

.9

';'2 (kg/s) -7

Operation Phase

\

a

r-- I><-

ivar

TABLE

If)

=

it =

a . 1 t.p

V1

xb

\ t\

a N

. 5fl

Fi g. 6.

\\

f-s- ~ i-A.

~ ~ ~'I'=

N

.flfl

~\

~ /tip= a

2

2.4UO kg/s 2

~

+ U.4UU Kg / s

2

+ 2.4UU kg/s

2

=

Id 1 e)