NONITERATIVE SOLUTION _
OF THE MINIMUM REGULATOR

AND MINIMUM IIMEWPROBLEMS. , "

For The Degree of Ph D
chhlgan State Un1vers1ty

GERALD E. BERNIER
1972  .‘

__M-..-4_ . m.”—

 

This is to certify that the

thesis entitled

NONITERATIVE SOLUTION
OF THE MINIMUM REGULATOR
AND MINIMUM TIME PROBLEMS

presented by
GERALD E. BERNIER

has been accepted towards fulfillment

of the requirements for
DOCTOR OF

PHILOSOPHY degree in ELECTRICAL ENGINEERING

 

 

Major professor

Date I“ /§'/~7/’

0-7 639

ABSTRACT

NONITERATIVE SOLUTION
OF THE MINIMMM REGULATOR
AND MINIMUM TIME PROBLEMS

BY

Gerald E. Bernier

This thesis considers the direct, noniterative computation of the
solutions of a set of optimal control problems. The specific control
problem considered is that of finding how close to the origin a given
initial state can be driven in a fixed run time. This problem is
studied for linear systems and for nonlinear systems which are linear
in the control vector. The set of optimal final states for all pos-
sible run times must also be continuous in the state space.

This locus of optimal final states is shown to be the state tra-
jectory of another system which is called the trajectory system. The
differential equation of the trajectory system is developed in the
dissertation. Computational solution requires that this differential
equation be solved simultaneously with a constraint condition, which
is also developed.

An algorithm is presented along with two different methods of
solution. One method is shown to be generally superior. It is used
in a special digital program to solve some example problems. The re-
sults of an example problem solved with an analog computer are also

included.

 

NONITERATIVE SOLUTION
OF THE MINIMUM REGULATOR
AND MINIMUM TIME PROBLEMS

By
I

l a
Gerald E: Bernier

A THESIS

Submitted to
Michigan State University
in partial fulfillment of the requirements
for the degree of

DOCTOR OF PHILOSOPHY

Department of Electrical Engineering and Systems Science

1972

List of

List of

1.

Tables . . . . . . . .

Figures . . . . . . .

Introduction . . . . . . . .

2. Preliminary Considerations .

2.1

2.2

2.3

2.4

2.5

2.6

2.7

2.8

Notation and Terminology

SYstem Definition . . .

TABLE OF CONTENTS

The Reachable Set and Controllability . . . . . . .

The Control Problem .

The Maximum Principle . .

Normality and the Optimal Control Vector . . . . . .

*
Dependence of u (T,t) on T . . . . . . . . .

Preliminary Mathematics .

O

3. The Time-invariant Linear System . . .

3.1
3.2

3.3

3.4
3.5
The

4.1

4.3

4.4

The Trajectory System .

Evaluation of ti(T)

Computational Solution of the Trajectory System

Differential Equation . .

Computational Results
Summary . . . . . .
General Linear System

The Trajectory System .

O O O O 0

Solution of the Trajectory System Equation . . . .

Computational Methods .

An Example Problem . . .

ii

 

10

12

13

23

26

26

33

34

45

48

52

52

55

57

59

4.5 Summary

5. A Class of Nonlinear Systems

0

 

5.1 The Trajectory System Differential Equation

5.2 Nonconvex Reachable Sets. . . . . . . .

5.3 Solution of the Trajectory System Equation

5.4 An Example Problem .

5.5 Summary

6. Conclusions and Extensions

6.1 Summary and Conclusions

6.2 Generalization and Future Research Possibilities -

Bibliography
Appendix I

Appendix II

111

60

61

62

64

66

68

68

70

7O

72

75

76

8O

3.4.1.

IIOIO

LIST OF TABLES

Error in the Solutions to Problem 2 for the
Double-Integrator Plant . . . . .

Fortran IV Program

iv

46

82

2.7.1.

2.7.2.

2.7.3.

2.7.4.

2.7.5.

2.7.6.

2.7.7.

3.4.1.

3.4.2.

3.4.3.

3.4.4.

4.4.1.

5.2.1.

5.4.1.

1.1.

1.2.

1.3.

LIST OF FIGURES

The Double Integrator Plant with x0 =l:i] .

~1

The Double Integrator Plant with x0 = [ 1] .

Switching Types According to the Location

of 50 in the Double Integrator Plant

The Set of Costate Trajectories . . . . . .
The Alternate Situation .

r Versus T For a Type I Switch. . . .

T Versus T For a Type II Switch .

Solutions of Problem 1 for the

Double-Integrator . . . . .

1.0
0.0

Plant with x_=

Analog Computer Soljtion of Double—Integrator

Solutions to Problems 1 and 2 for the
System of Equation (3.4.1.)

Solutions to Problems 1 and 2 for the
Third Order System of Equation (3.4.3.) .

Solutions of Problems 1 and 2 for the
System of Equation (4.4.1.)

A Nonconvex Reachable Set . . . . . .

Solutions of Problems 1 and 2 for the
System of Equation (5.4.1.) . . . . . . . .

Trajectories for the Double—Integrator Plant

Solutions of the Regulator Problem with a
TYpe I Switch C C O O O O O O O O O O O O 0

Solutions of the Regulator Problem with a
Type II Switch. . . . . . . . . . . . . .

l7

l8

l9

l9

19

46

49

50

60

65

69

77

79

79

1. INTRODUCTION

One of the branches of control theory which has received a great
deal of attention in the past twenty years is optimal control. In an
optimal control problem, one starts with a system and a cost functional,
J. The object of the optimal control problem is to apply inputs to
the system in such a way as to minimize the cost functional J. In this
dissertation, the optimal control problem of primary interest is to
apply bounded inputs to the system so as to drive its state from a
given initial state to the closest possible point to the origin within
a fixed run time. The solution is considered to be the final state.

An associated problem which is also considered is to apply bounded
inputs to the system in such a manner as to drive the state from a
given initial state to the origin in a minimum amount of time. The
solution is considered to be that minimum run time. These are referred
to as Problem 1 and Problem 2, reSpectively.

Probably the most significant event in the area of optimal control
theory was the publication of the Maximum Principle by L. S. Pontryagin
and his co—workers [12]. In using the Maximum Principle, a function
called the Hamiltonian which depends upon the cost functional J, the
derivative of the state vector, and a costate vector is defined. If
a control or input vector is optimal, i.e., minimizes the cost func—
tional, then it minimizes the Hamiltonian over the range of permissible
control vectors, and the costate vector at the end time is transversal
to the target set. The Maximum Principle describes a set of necessary
conditions for optimal controls, and usually these necessary conditions
do not provide sufficient information to calculate the Optimal controls

analytically. In general, there is no method for obtaining a closed

form solution to/Problem l or Problem 2 —-or any other realistic
optimal control problem-—-so the computation of optimal controls is
a very important branch of optimal control theory.

The methods for solving optimal control problems can be separated
into three broad classes. The first refers to those special cases where
it is possible to solve an optimal control problem analytically. Some
examples of this are given in Athans and Falb [l], but these methods
are not generally very useful because of the limited class of problems
to which they apply.

Secondly, consider a system which has discrete states and discrete
input levels which are defined at discrete points of time. In such a
system, all possible inputs could be exhaustively considered and the
results compared to find an optimum. This approach fails in most real-
istic problems because it requires too much time. The Dynamic Program-
ming method [4] utilizes the Principle of Optimality to eliminate most
of the permissible controls without first trying them, thus greatly
reducing the time necessary for the search procedure. Even though it
reduces the computation time, Dynamic Programming requires so much
storage space that the procedure is not often feasible.

The last group of computational methods —-which are of primary
significance —-are generally based on the Maximum Principle or the
Calculus of Variations. These methods can be further subdivided. Con-
sidering the solution of Problem 1, for an example, the primitive method
is to guess at a control function, run the system using that control
function, and in some way improve the control function — i.e., find
one that causes the final state to be closer to the origin. This pro-

cess is referred to as iterating over the control function. Such a

method was introduced by Ho [9]. It should be obvious that iteration
over an entire function would not normally be a very good method be-
cause of the large storage capacity that would be required. Most other
computational methods have improved this primitive approach by finding
single vectors or scalars to iterate over. For instance, one might
iterate over the initial or final value of the optimal costate vector
to solve Problem 1. This would be sufficient because knowledge of the
initial value of the costate would permit running the whole problem in
order to find the Optimum final state or the control function which
achieves it. Ho's method was improved by Fancher [6] who iterated
over the final value of the state vector rather than the whole control
function. A representative group of the other methods of this general
class are those due to (alphabetically) Barr [2], Bryson [5], Gilbert
[7], Neustadt [11] and Stratton [13].

The various methods in this third group are all iterative methods.
This means that an initial guess of some quantity representative of the
problem solution is made; then that and successive values of the given
quantity are evaluated and improved until an optimum solution is obtained.
Computational methods are called direct or indirect as well as itera—
tive or noniterative. A direct method is one which is primarily con-
cerned with the variable or variables which are considered to be the
problem solution. This may not always be a meaningful designation.
As an example, Ho's method is primarily concerned with the control
function so it is an indirect method for solving Problem 1. But if
the problem solution were defined to be the optimal control function,
Ho's method would be termed direct. Normally, however, practical con-

siderations force a Specific quantity to be designated the "solution”

of the problem. Fancher's method is a direct, iterative technique
for solving Problem 1.

This thesis is the reSult of research that was motivated by a
dual consideration:

1) Iterative methods are often undesirable because of convergence
problems and possible long computation times; and

2) A method that could be used on an analog computer alone could
prove to be advantageous.

These are termed a dual consideration because iterative methods cannot
be used with an analog computer alone —-there must be a digital or a
human interface. A computational procedure must also be direct if the
problem solution is to be presented on an analog computer. So the
research was aimed at developing a direct, noniterative computation
procedure for solving Problem 1. Such a procedure introduces a dif-
ferent philosophy to the solution of optimal control problems; it is
hoped that this philosophy will eventually lead to better methods of
solving a wide range of optimal control problems. The solution method
presented in this thesis has practical significance. Given an arbi—
trary initial state, there are an infinite number of meaningful cases
of Problem 1 and there is the additional case of Problem 2. The pro-
cedure presented in the sequel makes it possible to solve all of these
problems with a computation time which is of the same order of magnitude
as that which is required by other methods to solve just one case of
Problem 1.

After definitions, terminology and problem statements, some Lemmas
are presented in Chapter 2. The computation procedure is developed
for time invariant linear systems, time varying linear systems, and
nonlinear systems in Chapter 3, Chapter 4, and Chapter 5 respectively.

Then the results are summarized in Chapter 6.

2. PRELIMINARY CONSIDERATIONS

The purpose of this chapter is to define the problem which is
solved in Succeeding chapters and to present the background material
which is pertinent to solving that problem. In preparation for this
objective, section 2.1 is a presentation of the notation and termin-
ology to be used throughout. The systems to be studied are defined
in sections 2.2 and 2.3, and the control problem in section 2.4. The
results of the Maximum Principle, as it pertains to the control problem
at hand, are reviewed in section 2.5, and the normality condition is
attached to the control problem in section 2.6. The optimal control
vector and its switch times are discussed in section 2.7. Finally,

two Lemmas which are needed in later chapters are proved in section 2.8.

2.1 NOTATION AND TERMINOLOGY

All vectors are column vectors or n x 1 matrices -where n is
the dimension of the vector —-and are denoted by underlined small-case
Roman letters. Capital Roman letters are used to designate matrices,
the dimension of which is Specified unless clear from the context.

Superscripts are used to denote the components of a vector or a
matrix: a scalar component of a vector is denoted by a single super-
script (x1 is a component 0f.§): a scalar component of a matrix is
designated by a double superscript (bij is a component of B), and a
column of a matrix- 'being a vector —-is denoted by an underlined
small-case Roman letter with a single Superscript (b; is the jth column
or vector component of B). This allows the use of subscripts to dis-
tinguish between members of a set, such as a sequence. The transpose

of a vector or a matrix is denoted by a superscript T. A superscript

asterisk, xf, is used to denote an optimum vector-— a local optimum
or a global optimum-—-where optimum is defined according to the par-
ticular control problem at hand.

The scalar, t, is used to denote time. If a particular vector,
say x9 is a function of time, xfit) is its value at time t, and x(-)
is used to designate the whole function. The symbology xfi-) is asso-
ciated with a Specified time interval over which the function is de—
fined. Differentiation of a variable with respect to time is normally
written as that variable with a dot over it:

3(t) = dgét)

The notation “x“ is used for the Euclidean 2-norm of the vector
n
i 2 1/2
”35“ = x
i=1
The inner product of two vectors is
n
i i
(3.92.) = 2 :xy
i=1

The scalar signum function, sgn (a), is defined by

-1, a

sgn (a) = { 0, a

1, a

VIIA
OOO

 

and the vector signum function, sgn (x), is defined by

sgn<x1)

ass (>5) = . -
' n
Sgn(x )
2.2 SYSTEM DEFINITION
Because of the nature of the material to be presented in this
thesis, it is necessary to present that material as it applies to each
of three different systems. Each of these three systems is a special

case of the state equation.
1 = gag) (2.2.1)

8f f
where EiTand gfi-exist, and the m—dimensional control vector, u(-),

lies in the set
i
Q = {Elu : 1, i=1, ..., In} (2.2.2)

for all time. An admissible control is a measurable function whose
range is the set 9. The system state at t = 0 is designated 30.
The most general system to be discussed is described by the non-

