Cost approximation: A unified framework of descent algorithms for nonlinear programs^{∗}

Michael Patriksson^{†}

August 19, 1997

Abstract. The paper describes and analyzes the cost approximation algorithm. This class of iterative descent algorithms for nonlinear programs and variational inequalities places a large number of algorithms within a common framework and provides a means for analyzing relationships among seemingly unrelated methods. A common property of the methods included in the framework is that their subproblems may be characterized by monotone mappings, which replace an additive part of the original cost mapping in an iterative manner; alternately, a step is taken in the direction obtained in order to reduce the value of a merit function for the original problem. The generality of the framework is illustrated through examples, and the convergence characteristics of the algorithm are analyzed for applications to nondifferentiable optimization. The convergence results are applied to some example methods, demonstrating the strength of the analysis compared to existing results.

Key Words. Nondifferentiable optimization, cost approximation, partial linearization, descent algorithms, convergence analysis

AMS Subject Classification. 90C30, 49M37

Abbreviated Title. Cost approximation algorithms for nonlinear programs

1 Introduction

Let u : ℜ^{n} ↦ ℜ∪{+∞} be lower semicontinuous (l.s.c.), proper and convex and f : ℜ^{n} ↦ ℜ ∪ {+∞} continuously differentiable on an open neighbourhood of dom u. Consider the nondifferentiable optimization problem

[NDP]

min_{x∈ℜ^{n}} T (x) := f(x) + u(x).

This problem is generic in mathematical programming and encompasses the problem of minimizing a convex and/or continuously differentiable real-valued function over a nonempty, closed and convex set in ℜ^{n}.

The most common approach to solving [NDP] is to construct a sequence {x^{t}} of iterates in dom u such that the sequence {T (x^{t})} is strictly monotonically decreasing; typically, the sequence {x^{t}} converges to a solution x^{∗} to the generalized equation ([Rob79])

∗The author wishes to acknowledge the many constructive remarks made by the reviewers, the associate editor, and associate professor Torbjörn Larsson. The research leading to this paper was supported in part by grants from the Swedish Institute (303 GH/MLH) and the Swedish Research Council for Engineering Sciences (TFR) (282-93-1195)

†Department of Mathematics, Chalmers Institute of Technology, S-412 96 Gothenburg, Sweden (mipat@math.chalmers.se). Formerly with Department of Mathematics, Linköping Institute of Technology, S-581 83 Linköping, Sweden

[GE]

∇f(x^{∗}) + ∂u(x^{∗}) ∋ 0,

which, under the mild assumption below (cf. [Roc66]), constitutes the necessary optimality conditions of x^{∗} in [NDP], and further is sufficient for the global optimality of x^{∗} in [NDP] whenever either f is convex, or u is an indicator function of a nonempty, closed and convex set and f is pseudoconvex.

The directional derivative of u in the direction of d is

u′(x; d) := lim_{λ↓0} \frac{1}{λ}[u(x + λd) − u(x)].

Assumption 1.1 (Lower semicontinuity of u′(x; ·)) For any x ∈ dom u and d ∈ ℜ^{n},

u′(x; d) = lim inf_{e→d} u′(x; e). (1.1)

Note that (1.1) is fulfilled for any x ∈ rint (dom u) [Roc70a, Thm. 23.4]. A direct consequence of the assumption is that for any x ∈ dom u and d ∈ ℜ^{n} the directional derivative may be calculated as ([Roc66] and [Roc70a, Thm. 23.2])

u′(x; d) = sup{ ξ_{u}(x)^{T}d | ξ_{u}(x) ∈ ∂u(x) }. (1.2)

We let Ω denote the set of solutions to [GE], and assume that it is nonempty; a sufficient condition for this to be the case is that T is coercive, that is, dom T is bounded or

lim_{∥x∥→∞}{T (x)/∥x∥} = +∞. (1.3)

Moreover, Ω is a singleton set if it is nonempty and if T is strictly convex on its effective domain. Note that strong convexity of T, that is, there exists a constant m_{T} > 0 such that for any ξ_{T} (x) ∈ ∂T (x) and ξ_{T} (y) ∈ ∂T (y), x, y ∈ dom T,

[ξ_{T} (x) − ξ_{T} (y)]^{T}(x − y) ≥ m_{T} ∥x − y∥^{2}, (1.4)

implies both strict convexity and coercivity of T, and hence ensures the existence of a unique solution to [NDP]. (The expression (1.4) actually states that the subdifferential mapping ∂T is strongly monotone on dom T; we may use it as a definition of convexity of T by allowing m_{T} = 0.) Introducing some additional terms, a mapping Π is maximal monotone ([Min62, Bro68]) if it is monotone and its graph is not properly contained in the graph of any other monotone operator; its inverse mapping is denoted by Π^{−1}. Further, a mapping Π is Lipschitz continuous on dom Π (with modulus M_{Π} ≥ 0) if

∥Π(x) − Π(y)∥ ≤ M_{Π}∥x − y∥, ∀x, y ∈ dom Π.

An interesting application of [NDP] results from letting u be the indicator function of a nonempty, closed and convex subset X of ℜ^{n}, that is, u ≡ δ_{X}, where

