'1‘!" "3"“

 
 

”Into!” 1",.

H?"
L 1.175.111 J “'41.”:
\‘ '-‘:

 
 
     
 
  

  

   
 
   

   
   

   

  
      

"Pal-"3”

"at 3'1'I1I 1‘3"“. ;"""""‘A
39111

"W .
5" 1

   
 

 

v

  

1'1 "‘

      

"as“ .1
r"31"

 
 
  

 

I. ll

  
  
  
 

     
 
  
  

H'H'Inr £1?" ”*5";\L:‘V

1'3}- , 1

:‘f 3’ . 3}." v ‘. H 3' [Wu

l 1 a1 ‘5'! 1_ I ._ "MI ". K

"“37",!“ “If: 11 "1-3": 'g“:1"fl$\: fro.” "3-, :'fl:"fi"l 19"." 11‘" '15-‘31‘6

”‘3' ' ”pr‘fl ("“1“ "int“ 1,“ '2 "‘11V;'-"". ' In}; "21“"-

u r" "' ""23"": Dill-"4"";j‘i1—rvpflaw-‘.'."}I"‘"-""‘"3‘41' 17'1""‘juanAV‘KrV'

01'1” 1. L~15~11."1'1 1111.111, 11:»,

":33 )p‘é ¥"""I “I I I ""31"“ IQ”; E"I;\"3\33.T'1" :IliIfllI' ""1"“ $1312}, 3‘; ‘}r:‘
a- ,$ “"7” 0'1 (""5" “11““: NM

912331-11»? '0} 33.3151; "”3 " A“ 1 Ema

   
 
    
         

1 11 ‘l tn'I 1W‘3V'H
45:“... 3""‘3L'H'1I.

" “""“' '1‘? 111‘ Jr
1

1 .

- .1.'J'\2'1E.§§‘,'11'1'4.1111""'1'"‘1{(I
'3'91'1">"1".;'d'ai, 1
u . 1:11.

    

       
       
    
   
 

t- a
I
"1""; "7‘03““ p. :l‘",""1'"¢"f"lu >N".'.'_"fi;fi
V 1:3 r0133" "‘3" l‘ TV '
3‘1“: a?" 1." I .71! 1 “K1
9' ‘ ILL." I‘I :"'j}$."q:" ‘1‘: ”E:
"1«I\' .M ”fffitiuflu‘ If
v‘xI: "133 I :.:.'..:'
1‘ ‘35:“! -L .1 '1

     
  
   
 
 

.31 _
.1:

I9;

‘ )5: 1:1"

    

"3).: _.AL'\'

  

       
 

  
    

u ""3"!“ ‘1...” V:
3' "1% . 1113'“.
1I ivy-1L .
Lf

  
  

i'":"II I |.. x
V‘ '1'.) f‘“:“ “v.1!-

    

    

1 1 1' 1 _
' “1. '. "11 1 ' '1' I
‘ _ n I
9 5 1371' 1:.lw" 1'
1 n , ,1 1|
1 1 1 .
1111

 
 

31" 1.21111 "131.13 3:" 2‘

V1. 11',II\AIH ‘1)

(311-1 ll' iii-L1 "I "u“.

 

 

all i *2
LI}. "
EJCINA I ., f . "‘

3‘ h'
Y ‘ -. .
meCIlirl’r' . .
’ £1“

WWW?“- " 'WUM'h.‘

THESIS

   

This is to certify that the

thesis entitled

Existence and Computation of Static
Equilibria in Certain Economic Models
with Application to the PIES Model
of the Energy Sector

presented by

Paul Arthur Rubin

has been accepted towards fulfillment
of the requirements for

Ph. D. degree in Mathematics

MP. [W

Major professor

 

Date 5 980

0-7639

 

 

 

 

'ovenou rut: : '
25¢ per day por tu- -
RETHRNING LIBRARY MATERIALS:

Phco in book return to move
charge from circulation records

 

EXISTENCE AND COMPUTATION OF STATIC
EQUILIBRIA IN CERTAIN ECONOMIC MODELS,
WITH APPLICATION TO THE PIES MODEL

OF THE ENERGY SECTOR

BY

Paul Arthur Rubin

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

DOCTOR OF PHILOSOPHY

Department of Mathematics

1980

J; /

afloV

ABSTRACT
EXISTENCE AND COMPUTATION OF STATIC
EQUILIBRIA IN CERTAIN ECONOMIC MODELS,

WITH APPLICATION TO THE PIES MODEL
OF THE ENERGY SECTOR

BY

Paul Arthur Rubin

In this thesis we examine the Project Independence
Evaluation System (PIES) integrating model, a static
equilibrium model of the energy sector of the national
economy. After carefully formulating the equilibrium
problem, we establish conditions sufficient to ensure the
existence of an equilibrium, as well as conditions sufficient
to guarantee uniqueness of that equilibrium. We study the
PIES algorithm, prove that in certain cases it must converge
to an equilibrium, and exhibit an example on which it fails
to converge. The same is done for the PIES—VAR variant
algorithm. We pose a minimization problem related to the
task of locating equilibria, and propose a subgradient-
based algorithm for that problem. Finally, we describe the
implementation of the subgradient algorithm and discuss the

results of some computational trials.

To all who waited.

ii

ACKNOWLEDGMENTS

I would like to thank my advisor, Professor V. P.
Sreedharan, for his guidance in the research leading to
this thesis. I am also grateful to the members of my
guidance committee for their assistance, and in particular
to Professor P. K. Wong for his careful reading of the
dissertation. I further thank the faculty and staff of
the Department of Mathematics at Michigan State University
for their support of my efforts as a student. Finally,

I wish to express my gratitude to Mrs. Mary Reynolds and
Ms. Cindy Balzer for their excellent work in preparing the

manuscript.

iii

List of
Chapter
Chapter
Chapter
Chapter
Chapter
Chapter
Chapter

Chapter

List of

TABLE OF CONTENTS

Tables

0 Introduction

I The Integrating Model

II Existence of Equilibrium
III Uhiquness of Equilibria

IV The PIES Algorithm

V The PIES-VAR Algorithm
VI A Subgradient Projection Algorithm
VII Implementation of the Subgradient

Projection Algorithm

References

iv

page
page
page
Page
page
Page
Page

Page

Page

Page

21
33
42
65

76

101

110

LIST OF TABLES

Table 1 Page 51

CHAPTER O

INTRODUCTION

The Project Independence Evaluation System [5 ]
henceforth denoted PIES, is an aggregation of models
Which describes the energy sector of the national economy,
developed by the Federal Energy Administration (now
part of the Department of Energy) as a tool for the
evaluation of policy decisions. The various components
of PIES model the production, refinement, conversion,
transportation and consumption of a variety of energy
commodities. Of particular importance is the estimation
of a static partial equilibrium for the energy sector,

a vector of prices at which supply and demand will be in
agreement. The component models take as parameters
factors such as tax policies and the pricing of crude oil
by foreign producers. The static equilibrium predicted
by PIES based on specified values for these parameters

is taken as an indication of the expected market response

to those policies.

Central to the estimation of the static equilibrium

is the integrating model [ 5,I7], which computes the

2

equilibrium based on supply and demand functions generated
by other component models. The integrating model employs
an iterative technique which successively refines approxi-
mations to the equilibrium demand vector. Although
observed to converge rapidly in the examples reported by
Hogan and Wagner, the algorithm has evaded to date a
complete theoretical analysis regarding convergence, al—
though it has motivated a significant amount of research
and several variant algorithms [3 ,6 ]. In this paper

we attempt to shed some light on the existence and deter-
mination of equilibria. In particular, we propose a new

algorithm for computing the equilibria.

We begin with a mathematical formulation of the
integrating model, introducing the necessary concepts from
convex analysis. We then prove the existence of an equilib—
rium under fairly general conditions, conditions consistent
with the PIES models. We note that the authors of PIES,
in their publications, have assumed rather than demonstrated
the existence of equilibria. uniqueness of the equilibrium
has been demonstrated [15] under a number of hypotheses
made by the PIES modellers. We show that one of these
hypotheses appears to be inconsistent with the form of the
demand model. and that in the absence of that hypothesis

the equilibrium need not be unique.

We examine the PIES algorithm for computing the
equilibrium, exhibiting an example in which the algorithm

fails to converge to a solution. This explains the failure

of others to produce general convergence results, and
suggests the need for either additional hypotheses or a
different algorithm. We also examine one variant of the
PIES algorithm, the PIES-VAR algorithm [6»], and eXhibit

an example in which it fails to converge.

We suggest an algorithm for solving the equilibrium
problem in essentially the same form as that assumed by
the PIES algorithm. The algorithm we present is applicable
to a general class of problems. The proof of convergence
of the algorithm Which we present assumes that the demand
function has a potential. Implementation of the algorithm
does not require the existence of such a potential, which
suggests the distinct possibility of proving that the
algorithm converges even where the demand function has no

potential. We shall return to this problem at a later date.

CHAPTER I

THE INTEGRATING MODEL

We consider the problem of determining a static
market equilibrium in a context slightly more general
than that of the PIES algorithm [ 5]. Let us assume
that we have a finite collection of goods indexed by
the integers l,...,d. In the PIES models, goods repre-
sent various energy products, differentiated by type
of energy, region of production and region of consumption.
Wagner [17] suggests that d = 54 is typical for the

PIES model.

We represent demands, supplies and prices of these
goods as vectors of length d, i.e., elements of the
d-dimensional euclidean space I51. Let us pause to
establish some notation. We do not, in general, distinguish
between row and column vectors: the context determines
the shape of the vector. In particular, the usual inner-

product of two vectors will be denoted by the juxtaposition

of those two vectors; that is, for x,y 6 Rd ,
2
xy = '24 xiyi.
l=l

We denote by R3 the nonnegative orthant

[x e 151: xi 2,0, 1 = 1,...,d}.

When ordering vectors x,y E Eta, we adopt the following

notation:
X2y iff Xi2yi, i=l,...,d;

x > y iff x 2.y and x #’y;
loooopdo

x >> y lff Xi > yi. 1

For any subset S of I51, we denote by int S and bd S
the interior and boundary respectively of S. The affine
d

subspace of II of least dimension containing 8 is

the affine hull of S, denoted aff S. The interior and

 

boundary of S in the subspace topology of aff S are

the relative interior and relative boundary of S, de—

 

 

noted rel int S and rel bd S respectively. As much

of our interest will focus on sequences of vectors rather
than components of those vectors, we for convenience index
members of a sequence of vectors with subscripts. Thus
the symbol xn may represent the nth component of a
vector x or the nth member of a sequence of vectors,
depending on context. When such usage is potentially

ambiguous, we will exercise greater care.

To continue with our model, we assume that we are
given two vector-valued functions of a vector variable,

p,p'-1 : int IRS -* int IRE ,

which are C1 and are inverses of each other. The functions

6

p and p-1 represent respectively the vector of prices

at which a specified vector of quantities is demanded

(the indirect demand function) and the vector of quantities
demanded at a given vector of price levels (the direct

demand function). In the PIES models these functions are

 

predicted by econometric methods and the predictions then
approximated: the actual demand model is time-dependent,
whereas the input fed to the integrating model represents
a cross-section of the demand function with time fixed.
The form chosen by the PIES modellers is the log-linear

form, expressed by the equations

d m..
p. = k. H q.13

, i = 1,...,d (1.0.1)
1 1 3:1 3

where pi and qi are respectively the price of and
demand for the ith commodity. Adopting the notation

d

log x = (log x1.....log xd) for x 6 int Ig_,
we can rewrite equations (1.0.1) as
log p = Ki-M log q (1.0.2)
where K = (log k1,...,log kd) and M is the d)<d
matrix (mij)° The numbers mij represent the elasticity

of the price of good i with respect to demand for good j;
the significance of the log-linear form is that all price

elasticities are constant.

For the supply side of our model, we employ a linear
program. Suppose that there are s activities in which

producers may engage. Denote the levels of those activities

7

by a vector x E 353. These activities include production,
refining, conversion, transportation and storage. We de-
note by A the d><s matrix whose ith column represents
the output of the d commodities when the ith activity

is performed at unit level, so that Ax represents the
output from activities x. We assume that the only limita-
tions on the activities are that they be performed at
nonnegative levels and that their consumption of r given
resources not exceed the availability of these resources.
We denote by b E If: the vector of available resources
and by B the r><s matrix whose ith column represents
the consumption of resources when the ith activity is
performed at unit level. The vector c E 2le will repre-
sent the costs of performing the various activities at unit
levels. We assume that no commodity is both produced and
consumed, and that goods produced in excess of demand may
be disposed of at no cost. The significance of this last
assumption will be discussed later. The supply model is
predicated on the assumption that producers will elect to
meet a specified demand by adopting a linear program to
determine a vector of activities which meets or exceeds
that demand at least cost. Given demand q, the linear

program is

AX 2 q

Bx g_b
(1.0.3)

x 2'0

cx (min)

8
whose value (least cost) we denote by v(q). Note that
q is not required to be nonnegative, although negative
demands have no obvious economic interpretation in this
model. We extend the definition of v to all of JRd by

setting v(q) = +00 whenever the demand q is not pro-

ducible.

The supply model described by Hogan [4 , 5 ] and
Wagner [17] uses equality rather than inequality constraints.
In the case of the resource constraints, this is no problem.
Constraints on resource consumption are by nature inequali-
ties, but can be made equalities by the addition of slack
variables. The so-called "material balance" constraints,
such as requirement that material transported to a depot
exactly equal material transported from that depot, are by
nature equalities but can be written as pairs of inequalities.
In the case of output constraints, the distinction is some-
What more critical. In general, the most cost-efficient
way to produce at least q (allowing free disposal of
excesses) may be strictly more efficient than the least
expensive way to produce exactly q. We will return to the
matter of free disposal later, at which time we will show

that the PIES model tacitly assumes it.

To determine the level of prices at which producers
will meet demand q, we adopt marginal pricing, taking
the supplier's price for a unit of good i (at the current
demand level) to be the marginal cost of producing a unit

of good i. Other pricing mechanisms, such as using the

9

average cost of a unit of good i as the supply price,

have been mentioned in some of the literature on the PIES
model but will not be discussed here. To determine the
marginal costs in our formulation of the supply model,

we must introduce some concepts from convex analysis. These
concepts are treated in detail by Rockafellar [11], whose
notation we adopt where practical. We repeat some of the

definitions with the aim of establishing relevant notation.

