M V Solodov
email: solodov@impa.br

B F Svaiter

A TRULY GLOBALLY CONVERGENT NEWTON-TYPE METHOD FOR THE MONOTONE NONLINEAR COMPLEMENTARITY PROBLEM

Keywords: nonlinear complementarity problem, Newton method, proximal point method, projection method, global convergence, superlinear convergence

The Josephy-Newton method for solving a nonlinear complementarity problem consists of solving, possibly inexactly, a sequence of linear complementarity problems. Under appropriate regularity assumptions, this method is known to be locally (superlinearly) convergent. To enlarge the domain of convergence of the Newton method, some globalization strategy based on a chosen merit function is typically used. However, to ensure global convergence to a solution, some additional restrictive assumptions are needed. These assumptions imply boundedness of level sets of the merit function and often even (global) uniqueness of the solution. We present a new globalization strategy for monotone problems which is not based on any merit function. Our linesearch procedure utilizes the regularized Newton direction and the monotonicity structure of the problem to force global convergence by means of a (computationally explicit) projection step which reduces the distance to the solution set of the problem. The resulting algorithm is truly globally convergent in the sense that the subproblems are always solvable, and the whole sequence of iterates converges to a solution of the problem without any regularity assumptions. In fact, the solution set can even be unbounded. Each iteration of the new method has the same order of computational cost as an iteration of the damped Newton method. Under natural assumptions, the local superlinear rate of convergence is also achieved.

Introduction.

The classical nonlinear complementarity problem [START_REF] Pang | Complementarity problems, in Handbook of Global Optimization[END_REF][START_REF]Complementarity and Variational Problems: State of the Art[END_REF], NCP(F ), is to find a point x ∈ n such that x ≥ 0, F(x) ≥ 0,

x, F (x) = 0, (1.1) where F : n → n and •, • denotes the usual inner product in n . Throughout this paper, we shall assume that F (•) is continuous and monotone, i.e., F (x) -F (y), xy ≥ 0 for all x, y ∈ n .

Note that under this assumption, the solution set of (1.1) is convex.

While there exists a wide range of approaches to solving NCP(F ) (see [START_REF] Harker | Finite-dimensional variational inequality problems: A survey of theory, algorithms and applications[END_REF][START_REF] Pang | Complementarity problems, in Handbook of Global Optimization[END_REF][START_REF]Complementarity and Variational Problems: State of the Art[END_REF]), some of the most successful and widely used are Newton-type algorithms based on solving successive linearizations of the problem (see, for example, the more detailed discussion in [START_REF] Ferris | Feasible descent algorithms for mixed complementarity problems[END_REF]). Given a point x k , the (Josephy-) Newton method [START_REF] Josephy | Newton's Method for Generalized Equations[END_REF][START_REF] Pang | Iterative methods for variational and complementarity problems[END_REF][START_REF] Pang | Inexact Newton methods for the nonlinear complementarity problem[END_REF][START_REF] Harker | Finite-dimensional variational inequality problems: A survey of theory, algorithms and applications[END_REF][START_REF] Bonnans | Local analysis of Newton-type methods for variational inequalities and nonlinear programming[END_REF] generates the next iterate x k+1 by solving the linear complementarity problem x ≥ 0, F k (x) ≥ 0,

x, F k (x) = 0, (1.2) where F k (•) is the first-order approximation of F (•) at x k , F k (x) := F (x k ) + ∇F (x k )(xx k ), (1.3) assuming F (•) is differentiable. If the starting point is sufficiently close to some regular [START_REF] Robinson | Strongly regular generalized equations[END_REF] solution x of (1.1), the sequence generated by the Newton method is well defined and converges to x superlinearly or quadratically, depending on further assumptions. We note, in passing, that monotonicity of F is not needed for such local analysis.

There are two key difficulties with using the Newton method given by (1.2), (1.3) for solving NCP(F ). First, in the absence of regularity even local convergence cannot be ensured. Second, even if regularity holds at a solution, there are serious problems with ensuring global convergence. In particular, far from a regular solution of the problem, there is no guarantee that the linearization subproblems are solvable. And even if the subproblem solution exists, there is no guarantee that it actually constitutes some progress toward solving the original problem, NCP(F ). In this paper, we shall successfully address each of these difficulties in the context of monotone problems.

To enlarge the domain of convergence of the Newton method, some globalization strategy has to be used. However, as pointed out in [14, p. 182], "the global convergence remains a rare property for most of the modified methods." As mentioned in [14, p. 185], "the trouble with general variational inequalities and nonlinear complementarity problems is that valid merit functions which are relatively easy to compute are very difficult if not impossible to find." Although considerable progress has been made in the theory and numerical use of merit functions in recent years [START_REF] Fukushima | Equivalent differentiable optimization problems and descent methods for asymmetric variational inequality problems[END_REF][START_REF] Mangasarian | Nonlinear complementarity as unconstrained and constrained minimization[END_REF][START_REF] Fischer | An NCP function and its use for the solution of complementarity problems[END_REF][START_REF] Kanzow | Nonlinear complementarity as unconstrained optimization[END_REF][START_REF] Facchinei | A new merit function for nonlinear complementarity problems and a related algorithm[END_REF][START_REF] Solodov | Stationary points of bound constrained reformulations of complementarity problems[END_REF][START_REF] Luo | A new class of merit functions for the nonlinear complementarity problem[END_REF][START_REF] Mangasarian | A linearly convergent derivative-free descent method for strongly monotone complementarity problems[END_REF][START_REF] Solodov | Implicit Lagrangian[END_REF][START_REF] Solodov | Some optimization reformulations of the extended linear complementarity problem[END_REF] (see [START_REF] Fukushima | Merit functions for variational inequality and complementarity problems[END_REF] for a survey of merit functions for variational inequality and complementarity problems), most of the known merit functions do not appear to be useful for the specific task of globalizing the Newton method given by (1.2), (1.3). At this point, it is worth emphasizing that the Newton method under consideration should not be confused with other Newton-like methods which are not of the Josephy type-for example, Newton methods for minimizing a particular merit function or Newton methods for equation-based NCP reformulations [START_REF] Kanzow | Some equation-based methods for the nonlinear complementarity problem[END_REF][START_REF] Ralph | Global convergence of damped Newton's method for nonsmooth equations via the path search[END_REF][START_REF] Luca | A semismooth equation approach to the solution of nonlinear complementarity problems[END_REF][START_REF] Facchinei | Inexact Newton methods for semismooth equations with applications to variational inequality problems[END_REF][START_REF] Facchinei | A new merit function for nonlinear complementarity problems and a related algorithm[END_REF][START_REF] Jiang | A new nonsmooth equations approach to nonlinear complementarity problems[END_REF][START_REF] Sun | A computable generalized Hessian of the D-gap function and Newton-type methods for variational inequality problems[END_REF][START_REF] Yamashita | Modified Newton methods for solving a semismooth reformulation of monotone complementarity problems[END_REF][START_REF] Jiang | Global and local superlinear convergence analysis of Newtontype methods for semismooth equations with smooth least squares, in Reformulation-Nonsmooth, Piecewise Smooth, Semismooth and Smoothing Methods[END_REF][START_REF] Kanzow | Theoretical and numerical investigation of the D-gap function for box constrained variational inequalities[END_REF][START_REF] Han | Newton-Type Methods for Variational Inequalities[END_REF][START_REF] Sun | A regularization Newton method for solving nonlinear complementarity problems[END_REF].