linear state equation
3 = £(x) +B(t)u (2.2.3)

f
where gi-exists and E is admissible. A special case of this is the

time-varying linear system given by
x = A(t)§+B(t)2 (2.2.4)

where u_is admissible. Finally, the first system to be discussed is

 

the time-invariant linear system given by

a} = A_x_+B_u_ (2.2.5)

——

where u is admissible.

The material to be presented on optimal computing procedures is
much more specific for linear systems and thus each system (2.2.5),
(2.2.4), (2.2.3) is treated separately in Chapters 3, 4, 5 in increasing

order of complexity.

2.3 THE REACHABLE SET AND CONTROLLABILITY

Given a system S, its initial state £0, a set 9, and a run time
T, the set of all possible final states is called the reachable set
for £0 and T and is designated R(xb,T) or just R(T). Such a (S, O)
is called completely controllable if for any initial state ED and any
final state 5f’ xfeR(xo,T) for some finite T :_0. (8,9) is called 229:
_p1etely null controllable if for any x0, QFR(§0’T> for some finite T.Z 0.
Finally, for (5,9), the initial state x is said to be in the region

_0

of null controllability if QeR T) for some finite T_: 0. In the

(£09
rest of this dissertation, and unless otherwise specified, it is aSSumed

that go is in the region of null controllability.

2.4 THE CONTROL PROBLEM
The general control problem to be considered in this dissertation
is sometimes called the optimal regulator problem. The problem is to

transfer the system state from at t=0 to the point closest to the

£0
origin in the fixed run time T :_0 while using an admissible control.
This final state is called_x*(T) and the control vector used to attain

it is called uf(°). Such a point is called an optimal final state

and must be a member of the reachable set. Let T* be the shortesr

run time required to reach the origin. Then as long as T i T*, an
optimal final state must lie on the boundary of the reachable set.
This final state is called a local optimum if it has a neighborhood
which contains no points in the reachable set which are closer to the
origin; it is called a global optimum if no point of the reachable set
lies closer to the origin —-thus, a global optimum is also a local
optimum.

Determination of xf(T) is considered the primary problem and this
thesis is directed toward that end. The dissertation focuses on deter—
mining x*(T) for a variety of T's. For simplicity this is designated
Problem 1. The related problem of finding ué(~) is shown to be solv—
able immediately in terms of xf(T). The time optimal regulator problem,
called Problem 2, is the problem of determining T* —-the minimum run
time for which xf(T) =‘g.

Now there is an infinite set of control problems in category 1,
each one designated by its particular run time T. The solution of a
particular Problem 1 is then labelled x;(T) where the subscript denotes
an element from the set of run times. If the set {5;(T)} is known
for all T, the solution T* of Problem 2 can readily be found. The
method to be presented for solving Problem 1 generates x;(T) for all
T, 0 :_T :_T*. For T < T*, the solution of the control problem can
be used to find uf(T,'), the optimal control vector for run time T.

This is shown after the presentation of the Maximum Principle.

10

2.5 THE MAXIMUM PRINCIPLE

It is assumed that the reader is familiar with the Maximum Prin-
ciple, at least to the extent which it is covered in an introductory
text such as Athans and Falb [1] or Lee and Marcus [10]. This sec—
tion is presented only for the purpose of applying it to the problem
already defined.

For the control problem as Stated and the system defined in

(2.2.1), the cost functional is taken to be
ME) = 1/2 ”5(1)”2 . (2.5.1)

This can be rewritten as follows:

T

ME) = 1/2 ll£(0)||2 + 1/2] j—t “£(t)”2 at
o
T
= 1/2 “5(0)”2 +f gm x(t) dt . (2.5.2)

0

An equally good cost functional and the one which is used in the

sequel is

T
Jl(u) = f gangs) dt . (2.5.3)
0

since 1/2 “50“2 is a constant, the solution obtained using Jl(u) is

the same as that obtained using J(u), so the use of J is, indeed,

1
equivalent.

The Hamiltonian is formed as

H = (3+5, 1), (2.5.4)

 

11

where

é(t) = -:—x = —[g+AT(2+§)] . (2.5.5)

f
In (2.5.5) A is A(t) = g:-— this becomes A(t) or A in linear systems.

An adjoint vector is defined by

X. = 2.+ §,, (2.5.6)
so
0 O O T
y_ = B.+ x. = — A 2} (2.5.7)
This allows the Hamiltonian to be rewritten as
(2.5.8)

R(fislsfl) = < 193E >

where the dependence of H on the control vector u_is through 5. Since

the target set is all of R, the end—point transversality condition

states that (Athans and Falb, page 288)

*
31”) = o , (2.5.9)

so
(2.5.10)

gm = gm

The adjoint system is defined by (2.5.7) and (2.5.10).

In order for 2* to be an optimal control vector, it must minimize

the Hamiltonian according to
H(§*.z*.g*) : H(§*,1*.2) (2.5.11)
for all admissible B: If the system can be described by (2.2.3), this

can be combined with (2.5.8) to yield

 

12

EI(T,t) = ESE [EBT(t) y;(t)] (2.5.12)
for ts[0,T].

2.6 NORMALITY AND THE OPTIMAL CONTROL VECTOR

The control problem defined above is said to be normal if each
element of BT(t) y;(t) is zero only at a countable number of times in
the interval [0,T]. This is a necessary requirement if (2.5.12) is
to be a unique representation of the optimal control almost everywhere
throughout that time interval. A control problem which does not meet
the normality condition is said to be singular. Hereafter, all con-
trol problems of form (2.2.3) which are encountered are assumed to
be normal. In a normal control problem, each optimal control u1 is
either +1 or -1 almost everywhere. This is referred to as "bang-bang"
control.

The symbol Ef(T,°) is used to denote the control vector of a
local optimum. A nonlinear or a singular linear problem may have
more than one gf(T,-); however, in the case of a normal linear system,
“23(T,°) is unique and no ambiguity exists (see Athans and Falb page
404). For the normal control problems considered here, each component
of uf(T,°) is a bang-bang function.

The time t when the bang-bang function ui*(T,t) changes from +1

to -1 or from -1 to +1 is called a switch time. One of these switch

 

times occurs whenever an element of BT(t) y;(t) is zero. So the switch

times are defined by the equation

BT(t) z;(t) = 0 (2.6.1)

13

which is the equation of a moving hyperplane in the n-dimensional

costate Space. This hyperplane is called the switch hyperplane or

 

just the switch plane. Some element of 25(T,-) switches when y;(°)

 

intersects the switch plane. For some systems, there is a maximum
number of switches that can occur as stated in the following theorem

(see Athans and Falb, pages 402 and 403).

Theorem 2.6.1. Consider the normal control problem for the linear

 

time-invariant system

a} = A35+B_u , (2.6.2)

where the eigenvalues of A, Al, A2, ..., An, are all real. Let
uj*(T,0) be a component of the optimal control vector (if it exists).
Let ti be the switch times t of the bang-bang function uj*(T,t). Then
the maximum value of the number i is (n-l). In other words, each on-

off control uj*(T,°) can switch at most (n-l) times.

2.7 DEPENDENCE OF E*(T,t)0N T

In the previous sections, uf(T,t) has been discussed as a function
of t; it is also necessary to investigate uf(T,t) as a function of T.
The double-integrator plant, i = u, which is discussed in Appendix I,
is referred to extensively in the following material to illustrate
the concepts.

Attention is first focussed on a second-order normal linear system
with a one-dimensional control and an A matrix having real eigenvalues
—-the specific example being i = u. This means that only one switch

can occur .

 

14

To understand why it is necessary to examine the effect of run
time, consider a system such that for T = 1 second the optimal final
state is reached by letting u = -1 throughout the one second run. Con-
sider further that for T = 1.01 seconds the optimal final state is
reached if u = -l for 0 f_t < 1.009 and u = +1 for the rest of T. An
equally likely alternative is that u = +1 for 0 :_t < 0.001 and u = -l
for the rest of T. The interesting aspects of this situation are:

1) it is possible that for a sufficiently short run time T, no
control switching occurs, and

2) when T is large enough so that a switch occurs, that switch
could be either near the beginning or near the end of the
run.

Consider the double—integrator plant with x0 = [i]. If u = -l

is applied for one second, the final state is [165] and this is the
optimal final state for T = 1 second. Also for any run time less than

one second, the optimal final state is reached by applying u = —l with

no switching. Point 1 above has now been verified by an example. Now

consider that T is greater than one second, say T 1.2 seconds. In

-1 from zero to

this case the optimal final state is reached if u

something between 1.0 and 1.2 seconds, then u +1 for the remainder

of the 1.2 seconds. This example and several other cases of T > 1

are shown in Figure 2.7.1. The points on the dashed line are optimal

final states for this particular x0. This is an example of the switch

occuring near the end of the run time and is called Type I switching.
Taking again the double integrator plant, let x0 = [Di]. Now

for T small enough, the optimal final state is reached by letting

u = +1 with no switching. Let T0 be the maximum value of T for which

the optimal u does not switch and note that in this case it is some-

what less than 0.5 seconds (see Figure 2.7.2). Also, it should be

15

Switch
,Curve

Figure 2.7.1. The Double Integrator Plant withgg0 -[1]

Switch
Curve

x24 \

\
19 250

 

 

1

 

 

5 ‘ +4 "
--l I l
\ \ X
\ \ TO /
V
/
-1 qr 30
/
/
/
u - +1

Figure 2.7.2. The Double Integrator Plant with 3D -[:1]

-l

 

 

16

observed that a precise value for T0 is not obvious in this case

as it was in the previous one. For T > To, the optimal final state

is reached by letting u = -1 for a time and then switching to u = +1
for the rest of the run time. By taking T slightly greater than T0
note that this is an example of the switch occuring near the beginning

of the run time. This is called Type II switchigg.

 

Two special cases are now evident. The first is forgg.0 on the
x1 axis. In this case, the optimal control must be switched for any
T > O. For completeness, it is arbitrarily stated that a switch exists

for T = O and T0 = O. The other special case occurs when 50 is an

element of the switch curve. In that case, no switch is necessary be-

cause lies on a trajectory through the origin and this trajectory

50
can be followed for the whole run time. For completeness, it is arbi-

trarily stated that T0 = T*.

Whether Type I or Type II switching occurs for T > T0 in the
double integrator plant depends only on the location of_}_{O as shown
in Figure 2.7.3. The marginal or transition cases are also shown and
labeled according to the assumed values of To.
Consider a more general case such that for an arbitrary fixed

Ts[T0,T*] one and only one switch occurs in the function u*(T,-).

Let the time t when this switch occurs be denoted by r - r(T). The

 

notation suggests that the actual switch time 1 depends on the specific

run time T. When T = T the definition of T implies that the switch

0’ 0
must occur at one end of the interval [0,T0]; thus either 1(TO) = T0

or 1(T0) = 0.

17

Switch
Curve

 

 

    
   

TYPE I
[0 < To < 1*]
0 7’>
TO=0 xl
TYPE I s. [0 EYEE <HT*]
[0 < TO < T ] 0

Figure 2.7.3. Switching Types According to the Location of go
in the Double Integrator Plant

The above assertions can be shown by noting that u*(T,°) switches
when y;(-) crosses the switch plane. For T < To, y;(°) does not cross
Ithe switch plane, hence no switch occurs. The trajectory y;(') is
"attached" at one end to the set of optimal final states x;(T)I0:T
§T* , as shown by (2.5.10). So as T increases y;(~) can be pictured
as a curve which is defined over longer time intervals moving through
space-— this movement specified by the position of its endpoint y;(T).
A typical example is shown in Figure 2.7.4 where the state and costate
spaces are superimposed so x;(T) and y;(-) can be shown simultaneously
and where the switch hyperplane is assumed stationary for simplicity.
Assuming the set of optimal final states is continuous - this assump-

tion is discussed later —-and knowing that y§0(fi) does not intersect

 

l8

* .
the switch plane, ZTO(°) must indeed intersect the switch plane at an
endpoint. For the example shown in Figure 2.7.4 r(T0) = To, but this
is not necessarily the case as shown in Figure 2.7.5 where 1(Tb) - O.
2 2 -

x,y + /

Switch

/'
Hyperplane,’

/

 

 

 

 

Figure 2.7.4. The Set of Costate Trajectories

As already stated, the situation depicted in Figure 2.7.4 where

1(T0) = T is referred to as a Type I switch, and the other situation

0
where r(To) = 0 is referred to as a Type II switch. In either case

1(T) is not defined for T < To, so a value is arbitrarily assigned:

let r(T) - 1(T0) for 0 §_T < T ; note that for T in this range no

0

switch occurs. This arbitrary assignment makes r(T) continuous at T0,

and in no way affects the later results, because in neither case does

x203+

 

l9

 

350
* * 7* < >
x T = 1.
40(0) homo) Tl
*
y. ( )
T2
*
zT ( ) z
o . / ,
/' ’ /
* ff a’ Switch
1 ( ) , / Hyperpiane
T3 /
" * (0)
1T cg.» VTO
/ I
I I I
o I . 1’
X ,y1
Figure 2.7.5. The Alternate Situation

T(T) fall in the interval (0,T) until T :_TO.

versus T plots are shown in Figures 2.7.6 and 2.7.7.

 

 

 

r? {TI-1T
<:S§§§§§§

Figure 2.7.6.

I Versus T

For a Type I Switch

 

Sketches of typical T

In each case

:_T :_T*, r(T) must lie in the shaded trapezoidal region.

 

Figure 2.7.7.

[II—YA T2 *I
* 0 ‘T T

T Versus T

For a Type II Switch

 

20

Theorem 2.7.1. For a normal control problem Such that the set of
optimal final states {5;(T)IO_: T_: T*} forms a continuous curve in
state space, r(-) is a continuOus function in any sub-interval of

[0,T*] for which a switch occurs.

m: From (2.5.10) 13501) = gm, thus the set of optimal final
costates also forms a continuous curve in costate space. From (2.5.7),
the adjoint equation is a linear differential equation. The contin-
uity of y;(T)IO :_T :_T* means that for each Tls(0,T*) and for each

s > 0, there is a 6 > 0 such that
1L7; (Tl) — gm] < a (2.7.1)
1
whenever [T - Tll < 5. The switch time 1(T) is given by the equation
bTy* (r(r)) = 0 (2.7.2)
—- T

This means that if the costate equation is run in reverse time over
an interval of time whose length is (T — T(T)) from an initial value
equal to y;(T), the costate at the end of that time interval will be

an element of the switch hyperplane given by
<21) = 0 (2.7.3)

Since the costate equation is linear, trajectories with "nearly equal"
initial conditions will lie "close" to each other. More precisely
Since y; (T1) lies in an e—neighborhood of y;(T) when the initial time
T1 lies in a G—neighborhood of T, for each y > 0 there exists an e > 0
and thus a 6 > 0 such that

MAT”) + IT _ Tll> — z;(T(T)>n < V (2.7.4)

21

whenever [T-Tll < 6. Therefore y;(r(T) + IT-Tll) lies in a Y-neigh-

borhood of an element of the switch hyperplane. Furthermore, I; (r (T1))
1

is an element of the switch hyperplane by definition. So it

must also be true that for each m > 0 there exists a y > 0 and thus

a 6 > 0 such that
l(r(T) + IT - Tll) - r(Tl)l < w (2.7.5)
whenever [T - Tl] <25. This means that
Ir(T) - 1(r1)| < w + IT - T1] = 0 (2.7.6)

Tb summarize for each pl> 0 there must exist an m > 0 and thus a 6 > 0
such that |T(T) — r(Tl)I < p whenever IT - T1] < 6. This has been
shown for an arbitrary T1