1.1 Definition. A function f :]R(1 4 [-w,+w] is

 

convex iff for any x,y E 351 and any t 6 (0,1),
f(tx+(l-t)y) g tf(X)+(l-t)f(y) (1.1.1)

provided the right-hand side makes sense (more precisely,
if it does not involve expressions of the form co--°°). A
convex function is proper iff it never assumes the value

-w and is not identically +w.

Throughout this paper we deal only with proper convex
functions. It is immediate that the function v defined
in (1.0.3), being the value of a linear minimization pro-
gram, is convex; it will be shown later to be proper under

suitable hypotheses.

1.2 Definition. For f :]R(1 4 [-m,+w], convex, the

 

 

effective domain of f is the set

eff dom f = {x :f(x) < +w].

10

1.3 Definition. For f convex on IRd, u 6 Rd is

 

a subgradient of f at a iff f(a) is a finite real

 

number and

f(x) 2 f(a)+u(x—a) for all x 6 Rd. (1.3.1)

The collection of all subgradients of f at a is the

subdifferential of f at a and is denoted by Bf(a).

Clearly, af(a) is a closed convex (possibly empty)
set. It is known that af(a) is not empty if f is proper

and a E rel int (eff dom f).

We now return to the matter of marginal cost pricing
for our supply model. The function v defined by (1.0.3),
while convex and (as will be shown) pr0per, is in general,
not differentiable, the lack of smoothness occurring at
demand levels at which an activity begins or ceases to be
cost-effective. Our replacement for the marginal cost
vector is the subgradient, which provides us, through (1.3.1),
with a lower bound for the change in total cost as demand
changes. Unfortunately, at demand levels at which v is
nondifferentiable the subgradient is not unique, and so in
general we have an entire set of possible supply price

vectors for a given demand.

A static equilibrium between supply and demand is
defined as a vector of goods such that the price at Which
that vector is demanded is one possible price at Which that
vector is supplied. In terms of subdifferentials, this

becomes the following.

ll

*
1.4 Definition. q 6 If1 is an equilibrium for

our models iff

* *
p(q ) e av(q ).

where p is the indirect demand function and v is defined

by (1.0.3).

If our supply model is to be useful from a computational
standpoint, we must be able to compute Bv(q) at least
for those q which are producible, i.e., for q 6 eff dom v.
The following theorems from linear programming provide a

useful characterization of av(q).

1.5 Theorem. Let A be an arbitrary m><n matrix,
c an arbitrary n-vector. Denote by f(z) the value of

the linear program
Ax 2 z, x 2_0, cx (min), (1.5.1)

where the minimum over the empty set is taken to be +m;
then

. . . m
(1) f 18 a convex function on It ,

(ii) f is nondecreasing in the sense that if
w‘g 2
then
f(w) g f(z).
and
(iii) if f(b) is finite then f is subdifferentiable

at b.

12

In the event that f(b) is finite,
(iv) § 6 af(b) iff § is a maximal solution of

the dual problem

yA g_c, y 2.0, yb(max). (1.5.2)

.ggggf: Conclusions (i) and (ii) are immediate. Though
(iii) and (iv) can be derived from the general theory of
convex programming, we provide, for the reader's convenience,
a self-contained proof using only the duality theorems of

linear programming.

Assume f(b) is finite. By the duality theorems,

the problem (1.5.2) is feasible, and hence

Y={y€IR :y20,yAgc)
is nonempty. Note that Y is defined independently of z,
and so the dual to (1.5.1) is feasible for all z. By

the weak duality theorem of linear programming, then,

f(z) > -0. for all z 6 mm.

Since we have assumed that f(b) is finite, by the strong
duality theorem of linear programming the dual problem
(1.5.2) has a maximal solution, and any maximal solution §

of (1.5.2) satisfies
f(b) = 5113. (1.5.3)

We now show that every maximal solution y of (1.5.2)
belongs to af(b). Observe that for any 2, § is feasible

in the problem

l3
yA g_c, y‘Z O, yz(max) (1.5.4)

dual to (1.5.1). By the weak duality theorem, we have that

m

f(z) 2,§z for all z e Ii . (1.5.5)
In view of (1.5.3) and (1.5.5),
f(z) 2 f(b)+y(z-b) for all z e Rm, (1.5.6)

showing that y E Bf(b).

We now show that if § 6 Bf(b) then § is a maximal
solution of (1.5.2). We are given that (1.5.6) holds;
from (ii), it follows that for z g.b,

f(b) 2 f(z) 2 f(b)+y(z-b),
and so

y(z-—b) g_0 for all z g.b.

This shows that § 2 0. Taking z = Ax (x.2 0) in (1.5.6),

we have
f(Ax)-f(b ‘2 §(Ax-b).
Also, from (1.5.1)
f(Ax) g cx,
and so
cx-f(b) 2_y(Ax-b) for all x‘2 0.

Since f(b) is finite, there exists a minimal solution 2

of (1.5.1) when 2 = b, with f(b) = oi. Thus

cx-cx 2_y(Ax-b) for all x‘z O,

14
and so
(c-§A)x 2_c§-§b for all x 2.0. (1.5.7)

If the jth component of c-yA were negative, we could
violate (1.5.7) by taking x to be a sufficiently large
multiple of the jth standard basis vector in If}. Hence

(1.5.7) implies that

c-yA‘2 0,
which, together with the observation that § 2_O made
above, shows that i is feasible in (1.5.2). Taking x = o
in (1.5.7), we see that

Ci g_yb.
The opposite inequality also holds, by the weak duality

theorem, and so

ex = yb,

indicating that § is optimal in (1.5.2).

1.6 Theorem. Let A be an m><n matrix, B an r)(n

 

matrix, c an n-vector. Denote by f(z,w) the value of

the linear program

Ax‘z z, Bx 2_w, cx(min). (1.6.1)

Suppose that f(a,b) is a finite real number. For 2 6 191,

let
9(2) = f(z,b),

The function g is subdifferentiable at a, and for any
11 6 mm, u 6 ag(a) iff there exists a v 6 Br such that

(u.V) 6 af(a.b).

15

Proof: Assume that f(a,b) is finite; by Theorem

(1.5), of(a,b) is nonempty.

If (u,v) e of(a,b), then for all z e 112’“

9(2) = f(z,b) 2 f(a,b) +u(z-a)+v(b-b)
==gm)+um-a),
and so u E og(a). In particular, since of(a,b) is non—

empty, so is og(a), i.e., g is subdifferentiable at a.

Now let u E og(a) be given. In view of Theorem 1a5,
to show (u,v) E af(a,b) for some v it suffices to show
the existence of a v for which (u,v) is a maximal solu-

tion of
yA+ngc, yZO, $20, ya+sb(max), (1.6.2)

the problem dual to (1.6.1) when 2 and w are a and

b respectively. By precisely the same reasoning as in the
proof of Theorem 1.5, g is nondecreasing in the sense of
part (ii) of that theorem, and so u E Bg(a) implies u 2,0.

Since, by the weak duality theorem,
yA+ngc, yZO, $20, ya+sb2f(a,b)

implies (y,s) is Optimal in (1.6.2), it suffices to show

the existence of a vector v such that
uA+ngc, v20, ua+vb2f(a,b) (1.6.3)
Consider the auxiliary linear program

Bx 2_b, x‘Z O, cx-—qu(min) (1.6.4)

l6

and its dual
sB g_c-uA, 5'2 0, sb(max). (1.6.5)

Since f(a,b) is finite, (1.6.1) is feasible when
(z,w) = (a,b), and so (1.6.4) is feasible. Moreover, for
any x feasible in (1.6.4), we have from the definition

of g that
ex 2_g(AX).
This and the fact that u 6 og(a) imply that
cx 2 g(Ax) 2 g(a) +u(Ax-a)
and so
cx-qu 2_g(a)-ua. (1.6.6)

Thus the value of (1.6.4) is bounded below by g(a)-ua.
It follows that (1.6.4) and (1.6.5) both have Optimal solu-
tions, say i and v respectively. From the feasibility

of v in (1.6.5) we have

uA-thpg c, v 2.0,
leaving us only the task of showing that

uai-vb 2 f(a,b) = g(a). (1.6.7)

Since i and v are optimal, by the strong duality theorem

vb = cx-uA§;
applying (1.6.6) with x = i, we have

vb 2 g(a)-—ua

Which is (1.6.7).

17

We now introduce a special class of convex functions

needed in the sequel.

1.7 Definition. A.real-valued function on I53 is

said to be polyhedrally convex iff it is the upper envelope
of a finite collection of real affine functions defined on
all of Ifi. Thus f is polyhedrally convex iff there

exist affine functions fi (i = l,...,m) defined by

d

f (x) = a.x4—b. for all x 6 I! ,
1 1 1
with ai 6 Rd , bi 6 IR , such that
f(x) = max{fi(x) :i = l,...,m}

for all x E 351. Note that this is equivalent to saying

that f is piecewise-affine.

_l.8 Theorem. Let A be an m><n matrix, B an r)<n

 

matrix, c an n-vector, b an r-vector. Denote by g(z)

the value of the linear program
Ax 2_z, Bx 2.b, x‘2 O, cx(min) (1.8.1)

If there exists a E It“ such that g(a) is finite, then

there exists a polyhedrally convex function h on 19“
such that

g(z) = h(z) for all z 6 eff dom 9. (1.8.2)

Proof: Let a be such that g(a) is finite. Since
Ax 2_a, Bx‘z b, x'Z O, cx(min)
has an optimal solution, so must its dual

yAi—sB‘g_c, y 2.0, s 2_O, yai-sb(max). (1.8.3)

18

Thus the problem
yAi-sB g.c, y 2_O, 5‘2 0, yzi-sb(max), (1.8.4)

which is dual to (1.8.1), is feasible for all z, the
optimal solution to (1.8.3) in particular being feasible
in (1.8.4), and so the value of (1.8.1) is bounded below
for all z, i.e., g(z) > —w for all z. This means 9
is finite on eff dom g, and so it suffices to find h

polyhedrally convex such that
g(z) = h(z) whenever g(z) is finite.
Let F be the set of feasible solutions to (1.8.4), i.e.,
F = [(y.S) =Y 2,0: 5 2.0: ij'SB.S C}:

Which we have shown to be nonempty. Let 2 be any point
Where 9 is finite. Since F is polyhedrally convex and
(by virtue of the nonnegativity constraints) line-free,
problem (1.8.4) has Optimal solutions, at least one of
which must be an extreme point of F. F has finitely many
extreme points, which we may enumerate as {(ui,vi):

i=1,...,T]. Let
h (w) = u.w4-v.b, i = 1,...,T, w 6 IE“,
1 1 l
and let
h(w) = max{hi(w) :i = 1.....T), w 6 mm.

Since (1.8.4) is Optimized at one Of the extreme points Of

F, by the strong duality theorem
g(z) = max{uiz+vib : i = l,...,T} = h(z)

whenever g(z) is finite, which is the desired result.

19

We remark that equality (1.8.2) does not extend beyond
eff dom 9. Also, in the course of the proof we showed 9

to be prOper if it is ever finite.

1.9 Corollagy. The restriction of g to eff dom g

is continuous.

Proof: This follows from (1.8.2) and the observation that

h is continuous on all Of If“.

1.10 Proposition. Let g be as in Theorem 1.8, with

G = eff dom g. The set G is polyhedrally convex and closed.

Proof: Let
X: {xean:Bx2b,x20}.

X is polyhedrally convex and

G=A(X)‘]Rin I

so G is polyhedrally convex, and consequently closed.

We have now Obtained the desired description of Bv(g),
v the supply cost function defined by (1.0.3). Let
g 6 eff dom v. Replacing a,B,b and g in Theorem 1.6
with q,-B,-b and v respectively, we find that Bv(q)
is nonempty and u E av(q) iff there exists 5 E If:

such that
uqz 0, 5‘2 0, uA-sB g_c, uq-—sb = v(q),

i.e., u E av(q) iff (u,s) is an Optimal solution to

the problem dual to (1.0.3) for some 8. Also, by

20

Theorem 1.8, v is polyhedrally convex on eff dom v,
which by PrOposition 1.10 is polyhedrally convex and closed.

We thus have a procedure for computing the subgradients of v.

The time has come to face the issue of free disposal
of excess goods. In view of definition 1.4 and the hypothesis
that the indirect demand function p maps int Rf into
int I§?, we must have some reasonable hOpe that av(q)
contains at least one nonnegative vector for many (preferably
all) q in eff dom v. In our formulation, Theorem 1.6
guarantees that

d
av(q) C Ig..

In the absence of free disposal, however, the output con-
straints Of (1.0.3) become equalities. The corresponding
change in Theorem 1.6 would be to write Ax = z in (1.6.1)
and drOp the restriction Y.2 O in (1.6.2). In this event,
we no longer have restrictions on the signs of the marginal
prices. Economic arguments can be made for the nonnegativity
of the marginal costs, but close scrutiny shows that such
arguments require free disposal of excesses if they are to
be valid. We thus assume free disposal in order to ensure

nonnegative prices.

In the next chapter, we exhibit conditions under

which an equilibrium exists.

CHAPTER II

EXISTENCE OF EQ UI LI BRI UM

We turn to the task of establishing the existence

of an equilibrium demand, i.e., a vector q satisfying

under conditions sufficiently broad to include the supply
and demand models in PIES. We accomplish this in two
stages. Recall that the log-linear indirect demand func-
tion p of the PIES model, though continuous on

int IRE , is undefined on bd IRE . We begin with the

case in which p is continuous on all of RE .

Our original statement of the first theorem used
the conjugate function for v [11] to arrive at a formu-
lation under which the Kakutani fixed-point theorem could
be invoked. We present instead a version which combines

a more general statement with a more direct proof.

2.1 Lemma. Let f be a convex function on IRd;

 

then the set of minimizers of f is a convex (possibly

empty) subset of Ed .

21

22

Proof: Let
_ - d
m — 1nf{f(x) :x 6 1R 3.
The set Of minimizers of f is precisely the set
d
{x 6 1R : f(x) gm}.

Due to the convexity of f this set is convex.

2.2 Theorem. Let K C 351 be compact, convex and

nonempty, f : K 4 Rd continuous and g : K 4 IR con-

 

tinuous and convex. Define g(x) = +w for x E K. Then

there exists a point a 6 K such that

f(a) E Bg(a).

Proof: For each k E K, define fk : K 4 JR by

fk(x) = g(x)-f(k)x for all x E K.

Each f is continuous and convex on K. Since K is

k

compact, the set F(k) Of minimizers of fk over K is

nonempty; that is, if

= min f (x)
mk xEK k
then
F(k) = film) a! 91.

In view of Lemma 2.1, F(k) is also convex. Since fk is
continuous, F(k) must be closed. Thus for each k 6 K,
F(k) is a nonempty, compact, convex subset of K. We

now show that the point-tO-set map F satisfies the

hypotheses of Kakutani's fixed-point theorem [£3]. We

23

need only verify that F is a closed map, i.e., that if

k. E K, k. 4 k, s. 6 F(k.) and s. 4 s then s 6 F(k).
1 1 1 1 1

Since each si is a minimizer Of fk , we have
i

fk (x) 2 fk.(si) for all x 6 K,
1 1
g(x)-f(ki)x‘2 g(si)--f(ki)si for all x E K.
As i 4 m, we have by the continuity of f and g that
g(x)-f(k)x 2 g(s)-f(k)s for all x 6 K,
proving that s is a minimizer of fk’ In other words,

5 E F(k).

Hence by Kakutani's theorem, F has a fixed—point, so that

there exists a E K such that
a 6 F(a).
This says that
g(x)-f(a)x 2 g(a)-f(a)a for all x 6 K,
or equivalently
g(x) 2 g(a)4—f(a)(x-a) for all x E K.

Since 9 is infinite outside K, this characterizes f(a)

as a subgradient of g at a.

We introduce more assumptions regarding the supply
model, so that the hypotheses of the PIES models are met.
We then apply Theorem 2.2 to a sequence of approximations

to arrive at Theorem 2.6.

24

2.3 Definition. We define the set Q to be

 

eff dom v f) IRE

Q

[q20:qux for some x EX},
where X is given by
X = {x.2 0: Ex g b],

B and b as in (1.0.3). Q is the set of demand vectors

of interest to us from the economic standpoint. We have
shown as a consequence of Proposition 1.11 that eff dom v

is closed, and so Q is closed. We assume that X is
bounded, which is in accord with its economic interpretation.
It follows that Q is also bounded, and so is compact.

We further assume that Q is nonempty, and more specifically

that the linear program (1.0.3) is feasible for some q 22 0.

Under these assumptions, take K to be Q, g to
be v, and f to be p in Theorem 2.2. The theorem then
tells us that an equilibrium exists if p is continuous
on Q. As noted earlier, this condition is not satisfied
by the function p in the PIES model. In this case, p

is undefined on part of bd Q, namely Q (1 bd IRE .

The following lemmas will be used in the course of

the proof of the next theorem.

2.4 Lemma. Let f be a convex function on 151, x

a point at Which f is finite; then

2 6 af(x) iff f’(x;y) 2 zy for all y 6 Rd,

25

where

f’(x;y) = lim t-l{f(x4-ty)-f(x)}.
1110

Proof: See Rockafellar [11].

2.5 Lemma. Let v be as in (1.0.3) and Q as in
Definition 2.3; then as q ranges over Q. av(q) assumes

at most a finite number of subsets of I51.

‘nggg: As a consequence of Theorems 1.5 and 1.6,
av(q) is the set Of all y 6 Iii for which there exists
a vector 2 6 Hi: such that (y,z) is Optimal in the
dual to (1.0.3); that is, av(q) is the projection into
I51 of the set of Optimal solutions to the dual. The
set of optimal solutions to the dual is a face of the

polyhedrally convex set

{(y,z) :y20, z 20, yA-ngc}.

Polyhedrally convex sets have finitely many faces, proving

the lemma.

In the following theorem we use the subscripts i,j
and k to denote the components of a vector and the sub-

scripts m and n to denote the members Of a sequence.

2:6 Theorem. Let v be as in (1.0.3), 0 as in
Definition 2.3, and p continuous with

. d . d
p.1nt 1R+ 4 1nt1R+.

Assume that Q is compact and that p satisfies the

26

following conditions for i = l,...,d and q E Q:
if pi(q) 4 +co then qi 4 0. (2.6.1)

where pi and q1 are the ith components of p and

q respectively; and
if qi 4 0 then some pj(q) 4 +w. (2.6.2)

Under these assumptions there exists an equilibrium demand

q E Q, i.e., a vector q such that

p(q) E av(q).

Proof: Let e 6 3C1 be the vector (l,1,...,l).
_ . . d . d
For n — 1,2,... define pn..na_-4 1nt Ig_ by

pn(q) = p<q+ n'le).

Since each pn is continuous on It: and hence on O.
we have by Theorem 2.2 that for each n there exists a

qn E 0 such that

-1 _
p(qni-n e) — pn(qn) 6 av(qn).

Since 0 is compact, the sequence (qn) contains a con-

vergent subsequence, which we denote by Let

(qn)n€N'
q be the limit of this subsequence. We will show that

q )2 0, so that by the continuity Of p at q,
p(qn+n"le) 4 p(q) as n -o co, n e N.
Once this is shown, since
qn * Q. pn(qn) 4 p(q). pn(qn) 6 av(qn).

and the map 2'4 av(z) is closed, we see that

27

p(q) E av(q).

which is what we seek.

To show that q >> 0, it suffices to show that

the set

is empty. Suppose it is not. For i E I,

1 -1

e)i = qn,i + n 4 O as n 4 m, n 6 N.

