A Variable Metric Proximal Point Algorithm for Monotone Operators.
The proximal point algorithm (PPA) is a method for solving inclusions of the form $0\in T(z)$, where T is a monotone operator on a Hilbert space. The algorithm is one of the most powerful and versatile solution techniques for solving variational inequalities, convex programs, and convex-concave mini-max problems. It possesses a robust convergence theory for very general problem classes and is the basis for a wide variety of decomposition methods called splitting methods. Yet the classical PPA typically exhibits slow convergence in many applications. For this reason, acceleration methods for the PPA algorithm are of great practical importance. In this paper we propose a variable metric implementation of the proximal point algorithm. In essence, the method is a Newton-like scheme applied to the Moreau--Yosida resolvent of the operator T. In this article, we establish the global and linear convergence of the proposed method. In addition, we characterize the superlinear convergence of the method.  In a companion work, we establish the superlinear convergence of the method when implemented with Broyden updating (the nonsymmetric case) and BFGS updating (the symmetric case).
Introduction
The Proximal Point Algorithm (PPA) is one of the most powerful and versatile solution
techniques for problems of convex programming and mini-max convex-concave program-
ming. It possesses a robust convergence theory for very general problem classes in finite-
and infinite-dimensions (e.g. see [11, 16, 21, 22, 23, 28, 32, 41, 40]), and is the basis for a
wide variety of decomposition methods called splitting methods (e.g. see [4, 9, 12, 43, 44]).
Yet, the classical PPA typically exhibits slow convergence in many applications. For this
reason, acceleration methods for the PPA are of great practical importance. In this paper
we propose a variable metric implementation of the proximal point algorithm. Our approach
extends and refines results that originally appeared in [38] and is in the spirit of
several recent articles [3, 7, 10, 18, 20, 24, 25, 36]. However, there is a fundamental difference
between the method presented here and those studied in [3, 7, 10, 18, 20, 24, 25, 36].
This difference has a profound impact on the methodology applied in this article. All previous
work on this topic (except [38]) applies exclusively to monotone operators that arise
as the subdifferential of a finite-valued, finite-dimensional convex function. The results of
this article apply to general monotone operators on a Hilbert space. The resulting difference
in methodology roughly corresponding to the difference between methods for function
minimization and methods for solving systems of equations.
There are both advantages and disadvantages to the more general approach. The advantages
are that the method applies to a much broader class of problems. This is so
not only because the theory is developed in the Hilbert space setting, but, more impor-
tantly, because many monotone operators cannot be represented as the subdifferential of a
finite-valued, finite-dimensional convex function. General monotone operators do not possess
many of the rich structural properties associated with the subdifferential of a convex
function (e.g., subdifferentials of convex functions are the only maximal cyclically monotone
operators [33]). In addition, in the case where the operator is the subdifferential of a
convex function, we do not require the usual assumption that the underlying function be
finite-valued.
The disadvantages of our general approach arise from the fact that the method cannot
make use of the additional structure present when the operator is the subdifferential of a
convex function. This complicates both the structure of the method and its analysis. Of
particular note in this regard is the complexity of our global convergence result. If the
operator is the subdifferential of a convex function, then solving the inclusion
is equivalent to minimizing the underlying convex function. The global convergence of a
method is then typically driven by a line-search routine (e.g., see [3, 7, 10, 18, 20, 24,
25, 36]). In the general setting we do not have direct recourse to this strategy. This
complicates both the structure of the algorithm and its convergence theory. Nonetheless,
the proof technique developed in this paper can be refined in the convex programming
setting, thereby significantly simplifying both the global and local convergence results [5, 6].
Notwithstanding these differences in methodology, our approach is still nicely motivated
by recalling the behavior of the PPA in the context of convex programming:
min
where H is a Hilbert space and f : H 7! IR [ f+1g is a lower semi-continuous convex
function that is not identically +1. Define the Moreau-Yosida regularization of f to be
the function f - : H 7! IR given by
f -z) := min
The set of solutions to (1) corresponds precisely to the set of points at which f - attains
its minimum value. The function f - is continuously Fr'echet differentiable [28, Proposition
7.d]. The PPA applied to (1) is approximately the steepest decent algorithm applied to f -
[11]. This analogy immediately suggests that a variable metric approach could be applied
to the function f - to accelerate the method. This idea was first studied in [38] and is the
basis of the acceleration techniques described in [3, 7, 10, 18, 20, 24, 25, 36].
In [3], Bonnans, Gilbert, Lemar'echal, and Sagastiz'abal develop methods along an algorithmic
pattern originally suggested by Qian in [38]. This pattern circumvents many of the
difficulties associated with a variable metric approach applied directly to the function f - .
The key is to employ a matrix secant update based on the function f instead of f - . The
local convergence results in [3, Section 3] require some smoothness assumptions. In partic-
ular, linear convergence is established when the function f is differentiable with Lipschitz
continuous derivative, and super-linear convergence is established when f is twice strictly
Fr'echet differentiable at a unique solution -
z where the second derivative is positive definite
(we only speak of quotient or q-rate of convergence).
In [18, 20, 24, 25], the authors apply the bundle concept for nonsmooth convex minimization
[17] to approximate the Moreau-Yosida regularization f - and its derivative. Variable
metric updates, in particular, quasi-Newton updates, are then applied using these
approximate values. The super-linear convergence results in the papers [18, 20, 24] either
require strong smoothness assumptions on the function f (such as the Lipschitz continuity
of rf) or that the regularization parameter - diverges to +1. In [20], Lemar'echal, and
Sagastiz'abal propose a clever reversal quasi-Newton formula which uses the value of the
gradient of f - at a variety of points other than those strictly obtained by the iterates. This
promising idea deserves further theoretical and numerical study.
In [10] and [36], the authors develop an approach based on Newton's method for semi-
smooth functions as developed in [30, 31, 37, 34]. Properly speaking, these methods are
neither an adaptation of the PPA algorithm nor a variable metric method. Nonetheless, the
flavor of both these methodologies are present. In order to obtain super-linear convergence,
smoothness hypotheses are again required, however, these hypotheses are of a somewhat
more technical nature. Specifically, it is required that
(a) the function f be semi-smooth at a unique solution to (1) [37],
(b) every element of the set-valued mapping
exists for all
is nonsingular at the unique solution - z, and
(c) the sequence of Hessian approximates fV k g used to generate the iterates fz k g must
lim
dist
One can show that the semi-smoothness hypotheses is satisfied in many cases of interest
when f is finite-valued. Moreover, by Rademacher's theorem on the differentiability of
Lipschitz continuous functions, it follows that the set-valued mapping @ 2
B f(z) is always
well-defined and compact-valued in the finite dimensional, finite-valued case with the
non-singularity property being closely tied to the usual hypothesis of strong convexity.
Although the limiting hypotheses on the V k 's is a bit strong, it is not entirely unreasonable
in the absence of differentiability. In [36], Chen and Qi propose a very nice preconditioning
technique wherein an exact value for the gradient of a shifted Moreau-Yosida regularization
can be computed from inexact values for the gradient of f - . This technique is similar in
spirit to the reversal quasi-Newton formula found in [20]. Both of these techniques should
prove useful in numerical implementations.
The algorithm presented in this paper is most closely related to the methods proposed
by Chen and Fukushima [7] and Mifflin, Sun, and Qi [25]. However, there are several
fundamental distinctions. The foremost of which is that the methods in [7, 25] are restricted
to finite dimensional finite-valued convex programming problems. Within this
framework, these authors use bundle strategies to approximate f - and its gradient and establish
the global convergence of their methods with the aid of a line search routine. Chen
and Fukushima establish global and linear convergence results along with a generalization
of the Dennis-Mor'e characterization theorem for super-linear convergence [14]. One of
the most important features of the Chen-Fukushima algorithm is that the line search is
based on the function f rather than approximations to the function f - . This is very important
in practise since obtaining sufficiently accurate approximations to the function f -
is usually quite time consuming. Their linear and super-linear convergence results blend
bundle techniques with the theory of nonsmooth equations. Consequently, the convergence
hypotheses are reminiscent of those employed in [10] and [36], in particular, they
require semi-smoothness, CD-regularity, and the strong approximation property (2). In
[6], the methods of this paper are applied to the Chen-Fukushima algorithm to obtain the
super-linear convergence of the method when BFGS matrix secant updating is employed.
In [25], Mifflin, Sun, and Qi obtain the first super-linear convergence result for a variable
metric proximal point algorithm using the BFGS matrix secant update in the setting of
finite dimensional finite-valued convex programming. Their proposed algorithm uses a
line search based on approximations to the function f - and requires that the function f -
is strongly convex with rf - Fr'echet differentiable at the unique global solution to the
convex program. In addition it is assumed that the iterates satisfy a certain approximation
property involving the gradient rf - . In Section 4 of this paper, we discuss how these
hypotheses are related to those that are also required in our convergence analysis.
In this paper, we provide a general theory for a variable metric proximal point algorithm
applied to maximal monotone operators from a Hilbert space to itself. In the
important special case of convex programming, where T is taken to be the subdifferential of
the function f , we do not assume that f is finite-valued or differentiable on the whole space.
However, to obtain super-linear convergence, we do require certain smoothness hypotheses
at a unique global solution -
z. These smoothness hypotheses differ from those assumed in
[3, 18, 20, 24] since they are imposed on the operator T \Gamma1 rather than T . In this regard, they
are reminiscent of the hypotheses employed in [25]. The choice of smoothness hypotheses
has deep significance in the context of convex programming. Differentiability hypotheses on
imply the second-order differentiability of f , whereas differentiability hypotheses
are related to the standard strong second-order sufficiency conditions of
convex programming [40, Proposition 2] and thus reduce to the standard hypotheses used
in local analysis of convergence. In particular, the differentiability of (@f) \Gamma1 does not imply
that @f is single-valued or differentiable nor does it imply that f is finite-valued.
Our smoothness hypotheses also differ from those that appear in [7, 10, 36]. These
methods rely on the theory of nonsmooth equations and require hypotheses such as semi-
smoothness and non-singularity of the elements of @ 2
In addition, the proof theory for
these methods specifically requires that the underlying convex function be finite-valued in a
neighborhood of the unique solution to (1) (again, these methods assume that the function
is finite-valued on all of IR n ). This limits direct application to constrained problems since
in the constrained case solutions typically lie on the boundary of the constraint region (i.e.,
on the boundary of the domain of the essential objective function).
Throughout the paper we illustrate many of the ideas and results by applying them to
the case of convex programming. Our purpose here is not only to show how the results can
be applied, but also to ground them in the familiar surroundings of this concrete application.
Further details on the application of these results to the case of convex programming can
be found in [5].
The paper is structured as follows. We begin with a review of the classic proximal point
algorithm in x2. The VMPPA is introduced in x3. This section contains the approximation
criteria that must be satisfied at each iteration. Two criteria are presented. The first
is required to obtain global convergence and the second is required to accelerate the local
convergence of the method. This division into global and local criteria is one of the recurring
themes of the paper. On the global level the method behaves like a steepest descent method
while at the local level it becomes more Newton like. This feature is common to most general
purpose methods in nonlinear programming such as the non-monotone descent methods,
the dogleg method, and trust-region methods. In x4 we discuss the smoothness hypotheses
required for the local analysis. We also extend some of the differentiability results appearing
in [19, 35] to maximal monotone operators. In x5, we study the operators N k associated
with the Newton-like iteration proposed in x3. The focus of this section is to provide
conditions under which the operators N k are non-expansive at a solution to the inclusion
global convergence result paralleling Rockafellar's 1976 result [41] is given in x6.
In x7 we study local convergence rates. Linear convergence is established under a Lipschitz
continuity assumption on T \Gamma1 , and a characterization of super-linear convergence for the
VMPPA is also given. This characterization is modeled on the landmark characterization
of super-linear convergence of variable metric methods in nonlinear programming due to
Dennis and Mor'e [14]. In [6], we use this characterization result to establish the super-linear
convergence of the method when the derivatives are approximated using the BFGS
and Broyden updating strategies.
A word about our notation is in order. We denote the closed unit ball in the Hilbert
space H by IB. Then the ball with center a and radius r is denoted by a+ rIB. Given a set
and an element z 2 H, the distance of z to Z is dist (z;
be two Hilbert spaces. Given a multi-function (also referred to as a
mapping or an operator depending on the context) , the graph of T , gph T , is
the subset of the product space H 1 \Theta H 2 defined by gph (z)g.
The domain of T is the set domT := fz ;g. The identity mapping will be
denoted by I. The inverse of an operator T is defined by T
gph Tg.
Given a lower semi-continuous convex function f
S f+1g, the conjugate of f
is defined by f   (z   f(z)g.
Monotone Operators and the Classic Algorithm
Given a real Hilbert space H with inner product h\Delta; \Deltai, we say that the multi-function
for every z and z 0 in domT , and w 2 T (z) and w
have is said to be strongly
monotone with modulus -. The monotone operator T is said to be maximal if its graph
is not properly contained in the graph of any other monotone operator. An important
example of a monotone operator is the subgradient of a convex function (see Minty [27] and
[28]).
We are concerned with solving inclusions of the form
where T is a maximal monotone operator. In the case of the convex programming problem
(1), the operator T is the subdifferential of the convex function f and the inclusion
(3) characterizes the points z at which f attains its minimum value. A wide variety of
other problems can be cast in this framework, e.g. variational inequalities, complementary
problems, and mini-max problems. Existence results for inclusion (3) can be found in [41].
In 1962, Minty [27] showed that, when the operator T is maximal monotone, the Moreau-
Yosida Resolvent of T ,
is single-valued and non-expansive on H. This result suggests that a solution to the inclusion
can be iteratively approximated by the recursion z
can modify this scheme by varying the scalar - and by choosing the iterates z k+1 to be an
approximate solution to the equation . The proximal point algorithm applies
precisely these ideas. The algorithm, starting from any point z 0 , generates a sequence
in H by the approximation rule
z
The principle difficulty in applying the proximal point algorithm lies in executing the
operators In the case of convex programming, the iteration (4) reduces
to the iteration
z
Notice that executing the algorithm exactly (i.e., with "=" instead of "-" in the above
algorithm) can be as difficult as solving the original problem directly. Hence it is critical
that the convergence results are obtained under the assumption of approximation.
In [22] and [23], Martinet proved the convergence of the exact proximal point algorithm
for certain cases of the operator T with fixed c k j c. The first theorem on the convergence
of the general proximal point algorithm was proved by Rockafellar [41] in 1976. His theorem
not only insures the global convergence under an approximating rule, but also describes
the global behavior when the inclusion 0 2 T (z) has no solution.
The convergence rate of the proximal point algorithm depends on properties of the
operator T , the choice of the sequence fc k g, and the accuracy of the approximation in (4).
The first rate of convergence results were also obtained by Rockafellar [41] in 1976, under
the assumption that the solution set is a singleton f-zg. He proved that if the sequence
is bounded away from 0, and T \Gamma1 (w) is bounded by a linear function of kwk when w
is near 0, then the rate of convergence is at least linear. Luque [21] extended Rockafellar's
theorem to the case where T \Gamma1 (0) is not required to be a singleton, and showed that such
an estimate of the convergence rate is tight.
3 The Algorithm and Approximation Criteria
The algorithm proposed in this section is a Newton-like iteration for solving the resolvent
equation In the context of the convex programming problem, the iteration takes
the form
z
where the operator H k is used to approximate second-order properties of the function f - .
If f - is twice differentiable with [r method one sets
However, in general, f - is only known to be differentiable with Lipschitz
continuous gradient [28]. Thus, in the finite dimensional case, the Hessian r 2 f - (x) is only
guaranteed to exist on a dense subset by Rademacher's Theorem. Further results on the
second-order properties of f - can be found in [19, 35, 42].
It is well known that the negative gradient \Gammarf - (z k ) is the unique element w k solving
the problem
min
or equivalently, satisfying the inclusion
The proximal point algorithm for a general maximal monotone operator T can be formally
derived from equation (5) by replacing -, z k , and @f by c k , z k , and T respectively, to obtain
or equivalently,
where equality follows from the fact that w k is unique. This motivates us to define the
operator
This operator provides the analog of the direction of steepest descent in the operator setting.
The algorithm we propose for solving the inclusion 0 2 T (z) can be succinctly stated
as follows:
The Variable Metric Proximal Point Algorithm:
Having z k , set
z
and choose c k+1 - 1.
As mentioned in the previous section, it is critical that the convergence results are
obtained under the assumption that D k (z k ) can only be approximated. We use the following
approximation criteria:
and
The approximation criteria (G) is used to establish global convergence properties, while
criteria (L) is used to obtain local rates of convergence.
Although these criteria are used in the proof of convergence, they are impractical from
the perspective of implementation. In their stead, we provide criteria that are imple-
mentable. To obtain these criteria we recall the following result from Rockafellar [41].
Proposition 1 [41, Proposition 3] Let S k (w) := T (z k
we have the bound
dist (0; S k
Proposition 1 yields the following alternative approximation criteria for the w k 's. Since
this result is an immediate consequence of Proposition 1, its proof is omitted.
Proposition 2 Consider the following acceptance criteria for the w k 's:
dist (0; S k (w k
dist (0; S k (w k
We have
Remark Note that to satisfy either (G 0 ) or (L 0 ) it is not necessary to find an element of
least norm.
Before leaving this section we recall from [41] a few properties of the operators D k and
I that are essential in the analysis to follow.
Proposition 3 [41, Proposition 1]
a) The operator D k can be expressed as
and for any z
b) For any z; z
c) For any z; z
Remark An important consequence of Part c) above is that the operators P k and D k are
Lipschitz continuous with Lipschitz constant 1, that is, they are non-expansive. Henceforth,
we use of this fact.
4 On the Differentiability of T \Gamma1 and D k
Just as Newton's method for minimization locates roots of the gradient, one can view
the variable metric proximal point algorithm as a Newton-like method for locating roots
of the operator D k . This perspective motivates our approach to the local convergence
analysis. For this analysis, we require that the operator T \Gamma1 possesses certain smoothness
properties. These properties in turn imply the smoothness of the operators D k . Smoothness
hypotheses are used in the convergence analysis in much the same way as they are used
in the convergence analysis for Newton's method. For example, recall that to ensure the
quadratic convergence of Newton's method one requires the derivative at a solution to be
both locally Lipschitz and non-singular. Non-singularity insures that the iterates are well
defined and can be bounded, while the Lipschitzian hypothesis guarantees that the error
in the linearization is quadratically bounded (see [29, 3.2.12 and 10.2.2]). We make use of
similar properties in our analysis.
In order to discuss the smoothness of T \Gamma1 and D k , we recall various notions of differentiability
for multi-valued functions from the literature. For a more thorough treatment of
these ideas in the context of monotone operators, we refer the reader to [1, 19, 26, 35, 42].
Definition 4 We say that an operator continuous at a point -
(with modulus ff - 0 ) if the set \Psi( -
w) is nonempty and there is a - ? 0 such that
We say that \Psi is differentiable at a point -
w) consists of a single element -
z and there
is a continuous linear transformation J H such that for some ffi ? 0,
We then write
w).
Remarks 1) These definitions of Lipschitz continuity and differentiability for multifunction
are taken from [41, pages 885 and 887] (also see [2, page 41]). Note that these
notions of Lipschitz continuity and differentiability correspond to the usual notions
when \Psi is single-valued.
Rockafellar [41, Theorem 2] was the first to use Lipschitz continuity to establish rates
of convergence for the proximal point algorithm.
When the set \Psi( -
w) is restricted to be a singleton f-zg, the differentiability of \Psi at -
implies the Lipschitz continuity of \Psi at -
Moreover, one can take ff(-
This observation is verified in [41, Proposition 4].
It follows from the definition of monotonicity that if T is a maximal monotone operator,
then the operator rT (x) is positive semi-definite whenever it exists.
We now give a result that relates the differentiability of a multi-valued function to the
differentiability of its inverse. The proof is omitted since it parallels the proof of a similar
result for single-valued functions.
Lemma 5 Assume that \Psi : H \Gamma! \Gamma! H is differentiable at -
z with
wg and
with J \Gamma1 bounded. Also assume that \Psi \Gamma1 is Lipschitz continuous at -
w with \Psi
f-zg.
is differentiable at -
w with r\Psi
In the two examples that follow, we examine the concepts introduced in Definition 4
when the operator in question is the subdifferential of a convex function. The first example
illustrates that @f \Gamma1 can be Lipschitz continuous but not differentiable at the origin, while
in the second example @f \Gamma1 is differentiable at the origin but @f is not differentiable on
Example 6 Let
z if z - 0 and T (z) :=
T \Gamma1 is Lipschitz continuous at 0 but is not differentiable at 0.
Example 7 Let
\Gammaz if z ! 0
z 5=3 if z - 0 and T (z) :=
z 2=3 if z ? 0 .
y 3=2 if y - 0 .
T \Gamma1 is differentiable at 0 with but T is not differentiable on T \Gamma1 (0).
The super-linear convergence result of x7 requires the assumption that the operator T \Gamma1
be differentiable at the origin. Although this is a severe restriction on the applicability of
these results, it turns out that in the case of convex programming it is a consequence of the
standard second-order sufficiency conditions for constrained mathematical programs. This
and related results were established by Rockafellar in [40, Proposition 2]. In this context, it
is important to note that the second-order sufficiency condition is the standard hypothesis
used in the mathematical programming literature to insure the rapid local convergence of
numerical methods. So, at least in the context of constrained convex programming, such
a differentiability hypothesis is not as severe an assumption as one might at first suspect.
To the contrary, it is a bit weaker than the standard hypothesis employed for such results.
For the sake of completeness, we recall a portion of Rockafellar's result below.
Theorem 8 Consider the convex programming problem (1)
given by
otherwise,
with m. Suppose that the following conditions are
(i) The functions f i for are k - 2 times continuously differentiable in a
neighborhood of a point -
z 2 IR n .
(ii) There is a Kuhn-Tucker vector -
z such that -
(iii) The gradients frf are linearly independent.
(iv) The matrix
Then the operator @f \Gamma1 is continuously differentiable in a neighborhood of the
origin.
Remark Theorem 8 follows by applying the implicit function theorem to the Kuhn-Tucker
conditions for the parameterized problems minff(z) \Gamma hw; zig in a neighborhood of
The relationship to @f \Gamma1 comes from the fact that @f zig.
Rockafellar only establishes the result for 2. The extension to k ? 2 follows trivially
from the implicit function theorem.
We now examine the differentiability properties of the mapping D k . Two results in
this direction are given. The first uses equation (8) to relate the differentiability of the
operators T \Gamma1 and D k , while the second uses the definition of D k given in (6) to relate the
differentiability of the operators T and D k .
Proposition 9 Let
w. The operator T \Gamma1 is differentiable at -
y
with
bounded if and only if the operator D is differentiable at - z with
In either case, we have
Proof First assume that T \Gamma1 is differentiable at -
y with r(T \Gamma1 )(-y) bounded. The differentiability
of T \Gamma1 at -
y clearly implies that of D \Gamma1 at -
w with
Since D is Lipschitzian with
implies that D is differentiable at -
z with
derivative given by (10). Since r[D
we conclude that the latter is
bounded.
Conversely, assume that D is differentiable at -
z with (rD(-z)) \Gamma1 bounded. We show
that D \Gamma1 is single-valued and Lipschitzian at -
w. The result will then follow from Lemma
5.
be as in Definition 4 for rD(-z). Since D is single-valued and rD(-z) surjective
(it is invertible), we may apply a standard open mapping result from functional analysis
(e.g. [8, Theorem 15.5]) to obtain the existence of a ae ? 0 and a
Hence for each w 2 -
is bounded, there is a - ? 0 such that
Hence, by reducing ae and -
if necessary, we may assume that
wk
for
aeIB, where the second inequality follows since D is non-expansive. Therefore,
we can assume that o(kz \Gamma -
w+aeIB and z 2 D
ffiIB).
By substituting this into (12) and re-arranging, we obtain
w+ aeIB and z 2 D
We now show that (13) implies the existence of an ffl ? 0 such that D
ffiIB. Indeed, if this were not the case, then there would exist sequences fw i g and fz i g
such that z i
its images are convex, hence, by (11), there exists a sequence f-z i g with -
implies that
for all This contradicts the fact that w
w and k-z
and so such an ffl ? 0 must exist. This fact combined with (13) implies that
D \Gamma1 is Lipschitzian at -
w with D
now applies to yield the result.
be defined as in (9). Let -
D)(-z). The
operator T is differentiable at - y with [I bounded if and only if the operator D
is differentiable at -
z with [I +rD(-z)] \Gamma1 bounded. In either case we have the formula
Proof Replace D by P := I observe that D is differentiable at -
z
with [I +rD(-z)] \Gamma1 bounded if and only if P is differentiable at -
z with [rP (-z)] \Gamma1 bounded.
The proof now follows the same argument as in the proof of Proposition 9 with D replaced
by replaced by T and -
w replaced by - y.
Propositions 9 and 10 say quite different things about the differentiability of D k . To
illustrate this difference, observe that in Example 7, the operator T is not differentiable
at 0, while T \Gamma1 and D are differentiable at 0. On the other hand, if we take
with not differentiable at 0, while T and D are differentiable at
It is also important to note that even if neither T nor T \Gamma1 is differentiable, D may
be differentiable. But, in this case, we know from Propositions 9 and 10, that if D is
differentiable and neither T nor T \Gamma1 is differentiable, then both rD(-z) and rP (-z) have
to be singular or have unbounded inverses. For a further discussion of these issues in the
context of finite dimensional convex programming see [35].
When T is assumed to be the subdifferential of a convex function f , Propositions 9 and
can be refined by making use of the relation @f  is the convex conjugate
of f [39, Corollary 12A]. This allows us to extend [35, Theorem 1] and [35, Theorem 2] to the
Hilbert space setting (also see [19, Theorem 3.1]). However, some caution in terminology
is required since f   is not necessarily twice differentiable in the classical sense at points
where @f   is differentiable in the sense of Definition 4. Indeed, @f   may be multi-valued
arbitrarily close to a point of differentiability. The best way to interpret this result is
through Alexandrov's Theorem [1] which states that at almost every point -
z in the interior
of the domain of a convex function f : IR n 7! IR [f1g there is a quadratic function q - z such
that
In [19] and [35], the matrix r 2 q -
z is called a generalized
Hessian and is denoted Hf(x). Note that the existence of a generalized Hessian at the
point -
z guarantees that f is strictly differentiable at -
z. Moreover, if @f(x) is single-valued
in a neighborhood of a point -
z at which Hf(-z) exists, then r 2 f(-z) exists and equals Hf(-z).
We extend this terminology to the Hilbert space setting with the following definition.
Definition 11 Let OE: H 7! IR [ f1g be a function on the Hilbert space H. We say that
OE is twice differentiable in the generalized sense at a point -
there is a continuous
quadratic functional q - z such that
z is called
a generalized Hessian of OE at -
z and is denoted by HOE(-z).
With this terminology in hand, we apply Propositions 9 and 10 to the case of convex
programming. The proofs of these results are not required since they are a direct translation
of Propositions 9 and 10 into the terminology of convex programming.
Corollary
S f+1g be lower semi-continuous and convex. Let -
and set -
differentiable at -
z with
bounded if and only if f   has a generalized Hessian at -
y with [I
bounded. In either case we have
S f+1g be lower semi-continuous and convex. Let -
and set -
differentiable at -
z with [I +r 2 f -z)] \Gamma1
bounded if and only if f is twice differentiable in the generalized sense at -
y with [I
In either case we have
Remark As observed earlier, the generalized Hessian is necessarily positive semi-definite.
This observation can be used to further refine the statement of Corollaries 12 and 13.
5 Newton Operators
In this section we study the operators associated with the variable metric proximal point
iteration:
This notation emphasizes the fact that these operators produce Newton-like iterates. Just
as in the case of the classical Newton's method for equation solving [29, x12.6], one of
the keys to the convergence analysis is to show that these operators are contractive with
respect to the solution set T \Gamma1 (0). Clearly the operators N k are single-valued. Moreover,
fixed points of the operators N k are solutions to the inclusion 0 2 T (z) since
Thus, conditions that ensure that the operators N k are non-expansive with respect to
are important for the global analysis of the variable metric proximal point iteration.
To obtain this property, we impose the following conditions on the linear transformations
g.
Each H k is a continuous linear transformation with continuous inverse.
(H2) There is a nonempty closed bounded subset \Gamma of T \Gamma1 (0) such that
where
Remark The set \Gamma in (H2) is used to guarantee the boundedness of the sequence fz k g. By
taking one can show that every weak cluster point of the sequence fz k g is an
element of T \Gamma1 (0). It was observed by Iusem [13] that if T \Gamma1 (0) is bounded and one takes
then the sequence fz k g has a weak limit z Theorem 17 and
[41, Theorem 1]).
Hypothesis (H1) is standard and is automatically satisfied in the finite dimensional case.
On the other hand, hypothesis (H2) is quite technical and requires careful examination.
This hypothesis is problematic since it specifies that the matrices H k satisfy a condition
that depends on the unknown values oe k and kD k (z k )k. We will show that in certain cases
it is possible to satisfy (H2) without direct knowledge of these unknown values. This is
done in two steps. First it is shown in Lemma 14 that if T \Gamma1 is Lipschitz continuous or
differentiable at the origin, then fl k is bounded below by a positive constant (which can be
taken to be 1=6 as kD k (z k )k approaches zero). Then, in Lemma 15, it is shown that (H2)
is satisfied if a related condition in terms of H k and w k is satisfied. Taken together, these
results imply that at least locally (H2) can be satisfied by checking a condition based on
known quantities.
Further insight into hypothesis (H2) can be gained by considering the case in which T \Gamma1
is differentiable at the origin. In this case H k is intended to approximate \Gamma(rD k (0))
k J . Therefore, one can guarantee that (H2) is satisfied by choosing c k sufficiently
large and H k - I. This fact is used in [6] to establish the super-linear convergence
of the method when the H k 's are obtained via matrix secant updating techniques.
The purpose of hypothesis (H2) is to globalize what is essentially a local algorithm
(Newton's method). In the context of convex programming, one commonly obtains global
convergence properties with the aid of a line search routine applied to the objective function
f , or its regularization f - . However, in the operator setting there is no natural underlying
objective function to which a line search can be applied. This is a key difference between
the approach taken in this paper and those in [3, 7, 10, 18, 20, 24, 36]. In the convex
programming setting, the global convergence of the VMPPA is driven by a line search
routine applied to the objective function f (or its regularization f - ). In the operator
setting, hypothesis (H2) replaces the line search and the associated hypotheses needed
to make the line search strategy effective (such as the finite-valuedness of the objective
function f and the boundedness of the sequence fH k g). On the other hand, when it is
known that the operator T is the subdifferential of a finite-valued finite dimensional convex
function, then the algorithm of this paper can be modified to include the line search routine
of Chen and Fukushima [7] thereby avoiding the need for hypothesis (H2) [6].
We now show three cases where the fl k 's are bounded away from zero.
Lemma 14 Suppose T \Gamma1 (0) is nonempty.
(i) If the operator T is strongly monotone with modulus -, then T
and
5+2=-
for all k.
(ii) If the operator T \Gamma1 is Lipschitz continuous at the origin with modulus ff, then
dist
for all k such that kD k (z k )k - where - is given in Definition 4. Moreover, if
5+2ff
for all k such that kD k (z k )k - .
is differentiable at the origin with derivative J , then T there is
a such that for all k with kD k (z k )k - we have
and
for all k, where oe(-
Proof (i) If T is strongly monotone with modulus -, then kz \Gamma z
any z; z single-valued and
continuous. Let
z where
3 (a) we have
Hence
since I. By the definition of fl k ,
This establishes the result since c k - 1 for all k.
(ii) If kD k (z k )k - , Definition 4 implies that
or
hence (15) holds. If T then the lower bound on fl k follows as in Part (i).
(iii) This result follows as in Part (ii) using the second remark after Definition 4.
When w k - D k (z k ), one can establish the inequality in hypothesis (H2) from a related
condition on the vectors w k . A specific technique for accomplishing this is given in the
following lemma.
Lemma be such that
and let H k be a continuous linear transformation from H to itself. If z
then Therefore, if (H1) and criterion (L) are satisfied,
and if - and the sequence f(-fl
defined in (H2)), then hypothesis (H2) is satisfied.
Proof Now from (16) and (17), we have
hence
Again by (17),
since the inequality - fl k - 1=3 implies that -+ 6
- 1.
We conclude this section by showing that the operators N k are non-expansive with
respect to the set T \Gamma1 (0).
Proposition nonempty. If the linear transformations fH k g satisfies
hypotheses (H1) and (H2), then for all k we have kH k D k (z k )k - 3kD k (z k )k and
Proof Let -
z 2 \Gamma. From the definitions of P k and N k , we have
hence
From hypothesis (H2), we have
Hence
Then, again by hypothesis (H2),
Thus, from (20) and (21),
Letting z in Proposition 3 Part (c) yields
From (22) and (23) we have
We now consider ff
. If ff k - fl k, then (18) holds by (24). Suppose that
From (19), we have
Therefore, by (23),
Using the inequality
2a
for a ? b ? 0,
zk
But kD k (z k )k
From (25) we again obtain (18).
6 Global Convergence
The statement and proof of the global convergence result given below parallels the development
given by Rockafellar in [41, Theorem 1] for the classical proximal point algorithm.
Theorem 17 Let fz k g be any sequence generated by the variable metric proximal point
algorithm under criterion (G) (or (G 0 )). Suppose that the solution set T \Gamma1 (0) is nonempty
and the sequence of linear transformations fH k g satisfies the hypotheses (H1) and (H2).
Then the sequence fz k g is bounded, each weak cluster point of this sequence is an element of
0: If it is also assumed that T \Gamma1 (0) is bounded and
in (H2), then there is a - z converges weakly to -
z.
In order to establish this result we require the following technical lemma whose proof
is straightforward and so is omitted.
Lemma Suppose the nonnegative sequences fffl k g satisfies
is a
nonnegative sequence satisfying u is a Cauchy sequence.
Proof of Theorem 17 We begin by showing that the limit lim k kz
exists for
every -
z 2 \Gamma. To this end let -
observe that the definition of N k and Proposition
imply that
Therefore, Lemma implies that the sequence fkz
zkg is Cauchy and so -z) exists for
every -
z 2 \Gamma. An immediate consequence of the existence of these limits is the boundedness
of the sequences fz k g and oe k .
We now show that the sequence fD k (z k )g converges strongly to the origin. Indeed, if this
is not the case, then there is a subsequence J ae
This in turn implies that inf J due to the
boundedness of the sequence foe k g. Let - z 2 \Gamma. By Proposition 16,
z +N k (z k
with fC k g bounded, where the final inequality follows from criterion (G). Hence
whereby we obtain the contradiction
Therefore, lim k kD k (z k
Next let J ae be such that the subsequence fz k g J converges weakly to z 1 ,
i.e. z 1 is a weak cluster point of the sequence fz k g. We show that z 1 must be an element
of T \Gamma1 (0). From Proposition 3 (a), we have that \Gamma 1
or equivalently, hz \Gamma z
all k and z; w with w 2 T (z). Taking the limit over J yields the inequality
for all z; w with w 2 T (z). Since T is maximal monotone, we get
Under the assumption that (0), the argument showing that there is no more
than one weak cluster point of fz k g is identical to the one given by Rockafellar in ([41]
Theorem 1).
Remark To ensure the strong convergence of the sequence fz k g, one again requires a growth
condition on the inverse mapping T \Gamma1 in a neighborhood of the origin. Rockafellar has
shown that Lipschitz continuity of T \Gamma1 at the origin suffices for this purpose [41, Theorem
2]. Other conditions can be found in the work of Luque [21, Proposition 1.2]. The results of
Rockafellar and Luque are easily extended to the variable metric proximal point algorithm.
7 Convergence Rates
7.1 Linear Convergence
Just as in Rockafellar [41, Theorem 2], we require that the operator T \Gamma1 is Lipschitz
continuous at the origin in order to establish that the convergence rate is at least linear.
Theorem 19 Let fz k g be any sequence generated by the variable metric proximal point
algorithm satisfying both criterion (G) and (L) for all k. Assume that T \Gamma1 is Lipschitz
continuous at the origin with modulus ff and the solution set T \Gamma1 (0) is a singleton f-zg.
If the sequence fH k g satisfies the hypotheses (H1) and (H2) with then the
sequence fz k g strongly converges to the solution and there is an index -
k such that
where oe k satisfies limsup k!1 oe k ! 1. That is, the convergence rate is linear.
Proof By Theorem 17, we have kD k (z k )k ! 0. Hence, Part (ii) of Lemma 14 implies that
converges strongly to - z. We now establish the linear rate.
be as in Definition 4 and let ~ k be such that k 1
Proposition 3 (a) and the Lipschitz continuity of T \Gamma1 at 0, we have
Hence relation (14) and hypothesis (H2) yield
zk
Let a k := ff
Using (26) and (27),
Let
. By Proposition 16 and Lemma 14 we have, for k - ~
k, that
By (28) and (29), when k - ~ k
a 2
. From (30) we have
By (31), criterion (L) (or (L 0 )), and Proposition 3 (c),
zk
. Since there is a ~
for any k, and
sufficiently large. Moreover, we have limsup k!1 oe
limsup
ffi.
7.2 Super-linear Convergence
We now give an analog of Dennis and Mor'e 's [14] characterization theorem for the super-linear
convergence of variable metric methods in nonlinear programming that applies to the
VMPPA. This result is used in [6] to establish the super-linear convergence of the variable
metric proximal point algorithm when the Broyden (non-symmetric case) or the BFGS
(symmetric case) updating formulas are used to generate the matrices H k .
Theorem 20 Let fz k g be any sequence generated by the variable metric proximal point
algorithm satisfying criterion (L) for all k. Suppose that the operator T \Gamma1 is differentiable
at the origin with T
converges to the solution - z super-linearly if and only if
Remark By Proposition 9 we have
Consequently, condition (32) can
be recast in the more familiar form given in [15, Theorem 8.2.4]. Note that the assumption
in (32) on the sequence fH k g is much weaker than assuming that this sequence converges.
Specific choices of the linear transformations H k satisfying (32) are discussed in [6].
The proof of Theorem 20 requires the following lemma.
Lemma 21 Under the conditions in Theorem 20 we have
zk)IB,
for all k sufficiently large.
Proof For part (a), let ffi ? 0 be such that
z ae o(kwk)IB (33)
be such that whenever k ? - k 1 , kD k (z k )k - ffi. Then, by (33)
and Proposition 3 (c), when k ? -
We now prove (b). Note that N k (z k hence by criterion (L)
Therefore by (34) and Proposition 3 (c),
Proof of Theorem 20: Let ~
z k+1 := N k (z k
have ~ z
Hence
z
z
or equivalently,
z
z
z
z k+1 ))]
By Lemma 21 the first and third of the three terms appearing on the right hand side of this
inclusion can be bounded by an expression of the form o(kz
zk)IB. If (32) holds, then
Therefore there are positive sequences
fff 1k g and fff 2k g each converging to zero such that, for k ? -
zk
be such that ff 1k ! 1for all k ? - k 2 . Then, denoting ff 1k +ff 2k
by - k ,
zk
converges to -
z super-linearly.
Conversely, suppose that
lim
zk
zk
Divide (35) by kz From (36) and Lemma 21 we obtain
However, from (36) we have
zk
zk
zk
as k !1. Hence (32) holds.
Concluding Remarks
In this paper, we introduced a new proximal point algorithm for solving the inclusion
is an arbitrary maximal monotone operator. The global convergence of
the algorithm is demonstrated with an inexact solution at each step. This is important in
practice, since solving for the exact solution at each step is impractical, and may in fact
be almost as difficult as solving the original problem. If it is assumed that T \Gamma1 is Lipschitz
continuous at the origin, then the method is shown to be linearly convergent. If it is
further assume that T \Gamma1 is differentiable at the origin, then the classical characterization
of super-linear convergence due to Dennis and Mor'e also holds for the VMPPA. In [6],
this characterization of super-linear convergence is applied to establish the super-linear
convergence of the method when certain matrix secant updating strategies are employed
to generate the matrices H k . In [5], we give some of the implementation details in the case
of convex programming. We show how to apply the method to solve the associated primal,
dual, and Lagrangian saddle point problems. In particular, it is shown how the bundle
technique [17] can be applied to satisfy the approximation criteria (L) and (G) in both
the primal and saddle point solution techniques. Preliminary numerical results comparing
these three approaches are also presented.

Acknowledgments

The authors would like to thank the reviewers for their thorough work.
Their comments and suggestions have greatly contributed to our exposition. In particular,
we would like to thank Professor Alfredo Iusem of observing an error in an earlier version
of Theorem 17 and for his suggested revision of this result when the set T \Gamma1 (0) is assumed
to be bounded.



--R

The existence almost everywhere of the second differential of a convex function and some associated properties of convex surfaces.

A family of variable metric proximal point methods.
The method of successive projection for finding a common point of convex sets.
Application of a variable metric proximal point algorithm to convex programming.
On the super-linear convergence of the variable metric proximal point algorithm using Broyden and BFGS matrix secant updating
Proximal quasi-Newton methods for nondifferentiable convex optimization
Nonlinear Functional Analysis.
Splitting Methods for Monotone Operators with Application tp Parallel Optimization.
A globally and superlinearly convergent algorithm for non-smooth convex minimization
New proximal point algorithms for convex minimization.
A decomposition method and its application to convex programming.
Personal communication
A characterization of superlinear convergence and its application to quasi-Newton methods
Numerical Methods for Unconstrained Optimization and Nonlinear Equations.
The proximal points algorithm for reflexive Banach spaces.
Bundle methods in nonsmooth optimization.
An approach to variable metric bundle meth- ods
Practical aspects of the Moreau-Yosida regularization i: Theoretical preliminaries
Variable metric bundle methods: from conceptual to implementable forms.
Asymptotic convergence analysis of the proximal point algorithm.
Regularisation d'inequations variationelles par approximations successive.
Determination approach'ee d'un point fixe d'une application pseudo
A quasi-second-order proximal bundle algorithm

Control dan les inequations variationelles elliptiques.

Proximit'e et dualit'e dans un espace Hilbertien.
Iterative Solution of Nonlinear Equations in Several Variables.
Nonsmooth equations: Motivation and algorithms.
A globally convergent Newton method for SC 1 problems.
Weak convergence theorems for nonexpansive mappings in Banach spaces.
Convex Functions
Convergence analysis of some algorithms for solving nonsmooth equations.

A preconditioning proximal Newton method for nondifferentiable convex optimization.
A nonsmooth version of Newton's method.
The Variable Metric Proximal Point Algorithm: Theory and Application.
Conjugate Duality and Optimization.
Augmented Lagrangians and applications of the proximal point algorithm in convex programming.
Monotone operators and the proximal point algorithm.
Maximal monotone relations and the second derivatives of nonsmooth functions.
Partial inverse of a monotone operator.
Applications of the methods of partial inverses to convex programming: Decomposition.
--TR
variable metric;global convergence;maximal monotone operator;proximal point methods;convergence rates