Robust Optimal Service Analysis of Single-Server Re-Entrant Queues.
We generalize the analysis of J.A. Ball, M.V. Day, and P. Kachroo (Mathematics of Control, Signals, and Systems, vol. 12, pp. 307345, 1999) to a fluid model of a single server re-entrant queue. The approach is to solve the Hamilton-Jacobi-Isaacs equation associated with optimal robust control of the system. The method of staged characteristics is generalized from Ball et al. (1999) to construct the solution explicitly. Formulas are developed allowing explicit calculations for the Skorokhod problem involved in the system equations. Such formulas are particularly important for numerical verification of conditions on the boundary of the nonnegative orthant. The optimal control (server) strategy is shown to be of linear-index type. Dai-type stability properties are discussed. A modification of the model in which new customers are allowed only at a specified entry queue is considered in 2 dimensions. The same optimal strategy is found in that case as well.

Figure

1. Re-entrant Server
There is much current interest in developing optimal service strategies for queueing systems. The volume
by Kelly and Williams [9]includes several articles addressing this. Although queueing models are generally
integer-valued and stochastic, Dai [4]and others have developed connections between the stability of
Date: November 30, 2001.
Research supported by the ASPIRES program of Research and Graduate Studies, and the Millennium program of the College
of Artsand Sciences, both of Virginia Tech.
stochastic queuing systems and their deterministic fluid limits. Thus optimal strategies for fluid models
are recognized as significant for stochastic models. Fluid models for a large class of queueing systems can
be described by equations of the general form (5) introduced in Section 2 below. We pursue the same robust
control approach as in [2]for such models. Much of what we present here is a further development
of ideas from that paper. In particular, Section 2.1 gives explicit representations of the velocity projection
map (x, v) of the Skorokhod reflection mechanism which comes into play when one or more queues are
empty. Section 3 develops the construction of the value function for our control problem. Here, as in [2],
the construction proceeds without regard to the Skorokhod dynamics on the boundary of the nonnegative
orthant. (In more general multiple server examples the Skorokhod dynamics will play a more decisive role.)
Even so, the solution we construct must be shown to satisfy various inequalities associated with optimality
with respect to the Skorokhod dynamics on the boundary. We do not provide a deductive proof of these
inequalities, but rely instead on a system of numerical confirmation for individual test cases in Section 4. The
explicit representations of (x, v) are important for this, and for the optimality argument of Section 5. The
version of that argument given here improves on the one in [2]in that it applies to all admissible strategies,
not just those of state feedback form.
Our model allows new arrivals and unserved departures in the form of an exogenous load qi(t)foreach
xi; see (1) below. In some queueing applications this feature would be inappropriate. For instance in typical
re-entrant lines, new arrivals only occur at a specified entry queue and departures only as service is completed
at a designated final queue. In Section 6 we will look at the the 2-dimensional case of our model under the
more restrictive assumption that exogenous arrivals are only allowed in the entry queue x1. This requires a
number changes in our calculations. But we find that this change to the model does not eect the resulting
optimal service policy.
2. The Model and Approach of Optimal Control
We describe in this section the general model formulation and performance criteria that we will use. Fluid
models for a large class of queueing systems can be described by equations of the nominal form
(1) x(t)=q(t)  Gu(t).
The state variable is n-dimensional: x =(x1,.,xn)  Rn. For queueing models x(t) must remain in the
nonnegative orthant, K in (4) below. For that purpose we will couple (1) with Skorokhod problem dynamics,
resulting in (5) below. The term q(t)istheload on the system due to new arrivals (or unserved departures
if qi(t) < 0). The service allocation is specified by the control function u(t) whose values are taken from
a finite set U0  Rm of possible service control settings. For purposes of an adequate existence theory for
solutions to (5) we relax this to allow u(t) to be taken from the convex hull
(2) U =convU0.
For our single-server examples we willsimply take the standard unit vectors in Rn.Thus
1. The matrix G converts u(t) to the appropriate vector of
contributions to x. For the case to be considered here (Figure 1) G will be the lower triangular matrix
.
. sn1 0
The si > 0 are parameters which specify the service rates for the respective queues. Thus when u(t)=ek
(k<n), the eect of Gu(t) in (1) is to drain queue xk at rate sk with the served customers entering xk+1
at the same rate:
Multiple server examples are easily modeled by (1) and (5) as well. Consider for example the 2-server
re-entrant line in Figure 2. It would be natural to use
The first two columns in G correspond to the service allocation at server A and the second two at server
B. The u  U0 correspond to the 4 dierent combinations in which server A chooses between x1 or x3 and
server B chooses between x2 or x4.
x3

Figure