Perhaps the most natural globalization strategy for (1.2), (1.3) is a linesearch procedure in the obtained Newton direction (if it exists!) aimed at decreasing the value of some valid merit function. In this regard, we note the following. First, many of the known NCP merit functions are not differentiable, e.g., those based on the natural residual [START_REF] Pang | Inexact Newton methods for the nonlinear complementarity problem[END_REF] and the normal map [START_REF] Robinson | Normal maps induced by linear transformations[END_REF]. This makes linesearch difficult, although an alternative pathsearch approach is possible [START_REF] Ralph | Global convergence of damped Newton's method for nonsmooth equations via the path search[END_REF][START_REF] Dirkse | The PATH solver: A nonmonotone stabilization scheme for mixed complementarity problems[END_REF]. Second, for some differentiable merit functions, e.g., (the square of) the Fischer-Burmeister function [START_REF] Fischer | A special Newton-type optimization method[END_REF][START_REF] Kanzow | Nonlinear complementarity as unconstrained optimization[END_REF][START_REF] Facchinei | A new merit function for nonlinear complementarity problems and a related algorithm[END_REF], it appears that the Newton direction need not be a direction of descent, particularly far from a solution or in the absence of regularity. Thus the use of this function results in a rather indirect globalization of the Newton method [START_REF] Ferris | Feasible descent algorithms for mixed complementarity problems[END_REF]. (More on this later.)

To our knowledge, the only merit functions which have been used to globalize the Newton method of the form (1.2), (1.3) are the gap function [START_REF] Marcotte | A note on a globally convergent Newton method for solving monotone variational inequalities[END_REF], the regularized gap function [START_REF] Taji | A globally convergent Newton method for solving strongly monotone variational inequalities[END_REF][START_REF] Taji | Optimization based globally convergent methods for the nonlinear complementarity problems[END_REF], the D-gap function [START_REF] Peng | A hybrid Newton method for solving the variational inequality problem via the D-gap function[END_REF][START_REF] Peng | A hybrid Josephy-Newton method for solving box constrained variational inequality problem via the D-gap function[END_REF], and (the square of) the Fischer-Burmeister function [START_REF] Ferris | Feasible descent algorithms for mixed complementarity problems[END_REF]. (We note that some of these methods were developed for the more general mixed complementarity or variational inequality setting.) However, each of these globalizations has certain drawbacks. Using the gap function [START_REF] Marcotte | A note on a globally convergent Newton method for solving monotone variational inequalities[END_REF] requires exact minimization along the line to compute the stepsize, which is not implementable. Moreover, the gap function itself is not easy to compute in general.

In addition, compactness of the feasible set of the variational inequality problem is required, which precludes an application to complementarity problems. Using the regularized gap function [START_REF] Taji | A globally convergent Newton method for solving strongly monotone variational inequalities[END_REF][START_REF] Taji | Optimization based globally convergent methods for the nonlinear complementarity problems[END_REF] admits inexact Armijo-type linesearch but requires strong monotonicity of F for global convergence (see also [START_REF] Wu | A general descent framework for the monotone variational inequality problem[END_REF]). In addition, methods of [START_REF] Marcotte | A note on a globally convergent Newton method for solving monotone variational inequalities[END_REF][START_REF] Taji | A globally convergent Newton method for solving strongly monotone variational inequalities[END_REF][START_REF] Taji | Optimization based globally convergent methods for the nonlinear complementarity problems[END_REF] also need the (restrictive) strict complementarity assumption to establish superlinear/quadratic rate of convergence. Note that the subproblems in [START_REF] Marcotte | A note on a globally convergent Newton method for solving monotone variational inequalities[END_REF][START_REF] Taji | A globally convergent Newton method for solving strongly monotone variational inequalities[END_REF][START_REF] Taji | Optimization based globally convergent methods for the nonlinear complementarity problems[END_REF] are solvable due to the compactness of the feasible set and the strong monotonicity of F , respectively.

The methods based on the D-gap function [START_REF] Peng | A hybrid Newton method for solving the variational inequality problem via the D-gap function[END_REF][START_REF] Peng | A hybrid Josephy-Newton method for solving box constrained variational inequality problem via the D-gap function[END_REF] and the Fischer-Burmeister function [START_REF] Ferris | Feasible descent algorithms for mixed complementarity problems[END_REF] globalize the Newton method in a rather indirect way. In [START_REF] Peng | A hybrid Newton method for solving the variational inequality problem via the D-gap function[END_REF][START_REF] Peng | A hybrid Josephy-Newton method for solving box constrained variational inequality problem via the D-gap function[END_REF][START_REF] Ferris | Feasible descent algorithms for mixed complementarity problems[END_REF] there is no guarantee that the subproblems are solvable. If the Newton direction does not exist (which cannot be checked a priori, so some wasteful computing will inevitably be done), the method resorts to an "escape" mechanism of taking a gradient descent step for the merit function. Even if the Newton direction exists, it still may happen that it does not satisfy conditions needed to obtain sufficient descent of the merit function by means of a linesearch procedure. In this case, the Newton point will be discarded altogether, and again a gradient descent step will be taken. In a sense, this is a rather indirect globalization of the Newton method, because such a strategy does not correspond to the damped Newton methodology. As for the convergence results in [START_REF] Peng | A hybrid Newton method for solving the variational inequality problem via the D-gap function[END_REF][START_REF] Peng | A hybrid Josephy-Newton method for solving box constrained variational inequality problem via the D-gap function[END_REF][START_REF] Ferris | Feasible descent algorithms for mixed complementarity problems[END_REF], every accumulation point of the generated sequence of iterates is a stationary point of the employed merit function (even without the monotonicity assumption on F ). However, the existence of such accumulation points, and the equivalence of stationary points of the merit functions to solutions of NCP(F ), cannot be guaranteed without further assumptions. For example, in [START_REF] Peng | A hybrid Newton method for solving the variational inequality problem via the D-gap function[END_REF][START_REF] Peng | A hybrid Josephy-Newton method for solving box constrained variational inequality problem via the D-gap function[END_REF], F is assumed to be a uniform P -function, which implies that NCP(F ) has a (globally) unique solution.