(qni-n
In view of (2.6.2),

1

pn'j(qn) = pj(qn+n e) 4 +c=o (2.6.3)

for some j (not necessarily i).

For each k e {l,...,d}\I,

qn'kz‘ 0
and so by (2.6.1)
p (q ) = p (q +n'1e>/' +°°
n,k n k n °

Thus for such k, some subsequence of (pk(qn)) is bounded.

Hence we may extract a subsequence

(qn)n€L

from the subsequence (q ) such that for each

n nEN
k E {l,...,d}‘\I there exists a finite number xk such that

) 4 x as n 4 w, n E L. (2.6.4)

pn,k(qn k

By Lama 2.50

{6V(z) : z 6 Q}

28

is a finite collection of sets, and so we may pass to yet

another subsequence,

such that
ov(qn) = ov(qm) for all m,n 6 M.

We have postulated the existence of a vector y 6 Q such

that y 2) 0; since

qn,i 4 O for i 6 I, n 4 w, n 6 M,

surely we can find a number m 6 M such that

q < yi for all i 6 I and all n 6 M, n'2 m. (2.6.5)

n,i

Since v is convex,
v'(qm:y-qm).g V(y)-V(qm) < ”-

Since

pn(qn) E 6V(qn) Bv(qm) for all n 6 M, n > m,
by Lemma 2.4

pn(qn)(y-qm) g V’(qm:y-qm) < w for all n 6 M, n > m.

Hence
lim sup pn(qn)(y-qm) g V’(qm:y-qm) < w. (2.6.6)
n-Doo
n6M
Now

pn(qn) (y - qm) :21 pn'kmn) (yk - qm'k)

+ izejl pn'i(qn)(yi"q.m'i)r (2.6.7)

29

noting that a sum taken over an empty index set is zero.

For k E I,

Pn.k(qn)‘yk‘qm.k) " xk(Yk’qm,k) as n " °°' n 6 M'

and so the first sum in (2.6.7) converges to a finite limit

as n 4 m. For i 6 I,

yi"qm,i > O by (2.6.5)

and

pn'imn) > 0.

so the second sum in (2.6.7) consists exclusively of positive

terms. By (2.6.3), there exists a j such that
Pn j(qn) 4 +co as n 4 m, n 6 M;
I
in view Of (2.6.1), that j must belong to I, and so

pmmnnyj-qm) * +°° as n * °°'

Thus the second sum in (2.6.7) diverges to +m as n 4 m,

n 6 M, contradicting (2.6.6). The theorem is proved.

Conditions (2.6.1) and (2.6.2) are not satisfied by
every function p of the form (1.0.1), so they represent
nontrivial restrictions for the PIES model in particular.
Since q is allowed to move only in a bounded set, (2.6.2)
will hold if, for instance, the matrix M in (1.0.2) is
nonpositive with negative diagonal, a condition satisfied
at least by the example published in Nissen and Knapp [Si].

Condition (2.6.1) is an economic consequence if p

3O

realistically represents the behavior of actual consumers,
since violation Of (2.6.1) requires an infinite supply of

money in the economy.

The following example shows that (2.6.1), or some

similar condition, must be satisfied.

2.7 Example. Let d = 2, r = 1, s = 2,

A: (<1) ‘1’). B: (11).b= (1). C: (1.1) and
p(q) - (qilqgl.q;1) for q >> 0. (2.7.1)
Here
Q=£q20=q1+nglh
Observe that the inverse of p(-) is
q(p) = (pilp2.pgl):

if we let pl 4 +co and maintain p2 = p11, we see that

ql 4 l, violating (2.6.1). In searching for an equili-
brium, we need only consider q 2) 0, since p is
undefined for q1 = O or q2 = O. For q )2 O, we Obtain

by direct calculation:

{(1.1)} if q1+q2< l. ql >0.
qz > 0;

OV(q) = (2.7.2)
(12 > 0.

Now from ql > O and (2.7.1), we infer that

p1(q) 7‘ p2(q) for q 6 Q. q >> 0.

31

and so from (2.7.2) we see that

P(g) E 5V(Q) for all q E Q. q >> 0.

i.e., no equilibrium exists.

We close the chapter by noting that the use of equality
demand constraints in the PIES version Of (1.0.3), while
presenting difficulties mentioned previously, foreshadows
the following observation: at an equilibrium, the most
economical production program must meet all demands exactly.

We prove this below.

2.8 Proposition. Let v be given by (1.0.3), let

. d . d
p.1nt IR+ 4 1nt 1R+

be continuous, and suppose that there exists a vector q

such that

p(q) E av(q).
Let i be an Optimal solution of (1.0.3) for that q; then

Ax = q.

nggf: Note that for p(q) to be defined, we must
have q 2) 0. Also, the optimal solution i must exist,
since av(q) #’¢ implies v(q) is finite which in turn
implies that the linear program is feasible and bounded.
If the ith demand is satisfied with slack, i.e., if
(Ax)i > qi, by the complementary slackness principle of
linear programming the corresponding component of any optimal

solution to the dual of (1.0.3) must be zero. As was shown

32
in Chapter I, there exists a vector 2 such that

(p(q),z) is Optimal in the dual to (1.0.3), and so
Since the range of p lies by assumption in the interior
of the positive orthant, we must have

Ax=q.

Note that this result is not surprising, as we have
legislated it by requiring that both the argument and

the value of p(-) be strictly positive.

CHAPTER III

UNIQUENESS OF EQUILIBRIA

In this chapter, we address the question of whether
an equilibrium for our models, if one exists, must be
unique. We require a monotonicity condition on p, a
higher-dimensional generalization Of the idea that,
for one commodity, demand is a decreasing function of

price. We make the following definition.

3.1 Definition. Let X C an , f :X 4 Ian . We say

 

f is strictly monotonically decreasing iff

(f(x)-f(y))(x-y) < O for all x,y 6 X, x #y.

We now prove uniqueness Of the equilibrium When p

is strictly decreasing.

3.2 Theorem. If x,y 6 Q, p(x) 6 av(x), p(y) 6 av(y)
and p is strictly monotonically decreasing on Q, then

X = y.

Proof: Since p(x) 6 av(x),

v(y) 2v(x)+p(x)(y—X). (3.2.1)

33

34

Similarly, since p(y) 6 Bv(y),

V(x) 2V(y)+p(y)(x-y). (3.2.2)

Adding (3.2.1) and (3.2.2), we have
0.2 (P(X)-p(y))(y-x).
Since p is strictly monotonically decreasing,
O< (p(x)-p(y))(y-x) if x7‘Y.

and so x must equal y.

Sweeney [15], among others, has observed that under
the assumption that p' is globally negative definite,
the equilibrium, if one exists, is unique. Note that here
and in the sequel we apply the terms negative definite
and negative semidefinite to asymmetric as well as symme-
tric matrices. An arbitrary n><n matrix A is negative

definite iff for all y 6 If),
y # 0 implies yAy < 0.

This is equivalent to defining A to be negative definite
iff AntAt is negative definite, where At is the trans-
pose Of A. The following proposition shows that p’
negative definite implies p is strictly decreasing, so
that Theorem 3.2 applies. The proposition also includes

a partial converse.

3.3 PrOposition. Let X be a convex subset of If)

and let f:X4 Ian be Of class C1 on x. If f’(x)

is negative definite for all x 6 X, then f is strictly

35

monotonically decreasing on X. Conversely, if f is
. . 1
monotonically decrea51ng and C on an Open superset
G of X, then f’(y) is negative semidefinite for all

y 6 G.
Proof: Suppose first that f’ is negative definite
on X and x,y 6 X. Set h = y-x. we have
1 .
f(y)-f(x) =j f (x+th)h dt
0
and so

h(f(y)-f(x)) = fl hf’(x4—th)h dt.
0

The real-valued function
g(t) = hf’(x4—th)h
is continuous and negative on [0,1], and so

1
(y—x)(f<y) -f(x)) =f g(t)dt < 0.
O

proving strict monotonicity of f.

Now suppose that f is strictly decreasing and C1

on an Open superset G Of X, and suppose y 6 G. We

wish to show that for all h 6 If),
hf’(y)h g 0.

Suppose not, i.e., that there exists h 6 If) for which
hf’(y)h > 0.

Since G is Open and f’ is continuous at y, there

exists e > 0 such that

36
hf’(y4-th)h > O for all t 6 [-e,c].

Again using continuity Of f’, ‘we have

6
j hf’(y+th)h dt > o;
0

but

(6 hf’(y+th)h dt h(f(y+ch) -f(y))
o

%((y+eh) —y)(f(y+ch) -f(Y))

< 0.

a contradiction. This completes the proof.

Monotonicity Of p is central to the proof of unique-
ness of the equilibrium. In the PIES model, monotonicity
is a consequence Of the negative definiteness Of p', as
shown above. The authors of the PIES model made negative
definiteness of p’ a standing assumption.[ 5]. we shall
produce arguments to the effect that, for the log—linear
form of p, this assumption is exceedingly restrictive
when enforced globally, and will exhibit an example with
multiple equilibria in which p’ is not globally negative

definite.
We begin with the case d = 2.

3.4 Lemma. For p given by (1.0.1) with d = 2, if

p'(q) is negative definite for all q >2 0, either

m12 = m = 0 (3.4.1)

21

or

37

4-1, m = m 4-1 and

m = m 21 11

12 22

(3.4.2)

m i-m

11 224-1 < 0.

Proof: We remark first that a necessary condition

fer a 2)<2 matrix H to be negative definite is that
t

% det(H4-Ht) = det(§J%E—) > 0. (3.4.3)
We compute p’:
m -l m m m -l
11 12 11 12
k1m11q1 qz klm12ql q2
p’(q) = .
m —1 m m m -l
21 22 21 22
k2m21ql qz k2m22ql q2

Let q = (3,1), 5 > 0; then (3.4.3), with H

p’, becomes

m -m +1

1 l 2 2 ll 21
k1k2(m11m22"2 m12m21) ’ Z k1m12S
- -l
1 2 2 m21 m11
m11+m21—1
after division by s . For this inequality to hold

both as s 4 +w and as s 4 O, we require that either

(3.4.1) or

(1,5) shows that

hold. A similar analysis with q

either (3.4.1) or

must hold. Assuming that (3.4.1) fails to hold, we can

reduce (3.4.4) to

38

O < k k m m

l 2
1 2 11 22 ‘ Z [kl(m22+l)+k2(mll+l)]

_ - 1 _ 2
‘ "k1k2(m11+’“22+1) 4 [k1(m22+1) k2(m11+l)] '
and so p’(q) is negative definite only if

m + m + 1 < 0,

ll 22

completing the proof.

For higher dimensions, global negative definiteness

is even more restrictive.

3.5 Proposition. For p as in (1.0.1), d > 2, if
p’(q) is negative definite for all q >2 0 then for all

i # j, either

lj = mji = 0 (3.5.1)
or
m.. = m..4-1, m.. = m..4—l, and
13 33 31 11
(3.5.2)
m..4-m..+-1 < O.
11 3]

In addition, if (3.5.1) fails, then for all k such that
i 71 k 7‘ :1.

Proof: Assume that p'(q) is negative definite
for all q 2) O. This in particular requires that the
2 x 2 submatrix

.. d..
11 13 (3.5.4)

d.. d..
31 J]

39
Of

p’ = (dm).
formed by the intersection of rows i and j with columns
1 and j in p’ (i # j), must be negative definite,
and so (3.5.1) and (3.5.2) follow from Lemma 3.4, taking
qk = l for i #’k # j. Now suppose (3.5.1) fails to hold
for some i # j, and select h such that i #‘h #’j.
Let qk = 1 for all k #’h, qh = s > O, and examine

the submatrix (3.5.4), which reduces to

I

k.m..s 1h k.m..s
1 11 1 1]

m. m.

k.m..s 3h k.m..s 3h .
3 31 J 33

Condition (3.4.3) for this matrix becomes

m. +m. m

k.k.m..m..s 1h 3h - l [k.m..s lh-tk.m..s 3h]2 > 0,
131133 4 113 331
or equivalently
m. —m. m. -m.
4k.k.m..m.. > k?m.?s 1h 3h-I-2k.k.m. m .-tk?m?.s 3h 1h,
1 j 11 33 1 13 1 j 13 31 j 31

Since one of mij'mji is nonzero, this cannot hold for

all s 6 (O,+m) unless (3.5.3) holds.

In economic terms, mij is the elasticity of price
i with respect to demand j, a measure of consumer re-
action to prices. Proposition 3.5 places exceptionally
tight restrictions on what these elasticities could be.
In practice, the demand model is constructed to fit actual
consumer behavior, and it is quite unlikely that the

conclusions of PrOposition 3.5 would be satisfied, assuming

40

the specified commodities bear some relation to each
other. We are thus led to abandon the assumption that
p'(q) is negative definite for all q 22 O. we might
ask that p’(q) be negative definite at least for

q 6 Q, but if the proofs of Lemma 3.4 and Proposition
3.5 are revised to allow s to go to O (as it certainly
may in practice as q moves in Q) but not to go to
+m, we arrive at more complicated but still highly
restrictive necessary conditions. Hogan [5 ] describes
negative definiteness Of p' as "a weak economic
assumption and an Observed property of the PIES demand
functions over the relevant regions," but does not

explicitly describe the "relevant regions."

Proposition 3.3 indicates that we will be hard
pressed to relax the assumption Of negative definiteness
significantly while retaining monotonicity. Our next
example shows that without monotonicity, we may well

have multiple equilibria.

3.6 Example. Let d = 2, r = 1, s = 2,

A = (i i). B = (1 1). b = (2). c = (1,1) and

l 2 -2 -1

p(q) = (32q1 q; .32ql q2 ) for q >> 0.

We assert that (4,2), (2,4) and (3,3) are equilibrium

values for q. Observe that

p(4,2) (201):

ML!“ (1.2).

41

and

p(3.3) = (32/27,32/27).

The dual to problem (1.0.3) is

Zzl+z2

2 4-22 -w‘g 1
l 2 (3.6.1)

zl'ZZ'W‘Z 0

-w g_l

qlzl-tq222-2w (max)

The triples (2,1,4), (1,2,4) and (32/27,32/27,23/9) are Optimal
in (3.6.1) when q is (4,2), (2,4) and (3,3) respectively,

which, in light of Theorems 1.5 and 1.6, proves the assertion.

CHAPTER IV

THE PIES ALGORITHM

In the remainder of this paper we will assume the
existence of at least one equilibrium for our models, and
pursue methods for locating an equilibrium. The previous
chapter indicated that multiple equilibria are a possi-
bility When p is not strictly monotonic. From the
standpoint Of utilization of such models as tools in
policy making, we should be able to find all equilibria,
but as yet we cannot. Hogan and others sought equilibria
under the assumption that p’(q) is negative definite
for each q. We make the comparable assumption that p
is strictly decreasing, ensuring that the model has a

unique equilibrium.

With an eye toward simplifying future calculations,
and with Proposition 2.8 in hand, we revise (1.0.3) using

equality demand constraints, Obtaining

Ax = q

Bx g_b
(4.0.1)

x 2.0

cx (min).

42

43