2. 2-Server re-entrant line
2.1. Skorokhod Problem Dynamics. We denote by K the nonnegative orthant of Rn:
The faces of K are
The interior normal to iK is the standard unit vector We will use
to denote the set of all coordinate indices. For x  K,
I(x)={i
will denote the set of indices with zero coordinate values.
An essential feature of queueing models is that x(t) remains in K for all t. One could simply impose
this as a constraint the control functions u(t) and loads q(t) which are considered admissible. Although
some constraints on the load are reasonable, we find it much more natural in general to couple (1) with the
dynamics of a Skorokhod problem. On each face iK we specify a constraint vector di. If the solution of
(1) attempts to exit K through iK, then the idea is to add some positive multiple of di to the right side of
(1) to prevent the exit. A precise formulation is the following: given x(0)  K and q(t),u(t) (which we will
assume to be locally integrable), let
Problem is to find a continuous function x(t)  K, a measurable function r(t)  Rn and a
nondecreasing function (t)  0 which satisfy the following for t  0:

(0,t]
for each t, r(t)= iI(x(t)) idi for some i  0;
(t)= (0,t] 1x(s)K d(s).
By imposing a normalization there will be a unique solution, provided K and di satisfy certain
conditions. Dupuis and Ishii [5]and Dupuis and Ramanan [7]provide a substantial body of theory of
Skorokhod problems in general. In particular they show, using a velocity projection map (x, v), that a
Skorokhod problem can be expressed as a dierential system. The velocity projection map is of the form

for an appropriate choice of i  0. (See (6) below.) The result of coupling our (1) with the appropriate
Skorokhod problem is expressed as
holding almost surely.
The appropriate constraint vectors di are determined by the structure of the system in Figure 1. If
x  iK and server i is active (ui > 0) but the applied service rate siui exceeds the inflow qi to xi, then
according to (1) x would exit K through iK. In an actual network the system could not really use the
full service capacity siui allocated to xi. Instead, service would take place at a lower level which exactly
balances the inflow and outflow of queue xi. Mathematically this is achieved by adding a positive multiple of
the column Gei of G to the right side of the system (1), bringing xi to 0 and producing the correct reduction
of the throughput xi  xi+1 to the next queue. So we take di = s1 Gei (the normalization being so that
1). The same prescription is appropriate for the example of Figure 2: take di to be the unique
column of G having a positive entry in row i, normalized so that ni  di =1.
At this point we wish to highlight the fact that no restrictions on u(t)andq(t) are needed to keep x(t)
in the nonnegative orthant; (5) will determine a state trajectory with x(t)  K regardless. Thus we always

instruct the server to work at full capacity ( and the Skorokhod dynamics can be viewed as
automatically reducing the service rates to the levels that can actually be implemented. The model allows
a separate load term qi(t) for each queue. For re-entrant queues, one typically would only want to allow
for the entry queue in each re-entrant sequence. In Figure 1 for instance it would be natural to assume
nice feature of single servers with respect to the L2 performance criteria of Section 2.4
is that the the optimal strategy is the same for all loads, regardless of which coordinates might be zero. In
queueing applications one also naturally assumes that qi  0. However it was argued in [2]that, for purposes
of vehicular trac for instance, it is reasonable to consider qi < 0. This would correspond to customers that
leave the system without waiting to receive service all the way through. This is a reasonable consideration
in some applications. However it is hard to conceive of a realistic interpretation for qi < 0 when xi =0.
Even so, (5) will still yield a mathematical solution. The eect of the additional +idi terms in (x, q Gu)
might then be seen not as reductions to the service rates but as a transference of the reducing influence of
qi < 0 from the empty xi to the queues xj further along in sequence; fluid at xj would be drawn backwards
through the system to satisfy to external demand due to qi < 0atpreviousxi.
With di defined, we face the important technical issue of existence and regularity properties of the Skorokhod
problem. This issue is treated in detail in [5]and [7.] Those treatments consider a more general
convex polyhedron in place of our K. Our particular choice of the nonnegative orthant falls within the scope
of the earlier work [12]. Let D =[di]be the matrix with the constraint vectors as columns. In our case,
I  Q where Q is the subdiagonal matrix with entries

Assuming that Q has nonnegative entries and spectral radius less than 1, both clearly satisfied for us, [12]
provided a direct construction of the solution of the Skorokhod problem. In [7]it is shown that these
conditions from [12]fall within the scope of a more general set of sucient conditions for existence and
Lipschitz continuity of the Skorokhod map y()  x().

Drawing on the ideas of [12], we can give a direct construction of the (x, v) appearing in the dierential
formulation (5) of the Skorokhod problem. For any x  K and v  Rn, we will show that
be characterized using a linear complementarity problem:

subject to the following constraints for each i  I(x):
For i/I(x) there is no constraint on wi, and we consider to be implicit. Let  =[i]n .Using
I  Q, and rewriting (6) as
it is easy to see that the complementarity problem is equivalent to saying  is a fixed point  =x()of
the defined coordinate-wise by

The notation y+ refers to the usual positive part: fixed point representation is a
particular case of the general fixed point representation of variational inequalities in Chapter 1 of [10].) We
first observe the existence of a unique fixed point. The argument of [12]is to observe that (after a linear
change of variables)  is a contraction, under the nonnegativity and spectral radius assumption mentioned
above. However this is even simpler for our particular Q;  =x() reduces to

which determine the i sequentially. Iteration of x from any initial  will converge to the fixed point after
at most n steps. This makes it particularly simple to see that v   and thus v  (x, v) are continuous,

and to evaluate (x, v) numerically. Indeed (x, v) is Lipschitz in v for a fixed x, and is jointly continuous
in (x, v)ifx is restricted so that I(x) is constant.
We can easily check that (x, v) as defined by the above complementarity problem is indeed the velocity
projection map as identified in [5]. First, following an observation of [3], we can check that (y)=(0,y)is
the discrete projection map of Assumption 3.1 of [5]. Indeed, if y  K then clearly
the complementarity problem: (0,y)=y.Fory/K,

y  (0,y)=y
which is of the form  for some   0and  d(w), where d(x) is the set of reflection directions as defined
in [5, 3]. Next, suppose x  K, v  Rn, and let w,solve the complementarity problem for
above. We claim that
suciently small.
This will imply that
which is the characterization of (x, v)in[5,5.3]. Since w,solve (6), it follows that
di.We want to see that satify the complementarity conditions (7) -
associated with First consider i  I(x). Since xi =0,wehave
For the prodiuct (9) we have
Next consider i/I(x). Provided h>0 is suciently small we have
which also confirms the product condition. This verifies (11), as desired.
For purposes of our calculations, the characterizations of (x, v) in Lemma 1 below will be useful. If
J  N we will use
NJ =[nj]jJ and DJ =[dj]jJ
to denote the matrices whose columns are the normal vectors nj and constraint directions dj for the j  J.
Given v, and the corresponding i as described above, let
From the complementarity problem we know F0  L  I(x). Using any F0  F  L the values of i, i  F
are determined by setting wi =0,i  F in (6). In other words we can solve for F =[i]iF directly in
and consequently
where RF is the reflection matrix
For simply take
More precisely,
The fact that wi  0fori/F is equivalent to
Suppose that we dont know F0 or L at the outset, but just take an arbitrary F  I(x), calculate
v. By construction wi =0fori  F and i =0fori/F. So item 3 of the
complementarity problem is satisfied. If item 1 is satisfied of i  F, which is to say
and item 2 holds for i  I(x) \ F, which is to say
then we can say that is in fact (x, v). This discussion proves the following lemma.
Lemma 1. Given x  K, v  Rn,andF  I(x), the following are equivalent:
1. (x, v)=RF v;
2. both of the following hold:
(a) BF v  0 when

3. for some L with F  L  I(x) all of the following hold:
(a) BF v>0 when
(b) NL\F RF v =0when

(c) NIT(x)\LRF v>0 when

Notice that the strict inequality in 3 (a) simply identifies F as
2.2. The Optimal Control Policy. Our goal is to design a feedback control strategy (x), prescribing
a value in the extended control set U for each x  K, so that using u(t)=(x(t)) produces optimal
performance of the system. The criteria used to determine optimality is based on
where >0 is a parameter. Roughly speaking, the control should keep the integrated cost (14) small,
so that x(t) remains small compared to the load q(t) in a time-averaged sense. We will give this a more
precise formulation in Section 2.4 below. The running cost 1 x2  2 q2 of (14) has it roots in classical
H control, and is attractive for its broad familiarity and success in a wide range of control applications.
Other choices might be more appropriate for particular queueing applications, such as those associated with
optimal draining and time-to-empty criteria; see [18]and [1.] There are however considerations that favor
L2 in the trac setting. For a given total customer population, the L2 norm favors balanced queue lengths
over a situation in which some queues are empty and others are full. When each customer is a person who
has to wait in a queue, a cost structure that can be minimized by using excessive waits for some small class
of customers would be considered unacceptable.
The optimal policy itself is easy to identify by naive considerations at this point. In order to minimize
for a given q(t) one would minimize 1 x(t)2, for which one would naturally try to choose u(t)tominimize d 1 In the interior of K, where (1) applies, this suggests that the optimal u
is that which maximizes x  Gu over u U.OnK x is given by (5), which makes finding the u to minimize
potentially more dicult. However if we assume that all qi  0 then the Skorokhod dynamics
do not aect the minimizing u. To see this consider x  K with suppose u Umaximizes xGu.

It is easy to see that sup x  Gu > 0. Observe that for any i  I(x), since o-diagonal entries
of G are nonpositive, x  Gei  0. Thus the i-coordinate of u must be 0, from which we can conclude that
ni  Gu  0. Since qi  0 by hypothesis, we see that
ni  (q  Gu)  0, for all i  I(x).
This means that (x, q  Gu)=q  Gu. Also, for any i  I(x)wehavex  di  0 because
the i coordinate of di is positive. Therefore, for any u Uwe can say that

Thus the policy
U
is an obvious candidate for the optimal service policy. We will see below that it is indeed optimal in the
sense to be made precise in Section 2.4.
Several comments should be made at this point. First observe that (x) is set-valued. There are inherent
discontinuities in (x), as the optimum u jumps among the extreme points of U. When we replace u(t)
by (x(t)) in (5), we will want the resulting feedback system to have good existence properties. This is
addressed using the Filippov theory of dierential inclusions. For that it is important that (x) have closed
graph and be convex set-valued. (The notion of closed graph is often called lower semi-continuity for
set-valued functions.) It is easy to see that (15) has these properties.
Another important point to make is that the naive reasoning which suggests (x(t)) to us does not
actually imply that it achieves the smallest possible value of 1 x(T)2 for a given q() and target time T.Itis conceivable that it might be better to forgo pointwise minimization of d 1 x(t)2 in order to drive x(t)
into a dierent region (or section of the boundary) where larger reductions of x(t) could be achieved. Some
sort of dynamic programming argument, such as the Hamilton-Jacobi equation developed in Section 3, is
needed to adequately address such global optimality issues.
Although it may not have much practical import, one might ask whether allowing qi < 0 when xi =0
might aect the choice of u Uwhich minimizes x(x, q Gu). Indeed it can. If qi < 0 is large enough its
eect through the Skorokhod dynamics can produce cases in which a u/(x) minimizes x  (x, q  Gu).
This consideration would lead to an enhanced optimal policy which agrees with (x) on the interior of K,
but depends on both q and x when x  K. Although no longer state-feedback, this enhanced control would
(we expect) produce lower values of (14), but only for those negative loads q(t) which, as we described above,
have the eect of drawing customers backwards through the system. Even so, this enhanced control would
not improve the performance of the system in the worst case sense of the dierential game formulated in
Section 2.4, as Theorem 1 below will assert.
2.3. Minimum Performance Criteria. Our strategy  is expressed in state feedback form. Given a load
q(), the associated control function u(t) would be what results from solving the system
This system is a combination of a dierential inclusion, in the sense of Filippov, and a Skorokhod problem
as described above. The discussion in [2, Section 1.4]outlined how the arguments of [6]can be adapted to
establish the existence of a solution. A proof of uniqueness is more elusive. The usual Filippov uniqueness
condition would be that for some L
(xa  xb)  [(q  Gua)  (q  Gub)]= ( xa  xb)  (Gub  Gua)  Lxa  xb2.
This is immediate (using since by definition of (),
xa  Gub  xb  Gua and xb  Gua  xb  Gub
for all ua  (xa), ub  (xb). However, as noted in [2], when coupled with Skorokhod dynamics (16)
we are unable to conclude uniqueness based on existing results in the literature. Until that issue can be
addressed, we must allow the possibility of multiple solutions to (16). The uniqueness question is not
essential to our main result Theorem 1, however. We simply need to formulate its statement in such a way
that strategies are allowed to produce more than one control function u(t) for a given load q(t).
In general a service strategy () maps a pair x(0),q() to one or more control functions u(). We will
write u(t)=[x(0),q()](t), although this notation is not quite proper if there are actually more than one
u() associated with x(0),q()by. Rather than formulating a cumbersome notation to accommodate this,
we will simply use phrases like for any u(t)=[x(0),q()](t) to refer to all possible u(t). A strategy
should produce one or more control functions for any x(0)  K and load function q() which is locally
square-integrable. We insist that a strategy be nonanticipating, in the sense that if q(s)=q(s) for all s  t,
then for any u(t)=[x(0),q()](s) there is a u(t)=[x(0), q()](s) with u(s)=u(s) for all s  t.Given
any such x(0), q() and a resulting u(t), the general existence and uniqueness properties of the Skorokhod
problem (e.g. [5]) provide a unique state trajectory x(t)  K.
We will call a strategy () non-idling if for any nonnegative load qi(t)  0 for all i and all t  0, any
any u(t)=[x(0),q()](t), the resulting state trajectory x(t) has the property that ui(t) > 0
and occur simultaneously for some i only if In other words, all service eort is allocated
to nonempty queues, unless all queues are empty. In particular, our strategy (15) is non-idling, because if
x  K and is the index of the largest nonzero coordinate of x.
One of the features of single servers as in Figure 1 is that for nonnegative loads, a non-idling strategy
will never invoke the Skorokhod dynamics on K, until it reaches Indeed if x(t)  K \{0} but
the non-idling property means that 0, from which the structure of G implies that
ni  Gu(t)  0.
Since qi(t)  0, we conclude that
ni  (q(t)  Gu(t))  0.
Thus unless Multiple servers do not have this property. In the
case of Figure 2 for instance, if both then the service eort at B is wasted and the Skorokohd
dynamics will definitely come into play, regardless of x1 and x3. The Skorokhod dynamics will thus have a
stronger influence on the design of optimal strategies for multiple server models.
When considering those fluid models that arise as limits of discrete/stochastic queueing systems, the
stability criterion of Dai [4]is important for purposes of positive recurrence of the stochastic model. In that
setting the load q(t) is typically constant, with 1/qi equal to the mean time between new arrivals in queue
xi. The stability property of [4]is simply that for any x(0)  K, the state x(t) reaches at a finite
T  0. For our single server model, all non-idling strategies are equivalent in this respect. To see why,
consider the vector
Observe that T G =(1, 1,. ,1), so that for all u Uwe have   1. For any nonnegative load q(t)
and the u(t) resulting form any non-idling strategy, we have (on any interval prior to the first time T when
d
x(t)=  (x(t),q(t)  Gu(t))
dt
Said another way, W(x)=  x is a sort of universal Lyapunov function for all non-idling controls. Thus,
the first time T for which does not depend on the choice of non-idling control; it only depends on
the load q(t). For constant nonnegative loads q(t)  q, the Dai stability property simply boils down to
Moreover if q =(q1, 0,. ,0)T then this reduces to the familiar load condition of [4, (1.9)]:
q1 < 1.
Figure

3 illustrates this stability property. We have taken the optimal strategy (x) for our model with
subjected the system to the constant disturbance q(t)  (.4, 0). For these
parameters we find  =(2, 1) and q  so that the load condition (18) is indeed satisfied. The figure
illustrates the resulting trajectories of (16). When x(t) reaches the ray from the origin in the direction of ,
the solution of (16) in the Fillipov sense uses the averaged control value
which takes x(t) to the origin in finite time directly along the  ray. One may check that this ray consists
of those x for which (x)=U is multiple-valued.
Theorem 1 below considers the optimality of  with respect to all strategies  that satisfy the following
minimum performance criterion: given x(0)  K with
x(0) < 1,
there exists <1 so that whenever q() is a nonnegative load satisfying
q(t)  1 for all t,
and any u(t)=[x(0),q()](t) the resulting state trajectory satisfies
x(t) <for all t  0.
It is clear from our discussion above that every non-idling control satisfies the minimum performance criterion;
simply take
2.4. The Robust Control Problem. We now want to define more carefully the sense in which our service
strategy (x) is optimal. We follow the general approach of Soravia [14]to formulate a dierential game
based on (14). The focus is on a value function of the form

dt.
Here x  K, the outer infimum is over strategies (), the inner supremum is over locally square integrable
loads q(), all u(t)=[x, q()](t) and bounded time intervals [0,T], with x(t) the resulting solution of (1) for
x1

Figure

3. Controlled Trajectories for q(t)  (.4, 0).
The gain parameter >0 is customary in robust control formulations. However for the structure of our
problems  scales out of the game (19) in a natural way. To see this, consider a particular load q(t), control
u(t), and solution x(t) of (5). Make the change of variables
Then x(t)= d x(s), and because K is a cone, (x, q  Gu)=(x, q Gu). Thus x(s) solves (5) on the new
ds
time scale. With

ds.
If V ()=V1() is the value (19) for then the above implies that
V(x)=3V (1x).
From this point forward we simply take instead of V.
We can only expect V (x) <  to hold in a bounded region. To see why, imagine a load q(t) which is large
on some initial interval 0  t  s so as to drive the state out to a large value X, and then q(t)ischosenfor
t>sso as to maintain x(t)=X for t>s: q(t)=Gu(t). If X > supuU Gu, the integral in (19) grows
without bound as T , producing infinite value. We must exclude such scenarios from the definition of
the game. It turns out that the region  in which V (x) will be finite is described using the vector  of (17)
above:
We restrict the T in (19) to those for which x(t) remains in  for all 0  t  T.
This qualification on the state in turn requires us to place some limitations on the strategies () considered
as well. We need to exclude controls that cheat by encouraging x(t) to run quickly to the outer boundary
of  to force an early truncation of the integral in (19). Such controls could achieve an artificially low value
by having actually destabilized the system. To exclude such policies we insist that all control strategies ()
satisfy the minimum performance criterion stated at the end of Section 2.2. With these qualifications, we
can now state precisely the optimality properties of the feedback strategy (x).
Theorem 1. Let  and (x) be as defined above and suppose the boundary verifications of Section 4 have
been successfully completed. Using the control (x),forx  ,define

where the supremum is over all loads q(t), all resulting control functions u(t)=[x, q()](t), and those
<T< such that the controlled state from
other control strategy () satisfying the minimum performance criterion,

with the same qualifications on the supremum.
The proof of these assertions will be discussed in Section 5 below. The qualification regarding the boundary
verifications of Section 4 will be explained in the last paragraph before Section 3.1.
3. Construction of the Value Function by Staged Characteristics
The proof in Section 5 of Theorem 1 is based on showing that the function V (x) of (21) solves the
Hamilton-Jocobi-Issacs equation associated with the game (19):
The Hamiltonian function is complicated by the special reflection eects on K:

The essential property of our strategy (x) for the proof is that, given x
point for the supq infU defining H(x, DV (x)) in (24) is given by any u  (x). To be
specific, the minimum value of
over u Uis 0, achieved at and the maximum value of
over q  Rn is 0, achieved at q. Together these imply (23). Our primary task is to produce V (x)and
establish this property of .
In general (23) must be considered in the viscosity sense. Lions [13]considers the viscosity-sense formulation
of a general class of problems involving Skorokhod dynamics on K. Instead of working with H as in
(24), the viscosity sense solutions are described using only the interior form of the Hamiltonian (27), together
with special viscosity sense boundary conditions on K. In our case it will turn out that the solution V is
actually a classical one. We find the direct formulation in terms of H more natural for our development.
We will construct the desired solution V (x) by working in the interior K, where the complicating eects
of (x, v) are not present: (x, v)=v so

=infHu(x, p).
uU
Here Hu refers to the the individual Hamiltonian for u U:
The supremum is achieved for Also observe that for u Uto achieve the infimum in (27) means
simply that u maximizes p  Gu.Soforx  K, (23) and the saddle point conditions (25) and (26) simply
reduce to the statement that for any u  (x),
uU
We turn now to the construction of V (x) of by a generalized method of characteristics. We cover  with
a family of paths x(t) as described below. The idea is that at a point the gradient of DV (x(t))
should be given by the costate trajectory p(t) that accompanies x(t). Thus a simple covering of  by a
family of such paths will determine the values of DV (x) in . Knowing V
in the region. We itemize the essential features of this family of x(t),p(t) in (30)-(33), and then explain
their relation to the Hamilton-Jacobi-Isaacs equation and saddle point property above. To begin, the paths
must solve the system of ODEs
for some piecewise constant u(t) U. The value of u(t) may change from one time interval to another,
but at each time t we require the optimality condition
uU
Given an initial condition (depending on
x(0)) at which both x and p reach the origin:
Lastly 0  t<T we require
Observe that (30) is the Hamiltonian system p). This is intimately
connected with the propertyt that p(t)=DV (x(t)) for a solution of Hu (x, DV We will return to
this issue, near the end of Section 3.2 explaining why the manifold of (x, p) formed by our solution family is
truly the graph of a gradient also that for 0  t  T we have
The formula (28) for Hu shows that (34) is indeed satisfied at
(32). It is a general property of Hamiltonian systems as in (30) that the value of Hu (x(t),p(t)) is constant
with respect to t. Property (31) implies that the jumps in u(t) do not produce discontinuities with respect
to t in (34). Therefore (34) follows as a consequence of (30) - (32). Thus (34) and (31) give us (29) for
u(t) in particular. The construction of x(t),p(t) in Section 3.2 will show that u(t)  (x(t)) and that
extends to all u  (x(t)). This will provide the saddle point conditions (29) on the interior.
The equation H(x, DV has many solutions, if it has any at all. One property of the particular
solution we want is that V be associated with the stable manifold of (30), in accord with the general approach
of van der Schaft [15, 16, 17]to robust nonlinear control. We see this in the convergence to the origin of
above. Another important property is that

To this end, notice that the formula (28) for an individual Hamiltonian, together with (34), impliess that
So for p(t)=DV (x(t)), (33) is the same as saying
d
dt
If we stipulate that V which is (35). One
may wonder why we have insisted on (32). Observe that (33) (in the limit as t  0) implies that
necessary if
A family of x(t),p(t) as described above will give us a function V (x) which has the desired saddle point
properties at interior points. However for x  K both (25) and (26) are complicated by the nontrivial
structure of (x, v). We claim the V (x) so constructed does in fact satisfy the saddle point conditions (25)
and (26) at x  K as well. We do not give a mathematical proof of this. Instead we have developed a
scheme of numerical confirmation that can be applied to test this claim for any specification of si. This is
described in Section 4. We also note the requirement in (32) above that x(t)   for all 0  t<T, given
. This follows if we can verify that whenever x(t)  iK then
ni  (p  Gu)  0.
We rely on numerical tests for this fact as well. (See the discussion of (52) in Section 4.) Based on the
success of these tests for numerous examples, we conjecture that (25) and (26) are true in general. The
reference to the boundary verifications of Section 4 in Theorem 1 indicates that the validity of that result
depends on the success of those tests.
3.1. Identification and Properties of the Invariant Control Vectors. We will construct the family
x(t),p(t) as above by generalizing the development of [2]. The key is to look for solutions that approach
the origin as in (32) using a constant control u. The solution of (30) with constant
conditions x(T)=0=p(T)is
Observe that for 0  t  T  /2 the values of both  sin(t  T) and 1  cos(t  T) will be positive. Now
consider what (31) requires of (38):
uU
There are only a finite number of such they provide the key to the explicit representation of
the family of solutions x(t),p(t) that we desire.
We will call any   GU satisfying (39) an invariant control vector. To simplify our discussion here, let
denote the columns of G. (In more general models, gi would be the extreme points of GU.) To say
GU means that  is a convex combination of the gi: 1. For an  as in
consider the set of indicies
It follows from (39) that every j  J achieves the maximum value of   gi over i  N. Therefore

Our construction of V (x) depends on the fact that there is a unique such associated with every
nonempty subset J  N. The existence of J depends on properties of our particular set of gi, but the
uniqueness does not. So we present the uniqueness argument separately as the following lemma.
Lemma 2. Suppose gi, i =1,. ,m are nonzero vectors in Rn and J {1,. ,m} is nonempty. If there
exists a vector J as described in (40) then it is unique. Suppose J and J exist for both J  J.Then

Proof. We establish (41) first. Without assuming uniqueness, suppose exists as in (40). It
follows that

Now suppose both J and J exist for J  J. Then the same reasoning implies that


which is (41). Regarding uniqueness, suppose for the
same J. In that case (41) implies
But then (42) implies gj J. This means  is orthogonal to the span of {gj,
But since it is also in the span, we are forced to conclude that  =.
It is not dicult to determine whether or not J exists for a given J. If it does, the values
must be a nonnegative solution of the linear system


From such a nonnegative solution we can recover j from
check (40). With this observation we can prove that J exists for all J  N in our single-server model.
Theorem 2. Assume the specific G and U0 of our model (see Section 2). For every nonempty J  N there
exists a unique invariant control vector J . Moreover the j, j  J in (40) are strictly positive.
Proof. Let GJ =[gj]jJ be the matrix whose columns are the gj for just those j  J. Observe that (43)
simply says J =[j]jJ must solve
GTJ GJ J =1J .
For the existence of nonnegative j in (43) it is enough to show that GTJ GJ is invertible and that all entries
of its inverse are nonnegative. Consider the diagonal matrix
and let
Note that MJ is nothing but GJ for the particular case of all
it is enough to show that (MJT MJ )1 exists and has nonnegative entries. Now observe that MJT MJ is block
diagonal
. Ak1 0
where the A and B are tridiagonal of the form
121 .  121 .
One may check by explicit calculation that denoting the size of A,
c +1
Since all entries are positive in both cases, it follows that all entries of (MJT MJ )1 and hence (GTJ GJ )1
are nonnegative, as desired. Since no rows are identically 0, all i are positive in (43) and therefore the
respective j > 0 can always be found.
Next, we need to show that gj  J >gi  J for i/J, j  J. First observe that gj
constant over j  J. Also note that for
We know that j > 0. So if i/J,

Therefore, gi  J  0 for every i/J,andgj  J >gi  J . Lemma 2 gives the uniqueness.
Observe that  of (17) is a scalar multiple of N . Indeed, in the notation of the above proof,
It follows that  = N /N  N . In particular the  of (20) is alternately described as
Fundamental to our construction is the existence and uniqueness of the following representation of x
using a nested sequence of invariant control vectors.
Theorem 3. Assume the specific G and U0 of our model. Every nonzero x   has a unique representation
of the form

for some 0 <aj, aj < 1 and J1 . Jk  N. Moreover,
Proof. Consider any nonzero x  . We first solve
G.
The reader can check that G1 has all nonnegative entries, which implies that all i  0. Therefore every
x  K can be written as
Next, let and consider the invariant control vector
Now consider
(j  aj)gj
Our choice of a implies
for all j  J. However for one or more j  J,j  aj =0. By induction on the number of positive
coecients in (47) it is possible to write
J. Simply taking am = a and completes the induction argument.
Next notice that since Ji

2.

From the hypothesis that x   we conclude that aj < 1.
Now consider J1 in (45) and any j, j  J1. Then j, j  Ji for all i, which from (40) tells us that
However if j  J1 but j / J1 then
while for i>1, gj  Ji  gj  Ji (depending on whether j  Ji or not). We conclude that
gj  x>gj  x.
This proves that J1 is the set of j for which gj  x takes its largest possible value, as claimed.
Regarding uniqueness, since G is nonsingular the  in are uniquely determined, and then Jk from
the last term of (45) is necessarily the J above. Since Jk1 Jk,
k1
(j  akj)gj
still must have nonnegative coecients j  akj. But only those for j  Jk1 can be positive. Hence
J. This implies that ak = a as well. Thus Jk and ak are uniquely determined.
Repeating the argument on
k1
ajJjgives the uniqueness of the other aj,Jj .
The following lemma records two other facts that will be used below.
Lemma 3. Assuming the G and U0 of our model,
a) All coordinates of N are positive.
b) If x  K is as in Theorem 3, and if xi =0, then i/J1.
Proof. We have already observed that is as in (17). This proves a). Next suppose
It follows that x  gi  0. On the other hand, there does exist j with x  gj > 0 (take j to be largest
with xj > 0 for instance). Thus the set J1 of those j with
does not include i.
It is important to realize that the single server re-entrant queue being studied here is special in that a
unique J is defined for every J  N. This is not the case for more multiserver systems. For example,
consider the re-entrant line with two servers in Figure 2. The gi are Gui, ui  U0 as in (3). Simple test
calculations reveal many J for which no J exists. Moreover it turns out that the gi are linearly dependent
so the points representable as in (45) can account for at most a 3-dimensional subset of
K.
3.2. Construction of the Characteristic Family. We can now exhibit the desired family of solutions to
(30). Consider a nested sequence J1  J2 Jk, and parameters
coecient functions ai(t)andi(t)acordingtotheformulas
(48) .
.
.
(Here again we use y+ to denote the positive part: k. For all 0  t  T we
have
so that all ai(t)andi(t) are nonnegative. We claim that
provide a solution of of (30)-(33). The derivation of (48) is based on the calculations in [2]. Here we will
simply present as direct a calculation as possible. To that end, consider the partial sums appearing on the
left in (48):
Consider t in one of the intervals 1 <t<. Then for i<we have so that
For
Taking pairwise dierences we see that

Using this in (49), we find that for 1 <t<
Thus (30) is satisfied for 1 <t< using To confirm property (31) on that interval, observe
that since j(t)=0forj<,
we know from (46) that p(t)  Gu is maximized over u Uat any u
for which only the j  J coordinates are positive, in particular for
Implicit in this construction is a function defined by means of (49). Starting with x  ,
express x as in (45). Then determine p(x)using
iJiwhere the partial sums of the i are determined from those of the ai according to
for 0 <i  /2. This is the gradient map p(x)=DV (x) of our solution to (23). There are several facts to
record about p(x) before proceeding.
Theorem 4. The map p(x) described above is locally Lipschitz continuous in  and satisfies the strict
inequality
for all x  , x =0.

To see this, first consider what we will call a maximal sequence J1  J2  . Jn, i.e. and each Ji
has precisely i elements, with N. Consider the x representable as in (45) for this particular maximal
sequence. The maps x  ai  ai and i  i  p are linear. The maps ai  i are simply


These are Lipschitz so long as ai remains bounded below 1. Since ai  n in , we see that
p(x) is indeed Lipschitz in any compact subset of , constrained to those x associated with a fixed maximal
sequence of Ji. If in (45) we relax the positivity assumption to ai  0, then we can include additional Ji so
that every x   is associated with one or more maximal sequence of Ji. The Lipschitz continuity argument
extends to the x associated with any given maximal sequence in this way. To finish the continuity assertions
of the theorem we need to consider x on the boundary between the regions associated with distinct maximal
sequences:

The uniqueness assertion of Theorem 3 means that the nonzero terms of both representations agree, which
implies that the corresponding terms of the expressions for p also agree:

Thus p(x) is continuous across the boundaries of the regions associated with dierent maximal sequences.
From this it follows that the (local) Lipschitz continuity assertion of the theorem is valid in all of .
The argument that p(x)2 < x2 is the same as [2, pg. 334, 335]. The strict inequality comes from the
fact that

Equality occurs only for n 1). Notice that for x  , x
sin(n) < 1 implies that p(x)
The argument given in [2]that p(x) is indeed the gradient of a function V (x), x   also generalizes to
the present context. In brief, the standard reasoning from the method of characteristics can be applied to
each of the individual Hamiltonians Hu where to see that DV (x)=p(x)
in the region associated with a given maximal sequence Ji. Continuity across the boundaries between such
regions allows us to conclude that there is indeed a C1 function V in  with DV (x)=p(x). Taking V
implies that V (x) > by virtue of the discussion of (35) above.

Finally, we return to the connection of (31) with (x) and (29). Since U is convex,
and so (x) consists precisely of those u in the convex hull of ej, j  J1, J1 being as in (45) for x. Because
in p(x)= iJi the i are positive when the corresponding ai are positive, we see that (x) has the
alternate description
U
In particular, the u(t) of (31) belongs to (x(t)). Moreover p(x)Gu has the same value for all u  (x),
so
as desired in (29).
4. Verification of Conditions on the Boundary
We have completed the construction of V (x) satisfying (29) on the interior of . We now consider the
assertion of that for x  K the resulting V remains a solution when the H of (27) is replaced by H as in
(27), and that any u  (x) is a saddle point, as in (25) and (26). Specifically, we want
to confirm that for a given x    K, its associated any u  (x), the following hold:
Since we know Hu (x, imply (25), and (52) and (54) imply (26). Together these imply
(23).
Our validation of (52) - (54) consists of extensive numerical testing, as opposed to a deductive proof. We
will describe computational procedures below. Test calculations have been performed on numerous examples
(see Section 4.4), confirming (52) - (54) to within machine precision in each case. This gives us confidence
in the theoretical validity of (52) - (54), but until deductive arguments can be presented, their theoretical
validity must be considered conjectural.
4.1. Inactive Projection. Given x  , the corresponding any u  (x), (52) is
equivalent to the statement that
ni  (p  Gu)  0 for all i  I(x).
This would be easy to check by direct calculation at a given x. However the second part of the following
lemma provides an equivalent condition which is even easier to check.
Lemma 4. The following are equivalent
1. ni  (p(x)  Gu)  0 for all x    K with x =0,alu  (x),andi  I(x);

2. p(x)i  0 for all x    K and i  I(x);
3. p(x)i > 0 for all i and all x =0in the interior of .

Proof. Clearly (2) follows from (3) by continuity of p(x). To see that (2) implies (1), recall from our discussion
in Section 2.2 of the fact that  is a nonidling policy that that ni  Gu  0 for any u  (x). Therefore
(2) implies
ni  (p  Gu)  pi  0.
Finally, observe that (1) implies that the characteristic curves (49) do not exit K in forward time. From
any x(0) in the interior of , x(t) remains in K up to the time T at which
it follows that pi(x) > 0 for all i in x is in the interior of .
4.2. Control Optimality. Now we consider an approach to checking (53) at a given x  K with its
associated any u  (x). We want to check that u is the minimizer of
over u U. Observe that by virtue of (52)
Since pGu has the same value for all u  (x) it suces to consider any single u  (x) and to show that
it gives the minimum of p(x, pGu)overu U. Since this is a continuous function of u and U is compact,
we know that there does exist a minimizing u. Moreover for some F  I(x), p(s, pGu)=pRF (pGu),
according to (13). So given F we can identify u as a maximizer of p  RF Gu subject to the constraints of
Lemma 1 part 2). If (53) were false then an exception u would occur as a solution of such a constrained
minimization problem, for some F  I(x).
There is no exception to (53) for because in that case
solves a standard linear programming problem:

subject to u U
BF (p  Gu)  0, and
If u is an exception to (53) then so is any feasible maximizer uF to (55):
To verify (53) computationally we invoke a standard linear programming algorithm for (55) for each
nonempty subset F of I(x), and for each feasible maximizer so found, check that
We note that when I(x)={i} is a singleton we only need to check itself. In this case it is
sucient to check that
directly for each of j =1,. ,n. To see why, first observe that the last constraint in (55) is satisfied
vacuously. Since BF (p  Gu) > 0, the same must hold for all u Usuciently close to u. It follows that u
gives a local maximum of p  RF Gu over U. Since U is convex, it must be a global maximum. Therefore u
must be a convex combination of those ej for which p  RF observe that since
the constraint
BF (p  Gu)=ni  (p  Gu) > 0
is a scalar constraint. It must therefore be satisfied by one of the ej for which p  RF
This means that this ej also solves (55). Hence when I(x) is a singleton it suces to check just the ej as
candidates for u, rather than invoking the linear programming algorithm.
4.3. Load Optimality. Once (52) is confirmed we know that for any u  (x), (x, q Gu)=q Gu
and that
with respect to those q for which (x, q  Gu)=q  Gu, and that the maximal value is 0. To verify (54)
we need to be sure that there are not some other u  (x) and q with (x, q Gu) =q Gu and for which

Since (x, v) is continuous and piecewise linear, and since (x) is a compact set, it follows that there does
exist a u  (x) and q which maximizes (56) over q  Rn and u  (x). We derive necessary conditions
contingent on the specification of the subset F  I(x) for which (x, q Gu)=RF (q Gu). Using part 3
of Lemma 1 we know that u =uand q =qsatisfy
Consider the ane set of all q satisfying (58). Since the inequalities are strict in (57), all q near q and
satisfying (58) must also have (x, q  Gu)=RF (q  Gu). Thus q =qis a local maximum of(59) p  RF (q  Gu)  q2, subject to the constraint NLT\F RF (q  Gu)=0.A simple calculation shows that this implies
where PL\F is the orthogonal projection onto the kernel of NLT\F the constraint (58) is vacuous
and we take Substituting this back into (59) and considering the result as a function of u,it
follows that u =uis a local (and hence global by convexity) solution of the quadratic programming problem:
subject to u  (x),
BF PL\F (RFT p  Gu)  0, and
To verify (54) computationally, we consider all pairs of subsets F  L  I(x). For each, we invoke a
standard quadratic programming algorithm to find a feasible maximizer u, if any exists. If such a u is found,
we take
and then check by direct calculation whether this is an exception to (54), as in (56). If we consider all
F  L  I(x) but find no such exceptions, then (54) is confirmed for this x, p.
Again we note that the quadratic programming calculation can be skipped in some cases. If
then (x, q  Gu)=q Gu and we know there are no such exceptions to (54). Thus only need be

considered. Secondly, suppose I(x)={i} is a singleton. Then the only case to check is In
that case if there is an exception to (54), q must maximize
and satisfy
It follows that
But for the latter inequality simplifies to
ni  (p  Gu) <di  p.
Moreover since di = ni  ni+1 (with this is equivalent to
But i  I(x) means i/J1,soni  Gu  0 for all u  (x). So (61) would imply pj < 0 for some j.If
we have already checked that p  0 in accord with Lemma 4 and our confirmation of (52), then we can be
sure no exceptions to (54) occur when I(x) is a singleton. Thus we only need to appeal to the quadratic
programming calculations when two or more xi are zero.
4.4. Test Cases. We begin our test of (30)-(33) for a specific choice of parameters s1,. ,sn by calculating
all the invariant control vectors J . Then on each face iK a rectangular grid of points x  iK with
x N  N N is constructed. For each grid point x we then compute the representation (45) and then the
associated according to (51). We then check that all pi  0 in accord with Lemma 4
and carry out the constrained optimization calculations described above for all possible  F  L  I(x).
Obviously, the amount of computation involved will be prohibitive if the number of dimensions n is significant.
However, for modest n the calculations can be completed in a reasonable amount of time. We have carried
out these computations for numerous examples, including the following:
(s1,. ,sn)=(1, 1, 1)
No exceptions to (52)-(54) were found.
5. Proof of Optimality: Theorem 1
We turn now to the proof of the optimality assertions of Theorem 1. By hypothesis V (x) is as constructed
in Section 3.2, the saddle point conditions (25) and (26) have been confirmed, as well as the equivalent
conditions of Lemma 4. We know that () satisfies the minimum performance criterion of Section 2.4. As
explained above, V (x) > 0 for all x   with load q(). The argument of [2, Theorem

2.1]shows that with respect to (), on any interval [0,T]on which x(t) remains in , we have
dt.
For a given x(0)  , let x(t),p(t) be the particular path constructed according to (30), with x(T)=0.
We know that u(t)  (x(t)) so x(t) is the controlled path produced by  in response to the load q(t).
Along it we have from (36) that
and therefore

dt.
This establishes (21).
Next we consider an arbitrary strategy  satisfying the minimum performance criterion. We would like to
produce a load q(t) which is related to the resulting state trajectory x(t)byq(t)=DV (x(t)). In [2, Theorem
2.3]this was accomplished by limiting  to state-feedback strategies and appealing to an existence result for
Filippov solutions of the dierential inclusion [2, (2.23)]. Here we only approximate such a load. By taking
advantage of the properties of (x, ), our argument will not be limited to state-feedback strategies, and will
not need the Filippov existence result.
Given x   we will show that for any =>0 there exists a load q(t) satisfying qi(t)  0andq(t)    1
for all t>0, and such that (for some u(t)=[x, q()](t))
holds for all T. The dicult question of existence for the closed loop system
for an arbitrary strategy is easily resolved by introducing a small time lag:

The system can now be solved incrementally on a sequence of time intervals [tn1,tn]where
For t  [tn1,tn]the values of q(t) are determined by x(t) on the previous interval [tn2,tn1], so the basic
existence properties of the system under  subject to a prescribed q(t) insure the existence of x(t)andq(t)
as above. Let u(t)=[x, q()](t) be the associated control function. Since q(t) is always a value of DV (x)
at some x  , we know qi(t)  0andq(t)    1, and the minimum performance hypothesis insures that
x(t) remains in a compact subset of : x(t)    , <1. We must explain how the time lag leads to the
+= term in (62).
Observe that because (x, v)=RF v for one of only a finite number of possible matrices RF , and because
there is a uniform upper bound on x:
Consequently,
x(t)  x(t  =et)B=et.
Next, on the subset of x   with x    , DV (x) is Lipschitz; see Theorem 4. It follows that for some
constant C1 (independent of =) such that
We observed previously that (x, v) in Lipschitz in v. It follows that for some constant C2 and all =, t > 0
We know that
it follows that
integrating both sides over [0,T]and replacing = by =/C2 yields (62).
With this q(t) in hand the remainder of the argument proceeds as in [2]: if there exists a sequence
with x(Tn)  0, then V (x(Tn))  0 in (62) which implies
dt.
Suppose no such sequence exists. Then in addition to x(t)    <1 we know x(t) does not approach
0; it must remain in a compact subset M of  \{0}. From (63),

dt.
WealsoknowfromTheorem4that1 x2  1 DV (x(t))2 has a positive lower bound. Therefore the right
side in (64) is infinite. Thus (64) holds in either case. Since =>0 was arbitrary, (22) follows.
6. An Example with Restricted Entry
In this section we reconsider our model in modified so that the exogenous load
only applied to queue x1. This is illustrated in Figure 4. The system equations are now
where q(t) is a scalar and
,while the control matrix

the control values u Uand the constraint directions di all remain as before. We carry out the same general
approach to constructing V (x) as outlined at the beginning of Section 3. The details of the analysis are
dierent in several regards. This is significant because it shows that our general approach is not exclusive to
all the structural features of Sections 3.1 and 3.2. We will find that the optimal policy is the same (x)as
given in (15) above. In higher dimensions (n>2) it is interesting to speculate whether the optimal policy
would likewise remain unchanged if we removed the exogenous loads qi(t), i>1 . However, at present this
has only been explored in 2 dimensions.

Figure

4. Re-entrant Loop with Single Input Queue
The presence of M in (65) changes the individual Hamiltonian:

p1 being the first coordinate of p =(p1,p2). The supremum is achieved for p1. The corresponding
Hamiltonian system, for a given u U,is
We calculate the invariant control vector

as described in Section 3.1. Other than there are no additional J to consider.
To simplify notation we will drop the subscript N:

The first place we find a significant dierence from our previous analysis is in the calculation of a
solution to the Hamiltonian system associated with , analogous to (38). Previously, we did this using

being as determined by the construction of
because of the missing p2 term in the x2 equation of (66), we must use a Gu(t) which is both dierent from
and time dependent. We seek a solution x(t)=a(t),p(t)=(t) (both a(t)and(t) nonnegative) to

for some function 0  1. The overbar on x, p distinguishes this special solution from the others
encountered below. In light of the p equation and the terminal conditions (32), the solution we seek must
be of the form
for some function (t)  0. is a scalar multiple of , the right side of the x equation in
must also be a scalar multiple of . Since implies a relationship between (t)and(t),
which works out to be
Using this we can reduce (67) to a single second order dierential equation for (t):
(t)+A(t)=1,
where A is the constant
The solution (for initial conditions It will be convenient
for the rest of this discussion to fix so that
(One consequence of fixing is that for a given x the t<0 for which x(t)=x depends on x.) For
t  0 we confirm that 0  (t)  1, (t)  0anda(t)= (t)  0, as we wished.
We now have the desired solution:
This special solution provides the final stage of our family of paths x(t),p(t) as in (30)-(33) of Section 3,
but with some adjustment. In contrast to Section 3, our u(t)=[(t), 1(t)]T varies continuously, instead
of being piecewise constant. This means we have to pay closer attention to (34). Once again it is satisfied
at by virtue of the terminal conditions. When we calculate dt Hu(t)(x(t), p(t)), one term does not
automatically drop out:
d
dt
Since p is a scalar multiple of  and we know   g2, we do indeed find that Hu(t)(x(t), p(t))  0.
The analogue of (33) for this example is
This is because q2 x2 =(p1(t))2 < x(t)2 when q(t)=p(t). Also note that
taking advantage of the fact that Hu (x, to verify (33) along x, p in particular, simply observe
that
since both  and  are positive (excepting In the following xi(t), pi(t) will refer to the individual
coordinates of this particular solution.
Our special solution x(t), p(t) provides the final stage (t1 <t 0) for each of the solutions in the family
described at the beginning of Section 3. The initial stage (t<t1 < 0) will be a solution of (66), with
either e1 or e2, which joins x, p at some t1 < 0: x(t1)=x(t1), p(t1)=p(t1). In other words we solve (66)
backwards from x(t1), p(t1),for the appropriate choice of u. It turns out that using produces that
part of the family which covers a region below the line x() in the first quadrant, and using
the x(t) which cover a region above x(). This is illustrated in Figure 5, for parameter values s1 =4,s2 =1.
Note that the region covered by this family, and hence the domain  of V (x), is no longer the simple polygon
of (20).1.410.60.20

Figure

5. Characteristics for Restricted Entry Loop
We will need to verify that the resulting family indeed satisfies all the conditions outlined in Section 3.
These verifications are discussed below. Once confirmed, this implies that the optimal control (x) produces
e1 if x is below the line x(t), e2 if x is above the line, and any u Uif x is on the line. So although we
will not produce as explicit a construction for x  p as we did in Section 3, we still find the same optimal

control
uU
6.1. Interior Verifications. We have already discussed properties (30)-(33) of Section 3 for the final stage
of our family of solutions: x(t)=x(t), p(t)=p(t)fort1  t  0. However we still need to verify (31) and
(68) for the initial segment t<t1. In Section 3.2 this followed from properties of the J and the rather
explicit formulae for x(t)andp(t) in terms of them. Here we have not developed such an elaborate general
structure. Instead we resort to direct evaluation of the needed inequalities. By solving (66) for
x(t1)=x(t1), p(t1)=p(t1) we obtain the formulas for the lower half of our family: for t<t1 < 0,


For any   <t1 < 0, the above will be valid for t<t1 down to the first time at which x(1)(t) either reaches
the horizontal axis, or reaches the outer boundary of , curves appearing in
the figure, b(t1) <1(t1).) A formula for 1() is easily obtained from the expressions in (69). The value
can be identified as the point at which the determinant of the Jacobian of x(1) with respect to
vanishes. An explicit formula is possible for b() as well. (For brevity we omit both formulas.) Thus
is valid for
The points on the horizontal boundary 2K are x(1)(1(t1)) for those t1 with b(t1)  1(t1).
The analogous formulas for the upper half of our family are obtained by solving (66) for
x(t1)=x(t1), p(t1)=p(t1) to obtain the following expression for t<t1 < 0:


This time, for a given   <t1 < 0, the valid range of t<t1 is slightly dierent. It turns out that x(2)(t)
always reaches the outer boundary of , at a time prior to contacting the vertical boundary 1K.
(Once again, an explicit formula for b() is obtained by setting the Jacobian of (70) equal to 0.) Thus given
t1 < 0, (70) is valid for
The vertical boundary itself is traced out by the solution for t1 =0:

x (t)=
valid for b(0) <t 0.
The availability of these formulas makes it possible to check the inequalities we need for (31) and (33).
For (31) we want to verify
For (73), it turns out that
),which certainly is positive for t<t1. We resort to numerical calculation to confirm (72). We have already
noted that (33) should be replaced by (68):
x(i)(t)2  (p(i)(t))2 > 0,for both i =1, 2. It is a straightforward task to prepare a short computer program that, given values
for s1,s2, evaluates (72) and (68) for a large number of t<t1 pairs extending through the full range of
possibilities. In this way we have confirmed the above inequalities numerically.
6.2. The Horizontal Boundary. Finally we must consider the influence of the projection dynamics at
points x  K, confirming as we did in Section 4 that our remains a saddle point
when (x, v) is taken into account. This entails checking the same three facts, (52), (54), and (53) as before.
We consider the two faces of K separately.
The reflection matrix for 2K is

and
(x, v)= .
Now observe that
which is independent of q. Since it follows that (x, Mq  Gu)=
Mq  Gu, so that (52) reduces to (72). Moreover this independence of q also implies (54) since we know
is the saddle point in the absence of projection dynamics.
Next consider (53). The u Uare just

Note that corresponds to  = 1. So for (53) we want to show that the minimum of
occurs at little algebra shows that for 0    s1s+2s2 we have
so that
For s1s+2s2    1wehaven2  (Mp1  Gu)  0, so that
Thus the function of  in (75) is piecewise linear, in two segments. The slope of the right segment ( s2
1) is
which we already know to be negative, by virtue of our work in checking (72). Thus to establish (53) we
only need to check that the value for greater than that for  =0:
which is equivalent to p  g1  0. This we have confirmed numerically, by evaluating
for various choices of s1,s2 and t1 throughout its range.
6.3. The Vertical Boundary. Recall that along 1K we have and that from (71) we know
Thus (x, Mq  Gu)=(Mq  Gu), confirming (52).
The reflection matrix on 1K is

We already know that
over those q for which (x, Mq  g2)=Mq  g2. We need to consider the possibility of a global maximum
among those q with (x, Mq  g2)=R{1}(Mq g2), namely q with n1  (Mq g2)=q  0. However,
is maximized at its maximum over q  0 must occur at confirms (54).
Finally, we turn to (53). Since any u as in (74) we have
Therefore, after a little algebra,
which is minimized at corresponds to verifies (53).


--R


An Introduction to Variational Inequalities and their Applications
Queueing Systems: Theory and Applications
Reflected Brownian motion on an orthant
Neumann type boundary conditions for Hamilton-Jacobi equations


Nonlinear state space H1 control theory

On optimal draining of re-entrant fluid lines
Department of Mathematics

--TR
skorokhod problem;queueing;robust control