A Path-Following Infeasible-Interior-Point Algorithm for Linear Complementarity Problems^{∗}

Stephen Wright^{†}

November 18, 1992

PREPRINT MCS-P334-1192, MATHEMATICS AND COMPUTER SCIENCE DIVISION, ARGONNE NATIONAL LABORATORY

Abstract

We describe an infeasible-interior-point algorithm for monotone linear complementarity problems that has polynomial complexity, global linear convergence, and local superlinear convergence with a Q-order of 2. Only one matrix factorization is required per iteration, and the analysis assumes only that a strictly complementary solution exists.

1 Introduction

The monotone linear complementarity problem is to find a vector pair (x, y) ∈ IR^{n} ×IR^{n} such that

y = Mx + h,  (x, y) ≥ (0, 0),  x^{T} y = 0,  (1)

where h ∈ IR^{n} and M is an n × n positive semidefinite matrix. A vector pair (x^{∗}, y^{∗}) is called a strictly complementary solution of (1) if it satisfies the three conditions in (1) and, in addition, x_{i}^{∗} + y_{i}^{∗} > 0 for each component i = 1, 2, ⋯, n. We denote the solution set for (1) by S and the set of strictly complementary solutions by S^{c}.

A number of interior point methods have been proposed for (1). Among recent papers are the predictor-corrector algorithm of Ji, Potra, and Huang [1] which has polynomial complexity and two-step superlinear convergence, the predictor-corrector algorithm of Ye and Anstreicher [5] which is polynomial and quadratically convergent, the path-following infeasible-interior-point algorithm of Zhang [6] which has polynomial complexity, and the algorithm described by the present author in an earlier report [4], which modifies Zhang’s algorithm so that superlinear convergence with Q-order 2 is attained. All these algorithms assume existence of a strictly feasible point for (1), that is, a vector pair (x, y) such that y = Mx + h and (x, y) > 0. Recently, Potra [2, 3] has proposed infeasible-interior-point algorithms for linear programming that do not require this assumption but still attain (in the case of [3]) O(nL) iteration complexity and two-step Q-quadratic convergence.

∗This research was supported by the Applied Mathematical Sciences subprogram of the Office of Energy Research, U. S. Department of Energy, under Contract W-31-109-Eng-38.

†MCS Division, Argonne National Laboratory, 9700 South Cass Avenue, Argonne, Illinois 60439.

In this paper, we describe an algorithm that is similar to the algorithm in [4] in that it takes so-called safe steps in the early stages of the algorithm and fast steps in the later stages (with an intermediate stage in which both kinds of steps can occur). The main point of difference between the two approaches is in the calculation of the safe steps. Essentially, the algorithm in [4] requires that, for safe steps, the infeasibility ∥y − Mx − h∥ decreases at least as quickly as the complementarity gap x^{T} y, while safe steps taken by the algorithm in this paper reduce these two quantities by exactly the same rate. We are therefore able to show that the iterates generated by the new algorithm are uniformly bounded, without an a priori assumption that a strictly feasible initial point exists.

Each iteration requires the solution of one, two, or three linear systems with the same coefficient matrix, so only one matrix factorization need be performed. In the final “fast” phase of the algorithm, just one solve is required at each step. The convergence properties are identical to the algorithm of [4], namely, global linear convergence, polynomial complexity when correctly initialized, and local superlinear convergence with Q-order 2.

Throughout the remainder of the paper, we make just one assumption:

Assumption 1 S^{c} is nonempty.

When (1) is derived from a linear programming problem, all that is needed for Assumption 1 to be satisfied is that a solution exists! Moreover, Ye and Anstreicher [5, Proposition 5.1] give a simple example to show that Assumption 1 is, in a sense, necessary for superlinear convergence. In their example, the solution set S is nonempty while S^{c} is empty, and an affine scaling approach fails to give superlinear convergence even when started arbitrarily close to the solution. The steps generated by all superlinearly convergent interior-point methods proposed to date resemble affine scaling steps more and more closely as they approach the limit, so this example suggests that a completely different approach will be needed to design polynomial and superlinear algorithms that do not require Assumption 1.

The remainder of this paper is laid out as follows: In Section 2, the algorithm is specified. Some technical results that provide global bounds on the step components are proved in Section 3. Global linear convergence and polynomial complexity are proved in Section 4. Although we focus on a particular choice of starting point — a choice that guarantees convergence in polynomial time — the linear convergence result can be proved when any starting point (x^{0}, y^{0}) with (x^{0}, y^{0}) > 0 is used. We do not work through the analysis for the latter case, since only minor modifications to the proofs of this paper are required. In Section 5, we obtain bounds on the step components in terms of the complementarity gap. Finally, we prove the superlinear convergence result in Section 6.

Unless otherwise specified, ∥ · ∥ denotes the Euclidean norm of a vector. Iteration indices (usually k) appear as superscripts on vectors and matrices and as subscripts on scalars. Subscripts are used to indicate components of vectors and matrices.

In Sections 3 and 4, we adopt the following convention: The iteration index k appears explicitly in the statement of each lemma and theorem but is usually omitted in the proofs, unless doing so would cause confusion. For example, quantities such as β, τ, μ, (x, y), and D appearing in a proof should be taken as β_{k}, τ_{k}, μ_{k}, (x^{k}, y^{k}), and D^{k}, respectively.

2 The Algorithm

Given a starting point with (x^{0}, y^{0}) > (0, 0), the algorithm generates a sequence of iterates (x^{k}, y^{k}) > (0, 0), k = 1, 2, · · ·. To describe the step between successive iterates, we define

μ_{k} = (x^{k})^{T} y^{k}/n, e = (1, 1, ⋯, 1)^{T} ,
X^{k} = diag(x_{1}^{k}, x_{2}^{k}, ⋯, x_{n}^{k}), Y^{k} = diag(y_{1}^{k}, y_{2}^{k}, ⋯, y_{n}^{k}).

Given a centering parameter σ̃ ∈ [0, 1/2], vector pairs (u^{k}, v^{k}) and (ū^{k}, v̄^{k}) are defined as the solutions of the linear systems

[ \begin{matrix}M & −I \\ Y^{k} & X^{k} \end{matrix} ] [ \begin{matrix}u^{k} \\ v^{k} \end{matrix} ] = [ \begin{matrix}y^{k} − Mx^{k} − h \\ −X^{k}Y^{k}e + σ̃μ_{k}e \end{matrix} ]  (2)

and

[ \begin{matrix}M & −I \\ Y^{k} & X^{k} \end{matrix} ] [ \begin{matrix} & ū^{k} \\  & v̄^{k} \end{matrix} ] = [ \begin{matrix}0 \\ −X^{k}Y^{k}e + σ̃μ_{k}e \end{matrix} ] .  (3)

From (u^{k}, v^{k}) and (ū^{k}, v̄^{k}) we can calculate the scalar quantities

τ_{k} = (u^{k})^{T} v^{k}/(nμ_{k}),  η_{k} = [(u^{k})^{T} v̄^{k} + (v^{k})^{T} ū^{k}]/(nμ_{k}),  δ_{k} = (ū^{k})^{T} v̄^{k}/(nμ_{k}). (4)

Note from (3) that, since M is positive semidefinite, δ_{k} ≥ 0.

When a safe step is calculated, σ̃ is restricted to the interval [σ̄, 1/2], where σ̄ ∈ (0, 1/2) is a constant. Safe steps have the form

(x^{k}(φ), y^{k}(φ)) = (x^{k}, y^{k}) + φ(u^{k}, v^{k}) + θ(ū^{k}, v̄^{k}),  (5)

where, for any given φ ∈ [0, 1], θ is chosen so that

x^{k}(φ)^{T} y^{k}(φ) = (1 − φ)nμ_{k} = (1 − φ)(x^{k})^{T} y^{k}.  (6)

It is easy to verify by using (2) and (3) that

y^{k}(φ) − Mx^{k}(φ) − h = (1 − φ)(y^{k} − Mx^{k} − h),  (7)

so we see that the effect of (6) is to ensure that the reductions in infeasibility and complementarity gap are exactly the same.

By substituting (4) and (5) into (6), we see that θ and φ must satisfy the polynomial

θ^{2}δ_{k} − θ(1 − σ̃ − η_{k}φ) + φ(σ̃ + τ_{k}φ) = 0.  (8)

For safe steps, we constrain φ to lie in the interval

φ ∈ [0, min (1, \frac{σ̃}{|τ_{k}|}, \frac{1}{8|η_{k}|}, \frac{9}{512δ_{k}} )] ,  (9)

where, if any of τ_{k}, η_{k}, and δ_{k} are zero, the corresponding term is omitted from the formula for the upper bound. The effect of (9) is to ensure existence of a θ that satisfies (8), as can be seen from the following result:

Lemma 2.1 When σ̃ ∈ [0, 1/2] and φ satisfies (9), there is at least one nonnegative value of θ such that (8) is satisfied.

Proof. We write the polynomial (8) as

a(φ)θ^{2} − b(φ)θ + c(φ) = 0,

where, using (9),

a(φ) = δ_{k}
b(φ) = 1 − σ̃ − η_{k}φ ∈ [\frac{3}{8}, \frac{9}{8}] ,
c(φ) = φ[σ̃ + τ_{k}φ] ∈ [0, 2σ̃φ].

When δ_{k} = 0, we obtain

θ = \frac{c(φ)}{b(φ)} ≥ 0,

as required. Otherwise, δ_{k} > 0, and we have

θ = \frac{b(φ)}{2a(φ)} [1 ± \sqrt{1 − \frac{4a(φ)c(φ)}{b(φ)^{2}}} ] .

Since

\frac{b(φ)}{2a(φ)} > 0,

and

|\frac{4a(φ)c(φ)}{b(φ)^{2}} | = \frac{4a(φ)c(φ)}{b(φ)^{2}} ≤ \frac{4δ_{k}(2σ̃φ)}{(3/8)^{2}} ≤ \frac{(9/64)σ̃}{(9/64)} ≤ \frac{1}{2},

we find that both roots of (8) are real and nonnegative. If there are two roots for θ, we resolve the ambiguity by choosing the smaller one, and define (x^{k}(φ), y^{k}(φ)) accordingly.

In addition to satisfying (6), the actual step length φ_{k} should be chosen so that all vector pairs (x^{k}(φ), y^{k}(φ)), φ ∈ [0, φ_{k}], lie in the neighborhood N (γ̃) of the central path, where

N (γ̃) = {(x, y) ≥ 0 | x_{i}y_{i} ≥ γ̃(x^{T} y/n)},

and γ̃ ∈ (0, 1/2) is a given parameter.

The procedure for calculating a safe step can be stated formally as follows:

Safe step calculation: Given σ̃ ∈ [σ̄, 1/2] and γ̃ ∈ (0, 1/2],

Solve (2) and (3) to find (u^{k}, v^{k}), (ū^{k}, v̄^{k}); Choose φ̃ as the largest value in the interval (9) such that

(x^{k}(φ), y^{k}(φ)) ∈ N (γ̃)  (10)

for all φ ∈ [0, φ̃]; Return φ̃ and the corresponding θ̃.

Fast step calculations do not require (3) to be solved, and they replace the condition (6) by a relaxed condition, in which the decrease in complementarity gap is allowed to vary a little from the decrease in infeasibility. The procedure can be outlined formally as follows:

Fast step calculation: Given σ̃ ∈ [0, 1/2], γ̃ ∈ (0, 1/2], and β̃ ∈ [0, 1/2],

Solve (2) to find (u^{k}, v^{k}); Choose φ̃ as the largest value in [0, 1] such that

(x^{k}, y^{k}) + φ(u^{k}, v^{k}) ∈ N (γ̃),  (11)

and

(1 − β̃)(1 − φ)(x^{k})^{T} y^{k} ≤ (x^{k} + φu^{k})^{T} (y^{k} + φv^{k}) ≤ (1 + β̃)(1 − φ)(x^{k})^{T} y^{k},  (12)

for all φ ∈ [0, φ̃]. Return φ̃.

To obtain polynomial complexity, we define the starting point for our main algorithm as

(x^{0}, y^{0}) = (ξ_{x}e, ξ_{y}e),  (13)

where e is the vector whose components are all 1, and the scalars ξ_{x} and ξ_{y} are chosen to satisfy

ξ_{x} > 0,  ξ_{y} > 0,  ξ_{y} ≥ ∥h∥_{∞},  ξ_{y} ≥ ∥Me∥_{∞}ξ_{x}.  (14)

We can now state our algorithm.

Given γ̄ ∈ (0, 1/2), σ̄ ∈ (0, 1/2), ρ ∈ (0, γ̄), μ̄ > 0, and (x^{0}, y^{0}) defined by (13),(14)

t_{0} ← 1, γ_{0} ← 2γ̄;

for k = 0, 1, 2, ⋯

if μ_{k} ≤ μ̄ then Set σ̃ = μ_{k}, γ̃ = γ̄(1 + γ̄^{t_{k}}), and β̃ = γ̄^{t_{k}} and do a fast step calculation;

if  φ̃ ≥ (1 − ρ)

then φ_{k} ← φ̃, θ_{k} ← 0, γ_{k+1} ← γ̃, β_{k} = β̃, σ_{k} ← σ̃; t_{k+1} ← t_{k} + 1; (x^{k+1}, y^{k+1}) ← (x^{k}, y^{k}) + φ_{k}(u^{k}, v^{k}); go to next k;

end if end if

Choose σ̃ ∈ [σ̄, 1/2], set γ̃ = γ_{k}, and do a safe step calculation; φ_{k} ← φ̃, θ_{k} ← θ̃, γ_{k+1} ← γ̃, β_{k} = 0, σ_{k} ← σ̃; t_{k+1} ← t_{k}; (x^{k+1}, y^{k+1}) ← (x^{k}, y^{k}) + φ_{k}(u^{k}, v^{k}) + θ_{k}(ū^{k}, v̄^{k}); go to next k;

end for.

The algorithm takes safe steps until the complementarity gap μ_{k} falls below a threshold value μ̄. It then attempts to take fast steps in which the centering component of (u^{k}, v^{k}) is phased out by decreasing σ_{k} to zero in concert with μ_{k}. These fast steps are accepted only if the value φ̃ chosen by the step procedure exceeds 1 − ρ; otherwise, we revert to a safe step. The size of the central path neighborhood N (γ_{k}) is increased slightly on each fast iterate, to allow φ_{k} to approach 1.

3 Technical Results

In this section, we find bounds on the iterates (x^{k}, y^{k}), the step components (u^{k}, v^{k}) and (ū^{k}, v̄^{k}), and the scalars τ_{k}, η_{k}, and δ_{k} defined in (4). We define the infeasibility vector r^{k} as

r^{k} = y^{k} − Mx^{k} − h,  (15)

and note that, as observed in (7),

r^{k+1} = (1 − φ_{k})r^{k}, k = 0, 1, 2, · · · .

If we define the scalar ν_{k} by

ν_{k} = ∏_{j=0}^{k−1}(1 − φ_{j}),  (16)

it follows that

r^{k} = ν_{k}r^{0},  k = 0, 1, 2, · · · .  (17)

We can bound μ_{k} in terms of μ_{0} and ν_{k} by the following lemma:

Lemma 3.1 For all k = 0, 1, ⋯,

ν_{k}μ_{0} ∏_{j=0}^{k−1}(1 − β_{j}) ≤ μ_{k} ≤ ν_{k}μ_{0} ∏_{j=0}^{k−1}(1 + β_{j}).

Defining

β_{L} =^{△} ∏_{j=0}^{∞}(1 − β_{j}), β_{U} =^{△} ∏_{j=0}^{∞}(1 + β_{j}),

we have that

e^{−3/2} < β_{L} < 1 < β_{U} < e,

and

β_{L}ν_{k}μ_{0} ≤ μ_{k} ≤ β_{U}ν_{k}μ_{0}.  (18)

Proof. The first inequalities follow immediately from the algorithm.

Since we have either β_{k} = 0 (safe steps) or β_{k} = γ̄^{t}, t = 1, 2, ⋯ (fast steps), it follows that

β_{L} ≥ ∏_{t=1}^{∞}(1 − γ̄^{t}), β_{U} ≤ ∏_{t=1}^{∞}(1 + γ̄^{t}).

Now, since γ̄ ∈ (0, 1/2), and since

log(1 − α) ≥ (−3/2)α,  log(1 + α) ≤ α,  for α ∈ [0, 1/2],

we have

log ∏_{t=1}^{∞}(1 + γ̄^{t}) = ∑_{t=1}^{∞} log(1 + γ̄^{t}) ≤ ∑_{t=1}^{∞} γ̄^{t} = γ̄/(1 − γ̄).

Therefore

β_{U} ≤ e^{γ̄/(1−γ̄)} < e.

Similarly,

log ∏_{t=1}^{∞}(1 − γ̄^{t}) = ∑_{t=1}^{∞} log(1 − γ̄^{t}) ≥ (−3/2) ∑_{t=1}^{∞} γ̄^{t} = (−3/2)γ̄/(1 − γ̄).

Therefore

β_{L} ≥ e^{(−3/2)γ̄/(1−γ̄)} > e^{−3/2}.

The following theorem provides a bound on the iteration sequence. It is similar to Lemma 4.1 of Potra [2].

Lemma 3.2 If (x^{∗}, y^{∗}) ∈ S and ξ_{x} and ξ_{y} are defined as in (13),(14), the vector pair (x^{k}, y^{k}) satisfies