plane, so it must be true for all such run times.

such that 2; (°) intersects the switch hyper-
1

Q.E.D.

Since it was assumed in the proof that y;1(-) crossed the switch
hyperplane, this proof shows that r(T) is continuous on those inter-
vals of T for which switching occurs. This means that it is possible
that for short run times no control switching occurs, for longer run
times switching does occur, and for yet longer run times switching does
not occur again. But, over the interval for which switching does occur,
r(T) must be continuous. Also, if control switching occurs for the
run time T < T* but not for T +, then T(Tl) must be either 0 or T

1 1

This follows from the fact that costate trajectories which start close

1.

together must remain close together throughout.

22

Even with the condition of continuity required in Theorem 2.7.1,
it covers a wide range of problems; namely, all linear systems [3]
as well as many nonlinear systems. In particular, consider a nonlinear
system with the global optima such that {5;(T)IO :_T i T*} is discon-
tinuous. Even in this system it is quite possible that there are
(possibly several) local optima each of Which presents a continuous
set of "Optimal" final states. This point is discussed further in the
section on nonlinear systems.

The following Lemma, inserted without proof, is a Special case of

Theorem 7.1 found on page 22 of Hestenes [8] ‘where a proof is offered.

Lemma 2.7.2 Implicit Function Theorem
Given a function f(T,r) which is continuous and has a continuous
derivative with respect to r on an open set S in T, r-Space and such

that
f(T,T) = 0’

suppose

g§_

3T (TO’TO) # 0

f(TO,TO) = 0,

holds at a point (TO’TO) in S. Then there are a continuous function
r(T) on a neighborhood T' of T0 and a constant e > 0 such that
1(TO) = To, f(T,T(T)) = 0

and such that the relations

f(T,r) = o, lr - 1(I)| < e (T in T')

 

23

hold only in case I = T(T). If the function f is of class C(m) on

(m)

S, the function T(T) is of class C on T'.
This Lemma could have been used to prove Theorem 2.7.1, but the
approach used there was intended to be more instructive. In any event,

the Lemma is needed in the following.

 

Theorem 2.7.3. For any subinterval of [0,T*] over which 3x;(T)IO§t:T*£

is differentiable and switching occurs, 1(T) is also differentiable.
Proof: In Lemma 2.7.2 let f be defined as

f(T.T) = <12. gm» = blxiflr) + + bnzT“*(r) .

is differentiable, y;(T) is also. Therefore, When

 

 

Because [§;(T)

 

T represents the switch time, f(T,T) = O and-5; is continuous where S
is all of T,r. Also, as long as T is defined —-i.e., as long as a

switch occurs —-there is a set (TO’TO) where

_ E
f(T0,T0) - 0, Br (T0,TO) # O .

Then, from the Lemma, there is a neighborhood T' about TO on which

T(T) is differentiable. Since there is such a T0,r0 pair for each

time T in the subinterval of interest, T(T) is differentiable through—

out that interval.

2.8 PRELIMINARY MATHEMATICS
The purpose of this section is to present a few lemmas which are

needed in later proofs.

 

24

Lemma 2.8.1

 

lim .3;
ar+o 6T

(X(T + 6T) - X(T)) = A(T) X(T) (2.8.1)
where X(T) is a fundamental matrix defined by

.x(1) = A(T) X(T) (2.8.2)

with

X(O) = I (2.8.3)

Proof: This reSult follows trivially from the fact that the left hand

side of (2.8.1) is the definition of R(T).

 

Q.E.D.
Lemma 2.8.2
T+5T
lim' iL-x<r + 6T) X(-t) B(t) dt = B(T) (2.8 4)
5T+0 5T .
T

where X(T) is defined above and B(t) is any continuous matrix.

Proof: Since X(t) and B(t) are continuous in t, the mean value theorem

can be applied to yield:

T+6T
11m if X(T + (ST — t) B(t) dt

6T+0 6T
T
_ lim ;L_ _ _
_ 5T+0 [6T X(T + 5T r) B(T) (T + 6T 1)] (2.8.5)

where T :_T §_T + 6T.

25

Also,

T+6T

lim 1

6M Ff] X(T + 5T - t) B(t) dt
T

lim
GTWO X(T + 6T - T) B(T)
= X(T - T) B(T) = X(O) B(T) = B(T) (2.8.6)

because, as 6T+O, I must approach T.

Q.E.D.

A special case of this is:

Lemma 2.8.3

 

T+6T _At

lim .3; A(T+6T) dt =

(ST—>0 6T e I e I (2.8.7)
T

Proof: Since eAt satisfies (2.8.2) and (2.8.3), this follows immedi-

ately from Lemma 2.8.2 with B(T) = I.

Q.E.D.

3. THE TIME-INVARIANT LINEAR SYSTEM

This chapter is concerned entirely with the time—invariant linear
system as defined in Chapter 2. The reason for starting with a Special
case is two-fold: first, this problem is easier to handle mathematic-
ally and easier to understand physically than more general cases.
Second, the results are more extensive for this case than in the more
general cases. All of the results that appear in the more general
cases do appear in this case. Thus, this chapter aids the understanding
of later material.

The trajectory system is defined in section 3.1, and its differ-
ential equation is derived. A fact necessary for the solution of the
trajectory system differential equation is presented in section 3.2.
Three computational procedures are discussed in section 3.3, and the
results of some example problems are presented in section 3.4. The

results of the chapter are summarized in section 3.5.

3.1 THE TRAJECTORY SYSTEM
The state equation for the time-invariant linear system has already

been introduced in Chapter 2:
g_ = A§.+ 33 (2.2.5)

where g_is admissible and the initial state x0 is in the region of
null controllability.

The objective here is to develop a means of finding the set of
optimal final states x;(T)IOԤhT :_T* , which is called the solution

trajectory. This is done by developing a set of differential equations

 

whose solution is the solution trajectory. Thus, this set of equations

26

 

27

is called the trajectory gystem. The independent variable of the

 

trajectory system is T, the run time in the original system. It is
necessary to pay particular attention to the differences between T
and t in the following, because they both appear in the development
of the trajectory system equations.

The first step in developing the trajectory system equations is
to describe the location of a single point of the solution trajectory.
Such a point is the final point of an optimal trajectory of the orig—
inal system. An arbitrary point on that optimal trajectory, given

at the time T, satisfies the equation:
I
l§;(T) = eAT x + J{ e-At Buf(T,t) dt (3.1.1)
0

where the first argument of 3% denotes the run time and thus the par-
ticular problem, and the second argument denotes the particular point
of the trajectory. It is necessary that T :_T. The only point of
(3.1.1) known to lie on the solution trajectory —-besides x0, triv-
ially —-is the final point obtained by setting I = T. For convenience,

the primed notation
E (T) = 3:;(T) (3.1.2)

is used to designate a point of the solution trajectory. Thus, x'(T)

is the state variable of the trajectory system. Now,

T
gm = e 350+] e’At 33*(T,t)dt . (3.1.3)
o

 

28

Having obtained an equation for a point on the solution trajectory
as a function of T, it is necessary to put that equation in the form
most convenient for use in finding the whole solution trajectory.
That form is a differential equation obtained by differentiating
(3.1.3) with respect to T. This cannot be done by Leibnitz' rule be—
cause 2f(T,t) is not necessarily continuous, so the following approach
is used.

The differential equation which is sought is given in (3.1.15).
As the first step in the derivation of that differential equation x'

is written for a Slightly larger run time as

T+6T
5' (T + 5T) -_- eA(T+‘ST) £0 +f e"At 33*(1‘ + 6T,t) dt (3.1.4)
0
Now, let
(533' (T) = §'(T + (ST) -3{_'(T) , (3.1.5)
63*(T,t) = 3*(T + 6T,t) — _t_1_*(T,t) , (3.1.6)

and Subtract (3.1.3) from (3.1.4) to get
T

533'(T) = eA(T+<ST) _ eAT £0 +[ e-At 33*(T,t) dt

0

T+6T
+5 —
+ eA(T T) J{ e At BBI(T + 6T,t) dt
T

T
5 _
+ eA(T+ T) Jr e At Bdu*(T,t) dt
0

29

T+<ST
5 _

T

T
+ eA<T+5T)f e-At B63*(T,t) dt . (3.1.7)
o

If GT is chosen Small enough so that no control switching occurs in
the interval (T,T+6T), the control vector in the second term becomes

a constant —-call it gfe(T) —-and can be removed from the integral:

T+6T
eA<T+5T) I J“ 133*(1 + 5T,t) dt
T
T+6T
eA<T+6T) J[ e—Atdt Bu?(T) . (3.1.8)
-—e
T
From (2.5.10) and (2.5.12)
* T * _ T I
Ee(T) = SGNWg-B §¢(T)} - §§N_{-B x_(T)} . (3.1.9)

Substituting (3.1.8) into (3.1.7), dividing by GT, and taking

the limit as 6T vanishes yields

A6T
lim éx'gr) _ °, _ lim e ~ I ,
(ST->0]: (ST ] " 5 (T) ‘ 5m:[m5T 132(1)

Big: (T)

T+6T
+ "5'1— eA(T+6T) e-Atdt
T
T

T
+ If eA(T+6T) f e-At B<‘311_"‘(T.t:)dt }. (3.1.10)
0

30

Equation (3.1.10) is the differential equation of the trajectory
system, but each of the three terms on the right hand side must be

simplified. The first term can be reduced by a polynominal expansion

 

 

of eAGE:
2 2
A GT
11111 eAGT-I = lim (1+A5T’“ 2! + )-I ..., A
5T90 6T 6T+O 6T “
(3.1.11)
BY Lemma (2.8.3), the second term becomes
T+6T
11m _:_1_ A(T+OT) -At * = *
6T->O are I 9— dt BEeCT) Bgefl) (3.1.12)

T

To evaluate the final term in (3.1.10), first start by considering
a scalar control u* with only one switch. This switch occurs at T(T)
in u*(T,°) and at r(T+6T) in u*(T+5T,-). 6u*(T,t) is zero except be—
tween these two switch times, where it has a magnitude of two. Thus,

the term can be reduced using the above and the mean value theorem to:

T
lim lg; eA(T+6T) j“ e-At‘b§u*(T,t) dt

(ST-+0
o
(T+<ST)
_ lim _1_ A(T+6T) -At
’ arr—>0 are I 9 P112) dt
(T)
__ lim _1_ A(T+<ST-t') _
' (ST-+0 are E(:Z)(T(T+5T) 1(1)) (3.1.13)

where T(T) : t'_: T(T+6T) or T(T+6T) j_t' j_T(T), and the Sign of
5u*(T,t) is not yet known. Letting k take on one of the values :1,

(3.1.13) becomes

 

31

lim A(T+6T-t') T(T + 6T) - T(T)
6M [2kg 3 6T ]

A(T-T(T)) lim :(T + 6T) - T(T)
2ke 2 6T+O[ or ]. (3.1.14)

r(°) has been shown to be continuous in Theorem 2.7.1. On any interval
for Which it is also differentiable, (3.1.14) becomes 2kT(T)eA(T-T(T))b,

If u*(T,°) switches more than once, this becomes

A(T-t:2(T)) b + ...

2k1t1(T) eA(T‘tl(T)> 3+ 2k2t2(T) e

Also, if gf(T,-) is of order m.: 1 and each element switches once,

it becomes

m

E : ' A(T-t (1‘))
2k1t1(T) e 1 hi
i=1

*
where bi is the 1th column of B and ti(T) is the switch time of u1 (T,°).

*
In general, each element of u (T,-) can switch any number of

times, and each switch generates a separate term as above. Call the
total number of switches r: r does not depend on m. Now, (3.1.10)

reduces to

r
323(1) = 115(1) + 1322(1) + E Zkitifl) eA(T‘ti(T)) 131 (3.1.15)
i=1

This is the differential equation of the trajectory system.
AS in Figures 2.7.3 and 2.7.4, each input has some run time such
that for shorter run times it does not switch. For each of the terms

in the summation in (3.1.15), let T designate this threshold run

10

 

32

time. In (3.1.13) the jth term of the summation is zero until

T Z'T This fact is incorporated automatically by defining ti(T)

30'
as in section 2.7 thus making ti(T) = O for T < T10.
Equation (3.1.15) plus the initial condition given by §f(0) = go
defines the trajectory system. For clarity, the results of this sec-

tion are summarized in the following theorem.

Theorem 3.1.1. Given the system (2.2.5)

 

x = A__)_{_+B_u

where u is admissible and xfiO) = ED is in the region of null control-
lability, the set of optimal final states of the various optimal reg-
ulator problems.{x;(T)IO :_T :_T*} is coincident with the solution

of the trajectory system differential equation given by

3(1) = A_x_'(T)+Bu:(T)+E ti(1)11(1) (3.1.17)
i=1
where
x'(0) = x0 , (3.1.18)
3:31) = sgNg-BTg'ms (3.1.9)
and
31(1) = 2kieA(T'ti(T)) 31 . (3.1.16)

This theorem has been proven above. The main problem at this point
is the difficulty in solving the trajectory system differential equa-

tion due to the terms of the summation.

 

33

3.2 EVALUATION OF ti(T)

Equation (3.1.15) indicates that the trajectory system differen-

tial equation can be expressed in terms of ti(T) and ti(T). Before

proceeding to methods for solving this system for the locus of optimal

final states, a further fact which is useful for computing ti(T) is

presented in the following Theorem:

Theorem 3.2.1

 

<§'(1>.xi(1>> = 0.1=1. .
for all T 3_T10
Proof:
. AT<1-ti(1>) .

-Xt<t1(T)> = e zTc1)

from (2.5.7). And from (2.5.10)
T _
z;(ti(1>) = eA (T ti(T))Ԥ'(1)

But, by the definition of the switch time t1(T),
b eAT(T-ti(T)) *(T)> = 0
<-—i’ 11

Thus,

T
o = <1»... eA (Him) as»

A(1-t,(1)) b >

<35'(T), 9

Now, multiplying both sides of (3.2.5) by 2ki yields (3.2.1).

Q.E.D.

(3.2.1)

(3.2.2)

(3.2.3)

(3.2.4)

(3.2.5)

 

 

34

3.3 COMPUTATIONAL SOLUTION OF THE TRAJECTORY SYSTEM DIFFERENTIAL
EQUATION

Three methods are presented in this section for finding the locus
of Optimal final states by solving the differential equation of the
trajectory system. Method 1 is noniterative and well suited to use
on a digital computer. Method 2 is noniterative and can be used on
an analog computer as well as a digital computer. The last method is
called Method 2A because it is an extension of Method 2. Its use is
advantageous in certain limited situations which are discussed.

Method 1. The most straightforward method of solving (3.1.15)
is by discretizing (3.1.15) and finding the points of the solution
trajectory sequentially. Given §f(T) and each ti(T) for some Ts[O,T*]

and an increment AT, xf(T + AT) is given approximately by
gu+An =gwm+gwnAT was)

where xf(T) is given by (3.1.15) if the ti(T)'s are known. But the
t1(T)'s can be found by choosing values such that when each t1(T + AT)

is computed by
ti(T + AT) = ti(1) AT + ti(T) (3.3.2)

(3.2.1) will be satisifed at T + AT. This procedure does involve an
iteration on the set of ti(T)'s, but the method is straightforward
and even lends itself to hand computation with few enough switches.
An example is given later in the chapter.

In order to obtain the solution by Method 1, it is necessary to
have a means of computing a correction factor. That is, assume a

problem with only one switch time 1(T). At run time T,le(T), 1(T)

35

and KIT) are assumed to be known. Further it should be true that

