OPTIMAL CONTROL

Chapter_1 Static Optimization

1.1
1.2
1.3

Optimization Without Constraints, 1

Optimization With Equality Constraints, 5
Numerical Solution Methods, 21

Basics + 2,3D-Plots

The goals of this course are to expose the students to the mathematical tools of
parametric and dynamic optimization and their uses in designing optimally
behaving dynamic systems. The nonlinear systems course covers the following
topics:
1. Static Optimization
2.0 Optimal Control of Discrete-Time Systems
2.1 Linear Quadratic Regulator
2.2 Steady-State Closed-loop Control of Sub-Optimal Feedback
2.3 The Tracking Problem
2.4 Regulator with Function of Final State Fixed
3.0 Optimal Control of Continuous-Time Systems
3.1 Linear Quadratic Regulator
3.2 Steady-State Closed-loop Control of Sub-Optimal Feedback
3.3 The Tracking Problem
3.4 Regulator with Function of Final State Fixed
3.5 Final-Time-Free Problem
3.6 Constrained Input Problem
4.0 Dynamic Programming
4.1 Discrete-Time

Introduction
The state feedback and observer design approach is a fundamental tool in the
control of state equation systems. However, it is not always the most useful
method. Three obvious difficulties are:

The translation from design specifications (maximum desired

over and undershoot, settling time, etc.) to desired poles is not
always direct, particularly for complex systems; what is the best
pole configuration for the given specifications?
In MIMO systems the state feedback gains that achieve a given
pole configuration is not unique. What is the best K for a given
pole configuration?
The eigenvalues of the observer should be chosen faster than
those of the closed-loop system. Is there any other criterion
available to help decide one configuration over another?
The methods that we will now introduce give answers to these
questions. We will see how the state feedback and observer gains
can be found in an optimal way.

The Basic Optimal Control Problem

What does optimal mean? Optimal means doing a job in the best
possible way. However, before starting a search for an optimal
solution,

the job must be defined

a mathematical scale must be established to quantify what

we mean by best

the possible alternatives must be spelled out.

Unless these qualifiers are clear and consistent, a claim that a
system is optimal is really meaningless.
A crude, inaccurate system might be considered optimal because
it is inexpensive, is easy to fabricate, and gives adequate
performance.
Conversely, a very precise and elegant system could be rejected
as non-optimal because it is too expensive or it is too heavy or
would take too long to develop.

Further, the control as well as the state variables are generally

subject to constraints, which make many problems in optimal
control non-classical, since problems with path constraints can
hardly be handled in the classical calculus of variations.
That is, the problem of optimal control can then be stated as:
Determine the control signals that will cause a system to satisfy
the physical constraints and, at the same time, minimize (or
maximize) some performance criterion..

Despite its successes, however, optimal control theory is by no

means complete, especially when it comes to the question of
whether an optimal control exists for a given problems.

STATIC VERSUS DYNAMIC OPTIMIZATION?

Optimization is the process of minimizing or
maximizing the costs/benefits of some action.
Static Optimization refers to the process of
minimizing or maximizing the costs/benefits of some objective
function for one instant in time only.
Dynamic Optimization refers to the process of
minimizing or maximizing the costs/benefits of some objective
function over a period of time. Sometimes called optimal control.
In both cases equality and inequality constraints can be enforced.

1. Static Optimization is concerned with controlling a plant

under steady state conditions, i.e., the system variables are not
changing with respect to time. The plant is then described by
algebraic equations. Techniques used are ordinary calculus,
Lagrange multipliers, linear and nonlinear programming.

2. Dynamic Optimization concerns with the optimal control of

plants under dynamic conditions, i.e., the system variables are
changing with respect to time and thus the time is involved in
system description. Then the plant is described by differential (or
difference) equations. Techniques used are search techniques,
dynamic programming and variational calculus.

Static vs. Dynamic Optimization

Static Optimization

Optimal state, x*, and control

,u*, are fixed, i.e., they do not
change over time (for one

instant in time only)