ξ_{x}∥y^{k}∥_{1} + ξ_{y}∥x^{k}∥_{1} ≤ nμ_{0} + nμ_{k}/ν_{k} + ξ_{x}∥y^{∗}∥_{1} + ξ_{y}∥x^{∗}∥_{1} ≤ (β_{U} + 1)nμ_{0} + ξ_{x}∥y^{∗}∥_{1} + ξ_{y}∥x^{∗}∥_{1}.

Proof. Using (15), (17), and y^{∗} = Mx^{∗} + h, we have

M(νx^{0} + (1 − ν)x^{∗} − x)
= νMx^{0} + (1 − ν)Mx^{∗} − Mx
= ν(y^{0} − h − r^{0}) + (1 − ν)(y^{∗} − h) − (y − h − r)
= νy^{0} + (1 − ν)y^{∗} − y.

Therefore, by positive semidefiniteness of M,

0 ≤ (νx^{0} + (1 − ν)x^{∗} − x)^{T} (νy^{0} + (1 − ν)y^{∗} − y)
= ν^{2}(x^{0})^{T} y^{0} + (1 − ν)^{2}(x^{∗})^{T} y^{∗} + x^{T} y + ν(1 − ν)((x^{0})^{T} y^{∗} + (x^{∗})^{T} y^{0})
−ν((x^{0})^{T} y + x^{T} y^{0}) − (1 − ν)((x^{∗})^{T} y + x^{T} y^{∗}).

Now, since

x^{0} = ξ_{x}e,  y^{0} = ξ_{y}e,  (x^{∗})^{T} y^{∗} = 0,  (x^{∗})^{T} y + x^{T} y^{∗} ≥ 0,

we have

ν(ξ_{x}∥y∥_{1} + ξ_{y}∥x∥_{1}) ≤ ν^{2}(x^{0})^{T} y^{0} + x^{T} y + ν(1 − ν)(ξ_{x}∥y^{∗}∥_{1} + ξ_{y}∥x^{∗}∥_{1})
⇒ ξ_{x}∥y∥_{1} + ξ_{y}∥x∥_{1} ≤ nμ_{0} + nμ/ν + ξ_{x}∥y^{∗}∥_{1} + ξ_{y}∥x^{∗}∥_{1},

giving the first inequality. The second inequality follows directly from (18).

If we define

D^{k} = (X^{k})^{−1/2}(Y^{k})^{1/2},  (19)

then bounds on ∥D^{k}ū^{k}∥ and ∥(D^{k})^{−1}v̄^{k}∥ are derived in the following lemma.

Lemma 3.3

∥D^{k}ū^{k}∥ ≤ \frac{2}{γ̄^{1/2}} nμ_{k}^{1/2},
∥(D^{k})^{−1}v̄^{k}∥ ≤ \frac{2}{γ̄^{1/2}} nμ_{k}^{1/2}.

Proof. From the equations (3) and ū^{T} v̄ ≥ 0, we have

Dū + D^{−1}v̄ = (XY)^{−1/2}(σ̃μe − XY e)
⇒ ∥Dū∥^{2} + ∥D^{−1}v̄∥^{2} ≤ ∥(XY)^{−1/2}(σ̃μe − XY e)∥^{2}
⇒ ∥Dū∥ ≤ ∥(XY)^{−1/2}(σ̃μe − XY e)∥.

Now

∥(XY)^{−1/2}(σ̃μe − XY e)∥ ≤ ∥b∥_{∞}∥σ̃μe − XY e∥_{1},

where b is the vector whose components are (x_{i}y_{i})^{−1/2}. Since (x^{k}, y^{k}) ∈ N (γ_{k}) ⊂ N (γ̄), we have

∥Du∥ ≤ sup_{i}(x_{i}y_{i})^{−1/2}∥σ̃μe − XY e∥_{1} ≤ \frac{1}{(γ̄μ)^{1/2}} (1 + σ̃)nμ ≤ \frac{2}{(γ̄μ)^{1/2}} nμ,

giving the result.

Bounds on ∥D^{k}u^{k}∥ and ∥(D^{k})^{−1}v^{k}∥ are not much harder to find.

Lemma 3.4 For some (x^{∗}, y^{∗}) ∈ S, define

ξ_{x}^{∗} = \frac{1}{nξ_{x}} ∥x^{∗}∥_{1},  ξ_{y}^{∗} = \frac{1}{nξ_{y}} ∥y^{∗}∥_{1}.  (20)

Then