(35' (T), 1(T)) = 0. The first step was taken at T = 0 or from T = O
to T = AT. To do this an average value of T, say T(AT/2), would have
been used. At the second step, either T(3AT/2) could be computed and

used, or it could be computed but

, H921) + «13%
T(T) = 2 (3.3.3)

 

could be used for the next computation. The latter method is used
here, so at run time 13 it is also assumed that T(T - %;5 is known.

Now the value of x'(T + AT) is found by guessing a value of

T(T + 駧 = T and then using

1
§j(T + AT) 2 §f(T) + g}(T) AT
= x'(T) + AT Ax'(T) + bu:(T) + T1(T) 1(T) (3.3.4)
where
%1(T) — %-(%1 + T(T --%§)) (3.3.5)
also
1(T + AT) = 2keA[T + M ' HT) ’AT'ElflH 3 (3.3.6)

Now let the inner product of xf(T + AT) and vflT + AT) — which

is a function of T only —-be represented by the symbol 1(T1). If

1

12 can be found by using T2 = T1 + AT, it is true that

I(T2) = 0 (3.3.7)

36

 

 

 

 

 

but
I(T) = I(%)+3—I—A%+~- . (3.
2 l 8T1
Using only the terms in (3.3.8),
. I(% )
AT 2 - 311 (3.
ail
where
v v
LI = < 3’5 (T + AT) , v(T + AT)> +(x'(T + AT), EMT." AT) >
3% 8f - - 31
l l 1
(3.
From (3.3.4) and (3.3.6)
BK' + A AT
(g. T) = 7m) , <3.
1
and
32(T T AT) = --éI-AV(T + AT) . (3.
8T 2 -
1
Thus,
BI __ AT v
'éf—l - .2— {<X(T), X(T + AT)> -<_)£ (T + AT): AXCI‘ + AT»); (3-

By substituting (3.3.15) into (3.3.9) a correction factor can be

3.8)

3.9)

3.10)

3.11)

3.12)

3.13)

obtained

which leads to T2(T + 9—21), and the iteration is continued until the

inner product becomes as near zero as is desired.

An exception occurs where the first step is being made. In

that