Finally, we briefly comment on the interesting regularization approach proposed in [START_REF] Facchinei | Beyond monotonicity in regularization methods for nonlinear complementarity problems[END_REF], which converges globally when the solution set of the NCP is compact and F is a P 0 -function. First, it is important to note that the subproblems in the method of [START_REF] Facchinei | Beyond monotonicity in regularization methods for nonlinear complementarity problems[END_REF] are nonlinear complementarity problems, which are structurally just as difficult to solve as the original problem itself (although they are better behaved due to regularization). Therefore, this method is not of the Newton type. Second, the global convergence result of [START_REF] Facchinei | Beyond monotonicity in regularization methods for nonlinear complementarity problems[END_REF] states that the sequence of iterates remains bounded, and its every accumulation point solves the NCP. This is weaker than the full convergence of the whole sequence which we shall establish for our method. On the other hand, the P 0 assumption on F is, of course, weaker than our assumption of monotonicity. It should also be mentioned that it is often unknown whether the solution set is compact. When F is monotone, this holds if the NCP is strictly feasible, i.e., there exists an x such that x > 0 and F (x) > 0. However, this condition is not easy to verify in general. In any case, the key conceptual difference between our method and that of [START_REF] Facchinei | Beyond monotonicity in regularization methods for nonlinear complementarity problems[END_REF] is in the linear versus nonlinear structure of subproblems solved at each iteration. It is worth noting that the latter is also the important difference between our method and the proximal point algorithm [START_REF] Rockafellar | Monotone operators and the proximal point algorithm[END_REF], which does converge globally under the monotonicity assumption only.

We emphasize that for each of the cited Josephy-Newton algorithms, to ensure global convergence of the whole sequence of iterates to a solution of the problem, one needs assumptions which, among other things, imply that the solution is unique. In fact, relatively restrictive assumptions are required even to prove boundedness of the iterates and convergence to zero of the distance to the solution set.

In this paper, we present a new globalization strategy for monotone problems which overcomes the above mentioned drawbacks. In particular, our algorithm is truly globally convergent in the sense that from any starting point the whole sequence of iterates converges to a solution of NCP(F ) under no assumptions other than monotonicity and continuity of F (and of course, the existence of a solution); see Theorem 3.2. In particular, no regularity-type assumptions are needed. In fact, the solution set may even be unbounded. In addition, the linear complementarity subproblems are always solvable, and our algorithm allows for their inexact solution similar to the setting of [START_REF] Pang | Inexact Newton methods for the nonlinear complementarity problem[END_REF]. This feature of approximate subproblem solutions is of particular importance for large-scale problems. (We note that the effect of an inexact subproblem solution has not been analyzed in globalizations proposed in [START_REF] Taji | A globally convergent Newton method for solving strongly monotone variational inequalities[END_REF][START_REF] Taji | Optimization based globally convergent methods for the nonlinear complementarity problems[END_REF][START_REF] Peng | A hybrid Newton method for solving the variational inequality problem via the D-gap function[END_REF][START_REF] Peng | A hybrid Josephy-Newton method for solving box constrained variational inequality problem via the D-gap function[END_REF][START_REF] Ferris | Feasible descent algorithms for mixed complementarity problems[END_REF].) Under the assumptions of positive definiteness of the Jacobian of F at the solution, and its Hölder continuity in some neighborhood of it, the local superlinear rate of convergence is also established (see Theorem 4.3). Each iteration of our algorithm (see Algorithm 2.2) consists of an approximate solution of a (strongly monotone) linear complementarity problem, followed by a linesearch procedure and a (computationally trivial) projection step. Thus computational cost of each iteration is of the same order as that of the damped Newton method.

2. The algorithm. We start with some equivalent formulations of NCP(F ), each of which will be useful in the subsequent analysis. In particular, the following five statements are equivalent:

1. x solves NCP(F ).

2. x is a solution of the variational inequality problem over the nonnegative orthant n + :

x

∈ n + , F (x), x -x ≥ 0 for all x ∈ n + .
3. x yields a zero of the (maximal monotone) operator F + N :

0 ∈ F (x) + N (x),
where N (x) is the normal cone to n + at the point x. 4. x is a zero of the natural (projection) residual:

0 = r(x) := min{x; F (x)} = x -[x -F (x)] + ,
where the minimum is taken componentwise and [•] + stands for the orthogonal projection map onto n + . The approach presented in this paper is in some ways motivated by the hybrid projection-proximal point method of [START_REF] Solodov | A hybrid projection-proximal point algorithm[END_REF] and the projection method of [START_REF] Solodov | A new projection method for variational inequality problems[END_REF], which already proved to be useful for developing globally (and locally superlinearly) convergent Newton methods for systems of monotone equations [START_REF] Solodov | A globally convergent inexact Newton method for systems of monotone equations[END_REF]. In fact, in a sufficiently small neighborhood of a regular solution of NCP(F ) our Newton-type Algorithm 2.2 takes steps which can be viewed, in a certain sense, as iterations of the hybrid projectionproximal point method. This will prove to be the key to the local superlinear rate of convergence. We therefore first state the algorithm of [START_REF] Solodov | A hybrid projection-proximal point algorithm[END_REF] in the more general context of finding zeros of set-valued maximal monotone operators in a Hilbert space. Let T be a maximal monotone operator on a real Hilbert space H. And consider, for a moment, the problem of finding an x ∈ H such that 0 ∈ T (x). Note that NCP(F ) considered here is a particular instance of this problem with T (x) = (F + N )(x) and H = n . Algorithm 2.1 (hybrid projection-proximal point method [START_REF] Solodov | A hybrid projection-proximal point algorithm[END_REF]). Choose any

x 0 ∈ H and σ ∈ [0, 1); set k := 0. Inexact Proximal Step. Choose μ k > 0 and find y k ∈ H and v k ∈ T (y k ) such that 0 = v k + μ k (y k -x k ) + ε k , (2.1)
where

ε k ≤ σ max{ v k , μ k y k -x k }. (2.2) Stop if v k = 0 or y k = x k . Otherwise, Projection Step. Compute x k+1 = x k - v k , x k -y k v k 2 v k .
Set k := k + 1; and repeat. If problem 0 ∈ T (x) has a solution and the sequence {μ k } is bounded above, then the generated sequence {x k } either is finite and terminates at a solution or is infinite and converges (weakly) to a solution. The linear/superlinear rate of convergence is also achieved under standard assumptions. For complete properties of the method, see [START_REF] Solodov | A hybrid projection-proximal point algorithm[END_REF]. The idea of Algorithm 2.1 is to use an approximate proximal iteration to construct a hyperplane

H k := {x | v k , x -y k = 0}, (2.3)
which separates the current iterate x k from the solutions of 0 ∈ T (x). Let us make this more precise. If 0 ∈ T (x), then, by monotonicity of T , v k , xy k ≤ 0 for any y k and any v k ∈ T (y k ). Here we will consider a condition somewhat stronger than (2.2), namely,

ε k ≤ σμ k y k -x k , (2.4)
because this will be the only condition used in our Algorithm 2.2. Using (2.1) and the Cauchy-Schwarz inequality, we obtain