Dynamic Optimization

Optimal state and control vary

over time (over a period of

time). Some times called

optimal control.

J* = J(x*, u*)
Functional minimization (or
maximization)
Parameter optimization

J* = J[x*(t), u*(t)]
Optimal trajectory
Optimal feedback strategy

Equality and inequality

constraints can be enforced.

Equality and inequality

constraints can be enforced.

Optimized cost function, J*,

is a scalar, real number

Optimized cost function, J*,

is a scalar, real number

1.1 Optimization Without Constraint

Problem: A scalar performance index L(u) is given that is a
function of a control or decision vector u We want to:
Min L(u)
Solution: L(u + du ) = L(u) + dL (increment)
-Taylor expansion to 2nd order in du
with gradient vector Lu and curvature matrix Luu
-Stationary point u*?
dL = 0 to1st order term for arbitrary du

Taylor series
is a representation of a
function as an infinite
sum of terms calculated
from the values of its
derivatives at a single
point.

The exponential function (in blue), and the sum of the

first n+1 terms of its Taylor series at 0 (in red).

Taylor series

Plot of f(x) = sin(2x) from /4 to 5/4; note fs second derivative

is f(x) = 4sin(2x). Tangent is blue where curve is concave up
(above its own tangent), green where concave down (below its
tangent), and red at inflection points: 0, /2 and

dx
dt

 Graphical Interpretation of solution

First-order necessary condition for optimality

A point x* satisfying this condition is called a stationary point. The
condition is first-order because it is derived using the first-order
expansion
Second-order conditions for optimality
We now derive another necessary condition and also a
sufficient condition for optimality, under the stronger hypothesis
that f is a C2 function (twice continuously differentiable).

Note: Like the previous first-order necessary condition, this

second-order condition only applies to the unconstrained case.

 A necessary condition for x* to be a maximum is f (x*) = 0, called

such an x* a critical point.
A sufficient condition for x* to be a maximum is that it is a critical
point and f (x*) < 0 .
That is, if x* satisfies both equations, then it must be a maximum.
The condition is also called a second-order condition.

The height of the function at "a" is

greater than (or equal to) the height
anywhere else in that interval.
Or, more briefly: f(a) f(x) for all x in
the interval In other words, there is
no height greater than f(a).
Note: f(a) should be inside the
interval, not at one end or the other.

Comments:
f may have many critical points. The sufficient conditions helps
us pick the maximum out of the set of all critical points.
Both conditions only establish that x* is a local max. That is,
there is no point "close to" x* that has a higher f-value. But there
may be other points not "close to" x* that do. In that case there
must be multiple points fulfilling both conditions and we simply
calculate which of these has the highest f-value.
If f is strictly concave, there is only one critical point and this
must be the global max. The reason is that f (x) < 0 everywhere
(that is what it means to be strictly concave). Thus, there can only
be one point where f (x) = 0 and that satisfies the sufficient
condition automatically.

Examples: The function

http://en.wikipedia.org/wiki/Maxima_and_minima

x2 has a unique global minimum at x = 0.

x3 has no global minima or maxima. Although the first
derivative (3x2) is 0 at x = 0, this is an inflection point.
xx has a unique global maximum at x = e figure.
x-x has a unique global maximum over the positive real
numbers at x = 1/e.
x3/3 x has first derivative x2 1 and second derivative 2x.
Setting the first derivative to 0 and solving for x gives
stationary
points at 1 and +1. From the sign of the second
derivative we
can see that 1 is a local maximum and +1 is
a local minimum. Note that this function has no global
maximum or minimum.

Cont. examples: The Function

|x| has a global minimum at x=0 that cannot be found by
taking derivatives, because the derivative does not exist at x=0.
cos(x) has infinitely many global maxima at 0, 2, 4, ,
and infinitely many global minima at , 3, .
2 cos(x)x has infinitely many local maxima and minima, but
no global maximum or minimum.
cos(3x)/x with 0.1x1.1 has a global maximum at
x=0.1 (a boundary), a global minimum near x=0.3, a local
maximum near x = 0.6, and a local minimum near x=1.0.