case T(- %;5 does not exist, so T(T) = T(O) is computed, and later it

is called T(%§). In this case, instead of (3.2.13) the following is

used:

 

37

3‘:- = AT (X(T): X(T + AT» -<§,'(T + AT) . A10“ + AT)> (3.3.13a)

For the more general case where there is more than one switch,
values of t1(T) are guessed for each i. Then an iteration is carried

out over one ti(T) at a time with the others held constant. This makes

(3.3.4) become

a’ (T + AT) = x'(T) +AT{A§'(T) +_b_u:(T)

r
+ tkm 1km + 2 tin) 11m} (3.3.14)
i=1

i#k

Since the terms of the summation are constant while the iteration is

. 1
being carried out over tk(T), W is still given by (3.3.11).
k

So the correction factor for each tk(T) is found from (3.3.9) and

(3.3.13).

Since the solution trajectory is being computed sequentially from
T = O to T = T*, it is necessary that there be some method for detecting

a threshold switch time T10 when it occurs. Also the value !i(TiO)

must be determined. Considering the latter first: from section 2.7,

Tio , type I switch

ti(TiO) = (3.3.15)
0 , type II switch

This is used in (3.1.16) to find

2kihi , Type I switch

11(T10> = (3.3.16)

ATiO
2k e b Type II switch

1 —i ’

38

If the ith switch is a Type I switch, then it first occurs at

*
run time Tio where zTiO(T10) is an element of a switch hyperplane.

Thus, xf(TiO) is also a member of that switch hyperplane, so T10 can
be detected for a Type I switch by observing when u:(T) switches.
Also,by noting which element of u:(T) switches, Xi(TiO) is found to

be Zkihi where only ki must be determined.

If instead the ith switch is a Type II switch, then when it first

ATT

* *
occurs at T10, the starting point XTiO(O) e 10 XTiO(TiO) is an

element of a switch hyperplane. This T10 is detected by realizing that
AT
<§f(Tio), e 102i:>= 0. The second element of the inner product can

be found by defining

31(T) = eATgi (3.3.17)
Now,
éi(T) = AE1(T) , (3.3.18)
or
§i(T + AT) = §i(T) + AT(A§1(T)) (3.3.19)

Thus the solution trajectory is computed from‘xO at T = 0, and simul-

taneously (3.3.19) is solved from 2_‘t_>__i at T = 0, then T10 is the shortest

run time for whichm(xf(T), 51(T)> = 0. Also, when T

10
_ ATiO
-!i(TiO) - Zkie Bi is known except for the signed factor k1.

is reached,

At the time T = ti(T) 31(T) is added to the already existing

TiO’

equation for xf(T). At that time ki can be determined. The direction

of the additional vector term is known without k ki merely designates

19
the sense of that vector. Thus, at T = T10, k1 is chosen so that the

additional vector —-which is a part of xf(T) -has the same general

39

sense as those already existing terms of x'(T). In other words, k1

is chosen so that
0' - .
(35 (T10 ), ti(TiO) 11(Tio» > 0 (3.3.20)
The steps of Method 1 are shown in the following algorithm:

Algorithm 3.3.1

1. At run time T, xf(T), (T), j=1, ..., m, Xi(T)’ and ti

2

..j
AT _

(T - 7?), i=1, ..., r, are known. Let k—l.

AT

2. Let ti(T +—2—) = t1(

(T +AAT) using (3.3.6) and then find xf(T + AT) using (3.3.14).

AT
T -'3-), i=1, ..., r. Find eachli

Find each Zj(T) using (3.3.19).
___ v
3. Evaluate Ik <(vk(T + AT), x (T + AT)>.

4. If IR is not sufficiently small, improve tk(T + %§5 using

(3.3.9) and (3.3.13), and return to 2. When Ik is suffi—
ciently small, increase k by 1 for k<r and let k=l for k=r

then return to 2. Repeat until each I is sufficiently small.

k

5. Determine whether a threshold switch time has been reached
by checking for switches in 2:(-) or a sign change in<<§f(°),

§j(-)) from T to (T + AT). If no new switch, return to 1.

In case of a switch set X(TO) according to (3.3.16), deter-

mine ki from (3.3.20), and return to 1.

Method 2. Another method for computing 31(T) has been developed
in order to avoid the iteration necessary in Method 1. This method

is suitable for use on an analog computer. Differentiating both sides

40

of (3.1.16) with respect to T yields

31(T) = 2k1(1 - ti(T)) AeA(T‘ti(T))b

(1 - E1(T)) A11(T) (3.3.21)

where 3&(T10) is determined according to (3.3.16). If R1 and ti(T)
are known and if it is known whether the ith switch is a type I or a
type II switch, 31(T) can be found using (3.3.21) and (3.3.16).
Theorem 3.2.1 may be used to simultaneously remove the need to find
ti(T) in (3.3.21) and to determine t1(T) for use in (3.1.17). This
is accomplished by solving (3.3.21) for 21(T) where t1(T) is chosen
to be whatever value it must take on in order to satisfy (3.2.1). In
this way ti(T) is computed implicitly and can also be used in (3.1.17).
Using an analog computer, the factor (1 — ti(T)) would be generated
with a bootstrapping technique. This would become an iterative process
on a digital computer. The constant ki can be found using (3.3.20).

As long as ti(T) is continuous, (3.1.16) shows that 31(T) is
also continuous —-in fact, differentiable. Thus (3.1.17) shows that
.x'(T) is also differentiable where 2:(T) is continuous. Let the r

switch times be_arranged so that T < T < ... < T Each element

10- 20— - rO'
of 33(T) is piecewise constant, switching only at a threshold switch

* _
time T10. Thus for 0 :_T < T10"2e(T) is constant and each !i(T) — O,

, * .

so x_(T) is differentiable. For T10 < T < T20, ue(T) is also constant
and as long as ti(T) is continuous, 11(T) and thus x'(T) are differ-
entiable. This fact plus Theorem 2.7.3 show that 33(T) is continuous
if and only if t1(T) is continuous. Therefore, except at the times
T10, §f(T),'!1(T) and ti(T) are all differentiable. Furthermore,

they are all continuous for 0 §_T i T*.

 

41

As long as sufficient computer capacity is available, x'(T) can
readily be found by solving the equations of Method 2 on an analog
computer. The required analog capacity can be estimated by the number
of integrators and multipliers required. Equation (3.1.17) requires
n integrators and (3.3.21) requires (n x r) integrators. If T is
desired, an analog integrator is needed. And r integrators are needed
to find the t1(T) from ti(T). Thus, the total number of analog inte-
grators required is at least (n+1) x (r+l). Further, n multipliers
are needed for calculating each of the inner products in (3.2.1) and
(3.3.20). Thus, for example, in a second order system with one input
switch, six integrators and four multipliers are needed.

But, it was also noted in the previous section that in order to
locate the Type II switches it is necessary to solve (3.3.18) for each
11(T). Thus, in addition to the above number of integrators, (m x n)
more integrators are needed to search for additional Type II switches.
It is only if the number of Type II switches is known in advance that
these are not needed —-each can be put to use as soon as the last Type
II switch is found. An example of the use of an analog computer for
the solution of the problem is given later in the chapter.

Method 2 can be programmed on a digital computer using Algorithm

3.3.1 with the following changes: instead of (3.3.6) use
2(T + AT) = 3(T) + AT [(1 — T(T)) A3(T)] , (3.3.21a)
instead of (3.3.12) use

v + A
3—(E% T) = - gg'AXKT) , (3.3.21b)

l

 

42

and instead of (3.3.13) use

3 _AT

3, - —2-{<y_(T). X(T + AT>>

-<£'(T + AT). A103»; . (3.3.21c)

But this has no computational advantage over Method 1 in the case of

a linear time-invariant system.

Method 2A. In the case of a second-order system, Method 2 can
be simplified and the number of integrators reduced. This new method
is called Method 2A. In such a system, each Xi(T> is normal to x'(T)
and thus they are all parallel to one another. Define the vector

3(T) to be a unit vector normal to xf(T). Thus,

(B(T), 3(T)> = 1 (3.3.22)
and
<B.(T),£'(T>> = 0 (3.3.23)
define g(T). One possible vector is
xz'm
_ 1
R(T) - W 1' , (3.3.24)
-x (T)

the other is its negative. Now define a set of scalars 2ci(T)£ so

that

110:) = ci(T) B(T), i = l, r . (3.3.25)

For convenience, choose the polarity of 3(T) so that

43

Ci(T10) > 0 . (3.3.26)
Now from (3.2.1)
(x'(T), 11(T)> = 0 (3.3.27)
so it follows —-using (3.3.21) —-that

o = fiq'm, 11m» = <3'(T>,11(T)>
+<3_E'(T).i1(T)> = <:'(T>,11(T)>

+<§'<T). [1 - Elm] A1101»)

(gm + [1 — 131(1)] AT§'(T), 11m» (3.3.28)
But, using (3.1.17), this becomes

0 = (Elm. A35'(T) + 333T) + [1 - Elm] AT§'(T)

r .
+ Z tl(T) 11m) . (3.3.29)
i=1
This reduces to
0 = elm (3(1). A3<_'(T) + 323T)
o T r o
+ [1 - tl(T)] Ax'(T)>+ Z tl(’I‘) ci(T) . (3.3.30)
i=1
Since c1(T) # 0,
I‘ . *
2:31 t1(T) ci(T) = mum, A: (T) + Bge(T)

+ [l - t1(T)] AT_x_'(T)> . (3.3.31)

44
But (3.1.17) can be rewritten as

r 0
5y (T) = A_x_'(T) + 333T) + {:1 tl(T) ci(T) B(T) (3.3.32)

and this can be combined with (3.3.31) to yield
323m = 3(T) -<n(T>. 3m) + [1 - Elm] AT§'(T)> nu) (3.3.33)

where

_q(T) = Ax'(T)+Bu:(T) . (3.3.34)

If r is at least two, fewer integrators are needed because (3.3.21)
needs to be solved only for one set -31(T) —-and that only to deter-
mine [l-t1(T)]. But n additional multipliers are required to compute
the inner product in (3.3.33). The additional multipliers would not
be a problem on a digital computer. One potential difficulty is that
if 31(T) ever becomes zero, (3.3.31) does not follow from (3.3.30).
Another problem is that the switch times can not be computed — except
t1(T)- although this does result in a savings of an additional r
integrators. Also, the switch times are not part of the solution of
the problem defined in Chapter 2.

Theorem 2.6.1 gives the maximum number of switches for a system
in which A has real eigenvalues. If A has one or more complex con-
jugate pairs of eigenvalues, a large number of switches could exist.
Method 2A is useful as an efficient means of solving such a problem
for a second-order system. Method 2A can not be extended to higher
order systems because it depends on the fact that in a second-order

system there is a unique normal to x'(T).

45

3.4 COMPUTATIONAL RESULTS

Method 1 and Method 2 have been used to solve some problems. The
results are presented in this section. A computer program using
Algorithm 3.1.1 to solve Problem 1 and Problem 2 for any system with
a scalar control was written and is documented in Appendix II. It
was used to solve Problem 1 and Problem 2 for the double-integrator
plant by Method 1; the computed loci of optimal final states for four
different initial conditions are shown in Figure 3.4.1. The solution

to Problem 2, T*, as a function of is given in Athans and Falb -

x0
Equation (7-26), page 514. These values are compared with the com-
puted results in Table 3.4.1. The computation times on an IBM 1800
were about 20 seconds. It should be noted that only a moderate accur-
acy was required from the program and no attempt was made to minimize
either computation time or error. It might also be noted that the
initial conditions chosen represent Type I switches, Type II switches,
and a point on the switch line (hyperplane).

Method 2 was used to solve Problem 1 and Problem 2 for the double-
integrator plant on an analog computer. The initial condition con—
sidered was 30 =[3], which is one of the initial conditions shown in
Figure 3.4.1. The results of this analog solution are shown in Fig—
ure 3.4.2.

The program in Appendix II using Method 1 was also used to solve

Problems 1 and 2 for the system:

i“) = 3.0:) + u(t) (3.4.1)

46

 

 

 

0.4..

0.2"
0 . —+* : i ; p»
0.2 0.4 0.6 0.8 1.0 x1

-O.2«~

-O.4«»

 

Figure 3.4.1.

Table 3.4.1.

Solutions of Problem 1 for
the Double—Integrator Plant

Error in the Solutions to Problem 2 for
the Double—Integrator Plant

 

 

 

 

 

 

xOT Actual T* Computed T* Z Error
(1.00, 0.00) 2.000 2.019 0.95
(0.80, 0.50) 2.425 2.454 1.19
(0.60, -O.25) 1.340 1.329 0.82
(1.00, -0.50) 1.620 1.604 0.99

 

 

47

o. .I

flo.wgn x cues ucmam HoumquucH

IOHQSOQ mo coflusaom “mesmEou woamq< .N.q.m madman
comfiumano pow

coausaow Hmufiwfiv u x

15.0..

 

A
4

1b
4P
1P
4P

vn . .
H OH w 0 0.0 «.0 N.O

 

48

where

In! :_1 (3.4.2)

for several initial conditions. The results appear in Figure 3.4.3.
This system is discussed in Athans and Falb pages 526-536, but no
results useful for comparison are given there.

Method 1 was used to solve Problem 1 and Problem 2 for the system

 

 

 

 

given by
I- " F ...
-l
gt) = O _x_(t)+ 1 u(t) (3.4.3)
0 —l 1
b “ l- d
where
lul : 1 (3.4.4)
and
11.5
$0 = 1.3 . (3.4.5)
2.0
t J

 

 

This system is discussed in Athans and Falb, pages 536-551, but none
of the results given there are useful for comparison. The solution

trajectory computed for the problem above is shown in Figure 3.4.4.

3.5 SUMMARY

The trajectory system has been defined and a differential equa-
tion (3.1.15) describing it is derived. Three methods which use
(3.1.15) along with (3.2.1) to find the locus of optimal final states
for the problem defined in Theorem 3.1.1 are discussed in section 3.3.

Method 1 and Method 2 are equally suitable for use on a digital comr

49

3.133 coaumsvm mo Emumzm 9.3
now N new A wEmHnoum ou maOHuDHom .m.q.m shaman

 

 

long 1.0.01
tones... 3.0 u as
on; u e
um
11¢.OI
auN.Ol
34 u e.
um
A .1 v A o .«r
as 04 we o.o so To 0
K
4:

50

A.m.q.mV :ofiuwsvm mo Eoumzm poppo
usage was now N was H gemstone ou wcoauseom

.q.q.m ouswam

 

 

N H o a m .e. H o
4 141 .A1 s 4 1
AV H .74“
.ooo oHN.m u B
an
1.N 11m

 

 

 

') 5].

computer via Algorithm 3.1.1, but Method 2 can be used on an analog
computer as well. Method 2A is develOped for second-order systems
with several switches. The data resulting from solution of some

example problems is presented in section 3.4.

4. THE GENERAL LINEAR SYSTEM

Chapter 3 is concerned with the time-invariant linear system.
In this chapter, the results of Chapter 3 are extended to the case of
a general or time-varying linear system. The trajectory system dif-
ferential equation is develOped in section 4.1, and the normality of
Xi(.) is attached in section 4.2. Section 4.3 discusses the exten—
sions of Method 1 and Method 2 to the time-varying case. The solution
of a simple problem in section 4.4 demonstrates the computational pro-

cedure developed. The reSults are summarized in section 4.5.

4.1 THE TRAJECTORY SYSTEM

From Chapter 2, the state equation for a general linear system

can be given as
3'5 = A(t)3t_+l3(t)3 (2.2.4)

where 2_is admissible and the initial state £0 is in the region of
null controllability. Also, the matrices A(t) and B(t) are assumed
continuous. The object of this section is to develop the trajectory
system equations for this general linear system. The method used will

be similar to that used in Chapter 3. The results are stated in the

following theorem.

Theorem 4.1.1. For the linear system defined by

 

g = A(t)_x_+B(t)_u_ (2.2.4)

._0
null controllability, the set of optimal final states of the various

where_u is an admissible control and x(0) = x is in the region of

optimal regulator problems {3;(T)|O g T;g T*} is coincident with the

52

53

solution of the trajectory system given by:

for) = A(T) 5(1) + B(T) gm
r
+ E 2ki ti(T) X(T - ti(T)) bi(ti(T)) (4.1.1)
i=1
where
x'(O) = 350 (4.1.2)
and
2:”) = _S_C3_1‘I_{-BT(T) 32'”); (4.1.3)

Proof: As in Chapter 3, for the run time T, the point on the Optimal

trajectory corresponding to the time T is given by:
T
35:0) -.- x(T) 50 +f X(-t) B(t) 2*(T,t) dt (4.1.4)
0

where O :_T f_T and 0 §.T.i T*. X(T) is the fundamental matrix defined
in 2.8. Also, the arguments of E? are as in Chapter 3. The point of
that optimal trajectory which is also on the solution trajectory is

the final point:

T
33' (T) = X(T) 50 +f X(-t) B(t) 3*(T,t) dt . (4.1.5)
0

Also, for a slightly longer run time

T+0T
x' (T + OT) = X(T + 6T) 1‘0 +f X(-t) B(t) 3*(T + 6T,t) dt
0 (4.1.6)

54

Subtracting (4.1.5) from (4.1.6) yields

65' (T) = [X(T + (ST) - X(T)]

T
x + X(-t) B(t) E*(T,t)dt
{-0 f. }

T+ T
+ X(T + (SDI X(-t) B(t) 3*(T + 6T,t)dt
T
T
+ X(T + 6T)f X(-t) B(t) 63*(T,t)dt (4.1.7)
0

If GT is chosen small enough so that no control switching occurs in
the interval (T, T+6T), the control vector in the second term becomes

a constant and that term can be rewritten as

T+6T
X(T + OT) ] X(-t) B(t) 3*(T + 61:, t)dt
T

T+6T
= X(T = OT) f X(-t) B(t)dt 93(1)

T (4.1.8)

where 32(T) is as given in (4.1.3) by (2.5.10) and (2.5.12).
After substituting (4.1.8), dividing by OT, and taking the limit

as 6T vanishes, (4.1.7) becomes
111“ 935.01.) '1
6T->O [ 6T ] 35 (T)

_ lim x(T + 6T) - X(T) _ 1
" 390“: (ST ]X( T) i (T)

1 T+6T
+ E X(T + 5T) j; X(-t) B(t)dt 113T)

T
+ 311—, X(T + on] X(-t) B(t) 52*(T.t)dt} (4.1-9)
0

 

55

The first term of the right-hand side can be reduced by Lemma 2.8.1,
the second term by Lemma 2.8.2, and the third term as follows —-after

section 3.1 —-:

T
T 0 _-

6T+O
11m 1 r ti(T+<5T) J:
= — X(T + 6T) (2k ) X(-t) E (Udt
5TI° {ST [521 ti(T) i 1
_lim ifix(T+OT-T)b(‘r)t(T+5T)‘t(T)
' 6T+O i=1 6T 1 -i 1 i i
r .
= 12:31 Zki ti(T) X(T - ti(T)) 21(t1('r)) (4.1.10)

where ti(T)«: T .: ti(T + GT) or ti(T + 6T) :_T §_t1(T), i=1, ..., r

i i
and k1 = id. Now (4.1.9) reduces to (4.1.1). Equation (4.1.2) is

obvious.

Q.E.D.

By the same argument posed in section 3.3, xf(T) and ti(T) are contin-

uous except at the times T10.

4.2 SOLUTION OF THE TRAJECTORY SYSTEM EQUATION
The trajectory system differential equation (4.1.1) must be solved
in order to find the locus of Optimal final states. To simplify con-

sideration of the terms of the summation, a set of vectors are again

defined according to:

11m = 21.1 X(T - ti(T))1_)_i(ti(T)) , (4.2.1)

 

56

where 31(T10) is still given by (3.3.16). As in the case of the
time-invariant linear system, there are two basic ways of computing
31(T): either it can be computed directly from (4.2.1), or its dif-
ferential equation can be solved concurrently with (4.1.1). The dif-
ferential equation for 31(T) can be found by differentiating both

sides of (4.2.1) with reapect to T to get:
551(1) = [1 - £100] A(T — tim) 11(1)
+ ti(T) X(T - ti(T)) 21(ti(T)) (4.2.2)

Whichever method is used for computing Xi(T)’ an important fact

necessary for the solution of the trajectory system equations is given

by:

Theorem 4.2.1. For the system defined in 4.1,

 

<:xf(T),'gi(T)) = 0, i=1, ..., r (4.2.3)
for all T.: T10.
Proof:
* _ T *
main» - x (T — ticrn) 1T0)

XT(T - ti(T)) gm (4.2.4)
from (2.5.7) and from (2.5.10). By the definition of the switch time,
@1410», z;<ti<rr)>> = 0 (4.2.5)

Thus

0
l

<—b-i(ti(T))’ xT(T - tim) 1' Cr»

<§'(T). X(T - tim) hictim» (4.2.6)

 

57

Now (4.2.3) is obtained by multiplying both sides of (4.2.6) by 2R1.

Q.E.D.

4.3 COMPUTATIONAL METHODS
As stated above, there are basically two ways to compute 31(T).
Each leads to a different computational procedure. The two procedures

presented here are very similar to those presented in section 3.3.

Method 1. As in section 3.3, Method 1 involves the direction
computation of x(T) from its defining equation (4.2.1). If A(°) and
3(0) are known functions of time, then X(°) can be computed. Solu-
tion of (4.2.1) is then no more difficult with a digital computer than
was its counterpart in the time-invariant linear system. But the over-
all method still involves finding T(T) by an iterative technique; so,

a means for improving the estimate of T(T + €95 is needed.
The correction factor is formed in the same way as it was in

section 3.3. But since

1(T + AT) = Zki X[T + AT - T(T) - AT T1(T)] 241(1) + AT l°r1(T)].
(4.3.1)

there is a more complicated expression for

3% 2

M2. .. 2k1{- 93- A[T + AT - T(T) — ATT1(T)]
1

X[T + AT - T(T) - AT-E1(T)] EMT) + AT%1(T)]

+ X[T + AT - “[(T) - ATT1(T)] EEC?) + ATT1(T)]}

 

. - %1{A[T + AT - T(T) - ATT1(T)" x(T + AT)

- 2k1X[T + AT - T(T) - AT%1(T)] §[t(T) + ATTIUI)”
(4.3.2)

58

If b_is constant, (4.3.2) becomes the same as (3.3.13) except that
A(°) must be evaluated at the argument shown. In that case, the cor—

rection factor is changed very little. For the more general time-

varying b
31 AT .
F = 7 {<3flT) . x(T + AD)
1 .

- (at) (T + AT), A(T + AT - T(T) - ATT1(T)) 1(T + AT))
+ 2ki<3£'(T + AT), X(T + AT - T(T)

— AT T1(T)) _13_(t(T) + AT T1(T)))} (4.3.3)

must be used to determine the correction factor. This requires that
b(°) be differentiable. The solution of the trajectory system differ—
ential equation by Method 1 is very similar to the procedure outlined
in Algorithm 3.3.1. The similarities between the two versions of

Method 1 are even more pronounced if B is a constant matrix.

Method 2. As in Chapter 3, Method 2 is a solution method in
which (4.1.1), (4.2.1), and (4.2.2) are solved simultaneously —-even
on an analog computer. As in the case of Method 1, Method 2 is more
easily adapted to the time-varying linear system if B is a constant
matrix, because that removes the second term from the right-hand side
of (4.2.2). But even if B(-) is not constant, it must be differen-
tiable. Differentiability of B(°) is not necessary for the solution
of (4.1.1), only to use (4.2.2) in solving (4.1.1). Evaluation of
the first term of (4.2.2) requires evaluating A(°) at (T — ti(T))'
Doing this on an analog computer would probably be too involved, so
even Method 2 may as well be considered a digital computer method.

Again, this requires an iteration over ti(T).

59

If B is a constant matrix, Method 2 can be used with Algorithm

3.1.1. In this case

X(T + AT) = y_(T) + AT{[1 - an] A(T - T(T)) X(T)}. (4.3.4)

 

Thus,
32L§o+ AT) = - é-T— A(T - rm) x(T) . (4.3.5)
T1 2
making
3% = A}{<X(T + AT), X(T)) -<35' (T + AT). A(T - T(T)) x(T)): (4.3.6)
1

Now Method 2 is seen to be far superior to Method 1 because the fund-
amental matrix does not appear anywhere in the process. Instead, only
the known matrix A(o) must be evaluated in the system equations or in
the correction factor equations. When B is a function of time, the
second term on the right hand side of (4.2.2) shows that the fundamental
matrix must be used in the problem solution, so the superiority of

Method 2 is not extended to this case.

4.4 AN EXAMPLE PROBLEM
To illustrate the use of Method 2 in Algorithm 3.1.1 for a con—

stant B matrix, the double—integrator plant was changed to the system:

. 0 1 0
x(t) = x(t) + u(t) (4.4.1)
0 t 1
where
|u| : 1 (4.4.2)
and
1
x0 = . (4.4.3)

 

 

60

This problem was solved using the program in Appendix II. The results

are shown in Figure 4.4.1.

4.5 SUMMARY

The trajectory system is extended to the time-varying linear case
(4.1.1). The solution of (4.1.1) and (4.2.3) by Method 1 and Method 2
are discussed in section 4.3. In particular it is shown that Method 2
is far superior to Method 1 when B is a constant matrix. A simple
extension of the double-integrator plant is used in section 4.4 to

demonstrate the use of Method 2 in Algorithm 3.3.1 when B is constant.

A

 

1b
4
1r:
1
V

-0.2‘P

—O.4«»

 

x = a point from x=u

Figure 4.4.1. Solutions of Problems 1 and 2 for the
System of Equation (4.4.1)

 

5. A CLASS OF NONLINEAR SYSTEMS

In this chapter the results of Chapters 3 and 4 are extended to
a class of nonlinear systems. It should be noted at the outset that
the results of this chapter have not been proven, they are merely a
reasonable extension of the previous results. The problem encountered
is that a proof of these results requires that the set of optimal final
states be continuous. This is definitely not true for all nonlinear
systems. In fact, all systems can be put into one of three classifica-
tions: 1) some systems —-inc1uding all linear systems —-have no local
optima other than the global optimum and it is continuous; 2) some
systems have a discontinuous global Optimum which is made up of pieces
of various local optima each of which is continuous; and 3) some systems
have discontinuous local optima. For systems in the first group, the
methods presented in this thesis work as they are presented. For systems
in the second group, each local optimum can be found by the methods
presented in this paper, but it is then necessary to construct the
global optimum from these local optima. For systems in the third group,
this thesis provides no method of solution. At this time, no procedure
is available for determining a priori to which group a given nonlinear
system belongs. The only approach with a nonlinear system is to assume
it is in one of the first two groups. If it is in the third group,
it is not known that it would be evident that the solution method had
failed. In fact, this problem would be an interesting and probably a
fruitful area for further research.