'k
We define the set Q , corresponding to Q in Definition
2.3, to be

* d
Q = A(X) n 1R+ 1 (4.0.2)

where X is as in Definition 2.3. Moreover, we will
assume, unless otherwise stated, that p is strictly
monotonically decreasing, so that the equilibrium is

unique.

Lastly, we must face the fact that p(-) is not
defined on all of 0*. We have two options. We could
work on a compact subset of Q* bounded away from the
coordinate hyperplanes, or we could assume p(-) is
defined on all of 0*. The former alternative is plausible,
since equilibria are strictly positive vectors, but we
have no a priori knowledge Of where to truncate Q*, and
in addition, we would be adding constraints to (4.0.1).

We therefore choose the latter alternative.

The PIES algorithm is motivated by the observation
that if p has a potential f, i.e., if there exists
a function f on 0 such that p = Vf, and if all the
other assumptions above are met, then q is an equilibrium

demand if and only if q minimizes v-f.

4.1 Lemma. If f and g are prOper convex functions
and there exists a 6 rel int (eff dom f) D rel int (eff dom g),

then fi-g is convex and

44
a(f+ g) (X) = of(X) + 69(X)

for all x 6 eff dom f n eff dom 9.

Proof: See Rockafellar [11].

We remark that for S,T C I51, the expression S-tT

denotes the vector sum

{s-tt :s 6 S, t 6 T].

4.2 Proposition. If p is strictly monotonically

*
decreasing on Q and f : IR? 4 IR is a potential for p,

i.e., vf(x) = p(x) for all x 6 int R3, then q is

an equilibrium iff q minimizes v-f.

Proof: Assume that the hypotheses hold, and let
F = v-f. Since p is strictly monotonically decreasing,
* *
-f is convex, and so F is convex on Q . Since Q
d

*-
is a convex set and F = +w on IH_\ Q , F is convex on

JR:i . From Definition 1.4, q is an equilibrium iff

p(q) e av(q). (4.2.1)

Since f is finite on If?, F is proper, and so q
minimizes F iff

O 6 oF(q). (4.2.2)

We need only show the equivalence of (4.2.1) and (4.2.2).
We have assumed that there exists a vector y 6 Q such
that y 22 0; without loss of generality y 6 Q* and so
rel int (eff dom v) n rel int (eff dom -f) # ¢. Using

Lemma 4.1,

45
aF(q) = av(q)+a(-f)(q) = av(q) -p(q)

*
for all q 6 Q . Thus 0 6 aF(q) iff P(Q) 6 av(q).

proving the proposition.

Necessary and sufficient conditions for the potential

f Of a Cl function p on R3 to exist are that

api/qu = apj/aqi for all i,j = l,...,d.

Hogan [5 ] notes that these conditions are not met by the
PIES indirect demand function p, and so p is not
integrable independently of path, which implies that p
has no potential. We will nonetheless be guided at times

by Proposition 4.2.

The computationally easiest case to handle is when
each component function pi is a function only of the

corresponding variable qi, i.e.,

p(q) = (91(q1).....qd(qd)).

This occurs precisely when p’(q) is a diagonal matrix
for all q (with negative diagonal entries, since p is
strictly decreasing). Proposition 4.2 applies here, and

in fact the potential f is given by
d qi
f(q) = k + Z) [a qi(t)dt
i=1 i

for an arbitrary constant k and arbitrary a 6 int IRE .
Wagner [17] describes a method in which the problem Of
minimizing F, as in Proposition 4.2, can be approximately
solved by solving a linear program.which is an expansion

of (4.0.1). We will elucidate this now.

46

Let I I be closed intervals in [O,+m) such

1'...'d
that the rectangle I = Id x...xIa_ contains both a and

the unique equilibrium point q. Partition each interval

Ij into a finite number of subintervals with partition

points

w ...< w.

3'0 <...< w.

j.-n < am

such that wj O = aj. There is no need to use the same
I

number of points in each interval or even to use equal

number of points on either side of w, but it simpli-

 

 

3:0'
fies the notation to do so. Let
f . a
gl(wi,J) J — O'ooO'n
gi,j=( > I lzl’ooo'd
_ w. ' = _n'...'_l
L gl( 1'3) 3 J
and
Ailj = wirj'l'l-wigj' J = -n'...'n-l' l = 1I°°Oldo

The following result appears in Wagner's paper.

4.3 Lemma. The Optimal solution (x,y) to the linear

program
Bx glb
n-l -l
(AX)i = ai + J§OY1,j-jZ:J-n yl’j' l - l,...,d
x.2 0 (4.3.1)
0 g'yirj S-Aipj’ l = 1'.."d’ J = -nloooon-l
d n-l

cx — Z) . . . . min

/ i=1 jELn 91,3yi,3 ( )

exists and satisfies the following properties for i = l,...,d:

47

1f yi,m 2 O for some m 2 0, then yi,j = O

for all j < O and yi,j = Ai,j for

j: O'ooo'm—l; (403.2)
1f yi,m 2 O for some m < 0, then Yi,j = O

for all j 2_O and yi,j = Ai,j for

j = m+1,...,-1. (4.3.3)

In addition, the component ”1 corresponding tO demand
constraint i in any Optimal solution to the dual program

satisfies one of the following:

0 < Yi,j < Ai,jo j 2.0. vi = qi,j: (4.3.4)

0 < Yi,j < Ai,j' j < 0. Ii = -qi,j: (4.3.5)

yi,j = Ai,j' Yi,j+1 = 0. i 2.0. (4.3.6)
910' 2'”i 2-9i,j+1’

Yi.j = Ai,j' Yi,j_1 = 0: y < 0. (4.3.7)

'9i.j-1 2 "i 2 ”91.3”

Remarks. The key to the proof is that each component
function gi is monotonically decreasing and positive.
The partition essentially estimates the integral Of each
gi (i.e., the area under the graph Of gi) by a sum Of
rectangles, taking into account the sign reversal when

integrating from right to left. Letting

nil _ “g. _
z. = a. + y. . - y. .,
1 1 j=O 1,3 j=-n 1,3
the sum

n—l _

.2 91 . jyi .j

48

is just a Riemann-type sum for
z.

fa: gi(t)dt.
(4.3.2) and (4.3.3) state that the approximation by rectangles
is consistent with the integral, i.e., that rectangles do
not occur on both sides of ai and that no gaps appear
between rectangles. The proof consists of showing that
when this is not the case, a redistribution of weight among
the variables yi j can be made to reduce the value of

' n-l
the Objective function while preserving the value of Z) y. ..

. 1,3
=-n

(4.3.4)-(4.3.7) are proved by characterizing an Optimal

dual solution as a subgradient of the value Of (4.3.1) and

then observing the effect on the Objective function Of a

perturbation in the most extreme nonzero yi j for each i.
I

4.4 Notation. To connect problem (4.3.1) to the

problem

Ax = z

Bx g_b
(4.4.1)

x 2.0

d z.
v(z)-f(z) = cx-f(a) - Z} Ial gi(t)dt (min),
i=1 i

we establish a correspondence between vectors (x,y) feasible
in (4.3.1) and vectors (x,z) feasible in (4.4.1) with

z 6 I, via the following:

“’1 '21
Z. = a. + Z y. o - y. 0. (4.4.2)
1 1 j=0 1,3 j=-n 1,]

Each y clearly produces a z feasible for (4.4.1); each

2 produces a y feasible for (4.3.1) if we also require

49
that y satisfy (4.3.2) and (4.3.3) as well as the

constraints 0 g,yi j S'Ai j for all i,j.

4.5 Lemma. Let (x,y) be feasible in (4.3.1) and
let z be given by (4.4.2); then (x,z) is feasible in
(4.4.1) and the values of the Objective functions of
(4.3.1) and (4.4.1), evaluated at those x,y and z, differ

by no more than

nt4m

[max Ai,j][gi(wi,—n)'-gi(wi, )1. (4.5.1)

113 “

Proof: Feasibility of (x,z) in (4.4.1) is clear.

Since each gi is strictly decreasing,

w. .+y. .
1:] 1:]
91(Wi.j)yi.jwafL j gimdt 2 91(wi.j+1)yi.j