v k , x k -y k = μ k x k -y k 2 -ε k , x k -y k ≥ μ k x k -y k 2 -ε k x k -y k ≥ μ k (1 -σ) x k -y k 2 > 0, (2.5)
where the last inequality follows from (2.4). Thus whenever (2.4) holds, we have (2.5), and so H k given by (2.3) indeed separates x k from zeros of T . The last step of Algorithm 2.1 is equivalent to projecting x k onto this hyperplane. Separation arguments show that the distance to the solution set for a sequence thus constructed monotonically decreases, which essentially ensures global convergence of the algorithm. Algorithm 2.1 has certain advantages over the classical proximal point method [START_REF] Rockafellar | Monotone operators and the proximal point algorithm[END_REF] in the sense of less restrictive and more constructive tolerance requirements imposed on approximate solutions of proximal subproblems. See [START_REF] Solodov | A hybrid projection-proximal point algorithm[END_REF] for a more detailed discussion. Other related works are [START_REF] Solodov | Forcing strong convergence of proximal point iterations in a Hilbert space[END_REF][START_REF] Solodov | A hybrid approximate extragradient-proximal point algorithm using the enlargement of a maximal monotone operator[END_REF][START_REF] Solodov | A comparison of rates of convergence of two inexact proximal point algorithms[END_REF][START_REF] Solodov | Error bounds for proximal point subproblems and associated inexact proximal point algorithms[END_REF].

Of course, in the context of this paper straightforward application of Algorithm 2.1 to solving NCP(F ) is not practical, since this would involve solving a sequence of nonlinear subproblems which are in general just as difficult as the original NCP, even if they are better behaved due to the regularization. In this sense, the situation is similar to the method of [START_REF] Facchinei | Beyond monotonicity in regularization methods for nonlinear complementarity problems[END_REF] discussed in the introduction. However, the regularization and projection methodology of Algorithm 2.1 would prove to be crucial for devising a globally convergent method.

Given a current iterate x k and a regularization parameter μ k > 0, consider the regularized linear complementarity problem LCP(ϕ k )

x ≥ 0, ϕ k (x) ≥ 0, x, ϕ k (x) = 0, (2.6)
where

ϕ k (x) := F (x k ) + G k (x -x k ) + μ k (x -x k ) (2.7)
with G k being a positive semidefinite matrix (presumably, the Jacobian of F or its approximation, if F is differentiable at x k ). Suppose z k ≥ 0 is some approximate solution of this problem with e k being the associated natural residual [START_REF] Cottle | The Linear Complementarity Problem[END_REF][START_REF] Pang | Inexact Newton methods for the nonlinear complementarity problem[END_REF]:

min{z k ; ϕ k (z k )} = e k , (2.8)
where the minimum is taken componentwise (see (2.9), (2.10) in Algorithm 2.2 for conditions imposed on the error tolerance e k ). Note that since the matrix G k + μ k I is positive definite, LCP(ϕ k ) always has (unique) solution ẑk [START_REF] Cottle | The Linear Complementarity Problem[END_REF] for which the residual is zero:

0 = min{ẑ k ; ϕ k (ẑ k )}.
Hence, LCP(ϕ k ) always has inexact solutions such that e k ≤ δ k , whichever δ k ≥ 0 we choose. Therefore this inexact Newton step is well defined.

The next step of our algorithm is checking whether the inexact Newton point obtained by solving LCP(ϕ k ) provides an acceptable (in the sense of Algorithm 2.1) approximate solution for the proximal point subproblem

0 ∈ (F + N )(x) + μ k (x -x k ).
If this is the case, the previously described separation property holds, and our algorithm proceeds to make the projection step prescribed by Algorithm 2.1.

To make this more precise, we provide the following considerations. Since

e k = min{z k ; ϕ k (z k )} = z k -[z k -ϕ k (z k )] + ,
we have that

z k -e k = [z k -ϕ k (z k )] + ∈ n + .
By properties of the projection operator [33, p. 121],

z k -ϕ k (z k ) -(z k -e k ), x -(z k -e k ) ≤ 0 for all x ∈ n + . Therefore -ϕ k (z k ) + e k , x -(z k -e k ) ≤ 0 for all x ∈ n + ,
which implies that

h k := -ϕ k (z k ) + e k ∈ N (z k -e k ),
so that we have available an element in the normal cone N at the point z ke k ∈ n + . Consider now the pair

y k = z k -e k and v k = F (z k -e k ) -ϕ k (z k ) + e k ∈ (F + N )(z k -e k ) = (F + N )(y k ).
Let us analyze y k , v k as an approximate solution of the proximal subproblem (2.1) in Algorithm 2.1 with T = F + N . We have to check whether (2.4), the stronger version of condition (2.2), is satisfied with

ε k = -v k -μ k (y k -x k )
and y k , v k given above. If (2.4) is satisfied, then (2.5) holds and the hyperplane H k given by (2.3) separates the current iterate x k from zeros of T = F + N or, equivalently, from solutions of NCP(F ). Therefore we can make progress toward a solution of NCP(F ) by making the projection step of Algorithm 2.1 with y k and v k defined above (followed by projection onto n + to preserve feasibility). As we shall see, the test (2.4) and the resulting step will be crucial for obtaining fast local convergence when the regularity assumption holds.

However, far from the solution or if regularity does not hold, an approximate (and even exact) Newton point obtained from solving the linear model may not satisfy the tolerance requirements (2.4) for the nonlinear proximal subproblem. In this case, the preceding separation arguments are not valid, and the Newton point has to be refined. For this task, we employ a linesearch procedure in the approximate Newton direction z kx k (see (2.12) in Algorithm 2.2) which computes a point

y k = x k + α k (z k -x k ) such that 0 < F (y k ), x k -y k .
A similar linesearch technique was used in [START_REF] Solodov | A new projection method for variational inequality problems[END_REF][START_REF] Solodov | A globally convergent inexact Newton method for systems of monotone equations[END_REF]. Note that for any x which solves NCP(F ) we have

0 ≥ F (x), x -y k , where y k ∈ n + . Hence, by monotonicity of F , 0 ≥ F (y k ), x -y k .
Therefore in this case another hyperplane, namely,

H k := {x ∈ n | F (y k ), x -y k = 0},
strictly separates the current iterate x k from solutions of the problem. Once the separating hyperplane H k is obtained, the next iterate x k+1 is computed by projecting x k onto H k and then onto the nonnegative orthant n + . We now formally state the algorithm.

Algorithm 2.2. Choose any x 0 ∈ n , σ ∈ (0, 1), β ∈ (0, 1), and λ ∈ (0, 1). Set k := 0.

Inexact Newton

Step. Choose a positive semidefinite matrix G k and μ k > 0. Choose ρ k ∈ [0, 1) and compute z k ∈ n + , an inexact solution of LCP(ϕ k ) given by (2.6)-(2.8), such that