x3+3x22x+1 defined over the closed interval (segment)

[4,2] has a local maximum at x = 153, a local

minimum at x = 1+53, a global maximum at x = 2 and a global
minimum
at x = 4.

Contour
When the lines
close
together
magnitude
of
gradient is large:
variation is steep.

are
the
the
the

Contour

The bottom part of the diagram shows some contour lines with a
straight line running through the location of the maximum value.
The curve at the top represents the values along that straight line

singular point:
the singular points are those at which all the partial derivatives
simultaneously vanish.
A point of an algebraic variety which is not singular is said to be
regular. An algebraic variety which has no singular point is said to
be non singular or smooth. y2=x2(x + 1)
A plot of this curve is below with the singular point at the origin.
An example of singular point is when a graph crosses over itself:
The reason for that algebraically is that
both sides of the equation show powers
higher than 1 of the variables x and y. In
terms of differential calculus, if
F(x,y) = y2 x2(x + 1), so that the curve
has equation F(x,y) = 0,
then the partial derivatives of F with
respect to both x and y vanish at (0,0).

The linear quadratic regulator (LQR) is a special case of optimal control problem that
can be analytically solved, which arises when the objective (cost) function is a quadratic
function of state x and input and u, and a LTI state equation is considered.
The quadratic performance index to minimize is:

thus, the loss function and the terminal cost function can be expressed as quadratic
forms:

It is in general assumed that the weighting matrices Q, R and S are symmetric, with Q
and S positive semi definite and R positive definite.
Qualitative idea: if the quadratic term involving the state are large compared to those of
the input (Q>>R), then the optimal solution will exhibit good regulation performance,
that is, the response will go fast to its desired value, but at the expense of large
control energy.
A square matrix M is called symmetric if it verifies: M=MT
A real squared matrix M is called positive definite if xTMx>0, for every vector x0
A real squared matrix M is called positive definite if xTMx0, for every vector x0

Minimize an Absolute Criterion

Achieve a specific objective
Minimum time
Minimum fuel
Minimum financial cost
to achieve a goal

For a continuous-time linear system, defined on

with a quadratic cost function defined as

, described by

Example

-Slope becomes gradient vector:

- 2nd order derivative becomes Hessian

(or curvature matrix):

 Further Extension: vector function of vector

- Slope of f(u) with respect to u is a Jocabian of nm matrix

- 2nd order derivative becomes 3rd order tensor

(Appendix A.3)

-Meaning of matrix Luu>0 or <0?

Q positive definite (Q>0) if xTQx < 0 x 0
Q negative definite (Q<0) if
Q positive semidefinite (Q 0) if xTQx 0 x 0
Q indefinite if xTQx > 0 for some x 0 and xTQx < 0 for some x 0
Q negative semidefinite (Q 0) if xTQx 0 x 0

-Testing Positive/Negative Definiteness

Method -1. By Eigenvalues i of Q
Q positive definite (Q>0) if all > 0, i =1,L,n I
Q negative definite (Q<0) if ( 0) if all 0, i 1, , n Q positive
semidefinite Q i = L
0 and some 0, i 1, ,n i L = < > i Q indefinite if some
Q negative semidefinite xTQx > 0 x 0

Method -2. By Determinant

* Testing positive/negative-definiteness
Obtain leading minors of Q:

 Testing semi-positive/negative-definiteness : Obtain Principal

minor of Q:

* Q is indefinite if none of the above cases are satisfied

 Some useful Matrix Calculus formula (Appendix A.4)

 If Q is symmetric

 Some useful Hessians

If Q is symmetric

- Is u* max/min point of L(u) ?

dL at stationary point u* becomes

 Optimized value of performance index:

* Stationary point

Plug in the above equation

Example 1.1-2: Optimization by Scalar Manipulations

We have discussed optimization in terms of vectors and the gradient. As an alternative
approach, we could deal entirely in terms of scalar quantities.