for i = l,...,d and j 2,0, and so
wi .+yi .
(Din/3. '3 gi(t)dt-gi(wi
i,j

.j+1’Yi.j

A similar estimate holds for j < 0. If we add these

inequalities as j runs from -n to n-l and over-

estimate the right side by replacing y. . with max A. ,
1:3 k 10k

the right side teleSCOpes, and after summing over i we

arrive at (4.5.1).

We remark in (4.5.1) that the value Of gj at the
endpoints Of Ij are independent of the mesh size, and
so once the rectangle I is established, the accuracy
of the value Of (4.3.1) in estimating the value Of (4.4.1)
depends on the mesh sizes max A. . (i = l,...,d) with

1,3
3
(4.5.1) providing an explicit estimate.

50

A more critical issue is the accuracy Of the Optimal
§ Of (4.3.1) in estimating the Optimal z Of (4.4.1).
The following proposition shows that arbitrary accuracy
between § and 2 can be Obtained by taking a sufficiently
fine partition, but the crucial estimate relies on a

number which we cannot calculate in practice.

4.6 Proposition. Let (x,y) be Optimal in (4.3.1),
let 2 correspond to § as in (4.4.2), and let (§,E) be
optimal in (4.4.1). Given 6 > 0, there exists a 6 2 0

such that lE-§l < e if max Ai j < 6.
i,j '

Proof: We denote by (M) the problem Obtained from

(4.4.1) by adding the constraint
(2.212 e.

where e is assumed to be small enough such that (M) has
a feasible solution. Let the minimal values of (4.4.1)

and (M) be C and D respectively. Let G = D-C;

o > 0 since all vectors feasible in (M) are feasible

in (4.4.1) and any Optimal solution (x,z) of (4.4.1)

must have 2 = E and hence cannot be feasible in (M).

According to Lemma 4.5, there is a 6 2 0 such that
if max Ai . < 6, then for every (x,y) feasible in

i,j '3
(4.3.1), the value produced by (x,y) in the Objective
function of (4.3.1) and the value produced by the correspond-

ing (x,z) in the Objective function of (4.4.1) differ

by at most 0/3. It follows that the Optimal value E of

51
(4.3.1) differs from the Optimal value C of (4.4.1) by

at most 0/3, and in particular

C-E2-o/3. (4.6.1)

Now suppose we choose the mesh size to be less than
5 and still find that [24%) 26; then (32.2) is
feasible in (M). Problems (M) and (4.4.1) have the same
objective function. Let F be the value of that function

at (§,2). Since D is the optimal value of (M),
F 2 D = C+o. (4.6.2)

From our choice of 6 and the observation that (§,§)

produces the value E in (4.3.1),

E-F2-o/3. (4.6.3)
Adding (4.6.1)-(4.6.3).

C 2_C+-o/3,

which contradicts o 2 0.

Note that the choice of 6 depends, through 0, on
D, which is unknown in practice. Proposition 4.6 guarantees
however, that if we were tO repeatedly solve (4.3.1) with

mesh sizes tending to O, 2 would converge to 3.

Before moving to the more general case, we note that
if the partitioned rectangle I fails to contain the
Optimal solution of (4.4.1), the vector 2 corresponding
to the y Optimal in (4.3.1) will lie on bd I. In this

eventuality, we may partition a new rectangle about 2

and repeat the process.

52
For the more general case, in which p’ is not
even symmetric, let alone diagonal, the PIES algorithm
attempts to exploit problem (4.3.1) by approximating p(-)
by a function with diagonal derivative matrix. We now

present the algorithm.

4.7 PIES Algorithm. Start with a vector pO 6 R3.

 

k = Or e E [001]-

Step 1: Calculate q = p-1(pk).

Step 2: Define an approximation g(';q) to p(°)

by the following rule:

9i‘w’q) = pi(ql"'"qi-l'wi'qi+l'°"'qd)

for i = l,...,d.

Step 3: Solve approximately problem (4.4.1), taking
gi(°) to be gi(-;q), by solving (4.3.1), Obtaining a

solution (x,z), and set Wk = g(z;q)°

Step 4: If Wk = pk (to within predetermined
tolerances), stOp with (approximate) equilibrium demand q.

Otherwise, set

pk+l = 9Pk+ (1‘9” '

increment k by 1, and repeat from Step 1.

Both Hogan [5] and Wagner [17] report that the algorithm
typically converges within ten, and Often within six,
iterations. Wagner states that termination occurs when

the maximum deviation (presumably in p) is within 1%.

53

He further states, and HOgan concurs, that taking 9 =1/2
rather than 0 accelerates convergence. Despite these
findings, the algorithm has to date evaded complete
analysis, due in part to the fact that the function defined
in Step 2 does not approximate p(°) in a fashion that

we find tractable. The results we have on this algorithm
are fragmentary. We first introduce another concept

from convex analysis.

*
4.8 Definition. The conjugate v of a convex function
v : IRd4 [—co,+oo] is given by

* d d
v (x) = sup[xy-v(y) :y 6 It } for all x 6 ll .
We note that since v(y) = +co for y 6 eff dom v,

w
v (x) = sup[xy-v(y) :y 6 eff dom v}.

When v is the value of (1.0.3) or (4.0.1), we have

the following result.

*
4.9 Lemma. v is a proper convex function and

*‘k
v = (v ) .

Proof: Since v is both lower semicontinuous and

prOper, the lemma follows from a result in Rockafellar [11].

4.10 Lemma. y 6 av(x) iff x 6 av*(y).

 

Proof: See Rockafellar [11].

We note that ov*(p) is precisely the set Of demands

q for which p is a possible supply price.

54

In all subsequent analysis we will assume that in
Step 3 of the PIES algorithm, problem (4.4.1) is solved
exactly. We denote its Optimal solution 2 by z = Z(q)

to show the dependence on q.
:k * . '
4.11 Lemma. Z :Q 4 Q Is continuous.

Proof: We note first that gi(-;°) is jointly
continuous for each i. Let g = (gl....,gd). Since
g’(';q) is negative definite, there is a unique 2

Optimal in (4.4.1), and so
z = Z(q) iff g(z;q) 6 ov(z).

Let q 6 0* and suppose (qn) C 0* such that qn 4 q.
Let yn = Z(qn) 6 0*. Since 0* is compact, (yn)
clusters at some y 6 0*. Passing to a subsequence, we
may assume yn 4 y. By the closedness of the point-to-
set map av(-), we have g(y;q) 6 av(y) and so y = Z(q).
This holds for any cluster point of the original sequence

*
(yn), so by compactness of Q , Z(qn) 4 Z(q), i.e.,

Z is continuous.

4.12 Prgposition. Let g be an equilibrium demand,
_ _ w -
p = p(q). If v is differentiable at p, then the PIES
algorithm converges to equilibrium in one step when

started in a suitable neighborhood Of 5.

Proof: We begin by noting that since ov(q) is
closed for all q and only finitely many sets av(q)

exist [Lemma 2.5], there exists a neighborhood M of

55

5 such that for all q 6 0*.
5 (Z ov(q) implies av(q) n M = (5. (4.12.1)

Since g(-;-) is jointly continuous and g(g;g) = p(g)

= 5 6 M, there exists a neighborhood N Of a such
that y,z 6 N implies g(y;z) 6 M. Since Z is con-
tinuous and Z(g) = 6 [because g(g;g) = 5 6 av(g)], we
can find a neighborhood V of a such that V C N and
Z(V fl 0*) C N. Suppose that we start the algorithm with
pO 6 p(V), noting that p(V) is a neighborhood of B.
Let q0 = p-1(po). y = Z(qo). qO 6 V, so qO and y
belong to N and thus g(y;qo) 6 M. Since g(y;qo) 6 ov(y),
av(y) D M.# ¢, and so by (4.12.1) E 6 av(y), which
implies by Lemma 4.10 that y 6 av*(§). Since v* is
differentiable at 5 and g 6 av*(§). av*(§) = {5].

Thus y = g.

'1'. _. _
We remark without proof that av (p) = {q} when q

is an extreme point of Q and p 6 int av(g).

 

*
4.13 Lemma. For x,y 6 av (z),
v(x)-v(y) = z(x-y).
*
Proof: Using Lemma 4.10, x and y in 6v (2)
implies that z 6 av(x) and z 6 av(y), and so

v(y) 2 v(x)4—z(y-—x)
and

v(x) 2_v(y)-+z(x-y).

Combining these inequalities proves the lemma.

56

*
4.14 Lemma. 6v(x) = {y} for all x 6 int av (y).

 

*
Proof: If int av (y) = ¢, dbere is nothing to prove.
*
Suppose x 6 int av (y). Let u 6 Iii be arbitrary.
*
There exists t > 0 such that xi-tu 6 av (y), and so

by Lemma 4.13

v(xi—tu)-—v(x) = tyu. (4.14.1)
On the other hand, for any w 6 av(x)

v(xi-tu)-v(x) 2_twu. (4.14.2)

Combining (4.14.1) and (4.14.2),

(y—w)u 2 O for all u 6 Rd,

and so y = w. Since w 6 6v(x) was arbitrary, 6v(x) = {y}.

4.15 Proposition. If g is an equilibrium, p = p(g)
_ * -
and q 6 int av (p), then the PIES algorithm with e < 1
converges to 6 when started in a suitable neighborhood

of p.

_ *.
Proof: Let N = Z 1 (int av (p)), a neighborhood

of 5 since q 6 int 6v*(p), Z(q) = q and Z is
continuous by Lemma 4.11. Since p and p.1 are both
continuous on int Its, they are homeomorphisms there,
and so we can find a ball M about 5 such that
p-1(M) C N. Suppose that at some iteration k Of the
PIES algorithm, pk 6 M; then q = p-1(pk) 6 N and so

* _ _
y = Z(q) 6 int 5V (p), By Lemma 4.14. BV(Y) = (P)-

57

Since
”k = g(y:q) 6 av(y) = {5}.
we must have Wk = 6. Thus

)pk+l-§I = lepk+ (l- 9)7Tk"P)

 

elpk-is

Either pk = 5 (in which case the algorithm terminates)

 

or lpk+l-—p[ < [pk-p‘. In the latter case, pk+l 6 M

and the argument can be repeated, so that by induction

— m —
[pk-HH-p’ = 9 )pk-p) 4 0

We now present an example in which the indirect
demand function, while always possessing a negative
definite derivative matrix, contains a parameter we are
free to set. For large values of the parameter, the
PIES algorithm has failed to converge in over 100 itera-
tions, and is believed by us not to converge at all.

We shall give a heuristic argument for this and describe

the Observed behavior.

4.16 Example. Let d = 2, r = l, s = 2, T 2 l,

2 l

A = (1 2). B = (l l), b = (2). C = (1:1): With

p defined by (1.0.2) with

-T-1 -T 3(T4-l)log 2

M'= ( -T -T-1)' K = ((3T4-l)log 2)'

58

We note that the example is constructed to have the
unique equilibrium

q = (402): E) = (201):

*... -
and so that the hypothesis 6v (p) = [q] Of Proposition

4.12 holds.

Our heuristic argument against convergence of the
PIES algorithm on this example when T is large proceeds
as follows. Suppose that the sequence (pk) generated
by the algorithm converges to E, and so qk = p-1(pk)"&.
Recall from the proof of PrOposition 4.12 that when
ov*(§) = [a] and q is sufficiently near 5, Z(q) = a.

Thus for k large, Z(qk) = q and so

”k = 9(Z(qk):qk) = g(aqu).

Let
_ -T-l O _ O -T .
Md ‘ ) c) -T-1)' Mo ‘ (-T o)'
then
log wk = K-I-Md log q-I-MO log qk
and
log 5 = K+Md log g+MO log a,
and so

log Wk-log 5 = M0(log qk-log g).

Now M is invertible, so from (1.0.2)

log q = M-1 log p-M-1K

59
and hence
_ _ -1 -
log qk-log q — M (log pk-loq p).
Therefore

_ _ -1 _ -
log nk-log p - MOM (log pk log p).

. ' -1 T
The eigenvalues of MOM are -T and 2T+-1 .

write log pk-log p as a linear combination of the

If we

eigenvectors of M M-l, we arrive at log wk-log 5

O

by stretching one component by a factor of -T and

shrinking the other by a factor of '§é%fi, which approaches
1/2 as T grows large. If a were 0 (so that pk+l = wk),
this would prevent log pk-log E from converging to 0

unless by some chance log pk-1og 5 were an eigenvector

1

of M M- The situation when

._;£__
0 2T+1'

e 2 O is less clear; we present this only to motivate the

with eigenvalue

choice of the particular example.

A modification of the PIES algorithm was tried on
the example, using a program written and executed on a
Tektronix 4051 microcomputer, which holds numbers to
fourteen decimal digit accuracy. The modification occurred
in Step 3. Rather than using the linear program (4.3.1),
we solved (4.4.1) directly, using the method of steepest
descents. It was expected that, if anything, this would
be more accurate than the use of the linear program

approximation.

60
When run with e = 0.5 and T either 1 or 2,

the sequence appeared to be converging to 5, even from
fairly poor starts. When run with e = 0.5 and T = 3,
however, convergence did not occur. Table 1 shows the
distance (in the euclidean norm) between pk and 5

at various points during a number Of runs. In the first
run, the program was started close to 5. After a few
iterations, it began to eXhibit oscillatory behavior,
appearing to be approaching two distinct limit points.
To test this hypothesis, the program was restarted at
what appeared to be one of the two limit points. As

eXpected, the algorithm oscillated between

(1.9999995, 1.00000025)

P

and

(2.0000005, 0.99999975)

P

with the norm of the error constant to within machine
tolerances. To test the effect of roundoff errors, the
program was next started exactly at the solution. Ideally,
the sequence generated should be identically 5, but

in practice we would anticipate errors due to the compu-
tation of logarithms and eXponentials and the finite
accuracy of the machine. After ten iterations, the error
was less than 3><lO—13, which lends credence to the
claim that the oscillatory behavior observed earlier is

due to the algorithm and not to machine errors. Yet

another attempt was made with T = 3, this time with a

IE)

6 1
TABLE 1

Iteration

 

O
10
20

3O

10
20
3O
50
100
120

130

10

20

10

20

ML
1.4142 x lo‘6

5 .2628 x 10"7

5.5781 x 10‘7
5.5907 x 10'7
5.5902 x 10"7

5.5902 x 10'7

0

2.6671 x 10‘13

2.0616

6.8113 x 1010

2.6790 X 108

2.6524 x 105

3.6686 x 102

5.9964 x 10'2

1.0114 x 10"7

1.0288 x 10'7

1.0288 x 10"7

1.4142 x 10'6

3.2196 x 10"5

1.2403 x 10’3

1.4142 x 10‘6

3.2593 x 10"2

9.8588 x lo“1

lt-J

10

62
Iteration
4O
50

60

6.7905 x 10
.7

1

6.6315 x 10
1.2958 x 105
0

7.2116 x 10‘

O

9

63
mediocre start. The results were someWhat startling. In
two iterations the error went up ten orders of magnitude,
suggesting divergence. It then fell off for the next
one hundred iterations, suggesting convergence. Even
while apparently converging, the sequence (pk) appeared
to oscillate, but after one hundred twenty iterations it

appeared to be clustering at both
(1.99999990798, 1.00000004601)
and
(2.00000009202, 0.99999995399).

Attempts with 3 < T < 10 indicated neither convergence
nor a limit oscillation as in the case T = 3. Errors
did rise and fall, which is to be expected, since the
sequence (p-1(pk)), contained in the compact set 0*.
must cluster, and so (pk) must also cluster.

We did observe that when T was large, the algorithm
tended to wOrk better with e nearer to 1. In fact,
an attempt with T = 10, 6 = 0.999 and initial error
1.4142 x 10‘6 was still within 1.4254 x 10'6 of the correct
solution after ten iterations. This would tend to confirm

our heuristic reasoning, i.e., that Wk is further from

6 than pk is.

The results of our tests indicate that even for fairly
small values Of T, the PIES algorithm fails to converge

in a reasonable number of iterations (certainly not within

64
six to ten iterations, the published figures). Moreover,
when T is moderately large, the algorithm is numerically
unstable. The last entry in Table 1 shows that in a
benchmark run (9 = 0.5, T = 10, initial error 0), round-

off errors of 7x10-9

accumulated in just seven iterations.
Values of T above 10 are in no way unnatural; the data
published by Nissen and Knapp [59] include elasticities

as large as 25 in absolute value.

While we cannot prove analytically that the PIES
algorithm fails on Example 4.16, all computational evidence

indicates that it does.

CHAPTER V

THE PIES -VAR ALGORITHM

Irwin [6 ] has proposed a variant of the PIES
algorithm, which he has named PIES-VAR, that can be proved
to converge to the unique equilibrium when several hypotheses
are met. Among these hypotheses is the existence of a
"selection" of subgradients, a function f :Q* 4 151
such that f(q) 6 av(q) for all q E 0*, with
differentiable component functions. The following

result shows that for v defined by (1.0.3), this

hypothesis is generally not met.

5.1 Proposition. Let K C 151 be a convex set,

 

F : K 4 JR a convex function, x 6 int K and f : K 4 Rd

continuous at x. If
f(y) 6 BF(y) for all y E K , (5.1.1)

then F is differentiable at x and vF(x) = f(x).

Proof: It suffices to show that

1imly-xlu1 IF(y)-F(x)-f(x)(y-x)l = o . (5.1.2)
y4x

NOW (5.1.1) implies that for y 6 K,

F(y) -F(x) 2 f(X)(y-X)

65

66
and

F(X)*-F(Y) 2.f(Y)(X"Y)

and so

0g F(y)-F(x)-f(X)(y-X) _<. [f(y)-f(X)](y-X) .

Therefore

OglFW)-FM)-fm)w-XH ngWO-fhdlly-xl (SJWN

for all y 6 K. Since x 6 int K, (5.1.2) need only be
verified for y 6 K, in which case we have, due to

(5.1.3) and the continuity of f, that

lim [y—xl‘l [F(y) -F(X)-f(x)(y-X)I = lim [f(y)-f(x)| = o.
y4x Y‘X
Taking K = Q* and F = v in Proposition 5.1,
we see that Irwin's hypothesis that the component functions
of f be differentiable would require that v be
differentiable on int Q*, which in general is not the

case.

Despite the failure of our supply model to fit the
hypotheses of Irwin's convergence result, the PIES-VAR
algorithm will still converge When some stringent conditions
are met. This algorithm is similar to, though much more
straight forward than, the PIES algorithm, and so analysis
of PIES—VAR might suggest directions of study for the
PIES algorithm. With this in mind, we present the

PIES-VAR algorithm.

67
5.2 PIES-VAR.A;gorithm. Start with qO 6 Q*, k = 0.

Step 1: Construct g(-:qk) as in Step 2 of

algorithm 4.7.

Step 2: Solve problem (4.4.1) or its approximation

(4.3.1), obtaining the solution (X'qk+l)'

Step 3: If qk+l = qk (to within predetermined
tolerances), stOp with (approximate) equilibrium qk.

Otherwise, increase k by one and repeat from Step 10

From the standpoint of analysis, the major advantages
of PIES-VAR are that there is no shifting back and forth
between p and q, and no averaging. Even with these
advantages, we are able to prove convergence only in

special cases.

5.3 Proposition. If g. is an equilibrium,
5': p(g), and the conjugate v* of v is differentiable
at B, then the PIES-VAR algorithm converges to '5 in

one step when started sufficiently near a.

Proof: This is an immediate consequence of the proof
of PrOposition 4.12: for qO in the neighborhood N of
a; defined in the proof of 4.12, ql = Z(qo) = a and the

algorithm terminates.

The next result suggests that the structure of the
elasticity matrix M in (1.0.2) is relevant to convergence.

We decompose M as M.= Mdi-MO, Where M.d is a diagonal

68

matrix and M has zeros along its diagonal. We do not

0

require that p' be negative definite, but instead require

that Md have strictly negative diagonal entries, which

guarantees both that Md is invertible and that g(°;°)
is monotone in its first argument. This requirement
corresponds to the economic condition of negative own-
elasticities, which is a standing hypothesis of the PIES

demand model.

5.4 Proposition. Let p be as in (1.0.2),

M = Md-I-MO as above, 5' an equilibrium, and 5': p(g).
If 3'6 int av*(5) and Md-lMO has spectral radius

strictly less than unity, then the PIES-VAR sequence

(qk) converges to El when started sufficiently

near q.

Proof: The set int av*(§) is by assumption a

neighborhood Of a, and av(z) = (p) for all

z 6 int av*(E) by Lemma 4.14. Define E : int R34 Ed

by E(x) = log x-log a, where the notation log x is
as in (1.0.2). Let N = int av*(§). We assert that

for any member y = of a PIES-VAR sequence and

qk
any X 6 N,

MdE(X)4-MOE(y) = 0 implies X = (5.4.1)

qk+l '
Since Md log gA-MO log gA-K = log 5 and M.d log x4-
MO lOg y4-K = log g(x;y), MdE(x)4-MOE(y) = 0 implies

that g(x;y) =Ip. Moreover, x 6 N implies that

69

av(x) =‘E, and hence g(x:y) 6 av(x). This last relation
implies, by the definition of Z (preceding Lemma 4.11),

that x = Z(y) = Z(qk) = qk+l° Thus (5.4.1) is established.

Let H = —M - It is well known (cf. Varga [16])

0.
that Hn 4 0 as n 4 4 iff the spectral radius p(H) of

H is less than 1. Since we have assumed p(H) < l,

Hn 4 0 and hence there exists m such that

Hnnu < 1 for all n‘2 m , (5.4.2)
where H'H denotes the spectral norm (of. Varga EMS]
for the definition of the spectral norm). Since E is

a homeomorphism, the set E(N) is an Open neighborhood

d

of O in It . Let 80 be the largest Open ball with

center 0 contained in E(N) and let S1 = H-1(Si__l)
= {s :Hs 6 Si-—l}’ i = l,...,m. Each Si is an Open

neighborhood of 0, and hence so is S = (a Si‘
i=0

Finally, set U = E-1(S), a neighborhood of 5' contained
in N. Suppose that the sequence (qk) is generated by

PIES-VAR with qO 6 U. Let r0 = q0 and

r = E‘1(H*+ln(ro)). k = 0,1,... . (5.4.3)

k+1
We will show that (rk) and (qk) are in actuality

the same sequence. We will prove this by mathematical
induction. To begin, note that in view of (5.4.3) we

have

E(rk+_1) = nkl'l E(ro) = nkl'l E(qo) . (5.4.4)

7O

Observe that

MdE(rk+ l) +MOE(rk) = MdH E(rk) +MOE(rk)
(5.4.5)
= -MOE(rk)4-MOE(rk) = o
for k 2_0. ‘We will now show that rk = qk for k 2.0.
By definition, r0 = qo. we show that if rk = qk then
rk+1 = qk+-1° Taking rk = qk in (5.4.5) we have
Mdmrk+ 1) +MOE(qk) = o ; (5.4.6)

in view of (5.4.1) and (5.4.6) we can conclude that

if we show that 6 N. For k < m,

rk+1 qk+1 rk+1
) = E(qo) 6 S C'Ski-l' Since qo 6 U, and so

= Hk+1 E(ro) 6 HR"-1 S

r e N. For k2m, [E(r
k+1
HH

ki-l C S0 C E(N). Therefore

k+1” =1Hk+1 E(roH g [E(rOH.

H < 1- Since E(ro) 6 50' a ball centered

at 0, E(r also belongs to S and again we find

k4-1)
6 N.

oI

that rki—l

We can now show that qk 4 3’ when qO 6 U. Since
E(rk) = Hk E(ro) and Hk 4 O as k 4 m, E(rk) 4 O, and

so qk = rk 4 q.

The next prOposition guarantees convergence of
PIES-VAR from any start in 0*, when the indirect demand
function satisfies a very strong condition. We will use
the symbol g’(-;-) to denote the derivative of g with
respect to its first argument. We again require monotonicity

only of 9, not of p. For p as in (1.0.2), this reduces

71

to Md being negative definite.

5.5 Lemma. Let g be as in Step 2 of Algorithm

 

4.7, and let g’(x;y) be negative definite for all
x,y 6 Q*. Then for all x 6 0* there is a number

m(x) > 0 such that
b 2
g(a;X)(b-a) -f 9(srx) - ds 2 m(X) lb-al (5.5.1)
a
for all a,b 6 0*.
Proof: We note first that the integral on the
left hand side of (5.5.1) is well-defined, since g’(o;x)

diagonal implies that its integral is independent of

the path. Fix a,x 6 0*. Define G : [0,1] 4 ]R by

t
G(t) = tg(a;x)(b-a)-I g(Tb+ (l-T)a;x)(b-a)d'r. Clearly
O
G(O) = O and G'(O) = g(a;x)(b-—a)-g(a;x)(b-—a) = 00
Moreover, G”(t) = —(b-a)g’(tb4—(l-t)a;x)(b-a) for

0 < t < 1. Using the second Mean Value Theorem,
G(l) = -(b-a)g'(5b+ (l- 5)a:X)(b-a)

for some g 6 (0,1). Since g'(-;x) is negative definite
and continuous and Q* is compact, —g’(-:x) is uniformly
positive definite on 0*, i.e. there exists m(x) > 0

such that -ug'(w;x)u 2 m(x)]u]2 for all w 6 Q* and

all u 6 Rd. Thus

b
g(a;X)(b-a) -( g(S:X) - ds = 6(1)
a

= -(b—a)g’(gb+ (l- §)a;x)(b—a)

‘2 m(x)]b-—a|2 .

72
5.6 Proposition. Suppose that p 6 C1 such

that g(-;-) is monotone in the first argument and

satisfies the following Lipschitz condition:
lg(x;y)-—g(x;z)| g_k ly-—z] for all x,y,z 6 Q*.

Let ‘g be an equilibrium. If 1 < m(g), then the
PIES-VAR sequence converges to g' from any start

in q*.

Proof: Since solves (4.4.1) at iteration

qk+1

q _.
v(qk+ 1) - a“ 1 g(s:qk) 0as g v(q)

noting that the solution to (4.4.1) does not depend on

what we take to be the lower limit of integration.

Consequently
_ q _
v(qk+ l) -v(q) —I_k+ l g(srq) 'ds
q
q (5.6.1)
3 3‘” [g(s;qk) -g(s;§)] °ds.
q
since 5 = 8(3) 6 ME) and 95:6”) = mg).
v(qk+l)-v('g') 29(E:E)(qk+l-E) . (5.6.2)

In view of (5.6.2) and Lemma 5.5, we see that
._ q ._
V(q )--v(q)-f__k+l g(s;q) °ds
ki—l q

_._. _. q ._
2 g(QFQ)(qk+ 1-q)- __k+ l g(Siq) ' as

q

2m(§)|qk+l-§lz . (5.6.3)

73

Now
qk+1 __
J_ [g(Stqk)-g(s;q)] °ds
q
S -S—up lg(z;qk)—g(2;_q-)l iqk+l--q-l . (50604)

Since lg(z;qk)-g(z;g)l g_k lqk-g], we have by
combining (5.6.1), (5.6.3) and (5.6.4) that

ma) lqu-Elzg) )qk-g) lqu-El

