Parallel Variable Distribution for Constrained Optimization.
In the parallel variable distribution framework for solving optimization problems (PVD), the variables are distributed among parallel processors with each processor having the primary responsibility for updating its block of variables while allowing the remaining secondary variables to change in a restricted fashion along some easily computable directions. For constrained nonlinear programs convergence theory for PVD algorithms was previously available only for the case of convex feasible set. Additionally, one either had to assume that constraints are block-separable, or to use exact projected gradient directions for the change of secondary variables. In this paper, we propose two new variants of PVD for the constrained case. Without assuming convexity of constraints, but assuming block-separable structure, we show that PVD subproblems can be solved inexactly by solving their quadratic programming approximations. This extends PVD to nonconvex (separable) feasible sets, and provides a constructive practical way of solving the parallel subproblems. For inseparable constraints, but assuming convexity, we develop a PVD method based on suitable approximate projected gradient directions. The approximation criterion is based on a certain error bound result, and it is readily implementable. Using such approximate directions may be especially useful when the projection operation is computationally expensive.
AMS subject classications: 90C30, 49D27
1. Introduction and Motivation
We consider parallel algorithms for solving constrained optimization
problems
min
Research of the rst author is supported by CNPq and FAPERJ. The second author
is supported in part by CNPq Grant 300734/95-6, by PRONEX{Optimization,
and by FAPERJ.
c
2001 Kluwer Academic Publishers. Printed in the Netherlands.
Sagastizabal and Solodov
where C is a nonempty closed set in < n , and f : < continuously
dierentiable function. Our approach consists of partitioning the
problem variables x 2 < n into p blocks x , such that
and distributing them among p parallel processors. Note that in this
notation we shall not explicitly account for possible re-arranging of
variables.
In most parallel algorithms, an iteration usually consists of two
steps: parallelization and synchronization, [2]. Consider for the moment
the unconstrained case, corresponding to Typically, the
synchronization step aims at guaranteeing a su-cient decrease in the
objective function, while the parallel step produces candidate points (or
candidate directions) by simultaneously solving certain subproblems
Each (P l ) is dened on a subspace of dimension
smaller than n, in such a way that these p subspaces
\span" the whole variable space < n . Most methods, such as Block-
Jacobi [2], updated conjugate subspaces [10], coordinate descent [21],
and parallel gradient distribution [14], dene (P l ) as a minimization
problem on the l-th block of variables, i.e., in < n l . More recently,
Parallel Variable Distribution (PVD) algorithms, introduced in [6] and
further studied and extended in [19, 20, 7], advocated subproblems
l ) of slightly higher dimensions than n l . In the l-th subproblem, in
addition to the associated \primary" optimization variables x l , there
are representing in a condensed form all the
remaining n n l problem variables (see Algorithm 1 below). Since those
remaining n n l variables are allowed to change only in a restricted
fashion (in a subspace of dimension p 1), the additional computational
burden of solving such enlarged subproblems is not big. The idea is
that the \forget-me-not" terms add an extra degree of freedom yielding
algorithms with better robustness and faster convergence. We refer the
reader to [6, 22, 12] for numerical validation of PVD-type methods. We
note that not having any variables in parallel subproblems completely
xed can be especially important in the constrained case, i.e., when
this, the method can simply fail. We shall
return to this issue later on.
1.1. General PVD Framework
To formalize the notion of primary and secondary variables, we need to
introduce some notation. Let  l denote the complement of l in the index
set is the number of parallel processors. Given some
Parallel Variable Distribution 3
direction d i 2 < n , which we shall call PVD-direction, we denote by D i
l
the n l  (p 1) block-diagonal matrix formed by placing the blocks
of the chosen direction d i along
its block diagonal as follows:
d id id i
l 1
Note that if I denotes the identity matrix of appropriate dimension,
the linear transformation
l :=
I n l n l
l
l
maps
l  l In the unconstrained
case, more general transformations can be used, which give rise to a
fairly broad parallel variable transformation framework discussed in [7].
However, the theory of [7] does not appear to extend to the constrained
case, which is the focus of the present paper.
We describe next the basic PVD algorithm [6, 19, 20].
ALGORITHM 1. (PVD) Start with any x 0 2 C. Choose a PVD-
direction d
Having x i , check a stopping criterion. If it is not satised, compute
x i+1 as follows.
Parallelization. For each processor l 2 (possibly
solution (y i
l
l +D i
l  l )
l
l  l
Synchronization. Compute
l2f1;:::;pg
l (y i
choose new PVD-direction d i ; and repeat.
4 Sagastizabal and Solodov
The idea of PVD algorithm is to balance the reduction in the number
of variables for each (P l ) with allowing just enough freedom for
the change of other (secondary) variables. Due to this, the parallel
subproblems better approximate the original problem that has to be
resolved. Note that because the secondary variables can change only
along chosen xed directions, their inclusion does not signicantly increase
the dimensionality of the parallel subproblems (P l ). Indeed, the
number of variables in (P l ) is than if the
secondary variables were excluded. Of course, in the context of parallel
computing it is reasonable to assume that p is small relative to n l .
The synchronization step in Algorithm 1 may consist of minimizing
the objective function in the a-ne hull, subject to feasibility, of all the
points computed in parallel by the p processors. This would require
solving a p-dimensional problem, which is again small compared to
the original one. If (1.1) is a convex program, one can alternatively
dene x i+1 as a convex combination of the candidate points. In prin-
ciple, for convergence purposes, any point with the objective function
value at least as good as the smallest computed by all the processors
is acceptable. As for PVD-directions, they are typically some easily
computable feasible descent directions for the objective function f at
the current iterate x i , e.g., quasi-Newton or steepest descent directions
in the unconstrained case.
Algorithm 1 is a rather general framework, which has to be further
rened and specialized to obtain implementable/practical versions. In
this respect, the two principal questions are:
{ how to set up parallel subproblems; e.g., how to choose PVD-
directions, and
{ how to solve each parallel subproblem, including some criteria for
inexact resolution.
When (1.1) is unconstrained, some of these issues have been addressed
in [19], which contains improved convergence results (compared
to [6]), including linear rate of convergence. These results, as well
as useful generalizations such as algorithms with inexact subproblem
solution and a certain degree of asynchronization, were obtained by
imposing natural restrictions (of su-cient descent type) on the PVD-
directions. An even more general framework for the unconstrained case
was developed later in [7], where subproblems are obtained via certain
nondegenerate transformations of the original variable space. These
transformations can be very general, and there need not even be a
distinction between primary and secondary variables. However, this
approach does not seem to extend to the constrained case. It seems
Parallel Variable Distribution 5
also that intuitively justiable transformations do have some kind of
primary-secondary variable structure. For those reasons, in the present
paper we shall restrict our consideration to specic transformations of
the form (1.2).
When (1.1) is a constrained optimization problem, many questions
are still open, especially for a nonconvex feasible set C. When C is
convex with block-separable structure (i.e., C is a Cartesian product
of closed convex sets), it was shown in [6] that every accumulation
point of the PVD iterates satises the rst-order necessary optimality
conditions for problem (1.1). It was further stated that in the case
of inseparable convex constraints, the PVD approach may fail. This
conclusion was supported by a counter-example, which we reproduce
below:
This strongly convex quadratic program has the unique global solution
at 1). Consider any point ^
x 6=  x such that ^
observe that
Therefore, if we apply Algorithm 1 using ^ x as its starting point and
xing the secondary variables, then this PVD variant will stay at this
same nonoptimal point thus failing to solve the original problem.
This shows that in the constrained case, a special care should be taken
in setting up the parallel subproblems.
For a general convex feasible set, it was shown in [20] that using the
projected gradient direction d(x) := x P C [x rf(x)] for secondary
variables does the job (here P C [] stands for the orthogonal projection
map onto the closed convex set C). Specically, it was established that
setting d i := d(x i ) would ensure convergence of PVD methods for
problems with general (inseparable) convex constraints. Some criteria
for inexact subproblem solution were also given in [20]. When C is
a polyhedral set, computing the projected gradient direction requires
solving at every synchronization step of Algorithm 1 a single quadratic
programming problem. For this, a wealth of fast and reliable algorithms
is available [11, 5]. It should be noted that in the case of nonlinear
constraints, the task of computing the projected gradient directions
is considerably more computationally expensive. Actually, even in the
a-ne case computing those directions exactly (or very accurately) may
turn to be rather wasteful, especially when far from the solution of the
original problem. Therefore, improvements are necessary.
6 Sagastizabal and Solodov
In this paper, we propose two new versions of PVD for the nonlinearly
constrained case. The rst one applies to problems with block-
separable nonconvex feasible sets. It is based on the use of sequential
quadratic programming techniques (SQP). Our second proposal is for
inseparable convex feasible sets. We introduce a computable approximation
criterion which allows to employ inexact projected gradient
directions. This criterion is based on an error bound result, which is of
independent interest. We emphasize that its is readily implementable
and preserves global convergence of PVD methods based on exact
directions.
Our notation is fairly standard. The usual inner product of two vectors
denoted by hx; yi, and the associated norm is given by
using other norms, we shall specify them explicitly.
Analogous notation will be used for subspaces of any other dimensions;
for example, the reduced subspaces < n l . For the sake of simplicity,
we sometimes use a compact (transposed) notation when referring to
composite vectors. For instance, stands for the column vector
l
. For a dierentiable function will denote the
n-dimensional column vector of partial derivatives of f at the point
. For a dierentiable vector-function c : <
denote the m n Jacobian matrix whose rows are the transposed gradients
of the components of c. For a function h
we
its partial derivatives are Lipschitz-continuous
on the set X with modulus L > 0.
2. Nonconvex Separable Constraints
Suppose the feasible set C in (1.1) is described by a system of inequality
constraints:
. In this section, we assume that C has block-
separable structure. Specically,
Our proposal is to solve parallel subproblems (P l ) of Algorithm 1 in-
exactly, by making one step of the sequential quadratic programming
method (SQP) [3, Ch. 13]. Since SQP methods are local by nature,
Parallel Variable Distribution 7
we modify the synchronization step in Algorithm 1 by introducing a
suitable line-search based on the following exact penalty function [9]:
l
where  is a positive parameter, and y
taken componentwise. Given x i , the l-th block of rf(x i ) will be denoted
by
l := (I n l n l
l
Our algorithm is the following.
ALGORITHM 2. Start with any x 0 2 C. Choose parameters  > 0
and  ;  2 (0; 1), and positive denite n l n l matrices M 0
Having x i , check a stopping criterion. If it is not satised, compute
x i+1 as follows.
Parallelization. For each processor l 2
l
l as a KKT point of
h- l ; M i
l - l i
c l
l
l )- l
Synchronization. Dene
Line-search. Choose  i
l
l
Using the merit function (2.5), nd m i , the smallest nonnegative integer
m, such that
l ,
repeat.
The above proposal has to be compared to the original PVD algorithm
[6]. Two remarks are in order. First, Algorithm 2 can be viewed as
a PVD method with inexact solution of parallel subproblems. Indeed,
\hard" nonlinear PVD subproblems (P l ) of Algorithm 1 are approximated
here by \easy" quadratic programming subproblems (QP l ). And
secondly, the PVD approach is extended to the case of nonconvex
8 Sagastizabal and Solodov
constraints. Note that the inclusion of forget-me-not terms does not
seem to be crucial in the block-separable case (for example, such terms
were not specied in the analysis of [6]). Thus, we drop here secondary
variables in (P l ), by taking null PVD-directions. However, the use of
secondary variables deserves further study for feasible sets with inseparable
constraints. We also note that Algorithm 2 can be thought of as
a distributed parallel implementation of SQP, where a block-diagonal
matrix is chosen to generate quadratic approximations of the objective
function. We remark that the contribution of Algorithm 2 is meant
primarily to PVD framework rather than general SQP methods.
In Algorithm 2, we assume that subproblems (QP l ) have nonempty
feasible sets for every iteration (this is guaranteed, for example, if for
each l the Jacobian c 0
l maps < n l onto < m l , or if c l is convex). Alter-
natively, it is known that feasibility of subproblems can be forced by
introducing an extra \slack" variable [1, p. 377]. It is sometimes argued
that SQP methods are not convenient for solving large-scale (with n
and m large) nonlinear programs when there are inequality constraints
present. More precisely, it is argued that the combinatorial aspect introduced
by the complementarity condition in the KKT system associated
to each QP makes the subproblems resolution relatively costly. On the
other hand, SQP techniques are known to be very e-cient for small to
medium size problems, and they are often the method of choice in that
setting. In relation to Algorithm 2, it is important to note that each
subproblem (QP l ) is a quadratic programming problem of relatively
small dimension (n l and m l are presumed to be small compared to
n and m). There exist a number of very fast and reliable algorithms
for solving such problems (see, e.g., [5]). In Section 4, we report some
numerical results for Algorithm 2, obtained via simulation on a serial
computer.
For further reference, we state the optimality conditions for (QP l
the pair (- i
l solves the KKT system
l +M i
l
l
l
c l
l
l )- i
l
l  0 and h i
l
l )- i
Denoting by  0
(x; d) the usual directional derivative of the merit function
at x 2 < n in the direction d 2 < n , we next state that - i is
a descent direction for this function at the point x i . This in turn will
imply that the line-search procedure in Algorithm 2 is well-dened.
Parallel Variable Distribution 9
LEMMA 1. If f; c 2 C 1;1
h- i
l i
l
l
l are dened in Algorithm 2.
Proof. Dening the index-sets I
we have that (e.g., see [9, p. 301])
Using (2.7b) componentwise, we have
implying that
Hence,
l
l
To obtain the right-most inequality in (2.8) we proceed as follows:
l i
l
l
l
l i
l i
l
l )- i
l i
l
Sagastizabal and Solodov
where the second equality is by (2.7a), the fourth is by (2.7c), and the
inequality follows from (2.7c) and the fact that c +
l
l
l ).
Combining the latter relation with (2.9), we now have that
l
l
l
l )j 1
(j i
l
l )j 1
l i
l
l
where the second inequality is by the Cauchy-Schwarz inequality, and
the last is by the choice of  i . This completes the proof.
For Algorithm 2 to be globally convergent, we assume that the
penalization parameter  i is kept bounded above, which essentially
means that the multipliers  i
l stay bounded. From the practical point
of view, the latter assumption is natural. In what follows, we show
convergence of Algorithm 2 to Karush-Kuhn-Tucker (KKT) points of
(1.1), i.e., pairs (x;
THEOREM 1. Suppose that f; c 2 C 1;1
let the feasible
set C have block-separable structure given by (2.4). Let f(x be a
sequence generated by Algorithm 2. Assume there exists an iteration
and the matrices M i
l are
uniformly positive denite and bounded for all
iteration indices i. Then either   i
or every accumulation
point (x;
) of the sequence f(x i ;  i )g is a KKT point of the problem,
i.e., it satises (2.10).
Proof. If for some iteration index i it happens that
l
l
then (2.7a)-(2.7c) reduce to
l
l
l
c l
l
l  0 and h i
Parallel Variable Distribution 11
Now taking into account separability of constraints, it follows that
the KKT system (2.10).
Suppose now that (2.11) does not hold for any i. By (2.8) in Lemma
is then a direction of descent for the merit function   i
at x i .
By standard argument, the line-search procedure is well-dened and
terminates nitely with some stepsize t i > 0. The entire method is
then well-dened and generates an innite sequence of iterates.
We prove rst that the sequence of stepsizes ft i g is bounded away
from 0. Take any t 2 (0; 1]. Since f 2 C 1;1
we have that
Similarly, since c 2 C 1;1
using the equivalence of the norms in the
nite-dimensional setting, for some R 1 > 0 we have that
l
l
l
l
l )- i
l
l
l
l )- i
l
l )j 1
c l
l
l
l
l )- i
l
l
l
l
l
l
l )- i
l )+ (1 t)c l
l
where the equality is obtained by adding and subtracting tc l
l ), the
second inequality follows from jaj 1  ja
used in the last inequality.
Re-writing the above relation, we obtain
l
l
l )j 1  j
l
l )- i
l
tj
c l
l
l )- i
l
l
l
l
l
where the second inequality is by the convexity of j() and the
equality is by (2.7b).
Together with (2.12), the last inequality yields
l
l
l )j 1
l
l
12 Sagastizabal and Solodov
is a xed constant depending on R 1
. By a direct
comparison of the latter relation with (2.6), we conclude that (2.6) is
guaranteed to be satised once m is large enough so that
within the set of t satisfying In particular,
since the line-search procedure did not accept the stepsize
it follows that
either
By (2.8) in Lemma 1 and the uniform positive deniteness of matrices
l , there exists R 3 > 0 such that
l i
and using (2.13), we conclude that
By the assumption that
for all i  i 0 , from (2.6) it follows that
the sequence f
nonincreasing. Hence, it is either unbounded
below, or it converges. In the latter case, (2.6) implies that  t i
and since t i   t > 0, it holds that
By the denition of  i and the uniform positive deniteness of matrices
l , we conclude that
l
l
passing onto the limit in (2.7a)-(2.7c) as i !1, and taking into
account the boundedness of the matrices M i
l , we obtain the assertions
of the theorem.
3. Convex Inseparable Constraints
Suppose now that the feasible set C is dened by a system of convex
inequalities, in (2.3) has convex components
We emphasize that in this section we do not assume
separability of constraints. In this setting, it appears that the only currently
known way to ensure convergence of the PVD algorithm is to use
Parallel Variable Distribution 13
the projected gradient directions for the change of secondary variables
[20]. Omitting the iteration indices i, if x is the current iterate, to
compute these directions one has to solve a subproblem of the following
structure:
min
This is done by some iterative algorithm. As already discussed in Section
1, solving this problem in the general nonlinear case can be quite
costly. Moreover, when far from a solution of the original problem,
exact (or even very accurate) projection directions are perhaps unnec-
essary. This suggests developing a stopping rule for solving (3.14) or,
equivalently, an approximation criterion for projection directions to be
used in the PVD scheme. For algorithmic purposes, it is important to
make this criterion constructive and implementable.
Assuming some constraint qualication condition [13], we have that
z solves (3.14), i.e.,
if, and only if, the pair (z;
the KKT system
r z L(z;
u  0 and hu;
where
is the standard Lagrangian for problem (3.14).
Suppose z 2 C and u
is some current approximation to a
primal-dual optimal solution of (3.14), generated by an iterative algorithm
applied to solve this problem. Lemma 2 below establishes an error
bound for the distance from z to P C [x rf(x)] in terms of violations
of the KKT conditions (3.15) by the pair (z; u). A nice feature of this
estimate is that unlike some (perhaps, most) error bounds results in
the literature (see [18] for a survey), it does not involve any expressions
which are not readily computable or any other quantities which are
not observable. Furthermore, this error bound holds globally, i.e., not
just in some neighbourhood of the solution point. Thus it can be easily
employed for algorithmic purposes.
Lemma 2 is related to error bounds for strongly convex programs
obtained in [15]. However, Theorem 2.2 in [15] involves certain constants
which are in general not computable, while Corollary 2.4 in [15]
assumes not only that z is primal feasible, but also that (z; u) is dual
feasible, which here means that u  In Lemma 2
14 Sagastizabal and Solodov
we only assume that z is primal feasible and that the approximate
multiplier u is nonnegative. Our assumptions are not only weaker but
they also appear to be more suitable in an algorithmic framework, as
they will be satised at each iteration of many standard optimization
methods. Dual feasibility, in contrast, is unlikely to be satised along
iterations of typical algorithms, except in the limit.
LEMMA 2. Let c : < be convex and dierentiable, and suppose
that the set C given by (2.3) satises some constraint qualication.
Then for any z 2 C and u
holds that
where
"(z; u) :=2 jr z L(z; u)j +2
jr z L(z; u)j 2 4hu; c(z)i: (3.17)
Proof. Let
u be the associated multiplier, so
that the pair (z;
u) satises (3.15). Denote
where L(; ) is dened by (3.16). Take any z 2 C and u
. Denoting
further we have that
ui
jz
where the inequality follows from u  0,
and the facts that, by
the convexity of c(), c(z) c(z) c 0 (z)(z z)  0 and c(z) c(z)
On the other hand,
zi
jr z L(z; u)jjz
where the equality follows from the KKT conditions (3.15), and the
inequality follows from the facts that
c(z)  0, and the Cauchy-Schwarz inequality. Denoting t := jz  zj ;  :=
jr z L(z; u)j  0 and  := hu; c(z)i  0, and combining (3.18) with
(3.19), we obtain the following quadratic inequality in t:
Parallel Variable Distribution 15
Resolving this inequality, we obtain that
Recalling denitions of the quantities involved, we conclude that
jz  zj 2 jr z L(z; u)j +2
jr z L(z; u)j 2 4hu;
We propose the following PVD algorithm based approximations of
the projected gradient directions.
ALGORITHM 3. Choose parameters  1 2 (0; 1) and  2 2 (0; (1
Start with any x 0 2 C. Set i := 0.
Having x i , check a stopping criterion. If it is not satised, proceed as
follows.
PVD-direction choice. Compute z
and the associated approximate Lagrange multiplier u
for problem (3.14) satisfy
Compute
Parallelization. For each processor l 2 compute a solution
l
dened in Algorithm 1.
Synchronization. Compute
l2f1;:::;pg
l (y i
l
repeat.
Note that Algorithm 3.1 is a general-purpose method for problems
with no special structure. The example presented in the Introduction
shows that for such problems computing a meaningful direction
for secondary variables is indispensable. Compared to any alterna-
tive, computing an approximate projected gradient direction seems to
be quite favorable in terms of cost, and perhaps even the best one
Sagastizabal and Solodov
can do. Regarding the tolerance criterion (3.20), we note that if z i is
the exact projection point then "(z
the rst-order necessary optimality condition
From this, it is easy to see that the projection
problem (3.14) always has inexact solutions satisfying (3.20). Hence,
the method is well-dened.
As for the convergence properties of Algorithm 3, our result is the
following.
THEOREM 2. Let c : < be convex and dierentiable, and
L (C). Suppose fx i g is a sequence generated by Algorithm 3.
Then either f is unbounded from below on C, or the sequence ff(x i )g
converges, and every accumulation point of the sequence fx i g satises
the rst-order necessary optimality condition.
Proof. Consider any iteration processor l 2
and the point
l d i
We rst show that this point is feasible
for the corresponding subproblem (P l ) of minimizing the function
l
l
l  l ). Indeed, using the notation (1.2),
l d i
l D i
l e l
l d i
l d i
where the rst equality follows from the structure of D i
l
, and the inclusion
is by the facts that x 1), and the convexity of the
set C. As a result, we have that
l +D i
l  i
l
l (y i
l
l
l d i
l d i
l D i
l e l )
where the last inequality follows from f 2 C 1;1
A.24]).
Parallel Variable Distribution 17
By the construction of the algorithm and Lemma 2,
Hence, there exists some  i 2 < n such that
Dene the (continuous) residual function
With this notation,
By properties of the projection operator (e.g., see [1, Proposition B.11]),
since x i 2 C we have that
Hence,
Combining the latter relation with (3.21) and (3.22), and taking into
account the synchronization step of Algorithm 3, we further obtain
Observe that
where the rst relation is by the Cauchy-Schwarz inequality, the second
follows from (3.22), and the last from (3.20). Hence,
Using (3.23), the latter inequality, and again (3.20), we have that
By the choice of  and assumptions on parameters  1 and  2 , the
sequence ff(x i )g is nonincreasing. If f() is bounded below on C, then
bounded below, and hence it converges. In the latter case,
Sagastizabal and Solodov
(3.24). On the other
hand,
so that
We conclude that fR(x i )g also tends to zero. By the continuity of R(),
it then follows that for every accumulation point
x of the sequence
It is well known that the latter is equivalent to the the
minimum principle necessary optimality condition
xi  0 for all x 2 C.
4. Preliminary Numerical Experience
To get some insight into computational properties of our approach
in x 2, we considered test problems taken from the study of Sphere
Packing Problems [4]. In particular, the problems chosen are the same
as the ones used in [16]. Not having access to a parallel computer, we
have carried out a simulation, i.e., the subproblems are solved serially
on a serial machine. Even though this is admittedly a rather crude
experiment, it nevertheless gives some idea on what one might expect
in actual parallel implementation.
Given p spheres in <  with so that x 2 < p , the
problems are as follows.
Problem 1.
min
x (i 1)+t x (j 1)+t
s.t.
(i 1)+t
Problem 2. Given an integer   1,
min
s.t.
Parallel Variable Distribution 19
Problem 3.
min
s.t.
We have implemented in Matlab our Algorithm 2, the serial SQP
method, and the general PVD Algorithm 1. All codes were run under
Matlab version 6.0.0.88 (Release 12) on a Sun UltraSPARCstation. Details
of the implementation are as follows. Algorithm 1 is implemented
using the function constr.m from the Matlab Optimization Toolbox
for solving the subproblems (P l ). Algorithm 2 and the serial SQP
method are quite sophisticated quasi-Newton implementations with
line-search based on the merit function (2.5). Each (QP l ) is solved by
the Null-Space Method [8]. The penalty parameter  i is updated using a
modication of the Mayne and Polak rule [17] that allows  i to decrease,
if warranted. The stepsize t i is computed with an Armijo rule that uses
safeguarded quadratic interpolation. In order to prevent the Maratos'
eect (i.e., to ensure that the unit stepsize is asymptotically accepted
if possible), a second-order correction is added to the search direction
necessary [3, Ch. 13.4]. Finally, the quasi-Newton matrices
are updated by the BFGS formula, using also the Powell correction
and the Oren-Luenberger scaling. The methods stop when the relative
error, measured by the sum of the norm of the reduced gradient and
the norm of the constraints, is less than 10 5 .
We rst compare our Algorithm 2 with the general PVD Algorithm
1. Since this comparison turned out rather obvious and one-sided
(which is certainly not surprising), we do not report exhaustive testing
for this part. In Table II we report the number of iterations and the
running times for Algorithms 1 and 2 on Problem 1 (results for other
problems follow the same trend and do not yield any further insight).
Note that the running times reported are serial, without any regard to
the parallel nature of the algorithms. Because the iterative structure of
the two algorithms is the same, this seems to be a meaningful and fair
comparison. It is already clear from intuition that solving exactly the
general nonlinear subproblems (P l ) in Algorithm 1 is very costly, and
is unlikely to yield a competitive algorithm. This was easily conrmed
by our experience. Of course, one could use heuristic considerations
to (somehow) adjust dynamically the tolerance in solving the sub-problems
(in our experiment, all subproblems are solved to within the
20 Sagastizabal and Solodov

Table

I. Results for Algorithm 2 and the SQP method on sphere packing prob-
lems. Starting points are generated randomly, with coordinates in [ 10; 10]. Times of
solving QPs are measured (in seconds) using the intrinsic Matlab function cputime.
Problem p QP iter calls to QP \speedup" iter calls to
26 28
28 59
Prob.3
Parallel Variable Distribution 21

Table

II. Results for Algorithms 1 and 2 on
Problem 1. Starting points are feasible, generated
randomly. Times are measured (in seconds) using the
intrinsic Matlab function cputime.
Algorithm 1 Algorithm 2
Name p iter time iter time
tolerance of 10 5 , same as the original problem). However, there is no
convergence analysis to support such a strategy within Algorithm 1 in
the constrained case (strictly speaking, there is no even a convergence
proof for Algorithm 1 in the case under consideration, since the feasible
set is nonconvex!). As discussed above, one of the motivations for our
Algorithm 2 is precisely to provide a constructive implementable way of
solving subproblems (P l ) inexactly, by solving their quadratic approx-
imations. The results in Table II conrm that the proposed approach
certainly makes sense.
Our next set of experiments concerns with the comparison between
Algorithm 2 and serial SQP. It is not easy to come up with a meaningful
comparison of a serial method and a parallel method on a serial
machine. For example, the overall running time of a \parallel" method
obtained on a serial machine is considered notoriously unreliable to
predict time that would be required by its parallel implementation
(i.e., just dividing the time by the number of \processors" is not a
good measure). We therefore focus on indicators which we feel should
be still meaningful for the actual parallel implementation, as discussed
below. To evaluate the gain in computing the search directions, we
report the total time spent solving quadratic programming subproblems
within the SQP and within the serial implementation of Algorithm
2. Since no communication between the processors would be needed
within this phase, the \expected" time for solving QPs by the parallel
implementation of Algorithm 2 can indeed be obtained dividing by the
number of processors. It is somewhat surprising that even for smaller
problems, see Table I, the serial time for solving QPs in Algorithm 2 is
already considerably smaller than in standard SQP, with the dierence
22 Sagastizabal and Solodov
becoming drastic for larger problems. If we approximate the \speedup"
e-ciency for computing directions in parallel by
time spent solving QPs in SQP
time spent solving QPs in serial Algorithm 2  100 ;
then these e-ciencies vary from acceptable (around 80%) to very high
(over 10000%), and grow fast with the size of the problem. This conrms
our motivation as discussed in x 2, i.e., smaller QPs are signi-
cantly easier to solve. Obviously, the possible price to pay in Algorithm
2 is \deterioration" of directions compared to a good implementation
of the full SQP algorithm. This issue is addressed by reporting the
numbers of iterations and calls to the oracle evaluating the function
and derivatives values, see Table I. We note that the numbers of iterations
and function/derivatives evaluations are usually slightly higher
for Algorithm 2, but not always. Overall, those numbers for the two
methods are quite similar. Given the impressive gain that Algorithm
exhibits in solving QP subproblems, this indicates that the parallel
implementation of this algorithm should indeed be e-cient, at least
when computing the function and derivatives is cheap relative to solving
QPs. In particular, this should be the case for large-scale problems
where the functions and their derivatives are given explicitly.
5. Concluding remarks
Two new parallel constrained optimization algorithms based on the
variables distribution (PVD) have been presented. The rst one consists
of a parallel sequential quadratic programming approach for the
case of block-separable constraints. This is the rst PVD-type method
whose convergence has been established for nonconvex feasible sets. The
second proposed algorithm employs approximate projected gradient
directions for the case of general (inseparable) convex constraints. The
use of inexact directions is of particular relevance when the projection
operation is computationally costly.



--R

Nonlinear programming.
Parallel and Distributed Computation.

Sphere Packings
The Linear Complementarity Problem.







Nonlinear Programming.









--TR
Iterative methods for large convex quadratic programs: a survey
Sphere-packings, lattices, and groups
Parallel and distributed computation: numerical methods
Dual coordinate ascent methods for non-strictly convex minimization
Parallel Gradient Distribution in Unconstrained Optimization
New Inexact Parallel Variable Distribution Algorithms
bounds in mathematical programming
Two-phase model algorithm with global convergence for nonlinear programming
Testing Parallel Variable Transformation
Parallel Synchronous and Asynchronous Space-Decomposition Algorithms for Large-Scale Minimization Problems
On the Convergence of Constrained Parallel Variable Distribution Algorithms
Parallel Variable Transformation in Unconstrained Optimization
projected gradient;constrained optimization;sequential quadratic programming;parallel optimization;variable distribution