In section 5.1, the trajectory system differential equation is
developed for a class of nonlinear systems. The convexity of the reach—

able set is discussed in section 5.2. The normality of Xi(°) is

61

62

attached in section 5.3, and a solution method is discussed. A simple

example is presented in section 5.4, and section 5.5 is a summary.

5.1 THE TRAJECTORY SYSTEM DIFFERENTIAL EQUATION
The state equation for the class of systems to be discussed here

has been given in Chapter 2
3'5 = _f_(§) + B(t) 3 (2.2.3)

f
where gi-exists, u_is admissible and the initial states;O is in the

region of null controllability. The matrix B(t) is assumed to be con-
tinuous. The trajectory system equations for this class of nonlinear
systems is developed in this section subject to the assumption men-

tioned above.

Theorem 5.1.1. For the class of nonlinear systems defined by

 

g = :93) + B(t) 3 (2.2.3)

where 2.13 an admissible control, x(O) ='§0 is in the region of null
controllability and the set of optimal final states of the various
optimal regulator problems {3;(T)IO :_T j_T*} is continuous, this same

set of optimal final states is coincident with the solution of the

trajectory system given by:

i'm = gym» + B(T) 1120:)
r o
+ £1 2k1L t1<T> X(T — rim) hicicm (5.1.1)
where
x'(O) = (5.1.2)

To

63

and

gm §c_N{-BT(T) gym} . (5.1.3)

Also X(o) is the fundamental matrix of the system linearized along an

optimal trajectory —-i.e., X(°) is the fundamental matrix associated

3f

‘3‘:

with A(t) = .
§;(t)

Proof: As in the preceding chapters it is first necessary to obtain
6§f(T) = xf(T + 6T) - §f(T). Finding §f(T + GT) is done by first

extending uf(T,T) over the interval [T, T + 6T] to obtain x;(T + OT),
and then applying OET(T,°) = uf(T + 6T,-) - uf(T,-) to obtain x'(T-+ 6T).

The first step is simply given by

* * '* 2
'xT(T + OT) xT(T) + §T(T) OT + 0(6T )

= x*(T) + f(x*(T)) + B(T) u*(T)] ST + 0(6T2) (5.1.4)
—T —-—T —e
where 22(T) is as given in (5.1.3). Now the second step is accomplished
by linearizing the system about the trajectory x;(T) for ts[0, T + OT].
The validity of this linearization depends on the assumed continuity

of the set of optimal final states. From this linearization the fol-

lowing differential equation is obtained:
61:4) = A(t) age) + B(t) 6mm (5.1.5)

where

f
A(t) = -—- . (5.1.6)

 

 

 

64

Solving (5.1.5)

* + = +6
635T” OT) X(T T)

T+oT
535(0) +f X(-t) B(t) 53*(T,t) dt (5.1.7)
0

where X(t) is the fundamental matrix for A(t). Since 6§;(0) = 0, it

is easy to see that