Find L(1,-1) ?

Before start
section 1.2
How to draw
mesh, contours
plots
X and Y arrays for 3-D plots

X and Y arrays for 3-D plots

transforms the domain specified by vectors x and y into arrays X

and Y that can be used for the evaluation of functions of two
variables and 3-D surface plots
f (x,y) = x2 + y2. Then an easy calculation shows that f (2, 1) = 5.
Therefore, (2, 1, 5) is an ordered triple that lies on the graph of f .
This point is shown in Figure.
Figure Point (2, 1, 5) is on the graph of f .
To plot the graph of f (x) = x2 in the plane, we begin by making a
table of points that satisfy the equation, as shown in Table

Table A table of points satisfying f (x) = x2.

For creating a table of points that satisfy the equation f (x,y) = x 2 + y2. Matlab
accomplishes this with the meshgrid command.
>> [X,Y]=meshgrid([1,2,3,4,5])

Output is easily understood if one superimposes the matrix Y onto the matrix X to obtain
a grid of ordered pairs.
The rows of the output array X are copies of the vector x and the columns of the output
array Y are copies of the vector y.

Therefore, Table contains a set of points in the plane that we will substitute into the
function f (x,y) = x2 +y2. Matlabs array smart operators make this an easy proposition.

>> [X,Y]=meshgrid([1,2,3,4,5])
>> Z=X.2+Y.2

It is now an easy task to plot the surface to which these points belong.
>> mesh(X,Y,Z)

Figure : Plotting the surface

f (x,y) =x2 + y2.

>> x=-3:3
x = -3 -2 -1 0 1 2 3
>> y=-3:3
y = -3 -2 -1 0 1 2 3
>> [X,Y]=meshgrid(x,y)
>> Z=X.2+Y.2
>> mesh(X,Y,Z)
18

18

18

18

The surface z = x2 + y2 on the domain

The plot of z = x 2 + y 2.

20
15
z-axis

>> x=-3:.2:3;
>> y=-3:.2:3;
>> [X,Y]=meshgrid(x,y);
>> Z=X.^2+Y.^2;
>> mesh(X,Y,Z)
>> xlabel(x-axis)
>> ylabel(y-axis)
>> zlabel(z-axis)
>> title(The plot of z = x^2 + y^2.)

10
5
0
4
2
0
-2
y-axis

-4

-4

-2

0
x-axis