δ_{X}(x) := \{ \begin{matrix}0, & x ∈ X, \\ +∞, & x ∉ X. \end{matrix}

The subdifferential mapping of δ_{X} is the normal cone operator associated with the set X,

N_{X}(x) := \{ \begin{matrix}\{z ∈ ℜ^{n} | z^{T}(y − x) ≤ 0, ∀y ∈ X\}, & x ∈ X, \\ ∅, & x ∉ X. \end{matrix} (1.5)

(This is an example of a maximal monotone operator.) The problem [GE] then reduces to the problem of finding a vector x^{∗} such that

∇f(x^{∗}) + N_{X}(x^{∗}) ∋ 0,

which (without the need for Assumption 1.1) describes the first-order necessary conditions for the optimality of x^{∗} ∈ X in the constrained differentiable program

[CDP]

min_{x∈X} f(x).

We also note that the problem [NDP] encompasses dual formulations of convex programs (e.g., [FHN96]).

A large literature has been devoted to the study of computational methods for [NDP]. The present paper provides a framework of algorithms for [NDP] which includes many methods previously analyzed, and provides a convergence analysis which unifies and improves upon existing ones. The analysis includes investigations of inexact solutions of the auxiliary problems, different step length rules, and convergence rate results.

2 The cost approximation algorithm

The term cost approximation (CA for short) was coined in the author’s PhD thesis [Pat93a] to describe a framework of descent algorithms for nonlinear programs and variational inequality problems. In this section we present the algorithm framework for the solution of [NDP]. The main idea is to iteratively approximate the cost mapping ∇f, in order to obtain a problem which is more easily solved, in the sense that a (possibly) non-monotone mapping is replaced by a monotone one. (The term cost mapping is taken from applications of [GE] in equilibrium problems, where ∇f usually represents a cost vector.) The solution to this problem defines a search direction, in which a step is taken to decrease the value of T.

2.1 The subproblem phase

Let x ∈ dom u. We introduce a continuous and monotone cost approximating mapping Φ : dom u ↦ ℜ^{n}. If the mapping ∇f is replaced by Φ, then the error made in the approximation obviously is ∇f − Φ. This error is taken into account by adding to Φ the fixed error term ∇f(x) − Φ(x). The approximation may alternatively be interpreted as a fixation of the second term of the original cost [Φ + ∂u] + [∇f − Φ] at x.

Thus, we arrive at a generalized equation, in which a point y is sought such that

[GE_{Φ}]

Φ(y) + ∂u(y) + ∇f(x) − Φ(x) ∋ 0.  (2.1)

We let Y (x) denote the set of solutions y to [GE_{Φ}]. (Note that a fixed term may be added to Φ without affecting this subproblem.)

By construction, the approximation is exact at x. This fact is important, since it provides a termination criterion for the algorithm: if x solves the subproblem [GE_{Φ}], then it immediately follows that x also solves [GE]. The reverse is also true, and the subproblem thus provides a reformulation of [GE] as a fixed point problem in the (possibly point–to– set) mapping x ↦ Y (x).

Remark 2.1 The descent property of the search direction requires the mapping Φ to be monotone; this crucial property is, seemingly, not possible to relax within this framework.

In several known instances of CA algorithms, the corresponding mapping Φ may be identified as an approximation of ∇f (cf. Newton-type methods); the properties of ∇f in a given application may hence, implicitly through the requirements on Φ, restrict the possible choices of CA methods.

The convergence analysis provided in this paper concentrates on approximations made of the differentiable (and possibly non-convex) function f, and does not cover approximations of the convex function u, although such algorithms are easily devised within the framework; examples are subgradient optimization methods (e.g., [Sho85]) and the methods in [CoZ84]. Note that in such methods, however, the corresponding subproblems need not yield directions of descent with respect to T without proper modifications. □

If Φ is chosen as the gradient mapping of a function φ : dom u ↦ ℜ, which then is convex and continuously differentiable on dom u, then the subproblem [GE_{Φ}] reduces to solving the convex subproblem

[NDP_{φ}]

min_{y∈ℜ^{n}} T_{φ}(y, x) := φ(y) + u(y) + f(x) − φ(x) + [∇f(x) − ∇φ(x)]^{T}(y − x).

To give an example, assume that f is a convex function in C^{2} on dom u, and choose φ(y, x) := \frac{1}{2}(y − x)^{T}∇^{2}f(x)(y − x), whence one obtains T_{φ}(y, x) := u(y) + \{f(x) + ∇f(x)^{T}(y − x) + \frac{1}{2}(y − x)^{T}∇^{2}f(x)(y − x)\}; an extension of Newton’s method to nondifferentiable optimization is obtained, in which only the differentiable part of the objective is approximated. The choice φ(y, x) := \frac{1}{2γ} ∥y − x∥^{2}, γ > 0, leads to extensions of steepest descent, gradient projection and proximal point algorithms.

The construction of this subproblem may be given an alternative interpretation as a partial linearization ([Pat93a]) of T; expressing it as [φ + u] + [f − φ], the function T_{φ} is obtained from a first-order approximation (or, linearization) of the second term at x.

Thus, by applying the CA concept, the problem [NDP] can be solved through a sequence of convex optimization problems. Indeed, in most of the iterative methods for the solution of [NDP] that one may identify as special cases from the class of CA algorithms, subproblems of the form [NDP_{φ}] are solved. (Examples where the subproblem [GE_{Φ}] does not reduce to the convex optimization problem [NDP_{φ}] include some matrix splitting methods for quadratic programs [CPS92] and Gauss–Seidel/SOR methods [OrR70].)

2.2 The line search phase

In general, we can not expect the original problem to be solved by the subproblem solution y. An improved solution is therefore defined through a step taken in the direction of d := y − x such that the value of a merit function for [NDP] is reduced sufficiently. Although any merit function with sufficient continuity properties will do (such as the Euclidean distance to Ω), we choose T as the merit function. (In the solution of asymmetric variational inequality problems, gap functions associated with the CA subproblem [NDP_{φ}] are available; see [Pat93c, LaP94, Pat94a, Pat94b, Pat97] for further details.)

At the new point, the original mapping is again approximated—perhaps using another mapping Φ—and the algorithm proceeds until some stopping criterion is fulfilled.

2.3 The conceptual algorithm

A description of the CA algorithm is given in Table 2.1. A sequence {Φ^{t}} of monotone mappings is assumed to be given. (Note, however, that each mapping may also be constructed adaptively, given x^{s}, s = 1, 2, …, t, instead of being chosen a priori.)

Table 2.1: The cost approximation algorithm

0. (Initialization): Choose an initial point x^{0} ∈ dom u, and let t := 0.

1. (Search direction generation): Find a solution y^{t} to

[GE_{Φ^{t}}]

Φ^{t}(y^{t}) + ∂u(y^{t}) + ∇f(x^{t}) − Φ^{t}(x^{t}) ∋ 0. (2.2)

The resulting search direction is d^{t} := y^{t} − x^{t}.

2. (Termination criterion): If x^{t} solves [GE_{Φ^{t}}] → Stop (x^{t} ∈ Ω). Otherwise, continue.

3. (Line search): Choose a step length, ℓ_{t}, such that x^{t} + ℓ_{t}d^{t} ∈ dom u and the value of T is reduced sufficiently.

4. (Update): Let x^{t+1} := x^{t} + ℓ_{t}d^{t}, and t := t + 1.

5. (Termination criterion): If x^{t} is acceptable → Stop. Otherwise, go to Step 1.

It will be made clear in the sequel how accurate the steps 1 and 3 need to be performed.

2.4 Conventions and assumptions

In the convergence analysis of the CA algorithm, three forms of the cost approximating mapping will be used. In the most general case, a sequence {Φ^{t}} of continuous and monotone mappings on dom u is used; the mappings in this sequence are neither assumed to be dependent on each other nor to form a limit function (that is, equicontinuity is not assumed). In the second case, it is assumed that there exists a known continuous mapping Φ : dom u × dom u ↦ ℜ^{n} which is monotone on dom u in its first argument; each individual cost approximating mapping Φ^{t} in the sequence {Φ^{t}} then has the form Φ(·, x^{t}). (This form will be required whenever we do not wish to impose strong monotonicity conditions on Φ.) The third, and most simple, form of mapping to be used is the iteration independent form, that is, Φ^{t} ≡ Φ, for all t, for a given continuous and monotone mapping Φ : dom u ↦ ℜ^{n}. (The latter form will only be used in the linear convergence Theorem 4.6 and its Corollary 5.1.) Whenever the iteration dependency of Φ^{t} is insignificant in the context, the superscript will be suppressed. (This applies especially to the technical lemmas in Section 3.)

In general, the CA algorithm is only asymptotically convergent. (Finite convergence is obtained when the problem enjoys a sharpness property; see [Pat93c] for further details.) Therefore, we shall implicitly presume that the sequence of iterates is infinite. Further, we do not study the behaviour of the algorithm when a solution does not exist.

The remainder of the paper is organized as follows. In the next section, we provide some properties of the search directions and the step length rules considered, for different assumptions on the mappings Φ^{t}. In Section 4, a convergence analysis is made. Some consequences of this analysis are described in Section 5, demonstrating its strength compared to previous analyses of algorithms included in the framework.

3 Preliminaries

3.1 Analysis of the search directions

Lemma 3.1 (Characterizations of the search directions) Let x ∈ dom u, Y (x) be the (possibly empty) set of solutions to [GE_{Φ}], y ∈ Y (x), and d := y − x.

(a) (A fixed point characterization of Ω)

(1) x ∈ Ω ⟺ x ∈ Y (x).

(2) Let Φ ≡ ∇φ. Then, x ∈ Ω ⟺ x ∈ Y (x) ⟺ T_{φ}(y, x) = T_{φ}(x, x).

(b) (A characterization of Y (x)) Let Φ be maximal monotone. Then, Φ+∂u is maximal monotone and

Y (x) = [Φ + ∂u]^{−1}[Φ − ∇f](x). (3.1)

(c) (A well-posedness characterization)

(1) Let Φ be maximal monotone. Then, the set Y (x) is nonempty and locally bounded if

Φ(x) − ∇f(x) ∈ int (dom ([Φ + ∂u]^{−1})). (3.2)

(2) Let Φ ≡ ∇φ. Then, the set Y (x) is nonempty and bounded if and only if

∇φ(x) − ∇f(x) ∈ int (dom ([φ + u]^{∘})), (3.3)

where [φ + u]^{∘} is the conjugate function of φ + u.

(d) Let Φ : dom u × dom u ↦ ℜ^{n} be a continuous mapping on dom u × dom u of the form Φ(y, x), monotone on dom u with respect to y. Let S ⊆ dom u be an arbitrary, closed set in ℜ^{n}. Then, the mapping x ↦ D(x) := Y (x) − x is closed on S.

(e) Let u be strictly convex on dom u or Φ strictly monotone on dom u. If Y (x) is nonempty, then it is a singleton.

Proof.

(a) (1) Follows immediately from the construction of [GE_{Φ}].

(2) The second equivalence follows immediately from the construction of [NDP_{φ}].

(b) The first result follows from [Roc70b, Roc70c]; the second result then follows by direct calculation.

(c) (1) The set Y (x) is characterized by (3.1). Let Π := Φ + ∂u. Since, by definition, range Π = dom Π^{−1}, the set Y (x) is nonempty if and only if [Φ − ∇f](x) ∈ dom Π^{−1}; further, it is locally bounded if [Φ−∇f](x) ∈ int (dom Π^{−1}) ([Roc69]).

(2) From [Roc70a, Cor. 23.5.1 and Thm. 27.1.b], the set Y (x) is characterized by Y (x) = ∂([φ + u]^{∘})(∇φ(x) − ∇f(x)). By [Roc70a, Thms. 23.4 and 27.1.d], Y (x) is nonempty and bounded if and only if (3.3) holds.

(d) Let {x^{t}} ⊆ S and {y^{t}} be sequences in dom u with y^{t} ∈ Y (x^{t}), {x^{t}} → x and {y^{t}} → y. By the continuity of Φ and ∇f and the convexity of u, Theorem 24.4 of [Roc70a] yields that y ∈ Y (x). Hence, the mapping x ↦ Y (x) is closed on S. The graph of the mapping D is obtained from that of Y through an affine transformation; it then follows from [Ber63, p. 111] that x ↦ D(x) is closed on S, as desired.

(e) Follows from the strict monotonicity of Φ + ∂u. □

Remark 3.1 The fixed point property in (a) validates the termination criterion of Step 2 of the CA algorithm; the result (a)(2) shows that it is also easily checked whenever Φ is a gradient mapping.

The characterization (3.1) of the set Y (x) made in (b) bears resemblance to iterative formulas describing proximal point methods ([Mar70, Roc76]), and forward-backward splitting methods ([Bre73, EcB92]) for finding a zero of the sum of two maximal monotone mappings. The class of CA algorithms is actually far more general; for example, the mapping Φ is not necessarily affine or strongly monotone and can be chosen to adapt to problem structures (such as separability in u and f), and a (possibly inexact) line search is embedded in the algorithm. We also note that a monotone and single-valued mapping on ℜ^{n} automatically is maximal monotone.

The properties of Y (x) characterized in the result (c) are natural well-posedness requirements on the subproblem [GE_{Φ}]. A sufficient condition for (3.2) to hold is that Φ+∂u is weakly coercive ([Zei90, Def. 32.34.b]), since this condition implies that dom ([Φ+ ∂u]^{−1}) = ℜ^{n} (see [Zei90, Cor. 32.35] and [Bro68]). The corresponding sufficient condition in the case where Φ is a gradient is that φ + u is coercive [cf. (1.3)], since φ + u then is cofinite ([Roc70a, pp. 116 and 259]) and dom ([φ+u]^{∘}) = ℜ^{n} holds ([Roc70a, Cor. 13.3.1]). All these conditions are satisfied if dom u is bounded; an important special case is when u includes, as an additive term, an indicator function of a compact set.

A function φ : dom u × dom u ↦ ℜ of the form φ(y, x), continuous on dom u × dom u and convex and differentiable on dom u with respect to y, automatically has the continuity property of the mapping Φ = ∇_{y}φ required in the result (d). Note that dom u is not closed in general ([Roc70a, p. 52]), in which case it can not take the role of the set S. □

Lemma 3.2 (Descent properties) Let x ∈ dom u \ Ω, Y (x) be the (possibly empty) set of solutions to [GE_{Φ}], y ∈ Y (x), and d := y − x.

(a) If either u is strictly convex on dom u and ∂u(x) is bounded, or Φ is strictly monotone on dom u, then T′(x; d) < 0.

(b) Let Φ ≡ ∇φ, \overline{y} ∈ ℜ^{n} be any point such that T_{φ}(\overline{y}, x) < T_{φ}(x, x), and \overline{d} := \overline{y} − x. Then, T′(x; \overline{d}) < 0. Especially, T′(x; d) < 0.

(c) Let u be strongly convex on dom u or Φ strongly monotone on dom u. Let \overline{y} ∈ ℜ^{n} be an approximate solution to [GE_{Φ}] in the sense that for some vector r ∈ ℜ^{n},

Φ(\overline{y}) + ∂u(\overline{y}) + ∇f(x) − Φ(x) ∋ r, (3.4) and \overline{d} := \overline{y} − x. If ∥r∥ < (m_{u} + m_{Φ})∥\overline{d}∥, then T′(x; \overline{d}) < 0. Especially, Y (x) is nonempty and singleton, and

T′(x; d) ≤ −(m_{u} + m_{Φ})∥d∥^{2}. (3.5)

(d) Let Φ be strongly monotone on dom u and Lipschitz continuous on dom u. Then, there exists a ξ_{u}(y) ∈ ∂u(y) such that

−\frac{T′(x; d)}{∥∇f(x) + ξ_{u}(y)∥ · ∥d∥} ≥ ρ, (3.6)

where ρ := m_{Φ}/M_{Φ}.

Remark 3.2 The reader should note here that in (3.5) we make use of the convention that if a mapping Π is monotone, then it satisfies the inequality (1.4) with m_{Π} = 0. □

Proof.

(a) If u is strictly convex, then for all ξ_{u}(x) ∈ ∂u(x), [−Φ(y)−∇f(x)+Φ(x)−ξ_{u}(x)]^{T}d > 0 holds, which, from the monotonicity of Φ, yields

[∇f + ξ_{u}(x)]^{T}d < 0,  ∀ξ_{u}(x) ∈ ∂u(x).  (3.7)

If ∂u(x) is bounded, then the supremum of the left-hand side of the above inequality over ξ_{u}(x) ∈ ∂u(x) is attained, and the result follows.

If Φ is strictly monotone, then the monotonicity of ∂u ([Roc70a, Cor. 31.5.2]) yields [∇f + ξ_{u}(x)]^{T}d ≤ [Φ(x) − Φ(y)]^{T}d < 0, which implies the desired result.

(b) By the convexity of φ, T_{φ}(\overline{y}, x) < T_{φ}(x, x) implies that ∇f(x)^{T}\overline{d} < u(x) − u(\overline{y}). Then, from (1.2),

T′(x; \overline{d}) = ∇f(x)^{T}\overline{d} + sup\{ ξ_{u}(x)^{T}(\overline{y} − x) | ξ_{u}(x) ∈ ∂u(x) \}
< u(x) − u(\overline{y}) + sup\{ ξ_{u}(x)^{T}\overline{d} | ξ_{u}(x) ∈ ∂u(x) \},  (3.8)

and the right-hand side of (3.8) is nonpositive, since u is convex. The point y is just a special case of such a \overline{y}.

(c) We obtain from (3.4) that

[∇f(x) + ξ_{u}(x)]^{T}\overline{d} ≤ −(m_{u} + m_{Φ})∥\overline{d}∥^{2} + ∥r∥ · ∥\overline{d}∥

holds for all ξ_{u}(x) ∈ ∂u(x), from which the result easily follows.

(d) By (2.1) and the Lipschitz continuity of Φ, there exists a ξ_{u}(y) ∈ ∂u(y) such that ∥∇f(x) + ξ_{u}(y)∥ ≤ M_{Φ}∥d∥. Together with (3.5), with m_{u} = 0, (3.6) follows. □

Remark 3.3 The second result of (a) requires that Φ is strictly monotone if it is not a gradient mapping [cf. the result (b)]; that non-strict monotonicity is not sufficient is clear from the following example: take u := δ_{X}, where X := { x ∈ ℜ^{2} | 0 ≤ x ≤ 2 }, and f(x) := \frac{1}{2}(x_{1}^{2} + x_{2}^{2}). Let further x := (2, 2)^{T}, which clearly is non-optimal, and

Φ(x) := (\begin{matrix}0 & −1 \\ 1 & 0 \end{matrix}) x,

which is a non-strictly monotone and asymmetric mapping. By simple calculations, y = (0, 4)^{T} can be shown to solve [GE_{Φ}], but ∇f(x)^{T}(y − x) = 0.

The result (b) establishes that descent properties hold under very mild assumptions if Φ is a gradient mapping, and validates inexact solutions of [NDP_{φ}]. (In the above example, choose, for example, Φ := 0; then y = (0, 0)^{T} solves [GE_{Φ}], and ∇f(x)^{T}(y − x) = −8.)

The result (c) establishes the descent property under conditions on the accuracy of the solution of [GE_{Φ}]. While solving [GE_{Φ}] with some iterative procedure, the validity of (3.4) for a given value of r can often be checked easily.

The vector r appearing in (3.4) need not be associated only with a measure of the distance from \overline{y} to Y (x). The vector can include (or be equal to) an error in the calculation of ∇f(x); the result (c) thus provides conditions for obtaining descent when the value of ∇f(x) is observed in the presence of noise. See [Pol87, Chap. 4] for discussions on such errors and their effect on gradient methods.

The inequality (3.6) in (d) generalizes the property of gradient relatedness of search directions used in unconstrained differentiable optimization ([OrR70, Def. 14.3.1]). □

3.2 Step length rules

This subsection presents and validates the step length rules considered in the realization of Step 3 of the CA algorithm. Four different step length rules are used; the first is a conceptual one.

Definition 3.1 (Exact line search) Let x ∈ dom u, y ∈ Y (x), and d := y − x.

(E) Choose ℓ as a solution to min { T (x + ℓd) | ℓ ≥ 0 }.

We next define a new step length rule, which generalizes the Armijo [Arm66] rule from differentiable optimization. The original statement of the Armijo rule can not be applied here, due to the nondifferentiability of T; the fact that the generalized gradient can not be Lipschitz continuous makes it necessary to replace the directional derivative used in the original Armijo rule with a monotonicity measure in Φ. Also, the unit step is the largest for which the modified Armijo rule can be validated (see Lemma 3.3 below); we therefore choose a unit first trial step length. Further, Φ must be chosen strictly monotone.

Definition 3.2 (Armijo step length) Let Φ be strictly monotone on dom u. Let x ∈ dom u \ Ω, Y (x) := {y}, and d := y − x. Further, let α, β ∈ (0, 1).

(A’) Let ℓ := βī, where ī is the smallest nonnegative integer i such that

T (x + β^{i}d) − T (x) ≤ αβ^{i}[Φ(x) − Φ(y)]^{T}d. (3.9)

To see that this step length rule generalizes the Armijo rule, let u ≡ 0 in [NDP]. Then a subproblem solution y is characterized by the system of nonlinear equations ∇f(x) = Φ(x) − Φ(y) [cf. (2.2)]; using this relationship in (3.9) yields that Rule A’ reduces to the original Armijo rule, with a unit first trial step. (Note that Rule A’ uses an over-estimating model of T, since T′(x; d) ≤ [Φ(x) − Φ(y)]^{T}d holds.)

The Armijo Rule A’ is validated below.

Lemma 3.3 (Validity of the Armijo rule) Let Φ be strictly monotone on dom u. Let x ∈ dom u \ Ω, Y (x) := {y}, and d := y − x. Further, let α, β ∈ (0, 1).

(a) There exists a nonnegative integer ī satisfying (3.9).

(b) Assume further that ∇f is Lipschitz continuous on dom u. Then, either a unit step satisfies (3.9) or

βī > 2β(α − 1)\frac{[Φ(x) − Φ(y)]^{T}d}{M_{∇f}∥d∥^{2}}. (3.10)

Proof.

(a) Taylor’s formula yields, for ℓ ∈ [0, 1],

f(x + ℓd) − f(x) = ℓ∇f(x)^{T}d + o(ℓ). (3.11)

Further, for ℓ ∈ [0, 1] (cf. [Roc70a, Cor. 24.2.1 and Thm. 24.1], [Roc66]),

u(x + ℓd) − u(x) = ∫_{0}^{ℓ} u′(x + sd; d) ds ≤ −ℓu′(x + d; −d)
= ℓ inf{ ξ_{u}(y)^{T}d | ξ_{u}(y) ∈ ∂u(y) }
≤ ℓξ_{u}(y)^{T}d,  ∀ξ_{u}(y) ∈ ∂u(y).  (3.12)

In (3.12), let ξ_{u}(y) ∈ ∂u(y) be given by (2.1). Then, combining (3.11) and (3.12),

T (x + ℓd) − T (x) ≤ ℓ[Φ(x) − Φ(y)]^{T}d + o(ℓ).

For any value of α ∈ (0, 1), there must therefore be an \overline{ℓ} > 0 such that, for all ℓ ∈ (0, \overline{ℓ}), T (x + ℓd) − T (x) ≤ αℓ[Φ(x) − Φ(y)]^{T}d holds, whence the result.

(b) We begin by showing that for ℓ ∈ [0, 1],

T (x + ℓd) − T (x) ≤ ℓ[Φ(x) − Φ(y)]^{T}d + \frac{M_{∇f}}{2} ℓ^{2}∥d∥^{2} (3.13)

holds. Taylor’s formula and the Lipschitz continuity of ∇f yields, for ℓ ∈ [0, 1],

f(x + ℓd) − f(x) = ∫_{0}^{ℓ} ∇f(x + sd)^{T}d ds ≤ ℓ∇f(x)^{T}d + ∫_{0}^{ℓ} M_{∇f}∥d∥^{2}s ds
= ℓ∇f(x)^{T}d + \frac{M_{∇f}}{2} ℓ^{2}∥d∥^{2}.  (3.14)

Combining (3.14) with (3.12), the latter obtained as in (a) from (2.1), we obtain (3.13). Using (3.13) and the strict monotonicity of Φ, by choosing ℓ ≤ 2(α − 1)([Φ(x) − Φ(y)]^{T}d)/(M_{∇f}∥d∥^{2}), we obtain

T (x + ℓd) − T (x) ≤ αℓ[Φ(x) − Φ(y)]^{T}d. (3.15)

Replacing ℓ with βī, we see that ī is the smallest integer i to satisfy

β^{i} ≤ 2(α − 1)\frac{[Φ(x) − Φ(y)]^{T}d}{M_{∇f}∥d∥^{2}}. (3.16)

From (3.16) then follows the desired result. □

Remark 3.4 From (3.16) we may find conditions under which the unit step, that is, ī = 0, is accepted by Rule A’. Assume that Φ is strongly monotone on dom u. Then, a unit step is obtained if

\frac{2(1 − α)m_{Φ}}{M_{∇f}} ≥ 1. (3.17)

Note that, from (3.10), the step length given by Rule A’ always satisfies ℓ ≥ min{1, 2β(1− α)([Φ(y) − Φ(x)]^{T}d)/(M_{∇f}∥d∥^{2})}, which reduces to

ℓ ≥ min \{1, \frac{2β(1 − α)m_{Φ}}{M_{∇f}} \} (3.18)

whenever Φ is strongly monotone. □

If the calculation of T is expensive, then a step length formula which does not require its calculation is of interest. Such a rule is described next.

Definition 3.3 (Relaxation step length) Assume that ∇f is Lipschitz continuous on dom u, and let Φ be strongly monotone on dom u. Further, let 0 < ε_{i} < m_{Φ}/M_{∇f}, i = 1, 2.

(R) Choose ℓ in the closed interval [ε_{1}, min{1, 2m_{Φ}/M_{∇f} − ε_{2}}].

Remark 3.5 The inverse of the quotient m_{Φ}/M_{∇f} appearing in the description of the maximal step length is in some cases an approximate condition number of the Hessian of f; consider, for example, Newton’s method. □

The step length rule is validated below.

Lemma 3.4 (Validity of the relaxation step) Assume that ∇f is Lipschitz continuous on dom u, and let Φ be strongly monotone on dom u. Let x ∈ dom u \ Ω, Y (x) := {y}, and d := y − x. If ℓ is chosen according to Rule R, then T (x + ℓd) < T (x).

Proof. The inequality (3.13), together with the strong monotonicity of Φ, yields, for all ℓ ∈ [0, 1], that

T (x + ℓd) − T (x) ≤ ℓ (−m_{Φ} + \frac{M_{∇f}}{2} ℓ) ∥d∥^{2}, (3.19)

from which the result easily follows. □

Note that if

\frac{2m_{Φ}}{M_{∇f}} > 1 (3.20)

holds, then a unit step yields descent. [This condition is implied by (3.17).]

The relaxation step length has the advantage of being very simple, but it requires the knowledge of (an estimate of) the Lipschitz constant of ∇f. (In certain cases, such as when f is quadratic or separable, estimates of the Lipschitz constant M_{∇f} are available; see, for example, [WoZ96]. The maximal step length is usually determined through experiments.) The last step length rule introduced yields descent eventually, and does not require the knowledge (or estimation) of any problem or method parameters.

Definition 3.4 (Divergent series step length)

(D) Let the sequence {ℓ_{t}} satisfy

[0, 1] ⊃ {ℓ_{t}} → 0,  ∑_{t=0}^{∞} ℓ_{t} = ∞.  (3.21)

Remark 3.6 The requirements of the relaxation and divergent series step length rules need not be viewed as defining step length formulas per se, but can instead be viewed as conditions that some arbitrary step length rule must satisfy in order to yield convergence. The almost complete relaxation strategy described in [DeV85, Sec. 3.4], for example, allows for the utilization of an arbitrary step length selection rule through the (almost never used) upwards or downwards rounding of a tentative step length in order to fall between two step lengths (very small and very large, respectively) defined by preselected formulas both satisfying (3.21). □

4 Convergence results

4.1 The conceptual method

In our first convergence result, we will assume that the mappings Φ^{t} are monotone only. In order to ensure the convergence of the CA algorithm, this will require that the sequence {Φ^{t}} is constructed by a given function φ on dom u × dom u (that is, that Φ^{t} ≡ ∇φ(·, x^{t}) for each t), and that the line search is performed exactly.

Theorem 4.1 (Convergence under Rule E) Assume that u is continuous on dom u. Let the sequence {Φ^{t}} of cost approximating mappings be constructed such that, for each t, Φ^{t} ≡ ∇φ(·, x^{t}) for a given function φ : dom u × dom u ↦ ℜ of the form φ(y, x), continuous on dom u×dom u and convex and in C^{1} on dom u with respect to y. Assume that x^{0} ∈ dom u is such that the lower level set L(x^{0}) := { x ∈ dom u | T (x) ≤ T (x^{0}) } is bounded, and further that the problem [NDP_{φ}] is well defined, in the sense that (3.3) holds for every x ∈ L(x^{0}). Let Rule E be used. Then, \{T (x^{t})\} → T (\overline{x}) for some \overline{x} ∈ Ω, any accumulation point of the sequence {x^{t}} (at least one such point exists) lies in Ω, and

{ inf_{x∈Ω∩L(x^{0})} ∥x^{t} − x∥} → 0. (4.1)

Proof. In order to apply Zangwill’s Theorem A [Zan69, Sec. 4.5], we first identify the solution set with Ω, the descent function with T, and the algorithmic map with the composite mapping A := ED, where D is the direction finding mapping and E is the exact line search map. To fulfill the assumptions of Theorem A, we next show that (1) the sequence {x^{t}} lies in a compact set, (2) T (x^{t+1}) < T (x^{t}) if x^{t} ∉ Ω and the algorithm terminates if x^{t} ∈ Ω, and (3) the mapping A is closed at all points in L(x^{0}) \ Ω.

(1) From the boundedness assumption and [StW70, p. 134], the set L(x^{0}) is compact; the descent property then ensures that the sequence {x^{t}} lies in a compact set.

(2) The descent property follows from Lemma 3.2.b. The termination property (x^{t} ∈ Ω) follows from Lemma 3.1.a.2.

(3) The mapping D is closed on L(x^{0}) by Lemma 3.1.d [identifying S with L(x^{0})], and the same is true for the mapping E, since T is continuous on dom u ([Zan69, Le. 5.1]). The mapping ∂([φ(·, x) + u]^{∘}) is monotone, and therefore, by (3.3) and [Zei90, Prop. 32.33], locally bounded at ∇φ(x, x)−∇f(x). From (3.1), the mapping Y is therefore locally bounded on L(x^{0}). From (1), the sequence {x^{t}} is bounded. The local boundedness of Y implies that also the sequence {y^{t}} is bounded. This sequence then has an accumulation point. Then, [Zan69, Le. 4.2] ensures that the composite map A is closed on L(x^{0}).

The assumptions of Theorem A are thus fulfilled, and we conclude that {T (x^{t})} converges to T (\overline{x}) for some \overline{x} ∈ Ω, and any accumulation point of \{x^{t}\} lies in Ω. (The existence of at least one such point follows from the boundedness of {x^{t}}.) The last statement follows from [OrR70, Thm. 14.1.4]. □

Note that it is essential that u is continuous on L(x^{0}), since otherwise the line search mapping E need not be closed. (In the special case where u ≡ δ_{X}, the line search is made with respect to f, and the continuity assumption on u may be removed.)

4.2 A truncated algorithm

In this section we provide one example of the many possibilities for constructing realizations of truncated CA algorithms. The basis for the algorithm described and validated here is Lemma 3.2.b. (Another truncated CA algorithm, based on Lemma 3.2.c, is presented in Section 4.3.) The idea behind a truncated CA method is to reduce the work performed on [NDP_{φ}] by limiting the number of iterations performed when solving it using a descent algorithm. This strategy introduces a trade-off between the computational effort spent on solving the subproblem and the quality of the search direction obtained. (Examples of such methods are extensions of truncated Newton methods for systems of nonlinear equations, in which a limited number of steps of an iterative scheme is performed on the linear systems.)

The following assumption is made on the algorithm used for solving [NDP_{φ}].

Assumption 4.1 (Properties of a truncated algorithm) Given an x ∈ dom u, the problem [NDP_{φ}] is solved, starting from x, with an iterative algorithm where one iteration may be described by a mapping G : dom u × dom u ↦ 2^{dom u} × dom u having the following properties.

(1) (Fixed point) z ∈ Y (x) ⟺ (z, x) ∈ G(z, x), in which case the algorithm terminates.

(2) (Descent) If x ∉ Ω and z ∉ Y (x), then T_{φ}(y, x) < T_{φ}(z, x) for all y with (y, x) ∈ G(z, x).

(3) (Closedness) The mapping G is u.s.c., closed and compact valued on dom u × L(x^{0}).

(4) (Finite termination) The number of iterations performed is bounded from below by 1, and from above by a positive integer \overline{k}.

Remark 4.1 The condition (1) ensures that the fixed point property of the original algorithm is preserved, and that the algorithm terminates whenever a solution to [NDP_{φ}] has been obtained. The condition (2), together with the lower bound in (4) on the number of iterations performed and Lemma 3.2.b, ensures that the resulting approximate solution \overline{y} defines a direction of descent with respect to T. Condition (3) ensures the closedness of the mapping describing the overall algorithm. Whenever dom u is bounded, it can be replaced by the assumption that the mapping G is closed on dom u×dom u. Condition (4), finally, limits the work done in each main iteration of the CA algorithm, and, together with (iii), it ensures that any sequence \{\overline{y}\} of approximate solutions generated over a compact set is bounded.

If dom u is bounded, then the mapping G describing any truncated CA algorithm based on a function φ̂ : dom u×dom u ↦ ℜ of the form φ̂(y, x), together with any closed updating step which guarantees descent (such as Rule E or a fixed relaxation step), satisfies (1)–(3). □

The truncated CA algorithm is obtained from replacing Step 1 of the CA algorithm (see Table 2.1) by:

1’. (Search direction generation) Apply 1 ≤ k_{t} ≤ \overline{k} iterations of the algorithm described by G, starting from x^{t}, that is, let

(\overline{y}^{t}, x^{t}) ∈ \underbrace{G ∘ … ∘ G}_{k_{t} times G}(x^{t}, x^{t}).

The resulting search direction is d^{t} := \overline{y}^{t} − x^{t}.

Note that the value of k_{t} need not be determined a priori, and could instead be regarded as a product of the algorithm and the termination criteria chosen for [NDP_{φ}].

Theorem 4.2 (Convergence of a truncated algorithm) Replace Step 1 by Step 1’ in the CA algorithm. In Theorem 4.1, replace the assumption that (3.3) holds for every x ∈ L(x^{0}) with Assumption 4.1. Then, the conclusions of Theorem 4.1 hold for the truncated CA algorithm.

Proof. The proof utilizes Zangwill’s Theorem, and the Spacer Step Theorem [Lue84, p. 231]. We begin by identifying the algorithmic mapping A. We assume that {x^{t}} is infinite (otherwise, by Assumption 4.1.1, with z = x, the algorithm is terminated at a point in Ω). By Assumption 4.1.4, there must be at least one positive integer, say k, that occurs an infinite number of times in the sequence {k_{t}}. We define

\overline{Y} (x) := \{ \overline{y} ∈ ℜ^{n} | (\overline{y}, x) ∈ \underbrace{G ∘ … ∘ G}_{k times G}(x, x) \}.

The algorithmic mapping is A := E ∘ \overline{D}, where \overline{D} := \overline{Y} − · and E is the exact line search map.

We next establish that this mapping has the properties desired (boundedness, adaption, closedness). We first note that Assumption 4.1.2, together with Lemma 3.2.b, yields that the mapping \overline{D} has the descent property. From the use of Rule E, it then follows that T (\overline{x}) < T (x) for any \overline{x} ∈ E ∘ \overline{D}(x), and therefore A satisfies the second condition (adaption). The boundedness assumption on the set L(x^{0}) then implies that the sequence {x^{t}} is bounded, thereby satisfying the first condition of Zangwill (boundedness).

We now establish the third condition (closedness). The mapping G is closed on dom u×L(x^{0}) (Assumption 4.1.3). By the upper semicontinuity and compactness of G (Assumption 4.1.3 again), the image ∪_{x∈L(x^{0})}G(x, x) of G on the compact set L(x^{0}) is compact ([Nik68, Lem. 4.5]). Lemma 12 of [Mey79] then implies that if \overline{y}^{∞} is an accumulation point of a sequence \{\overline{y}^{t}\} given by (\overline{y}^{t+1}, x) ∈ G(\overline{y}^{t}, x), then the set \{ y | (y, x) ∈ G ∘ G(\overline{y}^{∞}, x) \} contains an accumulation point of \{\overline{y}^{t}\}. We are then in the position of using Lemma 4.2 of [Zan69] to conclude that the composite mapping G ∘ G is closed on dom u × L(x^{0}). This result used repeatedly then establishes that the mapping \overline{Y}, and therefore also \overline{D}, is closed on L(x^{0}). Applying Lemma 4.2 of [Zan69] again, we may conclude that the mapping E ∘ \overline{D} is closed on L(x^{0}). Moreover, the boundedness of L(x^{0}) and the descent property of E ∘ \overline{D} ensures that the range of this mapping over L(x^{0}) is a compact set. We invoke Corollary 4.2.1 of [Zan69] to conclude that the algorithmic mapping A is closed on L(x^{0}).

The algorithmic mapping identified by extracting the value k from the choices in the sequence {k_{t}} is clearly of the form

C(x) = { y ∈ dom u | T (y) ≤ T (x) }, x ∈ dom u.

We may therefore invoke the Spacer Step Theorem [Lue84, p. 231], which guarantees that the result holds, thanks to the properties of the mapping A. The rest of the proof follows that of Theorem 4.1. □

A natural sufficient condition for the boundedness of \{\overline{y}^{t}\} is that dom u is bounded; the special case u ≡ δ_{X} for a nonempty, compact and convex set X is covered in [Pat93b].

4.3 Implementable step length rules

Theorem 4.3 (Convergence under Rule A’) Let Φ : dom u × dom u ↦ ℜ^{n} be a continuous mapping on dom u × dom u of the form Φ(y, x), maximal and strictly monotone on dom u with respect to y. Assume that the point x^{0} is chosen so that the lower level set L(x^{0}) is bounded, and assume further that the problem [GE_{Φ}] is well defined, in the sense that (3.2) holds for every x ∈ L(x^{0}). Let Rule A’ be used. Then, any accumulation point of the sequence {x^{t}} (at least one such point exists) lies in Ω, and (4.1) holds.

Proof. Let x^{∞} be any accumulation point of the sequence {x^{t}}; the existence of such a point follows from the boundedness of L(x^{0}) and the descent property. Likewise, the sequence {y^{t}} must be bounded, by the well-posedness assumption (3.2); let y^{∞} be an arbitrary accumulation point. Hence, {d^{t}} is bounded.

We will next show that {d^{t}} → 0 must hold; the proof is by contradiction.

By Lemma 3.3.a,

T (x^{t+1}) − T (x^{t}) ≤ αℓ_{t}[Φ(x^{t}, x^{t}) − Φ(y^{t}, x^{t})]^{T}d^{t} < 0, t = 0, 1, ….

If {d^{t}} ↛ 0, then as {T (x^{t+1}) − T (x^{t})} → 0 holds, there must be a subsequence T such that \{ℓ_{t}\}_{T} → 0. There must then be an index \overline{t} such that for every t ≥ \overline{t} in T , the step length produced by Rule A’ is less than one, that is,

T (x^{t} + (ℓ_{t}/β)d^{t}) − T (x^{t}) > α(ℓ_{t}/β)[Φ(x^{t}, x^{t}) − Φ(y^{t}, x^{t})]^{T}d^{t}, t ≥ \overline{t}, t ∈ T . (4.2)

By the convexity of u, we obtain in the limit of T in the characterization of y^{t} that [Φ(x^{∞}, x^{∞}) − Φ(y^{∞}, x^{∞})]^{T}d^{∞} ≥ [∇f(x^{∞}) + ξ_{u}(x^{∞})]^{T}d^{∞} holds for any subgradient ξ_{u}(x^{∞}) of u at x^{∞}. Hence, [Φ(x^{∞}, x^{∞}) − Φ(y^{∞}, x^{∞})]^{T}d^{∞} ≥ T′(x^{∞}; d^{∞}) holds. Dividing the inequality (4.2) by ℓ_{t}/β, and taking the limit of T , we obtain, taking into account the above inequality, Assumption 1.1 and [Roc81, Prop. 3G], that T′(x^{∞}; d^{∞}) ≥ 0 holds.

A contradiction is now reached in the limit of the inequality (cf. Lemma 3.2.a) T′(x^{t}; d^{t}) ≤ [∇f(x^{t}) + ξ_{u}(y^{t})]^{T}d^{t} = [Φ(x^{t}, x^{t}) − Φ(y^{t}, x^{t})]^{T}d^{t} < 0. We conclude that {d^{t}} → 0. Using this result in the characterization of y^{∞} then yields that x^{∞} ∈ Ω.

The result (4.1) follows as in Theorem 4.1. □

In order to introduce simpler step length rules, we assume that the mappings Φ^{t} are strongly monotone; we then automatically reach a large freedom of choosing their form; for example, they are not necessarily determined by a continuous mapping on dom u × dom u.

Theorem 4.4 (Convergence under strongly monotone mappings Φ^{t}) Assume that T is bounded from below on ℜ^{n} and ∇f is Lipschitz continuous on dom u. For each t, let Φ^{t} be strongly monotone and Lipschitz continuous on dom u. Assume further that m_{Φ} := lim inf_{t→∞}{m_{Φ^{t}}} > 0 and M_{Φ} := lim sup_{t→∞}{M_{Φ^{t}}} < +∞.

(a) From any starting point x^{0} ∈ dom u and under Rules E, A’ or R, any accumulation point of {x^{t}} lies in Ω.

(b) If further {x^{t}} is bounded and Ω is nonempty, then under Rules E, A’ or R, (4.1) holds. Further, under Rule D, at least one accumulation point of {x^{t}} lies in Ω.

(c) If further Ω is finite, then under Rule A’ or R, {x^{t}} converges to one point in Ω. The same conclusion holds for Rule E whenever {ℓ_{t}} is bounded above.

Proof.

(a) We begin by showing that for each of the three step length rules, {d^{t}} → 0 must hold. From (3.19) and the use of Rule E follow that, for ℓ ∈ [0, 1],

T (x^{t+1}) ≤ T (x^{t} + ℓd^{t}) ≤ T (x^{t}) + ℓ (−m_{Φ^{t}} + \frac{M_{∇f}}{2} ℓ) ∥d^{t}∥^{2}, t = 0, 1, ….

Minimizing the right-hand side of this inequality, we obtain ℓ_{t}^{∗} = min{1, m_{Φ^{t}}/M_{∇f}}. Hence,

T (x^{t+1}) − T (x^{t}) ≤ −ρ∥d^{t}∥^{2},  t = 0, 1, …,  (4.3)

with ρ := min{m_{Φ}^{2}/(2M_{∇f}), M_{∇f}/2}. The Armijo rule satisfies (4.3) with ρ := (2βα(1 − α)m_{Φ}^{2})/M_{∇f}, so does the relaxation rule, with ρ := ε_{1}ε_{2}M_{∇f}/2.

Since T is lower bounded and, from the above, {T (x^{t})} is decreasing, {T (x^{t+1}) − T (x^{t})} → 0. Thus, by (4.3), in each of the three rules we obtain that {d^{t}} → 0.

From (2.1) and the upper bound on M_{Φ^{t}}, as {d^{t}} → 0,

{∇f(x^{t}) + ξ_{u}(y^{t})} = {Φ^{t}(x^{t}) − Φ^{t}(y^{t})} → 0. (4.4)

Assume now that {x^{t}} has an accumulation point, x^{∞}, corresponding to a convergent subsequence {x^{t}}_{t∈T}. Since {y^{t}−x^{t}} = {d^{t}} → 0, we have that {y^{t}}_{t∈T} → x^{∞}. From the closedness of ∂u ([Roc70a, Thm. 24.4]) and the continuity of ∇f, (4.4) yields that {ξ_{u}(y^{t})}_{t∈T} → ξ_{u}(x^{∞}) ∈ ∂u(x^{∞}), and, again appealing to (4.4), ∇f(x^{∞})+ξ_{u}(x^{∞}) = 0. Thus, x^{∞} ∈ Ω.

(b) The result (4.1) follows from (a) and [OrR70, Thm. 14.1.4]. Since {ℓ_{t}} → 0 in Rule D, it holds for all t ≥ \overline{t} [cf. (3.19)] that, for some σ > 0, T (x^{t+1}) − T (x^{t}) ≤ −ℓ_{t}σ∥d^{t}∥^{2}. Hence, eventually {T (x^{t})} is decreasing, and by the lower boundedness assumption, it converges to a finite value. We therefore must have that σ ∑_{t=\overline{t}}^{∞} ℓ_{t}∥d^{t}∥^{2} < ∞. But ∑_{t=0}^{∞} ℓ_{t} = ∞ holds, whence it follows that lim inf_{t→∞} ∥d^{t}∥ = 0 must hold. The corresponding accumulation point of {x^{t}} must then, by the construction of [GE_{Φ}], lie in Ω.

(c) The result follows from the equality ∥x^{t+1} − x^{t}∥ = ℓ_{t}∥d^{t}∥, the result of (a) that {d^{t}} → 0, the upper bound on ℓ_{t} and [OrR70, Thm. 14.1.5]. □

The proof for Rule A’ shows that convergence is guaranteed for any step length rule that results in a larger reduction of T at each iteration than the Armijo Rule. The above result can be improved to require Lipschitz continuity only for Rule R and D ([Pat98]).

It is to be observed that although the strong monotonicity property of Φ implies that (3.2) holds for every x ∈ dom u, neither {x^{t}} nor {y^{t}} are necessarily bounded. Since {T (x^{t})} converges and {d^{t}} → 0, however, boundedness of {x^{t}} and {y^{t}} is ensured if the lower boundedness of T is replaced by the stronger condition of weak coercivity (dom T is bounded or lim_{∥x∥→∞}{T (x)} = +∞).

The maximal allowed step length in Rule R is bounded by the constant 2m_{Φ^{t}}/M_{∇f}. Now, suppose that f may be written as f := f_{1} + f_{2}, where f_{1} and f_{2} are convex with Lipschitz continuous gradients, and that f_{1} is strongly convex. For the sake of this example, we further assume that f_{1} is quadratic. To enhance the speed of convergence of the CA algorithm, we consider redefining φ^{t} as φ^{t} := φ^{t} + f_{1}. (This operation corresponds to applying a cost approximation to the function [u + f_{1}] + f_{2}, instead of to u + f.) Indeed, this operation yields an increased maximal allowed step length, since

\frac{2m_{(φ^{t}+f_{1})}}{M_{∇f_{2}}} = \frac{2(m_{φ^{t}} + m_{f_{1}})}{M_{∇(f−f_{1})}} > \frac{2m_{φ^{t}}}{M_{∇f}} .

It also yields a faster convergence for the other step length rules; for Rule A’, observe the change in the linear convergence ratio in Theorem 4.6 below.

From this example, we are lead to conclude that, in order to achieve the best possible rate of convergence, the function φ^{t} := φ^{t}+f_{1} should be constructed such that f_{2} := f −f_{1} is not strongly convex. (The desire to obtain a high rate of convergence must of course be weighed against the computational difficulty of the subproblems.)

Chen and Rockafellar [ChR92] study splitting methods for the problem of finding a zero of the sum of two maximal monotone operators. They argue that all the strong monotonicity inherent in the problem mapping should be kept in the mapping defining the forward step; their result corresponds to the above for this special case of CA method.

We conclude this section by combining inexact solutions of [GE_{Φ^{t}}] with inexact line searches into an implementable algorithm.

Theorem 4.5 (An implementable algorithm) Let the assumptions of Theorem 4.4 hold. In the CA algorithm, let Step 1 be performed inexactly, in the sense that, for each t, the approximate solution \overline{y}^{t} satisfies (3.4), with ∥r^{t}∥ ≤ ω∥\overline{d}^{t}∥ for some ω < m_{Φ}. Then, all the conclusions of Theorem 4.4 hold, with the exception that in step length Rule R, the constant m_{Φ} should be replaced by m_{Φ} − ω.

The proof of this result is very similar to that of Theorem 4.4; its basis is a modification of the inequality (3.19) resulting from (3.4), in which m_{Φ} is replaced by m_{Φ} − ω.

A slightly weaker result can be established, wherein the condition ∥r^{t}∥ ≤ ω∥\overline{d}^{t}∥ is replaced by the requirement that {∥r^{t}∥} → 0 holds; the argument is, roughly, that if ∥r^{t}∥ ≥ m_{Φ}∥\overline{d}^{t}∥ holds for infinitely many t, then \{∥\overline{d}^{t}∥\} → 0 must hold; if not, then Lemma 3.2.c establishes that \overline{d}^{t} eventually is a direction of descent.

4.4 Linear convergence

For affine mappings Φ, it is possible to show that the CA algorithm converges Q-linearly to an element x^{∗} of Ω, that is, that

lim sup_{t→∞} \frac{∥x^{t+1} − x^{∗}∥}{∥x^{t} − x^{∗}∥} =: q < 1. (4.5)

Convergence rests on a nonexpansiveness result for a part of the mapping x ↦ Y (x). Let Φ(x) := Bx for some positive definite matrix B ∈ ℜ^{n×n}. Then, from the charac-

terization (2.1) of the set Y (x), for any x ∈ dom u,

y(x) = [B + ∂u]^{−1}[B − ∇f](x). (4.6)

We then rewrite (4.6) as y(x) := V W (x), where

V := [\frac{1}{m_{Φ}} (B + ∂u)]^{−1} , W := [\frac{1}{m_{Φ}} (B − ∇f)] ,

and m_{Φ} is the modulus of strong monotonicity of Φ.

This rewriting of the characterization of y(x) is done, because the operator V can be shown to be nonexpansive on dom u, that is,

∥V (x) − V (y)∥ ≤ ∥x − y∥, ∀x, y ∈ dom u.

Indeed, from the monotonicity of ∂u and the positive definiteness of B we have, for all x, y ∈ dom u, and ξ_{u}(x) ∈ ∂u(x), ξ_{u}(y) ∈ ∂u(y),

[ξ_{u}(x) − ξ_{u}(y)]^{T}(x − y) ≥ 0,
(x − y)^{T}B(x − y) ≥ m_{Φ}∥x − y∥^{2}.

Dividing both the above inequalities by m_{Φ}, and adding them, we obtain

([\frac{1}{m_{Φ}} (B + ξ_{u})] (x) − [\frac{1}{m_{Φ}} (B + ξ_{u})] (y))^{T} (x − y) ≥ ∥x − y∥^{2}.

But this inequality yields that the operator V is firmly nonexpansive, and therefore nonexpansive ([Roc76, EcB92]).

The importance of this result is that is that the convergence analysis is simplified, since, using Lemma 3.2.a.1, we immediately have that for any x^{∗} ∈ Ω,

∥y(x) − x^{∗}∥ = ∥V W (x) − V W (x^{∗})∥ ≤ ∥W (x) − W (x^{∗})∥,

and it suffices to study the operator norm of W in the vicinity of x^{∗}.

Theorem 4.6 (Linear convergence) Assume that ∇f is Lipschitz continuous on dom u. Let Φ(x) := Bx, where B ∈ ℜ^{n×n} is positive definite. Assume that, using Rule A’ or Rule R, {x^{t}} → x^{∗} ∈ Ω, where f is twice continuously differentiable and ∇^{2}f(x^{∗}) is positive definite. If

σ := \frac{M_{Φ}}{m_{Φ}} ∥I − B^{−1}∇^{2}f(x^{∗})∥ < 1, (4.7)

then {x^{t}} converges Q-linearly with the ratio

q ≤ max \{σ, 1 − \frac{2β(1−α)m_{Φ}}{M_{∇f}} (1 − σ)\} , for Rule A’, and
q ≤ max \{σ, 1 − \underline{ℓ}(1 − σ)\}  for Rule R,

where \underline{ℓ} := inf_{t}\{ℓ_{t}\} ≥ ε_{1}.

Proof. From Lemma 3.2.a.1 and the nonexpansiveness property of V, we have ∥y^{t}−x^{∗}∥ = ∥V W (x^{t})−V W (x^{∗})∥ ≤ ∥W (x^{t})−W (x^{∗})−∇W (x^{∗})(x^{t} −x^{∗})∥+∥∇W (x^{∗})(x^{t} −x^{∗})∥, where ∇W (x^{∗}) := (1/m_{Φ})(B − ∇^{2}f(x^{∗})).

Let ε > 0 be arbitrary. Since {x^{t}} → x^{∗}, we have for a sufficiently large t that ∥W (x^{t}) − W (x^{∗}) − ∇W (x^{∗})(x^{t} − x^{∗})∥ ≤ ε∥x^{t} − x^{∗}∥. We also have that ∥∇W (x^{∗})∥ ≤ (M_{Φ}/m_{Φ})∥I−B^{−1}∇^{2}f(x^{∗})∥ =: σ, so that ∥y^{t}−x^{∗}∥ ≤ (σ+ε)∥x^{t}−x^{∗}∥. Then, ∥x^{t+1}−x^{∗}∥ = ∥x^{t} + ℓ_{t}d^{t} − x^{∗}∥ ≤ ℓ_{t}∥y^{t} − x^{∗}∥ + (1 − ℓ_{t})∥x^{t} − x^{∗}∥ ≤ [1 − ℓ_{t}(1 − σ − ε)]∥x^{t} − x^{∗}∥.

For Rule A’, from (3.18), we then have, for a large enough t,

∥x^{t+1} − x^{∗}∥ ≤ q_{ε}∥x^{t} − x^{∗}∥, (4.8)

where q_{ε} := max{σ + ε, 1 − 2β(1 − α)m_{Φ}/M_{∇f}(1 − σ − ε)}, while, for Rule R, (4.8) holds for large enough t with q_{ε} := max\{σ + ε, 1 − \underline{ℓ}(1 − σ − ε)\}. Since ε was arbitrary we have asymptotically that

q := lim sup_{t→∞} \frac{∥x^{t+1} − x^{∗}∥}{∥x^{t} − x^{∗}∥} ≤ q_{0},

where q_{0} := max{σ, 1 − 2β(1 − α)m_{φ}/M_{∇f}(1 − σ)} in the case of Rule A’, and q_{0} := max\{σ, 1 − \underline{ℓ}(1 − σ)\} for Rule R. If σ < 1, then q_{0} < 1, and the theorem is proved. □

If (3.20) holds, then ℓ_{t} ≡ 1 is a valid step length in Rule R, and q ≤ σ holds. Hence, if B := ∇^{2}f(x^{∗}), then superlinear convergence is obtained, that is, q = 0 in (4.5). For Rule A’, the same conclusion may be drawn, provided that (3.17) holds.

The result of Theorem 4.4 suggests a means to choosing a matrix B to obtain linear convergence in practice; choosing Φ^{t} of the form x ↦ B^{t}x for some positive definite matrix B^{t} which, at least asymptotically, approach the Hessian ∇^{2}f(x^{t}), eventually the requirement (4.7) will be fulfilled. Any matrix B^{t} generated after this occurrence will be an appropriate choice for the fixed matrix B in the above theorem. Further, the longer this choice is postponed, the closer will the value of σ be to zero; therefore, a convergence rate close to superlinear may be obtained. (In fact, the above result is possible to extend to incorporate general iteration-dependent mappings Φ^{t} [Pat98].)

We have not been able to establish linear convergence for Rule E or D, since no lower bound on the respective step length is available.

5 Applications of the convergence analysis

We illustrate in this section further examples of the generality of the framework, and demonstrate the strength of the convergence analysis. We demonstrate that the results strengthen those previously obtained for special cases of CA methods, and introduce a flexibility in the realizations of existing methods under the same convergence conditions.

5.1 On the linear convergence of CA algorithms

Despite the fact that the linear convergence Theorem 4.6 is valid only for iteration independent and affine cost approximation mappings, it is strong enough to reproduce some well known linear convergence results in differentiable optimization, as we illustrate below.

Corollary 5.1 (Linear convergence) Let the assumptions of Theorem 4.6 hold, and let B := (1/γ)I, γ > 0.

(a) Let Rule R be applied, using unit steps. If the method converges, then it converges with a linear rate. Further, the convergence ratio is

q ≤ max {|1 − γm|, |1 − γM|} , (5.1)

where m and M are the smallest and largest eigenvalue of ∇^{2}f(x^{∗}), respectively.

(b) Let further u ≡ δ_{X}, where X is a nonempty, closed and convex set in ℜ^{n}. Then, the algorithm is equivalent to the Goldstein–Levitin–Polyak gradient projection algorithm ([Gol64, LeP66]). Moreover, the convergence criterion (3.20) and the linear convergence ratio (5.1) coincide with those of this algorithm.

(c) Let further X := ℜ^{n}. Then, the algorithm is equivalent to the steepest descent method with fixed step length γ ([Pol63]), with the corresponding convergence criterion and linear convergence ratio.

Proof.

(a) We first note that m_{Φ} = M_{Φ} = 1/γ. The condition (3.20) for convergence using unit steps then reduces to γ < 2/M_{∇f}, and since ∥∇^{2}f(x^{∗})∥ ≤ M_{∇f} and for any symmetric matrix A, ∥I − γA∥ < 1 if γ < 2/∥A∥, (4.7) follows.

With ℓ := 1 we have that q ≤ σ := ∥I − γ∇^{2}f(x^{∗})∥ = max{|1 − m|, |1 − M|}, where m and M is, respectively, the smallest and the largest eigenvalue of γ∇^{2}f(x^{∗}).

(b) Straightforward calculations yield y^{t} = P_{X}(x^{t} − γ∇f(x^{t})), where P_{X}(z) denotes the Euclidean projection of z onto X. Since ℓ_{t} ≡ 1, x^{t+1} = y^{t} follows, and we obtain the Goldstein–Levitin–Polyak gradient projection algorithm. The convergence criterion (γ < 2/M_{∇f}) and linear convergence ratio (5.1) can be found in [LeP66].

(c) Follows from (b), with X = ℜ^{n}. □

5.2 Improvements of existing results

We provide a few examples of previous convergence analyses of special cases of CA algorithms, and demonstrate that the analysis in this paper improves upon them.

Application 1

In the CA algorithm, choose the following: φ^{t} := 0; each subproblem [NDP_{φ}] is solved exactly; and the conceptual step length Rule E is used. We then obtain the method of Mine and Fukushima [MiF81]. Comparing their convergence results to that of Theorem 4.1, their convergence conditions include an unnecessary strict convexity assumption on u; further, they do not recognize that u needs to be continuous. Needless to say, our convergence analysis also allows for a larger flexibility in the realization of this method under the same convergence conditions; for example, using inexact solutions of [NDP_{φ}] is optional (see Theorem 4.2).

When u ≡ δ_{X} for a closed, convex set X in ℜ^{n}, this algorithm reduces to the Frank– Wolfe method. For this algorithm, the validation of truncated subproblems (Theorem 4.2) is new, and its successful use in several applications (where it is often referred to as an heuristic procedure; see, e.g., [LHB85]) is thus supported in theory. Further, in this special case of [NDP], the original Armijo step length rule can be used in place of Rule A’, and convergence established without the need for Φ to be strictly monotone (see [Pat93c]). As far as the author is aware, previous convergence results for Armijo step length rules within the Frank–Wolfe algorithm have always assumed that ∇f is Lipschitz continuous (e.g., [PsD78, Chap. III.3]), which, however, is unnecessary (cf. Theorem 4.3).

Application 2

In the CA algorithm, choose the following: α > 0; φ^{t}(x) := 1/(2γ_{t})∥x∥^{2}, where {γ_{t}} is chosen so that inf_{t}{γ_{t}} > 0 and sup_{t}{γ_{t}} < 1/α holds; each subproblem [NDP_{φ^{t}}] is solved exactly; and the Armijo step length rule A’ is used. We then obtain the method of Fukushima and Mine [FuM81]. The condition that sup_{t}{γ_{t}} < 1/α holds is unnecessary in the presence of the standard assumption that the acceptance parameter satisfies α ∈ (0, 1), as evidenced by Lemma 3.3 and Theorem 4.4. Further, the enforced dependency of α on the possible choices of γ_{t} in their method has detrimental effects both in terms of allowing less flexibility in the realization of the method, and in terms of the convergence rate, which becomes arbitrarily poor if γ_{t} is chosen close to 1/α. Optional realizations of this algorithm include the use of inexact solutions of [NDP_{φ^{t}}] as well as the use of the step length Rules R and D (cf. Theorems 4.4 and 4.5). Moreover, the larger flexibility with which the function φ can be chosen introduces possibilities to improve the convergence rate of the method, as the next application establishes.

Application 3

Assume that f ≡ 0; replacing Rule A’ with Rule R using unit steps reduces the above method to the proximal point algorithm ([Mar70, Roc76]). If also u is differentiable, then by the results of Rockafellar [Roc76] the proximal point algorithm can be made superlinearly convergent by letting {γ_{t}} → +∞. Superlinear convergence is not possible to obtain in the method of [FuM81] in the nondifferentiable case, since letting {γ_{t}} → +∞ would destroy the validity of their step length rule. The discussion following Theorem 4.6 however demonstrates that also in the nondifferentiable case, a convergence rate as close to superlinear as desired can be obtained in the CA algorithm by introducing secondorder information from f into Φ; we conjecture that the method which is obtained if the matrix B^{t} is never fixed actually is superlinearly convergent. We note that this result is consistent with the results for gradient related methods in unconstrained differentiable optimization, where superlinear convergence is characterized by the convergence of the search directions toward Newton directions.

Application 4

In the CA algorithm, choose the following: φ^{t} is strongly convex with a Lipschitz continuous gradient; and Rule R is used. We then obtain the auxiliary problem principle of Cohen [Coh78, Coh80] (see, e.g., [Pat93c, LaP94] for detailed comparisons). The convergence results of this scheme correspond essentially to that of Theorem 4.4.b for this special choice of CA algorithm. There are several optional realizations of the algorithm with the same convergence conditions, including the use of step length Rule A’, the use of mappings Φ^{t} that are not gradients, and the use of truncated subproblem solutions.

In the case where u ≡ δ_{X} for a closed, convex set in ℜ^{n}, we may identify several other algorithm frameworks that are included in that defined by the CA method, and the convergence analyses are, in all cases, more limited and weaker than ours; examples of such frameworks include those in [DFL86, LaM90, Tse91, Mig94, ZhM95]. Applications of the convergence analysis of the CA method to differentiable optimization can be found in [Pat93b, Pat93c].

5.3 Concluding remarks

The class of CA algorithms is very general; some examples of methods that are included have been mentioned. What we have not mentioned is that the generality present allows one to adapt the algorithm to problem structures (such as choosing Φ to adapt to a separability in u or f); further, the analysis provides convergence results for new methods.

In this paper a convergence analysis was made for the class of CA algorithms when applied to a nondifferentiable optimization program; consequences of the analysis were illustrated, which resulted in several improvements of previous analyses. A convergence analysis is performed for the special case of differentiable, constrained and unconstrained optimization problems in the author’s PhD thesis [Pat93c]; in several instances, in particular concerning the requirements on Φ^{t} and the convergence rate, the results reached in this paper are improved upon. These, and new developments, will be reported elsewhere in the future ([Pat98]).

References

[Arm66] L. Armijo, Minimization of functions having Lipschitz continuous first partial derivatives, Pacific J. Math., 16 (1966), pp. 1–3.

[Ber63] C. Berge, Topological Spaces, Oliver & Boyd, Edinburgh, 1963.

[Bre73] H. Brézis, Opérateurs Maximaux Monotones et Semi-Groupes de Contractions dans les Espaces de Hilbert, North-Holland, Amsterdam, 1973.

[Bro68] F. E. Browder, Nonlinear maximal monotone mappings in Banach spaces, Math. Ann., 175 (1968), pp. 81–113.

[ChR92] G. H.-G. Chen and R. T. Rockafellar, Forward-backward splitting methods in Lagrangian optimization, unpublished report, Department of Applied Mathematics, University of Washington, Seattle, WA, 1992.

[Coh78] G. Cohen, Optimization by decomposition and coordination: a unified approach, IEEE Trans. Automat. Control, AC-23 (1978), pp. 222–232.

[Coh80] G. Cohen, Auxiliary problem principle and decomposition of optimization problems, J. Optim. Theory Appl., 32 (1980), pp. 277–305.

[CoZ84] G. Cohen and D. L. Zhu, Decomposition coordination methods in large scale optimization problems: the nondifferentiable case and the use of augmented Lagrangians, in Advances in Large Scale Systems, Volume 1, J. B. Cruz, ed., JAI Press, Greenwich, CT, 1984, pp. 203–266.

[CPS92] R. W. Cottle, J.-S. Pang and R. E. Stone, The Linear Complementarity Problem, Academic Press, San Diego, CA, 1992.

[DeV85] V. F. Dem’yanov and L. V. Vasil’ev, Nondifferentiable Optimization, Optimization Software, New York, NY, 1985.

[DFL86] J.-P. Dussault, J. A. Ferland, and B. Lemaire, Convex quadratic programming with one constraint and bounded variables, Math. Programming, 36 (1986), pp. 90–104.

[EcB92] J. Eckstein and D. P. Bertsekas, On the Douglas-Rachford splitting method and the proximal point algorithm for maximal monotone operators, Math. Programming, 55 (1992), pp. 293– 318.

[FHN96] M. Fukushima, M. Haddou, V. H. Nguyen, J.-J. Strodiot, T. Sugimoto and E. Yamakawa, A parallel descent algorithm for convex programming, Comput. Optim. Appl. 5 (1996), pp. 5–37.

[FuM81] M. Fukushima and H. Mine, A generalized proximal point algorithm for certain non-convex minimization problems, Internat. J. Systems Sci., 12 (1981), pp. 989–1000.

[Gol64] A. A. Goldstein, Convex programming in Hilbert space, Bull. Amer. Math. Soc., 70 (1964), pp. 709–710.

[LaM90] T. Larsson and A. Migdalas, An algorithm for nonlinear programs over Cartesian product sets, Optimization, 21 (1990), pp. 535–542.

[LaP94] T. Larsson and M. Patriksson, A class of gap functions for variational inequalities, Math. Programming, 64 (1994), pp. 53–79.

[LHB85] L. J. LeBlanc, R. V. Helgason, and D. E. Boyce, Improved efficiency of the Frank–Wolfe algorithm for convex network programs, Transportation Science, 19 (1985), pp. 445–462.

[LeP66] E. S. Levitin and B. T. Polyak, Constrained minimization methods, U.S.S.R. Comput. Math. and Math. Phys., 6 (1966), pp. 1–50.

[Lue84] D. G. Luenberger, Linear and Nonlinear Programming, Addison-Wesley, Reading, MA, second ed., 1984.

[Mar70] B. Martinet, Regularisation d’inéquations variationnelles par approximations successives, Revue Française d’Informatique et de Recherche Opérationelle, R-3 (1970), pp. 154–158.

[Mey79] G. G. L. Meyer, Asymptotic properties of sequences iteratively generated by point–to–set maps, Math. Programming Study, 10 (1979), pp. 115–127.

[Mig94] A. Migdalas, A regularization of the Frank–Wolfe method and unification of certain nonlinear programming methods, Math. Programming, 65 (1994), pp. 331–345.

[MiF81] H. Mine and M. Fukushima, A minimization method for the sum of a convex function and a continuously differentiable function, J. Optim. Theory Appl., 33 (1981), pp. 9–23.

[Min62] G. J. Minty, Monotone (nonlinear) operators in Hilbert space, Duke Math. J., 29 (1962), pp. 341–346.

[Nik68] H. Nikaido, Convex Structures and Economic Theory, Academic Press, New York, NY, 1968.

[OrR70] J. M. Ortega and W. C. Rheinboldt, Iterative Solution of Nonlinear Equations in Several Variables, Academic Press, New York, NY, 1970.

[Pat93a] M. Patriksson, Partial linearization methods in nonlinear programming, J. Optim. Theory Appl., 78 (1993), pp. 227–246.

[Pat93b] M. Patriksson, A unified description of iterative algorithms for traffic equilibria, European J. Oper. Res., 71 (1993), pp. 154–176.

[Pat93c] M. Patriksson, A unified framework of descent algorithms for nonlinear programs and variational inequalities, Doctoral dissertation, Department of Mathematics, Linköping Institute of Technology, Linköping, Sweden, 1993.

[Pat94a] M. Patriksson, On the convergence of descent methods for monotone variational inequalities, Oper. Res. Lett., 16 (1994), pp. 265–269.

[Pat94b] M. Patriksson, The Traffic Assignment Problem: Models and Methods, VSP, Utrecht, 1994.

[Pat97] M. Patriksson, Merit functions and descent algorithms for a class of variational inequality problems, Optimization,41 (1997), pp. 37–55.

[Pat98] M. Patriksson, Nonlinear Programming and Variational Inequalities: A Unified Approach, Kluwer Academic Publishers, under preparation.

[Pol63] B. T. Polyak, Gradient methods for the minimisation of functionals, U.S.S.R. Comput. Math. and Math. Phys., 3 (1963), pp. 864–878.

[Pol87] B. T. Polyak, Introduction to Optimization, Optimization Software, New York, NY, 1987.

[PsD78] B. N. Pshenichny and Yu. M. Danilin, Numerical Methods in Extremal Problems, MIR Publishers, Moscow, 1978.

[Rob79] S. M. Robinson, Generalized equations and their solutions, part I: basic theory, Math. Programming Study, 10 (1979), pp. 128–141.

[Roc66] R. T. Rockafellar, Characterization of the subdifferentials of convex functions, Pacific J. Math., 17 (1966), pp. 497–510.

[Roc69] R. T. Rockafellar, Local boundedness of nonlinear, monotone operators, Michigan Math. J., 16 (1969), pp. 397–407.

[Roc70a] R. T. Rockafellar, Convex Analysis, Princeton University Press, Princeton, NJ, 1970.

[Roc70b] R. T. Rockafellar, On the maximal monotonicity of subdifferential mappings, Pacific J. Math., 33 (1970), pp. 209–216.

[Roc70c] R. T. Rockafellar, On the maximality of sums of nonlinear monotone operators, Trans. Amer. Math. Soc., 149 (1970), pp. 75–88.

[Roc76] R. T. Rockafellar, Monotone operators and the proximal point algorithm, SIAM J. Control Optim., 14 (1976), pp. 877–898.

[Roc81] R. T. Rockafellar, The Theory of Subgradients and its Applications to Problems of Optimization: Convex and Nonconvex Functions, Heldermann Verlag, Berlin, 1981.

[Sho85] N. Z. Shor, Minimization Methods for Non-Differentiable Functions, Springer-Verlag, Berlin, 1985.

[StW70] J. Stoer and C. Witzgall, Convexity and Optimization in Finite Dimensions I, Springer-Verlag, Berlin, 1970.

[Tse91] P. Tseng, Decomposition algorithm for convex differentiable minimization, J. Optim. Theory Appl., 70 (1991), pp. 109–135.

[Zan69] W. I. Zangwill, Nonlinear Programming: A Unified Approach, Prentice-Hall, Englewood Cliffs, NJ, 1969.

[WoZ96] G. R. Wood and B. P. Zhang, Estimation of the Lipschitz constant of a function, J. Global Optim., 8 (1996), pp. 91–103.

[Zei90] E. Zeidler, Nonlinear Functional Analysis and its Applications II/B: Nonlinear Monotone Operators, Springer-Verlag, New York, NY, 1990.

[ZhM95] D. Zhu and P. Marcotte, Coupling the auxiliary problem principle with descent methods of pseudoconvex programming, European J. Oper. Res., 83 (1995), pp. 670–685.