go: + 6T) - g (1) = [161' (1)) + B(T) 113(1)]

T+6T
(ST + 0(6T2) + X(T + 5T) f X(-t) B(t) 62*(T,t)dt. '
0 (5.1.8)

If the interval GT is chosen small enough so that no control switching
occurs in [T, T + 5T], then the integral is zero over that interval.
Taking that into account, and dividing both sides of (5.1.8) by OT and

taking the limit as 5T vanishes,

3'03) = gym) + B(T) 113m

r

+- E 2ki tl(T) X(T - ti(T)) 2i(Ti(T)) (5.1.1)
i=1
where k1 = :1. The limit in the last term is identical to (4.1.10).

Q.E.D.

5.2 NONCONVEX REACHABLE SETS
Note that the equations developed in the previous section are

valid regardless of the convexity of the reachable set. To see why,

65

assume the reachable set can have two lobes as shown in Figure 5.2.1.
Each lobe is a locally convex portion of the reachable set and each
contains an extremum. These extrema are connected by the curves A
and B. The trajectory system would follow one of these curves.

The fact that two curves are available means that at some point of
the solution of the trajectory system—3gs —-a choice of two possible
switches must be made. One choice results in curve A, the other choice
in curve B. Thus, all local optima can be found and the global optimum
can be constructed from them. One important point to be made here is
that the assumption in Theorem 4.1.1 is that each $2221_optimum is
continuous. As long as this is true, the global Optimum need not be

continuous.

.34 I,

A?

$9:

_....—-R(t)

 

 

w

Figure 5.2.1. A Nonconvex Reachable Set

66

5.3 SOLUTION OF THE TRAJECTORY SYSTEM EQUATION

As in the linear systems, define the vectors
31(T) = 2ki X(T - ti(T) bi(Ti(T)) , (5.3.1)

i=1, ..., r. It is still true that each 31(T) is normal to the tra-

jectory system state vector:

Theorem 5.3.1.

 

<1i(T)’ x'(T)) = 0, i=1, ..., r (5.3.2)
for all T Z-Ti0°
Proof: From (2.5.7)

i§<t) = -AT(t) z;(t) (5.3.3)

where the matrix A(t) is the same as that given in (5.1.6). Thus from

(2.5.10)
z;(ti(T)) = XT(T - t1(T)) y;(T) . (5.3.4)
But, by the definition of the switch time ti(T),
(simian), XT(T - ti(T)) y;(T)> = o . (5.3.5)

Thus,

0
ll

<p_i(t,(T>>. xT<T - tim) 35' (10>

<35} (T). X(T - tim) 111mm» . (5.3.6)

Q.E.D.

67

All of the results to this point are identical in form to those
of the time-varying linear system of Chapter 4; but it should be noted
that, although 21(T) in (5.3.1) looks identical to that in Chapter 4,
its computation is far different. Specifically X(T - ti(T)) must be
evaluated along the original system trajectory instead of the trajec-
tory system trajectory.

This problem can be circumvented as in Chapter 4 by letting B be
a constant matrix and using solution Method 2. Differentiating (5.3.1)

yields

gin) 21.1(1 - ti(T)) A(T — ti(T)) X(T - ti(T)) 311

(1 — t1(T)) A(T - ti(T)) 31m (5.3.7)

In Algorithm 3.1.1,

 

£1” + AT) = 11(T) + AT {(1 - ti(T)) A(T - ti(T)) Kim} (5.3.8)
where
8v (T + AT)
—4. AT
8,11 = - —2— A(T - ti(T)) 11(T) (5.3.9)
and
fi-I— = AZT—{<11(T + M). zi<T>>
i1

-<2_<.' (T + 4T), 40: - ti(T)) 11(1)); . (5.3.10)

It must be remembered that in the above equations A(°) is defined by
(5.1.6) and must be evaluated accordingly. This is the difficult part

of the solution.

68

5.4 AN EXAMPLE PROBLEM

In order to provide a simple example problem, the double-inte—

grator plant is extended as follows:

. 0 l 0
310:) = 2 _}_{_(t)+ u(t)
O -lx I 1
where
lul : 1
and
1
x =
‘0 O

This system is discussed in Athans and Falb pages 614-621.

system

0 1

A(t) 2
O -2[xT*(t)|

and in (5.3.8) and (5.3.10)
0 1

A(T - ti(T)) = - 2* _-
0 2|xT (T ti(T))l

(5.4.1)

(5.4.2)

(5.4.3)

In this

(5.4.4)

(5.4.5)

For a given run time T, this is evaluated by running the original

system backward from §;(T) = xf(T) by ti(T) seconds. At that point,

A(T - ti(T)) can be evaluated. Using a Special subroutine to perform

this calculation, the program in Appendix II was used to compute the

locus of Optimal final states which is shown in Figure 5.4.1.

5.5 SUMMARY

The trajectory system differential equation (5.1.1) has been

developed and the normality of Xi(T) attached (5.3.2).

This pair of

69

equations was solved by Method 2 for a simple extension of the

double-integrator plant.

-0.2<>

4.4+

 

i x = a point from §=u

Figure 5.4.1. Solutions of Problems 1 and 2 for the
System of Equation (5.4.1)

6. CONCLUSIONS AND EXTENSIONS

Section 6.1 summarizes the material presented in Chapters 2
through 5. The conclusions to be drawn from.that material are included
where they are pertinent. A possible extension of the trajectory
system method to other control problems is treated briefly in section

6.2.

6.1 SUMMARY AND CONCLUSIONS

The particular control problems treated in this thesis are defined
in Chapter 2. Briefly, Problem 1 is to determine how close to the
origin a given initial state can be driven within a fixed run time;
Problem 2 is to find the minimum run time necessary to drive a given
initial state to the origin. The method presented in this disserta-
tion solves all possible regulator problems for one arbitrary initial
state.

In this dissertation, the above problems are treated for a class
of linear and nonlinear systems. All of the systems in this class
have two things in common: each element of the optimal control vector
is a "bang-bang" function, and the locus of Optimal final states is
continuous for each local optimum. These properties are crucial to
the computational procedure. Instead of solving the optimal control
problem indirectly by solving the associated two—point boundary value
problem which comes from the application of the Maximum Principle,
this procedure solves the optimal control problem directly using the
general form of the results which the Maximum Principle makes avail-

able.

70

71

Considering the set of control problems which are collectively
labeled Problem 1, the specific Optimal control function depends on
the particular run time involved. The form of that dependence is in-
vestigated in section 2.7 and is one of the more important portions
of this dissertation. In particular, each switch time is shown to be
a function of the run time and each switch time ti(T) is shown to be
differentiable whenever the locus of optimal final states is differen-
tiable.

The locus of optimal final states can be computed by solving
simultaneously the differential equations of the trajectory system
(3.1.15), 4.1.1), (5.1.1) and the normality condition associated with
the set of auxiliary inputs {yi(T)} (3.2.1), (4.2.3), (5.3.2). This
locus of optimal final states is the complete set of solutions of
Problem 1 and readily results in the solution of Problem 2.

Method 1 and Method 2 presented in section 3.3 are two computa-
tional procedures which can be used in Algorithm 3.3.1 to solve the
trajectory system equations. In Method 1, each 31(T) is computed
directly from its defining equation (3.1.16), (4.2.1), (5.3.1). In
Method 2 each Xi<T> is generated from its differential equation instead.
These two methods are shown to be computationally equivalent for most
of the systems considered; however in a time-varying linear or a non-
linear system with a constant B matrix, Method 2 is far superior to
Method 1. Both methods describe noniterative computation procedures,
but an iteration is involved when using either procedure on a digital
computer. Only Method 2 can be utilized on an analog computer.

A simple Fortran program which uses Algorithm 3.3.1 is presented

in Appendix II. Either Method 1 or Method 2 can be used in this

72

program, but it only applies to systems with a scalar input that
switches at most twice. This program was used to solve the example
problems in sections 3.4, 4.4, 5.4. The use of an analog computer

to solve an example of Problem 1 in section 3.4 shows the feasibility

of analog solution.

6.2 GENERALIZATION AND FUTURE RESEARCH POSSIBILITIES

In this section the trajectory system approach to the solution
of Optimal control problems is generalized and then applied to the
minimum fuel regulator problem. Consider an optimal control problem
wherein the objective is to apply permissible inputs to a system so
as to drive its state from a given initial state to the closest pos—
sible point to the origin in a fixed run time and with minimum cost.
Assume the set of Optimal final states resulting from all possible
costs 0 §_J :_J* —-where J* is the minimum cost necessary to reach
the origin with an adequate run time —-forms a continuous curve in the
state space.

As in Chapter 3, consider the linear system

1': = Ax+Bu (6.2.1)

where 2_is admissible and x(0) = Let the fixed run time satisfy

£0.
T < T*. The equation for the Optimal final state is

3<_'(J) = $0) = eAT

T
150 “Ff e'Ath*(J.t) dt (6.2.2)
0

73

where the notation is similar to that in the preceding chapters.

Now (6.2.2) must be differentiated with respect to J.

Consider the minimum fuel regulator problem where

T m
J=F=f§|uildt
0

i=1

From the Maximum Principle, it can be shown that

3*(F,t) = -_QEZ_{B'y_;(t)}
where
—1, a_: -l
dez(a) = 0, -1 < a < l
, a 3_l
and
dez(al)
DEZ(a) = '
dez(an)

Now the trajectory system differential equation is

dx*(F) . r .
:g-f- = 35' (F) = E ti(F) kieMT'timnh
i=1

But since it is also true that

1:0“) = 8(5) 1:09

(6.

(6.

(6.

(6.

(6.

(6.

.3)

.4)

.5)

.6)

.7)

.8)

74

where u(F) is some unknown function, it can be shown that

l

' = — =
|<£ (F), li(F)>l 0‘07): 1 1: ..., r (6.2.9)
where
A(T-t (F))
= i
31(F) kie bi . (6.2.10)
With no fuel consumed and a run time T,
9ST) = eATlco ,’ (6.2.11)
so in the trajectory system
gm) = eATltO . (6.2.12)

Now the locus of optimal final states cannot be found because
u(F) is not known. This type of extension of the trajectory system

method would be a good area for further research.

 

 

 

 

 

 

BIBLIOGRAPHY

8.

10.

11.

12.

13.

BIBLIOGRAPHY

Athans, M. and P. L. Falb, Optimal Control, McGraw-Hill Book Co.,
New York, (1966).

Barr, R. 0., Computation ofygptimal Controls by quadratic Programming
on Convex Reachable Sets, Ph.D. thesis, UhiVersity of Michigan,

 

Barr, R. 0., and E. G. Gilbert, "Some Efficient Algorithms for a
Class of Abstract Optimization Problems Arising in Optimal Con-
trol," IEEE Trans. on Auto. Cont., AC-l4, pp. 640-652, (1969).

Bellman, R.,_Dynamic Proggamming, Princeton University Press,
Princeton, N. J., (1957).

 

Bryson, A. E., and W. F. Denham, "A Steepest-Ascent Method for
Solving Optimum Programming Problems," J. Appl. Mech. Ser. E,
V. 29, (1962), pp. 247-257.

Fancher, P. S., "Iterative Computation Procedures for an Optimum
Control Problem," IEEE Trans. on Automatic Control, AC-10,
pp. 346-348, (1965).

Gilbert, E. E., "An Iterative Procedure for Computing the Minimum
of a Quadratic Form on a Convex Set," SIAM J. on Control, Ser. A,
V. 4, No. 1, pp. 61-80, (1966).

Hestenes, M. R., Calculus of Variations and Optimal Control
Theory, John Wiley and Sons, Inc., New York, (1967).

Ho, Y. E., "A Successive Approximation Technique for Optimal
Control Systems Subject to Input Satruation," Trans. ASME, J.
Basic Eng., Series D, V. 82, (1960), pp. 33-40.

Lee, E. B., and L. Markus, Foundations of Optimal Control Theory,
John Wiley and Sons, Inc., New York (1967).

Neustadt, L. W., "Synthesizing Time Optimum Control Systems,"
J. Math. Anal. and Appl., V. 1, pp. 484-492, (1960).

Pontryagin, L. S., V. G. Boltyanskii, R. V. Gamkrelidze, and
E. F. Mischenko, The Mathematical Theory of Optimal Processes,
John Wiley and Sons, Inc., New York, (1962).

Stratton, R. B., Computation of Optimal Controls for Nonlinear
Systems Via Geometric Search, Ph.D. thesis, Michigan State

 

University, (1969).

75

APPENDIX I

 

 

 

 

APPENDIX I
THE DOUBLE-INTEGRATOR PLANT

The double-integrator plant is so named because it is represented

by the differential equation

x = u (I-l)

This equation may arise, for instance, in describing the reaction of
a unit mass when acted upon by.aforcecurthrust u. Such a system with

bounded input can be described by the state equations

where u is admissible.

For Optimal time or Optimal regulator control of this system

u = +1 or u = ‘1 are the only two input values of interest. For
u = +1, the state trajectories are parabolas opening along the positive
x-axis; for u = -1, they are parabolas opening along the negative

x—axis. Each of these two sets of trajectories has one member which
intersects the origin. These two trajectories are shown in Figure I.1
with the direction of motion shown by the arrows.

The half of each curve for which the state is moving toward the
origin is in solid line, while the other half is in dashed line. The
two solid half trajectories together form the minimum time switch
curve. Any initial state can be forced to the origin in minimum time
by letting u = :1 —-whichever is correct for that initial state —-
until the switch curve is reached, and then switching u and following

the switch curve in to the origin.

76

77

 

 

 

Figure 1.1. Trajectories for the Double-Integrator Plant

The problem of forcing the system state as close as possible to
the origin in a fixed run time is more complicated. To analyze this
problem, it is only necessary to consider initial states above the
switch curve in Figure I.1, because the trajectories for u = +1 and
u = -1 are symmetrical about the origin. For the rest of this dis-
cussion initial states are assumed to be above the switch curve. For
an initial state in the upper half-plane, the input is originally
u = -1. For short enough run times, the control never switches. When
there is sufficient run time to carry the state into the lower half-
plane, switching does occur. The longer the run time available, the
later the switch occurs and thus the farther the state penetrates into
the lower half-plane on the u = -1 trajectory before switching occurs.
The limit occurs at T* the optimal time when switching occurs at the

switch curve and the state reaches the origin. This is an example of

78

a Type I switch because switching first occurs at the end of the tra-
jectory. The situation is shown in Figure I.2 where the switch curve
and the set of optimal final states are shown as solid lines and the
trajectories as dashed lines.

For an initial state in the lower half-plane, the input is u = +1
throughout for short run times and it is u = -l originally and then
u = +1 for longer run times. This is demonstrated in Figure 1.3.
Again the switch curve and the set of optimal final states are shown
as solid lines while the trajectories are shown as dashed lines. This
is an example of a Type II switch because switching first occurs at

the beginning of the trajectory.

Switch
Curve

79