close all
y=-1:.1:1;
z=-1:.1:1;
[yy,zz]=meshgrid(y,z);
xx=zz.^3;
mesh(xx,yy,zz);
xlabel('x-axis')
ylabel('y-axis')
zlabel('z-axis')
title('The plot of X = Z^2
view([130,40])

Before start section 1.2

Solve section:1.1 Examples
and
End Chapter Problems of
section 1.1

Example (Not in book)

syms x y
f=x^3+y^3+3*x^2-3*y^2-8
fx=diff(f,x);
fy=diff(f,y);
s=solve(fx,fy);
[s.x,s.y]

Max & Min?

Classification of critical points

Luu=0 ID
Luu<0 +D
Luu>0 -D

Value of f(x,y) at critical points

f(0,0)=-8

f(-2,0)=(-8)+0+12+0-8=-4

f(-2,2)= -8+8+12-12-8=-8

f(0,2)= 0+8+0-12-8=12

Mesh Surface
50

[x,y]=meshgrid(-3:.2:3);
z=x.^3+y.^3+3*x.^2-3*y.^2-8;
mesh(x,y,z)
Difficult to
xlabel(x-axis)
view Mesh
ylabel(y-axis)
Surface

-50

-100
4
2
0
-2
-4

-2

-4

Contour map
[x,y]=meshgrid(-3:.2:3);
z=x.^3+y.^3+3*x.^2-3*y.^2-8;
contour(x,y,z)
xlabel( x-axis )
ylabel( y-axis )

-1

-2

-3
-3

-2

-1

[x,y]=meshgrid(-3:.2:3);
z=x.^3+y.^3+3*x.^2-3*y.^2-8;
contour(x,y,z,20)
xlabel( x-axis )
ylabel( y-axis )

-1

-2

Slop ??

-3
-3

-2

-1

1
y-axis

[x,y]=meshgrid(-3:.2:3);
z=x.^3+y.^3+3*x.^2-3*y.^2-8;
[c,h]=contour(x,y,z,-14:-4);
clabel(c,h)
xlabel( x-axis )
ylabel( y-axis )

-1

-2

-3
-3

-2

-1

0
x-axis

Approach from any

direction you are
walking downhill
crossing contour at
-8,-9,-11 and finally local
minimum at (0,2)

y-axis

If you move from point (0,0)

or (-2,2) in any one direction
the height of the contour
increases and decreases if
you move in opposite
direction

-1

Approach from any

direction you are walking
uphill crossing contour at
-7,-6,-5 and finally local
maximum at (-2,0)

-2

-3
-3

-2

-1

0
x-axis

Question:2 Assignment
Find the minimum value of L(x1,x2)=x12-x1x2+x22+3x1;

syms x y
f=x.^2-x*y+y.^2+3*x;
fx=diff(f,x);
fy=diff(f,y);
s=solve(fx,fy);
[s.x,s.y]
[x,y]=meshgrid(-3:.2:3);
f=x.^2-x*y+y.^2+3*x;
mesh(x,y,f)
xlabel(x-axis)
ylabel(y-axis)

let = L(x,y)=x2-xy+y2+3x

[x,y]=meshgrid(-3:.2:3);
f=x.^2-x*y+y.^2+3*x;
contour(x,y,f)
xlabel( x-axis )
ylabel( y-axis )
[x,y]=meshgrid(-3:.2:3);
f=x.^2-x*y+y.^2+3*x;;
contour(x,y,f,50)
xlabel( x-axis )
ylabel( y-axis )
[x,y]=meshgrid(-3:.2:3);
f=x.^2-x*y+y.^2+3*x;
[c,h]=contour(x,y,f,-14:-4);
clabel(c,h)
xlabel( x-axis )
ylabel( y-axis )

[x,y]=meshgrid(-3:.2:3);
f=x.^2-x*y+y.^2+3*x;;
contour(x,y,f,50)
xlabel( x-axis )
ylabel( y-axis )

ANIMATION _ JAVA
[x,y]=meshgrid(-3:.2:3);
z=x.^3+y.^3+3*x.^2-3*y.^2-8;
u = z(x,y);
quiver(x,y,u) ; grid on; axis equal
Check for quiver(x,y,u)

Problem 1.2-4 (Book)

Shortest
Distance
Parabola and a Point.

Between

A meteor is moving along the path

syms x y l
f=x^2+3*x-6-y;
2
Y=x
+3x-6. A space staion is at the point
J=(x-2)^2+(y-2)^2;
(x,y)=(2,2).
H=J+l*f;
H_x=simplify(diff(H,x));
a. Us Lagrange multiplier to find a cubic
H_y=simplify(diff(H,y));
equation for x at the closest point of
l=solve(H_y,'l');
approach.
y=solve(f,'y');
xeq=simplify(subs(subs(H_x,'l',l),'y',y));
b. Find the closest point of approach
xstar=roots(sym2poly(xeq));
(x,y), and the distance from this point
[dummy,index]=max(xstar>0);
to (2,2).
xx=xstar(index);
yy=subs(y,'x',xx);
Jstar=subs(subs(J,'x',xx),'y',yy);disp('** Final Answer');
x1=(-5:.1:5);y1=x1.^2+3*x1-6;
%[x2,y2]=circle(2,2,sqrt(Jstar));
figure(1);plot(x1,y1,2,2,'o',xx,yy,'x',xx,yy,':');
axis([-5,5,-10,10]);axis('square');
xlabel('x');ylabel('y');
legend('space station','closest approach','meteorite path');