∥D^{k}u^{k}∥ ≤ \frac{3n}{γ̄^{1/2}} μ_{k}^{1/2}(2 + β_{L}^{−1}(1 + ξ_{x}^{∗} + ξ_{y}^{∗}))  (21a)
∥(D^{k})^{−1}v^{k}∥ ≤ \frac{3n}{γ̄^{1/2}} μ_{k}^{1/2}(2 + β_{L}^{−1}(1 + ξ_{x}^{∗} + ξ_{y}^{∗}).  (21b)

Proof. We prove only (21a), since the proof of (21b) is nearly identical. We start by defining the vector pair (û, v̂) to be the solution of

[ \begin{matrix}M & −I \\ Y^{k} & X^{k} \end{matrix} ] [ \begin{matrix} & û^{k} \\  & v̂^{k} \end{matrix} ] = [ \begin{matrix}y^{k} − Mx^{k} − h \\ 0 \end{matrix} ] .  (22)

We have immediately from (22) that

Y û + Xv̂ = 0 ⇒ Dû + D^{−1}v̂ = 0 ⇒ Dû = −D^{−1}v̂,

and so ∥Dû∥ = ∥D^{−1}v̂∥. Substituting v̂ = −D^{2}û into the first block row of (22), we obtain

Mû − v̂ = r ⇒ (M + D^{2})û = r ⇒ û^{T} D^{2}û ≤ û^{T} r.

Therefore

∥Dû∥^{2} ≤ ∥Dû∥∥D^{−1}r∥ ⇒ ∥Dû∥ ≤ ∥D^{−1}r∥ = ν∥D^{−1}r^{0}∥.  (23)

Now

∥D^{−1}r^{0}∥ ≤ ∥d^{−1}∥_{1}∥y^{0} − Mx^{0} − h∥_{∞},

where d^{−1} is the vector with components (x_{i}/y_{i})^{1/2}. We have

∥d^{−1}∥_{1} = ∑_{i=1}^{n} (\frac{x_{i}}{y_{i}} )^{1/2} = ∑_{i=1}^{n} \frac{x_{i}}{(x_{i}y_{i})^{1/2}} ≤ \frac{1}{(γ̄μ)^{1/2}} ∥x∥_{1},

while

∥r^{0}∥_{∞} ≤ ∥y^{0}∥_{∞} + ∥Mx^{0}∥_{∞} + ∥h∥_{∞} ≤ ξ_{y} + ∥Me∥_{∞}ξ_{x} + ∥h∥_{∞}.

Because of the conditions (14), it follows that

∥r^{0}∥_{∞} ≤ 3ξ_{y}.

By substituting in (23), we have

∥Dû∥ ≤ ν \frac{3}{(γ̄μ)^{1/2}} ξ_{y}∥x∥_{1}.

From Lemma 3.2, we deduce that

∥Dû∥ ≤ \frac{3ν}{(γ̄μ)^{1/2}} (nμ_{0} + nμ/ν + ξ_{x}∥y^{∗}∥_{1} + ξ_{y}∥x^{∗}∥_{1}) ≤ \frac{3}{(γ̄μ)^{1/2}} nμ(1+β_{L}^{−1}(1+ξ_{y}^{∗}+ξ_{x}^{∗})), (24)

where the second inequality follows from (18) and the inequalities

νξ_{x}∥y^{∗}∥_{1} = nνξ_{x}ξ_{y}ξ_{y}^{∗} = nνμ_{0}ξ_{y}^{∗} ≤ nμξ_{y}^{∗}/β_{L}

and νξ_{y}∥x^{∗}∥_{1} ≤ nμξ_{x}^{∗}/β_{L}. Noting that

(u, v) = (û, v̂) + (ū, v̄),

we can combine (24) with the results of Lemma 3.3 to obtain

∥Du∥ ≤ ∥Dû∥ + ∥Dū∥ ≤ \frac{3n}{γ̄^{1/2}} μ^{1/2}(2 + β_{L}^{−1}(1 + ξ_{y}^{∗} + ξ_{x}^{∗})),

as required.

In the next theorem, we use Lemmas 3.3 and 3.4 to obtain bounds on a number of products involving components of the vector pairs (u^{k}, v^{k}) and (ū^{k}, v̄^{k}) that will be useful in the next section.

Theorem 3.5 Suppose we define a constant τ̃ independent of n by

τ̃ = \frac{9}{γ̄} (2 + β_{L}^{−1}(1 + ξ_{x}^{∗} + ξ_{y}^{∗}))^{2}.  (25)

Then for all k and for i = 1, 2, ⋯, n, we have

|u_{i}^{k}v_{i}^{k}| ≤ n^{2}τ̃μ_{k},
|ū_{i}^{k}v_{i}^{k} + u_{i}^{k}v̄_{i}^{k}| ≤ n^{2}τ̃μ_{k},
|ū_{i}^{k}v̄_{i}^{k}| ≤ n^{2}τ̃μ_{k}.

Also

|τ_{k}| = |\frac{(u^{k})^{T} v^{k}}{nμ_{k}} | ≤ nτ̃,
|η_{k}| = |\frac{(u^{k})^{T} v̄^{k} + (ū^{k})^{T} v^{k}}{nμ_{k}} | ≤ nτ̃,
0 ≤ δ_{k} = \frac{(ū^{k})^{T} v̄^{k}}{nμ_{k}} ≤ nτ̃.

Proof. These inequalities follow immediately from Lemmas 3.2, 3.3, and 3.4, with the help of the triangle inequality. For example, the second inequality is derived as follows.

|ū_{i}v_{i} + u_{i}v̄_{i}| ≤ |d_{i}ū_{i}||d_{i}^{−1}v_{i}| + |d_{i}u_{i}||d_{i}^{−1}v̄_{i}|,

where d_{i} = (y_{i}/x_{i})^{1/2} is the i-th diagonal element of D. Therefore

|ū_{i}v_{i} + u_{i}v̄_{i}| ≤ ∥Dū∥∥D^{−1}v∥ + ∥Du∥∥D^{−1}v̄∥
≤ 2\frac{3}{γ̄^{1/2}} (2 + β_{L}^{−1}(1 + ξ_{x}^{∗} + ξ_{y}^{∗}))nμ^{1/2} \frac{2}{γ̄^{1/2}} nμ^{1/2}.

Since ξ_{x}^{∗} ≥ 0, and ξ_{y}^{∗} ≥ 0, we have 4 < 3(2 + β_{L}^{−1}(1 + ξ_{x}^{∗} + ξ_{y}^{∗})), and so

|ū_{i}v_{i} + u_{i}v̄_{i}| ≤ \frac{9}{γ̄} (2 + β_{L}^{−1}(1 + ξ_{x}^{∗} + ξ_{y}^{∗}))^{2}n^{2}μ = n^{2}τ̃μ.

The other five inequalities are derived similarly.

Finally, we state a purely technical result that is used in Sections 4 and 6.

Lemma 3.6

\frac{1}{1 + α} ≥ 1 − α,  all α ≥ 0,  (26)
\sqrt{1 − α} ≤ 1 − \frac{1}{2}α,  all α ∈ [0, 1],  (27)
\sqrt{1 − α} ≥ 1 − \frac{1}{2}α − \frac{1}{2}α^{2},  all α ∈ [0, 1/2].  (28)

Proof. For (26), we have

\frac{1}{1 + α} = 1 − α + α^{2} − α^{3} + ⋯
= 1 − α + α^{2}(1 − α + α^{2} − ⋯)
= 1 − α + \frac{α^{2}}{1 + α}
≥ 1 − α.

For (27),

\sqrt{1 − α} = 1 − \frac{1}{2}α − \frac{1}{8}α^{2} − \frac{3}{48}α^{3} − ⋯ ≤ 1 − \frac{1}{2}α.

For (28), note that

\sqrt{1 − α} = 1 − ∑ c_{k}α^{k},

where

c_{k} = \frac{(\frac{1}{2})(\frac{1}{2})(\frac{3}{2})(\frac{5}{2}) ⋯ (k terms)}{k!} ≤ \frac{1}{4}.

Hence

\sqrt{1 − α} ≥ 1 − \frac{1}{2}α − \frac{1}{4}α^{2}(1 + α + α^{2} + ⋯) ≥ 1 − \frac{1}{2}α − \frac{1}{2}α^{2},

where the final inequality follows from α ∈ [0, 1/2].

4 Linear Convergence and Polynomial Complexity

In this section, we prove that the algorithm converges globally at a Q-linear rate and that, for sufficiently large choice of ξ_{x} and ξ_{y}, it has polynomial complexity.

The main result of this section shows that φ̄ defined by

φ̄ = \frac{.01σ̄(1 − 2γ̄)}{n^{2}τ̃}  (29)

is a lower bound on φ_{k} for all iterations k on which a safe step is taken. It is clear that φ̄ satisfies the bound (9) since, by Theorem 3.5 and the fact that σ̃ ≥ σ̄,

min (1, \frac{σ̃}{|τ_{k}|}, \frac{1}{8|η_{k}|}, \frac{9}{512|δ_{k}|}) ≥ \frac{9σ̄}{512nτ̃} ≥ φ̄.  (30)

It follows from Lemma 2.1 that there is a nonnegative θ satisfying (8) for each φ ∈ [0, φ̄]. In the case of δ_{k} = 0, we have the single root

θ = \frac{σ̃ + τ_{k}φ}{1 − σ̃ − η_{k}φ}φ,  (31)

while, when δ_{k} > 0, the root of interest to us (i.e., the smaller one) is

θ = \frac{1 − σ̃ − η_{k}φ}{2δ_{k}} [1 − \sqrt{1 − \frac{4δ_{k}φ(σ̄ + τ_{k}φ)}{(1 − σ̃ − η_{k}φ)^{2}}} ] .  (32)

In the following result, we define the range in which θ lies, in terms of φ.

Lemma 4.1 For φ ∈ [0, φ̄] and σ̃ ∈ [0, 1/2], the value of θ corresponding to φ satisfies

θ ≥ \frac{σ̃ + τ_{k}φ}{1 − σ̃ − η_{k}φ}φ  (33)

and

θ ≤ \frac{σ̃ + τ_{k}φ}{1 − σ̃ − η_{k}φ}φ [1 + \frac{4δ_{k}φ(σ̃ + τ_{k}φ)}{(1 − σ̃ − η_{k}φ)^{2}} ] ≤ 2φ.  (34)

Proof. When δ_{k} = 0, the inequalities (33) and (34) follow immediately from (31).

When δ_{k} > 0, we obtain the result by appealing to (32) and Lemma 3.6. From (30), we have that φ satisfies (9), so as in the proof of Lemma 2.1 we have

0 ≤ \frac{4δφ(σ̃ + τφ)}{(1 − σ̃ − ηφ)^{2}} ≤ \frac{1}{2}.  (35)

Therefore, from (27)

\sqrt{1 − \frac{4δφ(σ̃ + τφ)}{(1 − σ̃ − ηφ)^{2}}} ≤ 1 − \frac{2δφ(σ̃ + τφ)}{(1 − σ̃ − ηφ)^{2}} ,  (36)

while from (28)

\sqrt{1 − \frac{4δφ(σ̃ + τφ)}{(1 − σ̃ − ηφ)^{2}}} ≥ 1 − \frac{2δφ(σ̃ + τφ)}{(1 − σ̃ − ηφ)^{2}} − \frac{1}{2} [\frac{4δφ(σ̃ + τφ)}{(1 − σ̃ − ηφ)^{2}} ]^{2} .  (37)

It follows from (32) and (36) that

θ ≥ \frac{1 − σ̃ − ηφ}{2δ} \frac{2δφ(σ̃ + τφ)}{(1 − σ̃ − ηφ)^{2}} = \frac{σ̃ + τφ}{1 − σ̃ − ηφ}φ,

proving (33). From (32) and (37), we have

θ ≤ \frac{1 − σ̃ − ηφ}{2δ} [\frac{2δφ(σ̃ + τφ)}{(1 − σ̃ − ηφ)^{2}} + \frac{1}{2} [\frac{4δφ(σ̃ + τφ)}{(1 − σ̃ − ηφ)^{2}} ]^{2}]

≤ \frac{σ̃ + τφ}{1 − σ̃ − ηφ}φ [1 + \frac{4δφ(σ̃ + τφ)}{(1 − σ̃ − ηφ)^{2}} ] .

Now, by using σ̃ ≤ .5, (35), and the inequalities

|φτ| ≤ φ̄nτ̃ ≤ .01,  |φη| ≤ φ̄nτ̃ ≤ .01,  |φδ| ≤ φ̄nτ̃ ≤ .01,  (38)

we have

θ ≤ \frac{.51}{.49}φ [1 + \frac{1}{2}] ≤ 2φ,

proving (34).

To show that the central neighborhood condition (10) holds for φ ∈ [0, φ̄], we need to show that

(x_{i}^{k} + φu_{i}^{k} + θū_{i}^{k})(y_{i}^{k} + φv_{i}^{k} + θv̄_{i}^{k}) ≥ γ̃(1 − φ)μ_{k},  (39)

where, since this is a safe step, γ̃ = γ_{k}. From (2) and (3), we have

u_{i}^{k}y_{i}^{k} + x_{i}^{k}v_{i}^{k} = ū_{i}^{k}y_{i}^{k} + v̄_{i}^{k}x_{i}^{k} = −x_{i}^{k}y_{i}^{k} + σ̃μ_{k}.

Substituting these equations in (39) and dropping the iteration index k, we find that (10) is satisfied if

x_{i}y_{i}(1 − φ − θ) + (φ + θ)σ̃μ + φ^{2}u_{i}v_{i} + θ^{2}ū_{i}v̄_{i} + θφ(u_{i}v̄_{i} + ū_{i}v_{i}) ≥ γ(1 − φ)μ.  (40)

We can now use the fact that x_{i}y_{i} ≥ γμ, together with the bounds from Theorem 3.5, to deduce that (40) holds if

γμ(1 − φ − θ) + (φ + θ)σ̃μ − (θ^{2} + φ^{2} + θφ)n^{2}τ̃μ ≥ γ(1 − φ)μ,

or, equivalently,

− θγ + (φ + θ)σ̃ ≥ (θ^{2} + φ^{2} + θφ)n^{2}τ̃.  (41)

In the next two lemmas, we derive a lower bound for the left-hand side of this inequality and an upper bound for its right-hand side.

Lemma 4.2 If σ̃ ∈ [σ̄, 1/2] and φ ∈ [0, φ̄], then

− θγ_{k} + (φ + θ)σ̃ ≥ (.80)φσ̄(1 − 2γ̄),  (42)

where θ is the smaller root of (8).

Proof. From (33) and (34), we have

−θγ + (φ + θ)σ̃

≥ −γφ\frac{σ̃ + τφ}{1 − σ̃ − ηφ} [1 + \frac{4δφ(σ̃ + τφ)}{(1 − σ̃ − ηφ)^{2}} ] + φ [1 + \frac{σ̃ + τφ}{1 − σ̃ − ηφ}] σ̃
= φ\frac{σ̃(1 − γ)}{1 − σ̃ − ηφ} − φ^{2} \frac{σ̃η − σ̃τ + γτ}{1 − σ̃ − ηφ} − γφ^{2} \frac{4δ(σ̃ + τφ)^{2}}{(1 − σ̃ − ηφ)^{3}} .  (43)

Now, from Theorem 3.5, (38), σ̃ ∈ [σ̄, 1/2], and γ ≤ 2γ̄ < 1, we obtain

φ |\frac{σ̃η − σ̃τ + γτ}{1 − σ̃ − ηφ} | ≤ \frac{(1 + 2γ̄)nτ̃φ̄}{(.49)} ≤ \frac{(1 + 2γ̄)(.01)σ̄(1 − 2γ̄)}{(.49)} ≤ (.05)σ̄(1 − 2γ̄),  (44)

and

φ\frac{4δ(σ̃ + τφ)^{2}}{(1 − σ̃ − ηφ)^{3}} ≤ \frac{4nτ̃φ̄[.50 + .01]^{2}}{(.49)^{3}} ≤ 12σ̄(.01)(1 − 2γ̄) = .12σ̄(1 − 2γ̄),  (45)

while

\frac{σ̃(1 − γ)}{1 − σ̃ − ηφ} ≥ \frac{σ̄(1 − 2γ̄)}{1.01} ≥ (.99)σ̄(1 − 2γ̄).  (46)

By substituting (44), (45), and (46) into (43), we find that

−θγ + (φ + θ)σ̃ ≥ [(.99) − (.12) − (.05)] φσ̄(1 − 2γ̄) ≥ .80φσ̄(1 − 2γ̄),

as required. Lemma 4.3 For σ̃ ∈ [σ̄, 1/2] and φ ∈ [0, φ̄], we have

(θ^{2} + φθ + φ^{2})n^{2}τ̃ ≤ (.07)φσ̄(1 − 2γ̄),

where θ is the smaller root of (8).

Proof. From (34), 0 ≤ θ ≤ 2φ and so

θ^{2} + φθ + φ^{2} ≤ (4 + 2 + 1)φ^{2} ≤ (7φ̄)φ ≤ \frac{(.07)σ̄(1 − 2γ̄)}{n^{2}τ̃}φ.

We can now prove the main result of this section.

Theorem 4.4 If a safe step (with σ_{k} ∈ [σ̄, 1/2]) is taken at iteration k, then the value φ̃ = φ_{k} returned by the line search procedure satisfies

φ_{k} ≥ φ̄.

Proof. Lemmas 4.2 and 4.3 together imply that (41) holds for all φ ∈ [0, φ̄] and, therefore, that the central neighborhood condition (39) holds for these values of φ. Since we have already observed in (30) that φ̄ is smaller than the upper bound (9), it follows that φ_{k} lies somewhere between φ̄ and this upper bound. In particular, φ_{k} ≥ φ̄, as required.

Our linear convergence result follows.

Theorem 4.5 The sequence {μ_{k}} generated by the algorithm converges Q-linearly to zero with a rate constant of (1 − φ̄).

Proof. On safe steps, it follows from Theorem 4.4 that μ_{k+1} ≤ (1 − φ̄)μ_{k}. On fast steps, we have μ_{k+1} ≤ ρμ_{k}. However,

ρ < γ̄ < .50 < 1 − φ̄,

so we choose the more conservative value, namely, 1 − φ̄, as the global rate constant. The polynomiality result is almost immediate.

Theorem 4.6 Suppose that ξ_{x} and ξ_{y} are chosen large enough that for some (x^{∗}, y^{∗}) ∈ S, the scalars ξ_{x}^{∗} and ξ_{y}^{∗} defined by (20) are both O(1). Then, for any given ϵ > 0, the algorithm achieves μ_{K_{ϵ}} ≤ ϵ for

K_{ϵ} = O (n^{2} ln \frac{ξ_{x}ξ_{y}}{ϵ} ) .

Proof. When a safe step is taken at iteration k, it follows from (29) and the fact that τ̃ is independent of n that

μ_{k+1} ≤ (1 − \frac{C_{10}}{n^{2}} ) μ_{k},

where C_{10} = .01σ̄(1 − 2γ̄)/τ̃ is independent of n. When a fast step is taken, we have the even sharper decrease

μ_{k+1} ≤ ρμ_{k}.

Since C_{10}/n^{2} < .01 < 1 − ρ, we find that K_{ϵ} satisfies

ϵ = (1 − \frac{C_{10}}{n^{2}} )^{K_{ϵ}} μ_{0} = (1 − \frac{C_{10}}{n^{2}} )^{K_{ϵ}} ξ_{x}ξ_{y}.

Hence

K_{ϵ} = \frac{ln(ϵ/(ξ_{x}ξ_{y}))}{ln (1 − C_{10}/n^{2})} ≤ \frac{ln(ξ_{x}ξ_{y}/ϵ)}{C_{10}/n^{2}},

which gives the result.

5 Bounds for the Fast Step Components

In this section, we obtain a different set of bounds on the step (u^{k}, v^{k}) calculated in (2). These bounds are used to prove quadratic convergence in the next section.

Because of Assumption 1, we can choose any (x^{∗}, y^{∗}) ∈ S^{c} and define

B = {i | x_{i}^{∗} > 0}, N = {i | y_{i}^{∗} > 0},

where B ∩ N = ∅ and B ∪ N = {1, 2, ⋯, n}. (Note that the definition of B and N is independent of the particular choice of (x^{∗}, y^{∗}).)

For our main result, we take note of the similarity between (u^{k}, v^{k}) and the step (Δx^{k}, Δy^{k}) computed by the algorithm of [4]. Both steps are defined by identical formulae (2), and the relevant properties of the sequence (x^{k}, y^{k}) and the parameters σ̃ and γ̃ that give rise to the steps (Δx^{k}, Δy^{k}) and (u^{k}, v^{k}) are the same; namely

γ̃ ∈ [γ̄, 2γ̄],  σ̃ ∈ (0, 1/2],  (x^{k}, y^{k}) ∈ N (γ_{k}),
μ_{k} ≥ ν_{k}μ_{0}β_{L},  r^{k} = ν_{k}r^{0},  (x^{k}, y^{k}) is a bounded sequence.

(Note that we use β_{L} in place of β̂ in [4] and that, in the algorithm of this paper, we apply a restriction μ_{k} ≤ β_{U}ν_{k}μ_{0} that does not exist in [4].)

We also need several results like those in Section 3 of [4], which define upper and lower bounds for components of (x^{k}, y^{k}) in terms of μ_{k}. These results are given in the first three lemmas of this section. Our main result is Theorem 5.4, which uses the bounds on (x^{k}, y^{k}) to obtain bounds on (u^{k}, v^{k}) using a proof technique identical to one in [4, Section 5]. Finally, in Theorem 5.5, we restate these bounds in the form in which they will be used in Section 6.

We start by defining an “auxiliary sequence” similar to the one defined by Zhang [6, Section 4] and Wright [4, Section 3]. The first element of this sequence is any vector pair (w^{0}, z^{0}) ∈ IR^{n} × IR^{n} such that

z^{0} = Mw^{0} + h,

but not necessarily (w^{0}, z^{0}) ≥ 0. For instance, (w^{0}, z^{0}) could be taken as the least-norm solution of the underdetermined linear system z − Mw = h. Subsequent elements of the auxiliary sequence are defined as

w^{k+1} = w^{k} + φ_{k}(u^{k} + x^{k} − w^{k}) + θ_{k}ū^{k},  k = 0, 1, ⋯,  (47a)
z^{k+1} = z^{k} + φ_{k}(v^{k} + y^{k} − z^{k}) + θ_{k}v̄^{k},  k = 0, 1, · · · .  (47b)

We can state the following result.

Lemma 5.1 For k ≥ 0,

(i) z^{k} = Mw^{k} + h;

(ii) x^{k} − w^{k} = ν_{k}(x^{0} − w^{0}) and y^{k} − z^{k} = ν_{k}(y^{0} − z^{0}).

Proof. We prove (i) by induction. By definition, the result holds for k = 0. Assuming that z^{k} = Mw^{k} + h, we show that z^{k+1} = Mw^{k+1} + h. By (47a),

Mw^{k+1} + h = Mw^{k} + h + φ_{k}[Mu^{k} + Mx^{k} − Mw^{k}] + θMū^{k}.

From (2) and (3), we have

Mu^{k} = v^{k} + y^{k} − Mx^{k} − h, Mū^{k} = v̄^{k},

and so

Mw^{k+1} + h = Mw^{k} + h + φ_{k}[v^{k} + y^{k} − Mx^{k} − h + Mx^{k} − Mw^{k}] + θ_{k}v̄^{k}.

Since z^{k} = Mw^{k} + h, we have

Mw^{k+1} + h = z^{k} + φ_{k}[v^{k} + y^{k} − z^{k}] + θ_{k}v̄^{k} = z^{k+1},

by (47b).

For the first part of (ii), it suffices to show that

(x^{k+1} − w^{k+1}) = (1 − φ_{k})(x^{k} − w^{k}), k = 0, 1, · · · .

Substituting from (47a), and using x^{k+1} = x^{k} + φ_{k}u^{k} + θ_{k}ū^{k}, we have

x^{k+1} − w^{k+1} = (x^{k} + φ_{k}u^{k} + θ_{k}ū^{k}) − [w^{k} + φ_{k}(u^{k} + x^{k} − w^{k}) + θ_{k}ū^{k}]
= x^{k}(1 − φ_{k}) − w^{k}(1 − φ_{k}) = (1 − φ_{k})(x^{k} − w^{k}),

as required. The second part of (ii) is proved similarly.

There is no need to actually calculate the auxiliary sequence. It is used merely as a technical device to obtain the set of bounds in the following lemma.

Lemma 5.2 There is a constant C_{1} > 0 such that for all k ≥ 0, we have

i ∈ N ⇒ x_{i}^{k} ≤ C_{1}μ_{k},  (48a)
i ∈ B ⇒ y_{i}^{k} ≤ C_{1}μ_{k},  (48b)

and

i ∈ B ⇒ x_{i}^{k} ≥ γ̄/C_{1},  (49a)
i ∈ N ⇒ y_{i}^{k} ≥ γ̄/C_{1}.  (49b)

Proof. Let (x^{∗}, y^{∗}) ∈ S^{c}. Then

(x^{k} − x^{∗})^{T} (y^{k} − y^{∗})
= (x^{k} − w^{k} + w^{k} − x^{∗})^{T} (y^{k} − z^{k} + z^{k} − y^{∗})
= (x^{k} − w^{k})^{T} (y^{k} − z^{k}) + (w^{k} − x^{∗})^{T} (y^{k} − z^{k})
+(x^{k} − w^{k})^{T} (z^{k} − y^{∗}) + (w^{k} − x^{∗})^{T} (z^{k} − y^{∗}).

Now z^{k} = Mw^{k} + h and y^{∗} = Mx^{∗} + h, so (w^{k} − x^{∗})^{T} (z^{k} − y^{∗}) ≥ 0 and we have (x^{k} − x^{∗})^{T} (y^{k} − y^{∗}) ≥ (x^{k} − w^{k})^{T} (y^{k} − z^{k}) + (w^{k} − x^{∗})^{T} (y^{k} − z^{k}) + (x^{k} − w^{k})^{T} (z^{k} − y^{∗}). We can now use Lemma 5.1(ii) and (x^{∗})^{T} y^{∗} = 0 to write (x^{∗})^{T} y^{k} + (y^{∗})^{T} x^{k}

≤ (x^{k})^{T} y^{k} − (x^{k} − w^{k})^{T} (y^{k} − z^{k}) − (w^{k} − x^{∗})^{T} (y^{k} − z^{k}) − (x^{k} − w^{k})^{T} (z^{k} − y^{∗})
≤ nμ_{k} + ν_{k}^{2} |(x^{0} − w^{0})^{T} (y^{0} − z^{0})| + ν_{k} |(w^{k} − x^{∗})^{T} (y^{0} − z^{0})| + ν_{k} |(x^{0} − w^{0})^{T} (z^{k} − y^{∗})|
≤ nμ_{k} + ν_{k} [|(x^{0} − w^{0})^{T} (y^{0} − z^{0})| + |(w^{k} − x^{∗})^{T} (y^{0} − z^{0})| + |(x^{0} − w^{0})^{T} (z^{k} − y^{∗})|] .(50)

Now

∥w^{k} − x^{∗}∥ ≤ ∥w^{k} − x^{k}∥ + ∥x^{k} − x^{∗}∥ ≤ ν_{k}∥w^{0} − x^{0}∥ + ∥x^{k}∥ + ∥x^{∗}∥,
∥z^{k} − y^{∗}∥ ≤ ∥z^{k} − y^{k}∥ + ∥y^{k} − y^{∗}∥ ≤ ν_{k}∥z^{0} − y^{0}∥ + ∥y^{k}∥ + ∥y^{∗}∥.

Hence, by boundedness of (x^{k}, y^{k}), we can bound the bracketed term in (50) by a constant independent of k (C̄_{1}, say) and write

(x^{∗})^{T} y^{k} + (y^{∗})^{T} x^{k} ≤ nμ_{k} + C̄_{1}ν_{k}.

By (18), we have ν_{k} ≤ (β_{L}μ_{0})^{−1}μ_{k}, so

(x^{∗})^{T} y^{k} + (y^{∗})^{T} x^{k} ≤ μ_{k} [n + C̄_{1}(β_{L}μ_{0})^{−1}] .  (51)

Since (x^{∗}, y^{∗}) ≥ (0, 0) and (x^{k}, y^{k}) > (0, 0), we have from (51) that

i ∈ N ⇒ y_{i}^{∗}x_{i}^{k} ≤ μ_{k} [n + C̄_{1}(β_{L}μ_{0})^{−1}] ,  (52a)
i ∈ B ⇒ x_{i}^{∗}y_{i}^{k} ≤ μ_{k} [n + C̄_{1}(β_{L}μ_{0})^{−1}] .  (52b)

Since x_{i}^{∗} > 0 for i ∈ B and y_{i}^{∗} > 0 for i ∈ N, we can define C_{1} < ∞ by

C_{1} = [n + C̄_{1}(β_{L}μ_{0})^{−1}] max (sup_{i∈B} \frac{1}{x_{i}^{∗}} , sup_{i∈N} \frac{1}{y_{i}^{∗}} ) .  (53)

The inequalities (48) follow immediately from (52) and (53).

For (49a), we have using γ_{k} ≥ γ̄ that

x_{i}^{k}y_{i}^{k} ≥ γ_{k}μ_{k} ≥ γ̄μ_{k},

so, using (48b), we have

i ∈ B ⇒ x_{i}^{k} ≥ \frac{γ̄μ_{k}}{y_{i}^{k}} ≥ \frac{γ̄μ_{k}}{C_{1}μ_{k}} = \frac{γ̄}{C_{1}} ,

as required. The remaining inequality (49b) is proved similarly.

Note that, in contrast to Zhang [6] and Wright [4], we do not require

(x^{0}, y^{0}) ≥ (w^{0}, z^{0}).  (54)

In [4], we needed (54) to prove (48) and (49) but, in this paper, our a priori knowledge of the boundedness of {(x^{k}, y^{k})} makes this assumption unnecessary.

We now prove a simple result that establishes lower bounds on y_{i}^{k}, i ∈ B, and x_{i}^{k}, i ∈ N.

Lemma 5.3 For all k ≥ 0, there is a positive constant C_{2} such that

i ∈ B ⇒ y_{i}^{k} ≥ \frac{1}{C_{2}} γ̄μ_{k},  (55a)
i ∈ N ⇒ x_{i}^{k} ≥ \frac{1}{C_{2}} γ̄μ_{k}.  (55b)

Proof. Since, by Lemma 3.2, {(x^{k}, y^{k})} is bounded, we can choose C_{2} < ∞ such that, for all k,

0 < x_{i}^{k} ≤ C_{2},  0 < y_{i}^{k} ≤ C_{2},  i = 1, 2, ⋯, n.

For (55a) we have

i ∈ B, x_{i}^{k}y_{i}^{k} ≥ γ_{k}μ_{k} ⇒ y_{i}^{k} ≥ \frac{γ̄μ_{k}}{x_{i}^{k}} ≥ \frac{γ̄μ_{k}}{C_{2}} ,

as required. The proof of (55b) is similar.

We now have all the inequalities needed to prove the main results for boundedness of (u^{k}, v^{k}), using identical techniques to those in [4, Section 5].

Theorem 5.4 There are positive constants C_{4} and C_{5} such that the components of the step (u^{k}, v^{k}) satisfy the bounds

i ∈ N ⇒ |u_{i}^{k}| ≤ C_{4}μ_{k},  (56a)
i ∈ B ⇒ |v_{i}^{k}| ≤ C_{4}μ_{k},  (56b)

and

i ∈ B ⇒ |u_{i}^{k}| ≤ \frac{C_{5}}{2} (μ_{k} + σ̃),  (57a)
i ∈ N ⇒ |v_{i}^{k}| ≤ \frac{C_{5}}{2} (μ_{k} + σ̃).  (57b)

Proof. From (21a), we have

|(\frac{y_{i}^{k}}{x_{i}^{k}} )^{1/2} u_{i}^{k}| ≤ ∥D^{k}u^{k}∥ ≤ \frac{3n}{γ̄^{1/2}} (2 + β_{L}^{−1}(1 + ξ_{x}^{∗} + ξ_{y}^{∗}))μ_{k}^{1/2} = C̄_{4}μ_{k}^{1/2},

for C̄_{4} defined in an obvious way. Hence, using (48a) and x_{i}^{k}y_{i}^{k} ≥ γ̄μ_{k}, we have for i ∈ N that

|u_{i}^{k}| ≤ (\frac{x_{i}^{k}}{y_{i}^{k}} )^{1/2} C̄_{4}μ_{k}^{1/2} = \frac{x_{i}^{k}}{(x_{i}^{k}y_{i}^{k})^{1/2}} C̄_{4}μ_{k}^{1/2} ≤ \frac{x_{i}^{k}}{(γ̄μ_{k})^{1/2}} C̄_{4}μ_{k}^{1/2} ≤ \frac{C_{1}μ_{k}}{γ̄^{1/2}} C̄_{4}.

The inequality (56a) is obtained by setting C_{4} = C_{1}C̄_{4}/γ̄^{1/2}.

The proof of (56b) is similar.

The remaining inequalities (57) can be proved as in Lemmas 5.2 and 5.3 of [4]. Lemma 5.2 is a technical result that is an extension of an earlier result of Ye and Anstreicher [5, Lemma 3.5]. In the proof of [4, Lemma 5.3], we use the inequalities (18), (48), (49), (55), (56), and boundedness of the sequence {(x^{k}, y^{k})}, all of which have been proved above. We omit further details and refer the interested reader to [4].

Theorem 5.5 Suppose we define the constant C_{6} as

C_{6} =^{△} 2 max(1, C_{4}C_{5}).

Then if σ̃ = σ_{k} = μ_{k}, we have

|u_{i}^{k}v_{i}^{k}| ≤ (C_{6}/2)μ_{k}^{2}, i = 1, ⋯, n,

and

|(u^{k})^{T} v^{k}| ≤ n(C_{6}/2)μ_{k}^{2}.

Proof. The result follows immediately from Theorem 5.4.

6 Superlinear Convergence

In this final section, we show that the sequence {μ_{k}} converges superlinearly to zero with Qorder 2. The development follows that in [4, Section 6] quite closely. We write out the proofs of our results where there is enough difference from the proofs of [4] to cause confusion, and omit them otherwise.

We start by defining a threshold condition involving μ_{k} and γ_{k}, and finding bounds on the step length φ̃ given by the fast step procedure when this condition is satisfied.

Lemma 6.1 Suppose at iteration k that

\frac{μ_{k}}{(γ_{k} − γ̄)(1 − γ̄)} ≤ \frac{ρ}{3C_{6}}  (58)

and that a fast step is calculated. Then the step length φ̃ will satisfy

1 ≥ φ̃ ≥ 1 − C_{6} \frac{μ_{k}}{(γ_{k} − γ̄)(1 − γ̄)}.

Proof. Before proceeding, note that if the fast step is successful, the algorithm sets γ_{k+1} to γ̄(1 + γ̄^{t_{k}}), and so

γ_{k} − γ_{k+1} = γ̄(1 + γ̄^{t_{k}−1}) − γ̄(1 + γ̄^{t_{k}}) = γ̄^{t_{k}}(1 − γ̄) = (γ_{k} − γ̄)(1 − γ̄).  (59)

Under these circumstances, condition (58) is equivalent to

\frac{μ_{k}}{γ_{k} − γ_{k+1}} ≤ \frac{ρ}{3C_{6}} .  (60)

We use (58) and (60) interchangeably for the rest of this section.

The remainder of the proof is similar to the proof of Lemma 6.1 in [4]. First, we show that the centrality condition (11) is satisfied for all φ in the interval

φ ∈ [0, 1 − C_{6} \frac{μ_{k}}{γ_{k} − γ_{k+1}} ] .  (61)

Second, we show that the complementarity reduction condition (12) is also satisfied for all φ in the interval (61). The result of the lemma then follows immediately from the way in which φ̃ is chosen by the fast step procedure.

We start with (11). Using the parameter settings γ̃ = γ_{k+1} and σ̃ = σ_{k} = μ_{k}, we have from (2) and (x^{k}, y^{k}) ∈ N (γ_{k}) that

x_{i}^{k}y_{i}^{k} ≥ γ_{k}μ_{k},
x_{i}^{k}v_{i}^{k} + y_{i}^{k}u_{i}^{k} = −x_{i}^{k}y_{i}^{k} + σ̃μ_{k} = −x_{i}^{k}y_{i}^{k} + μ_{k}^{2}.

Hence, using Theorem 5.5, we obtain

(x_{i}^{k} + φu_{i}^{k})(y_{i}^{k} + φv_{i}^{k})
= x_{i}^{k}y_{i}^{k} + φ(x_{i}^{k}v_{i}^{k} + y_{i}^{k}u_{i}^{k}) + φ^{2}u_{i}^{k}v_{i}^{k}
≥ γ_{k}μ_{k}(1 − φ) + φμ_{k}^{2} − φ^{2}(C_{6}/2)μ_{k}^{2}.

On the other hand,

\frac{1}{n}γ_{k+1}(x^{k} + φu_{k})^{T} (y^{k} + φv_{k})

≤ γ_{k+1} [μ_{k}(1 − φ + φμ_{k}) + φ^{2}(C_{6}/2)μ_{k}^{2}] .

Therefore condition (11) will be satisfied provided that

γ_{k}μ_{k}(1 − φ) + φμ_{k}^{2} − φ^{2}(C_{6}/2)μ_{k}^{2}
≥ γ_{k+1} [μ_{k}(1 − φ + φμ_{k}) + φ^{2}(C_{6}/2)μ_{k}^{2}] ,

or, equivalently,

(γ_{k} − γ_{k+1})μ_{k}(1 − φ) + φμ_{k}^{2}(1 − γ_{k+1}) − φ^{2}μ_{k}^{2}(C_{6}/2)(1 + γ_{k+1}) ≥ 0.

Since γ_{k+1} ∈ [γ̄, 1), α ∈ (0, 1], and φ ∈ (0, 1], this last inequality will hold if

(γ_{k} − γ_{k+1})μ_{k}(1 − φ) − φC_{6}μ_{k}^{2} ≥ 0 ⇒ 0 ≤ φ ≤ φ̄_{I} =^{△} \frac{1}{1 + C_{6}μ_{k}/(γ_{k} − γ_{k+1})}.

We can use (26) to obtain

φ̄_{I} ≥ 1 − C_{6} \frac{μ_{k}}{γ_{k} − γ_{k+1}} ,

and so (11) is satisfied for φ in the range (61).

It suffices to show for the second part of the proof that

|\frac{1}{n}(x^{k} + φu^{k})^{T} (y^{k} + φv^{k}) − (1 − φ)μ_{k}| ≤ β̃(1 − φ)μ_{k}  (62)

for all φ satisfying (61). Since

\frac{1}{n}(x^{k} + φu^{k})^{T} (y^{k} + φv^{k}) = μ_{k}(1 − φ + φμ_{k}) + \frac{1}{n}φ^{2}(u^{k})^{T} v^{k},

we have

|\frac{1}{n}(x^{k} + φu^{k})^{T} (y^{k} + φv^{k}) − (1 − φ)μ_{k}| ≤ φμ_{k}^{2} + \frac{1}{n}φ^{2} |(u^{k})^{T} v^{k}| ≤ φμ_{k}^{2} + φ^{2}(C_{6}/2)μ_{k}^{2}.

Therefore, (62) holds if

φμ_{k} + φ^{2}(C_{6}/2)μ_{k} ≤ β̃(1 − φ),

which in turn is true provided that

φ [(1 + C_{6}/2)μ_{k} + β̃] ≤ β̃.  (63)

Now, using C_{6} ≥ 2, (63) is satisfied if

φ ≤ φ̄_{II} =^{△} \frac{β̃}{β̃ + C_{6}μ_{k}} .

Since β̃ = γ̄^{t_{k}} ≥ γ_{k} − γ_{k+1} by (59), we have

φ̄_{II} = \frac{1}{1 + C_{6}μ_{k}/γ̄^{t_{k}}} ≥ \frac{1}{1 + C_{6}μ_{k}/(γ_{k} − γ_{k+1})}.

Using (26) again, we obtain

φ̄_{II} ≥ 1 − \frac{C_{6}μ_{k}}{γ_{k} − γ_{k+1}} ,

which implies that (12) holds for all φ satisfying (61).

We can now show that, when a fast step calculation is performed and the condition (58) holds, the fast step will be accepted by the algorithm, since μ_{k+1} ≤ ρμ_{k}.

Lemma 6.2 If (58) is satisfied at iteration k and a fast step is calculated, then

μ_{k+1} ≤ \frac{3C_{6}}{γ_{k} − γ_{k+1}} μ_{k}^{2} ≤ ρμ_{k},  (64)

and

\frac{μ_{k+1}}{(γ_{k+1} − γ̄)(1 − γ̄)} ≤ \frac{μ_{k}}{(γ_{k} − γ̄)(1 − γ̄)}.  (65)

Proof. We start by showing that μ_{k}(φ) defined by

μ_{k}(φ) =^{△} \frac{1}{n}(x^{k} + φu^{k})^{T} (y^{k} + φv^{k}) = μ_{k}(1 − φ + φμ_{k}) + \frac{1}{n}φ^{2}(u^{k})^{T} v^{k}

is monotonically decreasing for φ ∈ [0, 1] when μ_{k} > 0. Using Theorem 5.5 and the fact that C_{6} ≥ 1, we have

μ′_{k}(φ) = −μ_{k}(1 − μ_{k}) + φ\frac{2}{n}(u^{k})^{T} v^{k} ≤ −μ_{k}(1 − μ_{k}) + φC_{6}μ_{k}^{2} ≤ −μ_{k}(1 − 2C_{6}μ_{k}).

It follows from (58) that 2C_{6}μ_{k} < 1, and so μ′_{k}(φ) < 0 as required.

Since a successful fast step will choose φ_{k} in the interval (61), we have

μ_{k+1} = μ_{k}(1 − φ_{k} + φ_{k}μ_{k}) + \frac{1}{n}φ_{k}^{2}(u^{k})^{T} v^{k}
≤ μ_{k}(1 − φ_{k}) + μ_{k}^{2} + \frac{1}{n}(u^{k})^{T} v^{k}
≤ \frac{C_{6}μ_{k}^{2}}{γ_{k} − γ_{k+1}} + μ_{k}^{2} + \frac{C_{6}}{2} μ_{k}^{2}.

Since 1 ≤ C_{6}/2 ≤ C_{6}/(γ_{k} − γ_{k+1}), we can bound all three terms in the last expression by C_{6}μ_{k}^{2}/(γ_{k} − γ_{k+1}) and obtain

μ_{k+1} ≤ \frac{3C_{6}}{γ_{k} − γ_{k+1}} μ_{k}^{2},

giving the first inequality in (64). The second inequality in (64) follows immediately from (60).

For (65), we can use (64) to write

\frac{μ_{k+1}}{(γ_{k+1} − γ̄)(1 − γ̄)} = \frac{μ_{k+1}}{γ̄^{t_{k}+1}(1 − γ̄)} ≤ \frac{1}{γ̄^{t_{k}+1}(1 − γ̄)} \frac{3C_{6}}{γ̄^{t_{k}}(1 − γ̄)}μ_{k}^{2} = \frac{3C_{6}}{γ̄} [\frac{μ_{k}}{γ̄^{t_{k}}(1 − γ̄)}]^{2} . (66)

From (58) and (59), we have

\frac{3C_{6}}{γ̄} \frac{μ_{k}}{γ̄^{t_{k}}(1 − γ̄)} = \frac{3C_{6}}{γ̄} \frac{μ_{k}}{(γ_{k} − γ̄)(1 − γ̄)} ≤ \frac{ρ}{γ̄} < 1,

where the last inequality follows from the definition of ρ and γ̄ in the algorithm. Substitution of this inequality into (66) yields the desired result.

Theorem 6.3 Suppose that condition (58) is satisfied at iteration K and that

μ_{K} ≤ μ̄.

Then

(i) the algorithm takes fast steps at iteration K and at all subsequent iterations, and

(ii) the sequence {μ_{k}} converges superlinearly to zero with a Q-order of 2.

Proof. The second inequality in (64) guarantees that the fast step is accepted at iteration K, while (65) ensures that the threshold condition (58) still holds at iteration K + 1. By induction, it follows that (i) is true.

For (ii), we apply an argument from [4] to the first inequality in (64). See [4, Theorem 6.3] for the details.

Finally, we show that the algorithm will eventually reach an iterate K at which both (58) and μ_{K} ≤ μ̄ are satisfied, and so superlinear convergence is guaranteed to occur.

Theorem 6.4 Suppose we define constants f and ϵ̂ by

f =^{△} \frac{log γ̄}{log ρ} ∈ (0, 1)

and

ϵ̂ = \{ \begin{matrix}μ̄ & if \frac{μ̄}{γ̄(1−γ̄)} ≤ \frac{ρ}{3C_{6}} , \\ [\frac{ρ}{3C_{6}} \frac{γ̄^{2}(1−γ̄)}{μ̄}]^{1/(1−f)} μ̄ & otherwise. \end{matrix}  (67)

Then if K is the smallest positive integer for which

μ_{K} ≤ ϵ̂,  (68)

we have that

\frac{μ_{K}}{γ̄^{t_{K}}(1 − γ̄)} ≤ \frac{ρ}{3C_{6}}  (69)

and

μ_{K} ≤ μ̄,  (70)

and hence the conditions of Theorem 6.3 are satisfied.

Proof. The proof is almost identical to that of Theorem 6.4 in [4], so we omit it.

This section culminates in the following result, which is immediate from Theorems 6.3 and 6.4.

Corollary 6.5 The sequence {μ_{k}} converges to zero superlinearly with a Q-order of 2.

Acknowledgments

I am grateful to Florian Potra, who motivated this investigation by remarking that there should be a way to drop the assumption that a strictly feasible point exists.

References

[1] J. Ji, F. Potra, and S. Huang, A predictor-corrector method for linear complementarity problems with polynomial complexity and superlinear convergence, Tech. Rep. 18, Department of Mathematics, University of Iowa, Iowa City, Iowa, August 1991.

[2] F. A. Potra, An infeasible interior-point predictor-corrector algorithm for linear programming, Tech. Rep. 26, Department of Mathematics, University of Iowa, Iowa City, Iowa, June 1992.

[3] , A quadratically convergent infeasible interior-point algorithm for linear programming, Tech. Rep. 28, Department of Mathematics, University of Iowa, Iowa City, Iowa, July 1992.

[4] S. J. Wright, An infeasible-interior-point algorithm for linear complementarity problems, Preprint MCS–P331–0991, Mathematics and Computer Science Division, Argonne National Laboratory, Argonne, Illinois, October 1992.

[5] Y. Ye and K. Anstreicher, On quadratic and O(\sqrt{n}L) convergence of a predictorcorrector algorithm for LCP, Tech. Rep. 91–20, Department of Management Sciences, University of Iowa, Iowa City, Iowa, November 1991.

[6] Y. Zhang, On the convergence of an infeasible interior-point algorithm for linear programming and other problems, Tech. Rep. 92–07, Department of Mathematics and Statistics, University of Maryland, Baltimore County, Maryland, April 1992.