and so

lqk+l-'<il smarl x lqk-El .

If 1 < m(g), we have qk 4 g. at least as fast as a

geometric sequence.

Turning back now to PrOposition 5.4, we might

ask what happens if all other hypotheses are satisfied

but Md-l MO has Spectral radius greater than unity.

This question, which leads us to the following example,
also prompted example 4.16 and the corresponding remarks

there.

5.7 Example. Let d = 2, r = l, s = 2,

2 1
. B=(l 1). b=(2). c=(1.1). with
l 2

>
II

p as in (1.0.2) with

K 1.1 log 2-log 3 M -2 0-9
= 1.1 lo 2-1o 3 '
g g -1 —0.1

= (2 . 2). E: (1/3 . 1/3).

.QI

74
This example is constructed with p’ negative definite
throughout 0* (guaranteeing a unique equilibrium,
namely a) and with '5 6 int av*(§) and av(g) = {5].
To verify this last claim, we note that for all q 6 Q*.
the unique solution to (1.0.3) is x(q) = ((2q1-q2)/3,
(-q14-2q2)/3) With V(q) = cx(q) = (qli-q2)/3. Thus

{(l/3,l/3)} = {5). Moreover, for all p 6 R2 we have

EVE)

SUP (pq-V(q) Hi 6 R2} = sup (pq-V(q) :q 6 0*}

v*(p)
= SUP {pq- (ql+q2)/3 :q 6 0*} = sup ((p-Em =q 6 0*}.

0. Since ‘6 6 int Q*.

and in particular v*(p)

V*(p) -v*('§) = V*(p) .>_ (p-B)q

for all p 6 JR2 and any q in a neighborhood of g
in Q*, from which it follows that q 6 Bv*(§) for all

q near 5' and hence 5'6 int av*(§).

5.8 Lemma. In example 5.7, g(§:q) 6 ov(g)

implies that q =‘g.

Proof: Since ov(g) = [E]. g(g;q) 6 Bv(g)

implies that g(grq) = 5. NOW p(x) = (kxl—zxzo’g,
kxl-lxz-O'l) for all x 22 0, where k = 3-121'1.
Hence g(EIQ) = (k2-2q20'9,k2-O'lq1-1). since 3': (2,2).
Therefore g(g;q) = 5' iff qzo'9 = 20'9 and q1-1 = 2-1.

i.e. iff q=‘q'.

5.9 Proposition. If qO #“g, the PIES-VAR

sequence (qk) cannot converge to q.

75
_P_r_o_c_>_f_: By Lemma 5.8, if qk+1 = (i then
g(aqu) 6 OV(E) and so qk = 6. Since qO 316, it
follows that qk 7! g for all k. Now suppose
qk 4 6'6 int av*(p). By Lemma 4.14, av(q) = (P)
for all q 6 int ov*(p). Since gk 4 a, there exists

m > 0 such that qk 6 int av*(§) for all k 2 m.

Using the decomposition M = Md+ MO'

have g(qk+ l:qk) = p and so K+Md log qk+l+M0 log qk =

for k2m we

log 5. Also K+ Md log g-I-MO log 3 = log 5, and so

for k2m

_ _1 _
log qk+ 1- log q = -Md MO[log qk- 109 q] .

It follows that for k 2 m

- __, -1 2 —
log qk+ 2 -log q — (Md MO) [log qk- log q] .
Now
-4.5 O
(Md-l M0)2 =
O -4.5
and so for k2m I log qk+2-logg I = 4.5 I log qk-log El.
From this we deduce that for k = 1,2,...
lloqq -1093) =4s" (1.. -10951
m+2k ' qm '
Since qk 4 q by assumption, I log qm+ 2k- log q l 4 0
as k 4 co, and so we must have q = a, a contradiction.

m

CHAPTER VI

A SUBGRADIENT PROJECTION ALGORITHM

The lack of theoretical justification for convergence
of the PIES algorithm, and the inherent difficulties in
analyzing it, have led a number of researchers, including
Eaves [3] and Irwin [6], to prOpose other algorithms for
locating equilibria in models of the PIES form. The failure
of the PIES algorithm on example 4.16 underscores the need
for a more generally convergent method. In this chapter we
state an algorithm for minimizing a convex function of a
particular type over a polytoPe. The algorithm and a proof
of convergence were prOposed by Sreedharan. When the indirect
demand function p has a concave potential f, the function
v-f is of the specified type, and we will prove that our
algorithm produces a sequence converging to the minimizer
of v-f over 0*. The extremal point located is, by
virtue of prOposition 4.2, an equilibrium. In addition, the
function f is used as a tool to prove convergence but is
never specifically evaluated by the algorithm, leaving Open
the possibility that the algorithm will prove efficacious

on some problems in Which p does not have a potential f.

The problem of minimizing a convex function over

a linearly constrained set has been the subject of much

76

77
study. We single out the approach taken by Rosen [12]
where the objective function is smooth. Rosen attempts
to exploit the known convergence of the method of steepest
descent in the unconstrained case. As is typical of
so-called "feasible direction" methods, Rosen computes
at each iteration a direction of descent Which points
into the constrained region from the current point.
He searches in that direction until he reaches either
a relative minimum along the ray of search or the
boundary of the constrained region. The process then
repeats. Rosen's contribution is the choice of direction.
When possible, he uses the negative gradient as the
direction; When this direction points out of the set,
he projects it onto a face of the set. Rosen's method
is susceptible to a phenomenon known variously as
"jamming" or "zigzagging", in Which the sequence generated
clusters at, or even converges to, nonOptimal points.
The trouble lies in the possibility that the sequence is
alternating among two or more faces of the constrained
region in such a way that the distance along the direction
of search from the current point to the boundary is going
to zero. Various modifications have been prOposed to avoid
this. In particular, Polak [10] has adOpted a technique
which prohibits the sequence generated from approaching
arbitrarily close to a face when conditions for a

constrained minimum are not being met at the limit.

78

Rosen's method, and much of the other work in the
area, requires that the objective function be differentiable.
Even when p has a potential f, our objective function
v-f is not differentiable everywhere because v is not.
Attention has recently been focused on algorithms for
Optimizing nondifferentiable convex functions. The
algorithms of Wolfe [l9] and Lemarechal [ 7], which
generalize classical methods for unconstrained Opti-
mization by replacing the gradient with a carefully
chosen subgradient, do not treat the constrained case.
The algorithm of Bertsekas and Mitter [1.], which
does handle constrained problems, requires computation
of the "c-subdifferential" of the objective function,
Which is prohibitive in the problem we consider here.
The algorithnt proposed here is prompted by those of
Sreedharan [13,14] , Rosen [12] and Polak [10]. It resembles the
Bertsekas - Mitter algorithm but requires the computation
of only a manageable portion of the c-subdifferential.
In this chapter we pose the algorithm and prove its
convergence; in the next chapter, we discuss the actual

implementation and report on some trial applications.

6.1 Problem. Let X C IVE be a nonempty convex

 

polyt0pe given by X = {x 6 Rd] aix g bi' i = l,...,m].

NOte that.polytOpes are, by definition, bounded.

Let v.:]Rd 41R (j=l,...,r) be given by

3
d
v.x= .x+c., .6]R,c.61R andletthe
J() 93 3 93 J .

79

polyhedrally convex function v : Rd 4 R be given by

v(x) = max{vj(x)l j = l,...,r}, x 6 Eta. (6.1.1)

Let f : Rci 4 R be strictly convex and of class C1 on

a neighborhood of X. Finally, denote by 6

the indicator function (cf. Rockafellar [11]) of X,

 

i.e. 6(x) = 0 if x 6 X and 6(x) = +co if x K'X.

We consider the problem

ll
H
§
3

aix g b1' i
. (6.1.2)

f(x)-tv(x) (min)

which is equivalent to the problem of locating an

. . . . d
unconstrained minimizer of F = f4-v4-6 over 11 .

We note before proceeding that the problem of
locating the equilibrium of the PIES model fits the
form of (6.1.2) when the indirect demand function p
has a strictly concave potential, for if we take
X = 0* and take f to be the negative of the potential
of p, then by theorem 1.8 the value v Of the linear
program (4.0.1) is pothedrally convex on Q*, and by

proposition 4.2 q is an equilibrium iff q solves (6.1.2).
We introduce some needed notation.

6.2 Notation. For x 6 X and 6‘2 0, we define

the following sets of indices:

80

I€(x) = (l g,i g_m) aiX-2 bi-e): (6.2.1)
J€(x) = {l ng g_r| vj(x) 2_v(x)-c}. (6.2.2)
Note that
10(x) = {1 g,i g_ml aix = bi} (6.2.3)
and
JO(X) = {l g_j g,r| vj(x) = v(x)}. (6.2.4)

We also define two convex subsets of 151, namely

Cc(x) = cone{ail i e I€(x)} (6.2.5)

v

and

K€(x) = conv{gjl j 6 J€(x)}, (6.2.6)

where for any set S we use cone S and conv S to denote
respectively the convex cone, with apex at the origin,

generated by S and the convex hull of S.

For any nonempty closed convex set S C 151 there is
a unique point x 6 S nearest to the origin, which we
denote by N[S]. The point a = N[S] is characterized

by the following:

a = N[S] iff a(x-a) 2.0 for all x 6 S. (6.2.7)

We now present a subgradient projection algorithm for

problem (6.1.2).

81

6.3 Algorithm. Begin with arbitrary x0 6 X,

 

co 2 0 and with k = 0.

Step 1: Compute y0 = N[vf(xk)i—Ko(xk)-tco(xk)].

If y0 = 0, stop: x solves (6.1.2). If yO #’0,

k

set 8 80.

Step 2: Compute y8 = N[vf(xk)4-K€(xk)-tC€(xk)].

2 _ _
Step 3. If lye) 2 8, set 6k — 8, 5k — y€

and go to step 5.

Step 4: Replace c with 8/2 and go to step 2.

Step 5: Compute Gk = max{o 2 0| xkndsk 6 X].
(It will be shown that 6% 2 0). Find Gk 6 [0.3%] such
that there exists 2k 6 vf(xk-dksk)4-K0(xk-dksk)

satisfying zksk = 0; if no such Gk exists, set
ak=ak°

Step 6: Set xk+1 = xk-dksk, increase k by

l, and go to step 1.

The implementation of steps 1 and 2, which can
be treated as quadratic programs, and of step 5, Which
requires a special line search procedure, are discussed

in the next chapter.

We next state a sequence of lemmas leading to a

proof that algorithm 6.3 solves problem (6.1.2), or

82

equivalently locates an unconstrained minimizer of

F = f4-v4-6. ‘We must first define the "c-subdifferential."

6.4 Definition. Let G : Rn 4 [49,00] be a convex

 

function. The e-subdifferential of G at the point x,

 

denoted a€G(x), is defined by

a€G(x) = (u e n“ [G(y) 2 G(x)+u(y-x) -e

for all y e Rn}.

The usual subdifferential 5G(x) of G at x is just

BOG(x).

We now state a sequence of lemmas, using the earlier

notation.

\

6.5 Lemma. For all e 2_0 and all x 6 351,

 

K€(x) C o€v(x).

Proof: If u 6 K€(x), then by (6.2.6) there exist

Xj-Z 0, j 6 J€(x) such that Z:lj = l and

u = :3 1.9..
j6J€(X) 3 3

For j 6 J€(X), we have

vj(y) = vj(X)4-gj(y-X) 2_v(x)-—c-+gj(y-x). (6.5.1)

Therefore for every j 6 J€(x),

v(y) = max vi(y) 2 v(x)-c-+g.(y-x),
i=1,...,r 3

and so

83

v(y) = 2'3 X.V(y) 2 V(X)-e+ '23 Rely-X)
j6J€(x) 3 j6J€(x) 3 3
(6.5.2)
=V(X)-e+u(y-X).

Since (6.5.2) holds for all y 6 I51, the lemma is proved.

6.6 Lemma. Given x 6 151, there exists a neighborhood

 

V of x such that JO(y) C JO(x) for all y 6 V.

gpppfi: The functions wj = v-vj, j = l,...,r
are continuous with wj(x) 2 0 iff j z JO(X). Thus
there exists a neighborhood V Of x such that wj is
positive throughout V for each j g JO(X). If
j g J0(x) and y 6 V, wj(y) 2 0, and so j z J0(y),

proving the lemma.

6.7 Definition. Let S C 191 be a nonempty set.

 

The support function P of S is defined by

 

Cp(X) = SUP(XYI y 6 S}. x 6 Rn.

6.8 Lemma. Two closed convex subsets of 151 are

 

identical iff their support functions are identical.
Proof: See Rockafellar [11].

6.9 Lemma. 6v(x) = KO(x) for every x 6 X.

Proof: Let x 6 X. Since both av(x) and K0(x)
are closed and convex, it suffices, in view of lemma 6.8,

to ShOW’that they have the same support function, namely

84

I I

v (x;-) [v defined as in lemma 2.4]. It is well-known
that since v is everywhere finite—valued,

v'(x;y) = sup {yul u e av(x)} for all y 6 Rd, (6.9.1)

i.e. v’(x;-) is the support function of av(x). Given

y 6 151, by lemma 6.6 there exists an e 2 0 such that

JO(x+-dy) C J0(x) for all a 6 [0,6]. Since for 0 g_d g_e

v(x4—oy) = max v.(x4-ay)
j6JO(x+ay)
and
V(X) = max V-(X)r
j6J0(x) 3

we see that

v(xi-dy)-v(x) = max [v.(x4—oy)-—v.(x)]
' J
36J0(x)
= _max dgjy.
36J0(x)
This shows that
v’(x;y) = rmnc g.y = max [uyl u 6 Ko(x)}, (6.9.2)

- J
36J0(X)
and so v'(x;-) is also the support function of KO(X).
We note that (6.9.2) proves the following statement:

6.10 Corollary. v’(x;s) = max {sul u 6 K0(x)}.

6.11 Lemma. For each x 6 X and e 2 0 there

exists a y > 0 such that

JO(x) C J€(y) whenever lx-y| < y.

85

Proof. Choose y 2 0 such that lgjly < %- for
j = l,...,r and [v(x)-v(y)] <-% if Ix-yl < y. Now

if j6 J0(x) and [x-yl < y, then

V(y) - vj(y) = My) -V(X) +vj(x) -Vj(y)

(X‘Y)l < 60

<-§-+lqJ

and so j 6 J€(y).

6.12 Lemma. oF(x) = vf(x)4-Ko(x)-tC (x) for all x 6 X.

0

 

ugpppg: The indicator function 6 is clearly prOper
and convex, while f and v are everywhere finite valued.
It is well-known that for x 6 X, 66(x) = CO(x). Moreover,
any a 6 rel int X belongs to rel int (eff dom f) n

rel int (eff dom v) n rel int (eff dom 6). The result now

follows from lemma 4.1.

The next lemma shows that the stOpping criterion

in step 1 of the algorithm is well chosen.

6.13 Lemma. If y0 = O in step 1 of algorithm

 

6.3, then xk is the minimizer of F.

Proof: y0 = 0 implies that 0 6 aF(xk), a
necessary and sufficient condition for xk to minimize F.

The strict convexity of f ensures that the minimizer of

F is unique.

6.14 Lemma. Step 4 of algorithm 6.3 is not executed

 

infinitehy often in any one iteration.

86

Proof: If step 4 is executed infinitely often,

then 3 4 O and ye 4 0. Now

61 2 62.2 0 implies that K€ (xk) C K€ (x and

)
2 l k

C€2(Xk) C C€1(Xk)

and hence that
[Y lily l.
61 ‘32
so that y€ 4 O as c 4 0 implies that y6 = 0 for

every 6.2 0; but then y0 = O, and we cannot have

reached step 4, a contradiction.

We now show the practicability of step 5 of the algorithm.

 

6.15 Lemma. If sk # 0, then -sk is a feasible

direction of strict descent at the point xk.

Proof: From the definition of s in step 3 of

k
the algorithm,

sk = N[vf(xk) + K€k(xk) +C€k(xk)] .

Let 1 6 10(xk) C I€ (Xk); then ai 6 C (x

k 6k

k) and so

5 -+ai 6 vf(xk)+K€

k (xk)+C€

(xk).

k k

using the fact that Ce (Xk) is a convex cone. By
k

(6.2.7) we have sk(sk4-ai-sk)‘2 0. Thus aisk 2 0

for every 1 6 Io(xk). Since aixk < bi for i z 10(xk),

87

there eXists d 2 0 such that ai(xk-dsk) S-bi for
all i = l,...,r. Hence -sk is a feasible direction at

Xk.

To show that -s

k is a direction of strict descent,

we show that

F'(xk;-sk) = :13 [F(xk-dsk)-F(xk)]/o < 0. (6.15.1)

From the first part of the proof, there exists ‘5 2 0 such
6 X for 0 g_d g/d. For a in this range,

that x -as

k k
F(Xk-Gsk) = f(xk-csk)4-v(xk-dsk) and so by corollary

6.10

F’(Xk:-sk) = f’(X—k;—sk) +V’(xk:-sk)

-vf(xk)sk4-max{-skyl y 6 KO(xk)}

= -min((vf(xk)4-y)skl y 6 KO(Xk)]. (6.15.2)
When y 6 Ko(xk) C Kc (Xk) we have
k

vf(xk) + y 6 Vf(Xk) + K€k(xk) +C€k(xk)

and so by (6.2.7) sk(vf(xk)4-y-sk).2 0 and consequently
(vf(xk)+y)sk 2|sk|2 > 0. Combining this with (6.15.2),
we have

F’(Xk7_sk)-§ -lsk|2 < 0, (6.15.3)

completing the proof.

From the first half of lemma 6.15 we have the

following corollary.

88

6.16 Corollary. The number a; defined in step 5

of algorithm 6.3 is positive.

The next lemma shows that in the relevant case the

vector 2k in step 5 of the algorithm exists.

6.17 Lemma. Let sk # 0 and define m on (0,6k]
by m(d) = F(Xk-Gsk). If 5% is not a minimizer of T

on [0, k], then 2k satisfying step 5 of algorithm 6.3

exists.

ggppfi: By lemma 6.15, m'(0) = F’(xk;-sk) < 0,
so that there is some a 6 (0,5k] such that p(q) < m(O).
Since we have hypothesized that 6k does not minimize cp,
there exists Gk 6 (0,5k) minimizing m over [0,ak].
Set y = xk-dksk. There exists 8 2 0 such that

F(y) g F(y+lsk) for [11 g c. It follows that

[F(y+).sk)-F(y)]/). 20 (6.17.1)

and

[F(y-lsk)—F(y)]A 20. (6.17.2)

0 < l g.c. Since F is convex, the directional derivatives
F’(y;sk) and F’(y;-sk) both exist, and from (6.17.1) and
(6.17.2) we conclude that F'(y;sk) 2 0 and F'(y;—sk) 2_0.

Using corollary 6.10,

F'(y;:tsk) f’(y;i:sk)+v'(y;:ts

k)

ivf(y)sk+max{ :tuskl u 6 KO(y)].

89

Since K0(y) is compact, there exist u,w 6 Ko(y) such

that
._ ’ .
Vf(y)sk-I-usk — F (y,sk) 2 O
and
__’ .-
vf(y)sk+wsk — F (y, sk) g_0,

and so for an apprOpriately chosen convex combination h of

u and w we have h 6 K0(y) and

vf(y)sk'+hsk = 0.

Taking zk = vf(y)-th 6 vf(y)+-Ko(y) satisfies the requirement

in step 5 Of the algorithm.

The number 0k determined in step 5 of algorithm

6.3 has the following prOperty.

6.18 Lemma. Let sk # 0 and p be as in the
previous lemma. Then Gk is the unique minimizer of

p on [0,Ei]. Moreover, 0k is positive.

gpppg: Since F’(xk;-sk) < 0 by lemma 6.15, the
conclusion that Gk 2 0 follows immediately once we
show that Gk minimizes m over [0.6%]. Uniqueness
of this minimizer follows from the strict convexity of

F.

If zk satisfying step 5 of the algoirthm cannot

be found, then by lemma 6.17 5% minimizes p over

90
[0,Ek], and in step 5 we set Gk = 3%. Suppose then
that Gk 6 (0,3k] is located such that an apprOpriate
vector zk exists. Set y = xk-dksk. Since
zk 6 vf(y)+-Ko(y) C BF(y), for any a 6 [0.Ek] we

have by the subgradient inequality that
cp(c) = F(xk-csk) 2 F(y)+ (<Jlk-Ol)zksk = F(y) = m(ok).
so that Gk minimizes T over [0.3%].

6.19 Corollapy. Let 3k #’0 and xk+1 = xkncxksk

as in step 6 of algorithm 6.3; then F(xk+1) < F(xk).

Proof: This follows from lemma 6.18 and the