e k ≤ ρ k μ k z k -x k (2.9) and e k , ϕ k (z k ) + z k -x k ≤ ρ k μ k z k -x k 2 .
(2.10)

Stop if z k = x k . Otherwise, Linesearch Step. Set v k := F (z k -e k ) -ϕ k (z k ) + e k and y k := z k -e k .
Let

ε k = -v k -μ k (y k -x k ). If ε k ≤ σμ k y k -x k , (2.11)
then go to the Projection Step.

Otherwise, find y k = x k + α k (z kx k ), where α k = β m k with m k being the smallest nonnegative integer m such that

F (x k + β m (z k -x k )), x k -z k ≥ λ(1 -ρ k )μ k z k -x k 2 .
(2.12) Set v k := F (y k ) and go to the Projection Step.

Projection

Step. Compute

x k+1 = max 0; x k - v k , x k -y k v k 2 v k . (2.

13)

Set k := k + 1; and repeat.

To compute z k , an approximate solution of LCP(ϕ k ) satisfying (2.6)-(2.10), one can employ any appropriate algorithm known to converge for the strongly monotone linear complementarity problem. There are many algorithms which would generate a sequence converging to the unique solution of LCP(ϕ k ) with a (global) quadratic rate. This guarantees that after finitely many iterations (just a few, one hopes), the LCP(ϕ k ) residual e k would be small enough, so that (2.9)-(2.10) holds.

Note that Algorithm 2.2 has computational cost per iteration of the same order as any other damped Newton method: solving, possibly inexactly, a linear complementarity problem followed by a simple linesearch procedure (if (2.11) is not satisfied). The projection step is explicit, and therefore its computational cost is negligible. The advantage of Algorithm 2.2 over other Newton-type methods is that it is truly globally convergent under minimal assumptions.

Global convergence analysis.

We start with the global convergence analysis. Throughout we assume that the solution set of the problem is nonempty. First note that if Algorithm 2.2 terminates with z k = x k , then e k = 0 by (2.9), and we have that z k is the exact solution of LCP(ϕ k ). Therefore ϕ k (z k ), xz k ≥ 0 for all x ∈ n + . Because z k = x k implies that ϕ k (z k ) = F (z k ), it follows from the latter inequality that z k solves NCP(F ). From now on, we assume that z k = x k for all k, and an infinite sequence {x k } is generated.

We first state a preliminary result [START_REF] Solodov | A hybrid projection-proximal point algorithm[END_REF] whose simple proof we include for completeness.

Lemma 3.1. Let x, y, v, x be any elements of n such that v, xy > 0 and v, xy ≤ 0.

Let x = x - v, x -y v 2 v. Then x -x 2 ≤ x -x 2 -x -x 2 .
Proof. It follows from the hypothesis that the hyperplane H = {s | v, sy = 0} separates x from x. Moreover, x is the projection of x onto the half-space {s | v, sy ≤ 0}. Since x belongs to this half-space, it follows from the basic properties of the projection operator (see [33, p. 121

]) that x -x, x -x ≥ 0. Therefore x -x 2 = x -x 2 + x -x 2 + 2 x -x, x -x ≥ x -x 2 + x -x 2 ,
which completes the proof.

We are now ready to prove our main global convergence result. The remarkable property established in Theorem 3.2 is that the whole sequence of iterates always converges to some solution of NCP(F ) under the assumptions of merely monotonicity and continuity of F . NCP(F ) need not have a unique solution; in fact, the solution set may even be unbounded. In the latter case the level sets of any merit function for NCP(F ) are also unbounded, so that for typical Newton methods with globalization strategies based on merit functions, even boundedness of iterates cannot be established.

Theorem 3.2. Suppose that F is continuous and monotone. Then any sequence {x k } generated by Algorithm 2.2 is bounded.

Suppose further that there exist constants C 1 , C 2 , C 3 > 0 and t > 0 such that G k ≤ C 1 for all k and, starting with some index k 0 ,

C 2 ≥ μ k ≥ C 3 r(x k ) t . Suppose that min{1; 1/C 2 } > lim sup k→∞ ρ k .
Then {x k } converges to some x, which is a solution of NCP(F ).

Proof. First note that because

z k -e k = [z k -ϕ k (z k )] + and x k ∈ n
+ , by properties of the projection operator [33, p. 121] it follows that

z k -ϕ k (z k ) -(z k -e k ), x k -(z k -e k ) ≤ 0. Therefore -ϕ k (z k ) + e k , x k -z k + e k ≤ 0.
Making use of the latter inequality, we further obtain

F (x k ), x k -z k ≥ F (x k ) -ϕ k (z k ) + e k , x k -z k + e k -F (x k ), e k = (G k + μ k I)(x k -z k ), x k -z k + e k + e k , x k -z k -F (x k ) + e k 2 ≥ μ k z k -x k 2 -e k , (G k + μ k I)(z k -x k ) + F (x k ) + z k -x k = μ k z k -x k 2 -e k , ϕ k (z k ) + z k -x k ≥ μ k (1 -ρ k ) z k -x k 2 , (3.1)
where the last inequality follows from (2.10) and the next to last follows from positive semidefiniteness of G k .

We next show that the linesearch procedure (2.12), if activated, always terminates with a positive stepsize α k . Suppose that for some iteration index k this is not the case. That is, for all integers m we have

F (x k + β m (z k -x k )), x k -z k < λ(1 -ρ k )μ k z k -x k 2 .
Since F is continuous, passing onto the limit as m → ∞, we obtain

F (x k ), x k -z k ≤ λ(1 -ρ k )μ k z k -x k 2 . (3.2)
Now since λ ∈ (0, 1), ρ k ∈ [0, 1), and z k = x k , (3.2) contradicts (3.1). Therefore the linesearch step is well defined. Denote

xk := x k - v k , x k -y k v k 2 v k . (3.3)
Observe that this is well defined because v k = 0. Indeed, if v k = 0 in the case when the linesearch is not activated, then (2.11) must hold. Furthermore, in that case v k = 0 implies ε k = -μ k (y kx k ), which together with (2.11) implies that ε k = y k -x k = 0. Therefore, using the definition of y k , we conclude that

x k = z k -e k . If μ k ≤ 1, then ρ k μ k < 1, so (2.9) implies that z k = x k , in contradiction with the stopping test. Suppose now that μ k > 1. Note that 0 = v k = F (x k ) -ϕ k (z k ) + e k = (-G k +(1-μ k )I)e k . It then follows that 0 = (-G k +(1-μ k )I)e k , e k ≤ (1-μ k ) e k 2 .
Since 1 < μ k , we again have that e k = 0, i.e., z k = x k , which contradicts the stopping rule. If the linesearch is used to compute y k , then v k = F (y k ), which again cannot be zero by (2.12). The observation that v k = 0 concludes the proof that the whole algorithm is well defined.