([34

 

Figure I.2.

Switch
Curve

'01

 

Solutions of the Regulator Problem

with a Type I Switch

x2+

 

Figure I.3.

 

 

Solutions of the Regulator Problem

with a Type II Switch

APPENDIX II

 

 

APPENDIX II

The Fortran IV program presented at the end of this appendix
uses Algorithm 3.3.1 and either Method 1 or Method 2 to solve both
Problem 1 and Problem 2. The program given here can only be used for
a system with a scalar input that switches at most twice. It would
not be a very difficult matter to change the program so that it can
handle any size control vector and any number of switches. It is also
necessary that the B matrix be constant.

The name of the program is OPTML; the various Subroutines are
names PROB, REV, RENEW, SCALE, SWPLN, CHANG, END, DOT and SGN. The
user‘Supplies data cards and parts of PROB and REV. The first data
card must have an integer number in column 1 which is the order of
the system, N. The next N cards contain the elements of the vector
—22; after that, eaCh set of N cards is taken to be an initial state
50. These entries must eaCh be in E format and right justified to
column 15.

In PROB the user supplies the system equations.' There are six
groups of equations which are, in order: 3;: A§.+ by, 3(T + AT) = 5(T)
+ AT Ag, 2(T + AT), Hyp = IB?§A w) and alt. Alt is an expression to
be used in computing y = u:(T) when §(T) is an element of the switch
hyperplane. The vector w_is an element of the switch hyperplane.

The vectOr w_is used in computing the correction factor for T(T).
Using Method 2, XL=.Z.+ AT(l — T)Aw_and w;=‘%-A_. Using Method 1,
1 =1 2e—Asb and w = A 1733.

The subroutine REV is only necessary when the system is nonlinear.

Using the vector A as the final value of the state and the costate,

both sets of equations are run backwards in time an amount TINT. The

80

81

step size AT = DELT is available. At the end, the value of‘x is put
in the location B. The example of PROB shown is for the double-
integrator plant, and the example of REV is for the nonlinear system

of the example in section 4.4.

82

.Ouxq
~U02m~
HUFOZ
o 4

IIOCHDO
Fit—o"?!

m h
accoumdm
mooonhmmo
doonwzoo
AQHV x Amofiemv o<wa CON
.0~.O>u.o~.>
AQ~.O>UADH.N
zeauom com 00 on
Aom.o> .modemv 04mm Mom
zeauom HON Do 00H
oodefideha szmm
2 A¢OH.N. Q<wm
..mozoumm IO.M.>¢. mH ZOHHDLOm mZHP znzmzwz thIomex~v h<zmom mod
AHH.h<sz¢ ¢OH
Amomdm vp<zmom moH
A.xm.m.om.eexflvb<zmom NOH
.\.D<hoz¢eon H
.Amvx1¢.xo~..N.xI¢eon.AHvxr¢exHH.thrmexmfiehrdexmexdv h<zaom Hod
AmmemomdwewomawexHvh<zmom am
x<.Z.OD<Ho.D<hosoaqpez<behmwoeh 202200
A©.ON.AOVN A
aonze>..o.zoxeaovxo.Amva>.Aw.>..ovoxaaov3>eao.x ZOHmZmz~oo
mdmwmbzm 3103 mzow
zqawoxa wumDDm hmumw
”1mHoemwhzmxam¢¢~.omquvmuom*
zqawoma mmwuomazozw
mzhao «Om \\

Emuwoum >H :muuuom .H.HH manna

 

 

 

83

beam.>.x0qu.mH ¢H
ome¢H.@m Ahmqhmva
.zezc>.zoxvhoouxuwzu
on H

eo<zoxeo<za>ehm<ed>xe>emezaxezo>eNeoNexexoeo>.>. zwzwm mmdu

szHHzoo m
A; H

.o.eoz.ea<em.es<.a>:.>.m.zox.zo>.N.oN.x.xo. >.4L. zsczm unao

Ahm<e

tzszou m
memem H>va
HQ>I.ZomU>

a>I.>.m.zox.zo>.~.o~.x.xo. >.:. moms anao

¢uz
wDZthou wH
mHeNHemH .lezamH
szemvzmhdo 44<U H
NH OH 00
h HmoHem. wbHaz ¢
¢eHeH Hhxwzva
bewzemeH .oszex.QZw 4440
tzthou N
Hem em Abmwhva
wzoolzdozxnbmmh
Hazexexvhoovhmomuzmozx
D<hQNOD<ho
.ZeNexvhoouazou
H+Hoznhoz

zmozxeamvxeAm.xeAHexea<pek HmoHem. whHmz

mDZHhZOu o
AHoHamvwkaz
moouaqho
.ou34k

.ouw HF
Hh.Hh.mH .mzoov uH
.thaqhm

H.o.oso6v .H.HH wanna

84

02m

Hme mmqu
wDZHhZOu

OwehH.C> AHIZiva
AzzeoHvxmhqo mmqu

mDZHHZOu
H OH 00
Ahmqeo<zox
eo<za>ea2wpe Q>Ie>emszoxezo>aN.ON.X.XQ. >eboz. ozqru mmqu
mDZHHZOU
0 Oh Do

i12.xo.xoceoe.emamhenmonmzoo

iomeoNuioHcN

iosco>nioic>

Aomioxuiofiex

z.HnoH new as

gaseoussso

.m\io:<so+:<so.um

emaemwe + ixaczumsxswoaqsun<s

snmo+suh

Hm.om.s .exmzeum

ioN.o>.xo.>.eoz .n.es.w.e¢<em.sxmz.H.>.oa>:.mn<um undo
ousxmz

ien<.on>:.>.m.zox.za>.N.o~.ex.xc.>.z. moan ns<u
an:

maziezou

om.om.om licum

Hm.os.s iexinuH

.a~.o>.xo.>.eoz .L.HH.a.eaaem.exmz.o.u::o.n2me.mnqum Esau
Husxmz

o¢.o¢.os locus

iz.oN.ox.eooun2me

LH.om.omiemLe.iH

mawlmmowmnhmmh
xumluw AxuwluvzwmnmeHm

H.v.usoov .H.HH OHAmH

NH

mH
NH

H

HH

>¢h

Hm
om

ON
OH

O¢

0m

5

85

 

02w
zmshma
AH.x Uhm<
mDZHhZGu o
Zmnhwm
.OHANVZ
m\AN.>u.H.3
szszou m
Zxahwx
.mvxlnd>I
wzthzou ¢
zmahma
.m.3u.~.>
.mvzwA..N\.OD<HQ+D<FQVVI.H.whmmc+HH.3nHH.>
szthou m
zmshwm
HNVNUHN.ON
.NVNwhmmD+HH_NnAHVON
tzHhZOQ N
szhwm
x<n ..N\.OD<HD+D<PQ..*AN.>+> "vaxo
x<* A.N\AOD<HQ+D<hQ..%HH.>+ANVXUHH.xo

tzthou H
H..oem.¢.memeH. OH ow
.OHAHVm
14.2.334Ho.D<HQ.OD<#.D<b.hmwoep zozzou
Hovza>sHOVJN..OVX.HO.N.HG.XQ.Hov> .Ao.m..@.3 onmezHo

Ap4<ea>Ie>eme 3azo>eNeoNexexoe>eH.mOaa wthDDamDm
mmmowth 0103 m2:*
deosmd mumDOm hmHmw
z<mooma mmwuomazozw
«Om \\

H.v.u:oov .H.HH OHQMH

86

 

ozm

zmzhmm

.Hvquva

zeHnH o 00

N Ch 00

hqmo+mum

.H.oauHH.a

.H.OquH.x

zeHuH ¢ 09
.HN.Q*.HN.x.zwm*HN.x*oN+HH.dlvnhmwolHNvauHNvoa
HH.auHH.OQ
..NvavzommewolH.N.x.szNHN.waN.x*h4wol.N.quNvax
.viwhmwQIHHvquHvox

memem Hszhtmva

wDZthou

.oum

.H.<UHH.Q

HHL<UHHvx

zeHuH H 00

x<.Z.OD<HoeD<baeoa<hez<hehgmoeh zozzou

.0.0d..0.0xeHovae.o.ero.meHov< ZOHmzwIHo
.szhem.<. >wm wthDOmmDm

wmmwwth Odo

zqawomd muxsom

z<mwoaa mmwuo

i.e.oeooc

3 mzow
hmHLN
mazozw
mom \\

. H . HH anmH.

87

 

02w
sthwm
w4zo>\.oHvO>u.oH.zo> OON
w4ZOX\.oH.Oxu.QHvzox
zeHnoH OON Do
HHzeo>eo>vhoovhmOmuw4ZO>
4H44. Q>Ie>eme >.Z®>.N.0N.xexoe0>ez. mama 4440
mu:
h4molhl034hum
.N\AOD4FQ+D4FQ4%H4wo+b4bn034h
.Hzeoxeaxvhoovhmomno4zox
”~44. a>Ie>emezoxezo>eN.ONexexo. >.Iv moan 444v
Nu:
A HOH.xo.*h4wo+.oH.xn.oH.ax mON
zeHucH mON OD
4H44. Q>I.>.mezoxezo>.Nsowexexo. >.z. moan 4440

H":

14e2.034ho.D4hoeon4pea4heh4moeh zozzou
AHvox. .H.zaxe.Hvzo>.4H.N.4H.ON.AHverHvxo.AHvo>.AHv> onmzwzmo
40x H

.o4zox.04z3>eh44ea>re>emezox.zo>.N.ONexexoea>e>.szwm wZHHDOmmDm
mdmwmth omoz mzow
Z4m00dd wumDOm emu *

moi \\

A.p.ucouV .H.HH mHan

 

88

 

m 34 oo
sumo
muse

.auwaemm

.ouaaemm
ma.~m.s amine
e.mm.m4 amine

o.w.m issue

as as on

41"“

muss

sue

Haemm\~<emmwenmouesmo

mouse

maemmumaemm

Naemmuaaemm

i.monao.\so.mm<w.Haemmnmsemm.+4<emnum<emm

cums

~.m.m iaoooo.1.<.mm<.ua

ms.4.e lessens

i\..0220u 2mmn we: :o4.4_. smmzaz Inseam nosezou1-.\.xH.e<zsou
.~_.m.m4m.m.m4m.x4.eszaou

i.m4m4mmoaz_ mt» QZDOu mat mzspzosmnm mtetos.XH.e<zsod
¥<.z.o:<eo.:aeo.o:<e.:<a.pamo.e zozzoo
laceN..4.o>.iacxo.iae> ZOHmZLzao

.o~.o>.xo.>.eoz.s .44.i.eaaem.exmz.a.n.<.m4aom mziesoemzm

msmoLezs one: mzos
zssooma momDOW smash
zasoosa mmmuoaczozn
sou \\

A.c.uooov .H.HH oases

¢m~or~

non
Hm
om

 

89

 

02m

zmahwa

w.m.4 .Hmem. whHmz
Homem. whHm:

zmahma

OUhxwz

zmahwa

Hluhxwz

szhwa

thxmz
4N4mel.H4\P4wOUH4wo
h4wo+huh

ouH

Cum

0H Oh Ow

Nun
N4hwm\AN4pwml.H.Nb4mOUh4wo
h4wo+hup
p4mo+huon4h

NH Oh Do
AxvoNqu40>

zeHux OH Do

.OUOD4P

oeHHeNm .hxwzva

4 HooHem. whHmz
H+4N4

m.onos4ho

HIUFOZ
AAZ.XQ.>.HOQ.Zwmux4
.oupm4hm

A.s.osoov .H.HH sassy

Nm

mH

4H

0H

mH

NH
HH

OH

w

 

 

 

m.ou:aeo
emcee
oueoz
.oueaaem
4 4oo4.m. mess:
4+4u4
442.xo.>4eooezumu¥<
4o4.3 "4o44x o4m
2.4u04 o4m oo
>11»
444s. a>t.>.m.zox.zo>.N.oN.x.xo. >.z. moms 44cc
4":
4o44xoas4mouio44z uioaex 34w
4444x1404.z
2.4u04 s4~ ca
4444. a>x.>.m.zax.zo>.N.o~.x.xo. >.z. some 4444
4a:
>1u>
mazaezoo m4

% 44.4m .m4. oe oo o4
o4.m .o44>4u4
444aizomu>
444a. a>I.>.m.zox.zo>.N.oN.x.xo. >.z. moan 44<o
on:

4\..oz:ou 244m ma: Io4.44. mmmzsz toeszm 4oaezouzmm.\.x4.e<:mOu oo4
44.2.osaeo.zseo.a:<e.:<e.e4mc.4 . zozzoo
.o43.44.zox.444zo>..44N.44.oN.44.x..44xo.444> ZQLWZLZHQ
44 4

epahozehm4hmeb44ea>Ia>emezcxezo>emeawexexae >._.:Z4n43m wzHHDDmmsm
mmmowtaH omox mzow
Z4mooma womaom ...mH4ow
do”. \\

A.o.osoov .H.HH manna

 

 

 

 

02m
Zmahma

aou+034hou034ko Ho
zmahmm
OD4HQHD4HQ

.N\moo+034bon034ho oo
Hoeooeoo AHOZVmH
AAnimhwb4m04\HZ.ZG>.ZOx4hDQvwoNumou
14%404ZOX\AZ.ZO>.>VHOOVlazezezoxthDUQth
Ab44ea>Ie>eo4zo>ezezo>emeoNexexoa >ezv meme 444D

mu:

X4.Zeoz4hoea4hoeos4hea4h.H4moah 202200
.043. 4HvzaxeAH420>e.H4N.AHVONefiHvxaaHvxoeHH4> ZOHmzmzHo
4H44ew4zox H

.04zo>.a2mhe a>I.>.m.zox.ZO>.N.ONexexoe >ehozvoZ4Iu mthDOmmDm
mmmothH Odo: mzow
Z4mwoma mumaom bmH4N
Z4mooma mmmuomazozw
mom \\

91

OZw

Zmtha
mDZHhZDu m

NUHH

HIUfi

D4honos4ho

A.m.uaoov .H.HH OHQMH

 

92

 

x4ezeas4ho.L4hae334pez4heb4mceh

4hxmzemeH

02m
zsaems
4uexmz
.o4\e4moueu e
41ozwa
ouw
one
.o4we4maue4mo o
o.om.o .xzaxz. 14
4444x42omuxz o4
o.o4.o im><znemmz4 u4
aname.zomuemmz s
242444
ouexmz om
~u02m_
.o4\e4moup4mo
sins
anus
4444x4zomuxz
4czmeezomuw><z s
o2m4..m.s. as oo
losexwazmeuazme ~4~
2.4u04 ~4~ oo
.4uatme
zozzoo
14.x 204m24240
.c2m4.x.ozL mziesoanam
memomsz4 case 4244
£4mw3xa mumDOm hmH4N
sou \\

A.e.osoov .H.HH oases

93'

 

02w
szhwa 4
.Huzom m
4 OH 00
.ouzom N
4 OF 06
.Hluzom H
meNeH Hxva
Hxvzow ZOHHUZDm
mmmowth Odo: mzow
Z4mooma wUdDOm bmH4w
z4mwoaa mmmuomazozw
mow \\

02w
szhwm
HHHmeH44+ hOOMHDQ H

ZeHuH H De

.OUHDQ

HHHmeHH44 ZOHmZmzHC

Hzeme4vbao ZOHhUZDm
mamothH Odo: m2©%
Z4aooma modDOm hmH4w
Z4mooma mmmuOmazozw
mod \\

A.U.ucouv .H.HH OHHMH

 

 

”'111111111111111111111111111“