. I .-
observation that F (Xk' sk) < 0.

The lemmas stated up to this point prove that the
algorithm is feasible and that F decreases at each
iteration. We now turn to lemmas leading to a convergence

proof.

6.20 Lemma. Let E's X be the minimizer of F and

 

E' be a cluster point of (Xk)’ Then xk converges to E;

Proof: Let 9 be any cluster point of (xk). We
have F(Q) = F(x). Since F is strictly convex, E" is
the unique minimizer, so that Q = I; Thus (xk) is a
sequence, from a compact set X, having only one cluster

point I. and so xk 4 I.

Then

and
show
Xk".

from

91

6.21 Lemma. Let 0 be a cluster point of (s

 

k).

the sequence (xk) converges to ‘2. the minimizer of F.

Proof: We pass to corresponding subsequences (sk.)

A
(Xkl) Asuch that sk, 4 0 and xk, 4 x 6 x. We shall
that x minimizes F, so that by the previous lemma

x. Since the restriction of F to X is continuous

A
within X, to prove that x is a minimizer of F,

A
it suffices to show that F(y) 2 F(x) for all y 6 rel int X.

A A
Let y 6 rel int X. For all i 6 Io(x), aiy < bi = aix,
and so for k’ sufficiently large
. A
ai(y-xk,) < 0 for all 1 6 Io(x). (6.21.1)
. 2
Since Sk’ 4 0 and 6k’ g_lsk.l , ek’ 4 0, so
that for k’ sufficiently large
A
Iek’(xk.) C IO(x)° (6.21.2)
Now there ex1st uk. 6 ch,(xk’) and w . 6 Cck.(xk')
such that Sk’ = vf(xk.)4-uk,4-wk.. By lemma 6.5,
KB (xk) C as v(xk), so that
k k
V(y)-V(x.k.) 2uk.(y-xk.)-ck.. (6.21.3)

Since

f is convex, it follows that

F(Y) -F(xk') 2 Vf(X-kl)(Y-Xkl)+ukl(y-Jfi<l) "Ck!

= Ska(Y-)fi(I)-Wkr(y-Xkl)-€klo (6.2.1.4)

92

Assume that k' is large enough that (6.21.1) and
(6.21.2) hold. Since wk: belongs to the convex cone

generated by {ail 1 6 Iek.(xk')}' in View of (6.21.1)

and (6.21.2) we have wk.(y-xkl) g_0, and so from

(6.21.4)
IMy)—Fm A 2sk4y-xkd-ek. (6215)
when k} is sufficiently large. In the limit (6.21.5) gives
A
F(y)-F(X) 2 O,
proving the lemma.
6.22 Lemma. If 0 is a cluster point of the sequence

(ck) defined in algorithm 6.3, then (Xk) converges to I,

the ndnimizer of F.

.gpppg: Passing to corresponding subsequences
(€kl) and (Xk’)' we may assume that ck: 4 O and
Xk' 4 Q 6 X. By lemma 6.14, step 4 of the algorithm is
executed finitely Often per iteration, and hence the

subsequence (6k!) can be chosen such that
2 _
lye] g|c, y€ — N[vf(xk.)+-K€(xk.)4-C€(xk.)]
and
2 c _ _
lye/2| > 2 — 6km y‘s/2 — N[Vf(XkI)+K£(Xkl)+C£(xka)].
2 2
From these we see that Iy2 c [123 28k), showing that
k

YZE: I 4 0. We now repeat the proof of lemma 6.21,
k

93

replacing s with yze: 7' concluding that x 4 x.

I
k k

6.23 Lemma. The sequence (sk) is bounded.

Proof: Note that

Vf(xk) + K0(xk) C vf(xk) + K€k(xk) + C€k(xk) .

so that
(SR) 3 lvf(xk)l + 1N1K0(xk)11. (6.23.1)

K0(Xk) is one of a finite number of possible polytOpes,
so that there is an upper bound on [N[K0(xk)]l inde-
pendent of k. As f is of class C1 on the compact
set X, the right hand side of (6.23.1) is bounded,

proving the lemma.

 

6.24 Lemma. If the sequence (sk) is bounded

away from 0, then (0k) converges to 0.

Proof: Suppose that (s is bounded away

k)
from O and that Gk f 0. Since dklskl is bounded

above by the diameter of X and (sk) is bounded away
from 0, (0k) is bounded. Given this and the compactness
of X, we can pass to corresponding subsequences (Sk')'
(0k.) and (xk.) such that Sk' 4 S # 0, 0k. 4 d 2 0
and Xk' 4 x 6 X. By lemma 6.19, the sequence (F(Xk'))

is monotone decreasing, so that all of its subsequences have

the same limit, namely F(x). In particular, F(xk'+1) 4 F(x);

94
but Xk’+l = Xk’ ”ak’sk' 4 x-ds , so that
Fm-os)==FM). (624A)

Since F is convex and F(xka-ak.sk.) g F(xk:-Ask:)

for all A 6 [0.6%.], we have
F(Xk' was.” s F(xx' - skew/2) s “X.”
and so in the limit
F(x-ds) gF(x-o(s/2) gF(x). (6.24.2)
Since a 2 O and s #’0, (6.24.1) and (6.24.2) taken together

contradict the strict convexity of F.

6.25 Lemma. Let the sequence ) defined in algorithm

 

(ck
6.3 be such that there exists 6.2 0 satisfying ck.2 c

for every k. For any index i, the inequality

binaixk g_c (6.25.1)

implies the inequality

bimaixk ng._--aixk_*_l . (6.25.2)

Proof: If (6.25.1) holds, then i 6 I (xk), and
€k
so a. 6 C (xk). As was noted in the proof of lemma 6.15,
it follows that aisk 2_0. Since xk+1 = xk-dksk,

(6.25.2) must hold.

6.26 Lemma. Assume that the following hold.
(1) The sequence (ck) in algorithm 6.3 is

such that there exists 6 2 0 with ek‘2 c for all k.

95

(ii) The sequence (0k) converges to 0.

(iii) Some subsequence (Xk') of (xk) converges
to the point x.
Then there exists a subsequence of (xk.), again denoted

(xk.), such that Io(xk.) = IO(X) for every index k'.

ggppfi: Assume that (i), (ii) and (iii) hold.
Since the index sets 10(xk.) are subsets of the
finite set {l,...,m}, we can pass to a subsequence
of (Xk’)' again denoted (xk.), such that for some
subset I of {l,...,m] we have 10(xk.) = I for all
k'. We must show that IO(x) = I. If i 6 I, then

I = bi for all k’, so that in the limit aix = b..

a.x
i k 1

Therefore I C Io(x). Now suppose that i 6 Io(x)\\I.

We derive a contradiction. Since xk+1 = xk-dksk,

with (s ) shown bounded in lemma 6.23 and Gk 4 0,

k
we see that lxk+l-xkl 4 0 as k 4 m. Hence there

exists k such that

O
ai(xk+l-xk) < %- for all k 2_k0. (6.26.1)
Choose p 2kO such that 10(xp) = I and
6* = bi"aixp < 8. (6.26.2)

Such an index p exists because i 6 IO(x) implies
that bi-aixk. 4 bi-aix = 0. Also 6* 2 0, Since

i f I. Let q be the first index such that q 2 p and

bi-aixq g .*/2 . (6.26.3)

96

Now by (6.26.1). (6.26.2) and (6.26.3)

bi"ain-l = bi-aixq4-ai(xq-xq_1)

<e*/2+€/2<€,

and so by lemma 6.25

bi-aixq_lgbi-aixqg €*/2 . (6.26.4)

Note that q- ‘2 p. If q-l = p, then (6.26.4) contradicts
(6.26.2). If q-—1 2 p, then (6.26.4) contradicts the
choice of q as the smallest index greater than p such

that (6.26.3) holds.

6.27 Corollary. Suppose that the following hold.
(i) There exists 8 2 0 such that 5k 2_e for all k.
(ii) There exists n 2 0 such that Iskl 2 n for
all k.
(iii) Some subsequence (xk.) of (xk) converges
to x.
Then there is a subsequence of (xk.), again denoted

(in), such that I0(xk)) = Io(x) for all k'.

Proof: Hypothesis (ii) of this corollary implies

hypothesis (ii) of lemma 6.26 by lemma 6.24.

We are at last prepared to prove the convergence

of our algorithm.

97

6.28 Theorem. Algorithm 6.3 generates either a
terminating sequence whose last term solves problem (6.1.2)
or an infinite sequence converging to the solution of

(6.1.2).

gpppg: In view of lemma 6.13, we need only consider
the case in which algorithm 6.3 generates an infinite
sequence (Xk)‘ In this case, Sk #’0 for every k.
We assume that (xk) fails to converge to the solution

of (6.1.2), and derive a contradiction.

By lemma 6.21 we may suppose: that there exists
n 2 0 such that [ski 2.n for all k. Similarly,
by lemma 6.22 we may assume that there exists 6 2 0
such that [ck] 2 c for all k. Since X is compact

and, by lemma 6.23, (s is bounded, we may pass to

k)

a subsequence (k') of positive integers such that

Xk' 4 x 6 X and s . 4 S y'o. (6.28.1)

k
From step 6 of the algorithm, ku+l = Xk.-Gk:sk..

Since Iskl 2 n for all k, lemma 6.24 ensures that

Gk 4 0, and so Xk’+l 4 x. Passing to a subsequence

of (k’), again denoted (k'), we may suppose that there

exist sets I, J and J’ of indices such that

Iek.(xk’) = I, Jek’(xk:) = J, Jo(xk.+l) = J’ (6.28.2)

I I .
for all k . We assert that J C J. Since xk.+1 4 x,

98

I - . .
J0(xk’+l’ C J0(x) for k; suffiCiently large. Moreover, Since

xk, 4 x and e 2 O, by lemma 6.11 JO(X) C J€(Xkl)

for k’ large enough. AS ck. 2_e. we must have J€(xk,) C

J (x I). Thus for k’ sufficiently large
Ck: k
I _ -
J - JO(Xk'+l) C J0(X’ C Jc(xk” C Jckz(xk" — J’

and so J’ C J.

Using corollary 6.27, since 4 x we can pass

xk’+1
to yet another subsequence, again denoted (k’) such that,

10(Xkr) = IO(X) = IO(Xk’+l) (6.28.3)
for all k'. Now set
K = conv{gj| j 6 J} and C = cone{ai] i 6 I]. (6.28.4)

Due to (6.28.2) and (6.28.4), we see that for all k’

Kek’(xk,) = K and Cck)(xk’) = C. (6.28.5)

From (6.28.3) we deduce that ck. < 6%. for each k’;
for if 0k. = Ck), some constraint inactive at Xk'
becomes active at xk’+l' and so 10(xk.) # 10(xk’+1"

Thus for each k’, the vector 2k. Specified in step

5 of algorithm 6.3 must exist, i.e.
Zkl E Vf(XkI+l)4'KO(XkI+l)
and

I = 0. (6.28.6)

99
Taking into account (6.28.2) and (6.28.4), we have
KO(Xk’+l’ = conv(9j) j 6 J’} C conv{qj)j E J}
and so
zk. 6 vf(xk.+l)4-K. (6.28.7)

Since vf(xk.+1) 4 vf(x) and K is compact, by passing
to still another subsequence (k’) and applying (6.28.7)
we may assume that there exists 2 6 vf(x)4-K such that

Z I 4 Z.

k
From steps 2 and 3 of the algorithm, we have that

Ski = N[Vf(XkI) +K€kl(xk’) +C€k

and so, in view of (6.28.2) and (6.28.4),

’(Xkl)]

sk’ = N[Vf(X_kl) +K+ C] .
Since Xk' 4 x and sk. 4 s, it follows easily that

s = N[vf(x)+K+C]. (6.28.8)

Now 2 6 vf(x)4-K.C vf(x)+-K4—C, and so by (6.2.7) we
have s(z-S)20, i.e.
2
zs 2 [SI . (6.28.9)

As Iskzl 2.0 for all k’, clearly (6.28.9) implies
that

25 2,02 > O. (6.28.10)

100
On the other hand, letting k’ 4 w in (6.28.6) yields

25 = 0, (6.28.11)

contradicting (6.28.10). Thus our assumption that (xk)
fails to converge to the solution of (6.1.2) cannot be
valid. The proof that the algorithm generates a sequence

converging to the optimal solution is now complete.

CHAPTER VII

IMPLEMENTATION OF THE SUBGRADIENT PROJECTION ALGORITHM

In this chapter we propose methods for implementing
algorithm 6.4 and describe computational experiences. The
algorithm was coded in FORTRAN on a CDC 6500 computer and
tested on several problems of the form described in section
6.1. We present here the results of those tests.

The program used a value of 10'.5 for the parameter
co in algorithm 6.4. For computational efficiency, we
divided 8 by 10 rather than by 2 in step 4 of the
algorithm. Step 4 can be rewritten so that e is replaced
by e/a for any fixed a 2 l. The algorithm terminated
when the euclidean length of the projection yO in step 1
was less than a specified figure, usually 10-10.

The two major problems in implementation were the
nearest-point projection subroutine, required in steps 1
and 2 of the algorithm, and the line search procedure used
in step 5. We first describe the line search.

In step 5, it is required to determine a 6 [0,5]
such that 25 = O for some 2 6 vf(x-as)4-KO(x-ds).
Using the notation established in section 6.3, the
orthogonality condition becomes

23 Bigis + vf(x-ds)s = O , (7.0.1)
i6JO(x-ds)

8120 forall i, 23 B.=l.
i6JO(x-ds)

lOl

102

Rewriting (7.0.1) as

Z) Bi[g.s + vf(x-ds)s] = O , (7.0.2)
i6JO(x-ds) l

we note that a solution (Bi) exists iff for

i6JO(x-ds)

some i,j 6 J x-—ds) either gis-tvf(x-ds)s = O or

O(
gis+-vf(x-—ds)s and gjs4-vf(x-ds)s have Opposite
signs. We can summarize the characterization of a as

follows.

7.1 Lemma. 0 6 [0,5] is the desired solution in

 

step 5 of algorithm 6.4 if and only if either

(a) there exist i,j 6 JO(X-ds) such that

[915+vf(X-qs)S][qu+vf(x-qs)S] g0 (7.1.1)

or

(b) d = a and no triple (0,1,j) with O g_d £23

and i,j 6 J (x-as) satisfies (7.1.1).

O

Our line search proceeds as follows. We begin with
d = 0 and locate the first value 3 2 a at which
JO(x-as)\JO(x-ds) # q. We check whether there exist
i,j 6 JO(x-3s) satisfying (7.1.1). If not, we check

whether there exists i 6 J (x-as) n J0(x-ds) such

0
~ A

that gis4-vf(x-ds)s and gis-tvf(x-os)s have Opposite

signs. If so, we locate the value of 0 between 6 and

A
a for which gis4-vf(x-ds) = O and terminate the search.

This value of a can be found by Newton's method When f is

103

A ._
sufficiently smooth. If no such 1 exists and 0 < 0,

~ A
we replace a with d and repeat the process.

7.2 Line searchAprocedure. Begin with x,s,a

 

given.

Step 0: Set a = O and consider all indices in

Jo(x) as untested.

Step 1: Choose an untested index i 6 JO(x-ds).

A 1..
Step 2: Compute d = max {0 6 [0,0] :vi(x-os) =

A a.
v(x-ds)]. If a = 0, go to step 1.

A
Step 3: For each j 6 JO(x-ds) do the following:
2;: If

A
[gis4-vf(x-ds)s][gjs-tvf(x-as)s] g_0 ,

A
terminate the search with a = d;

.gp: If j 6 JO(x-ds) and
.. A
[gjs+-vf(x-ds)s][ngA-vf(x-ds)s] < O ,
go to step 5.
A — ~ , A
Step 4: If a < d, replace 0 With a and go to

step 1, treating all indices as untested. Otherwise, terminate

the search with 0 =‘E.

~ A
Step 5: Locate the value a 6 [0,0] for which

gjs4-vf(x-ds)s = O and terminate the search.

104

7.3 Projection procedure. The other major problem
in implementing the algorithm is the calculation of the
nearest point to the origin in a convex set. The convex
set in question is vf(x)4-K€(x)-tC€(x) with c 2_0.
For convenience, let uj = gji-vf(x), j = l,...,r.

Since vf(x)4-K€(x) = conv {uj :j 6 J€(x)}, the

projection problem can be posed as:

 

 

( h k W
061R+.fi361R+
k
v
4.15:]. 5
< v=l v (7.3.1)
1 h k 2
3 2: 0a:.L + Z Bu.l (min)
0:1 H p v=l V 3v
J

where I€(x) = {il""'ih} and J€(x) = [j1,...,jk}. The

objective function is a quadratic form in 0 and B, and so
the projection problem is simply a quadratic programming
problem. Letting g = (0,6) and H = (a. ...a. u. ...u. )

1l 1h 31 3k '

we can write the objective function as 1/2 5 HtH 6, where

(using superscript t to denote transposes)

GtG GtU

!
I
UtG UtU /

I

G = (a. ...a. ) and U = (u. ...u. ).
1 1h 11

HtH is positive semidefinite but not positive definite,
and problem (7.3.1) can in general have multiple optimal

solutions. Despite this, there exists a unique nearest

105

point N[vf(x)4-K€(x)-tC€(x)], since vf(x)4-K€(x)4-C€(x)
is a convex set and the euclidean norm is a strictly convex
norm.

Our approach is to locate a vector 6 = (0,5)
satisfying the Kuhn-Tucker necessary and sufficient

conditions for Optimality in the quadratic programming

 

 

problem (7.3.1). These conditions are:
(
§.y6m:'+k.w613
e§ = l g
( t (7.3.2)
H H§-y4-we = 0
L 6y = 0 2
h k .
where e= (0,...,0)x(l,...,1) 6R xR and w is

a sign-unrestricted scalar. We initially employed the

first phase of the Dantzig Two-Phase linear programming
procedure [12], with a restricted basis-entry rule to
maintain the complementarity condition gy = O. This

method can encounter nonOptimal tableaus in which further
pivoting is blocked by the complementarity restriction,

even when HtH is positive definite. A Simple modification,
however, solves this. Since all equality constraints,

with one exception, are homogeneous, we can construct an
initial Phase I - feasible solution by setting any one

of the BV equal to unity and completing the basis

106

with artificial variables (which may require multiplication
of some of the homogeneous equations by -l). From this
point onward, the modified Phase I procedure reduces to

an algorithm of Wolfe [18]. Our elimination of all

linear terms in the objective function, by means of the
substitution uj = vf(x)4-gj, causes our problem to
satisfy conditions under which Wolfe's proof of convergence
of the algorithm applies with only minor alterations, even
through HtH is only semidefinite. Wolfe assumes that

the system of equality constraints is nondegenerate,

i.e. that the constraint equations are linearly independent
and that no basic feasible solutions with a basic variable
equal to zero occur. The problem before us clearly satisfies
the assumption of independent constraints, while Wolfe's
proof can be modified to obviate the assumption that no

basic variable vanishes.

7.4 Computational results: the PIES counterexample.
Algorithm 6.4 was first tested on example 4.16,
although convergence in this case has not been proved due

to the lack of a potential for -p. When tested with
the parameter T = 3 and a convergence criterion of less

than 10-10

error in the euclidean norm of q, the
algorithm converged in one to three iterations from a
variety of starts. The general pattern was one step

from the starting point to an edge of Q* containing

the equilibrium vector (4,2), and then another step to

107

the solution. Where a third step was required, it
appeared to be a very short step, perhaps correcting
rounding errors. The only difficulties occurred when
the initial point was too close to the origin, where we
believe the size of the components of -p(Q) is SO much
greater than that of the generators ai and gj as to
cause severe roundoff errors. The algorithm was also
tested with parameter values T = 7 and T = 10, with

similar convergence results.

7.5 Computational results: the PIES-VAR counterexample.

We next tested algorithm 6.4 on example 5.7. AS with
the previous example, the "gradient" -p does not actually
have a potential, and so convergence of the algorithm is
not guaranteed. Example 5.7 was attempted with three
different starts. When started at q = (1,1), the
algorithm reached the equilibrium (2,2), to within 10”8
euclidean norm, in one iteration. From a start at (4,2),
convergence was oscillatory: the error was approximately

10 after

2>(lO-7 after 40 iterations and 2.2)(10-
60 iterations. Started from (1,2), however, the algorithm
appeared to fall into a four-iteration cycle. unlike the

previous attempt, in which the norm of the error decreased

monotonically, the norm of the error in the four step

cycle fluctuated between .511 and .961.

108

7.6 Computational results: Wolfe's example.
Wolfe [19] considers an example in two variables with

f identically zero:
v(x) = max {V1,v2,v3] (min)

where vl(x) = -x1, v2(x) = xl-I—x2 and v3(x) = xl-2x2.
The level sets of v are triangles nested about the origin.
The example was solved under the added constraints
-10 S.le x2 g 10. The global minimizer (0,0), interior
to the constraint region, was reached in at most two
iterations from a variety of starts. Similar results were
obtained using constraints which placed (0,0) on the
boundary of the feasible region. The method also reached
the constrained Optimum in at most two iterations using
constraints which put (0,0) exterior to the feasible
region.

We note that Wolfe's example (and Powell's example,
which follows) can be posed as linear programs, but with no

gain in Speed of convergence.

7.7 Computational results: Powell's example.

 

Wolfe [19] reports the following example due to Powell,

on which the conjugate gradient method converges only linearly:

v(x) max {vj(x) :j = O,...,4} (min)

where vj(x) cos(2wj/5)4—xzsrn(2wj/5). Linear convergence

x1
of the conjugate gradient method is Observed when started

at any point of the form (p cos(Wj/5). P sin (Wj/5)). The

109

contours of v are regular pentagons centered at the
minimizer (0,0). Algorithm 6.4 converged to the solution
in two iterations from any feasible starting point using
constraints which placed the origin in the interior of the
constraint region. Similar results occurred when constraints
were used which placed the origin on the boundary of the
feasible region. USing the constraints x1 g_O, x2 g_0,
xl-i-x2 g -l, for which the origin is infeasible, the
algorithm typically took approximately five to ten iterations

to reach the constrained minimum.

7.8 Computational results: a larger example.

The algorithm was tested on an example having nine
variables. The feasible region was a hypercube with the
Optimal solution at one vertex. The smooth part f of the
objective was a strictly convex quadratic, and there were
three affine functions Vi' The algorithm converged to

13

within 10- of the solution (which was a unit vector),

in at most ten iterations, from a variety of starts.

LIST OF REFERENCES

[l]

[2]

[3]

[4]

[5]

[6]

[7]

[8]

[9]

[10]

[ll]

[12]

LIST OF REFERENCES

D.P. Bertsekas and S.K. Mitter, "A Descent Numerical
Method for Optimization Problems with Nondifferentiable
Cost Functionals," SIAM J. CONTROL, 11, 4 (1973).

G. Dantzig, Linear Programmipg and Extensions,
Princeton University Press, Princeton, NJ (1963).

B. Curtis Eaves, "A Locally Quadratically Convergent
Algorithm for Computing Stationary Points," TR SOL
78-13, Department of Operations Research, Stanford
University (1978).

W.W. Hogan, "Energy Policy Models for Project
Independence," Comput. and Ops. Res., 2 (1975).

 

W.W. Hogan, "Project Independence Evaluation System:
Structure and Algorithms," Proc. p£_§ymposium.ip Appl.
Math. 9; the A.M.S., 21 (1976).

 

C.L. Irwin, "Analysis of a PIES-type Algorithm,"
presented at the Symposium on Energy Modeling and
Net Energy Analysis, Colorado Springs (1978).

C. Lemarechal, "An Extension of Davidon Methods to
Nondifferentiable Problems," Math. Prog. Study 3
(1975).

 

H. Nikaido, Convex Structures and Economic Theory,
Academic Press, New York, NY (1968).

D. Nissen and D. Knapp, "A Regional Model of Interfuel
Substitution," Energy: Mathematics and Models, Proc.
of a SIMS Conference on Energy. SIAM (1975).

E. Polak, Computational Methods ip_0ptimization,
Academic Press, New York, NY (1971).

 

R.T. Rockafellar, Convex Analysis, Princeton University
Press, Princeton, NJ (1970).

 

J.B. Rosen, "The Gradient Projection Method for
Nonlinear Programming, Part I," J. SIAM, 8, 1 (1960).

110

[l3]

[14]

[15]

[l6]

[17]

[18]

[19]

111

V.P. Sreedharan, "Least Squares Algorithms for Finding
Solutions of Overdetermined Linear Equations which
Minimize Error in an Abstract NOrm," Numer. Math..

17 (1971).

 

V.P. Sreedharan, "Least Squares Algorithms for Finding
Solutions of Overdetermined Systems of Linear Equations
which Minimize Error in a Smooth Strictly Convex Norm,"

.q. Approx. Thry.. 8, l (1973).

J. Sweeney, "On the Uniqueness of PIES Solutions,"
Department of Engineering Economic Systems, Stanford
University (1976).

R.S. Varga, Matrix Iterative Analysis, Prentice-Hall,
Inc., Englewood Cliffs, NJ (1962).

M. Wagner, "Project Independence Evaluation System
Integrating Model," Energy: Mathematics and Models,
Proc. of a SIMS Conference on Energy, SIAM (1975).

P. Wolfe, "The Simplex Method for Quadratic
Programming," Econometrica, 27, 3 (1959).

P. Wolfe, "A Method of Conjugate Subgradients for
Minimizing NOndifferentiable Functions," Math. Proq.
Study 3 (1975).

 

"I(l((l'lllfllllttlfliis