With the notation (3.3), we have that x k+1 = [x k ] + . Let x be any solution of NCP(F ). Since x ∈ n + , it is easy to see that (2.12). In this respect, the only difference is the choice of y k and v k . As discussed before, in either case v k , xy k ≤ 0, so that by Lemma 3.1, it follows that xk -

x k+1 -x ≤ xk -x . (3.4) Note that if (2.11) is satisfied, then v k , x k -y k > 0 (see (2.5)). If (2.11) is not satisfied, then v k , x k -y k = α k F (y k ), x k -z k > 0 still holds by
x 2 ≤ x k -x 2 -xk -x k 2 .
Combining the latter relation with (3.4), we obtain

x k+1 -x 2 ≤ x k -x 2 -xk -x k 2 . (3.5)
It immediately follows that the sequence { x kx } is monotone, so it converges. Therefore, {x k } is bounded.

We next consider the two possible cases:

0 = lim inf k→∞ r(x k ) (3.6) and 0 < lim inf k→∞ r(x k ) . (3.7)
In the first case, by continuity of r(•) and boundedness of {x k }, there exists x, an accumulation point of {x k }, such that r(x) = 0. Therefore x is a solution of NCP(F ). Since x was an arbitrary solution, we can now choose x = x in (3.5). Because the sequence { x kx } converges and x is an accumulation point of {x k }, it must be the case that {x k } converges to x, which is a solution of NCP(F ).

Consider now the second case. By (3.5), it follows that

0 = lim k→∞ xk -x k or, equivalently, 0 = lim k→∞ v k , x k -y k v k . (3.8)
By the choice of μ k and (3.7), it then follows that μ k ≥ C 3 r(x k ) t ≥ C 4 > 0 for all k. By (3.1) and the Cauchy-Schwarz inequality, we obtain

F (x k ) x k -z k ≥ F (x k ), x k -z k ≥ μ k (1 -ρ k ) x k -z k 2 .
Hence,

F (x k ) ≥ C 4 (1 -ρ k ) x k -z k .
Taking into account boundedness of {x k } and continuity of F , and the fact that 1 > lim sup k→∞ ρ k , we conclude that the sequence {x kz k } is bounded. It now easily follows that the sequences {z k }, {e k }, and {y k } are all bounded. By the triangle and Cauchy-Schwarz inequalities, and the nonexpansiveness of the projection operator, we have

x k -z k ≥ x k -(z k -e k ) -e k = x k -[z k -ϕ k (z k )] + -e k ≥ x k -[x k -F (x k )] + -[x k -F (x k )] + -[z k -ϕ k (z k )] + -e k ≥ r(x k ) -x k -z k -F (x k ) + ϕ k (z k ) -e k ≥ r(x k ) -(1 + ρ k μ k ) x k -z k -(G k + μ k I)(x k -z k ) ≥ r(x k ) -(1 + C 1 + 2C 2 ) x k -z k .
Hence,

(2 + C 1 + 2C 2 ) x k -z k ≥ r(x k ) . It follows from (3.7) that 0 < lim inf k→∞ x k -z k . (3.9)
Suppose that condition (2.11) in Algorithm 2.2 holds an infinite number of times. For such iterations k, by (2.5), we have (recall also that v k = 0)

v k , x k -y k v k ≥ μ k (1 -σ) x k -y k 2 v k ≥ C 4 (1 -σ) x k -z k + e k 2 F (z k -e k ) -ϕ k (z k ) + e k . (3.10) Since G k ≤ C 1 , μ k ≤ C 2 and {z k }, {e k } are bounded, it follows that {F (z k -e k )} and {ϕ k (z k )} are bounded. Therefore, for some C 5 > 0, C 4 (1 -σ) F (z k -e k ) -ϕ k (z k ) + e k ≥ C 5 .
Passing onto the limit in (3.10) (along the indices k for which (2.11) holds) and taking into account (3.8), we obtain that

0 = lim inf k→∞ x k -z k + e k .
By the triangle inequality and (2.9), we have

x k -z k + e k ≥ x k -z k -e k ≥ (1 -ρ k μ k ) x k -z k ≥ (1 -ρ k C 2 ) x k -z k . Because 1/C 2 > lim sup k→∞ ρ k , we further conclude that 0 = lim inf k→∞ x k -z k
which contradicts (3.9). We conclude that if 0 < lim inf k→∞ r(x k ) , then condition (2.11) in Algorithm 2.2 may hold no more than a finite number of times. Hence, we can assume that for all k sufficiently large, y k and v k are obtained through the linesearch step (2.12), in which case

v k , x k -y k v k = α k F (y k ), x k -z k F (y k ) ≥ α k λ(1 -ρ k )μ k x k -z k 2 F (y k ) .
Using (3.8) and taking into account the boundedness of {F (y k )} and the fact that μ k ≥ C 4 and 1 > lim sup k→∞ ρ k , we obtain that

0 = lim k→∞ α k x k -z k .
Now, because of (3.9), we conclude that it must be the case that

0 = lim k→∞ α k .
The latter is equivalent to saying that m k → ∞. It follows that for every sufficiently large k, the stepsize is decreased at least twice, i.e., m k ≥ 2. Hence, the stepsize rule (2.12) is not satisfied for the value of β m k -1 , i.e.,

F (x k + β m k -1 (z k -x k )), x k -z k < λ(1 -ρ k )μ k z k -x k 2 .
Taking into account boundedness of the sequences {x k }, {μ k }, {ρ k }, and {z k }, and passing onto a subsequence if necessary, as k → ∞ we obtain

F (x), x -ẑ ≤ λ(1 -ρ)μ ẑ -x 2 ,
where x, μ, ρ, and ẑ are limits of corresponding subsequences. On the other hand, passing onto the limit in (3.1), we have that

F (x), x -ẑ ≥ (1 -ρ)μ ẑ -x 2 .
Taking into account that μ > 0 and ẑx > 0 (by (3.9)), and ρ ≤ lim sup k→∞ ρ k < 1, the last two relations are a contradiction because λ ∈ (0, 1). Hence the case 0 < lim inf k→∞ r(x k ) is not possible. This completes the proof.

Remark. Since boundedness of {x k } was established without any boundedness assumptions on {G k }, for the special choice of G k = ∇F (x k ) the condition G k ≤ C 1 can be removed, due to the continuity of ∇F (•).

Local convergence analysis.

The following error bound [START_REF] Pang | Error bounds in mathematical programming[END_REF] result will be crucial for establishing the superlinear rate of convergence of our algorithm. Note that this error bound actually holds under the more general assumption of x being a regular solution [START_REF] Robinson | Strongly regular generalized equations[END_REF][START_REF] Pang | Inexact Newton methods for the nonlinear complementarity problem[END_REF][START_REF] Harker | Finite-dimensional variational inequality problems: A survey of theory, algorithms and applications[END_REF]. Positive definiteness of ∇F (x) is a simple sufficient condition for x to be regular [START_REF] Peng | A hybrid Newton method for solving the variational inequality problem via the D-gap function[END_REF]. Here we state only the simplified result.

Lemma 4.1 ([14, Proposition 4.4]). If x is a solution of NCP(F ) where ∇F (x) is positive definite, then there exist a constant θ > 0 and a neighborhood B of x such that xx ≤ θ r(x) for all x ∈ B.

We will also need the following error bound result for strongly monotone linear complementarity problems [START_REF] Pang | A posteriori error bounds for the linearly-constrained variational inequality problem[END_REF]. This result is actually related to Lemma 4.1, but note that the constant θ can be estimated explicitly and the error bound holds globally.

Lemma 4.2. Let ẑ be the solution of the linear complementarity problem

z ≥ 0, Mz + q ≥ 0, Mz + q, z = 0,
where M is a positive definite matrix and q is an arbitrary vector. For all z ∈ n it holds that

z -ẑ ≤ 1 + M C(M ) min{z; Mz + q} ,
where C(M ) > 0 is the smallest eigenvalue of (M + M )/2.

We are now ready to establish the superlinear convergence of our algorithm for solving the monotone nonlinear complementarity problems under the assumptions of positive definiteness of ∇F at the solution and its local Hölder continuity. The proof relies on the fact that in a sufficiently small neighborhood of such a solution, the approximate Newton point computed by Algorithm 2.2 satisfies the error tolerance requirements of Algorithm 2.1 for solving the corresponding nonlinear proximal point subproblem. Therefore the corresponding projection step in Algorithm 2.2 is taken immediately, and the linesearch procedure is not used.

Theorem 4.3. Let F be monotone and continuous on n . Let x be the (unique) solution of NCP(F ) at which F is differentiable with ∇F (x) positive definite. Let ∇F be locally Hölder continuous around x with degree p ∈ (0, 1]. Suppose that the assumptions of Theorem 3.2 are satisfied. In addition, suppose that

μ k = r(x k ) t , t ∈ (0, p), 0 = lim k→∞ ρ k ,
and starting with some index k 0 , G k = ∇F (x k ).

Then the sequence {x k } converges to x Q-superlinearly. Proof. By Theorem 3.2, we already know that the sequence {x k } converges to x. (Note that monotonicity of F and nonsingularity of ∇F (x) imply that x is the unique solution.) Note that by the choice of μ k , we have that 0 = lim k→∞ μ k .

By Hölder continuity of ∇F in the neighborhood of x, it follows that there exists some L > 0 such that for all u ∈ n sufficiently small and all indices k sufficiently large

∇F (x k + u) -∇F (x k ) ≤ L u p , p ∈ (0, 1]. Therefore F (x k + u) -F (x k ) = 1 0 ∇F (x k + su)u ds = ∇F (x k )u + 1 0 (∇F (x k + su) -∇F (x k ))u ds = ∇F (x k )u + R k (u), where R k (u) ≤ 1 0 ∇F (x k + su) -∇F (x k ) u ds ≤ L 1 0 s p u 1+p ds = L 1 + p u 1+p .
Hence, defining C 6 := L/(1 + p), we have that

F (x k + u) -F (x k ) -∇F (x k )u = R k (u), (4.1) R k (u) ≤ C 6 u 1+p .
We next show that (2.11) is always satisfied for all k sufficiently large, so that the linesearch step is never used. Thus let

y k = z k -e k and v k = F (z k -e k )-ϕ k (z k )+e k .
Assume that k is sufficiently large so that (4.1) holds and y kx k is sufficiently small (later we shall verify that {y kx k } → 0). Then we have that

-ε k = v k + μ k (y k -x k ) = F (z k -e k ) -ϕ k (z k ) + e k + μ k (z k -e k -x k ) = F (z k -e k ) -F (x k ) -∇F (x k )(z k -e k -x k ) -∇F (x k )e k + (1 -μ k )e k = R k (z k -e k -x k ) -∇F (x k )e k + (1 -μ k )e k . (4.2)
By (4.1) and the Cauchy-Schwarz and triangle inequalities, it follows that

ε k ≤ C 6 y k -x k 1+p + (C 1 + 1) e k ≤ C 6 y k -x k 1+p + (C 1 + 1)ρ k μ k z k -x k . (4.3) Furthermore, z k -x k ≤ y k -x k + e k ≤ y k -x k + ρ k μ k z k -x k .
Hence,

z k -x k ≤ C 7 y k -x k , (4.4) where C 7 ≥ 1/(1 -ρ k μ k ) (recall that ρ k μ k → 0). Defining C 8 := C 7 (C 1 + 1)
, and combining (4.4) with (4.3), we obtain

ε k ≤ (C 6 y k -x k p + C 8 ρ k μ k ) y k -x k = C 6 y k -x k p r(x k ) t + C 8 ρ k μ k y k -x k .
Clearly, condition (2.11) is satisfied whenever Thus for establishing that the linesearch procedure is never used for indices k sufficiently large, it is left to prove (4.5). (This would also imply that {y kx k } → 0.)

C 6 y k -x k p r(x k ) t + C 8 ρ k ≤ σ.
Let ẑk be the exact solution of LCP(ϕ k ). It holds that

F (x) -ϕ k (ẑ k ), x -ẑk = -ẑk , F (x) -x, ϕ k (ẑ k ) ≤ 0, because x, F (x), ẑk , ϕ k (ẑ k
) are all nonnegative and x, F (x) = ẑk , ϕ k (ẑ k ) = 0. Furthermore, by (4.1) we obtain

F (x) -ϕ k (ẑ k ) = F (x) -F (x k ) -(∇F (x k ) + μ k I)(ẑ k -x k ) = F (x) -F (x k ) -∇F (x k )(x -x k ) -∇F (x k )(ẑ k -x) -μ k (ẑ k -x k ) = R k (x -x k ) -∇F (x k )(ẑ k -x) -μ k (ẑ k -x k ). Hence, 0 ≥ F (x) -ϕ k (ẑ k ), x -ẑk = R k (x -x k ), x -ẑk -∇F (x k )(ẑ k -x), x -ẑk -μ k ẑk -x k , x -ẑk . (4.6)
By positive definiteness of ∇F (x), it follows that for some C 9 > 0 and all indices k sufficiently large

∇F (x k )d, d ≥ C 9 d 2 for all d ∈ n . (4.7)
By (4.6) and the latter relation we obtain

C 9 ẑk -x 2 ≤ ∇F (x k )(ẑ k -x), ẑk -x ≤ R k (x -x k ), ẑk -x + μ k ẑk -x k , x -ẑk ≤ C 6 x k -x 1+p ẑk -x + μ k ẑk -x k ẑk -x ,
where the last inequality follows from (4.1) and the Cauchy-Schwarz inequality. Therefore,

C 9 ẑk -x ≤ C 6 x k -x 1+p + μ k ẑk -x k . (4.8) By the triangle inequality, ẑk -x k ≤ ẑk -x + x k -x .
Combining the latter inequality with (4.8), we obtain ẑk -

x k ≤ C -1 9 C 6 x k -x 1+p + μ k ẑk -x k + x k -x . Hence, (1 -μ k /C 9 ) ẑk -x k ≤ C 6 /C 9 x k -x 1+p + x k -x .
And for

C 11 ≥ 1/(1 -μ k /C 9 ) (recall that μ k → 0), we obtain ẑk -x k ≤ C 11 1 + C 6 /C 9 x k -x p x k -x . (4.9)
By the triangle inequality,

z k -x k ≤ ẑk -z k + ẑk -x k . (4.10)
Furthermore by Lemma 4.2, we have that

z k -ẑk ≤ 1 + ∇F (x k ) + μ k I c k e k ,
where c k is the smallest eigenvalue of (∇F (x k ) + ∇F (x k ) )/2 + μ k I. In view of (4.7), clearly c k ≥ C 9 . Therefore

z k -ẑk ≤ 2 + C 1 C 9 e k .
Furthermore, defining C 10 := (2 + C 1 )/C 9 , we obtain

z k -ẑk ≤ C 10 ρ k μ k z k -x k . (4.11)
Combining the latter relation with (4.10), we have that

z k -x k ≤ C 10 ρ k μ k z k -x k + ẑk -x k .
Since ρ k μ k → 0, taking into account (4.9), it is now clear that there exists C 12 > 0 such that

z k -x k ≤ C 12 x k -x . (4.12) Note that y k -x k ≤ z k -x k + e k ≤ (1 + ρ k μ k ) z k -x k ≤ 2 z k -x k ≤ 2C 12 x k -x , (4.13)
where the last inequality follows from (4.12). Hence, by Lemma 4.1,

y k -x k ≤ 2C 12 θ r(x k ) . (4.14) Therefore y k -x k p r(x k ) t ≤ 2C 12 θ r(x k ) p-t , (4.15)
which establishes (4.5) by the choice of t ∈ (0, p).

Hence, from now on, we can assume that the linesearch is never used, so that

y k = z k -e k and v k = F (z k -e k ) -ϕ k (z k ) + e k . Recalling the definition xk := x k - v k , x k -y k v k 2 v k ,
it is easy to see that x k+1x ≤ xkx because x ∈ n + and x k+1 is the orthogonal projection of xk onto n + . Applying further the triangle inequality, we obtain

x k+1 -x ≤ xk -y k + y k -z k + z k -ẑk + ẑk -x . (4.16)
We proceed to analyze the four terms in the right-hand side of (4.16).

For the second term in (4.16) we have

y k -z k = e k ≤ ρ k μ k z k -x k ,
and for the third term (see (4.11)) we have

z k -ẑk ≤ C 10 e k ≤ C 10 ρ k μ k z k -x k .
Using (4.12) we obtain

y k -z k + z k -ẑk ≤ (1 + C 10 )ρ k μ k z k -x k ≤ (1 + C 10 )C 12 ρ k r(x k ) t x k -x . (4.17)
Next, we consider the last term in (4.16). By (4.8), we have

C 9 ẑk -x ≤ C 6 x k -x 1+p + μ k ẑk -x k ≤ C 6 x k -x 1+p + μ k ( ẑk -x + x k -x ),
where the last inequality follows from the triangle inequality. Therefore,

(C 9 -μ k ) ẑk -x ≤ C 6 x k -x 1+p + μ k x k -x .
Since μ k → 0, for C 13 ≥ max{1; C 6 }/(C 9μ k ), we have ẑkx ≤ C 13 x kx p + r(x k ) t x kx . where ξ k is the angle between xkx k and y kx k . Because xkx k = -s k v k for a certain s k > 0, the angle between the vectors v k and -μ k (y kx k ) is also ξ k . Observe from (4.2) that

v k = -μ k (y k -x k ) + R k (y k -x k ) -∇F (x k )e k + (1 -μ k )e k .
Given the above relation, sin of the angle between v k and -μ k (y kx k ) can be easily bounded:

sin ξ k ≤ R k (y k -x k ) + ((1 -μ k )I -∇F (x k ))e k μ k y k -x k ≤ C 6 y k -x k 1+p + (1 + C 1 ) e k μ k y k -x k ≤ C 6 y k -x k p r(x k ) t + (1 + C 1 ) ρ k z k -x k y k -x k ≤ 2C 6 C 12 θ r(x k ) p-t + C 7 (1 + C 1 )ρ k ,
where the second inequality follows from (4.1) and the last inequality follows from (4.15) and (4.4). Hence, for some C 14 > 0, sin ξ k ≤ C 14 ( r(x k ) p-t + ρ k ). By (4.16), combining (4.21), (4.17), and (4.18), we conclude that there exists C 15 > 0 such that x k+1x ≤ C 15 r(x k ) p-t + ρ k + (1 + ρ k ) r(x k ) t + x kx p × x kx , (4.22) which means that {x k } converges to x at least superlinearly.

Remark. Finally, we note that

r(x) = x -[x -F (x)] + -x + [x -F (x)] + ≤ x -x + [x -F (x)] + -[x -F (x)] + ≤ (2 + C 16 ) x -x ,
where C 16 is the constant of local Lipschitz continuity of F . Therefore if one further chooses relation (4.22) implies that the order of superlinear convergence is at least 1 + pt, where t ∈ (0, p).

ρ k = r(x k ) p ,

Concluding remarks.

We presented a new globalization strategy for the Newton method applied to the monotone nonlinear complementarity problem. Our strategy is based on the projection-proximal point methodology and makes full use of the monotonicity structure of the problem. The resulting hybrid algorithm is truly globally convergent to a solution without any additional assumptions, even if the solution set is unbounded. This is an important property which is not possessed by globalization approaches based on merit functions. Under natural assumptions, locally superlinear rate of convergence was also established.

  is satisfied for all indices k sufficiently large if it holds that 0 = lim k→∞ y kx k p r(x k ) t . (4.5)

  (4.18) Finally, we consider the first term in the right-hand side of (4.16). Since the point xk is the projection of x k onto the hyperplaneH k = {x | v k , xy k = 0},and y k ∈ H k , the vectors xkx k and xky k are orthogonal. Hence, xky k = y kx k sin ξ k , (4.19)

(4. 20 )

 20 By(4.19), using (4.13) and (4.20), we obtainxky k = y kx k sin ξ k ≤ 2C 12 C 14 r(x k ) p-t + ρ k x kx . (4.21)

Acknowledgments. We thank the two anonymous referees and the associate editor for constructive suggestions which helped us to improve the paper.

The research of the first author was supported by CNPq grant 300734/95-6 and by PRONEX-Optimization. The research of the second author was supported by CNPq grant 301200/93-9(RN) and by PRONEX-Optimization.