Mounir Haddou
email: mounir.haddou@insa-rennes.fr

Tangi Migot
email: tangi.migot@insa-rennes.fr

Jérémy Omer
email: jeremy.omer@insa-rennes.fr

A Generalized Direction in Interior Point Method for Monotone Linear Complementarity Problems

Keywords: concave functions, interior-point methods, linear programming, linear complementarity problems, polynomial time complexity Mathematics Subject Classification

. We prove that this new method possesses the best known upper bound complexity for these methods. Moreover, we extend results known in the literature since we consider a general family of smooth concave functions in the Newton system instead of the square root.

Introduction

In this paper we focus on the simplest interior point methods (IPMs) : the full Newton step interior point method [START_REF] Roos | Interior point methods for linear optimization[END_REF], this method is one among the IPMs with the best worst-case complexity. Therefore, naturally any new try in the IPM framework must be validated on this method. In particular, we focus on this method applied to the monotone complementarity problem [START_REF] Kojima | A unified approach to interior point algorithms for linear complementarity problems[END_REF]. This problem is a generalization of the linear optimization problem and has been very popular in the literature due to its numerous applications [START_REF] Facchinei | Finite-dimensional variational inequalities and complementarity problems[END_REF]. In the literature of IPM, this problem has also been an active subject [START_REF] Kojima | A unified approach to interior point algorithms for linear complementarity problems[END_REF].

Darvay [START_REF] Zs | A new algorithm for solving self-dual linear optimization problems[END_REF][START_REF] Zs | New interior-point algorithms in linear programming[END_REF], introduces a modification in the interior point method for finding search directions for linear optimization problems, based on an algebraic transformation of the central path. In particular, he applied the square root function to both sides of the centering equation, and he used Newton's method to obtain the new direction. He proved that the corresponding full Newton step algorithm has O( √ nL) iteration complexity. This new direction using the square root has become an active subject in the past few years [START_REF] Achache | A weighted path-following method for the linear complementarity problem[END_REF][START_REF] Achache | A new primal-dual path-following method for convex quadratic programming[END_REF][START_REF] Achache | Complexity analysis and numerical implementation of a short-step primaldual algorithm for linear complementarity problems[END_REF][START_REF] Asadi | Polynomial interior-point algorithm for P * (κ) horizontal linear complementarity problems[END_REF][START_REF] Asadi | A new full-Newton step O(n) infeasible interior-point algorithm for P * (κ)-horizontal linear complementarity problems[END_REF][START_REF] Asadi | An infeasible interior-point method with improved centering steps for monotone linear complementarity problems[END_REF][START_REF] Asadi | On the P * (κ) horizontal linear complementarity problems over Cartesian product of symmetric cones[END_REF][START_REF] Asadi | An infeasible full-Newton step interior point method for horizontal linear complementarity problem over Cartesian product of symmetric cones[END_REF][START_REF] Asadi | A corrector-predictor interior-point algorithm for P * (κ)-HLCPs over Cartesian product of symmetric cones[END_REF][START_REF] Wang | A new primal-dual path-following interior-point algorithm for semidefinite optimization[END_REF][START_REF] Wang | A primal-dual interior-point algorithm for second-order cone optimization with full Nesterov-Todd step[END_REF][START_REF] Wang | A new full Nesterov-Todd step primal-dual path-following interior-point algorithm for symmetric optimization[END_REF][23]. Several authors generalized this approach to a wide class of optimization problems, for example for linear complementarity problems [1, 3-6, 8, 23], convex quadratic programming [START_REF] Achache | A new primal-dual path-following method for convex quadratic programming[END_REF], second-order cone optimization [START_REF] Wang | A primal-dual interior-point algorithm for second-order cone optimization with full Nesterov-Todd step[END_REF], semidefinite optimization [START_REF] Wang | A new primal-dual path-following interior-point algorithm for semidefinite optimization[END_REF] and symmetric cone optimization [START_REF] Asadi | On the P * (κ) horizontal linear complementarity problems over Cartesian product of symmetric cones[END_REF][START_REF] Asadi | An infeasible full-Newton step interior point method for horizontal linear complementarity problem over Cartesian product of symmetric cones[END_REF][START_REF] Asadi | A corrector-predictor interior-point algorithm for P * (κ)-HLCPs over Cartesian product of symmetric cones[END_REF][START_REF] Wang | A new full Nesterov-Todd step primal-dual path-following interior-point algorithm for symmetric optimization[END_REF].

Inspired by Darvay's new approach we introduce here a new class of IPMs by considering a large family of smooth concave functions instead of the square root. This new class of methods generalizes the classical path-following IPMs, since we obtain them as a special case. The technique presented here does not include Darvay's algorithm, but we can consider a smoother version that belongs to our family of methods. Our main contribution is that we prove that the algorithm with the new directions converges to a solution with the best known complexity for this family of methods.

In Section 2, we introduce the problem and our new directions. In Section 3, we show the polynomial complexity of our new class of methods.

Notations

Through this paper we will use the following notations: R n

+ = {x ∈ R n | x ≥ 0}, R n ++ = {x ∈ R n | x >
0} and e denotes the vector with all entries equal to one and whose dimension can be deduced from the context. Given two vectors z, s ∈ R n , we denote by z T s the usual scalar product and by zs the Hadamard product of two vectors, that is zs = (z i s i ) 1≤i≤n . Moreover, we extend this component-wise operation to the division of two vectors and to the square root, that is

∀z ∈ R n , s ∈ (R \ {0}) n , z/s = (z i /s i ) 1≤i≤n and ∀z ∈ R n + , √ z = ( √ z i ) 1≤i≤n .

Preliminaries and Problem Settings

An LCP consists in finding z, s ∈ R n such that for a square matrix M ∈ R n×n and a vector q ∈ R n s = M z + q, z, s ≥ 0, zs = 0. (LCP)

A couple (z, s) such that s = M z + q is said to be feasible for the LCP if we have z, s ≥ 0 and strictly feasible if z, s > 0. From now on, we use standard notation

F + := {(z, s) ∈ R n ++ × R n ++ | s = M z + q}
for the set of strictly feasible points of (LCP). In this paper we consider a monotone linear complementarity problem, i.e. an LCP where M is positive semi-definite. In this case the set of solutions of (LCP) is a convex set.

The main strategy of IPMs is to follow the central path (z(µ), s(µ)) for µ ∈ R n ++ , defined by

s = M z + q, z, s ≥ 0, zs = µ . (LCP µ )
The couples (z(µ), s(µ)) are also called µ-centers and define an homotopic path. The limit when µ → 0 satisfies the complementarity condition, and hence yields optimal solutions whenever the limit exists. IPMs follow the central path approximately by solving an approached version of the non-linear system in (LCP µ ) for several values of µ. The main tool to solve such a system is the Newton method. A Newton step (∆z, ∆s) is given as the solution of the following linear system

M ∆z = ∆s z∆s + s∆z = µ -zs (1)
There exists a wide variety of different IPMs that are based on this principle. In this paper we focus on the simplest IPM (see Algorithm 1): the full Newton step interior point method (FN-IPM).

Input: an accuracy parameter ϵ > 0 ; a sequence of update parameters In [START_REF] Zs | New interior-point algorithms in linear programming[END_REF], Darvay introduces a modification in (LCP µ ) by considering

{θ k }, 0 < θ k < 1 ∀k ∈ N ; initial values (z 0 , s 0 ) ∈ F + , µ 0 = z 0 s 0 ; 1 z := z 0 , s := s 0 , µ := µ 0 , k := 0 ; 2 while z T s ≥ nϵ do
s = M z + q, z, s ≥ 0, φ(zs) = µ , (2)
where φ : R n + → R n + is assumed to be a smooth function such that φ(0) = 0 and φ is defined by a component-wise extension of a real-valued function

φ to R n , i.e. for t ∈ R n φ(t) = (φ(t i )) 1≤i≤n .
Darvay's method modifies the Newton steps. More precisely, the Newton step (∆z, ∆s) is given by the linear system

M ∆z = ∆s φ ′ (zs)(z∆s + s∆z) = µ -φ(zs) (3)
In this paper we consider functions φ ∈ C 3 (R + ) such that φ(0) = 0, φ is increasing and concave, and φ ′′′ (t) ≥ 0 ∀t ∈ R + . This class of functions allows for a generalization of the classical IPMs, since we obtain the classical central path system (LCP µ ) for φ(t) = t. The square root function does not belong to this family since it is not differentiable at 0, but we can build a smooth version with φ ϵ>0 : t → √ t + ϵ-√ ϵ. We modify Algorithm 1 to solve (3) instead of (1) at step 4, and call the resulting algorithm φ-FN-IPM. The main result of this article is that φ-FN-IPM, converges to an ϵ-solution in at most O √ n log( n ϵ ) iterations. This upper bound is the best one known for the FN-IPM.

Polynomial Complexity

In this section, we consider the worst-case complexity of the φ-FN-IPM described in Algorithm 2 with φ :

R + → R + verifying (i) φ(0) = 0; (ii) φ ∈ C 3 ([0, +∞)); (iii) φ ′ (t) > 0, ∀t ≥ 0, i.e. φ is increasing; (iv) φ ′′ (t) ≤ 0, ∀t ≥ 0, i.e. φ is concave; (v) φ ′′′ (t) ≥ 0, ∀t ≥ 0.
These functions are invertible and can be extended in a smooth way for negative t by considering :

φ(t) = tφ ′ (0)+ t 2 2 φ ′′ (0)+ t 3 6 φ ′′′ (0). Function φ : R n → R n is then defined component-wisely: ∀t ∈ R n , φ(t) = (φ(t i )) 1≤i≤n .
Input: an accuracy parameter ϵ > 0 ; a sequence of update parameters {θ k }, 0 < θ k < 1 ∀k ∈ N ; initial values (z 0 , s 0 ) ∈ F + , µ 0 = z 0 s 0 ; 1 z := z 0 , s := s 0 , µ := µ 0 , k := 0 ; 2 while n i=1 φ(z i s i ) ≥ nϵ do Algorithm 2: φ-Full Newton step IPM (φ-FN-IPM)

One important characteristic of φ is the existence of a constant T defined by

-φ ′′ (0) = T (φ ′ (0)) 2 . ( 4
)
Noticing the conditions (iii) and (iv), we have T ≥ 0 and T = 0 for φ(t) = t.

It should be noted that arbitrary values of T > 0 can be achieved by scaling φ.

Note that t → √ t, in the same way as any function t → t q , 0 < q < 1, does not satisfy these hypotheses since it is not differentiable in 0. However we can consider a smooth version for ϵ > 0 with t → (t + ϵ) q -(ϵ) q . As said in the introduction the classical method is given by φ(t) = t. Other examples are φ : t → log(1 + t) and functions constructed as in [START_REF] Haddou | Smoothing methods for nonlinear complementarity problems[END_REF][START_REF] Migot | A smoothing method for sparse optimization over polyhedral sets[END_REF]: for instance φ : t → t t+1 and φ : t → 1 -exp(-t). It is interesting to note two latter functions φ are bounded. Moreover, notice that if a function φ satisfies all these hypotheses, then t → αφ(Ct) with α, C ∈ R ++ also satisfies these hypotheses.

The main result of this section (Theorem 3.5) states the polynomial worstcase complexity of the φ-FN-IPM described in Algorithm 2. In order to achieve this result we define a measure of the proximity to the central path in Section 3.1. Then, in Section 3.2, Theorem 3.1 estimates the error made at each Newton step. Sections 3.3 and 3.4 present conditions to ensure the correct behaviour of the algorithm: strict feasibility of the iterates (Theorem 3.2) and quadratic convergence of the Newton process (Theorem 3.3). Section 3.5 provides the sequence of update parameters (Theorem 3.4).

Proximity Measure

At each iteration, after updating parameter µ, we compute the Newton direction (∆z, ∆s) as a solution of system (3). Then we update the iterates with z + = z + α∆z and s + = s + α∆s .

(

) 5
Note that here we consider a damping factor α ∈ [0, 1] to be more general. In this case we denote by (α∆z, α∆s) the Newton step with length α and call it the α-Newton step. Then the full Newton step is given for α = 1.

In order to measure the distance to the target on the central path we consider a proximity measure δ φ (z, s, µ) defined by

δ φ (z, s, µ) := 1 2 φ ′ (0) φ ′ (zs) (v φ (z, s, µ) -1 -v φ (z, s, µ) 2 , with v φ (z, s, µ) := φ(zs) µ and v φ (z, s, µ) -1 := µ φ(zs)
.

We may omit the arguments of v φ (z, s, µ) and δ φ (z, s, µ), when it is clear from the context. Notice that this proximity measure is a generalization of the one presented in [START_REF] Roos | Interior point methods for linear optimization[END_REF], where the authors consider

δ t →t := 1 2 µ zs - zs µ 2 .
Both proximity measures are equal for φ(t) = t. Moreover for any function φ, the two proximity measures are asymptotically similar (for zs ↓ 0).

The following two lemmas link the iterates and the proximity measure.

Lemma 3.1 If (z, s) ∈ F + and δ φ ≤ 1, then φ(zs) ≤ 6µ.
Proof Assume by contradiction that there exists i ∈ {1, . . . , n} such that φ(z i s i ) > 6µ i . Since φ ′ is decreasing and zs > 0:

2δ φ ≥ v φ -1 -v φ 2 = √ µ φ(zs) - φ(zs) √ µ 2 ≥ √ µ i φ(z i s i ) - φ(z i s i ) √ µ i = φ(z i s i ) √ µ i - √ µ i φ(z i s i ) > √ 6 - 1 √ 6 ≈ 2, 04 ,
where the penultimate step comes from the increasing property of the function

x → x -1/x on R ++ . This is in contradiction with δ φ ≤ 1. ⊓ ⊔
The following lemma gives bounds on δ φ that depend on some constant Γ (∥µ∥ ∞ ) defined as

Γ (∥µ∥ ∞ ) := 1 - φ ′′ (0)φ -1 (6∥µ∥ ∞ ) φ ′ (φ -1 (6∥µ∥ ∞ )) . Lemma 3.2 Let δ = 1 2 v φ -1 -v φ 2 , then δ ≤ δ φ ≤ 1 - φ ′′ (0)∥zs∥ ∞ φ ′ (∥zs∥ ∞ ) δ .
Furthermore in a close neighbourhood of the central path, i.e. δ φ ≤ 1, we have

δ ≤ δ φ ≤ Γ (∥µ∥ ∞ )δ . (6)
Proof By convexity of function φ ′ for all i ∈ {1, . . . , n}

φ ′ (z i s i ) ≥ φ ′ (0) + φ ′′ (0)z i s i .
Then, for all i ∈ {1, . . . , n}

1 ≤ φ ′ (0) φ ′ (z i s i ) ≤ 1 - φ ′′ (0)z i s i φ ′ (z i s i ) ≤ max i 1 - φ ′′ (0)z i s i φ ′ (z i s i ) .
Hence, by definition of δ and

δ φ δ ≤ δ φ ≤ 1 - φ ′′ (0)∥zs∥ ∞ φ ′ (∥zs∥ ∞ ) δ .
The sharpest result when δ φ ≤ 1 is deduced from Lemma 3.1.

⊓ ⊔

In the previous lemma, equation [START_REF] Asadi | An infeasible interior-point method with improved centering steps for monotone linear complementarity problems[END_REF], we define Γ as a function of ∥µ∥ ∞ which depends on the choice of φ. For φ(t) = t, we get Γ (∥µ∥ ∞ ) = 1 for all µ. Moreover, for any function φ, Γ is increasing with respect to ∥µ∥ ∞ , and converges to 1 as ∥µ∥ ∞ ↓ 0. Moreover, in the course of the proof we showed that if

δ φ ≤ 1 φ ′ (0) φ ′ (z i s i ) ≤ Γ (∥µ∥ ∞ ), ∀i ∈ {1, . . . , n}. (7)
This result will be useful in a future proof.

Error Bound of the Newton Step

We use the first order Taylor-Lagrange formula applied to φ in zs. There exists

ξ ∈ [z + s + , zs] (or ξ ∈ [zs, z + s + ] if zs < z + s + ) such that φ(z + s + ) = φ(zs) + αφ ′ (zs)(z∆s + s∆z) + α 2 ∆z∆sφ ′ (zs) + φ ′′ (ξ) 2 z + s + -zs 2 , (8)
with (∆z, ∆s) solution of (3). The update of µ will be chosen such that φ(z + s + ) < φ(zs) and thus 0 < z + s + < zs. Therefore the error we make when we say that φ(z

+ s + ) is the µ-center is η(α) := α 2 ∆z∆sφ ′ (zs) + φ ′′ (ξ) 2 z + s + -zs 2 . ( 9
)
The following sequence of lemmas aims to bound this error in terms of the proximity measure. Before doing so, we recall a useful lemma from [START_REF] Illés | Polynomial affine-scaling algorithms for P * (κ) linear complementary problems[END_REF].

Lemma 3.3 (Lemma 5.1, [START_REF] Illés | Polynomial affine-scaling algorithms for P * (κ) linear complementary problems[END_REF]) Let (z, s) ∈ F + and a ∈ R n . Assume that matrix M is a positive semidefinite matrix. Let (∆z, ∆s) be the solution of

-M ∆z + ∆s = 0 s∆z + z∆s = a .
Then

||∆z∆s|| 1 ≤ C 1 a √ zs 2 2 , ||∆z∆s|| 2 ≤ C 2 a √ zs 2 2 , ||∆z∆s|| ∞ ≤ C ∞ a √ zs 2 2 , with C 1 = 1/2, C 2 = 1/(2 √ 2) and C ∞ = 1/4.
Through the rest of this article we denote by C p for p ∈ {1, 2, ∞} the constants defined as

C 1 = 1 2 , C 2 = 1 2 √ 2 and C ∞ = 1 4 .
Straightforward application of this lemma for a= µ-φ(zs) φ ′ (zs) and v φ = φ(zs) µ

gives the following lemma.

Lemma 3.4 Let (z, s) ∈ F + , (∆z, ∆s) be the solution of (3) and p ∈ {1, 2, ∞}, then

∥∆z∆s∥ p ≤ C p φ(zs) √ µ √ zs φ ′ (zs) (v φ -1 -v φ ) 2 2
.

The next lemma will bound ||∆z∆sφ ′ (zs)|| p for p ∈ {1, 2, ∞}.

Lemma 3.5 Let (z, s) ∈ F + , (∆z, ∆s) be the solution of (3) and p ∈ {1, 2, ∞},

then ||∆z∆sφ ′ (zs)|| p ≤ ||∆z∆sφ ′ (0)|| p ≤ 4C p ∥µ∥ ∞ δ 2 φ .
Proof By concavity of φ we have that φ(z i s i ) ≤ φ ′ (0)z i s i , ∀i, so

φ(zs)φ ′ (0) zs(φ ′ (zs)) 2 ≤ (φ ′ (0)) 2 (φ ′ (zs)) 2 . ( 10
)
Furthermore for p ∈ {1, 2, ∞} and using Lemma 3.4, followed by ( 10)

||∆z∆sφ ′ (zs)|| p ≤ ||∆z∆s|| p φ ′ (0) ≤ C p φ(zs) √ µ √ zs φ ′ (zs) ((v φ -1 -v φ )) 2 2 φ ′ (0) ≤ n i=1 C p φ(z i s i )µ i z i s i (φ ′ (z i s i )) 2 (v φ -1 i -v φ i ) 2 φ ′ (0) ≤ n i=1 C p µ i φ ′ (0) φ ′ (z i s i ) 2 (v φ -1 i -v φ i ) 2 ≤ 4C p ∥µ∥ ∞ δ 2 φ .
⊓ ⊔

Now we move to the main result which gives a bound for the complete error.

Lemma 3.6 Let (z, s) ∈ F + and (∆z, ∆s) be the solution of (3). For p ∈ {2, ∞} we have

∥η∥ p ≤ (4C p + 2T φ(∥zs∥ ∞ )) α 2 ∥µ∥ ∞ δ 2 φ + C ∞ φ(∥zs∥ ∞ ) + αC 2 p ∥µ∥ ∞ δ φ 8α 3 T ∥µ∥ ∞ 3/2 δ 3 φ
Proof By equation ( 9) we have

η =α 2 ∆z∆sφ ′ (zs) + α 2 φ ′′ (ξ) 2 (z∆s + s∆z) 2 + α 4 φ ′′ (ξ) 2 (∆s∆z) 2 + α 3 φ ′′ (ξ)(∆s∆z)(z∆s + s∆z)
Taking the p-norm for p ∈ {2, ∞} and using triangle inequalities

∥η∥ p ≤α 2 ∥∆z∆sφ ′ (zs)∥ p +α 2 φ ′′ (ξ) 2 (z∆s+s∆z) 2 p + α 4 φ ′′ (ξ) 2 (∆s∆z) 2 p + α 3 ∥φ ′′ (ξ)(∆s∆z)(z∆s + s∆z)∥ p .
We now bound each term of the above right-hand side. First, Lemma 3.5 gives

∥∆z∆sφ ′ (zs)∥ p ≤ 4C p ∥µ∥ ∞ δ 2 φ .
Using successively ( 4), (µ-φ(zs)

) 2 = (µv φ (v φ -1 -v φ )) 2 = µφ(zs)(v φ -1 -v φ ) 2 and φ(z i s i ) ≤ φ(∥zs∥ ∞ ), ∀i ∈ {1, ..., n}, we obtain step by step φ ′′ (ξ) 2 (z∆s + s∆z) 2 p ≤ φ ′′ (0) 2 (z∆s + s∆z) 2 p = T 2 φ ′ (0) φ ′ (zs) 2 (µ -φ(zs)) 2 p ≤ T φ(∥zs∥ ∞ )∥µ∥ ∞ 2 φ ′ (0) φ ′ (zs) (v φ -1 -v φ ) 2 p ≤ 2T φ(∥zs∥ ∞ )∥µ∥ ∞ δ 2 φ .
To bound the third term, we use Lemma 3.5, equality (4) and the inequality ∥u 2 ∥ p ≤ ∥u∥ 2 p for all u ∈ R n , and we obtain

φ ′′ (ξ) 2 (∆s∆z) 2 p ≤ T 2 (φ ′ (0)∆s∆z) 2 p ≤ T 2 ∥φ ′ (0)∆s∆z∥ 2 p ≤ T 2 (∥µ∥ ∞ 4C p δ 2 φ ) 2 = 8T C 2 p ∥µ∥ ∞ 2 δ 4 φ .
Finally, the definition of v φ implies that

∥φ ′ (0)(z∆s + s∆z)∥ p = φ ′ (0) φ ′ (zs) µv φ (v φ -1 -v φ ) p ≤ φ ′ (0) φ ′ (zs) ∥µ∥ ∞ φ(∥zs∥ ∞ )(v φ -1 -v φ ) p
Using the above inequality, as well as ∥uw∥ p ≤ ∥u∥ ∞ ∥w∥ p , ∀(u, w) ∈ R n ×R n , and (4), we get

∥φ ′′ (ξ)(∆s∆z)(z∆s + s∆z)∥ p ≤ T ∥φ ′ (0)(∆s∆z)∥ ∞ ∥φ ′ (0)(z∆s + s∆z)∥ p ≤ 8T C ∞ φ(∥zs∥ ∞ )∥µ∥ ∞ 3/2 δ 3 φ ,
which completes the proof.

⊓ ⊔

In the special case where we are in a close neighbourhood of the central path we get an improved version of the result: Theorem 3.1 Let (z, s) ∈ F + , δ φ ≤ 1 and (∆z, ∆s) be the solution of (3) and p ∈ {2, ∞}, we have

∥η∥ p ≤ (4C p + 12T ∥µ∥ ∞ )α 2 ∥µ∥ ∞ δ 2 φ + ( √ 6C ∞ + αC 2 p δ φ )α 3 8T ∥µ∥ ∞ 2 δ 3 φ (11)
Proof The proof is similar to the proof of Lemma 3.6, but we use φ(zs) ≤ 6µ from Lemma 3.1 instead of φ(z i s i ) ≤ φ(∥zs∥ ∞ ), ∀i. ⊓ ⊔

For instance, with α = 1, using δ q φ ≤ δ φ for q ≥ 1, (11) becomes

∥η∥ ∞ ≤ ∥µ∥ ∞ + 25 + 4 √ 6 2 T ∥µ∥ ∞ 2 δ 2 φ .
For φ(t) = t we get the same result as in [START_REF] Roos | A full-Newton step O(n) infeasible interior-point algorithm for linear optimization[END_REF]: ∥η∥ ∞ ≤ ∥µ∥ ∞ δ 2 φ .

Feasibility of the Newton

Step A Newton step is feasible (strictly feasible) if the couple (z + , s + ) defined by ( 5) is feasible (strictly feasible).

Theorem 3.2 Let α be in [0,1] and δ φ ≤ 1. The α-Newton step is strictly feasible for (z, s) if

αδ 2 φ < 1 1 + 25+4 √ 6 2 T ∥µ∥ ∞ . ( 12
)
Note that for α = 1, the above condition makes the full Newton step be strictly feasible.

For α = 1 condition [START_REF] Facchinei | Finite-dimensional variational inequalities and complementarity problems[END_REF] implies that the proximity measure, δ φ , must be less than 1, which will not be a restrictive assumption.

Proof For α = 0 the result trivially holds. Therefore, we can assume without loss of generality that α ∈]0, 1]. Let β ∈]0, α], and define z β := z + β∆z and s β := s + β∆s, where (∆z, ∆s) is the solution of (3). The proof rests upon a continuity argument. Indeed, we assume z, s > 0 so φ(zs) > 0 and, using equations ( 3) and (8):

φ(z β s β ) = φ(zs) + βφ ′ (zs)(z∆s + s∆z) + η(β) = φ(zs)(1 -β) + β µ + η(β) β ≥ φ(zs)(1 -β) + β µ -e ∥η(β)∥ ∞ β which is positive for all z, s > 0 if ∥µ∥ ∞ > ∥η(β)∥ ∞ /β. Using Lemma 3.6 this condition holds if ∥µ∥ ∞ >(4C ∞ + 2T φ(∥zs∥ ∞ ))β∥µ∥ ∞ δ 2 φ + C ∞ φ(∥zs∥ ∞ ) + βC 2 ∞ ∥µ∥ ∞ δ φ 8β 2 T ∥µ∥ ∞ 3/2 δ 3 φ .
The right-hand side is increasing with respect to β, so it is sufficient to verify

∥µ∥ ∞ >(4C ∞ + 2T φ(∥zs∥ ∞ ))α∥µ∥ ∞ δ 2 φ + C ∞ φ(∥zs∥ ∞ ) + αC 2 ∞ ∥µ∥ ∞ δ φ 8α 2 T ∥µ∥ ∞ 3/2 δ 3 φ .
Therefore, since δ φ ≤ 1, using Lemma 3.1, it suffices to have

1 -δ 2 φ α 4C ∞ + 12 + 8 √ 6C ∞ + 8C 2 ∞ T ∥µ∥ ∞ > 0 ,
which corresponds to our assumption. It follows that for all β ∈ [0, α], φ(z β s β ) > 0. By continuity of φ, this implies that none of z β or s β vanish for β ∈ [0, α], so the result follows. ⊓ ⊔

Quadratic Decrease of the Proximity Measure

The Newton method is known to behave well in a close neighbourhood of the solution. The following theorem states a condition on the proximity measure, δ φ := δ φ (z, s, µ), that ensures a quadratic convergence of the Newton step. We denote by δ φ := δ φ (z + , s + , µ) and v φ := v φ (z + , s + , µ) the proximity measure and the function v φ after the Newton step.

Theorem 3.3 Let (z, s) ∈ F + and (z + , s + ) be defined as in (5) for α = 1.

Let Γ (∥µ∥ ∞ ) defined in [START_REF] Asadi | An infeasible interior-point method with improved centering steps for monotone linear complementarity problems[END_REF] and

Q(∥µ∥ ∞ ) := 1 -Γ (∥µ∥ ∞ )( √ 2 + (13 + 2 √ 6)T ∥µ∥ ∞ ) 2 /4 1 + 25+4 √ 6 2 T ∥µ∥ ∞ . ( 13
)
If δ 2 φ ≤ Q(∥µ∥ ∞ ), then δ φ ≤ δ 2 φ .
Proof Let Γ (∥µ∥ ∞ ) be defined as in Section 3.1:

δ φ ≤ Γ (∥µ∥ ∞ ) 2 ∥(v φ ) -1 -(v φ )∥ 2 = Γ (∥µ∥ ∞ ) 2 ∥(v φ ) -1 (e -(v φ ) 2 )∥ 2 (14) ≤ Γ (∥µ∥ ∞ ) 2 η(α) µ e + η(α) µ 2 (15) ≤ Γ (∥µ∥ ∞ ) η(α) µ 2 2 1 -η(α) µ ∞ . ( 16
)
Let δ 2 φ ≤ Q(∥µ∥ ∞ ) ≤ 1, using Theorem 3.1 with a full step, i.e. α = 1, it follows

δ φ ≤ δ 2 φ Γ (∥µ∥ ∞ ) 4C 2 + (12 + 8 √ 6C ∞ + 8C 2 2 )T ∥µ∥ ∞ 2 1 -δ 2 φ (4C ∞ + 12 + 8 √ 6C ∞ + 8C 2 ∞ T ∥µ∥ ∞ ) . So, δ φ ≤ δ 2 φ if δ 2 φ ≤ Q(∥µ∥ ∞ ). ⊓ ⊔
Considering φ(t) = t, the condition of Theorem 3.3 becomes the same as in [START_REF] Roos | A full-Newton step O(n) infeasible interior-point algorithm for linear optimization[END_REF]:

δ 2 φ ≤ 1/2.
Remark 1 The condition in Theorem 3.3 implies the condition in Theorem 3.2. So, if the iterates locate in the neighbourhood of quadratic convergence, the full Newton step will provide strictly feasible iterates.

Remark 2 Notice that since the proximity measure is always non-negative, the condition from Theorem 3.3 can hold only when µ is sufficiently small, i.e. when Γ (∥µ∥ ∞ ) √ 2 + (13 + 2 √ 6)T ∥µ∥ ∞ ≤ 2. This is not a restrictive assumption, because we can always scale a given initial point so that it satisfies this condition.

Updating Parameter Strategy

The sequence of parameters {θ k } must be chosen such that the iterates remain strictly feasible and satisfy the condition of Theorem 3.3. In this section, Proposition 3.1 gives an upper bound on the proximity measure after an update on µ, that is µ + = µ(1 -θ k ), and then Theorem 3.4 describes how to build the sequence {θ k }.

First, we provide an upper bound of the proximity measure after an update of the parameter, denoted δ + φ := δ φ (z + , s + , µ(1 -θ k )), in terms of the update θ k , the proximity measure before this update, denoted δ φ := δ φ (z + , s + , µ), and the proximity measure before the Newton step, denoted δ φ := δ φ (z, s, µ). The computation of this upper bound is based on the following lemma. Lemma 3.7 Let (z, s) ∈ F + and (∆z, ∆s) be the solution of (3). Then,

n i=1 φ(z + i s + i ) µ i ≤ n + 2δ 2 φ .
Proof We first show that

∆z∆sφ ′ (zs) µ 1 ≤ 2δ 2 φ . (17)
For this, we apply Lemma 3.3 with ∆z/ √ µ, ∆s/ √ µ and a= µ-φ(zs) √ µφ ′ (zs) to obtain

∆z∆s µ 1 ≤ 1 2 µ -φ(zs) √ µφ ′ (zs) √ zs 2 2
.

Then, we get ( 17) by following the same steps as in the proof of Lemma 3.5.

From equation [START_REF] Asadi | An infeasible full-Newton step interior point method for horizontal linear complementarity problem over Cartesian product of symmetric cones[END_REF] and by concavity of φ, φ(z

+ s + )/µ ≤ e+∆z∆sφ ′ (zs)/µ. So, n i=1 φ(z + i s + i )/µ i ≤ n i=1 1 + |∆z i ∆s i φ ′ (z i s i )/µ i | ≤ n + 2δ 2 φ . ⊓ ⊔ Proposition 3.1 Let v φ := v φ (z + , s + , µ) and v + φ := (z + , s + , µ + ), where µ + := (1 -θ k )µ. Then, (δ + φ ) 2 ≤ (1 -θ k )(δ φ ) 2 + Γ (∥µ∥ ∞ ) 2 4(1 -θ k ) n(θ k ) 2 + (4θ k -2(θ k ) 2 )δ 2 φ .
Furthermore, assuming that

δ 2 φ ≤ Q(∥µ∥ ∞ ) yields (δ + φ ) 2 ≤ (1 -θ k )Q(∥µ∥ ∞ ) 2 + Γ (∥µ∥ ∞ ) 2 4(1 -θ k ) n(θ k ) 2 + (4θ k -2(θ k ) 2 )Q(∥µ∥ ∞ ) . Proof Noticing that v + φ = v φ / √ 1 -θ k , it follows that (δ + φ ) 2 = 1 4 n i=1 φ ′ (0) φ ′ (z + i s + i ) 2 ((v + φ i ) -2 + (v + φ i ) 2 -2) = 1 4 n i=1 φ ′ (0) φ ′ (z + i s + i ) 2 (1 -θ k )(v φ i ) -2 + v φ 2 i (1 -θ k ) -2 = (1 -θ k )( 1 4 n i=1 φ ′ (0) φ ′ (z + i s + i ) 2 (v φ i ) -2 + v φ 2 i -2) + 1 4 n i=1 φ ′ (0) φ ′ (z + i s + i ) 2 -2θ k + 2θ k -(θ k ) 2 1 -θ k v φ 2 i = (1 -θ k )(δ φ ) 2 + 1 4 n i=1 φ ′ (0) φ ′ (z + i s + i ) 2 -2θ k + 2θ k -(θ k ) 2 1 -θ k v φ 2 i
Using successively equation ( 7), Lemma 3.7 and Γ (∥µ

+ ∥ ∞ ) ≤ Γ (∥µ∥ ∞ ), we obtain (δ + φ ) 2 ≤ (1 -θ k )(δ φ ) 2 + Γ (∥µ + ∥ ∞ ) 2 - 2θ k 4 n + 2θ k -(θ k ) 2 4(1 -θ k ) (n + 2δ 2 φ ) ≤ (1 -θ k )(δ φ ) 2 + Γ (∥µ∥ ∞ ) 2 n(θ k ) 2 4(1 -θ k ) + 2θ k -(θ k ) 2 2(1 -θ k ) δ 2 φ .
This proves the first part of the proposition. Now, assuming that δ φ ≤ Q(∥µ∥ ∞ ) allows us to use Theorem 3.3 and so

δ φ ≤ δ 2 φ ≤ Q(∥µ∥ ∞ ) 2 gives the result. ⊓ ⊔
We conclude this section by a description of the choice of the update parameters θ k .

Theorem 3.4 Let µ be such that Q(∥µ∥ ∞ ) > 0 and (z, s) ∈ F + such that δ 2 φ ≤ Q(∥µ∥ ∞ ). Define θ k as θ k = -b + √ b 2 -4ac 2a , with a = Γ (∥µ∥ ∞ ) 2 n -2Γ (∥µ∥ ∞ ) 2 Q(∥µ∥ ∞ ) + 4Q(∥µ∥ ∞ ) 2 , b = 4Γ (∥µ∥ ∞ ) 2 Q(∥µ∥ ∞ ) -8Q(∥µ∥ ∞ ) 2 + 4Q(∥µ∥ ∞ ) , c = 4Q(∥µ∥ ∞ ) 2 -4Q(∥µ∥ ∞ ) .
The proximity measure δ + φ := δ φ (z + , s + , µ(1 -θ k )) satisfies the conditions of feasibility in Theorem 3.2 and quadratic convergence of the Newton step in Theorem 3.3.

By definition Q(∥µ∥ ∞ ) < 1, thus c is negative. Furthermore, for n sufficiently large a is positive and so b 2 -4ac and θ k are positive.

We would also like to point out that the value of θ k is of order 1/ √ n for n large. This observation is fundamental considering the complexity of the algorithm.

Proof As pointed out earlier in Remark 1 the condition of Theorem 3.2 is weaker than the condition of Theorem 3.3. Thus, it is sufficient to satisfy the latter condition to ensure strict feasibility of the iterates.

According to the condition of Theorem 3.3, after an update of µ, i.e. µ + = (1 -θ k )µ, the proximity measure δ + φ must satisfy

(δ + φ ) 2 ≤ Q(∥µ + ∥ ∞ ) .
As Q is decreasing with respect to ∥µ∥ ∞ , it is sufficient to ensure that

(δ + φ ) 2 ≤ Q(∥µ∥ ∞ ) .
By Proposition 3.1 in the case δ 2 φ ≤ Q(∥µ∥ ∞ ), we can choose any θ k satisfying

(1-θ k )Q(∥µ∥ ∞ ) 2 + Γ (∥µ∥ ∞ ) 2 4(1 -θ k ) n(θ k ) 2 + (4θ k -2(θ k ) 2 )Q(∥µ∥ ∞ ) ≤ Q(∥µ∥ ∞ ) .
Therefore, it is sufficient to choose θ k > 0 such that

θ k ≤ -b + √ b 2 -4ac 2a ,
with a, b and c defined as in the statement of the theorem.

⊓ ⊔

Remark 3 For a more explicit characterization of θ k , we can study its asymptotical behaviour. By definition of Γ and Q,

lim µ→0 Q(∥µ∥ ∞ ) = 1/2 and lim µ→0 Γ (∥µ∥ ∞ ) = 1, so lim µ→0 θ k = -1+ √ 1+4n 2n ≤ 1 √ n .

Complexity Analysis of the Full Newton Step IPM

The complexity of this algorithm is obtained by the extension of a classical lemma, whose proof can be found for instance in [START_REF] Roos | Interior point methods for linear optimization[END_REF]. Lemma 3.8 Let θ be such that 0 < θ ≤ θ k ∀k ∈ N. The φ-FN-IPM for monotone LCP described in Algorithm 2 provides an ϵ-solution (z, s), which satisfies n i=1 φ(z i s i ) ≤ nϵ after at most log n i=1 φ(z 0 i s 0 i )/ϵ / θ iterations.

The sequence {θ k } is given by Theorem 3.4. As already stated, θ k is of order 1/ √ n for n large, which justify the existence of θ. Moreover, we can choose θ of order 1/ √ n without loss of generality. As a result, we can now state our main theorem. Theorem 3.5 Let µ 0 = z 0 s 0 . Algorithm 2, with the sequence of update parameters θ k described above, guarantees strict feasibility of the iterates and quadratic convergence of the proximity measure. Moreover, it provides an ϵsolution (z, s), which satisfies

n i=1 φ(z i s i ) ≤ nϵ after at most O √ n log n ϵ iterations.
Preliminary computational experiments give similar results to the classical method on a small selection of LPs from the NETLIB repository1 . This confirms the validity of our approach. Some informations regarding these results are presented on Appendix A. Further investigations on more sophisticated methods may get the best out of this new direction.

Conclusions

The method presented in this article shows a generalization of the FN-IPM with polynomial upper bound complexity for monotone LCPs considering a general family of smooth increasing concave functions. The main contributions of this article are that we extend the classical path-following method and Darvay's method with φ(t) = √ t and prove that these new methods have the best known worst-case complexity.

Further research may extend this result to a more general family of LCPs such that P * (κ)-LCPs or P 0 -LCPs with bounded level sets as in [START_REF] Asadi | Polynomial interior-point algorithm for P * (κ) horizontal linear complementarity problems[END_REF][START_REF] Asadi | A new full-Newton step O(n) infeasible interior-point algorithm for P * (κ)-horizontal linear complementarity problems[END_REF]23].

Recent developments by Roos in [START_REF] Roos | An improved and simplified full-Newton step O(n) infeasible interior-point method for linear optimization[END_REF] consider an infeasible IPM with full Newton step using only one feasibility step. More investigations regarding the method presented in this paper could extend the results in [START_REF] Roos | An improved and simplified full-Newton step O(n) infeasible interior-point method for linear optimization[END_REF].

Despite having the best worst case upper-bound for IPMs the full Newton step is not the most used approach for numerics, since this upper-bound is attained in general [START_REF] Roos | Interior point methods for linear optimization[END_REF]. We believe that the philosophy applied in this paper can be generalized to other IPMs approaches. In particular, we are planning to study a predictor-corrector implementation of this approach.

A Numerics

To validate the theoretical results, we first implemented the φ-FN-IPM described in Algorithm 2 with the sequence of update parameters given by Theorem 3.4 and with the functions φ given in Table 1. The datasets and the residuals after convergence are detailed Name Id. θ 1 log ( ) α Fct. t

t t+1
log(1 + t) (t + 10 -3 ) α -(10 -3 ) α with α ∈ (0, 1)

Table 1 φ-functions used in the computational tests in Table 2. The stopping criterion is z T s ≤ nϵ, where ϵ = 10 -6 , and we compute the residuals as ∥zs∥∞ × 10 6 . For every φ function that appears in Table 1, the method converged in the same number of iterations, so we only display the number of iterations once. This phenomenon is not surprising, since the φ-FN-IPM used here stays very close to the central path, as shown by Theorem 3.3. Therefore, this implementation does not exploit fully the new directions. Nonetheless, some small differences remain in the residuals as illustrated in To illustrate the possible differences between the functions φ, we run another experiment that considers a fixed value for θ k = 1 √ n for all k ∈ N. Table 3 illustrates the different behaviours observed for different choices of φ. We notice that a smaller number of iterations seems to be required when the derivative in zero is larger for the methods with α = 0.25, 0.5 and 0.75. Thus, further research exploring interior-point methods in large neighbourhood of the central path may get the best out of these differences.

3 µ 6 k := k + 1 ; 1 :

 3611 := (1 -θ k )µ; 4 solve system (1) to find (∆z, ∆s) ; 5 (z, s) := (z, s) + (∆z, ∆s); Algorithm Full Newton step IPM (FN-IPM)

3 µ := ( 1 -θ k )µ; 4 solve system ( 3 ) 5 ( 6 k

 314356 to find (∆z, ∆s) ; z, s) := (z, s) + (∆z, ∆s); := k + 1 ;

  23. G.Q. Wang, X.J. Fan, D.T. Zhu, and D.Z. Wang. New complexity analysis of a full-Newton step feasible interior-point algorithm for P * (κ)-LCP. Optimization Letters, 9(6):1105-1119, 2015.

Table 2 .

 2

	Instance	Iter.	Id.	θ 1	log	α = 0.25	α = 0.5	α = 0.75
	ADLITTLE	268	0.995583 0.995853 0.995856 0.995851 0.995855 0.995854
	AFIRO	170	0.943919 0.944204 0.944207 0.944200 0.944205 0.944205
	BEACONFD	418	0.995479 0.995728 0.995730 0.995724 0.995727 0.995728
	BOEING2	420	0.974064 0.974307 0.974309 0.974309 0.974309 0.974307
	BLEND	264	0.984332 0.984599 0.984602 0.984600 0.984599 0.984599
	GROW7	406	0.983538 0.983785 0.983788 0.983780 0.983787 0.983785
	ISRAEL	428	0.988911 0.989157 0.989155 0.989148 0.989154 0.989156
	KB2	203	0.944257 0.944530 0.944532 0.944530 0.944531 0.944530
	RECIPELP	332	0.960634 0.960884 0.960887 0.960884 0.960885 0.960885
	SC50A	218	0.951084 0.951353 0.951356 0.951350 0.951354 0.951354
	SC50B	218	0.951084 0.951353 0.951356 0.951350 0.951354 0.951354
	SC105	316	0.973217 0.973472 0.973475 0.973473 0.973474 0.973473
	SCAGR7	342	0.985376 0.985631 0.985633 0.985632 0.985632 0.985632
	SHARE1B	372	0.968387 0.968633 0.968635 0.968634 0.968632 0.968634
	SHARE2B	310	0.972849 0.973105 0.973107 0.973102 0.973107 0.973106
	STOCFOR1	324	0.982774 0.983031 0.983033 0.983033 0.983031 0.983031

Table 2

 2 Value of Res.=∥zs∥∞ × 10 6 , after the algorithm reaches z T s ≤ nϵ = 10 -6 n

http://www.netlib.org/lp/

Instance

Id. θ 1 log α = 0.25 α = 0.


Jean-Pierre • Dussault
email: jean-pierre.dussault@usherbrooke.ca

Mounir • Haddou
email: mounir.haddou@insa-rennes.fr

Abdeslam • Kadrani
email: kadrania@gmail.com

Migot

Tangi

Mounir Haddou

Tangi Migot
email: tangi.migot@gmail.com

On Approximate Stationary Points of the Regularized Mathematical Program with Complementarity Constraints

Keywords: non-linear programming, MPCC, regularization methods, constraint qualification, optimization model with complementarity constraints Mathematics Subject Classification (2000) 90C30, 90C33, 49M37, 65K05

We discuss the convergence of regularization methods for the mathematical program with complementarity constraints with approximate sequence of stationary points. It is now well accepted in the literature that, under some tailored constraint qualification, the genuine necessary optimality condition for this problem is the M-stationarity condition. It has been pointed out, in Kanzow, C., and Schwartz, A., Mathematics of Operations Research (2015), that relaxation methods with approximate stationary points fail to ensure convergence to M-stationary points. We define a new strong approximate stationarity concept and we prove that a sequence of strong approximate stationary points always converges to an M-stationary solution. We also prove under weak assumptions the existence of strong approximate stationary points in the neighborhood of an M-stationary solution.

Introduction

The Mathematical Program with Complementarity Constraints (MPCC) is a non-linear optimization problem including a complementarity problem in the constraints. It plays a very important role in many fields, such as engineering design, economic equilibria, transportation science, multilevel game, and mathematical programming itself. However, this kind of problems is generally difficult to deal with because it fails to satisfy the standard Mangasarian-Fromovitz constraint qualification (MFCQ) at any feasible point [START_REF] Ye | Exact penalization and necessary optimality conditions for generalized bilevel programming problems[END_REF]. See the monographs [START_REF] Luo | Mathematical programs with equilibrium constraints[END_REF][START_REF] Outrata | Nonsmooth approach to optimization problems with equilibrium constraints: theory, applications and numerical results[END_REF] for details about the basic theory, effective algorithms, and various applications of the MPCC.

During the past two decades, many researchers introduced necessary optimality conditions such as the Clarke (C-), Mordukhovich (M-), strong (S-), and Bouligand (B-) stationarity conditions for the MPCC; see, e.g., [START_REF] Ye | Exact penalization and necessary optimality conditions for generalized bilevel programming problems[END_REF][START_REF] Flegel | A direct proof for M-stationarity under MPEC-GCQ for mathematical programs with equilibrium constraints[END_REF][START_REF] Kanzow | Mathematical programs with equilibrium constraints: enhanced fritz john-conditions, new constraint qualifications, and improved exact penalty results[END_REF][START_REF] Scheel | Mathematical programs with complementarity constraints: Stationarity, optimality, and sensitivity[END_REF][START_REF] Ye | Constraint qualifications and necessary optimality conditions for optimization problems with variational inequality constraints[END_REF][START_REF] Ye | Necessary and sufficient optimality conditions for mathematical programs with equilibrium constraints[END_REF][START_REF] Ye | Necessary optimality conditions for optimization problems with variational inequality constraints[END_REF]. Among these stationarities, the B-stationarity is known to be a good candidate for optimality, but since it is computationally difficult, it is rarely used in algorithmic analysis; the S-stationarity is the strongest and equivalent to the KKT conditions (see, e.g., [START_REF] Flegel | On the guignard constraint qualification for mathematical programs with equilibrium constraints[END_REF][START_REF] Guo | Solving mathematical programs with equilibrium constraints[END_REF]), but its interest is reduced since it does not always hold for the MPCC. The M-stationarity, which has already widely been investigated (see, e.g., [START_REF] Ye | Exact penalization and necessary optimality conditions for generalized bilevel programming problems[END_REF][START_REF] Flegel | A direct proof for M-stationarity under MPEC-GCQ for mathematical programs with equilibrium constraints[END_REF][START_REF] Kanzow | Mathematical programs with equilibrium constraints: enhanced fritz john-conditions, new constraint qualifications, and improved exact penalty results[END_REF][START_REF] Ye | Constraint qualifications and necessary optimality conditions for optimization problems with variational inequality constraints[END_REF][START_REF] Ye | Necessary and sufficient optimality conditions for mathematical programs with equilibrium constraints[END_REF][START_REF] Ye | Necessary optimality conditions for optimization problems with variational inequality constraints[END_REF][START_REF] Outrata | Optimality conditions for a class of mathematical programs with equilibrium constraints[END_REF]), is the most relevant concept since it is the weakest necessary condition holding, under suitable constraint qualifications, at any local minimizer of the MPCC and is computationaly tractable.

The feasible set of the MPCC involves a complementarity constraint equivalent to G(x)=0 OR H(x) = 0. This is a thin set exhibiting some irregularity when G(x)=0 AND H(x)=0. It is this thinness that makes constraint qualifications fail at any feasible point. In view of the constraint qualifications issues that plague the MPCC the regularization methods provide an intuitive answer. The complementarity constraint is regularized using a parameter so that the new feasible domain is not thin anymore. Regularization methods drive this parameter to zero so that the feasible set of the regularized problem converges to the feasible set of (3) similar to an homotopy technique; see [START_REF] Scholtes | Convergence properties of a regularization scheme for mathematical programs with complementarity constraints[END_REF][START_REF] Demiguel | A two-sided relaxation scheme for mathematical programs with equilibrium constraints[END_REF][START_REF] Lin | A modified relaxation scheme for mathematical programs with complementarity constraints[END_REF][START_REF] Steffensen | A new relaxation scheme for mathematical programs with equilibrium constraints[END_REF][START_REF] Kadrani | A new regularization scheme for mathematical programs with complementarity constraints[END_REF][START_REF] Kanzow | A new regularization method for mathematical programs with complementarity constraints with strong convergence properties[END_REF]. In the corresponding papers it has been shown that, under suitable conditions providing convergence of the methods, they still might converge to some spurious points, called C-stationary points, the convergence to M-stationary points being guaranteed only under some additional second-order conditions, strict complementarity-type condition or conditions on the sequence of iterates that cannot be guaranteed a priori. The only exception are the methods KDB and KS proposed respectively in [START_REF] Kadrani | A new regularization scheme for mathematical programs with complementarity constraints[END_REF] and [START_REF] Kanzow | A new regularization method for mathematical programs with complementarity constraints with strong convergence properties[END_REF]. In those references, the authors prove convergence of their regularization scheme to an M-stationary point under a classical constraint qualification. However, in [START_REF] Kanzow | The price of inexactness: Convergence properties of relaxation methods for mathematical programs with complementarity constraints revisited[END_REF], Kanzow and Schwartz discuss convergence of the methods considering a sequence of approximate stationary points, that is points that satisfy approximately the KKT conditions. They illustrate the fact that, for such sequences of approximate stationary points of the regularized subproblems, the method may converge to spurious weak-stationary points. Our motivation is to deal with this issue.

We focus in this paper on the KDB regularization method introduced in [START_REF] Kadrani | A new regularization scheme for mathematical programs with complementarity constraints[END_REF]. By defining a new kind of approximate stationary point called strong epsilon-stationary point we prove that any cluster point of a sequence of strong epsilon-stationary points of the KDB regularization is an M-stationary points. Moreover, we analyze existence of the strong epsilon-stationary points in the neighborhood of a solution. Previous studies of regularization methods require strong assumption to guarantee the existence of stationary point [START_REF] Kadrani | A new regularization scheme for mathematical programs with complementarity constraints[END_REF][START_REF] Kanzow | A new regularization method for mathematical programs with complementarity constraints with strong convergence properties[END_REF]. We show that strong epsilon-stationary points have a good behavior with this respect by proving the existence in the neighborhood of any M-stationary point without the need of constraint qualification.

The results presented in this paper are theoretical and present the advantages of strong epsilon-stationary points. In particular, we cover two important questions in the study of regularization methods for the MPCC namely the global convergence of strongly epsilon-stationary points to M-stationary solutions and weak assumptions allowing to prove the existence of such approximate stationary points. This study paves the way to further studies on an algorithm able to compute strong epsilon-stationary points.

The rest of the paper is organized as follows. In the next section, we sum up some stationarity concepts and constraint qualifications for the standard non-linear program and the MPCC. In Section 3, we recall the regularization method from [START_REF] Kadrani | A new regularization scheme for mathematical programs with complementarity constraints[END_REF] and motivate the difficulty arising when dealing with approximate stationary points of the regularized subproblems. We introduce our new strong epsilon-stationary point and prove convergence of the method to an M-stationary point in Section 4. In Section 5, after introducing the regularization method with slack variables, we prove our main existence result of strong epsilon-stationary point. Finally, we present some perspectives of this work in Section 6.

Preliminaries

Non-Linear Programming

Let a general non-linear program be

min x∈R n f (x) s.t. g(x) ≤ 0, h(x) = 0, (1)
with g : R n → R m , h : R n → R p and f : R n → R. Denote F the feasible region of (1), the set of active indices I g (x) := {i ∈ {1, ..., m} : g i (x) = 0} at x and the Lagrangian L(x, λ) := f (x) + g(x) T λ g + h(x) T λ h with λ = (λ g , λ h ). A fundamental tool to address such problems is a necessary optimality condition, the best known being the Karush-Kuhn-Tucker (KKT) condition. Let x * ∈ F be a local minimum of (1) that satisfies a constraint qualification, then there exists λ ∈ R m × R p such that (x * , λ) satisfy the equation ∇ x L(x * , λ) = 0, with the conditions λ g ≥ 0, supp(λ g ) := {i : λ g i = 0} ⊆ I g (x * ). A couple (x * , λ) satisfying the above conditions is called a KKT point; it is a stationary point of the Lagrangian function, its gradient vanishes. The stationary condition is often expressed as follows.

∇ x L(x * , λ) = 0, g i (x * ) ≤ 0, λ g i ≥ 0, λ g i g i (x * ) = 0, ∀i ∈ {1, . . . , m}, h i (x * ) = 0, ∀i ∈ {1, . . . , p}. (2)
Algorithms to solve (1) aim at computing a stationary point. In a practical context, it can be difficult, even impossible, to compute such stationary points, therefore approximate, or epsilon-stationary points are usually considered. We refer the reader for instance to [START_REF] Kanzow | The price of inexactness: Convergence properties of relaxation methods for mathematical programs with complementarity constraints revisited[END_REF] or [START_REF] Andreani | On sequential optimality conditions for smooth constrained optimization[END_REF] for a deeper motivation on this condition.

Definition 2.1 Let ≥ 0. We say that x * is an -stationary point of (1), if there exists λ ∈ R m × R p satisfying

∇ x L(x * , λ) ∞ ≤ , g i (x * ) ≤ , λ g i ≥ 0, |λ g i g i (x * )| ≤ , ∀i ∈ {1, . . . , m}, |h i (x * )| ≤ , ∀i ∈ {1, . . . , p}.
At = 0 we get the usual definition (2).

The Mathematical Program with Complementarity Constraints

We consider the Mathematical Program with Complementarity Constraints (MPCC) min

x∈R n f (x) s.t. g(x) ≤ 0, h(x) = 0, 0 ≤ G(x) ⊥ H(x) ≥ 0, (3)
with f : R n → R, g : R n → R m , h : R n → R p and G, H : R n → R q . All these functions are assumed to be continuously differentiable through this paper.

The notation 0 ≤ u ⊥ v ≥ 0 for two vectors u and v in R q is a shortcut for u ≥ 0, v ≥ 0 and u T v = 0. Let Z be the set of feasible points of (3). Given x ∈ Z, we denote I +0 (x) := {i ∈ {1, . . . , q} : G i (x) > 0 and H i (x) = 0},

I 0+ (x) := {i ∈ {1, . . . , q} : G i (x) = 0 and H i (x) > 0},
I 00 (x) := {i ∈ {1, . . . , q} : G i (x) = 0 and H i (x) = 0}.

We define the MPCC-Lagrangian function of (3) as

L M P CC (x, λ) := f (x) + g(x) T λ g + h(x) T λ h -G(x) T λ G -H(x) T λ H , with λ := (λ g , λ h , λ G , λ H ).
We introduce more stationary concepts as in [START_REF] Ye | Necessary and sufficient optimality conditions for mathematical programs with equilibrium constraints[END_REF][START_REF] Ye | Necessary optimality conditions for optimization problems with variational inequality constraints[END_REF][START_REF] Ye | Constraint qualifications and necessary optimality conditions for optimization problems with variational inequality constraints[END_REF][START_REF] Outrata | Optimality conditions for a class of mathematical programs with equilibrium constraints[END_REF][START_REF] Scheel | Mathematical programs with complementarity constraints: Stationarity, optimality, and sensitivity[END_REF][START_REF] Outrata | A generalized mathematical program with equilibrium constraints[END_REF][START_REF] Luo | Mathematical programs with equilibrium constraints[END_REF]. Those concepts are needed for two reasons:

unless assuming a restrictive constraint qualification, a local minimizer x * may fail to be a KKT point, so that (2) needs to be weakened in order to obtain a necessary condition;

when analyzing cluster points of algorithms, other weak stationarity conditions appear naturally.

Definition 2.2 x * ∈ Z is said -W-stationary if there exists λ ∈ R m × R p × R q × R q such that ∇ x L M P CC (x * , λ) = 0, λ g ≥ 0, λ g i = 0, ∀i / ∈ I g (x * ), λ G i = 0, ∀i ∈ I +0 (x * ), and, λ H i = 0, ∀i ∈ I 0+ (x * ); -C-stationary, if it is W-stationary and λ G i λ H i ≥ 0 for all i ∈ I 00 (x * ); -M-stationary, if it is W-stationary and either λ G i > 0, λ H i > 0 or λ G i λ H i = 0 for all i ∈ I 00 (x * ); -S-stationary, if it is W-stationary and λ G i ≥ 0, λ H i ≥ 0 for all i ∈ I 00 (x * ).
As pointed out in [START_REF] Flegel | On the guignard constraint qualification for mathematical programs with equilibrium constraints[END_REF], strong stationarity is equivalent to the standard KKT conditions of an MPCC. In order to guarantee that a local minimum x * of (1) is a stationary point in any of the previous senses, one needs to assume that a constraint qualification (CQ) is satisfied in x * . Since most standard CQs are violated at any feasible points of (3), many MPCC-analogues of these CQs have been developed. Here, we mention only those needed later.

Definition 2.3 A feasible point x * of (3), i.e. x * ∈ Z, is said to satisfy 1. MPCC-linear independence CQ (MPCC-LICQ), if the gradients ∇g i (x * ) (i ∈ I g (x * )), ∇h i (x * ) (i = 1, . . . , p), ∇G i (x * ) (i ∈ I 00 (x * ) ∪ I 0+ (x * )), ∇H i (x * ) (i ∈ I 00 (x * ) ∪ I +0 (x * )),
are linearly independent.

MPCC-Mangasarian Fromovitz CQ (MPCC-MFCQ), if the only solution

to the equation

i∈Ig λ g i ∇g i (x * ) + p i=1 λ h i ∇h i (x * ) + i∈I 00 (x * )∪I 0+ (x * ) λ G i ∇G i (x * ) + i∈I 00 (x * )∪I +0 (x * ) λ H i ∇H i (x * ) = 0 with λ g i ≥ 0 is the trivial solution, i.e. λ g = λ h = λ G = λ H = 0.
Note that each of these CQs implies that a local minimum is M-stationary, see [START_REF] Flegel | Abadie-type constraint qualification for mathematical programs with equilibrium constraints[END_REF][START_REF] Ye | Necessary and sufficient optimality conditions for mathematical programs with equilibrium constraints[END_REF], but only MPCC-LICQ is sufficient to guarantee strong stationarity of a local minimum, see [START_REF] Flegel | On the guignard constraint qualification for mathematical programs with equilibrium constraints[END_REF][START_REF] Luo | Mathematical programs with equilibrium constraints[END_REF][START_REF] Pang | Complementarity constraint qualifications and simplified bstationarity conditions for mathematical programs with equilibrium constraints[END_REF]. The MPCC-LICQ is among the first MPCCtailored constraint qualifications and may already be found in [START_REF] Luo | Mathematical programs with equilibrium constraints[END_REF] and [START_REF] Scheel | Mathematical programs with complementarity constraints: Stationarity, optimality, and sensitivity[END_REF], the MPCC-MFCQ was introduced in [START_REF] Scheel | Mathematical programs with complementarity constraints: Stationarity, optimality, and sensitivity[END_REF] and presented in the form above in [START_REF] Hoheisel | Theoretical and numerical comparison of relaxation methods for mathematical programs with complementarity constraints[END_REF].

The following example due to Scheel and Scholtes [START_REF] Scheel | Mathematical programs with complementarity constraints: Stationarity, optimality, and sensitivity[END_REF] exhibits a situation where the global minimizer is not an S-stationary point but an M-stationary point. We will return to this example later on.

Example 2.1

min x∈R 3 x 1 + x 2 -x 3 s.t. -4x 1 + x 3 ≤ 0, -4x 2 + x 3 ≤ 0, 0 ≤ x 1 ⊥ x 2 ≥ 0.
The global solution is (0, 0, 0) T but is not an S-stationary point. Indeed, the gradient of the Lagrangian equal to zero yields

0 =   1 1 -1   + λ g 1   -4 0 1   + λ g 2   0 -4 1   -λ G   1 0 0   -λ H   0 1 0   ,
and since λ g 1 + λ g 2 = 1(third line), summing the first two lines yields 2 -4(λ g 1 + λ g

2 ) -λ G -λ H = 0 and therefore λ G + λ H = -2; both cannot be non-negative.

All in all, unless MPCC-LICQ is assumed, devising algorithms to provably reach stationary points (equivalent to S-stationary point) is not possible in general, and we must satisfy ourselves in devising algorithms reaching Mstationary points.

The KDB Regularization Method and Its Convergence

Any complementarity constraint represents a thin geometric locus troublesome with respect to CQ. A regularization introduces a parametric way to thicken these thin constraints and hopefully the regularized subproblem will satisfy a CQ.

Given a map Φ : R q × R q × R → R q regularization methods consider the following parametric non-linear program R t, t parametrized by two positive parameters t and t:

min x∈R n f (x) s.t. g(x) ≤ 0, h(x) = 0, G(x) ≥ -t1l, H(x) ≥ -t1l, Φ(G(x), H(x); t) ≤ 0. (4)
In [START_REF] Kadrani | A new regularization scheme for mathematical programs with complementarity constraints[END_REF], Kadrani, Dussault and Benchakroun introduced a method which enjoys the desired goal to converge to an M-stationary point without additional assumptions such as second-order conditions. Their method considers an approximation of the complementarity constraints defined for all i ∈ {1, . . . , q}

as Φ i (G(x), H(x); t) := (G i (x) -t)(H i (x) -t). (5)
The only difference here compared to the initial presentation of the method in [START_REF] Kadrani | A new regularization scheme for mathematical programs with complementarity constraints[END_REF] is that t might be different from t, i.e. thicken asymetrically the complementarity constraints. In the sequel we use the following notations. The Lagrangian function of ( 4) is defined as

L R t,t (x, ν) := f (x) + g(x) T ν g + h(x) T ν h -(G(x) + t1l) T ν G -(H(x) + t1l) T ν H + Φ(G(x), H(x); t) T ν Φ , with ν := (ν g , ν h , ν G , ν H , ν Φ ) ∈ R m × R p × R q × R q × R q .
Denote by I Φ (x; t) := {i : Φ i (G(x), H(x); t) = 0} the set of active indices of the regularization constraint.

The method has been shown to converge to an M-stationary point as recalled in the following theorem from [START_REF] Hoheisel | Theoretical and numerical comparison of relaxation methods for mathematical programs with complementarity constraints[END_REF] (the original convergence result from [START_REF] Kadrani | A new regularization scheme for mathematical programs with complementarity constraints[END_REF] being stated with MPCC-LICQ). The extension of [START_REF] Hoheisel | Theoretical and numerical comparison of relaxation methods for mathematical programs with complementarity constraints[END_REF]Theorem 3.5] to the case t = t is straightforward. Theorem 3.1 Let {t k } ↓ 0, { tk } ↓ 0, and {x k } be a sequence of stationary points of (4). If x k → x * and MPCC -MFCQ holds in x * , then x * is an M-stationary point of (3). This result was a clear improvement over other existing methods such as [START_REF] Scholtes | Convergence properties of a regularization scheme for mathematical programs with complementarity constraints[END_REF][START_REF] Demiguel | A two-sided relaxation scheme for mathematical programs with equilibrium constraints[END_REF][START_REF] Lin | A modified relaxation scheme for mathematical programs with complementarity constraints[END_REF][START_REF] Steffensen | A new relaxation scheme for mathematical programs with equilibrium constraints[END_REF], which only converge to a C-stationary point with similar assumptions. An additional necessary second order condition is often required to ensure convergence to an M-stationary point. Only the method developed in [START_REF] Kanzow | A new regularization method for mathematical programs with complementarity constraints with strong convergence properties[END_REF], extending the KDB method as a regularization, shares the same strong convergence property.

However, as shown in [START_REF] Kanzow | The price of inexactness: Convergence properties of relaxation methods for mathematical programs with complementarity constraints revisited[END_REF], the situation is different when considering sequence of approximate stationary points. We recall here the corresponding result, where the extension to the case t = t is again straightforward.

) Let {t k } ↓ 0, { tk } ↓ 0,{ k } ↓ 0, {x k } be a sequence of k -stationary points of (4). If x k → x * and MPCC -MFCQ holds in x * , then x * is a W-stationary point of (3).
Suppose further that k = o( tk ) and there is a constant c > 0 such that, for all i ∈ I 00 (x * ) and all k sufficiently large, the iterates

(G i (x k ), H i (x k )) satisfy (G i (x k ), H i (x k )) / ∈ ](1 -c)t k , (1 + c)t k [ 2 \ ](1 -c)t k , t k [ 2 . (6)
Then x * is an M-stationary point.

Note that the authors also gave similar type of conditions to get convergence to a C-or an S-stationary point. The drawback of condition (6), however, is that it cannot be checked a priori and the strong convergence result of the method no longer holds.

In [START_REF] Kanzow | A new regularization method for mathematical programs with complementarity constraints with strong convergence properties[END_REF]Example 3], the authors provide an example satisfying MPCC-LICQ and exhibiting convergence to a W-stationary point. One of the main reason for this behavior being that if we perturb the relation |ν Φ i Φ i (x 1 , x 2 ; t)| ≤ (leaving the other conditions ν Φ ≥ 0, Φ(x 1 , x 2 ; t) ≤ 0), ν Φ may be positive when the constraint Φ(x 1 , x 2 ; t) is not active. Now, if the complementarity constraint is satisfied approximately, but an approximate complementarity condition of the form ν Φ i (Φ i (x 1 , x 2 ; t) -) = 0 is guaranteed, convergence may occur to C-stationary points as shown in the following example.

Example 3.1 Consider the problem

min x∈R 2 1 2 ((x 1 -1) 2 + (x 2 -1) 2 ) s.t. 0 ≤ x 1 ⊥ x 2 ≥ 0.
We specialize the relations in Definition 2.1 for (4) as the following

x 1 -1 x 2 -1 -ν G 1 0 -ν H 0 1 + ν Φ x 2 -t x 1 -t ∞ ≤ , 0 ≤ ν G , (x 1 + t) ≥ 0, ν G (x 1 + t) ≤ , 0 ≤ ν H , (x 2 + t) ≥ 0, ν H (x 2 + t) ≤ , 0 ≤ ν Φ , (x 1 -t)(x 2 -t) ≤ , ν Φ [(x 1 -t)(x 2 -t) -] = 0.
The points

(t + √ , t + √ ) T together with ν G = ν H = 0 and ν Φ = 1-t- √ √
satisfy the above relations. The limit point when t, -→ 0 is the origin, which is a C-stationary point with λ G = λ H = -1.

Epsilon-stationary points have two main drawbacks considering our goal of computing an M-stationary point. The difficulties may come from the approximation of the complementarity condition as shown by [START_REF] Kanzow | A new regularization method for mathematical programs with complementarity constraints with strong convergence properties[END_REF]Example 3] or from the approximation of the feasibility of the regularized constraint as illustrated in Example 3.1.

We focus in this paper on studying approximate stationary points. Note that an alternative but complementary approach would be to study sequential optimality conditions. This approach first designed for nonlinear programming [START_REF] Andreani | On sequential optimality conditions for smooth constrained optimization[END_REF] has recently been extended to the MPCC in [START_REF] Andreani | New sequential optimality conditions for mathematical programs with complementarity constraints and algorithmic consequences[END_REF] or in [START_REF] Ramos | Mathematical programs with equilibrium constraints: a sequential optimality condition, new constraint qualifications and algorithmic consequences[END_REF].

Convergence of Strong Epsilon-Stationary Sequences

We now address the convergence of sequences of epsilon-stationary points. The previous section motivates the definition of a new kind of epsilon-stationary point called strong epsilon-stationary point, which is more stringent regarding the complementarity constraint. Definition 4.1 Let ≥ 0. We say that x is a strong -stationary point for (4) if there exists

ν ∈ R m × R p × R q × R q × R q such that ∇L R t,t (x, ν) ∞ ≤ , g i (x) ≤ , ν g i ≥ 0, |g i (x)ν g i | ≤ , ∀i ∈ {1, . . . , m}, |h i (x)| ≤ t + 3 , ∀i ∈ {1, . . . , p}, G i (x) + t ≥ -, ν G i ≥ 0, ν G i (G i (x) + t) ≤ , ∀i ∈ {1, . . . , q}, H i (x) + t ≥ -, ν H i ≥ 0, ν H i (H i (x) + t) ≤ , ∀i ∈ {1, . . . , q}, Φ i (G(x), H(x); t) ≤ 0, ν Φ i ≥ 0, ν Φ i Φ i (G(x)
, H(x); t) = 0, ∀i ∈ {1, . . . , q}. Note that we are less demanding on the equality constraint. This will become clearer when considering slack variables to the MPCC to get our result on existence of strong epsilon-stationary point in Section 5.

Since a strong epsilon-stationary point is, in particular, an epsilon-stationary point, we can deduce from the first part of Theorem 3.2 the following result.

Lemma 4.1 Let {t k } ↓ 0, { tk } ↓ 0,{ k } ↓ 0, {x k } be a sequence of strong k -stationary points of (4). If x k → x * and MPCC -MFCQ holds in x * , then x * is a W-stationary point of (3).
We can now move to our main convergence result. Theorem 4.1 Let {t k } ↓ 0, { tk } ↓ 0, k = o( tk ), and {x k } be a sequence of strong k -stationary points of (4). If x k → x * and MPCC -MFCQ holds in x * , then x * is an M-stationary point of (3).

Proof The proof of [START_REF] Kanzow | The price of inexactness: Convergence properties of relaxation methods for mathematical programs with complementarity constraints revisited[END_REF]Theorem 9] shows that the sequence {η g,k , η h,k , η G,k , η H,k } is bounded with η g,k := ν g,k , η h,k := ν h,k and for all i

η G,k i := ν G,k i -ν Φ,k i (H i (x k ) -t k ), η H,k i := ν H,k i -ν Φ,k i (G i (x k ) -t k ). (7)
Up to a subsequence, the limit point of this sequence is a mutliplier, denoted (η g, * , η h, * , η G, * , η H, * ), such that x * is a W-stationary point. We now consider indices i ∈ I 00 (x * ). Our aim here is to prove that either η

G, * i > 0, η H, * i > 0 or η G, * i η H, * i = 0 by showing that η G, * i < 0 =⇒ η H, * i = 0. Assume that η G, * i < 0, then for k sufficiently large it holds that η G,k i < 0 ⇐⇒ ν G,k i -ν Φ,k i (H i (x k ) -t k ) < 0 =⇒ ν Φ,k i > 0 and H i (x k ) > t k , since ν G,k and ν Φ,k are non-negative. Using the complementarity condition ν Φ,k i (G i (x k ) -t k )(H i (x k ) -t k ) = 0 it follows that G i (x k ) = t k . Moreover, using that H i (x k ) > t k , we obtain k ≥ ν H,k i (H i (x k ) + tk ) ≥ ν H,k i tk , therefore ν H,k i ↓ 0 as k = o( tk ). Finally, since G i (x k ) = t k and ν H,k i ↓ 0 it holds that η H, * i = 0.
Further work may focus on weakening the constraint qualification necessary to prove this convergence theorem. In the literature, MPCC-RCPLD [START_REF] Andreani | Convergence properties of a second order augmented lagrangian method for mathematical programs with complementarity constraints[END_REF] or MPCC-CCP [START_REF] Ramos | Mathematical programs with equilibrium constraints: a sequential optimality condition, new constraint qualifications and algorithmic consequences[END_REF] have been used to prove convergence results involving Mstationary points.

Existence of Strong Epsilon-Stationary Points for the Relaxation with Slack Variables

On Lagrange Multipliers of the Relaxation

Beyond the convergence of a sequence of (approximate) stationary points of the regularization, another important question is wether such points exist in a neighborhood of an M-stationary point.

Following Example 2.1, we can encounter a situation where the regularized subproblems have no stationary point.

Example 5.1 The KDB regularized problem for Example 2.1 is

min x∈R 3 x 1 + x 2 -x 3 s.t. -4x 1 + x 3 ≤ 0, -4x 2 + x 3 ≤ 0, x 1 ≥ -t, x 2 ≥ -t, (x 1 -t)(x 2 -t) ≤ 0.
The point (t, t, 4t) T is feasible so that the minimum value of this program is ≤ -2t. Moreover, whenever x 1 > t, we must have x 2 ≤ t to satisfy (x 1t)(x 2 -t) ≤ 0. This allows to conclude that (t, t, 4t) T is the global minimum of the regularized problem. ν G = ν H = 0 and the gradient of the Lagrangian equal to zero yields

0 =   1 1 -1   + ν g 1   -4 0 1   + ν g 2   0 -4 1   + ν Φ   0 0 0   ,
which cannot be satisfied.

This last example seems to contradict Theorem 4.6 in [START_REF] Kanzow | A new regularization method for mathematical programs with complementarity constraints with strong convergence properties[END_REF], but MPCC-LICQ is not satisfied at the origin. It has been pointed out earlier that a practical algorithm may not be able to compute stationary points of the regularized subproblem, but only some approximate stationary point. In the previous section, in particular, we underline the interest for strong epsilon-stationary points. An intuitive idea would be that a weaker constraint qualification (like MPCC-MFCQ) would guarantee existence of such points in the neighborhood of an M-stationary point. However, the following one-dimensional example shows that things are not that simple. The constraint set involves only linear functions, and the objective is linear. There is only one feasible point, x * = 0. Moreover, it is not hard to verify that this point is an M-stationary point.

Given 1 > > 0, let us now verify that there is no strong -stationary point for this problem. The condition on the gradient of the Lagrangian yields

| -1 -ν G -ν H + 2ν Φ (x -t)| ≤ , (8)
since in this case ∇Φ(G(x), H(x); t) = 2(x -t). Since x must be feasible, we get (x -t) 2 ≤ 0, therefore x = t. Hence, there is no x satisfying (8) as ν G and ν H are non-negative.

The previous example shows that existence of a strong epsilon-stationary point is not ensured in the neighborhood of an M-stationary point even so a constraint qualification holds. However, by reformulating (3) with slack variables things could be slightly different as shown by the continuation of the example.

Example 5.3 (Example 5.2 continued) Given > 0 and t ≤ t. We now verify that there is a strong epsilon-stationary point of the approximation KDB written with slack variables for this example,

min x∈R -x s.t. s G = x, s H = x , s G ≥ -t, s H ≥ -t, (s G -t)(s H -t) ≤ 0.
Let δ be a constant satisfying 2 ≥ δ > 0. Consider the point x = 0 and (s G , s H ) = (t, t+δ) with Lagrange multiplier (ν s G , ν s H , ν G , ν H , ν Φ ) = (-1, 0, 0, 0, 1 δ ). Let us verify that this point is a strong epsilon-stationary point of the regularization:

-Condition on the gradient of the Lagrangian:

|∇ x L(x, s, ν)| = | -1 -ν s G -ν s H | = 0, |∇ s G L(x, s, ν)| = |ν s G -ν G + ν Φ (s H -t)| = 0, |∇ s H L(x, s, ν)| = |ν s H -ν H + ν Φ (s G -t)| = 0.
-Condition on the feasibility:

|x -s G | = t ≤ t + 2 , |x -s H | = t + δ ≤ t + 2 , s G + t = t + t ≥ -, s H + t = t + t + δ ≥ -, and (s G -t)(s H -t) = 0. -Condition on the complementarity: |(s G + t)ν G | = 0, |(s H + t)ν H | = 0, and |(s G -t)(s H -t)ν Φ | = 0.
This example shows that the formulation with slack variables has a strong epsilon-stationary point.

The Regularization with Slack Variables

Consider the following MPCC reformulated with slacks variables:

min (x,s)∈R n ×R 2q f (x) s.t. g(x) ≤ 0, h(x) = 0, s G = G(x), s H = H(x), 0 ≤ s G ⊥ s H ≥ 0. ( 9
)
Consider now the following parametric non-linear program with t, t ≥ 0:

min (x,s)∈R n ×R 2q f (x) s.t. g(x) ≤ 0, h(x) = 0, s G = G(x), s H = H(x), s G ≥ -t1l, s H ≥ -t1l, Φ(s G , s H ; t) ≤ 0. ( 10
)
The Lagrangian function of ( 10) is defined as

L R s t,t (x, s, ν) := f (x)+g(x) T ν g +h(x) T ν h -(G(x)-s G ) T ν s G -(H(x)-s H ) T ν s H -(s G + t1l) T ν G -(s H + t1l) T ν H + Φ(s G , s H ; t) T ν Φ .
The following result is a direct corollary of Theorem 4.1 stating that the reformulation with slack variables does not alter the convergence result.

Corollary 5.1 Let {t k } ↓ 0, { tk } ↓ 0, k = o( tk ), {x k , s k } be a sequence of strong epsilon-stationary points of [START_REF] Flegel | On the guignard constraint qualification for mathematical programs with equilibrium constraints[END_REF]. If x k → x * and MPCC -MFCQ holds at x * , then x * is an M-stationary point of (3).

Proof Let h(x) : R n → R p ×R q ×R q be such that h(x) := (h(x), s G -G(x), s H -H(x)) and x := (x, s G , s H ). It is clear that the non-linear program [START_REF] Flegel | On the guignard constraint qualification for mathematical programs with equilibrium constraints[END_REF] fall under the formulation (4). Moreover, MPCC-MFCQ for ( 9) is equivalent to MPCC-MFCQ applied to (3) with slack variables. Therefore, we can apply Theorem 4.1 to conclude this proof.

The following lemma giving an explicit form of the gradient of the Lagrangian function of (10) can be deduced through direct computations.

Lemma 5.1 The gradient of L R s t,t (x, s, ν) is given by

∇ x L R s t,t (x, s, ν) = ∇L M P CC (x, ν g , ν h , ν s G , ν s H ), ∇ s G L R s t,t (x, s, ν) = ν s G -ν G + ∇ s G Φ(s G , s H ; t) T ν Φ , ∇ s H L R s t,t (x, s, ν) = ν s H -ν H + ∇ s H Φ(s G , s H ; t) T ν Φ .
There are two direct consequences of this result. First, it is easy to see from this lemma that computing a stationary point of L R s t,t (x, s, ν) is equivalent to computing a stationary point of L R t,t (x, ν). Secondly, a stationary point of L R s t,t (x, s, ν) satisfies ∇L M P CC (x, ν g , ν h , ν s G , ν s H ) = 0, which is one of the conditions of W-stationarity for (3).

Existence

The following lemma proves that given any M-stationary point x * there is a neighborhood such that for all x in this neighborhood there exists s such that (x, s) is a strong epsilon-stationary point for the regularized problem. Lemma 5.2 Assume that f, g, h, G, H have Lipschitz continuous gradients. Let x * be an M-stationary point of (3) and > 0 be small enough. Then, there exist positive constants c, δ, t * ≥ c , t * such that, for all t ∈]c , t * [ and t ∈]0, min(t * , t + c )[, it holds that ∀x ∈ B δ (x * ), ∃ŝ ∈ B ((G(x * ), H(x * )) T ), (x, ŝ) is a strong epsilon-stationary point of [START_REF] Flegel | On the guignard constraint qualification for mathematical programs with equilibrium constraints[END_REF] .

Proof The proof is conducted in two steps. First, taking any x ∈ B δ (x * ) we construct a point that is a candidate to be strong epsilon-stationary. Then, we verify that this candidate is actually a strong epsilon-stationary point.

Let x * be an M-stationary point. Therefore, x * is feasible for (3) and there exists λ = (λ g , λ h , λ G , λ H ) such that:

∇ x L M P CC (x * , λ) = 0, λ g ≥ 0, λ g i = 0 (i / ∈ I g (x * )), λ G i = 0 (i ∈ I +0 (x * )), λ H i = 0 (i ∈ I 0+ (x * )), either λ G i > 0, λ H i > 0 or λ G i λ H i = 0 (i ∈ I 00 (x * )).
Let c be the positive constant defined as follows

c := sup i∈supp(λ G ),j∈supp(λ H ),k∈supp(λ g ) 1 |λ G i | + 1 |λ H j | + 1 |λ g k | , with the convention that c := 0 if supp(λ G ) ∪ supp(λ H ) ∪ supp(λ g ) = ∅.
Besides, all of the functions involved here are assumed to have Lipschitz continuous gradients, thus there exist constants L 1 and L 2 such that

L 1 := sup i,j,k,l ∇f (x * ) ∞ + ∇g i (x * ) ∞ + ∇h j (x * ) ∞ + ∇G k (x * ) ∞ + ∇H l (x * ) ∞ ,
and for ν depending of x * and λ defined below

L 2 := ∇ 2 xx L R s t,t (x * , G(x * ), H(x * ), ν) ∞ . Let δ > 0 be a constant satisfying max L 1 , L 1 λ ∞ , L 2 δ + o(δ) ≤ . (11)
In other words, δ = O( ).

Construction of the point (x, ŝ, ν) For any x ∈ B δ (x * ), let us construct a point (x, ŝ, ν) that satisfies Definition 4.1. Set

ξ := x -x * , νg := λ g , νh := λ h , νs G := λ G , νs H := λ H ,
where ξ 2 < δ. Using the following notations

I 00 -0 := {i ∈ I 00 : λ G i ≤ 0, λ H i = 0}, I 00 0-:= {i ∈ I 00 : λ G i = 0, λ H i ≤ 0}, I +0 -:= {i ∈ I +0 : λ G i = 0, λ H i < 0}, I +0 + := {i ∈ I +0 : λ G i = 0, λ H i ≥ 0}, I 0+ -:= {i ∈ I 0+ : λ G i < 0, λ H i = 0}, I 0+ + := {i ∈ I 0+ : λ G i ≥ 0, λ H i = 0}, I 00 ++ := {i ∈ I 00 : λ G i > 0, λ H i > 0}
, let ŝG , ŝH , νG , νH and νΦ be as follows:

ŝG :=                      t, i ∈ I 00 -0 , t + c , i ∈ I 00 0-, /2-tλ G i λ G i , i ∈ I 00 ++ , G i (x), i ∈ I +0 , t, i ∈ I 0+ -, -t, i ∈ I 0+ + . ŝH :=                    t + c , i ∈ I 00 -0 , t, i ∈ I 00 0-, t, i ∈ I 00 ++ , t, i ∈ I +0 -, -t, i ∈ I +0 + , H i (x), i ∈ I 0+ , (12) νG
:= λ G i , i ∈ I 00 ++ ∪ I 0+ + , 0, otherwise, νH := λ H i , i ∈ I +0 + , 0, otherwise, νΦ :=        -λ G i ŝH,i -t , i ∈ I 00 -0 ∪ I 0+ -, -λ H i ŝG,i -t , i ∈ I 00 0-∪ I 00 ++ ∪ I +0 -, 0, otherwise.
(x, ŝ, ν) satisfy the stationarity: We verify that (x, ŝ, ν) satisfies the strong epsilon-stationarity conditions, that we divide into the next 6 parts: Let us now run through these 6 conditions. Note that the cases i ∈ I 00 -0 (resp. i ∈ I +0 ) and i ∈ I 00 0-(resp. i ∈ I 0+ ) are symmetric, hence it is sufficient to verify these conditions for i ∈ I 00 -0 ∪ I 00 ++ ∪ I +0 . I. Since (ν g , νh , νs G , νs H ) = (λ g , λ h , λ G , λ H ) it holds that

I. ∇ x L R s t,t (x, ŝ, ν) ∞ ≤ , g i (x) ≤ , νg i ≥ 0, |g i (x)ν g i | ≤ ∀i
∇ x L R s t,t (x * , ŝ, ν) ∞ = 0, g i (x * ) ≤ 0, νg i ≥ 0, |g i (x * )ν g i | = 0 ∀i ∈ {1, . . . , m}, |h i (x * )| = 0 ∀i ∈ {1, . . . , p}.
Now, since x = x * + ξ and ξ ≤ δ, using Taylor theorem at x * yields

g i (x) = g i (x * ) + ξ i g i (x * ) + o(ξ i ) ≤ δL 1 + o(δ) ≤ ,
where the last inequality follows by definition of δ in [START_REF] Guo | Solving mathematical programs with equilibrium constraints[END_REF]. We can proceed in the same way to prove that |h i (x)| ≤ , |g i (x)ν g i | ≤ and ∇ x L R s t,t (x, ŝ, ν) ∞ ≤ .

II. As pointed in Lemma 5.1, the gradient of the Lagrangian with respect to s is given by

∇ s G,i L R s t,t (x, ŝ, ν) = νs G i -νG i + νΦ i (ŝ H,i -t), ∇ s H,i L R s t,t (x, ŝ, ν) = νs H i -νH i + νΦ i (ŝ G,i -t). (13)
It is sufficient to verify these conditions for i ∈ I 00 -0 ∪ I 00 ++ ∪ I +0 + ∪ I +0 -. The following tables summarize the relevant quantities to prove that the gradient (13) vanishes.

∇ s G,i L R s t,t νs G i νG i νΦ i (ŝ H,i -t) ∇ s G,i L R s t,t I 00 ++ λ G i λ G i -λ H i ŝG,i -t 0 0 I 00 -0 λ G i 0 -λ G i ŝH,i -t ŝH,i -t 0 I +0 + λ G i = 0 0 0 ŝH,i -t 0 I +0 - λ G i = 0 0 -λ H i ŝG,i -t 0 0 set νs H i νH i νΦ i (ŝ G,i -t) ∇ s H,i L R s t,t I 00 ++ λ H i 0 -λ H i ŝG,i -t ŝG,i -t 0 I 00 -0 λ H i = 0 0 -λ G i ŝH,i -t 0 0 I +0 + λ H i λ H i 0 ŝG,i -t 0 I +0 - λ H i 0 -λ H i ŝG,i -t
ŝG,i -t 0 III. Using the same technique as in I. it holds that

|G i (x) -ŝG,i | ≤ + |ŝ G,i |, ∀i ∈ I 0+ ∪ I 00 , |H i (x) -ŝH,i | ≤ + |ŝ H,i |, ∀i ∈ I +0 ∪ I 00 , while |G i (x) -ŝG,i | = 0 for i ∈ I +0 and |H i (x) -ŝH,i | = 0 for i ∈ I 0+ . For i ∈ I 00 -0 ∪ I 00 0-, |ŝ G,i | ≤ t + c follows from ŝG,i = t or ŝG,i = t + c . For i ∈ I 00 ++ ,| /2-tλ G i λ G i | ≤ t + c . IV.
By construction, ŝG,i + t ≥ 0 and ŝH,i + t ≥ 0. Now, we verify the complementarity condition for i ∈ I 00 ++ ∪ I 0+ + as otherwise νG i = 0 (the other case follows by symmetry).

νG i (ŝ G,i + t) = 0, for i ∈ I 0+ + , λ G i 2λ G i = 2 , for i ∈ I 00 ++ .
V. By construction of ŝG and ŝH , we have Φ i (ŝ G , ŝH ; t) = 0 for i ∈ I 00 ∪ I 0+ -∪ I +0

-and Φ i (ŝ G , ŝH ; t) < 0 for the remaining indices taking t sufficiently small such that G i (x) > t for i ∈ I +0 + and H i (x) > t for i ∈ I 0+ + . For indices I +0 + ∪ I 0+ + , where Φ i (ŝ G , ŝH ; t) = 0 the complementarity condition is satisfied as νΦ i = 0. VI. By construction it is clear that νG i ≥ 0, νH i ≥ 0. Consider indices i ∈ I 00 0-∪ I 00 ++ ∪ I +0 -, i.e. νΦ i =

-λ H i ŝG,i -t , as i ∈ I 00 -0 ∪ I 0+ -will follow by symmetry.

-λ H i ŝG,i -t is non-negative since ŝG,i > t when λ H i < 0 and ŝG,i < t when λ H i ≥ 0. For the case I 00 0-, λ H i ≤ 0 and ŝG,i = t + c > t. For the case I +0 -, λ H i < 0 and ŝG

,i = G i (x) > G i (x * ) -L 1 δ -o(δ) > G i (x * ) -> t as t is chosen sufficiently small. For the case I 00 ++ , λ H i > 0 and ŝG,i = 2λ G i -t < 0 < t.
This result cannot be extended to a general neighborhood for s due to the definition of strong epsilon-stationary point. However, we prove that s can live in a neighborhood of ŝ (defined in Lemma 5.2) on the boundary of the feasible set.

Theorem 5.1 Assume that f, g, h, G, H have Lipschitz continuous gradients. Let x * be an M-stationary point of (3) and > 0 be small enough. Then, there exist positive constants c, δ, δ , t * ≥ c , t * such that, for all t ∈]c , t * [ and t ∈]0, min(t * , t + c )[, any point in the following set

(x, s) ∈ B δ (x * ) × B δ (ŝ) ∩ {(s G , s H ) : Φ i (s G , s H ; t) = 0, ∀i ∈ I Φ (x * ; t)},
with ŝ defined in [START_REF] Outrata | Optimality conditions for a class of mathematical programs with equilibrium constraints[END_REF], is a strong epsilon-stationary point of the regularization [START_REF] Flegel | On the guignard constraint qualification for mathematical programs with equilibrium constraints[END_REF].

Proof First we stress that we keep the multipliers ν as defined in lemma 5.2

The proof is an adaptation of the one of Lemma 5.2, where are unchanged the parts I, V, VI. Part I remains true since s only appears in ∇ x L R s t,t which is actually independent on s by Lemma 5.1. Part V remains true as the active constraints remain active for all s ∈ B δ ((G(x * ), H(x * )) T ) ∩ {(s G , s H ) : Φ i (s G , s H ; t) = 0 for all i ∈ I Φ (x * ; t)}. Part VI remains true as we keep the definition of multipliers νΦ using ŝ.

Hence, we discuss only the parts II, III and IV.

Let δ be such that max 1, 2 λ ∞ , λ ∞ min(min i∈I +0 {G i (x * )} -t -, t + t, t + c -t) δ < , and s = ŝ + α with α < δ . II. In Lemma 5.2, we proved stronger than strong-epsilon stationarity since we obtained vanishing gradient of the Lagrangian with respect to s, stronger than normwise bounded by . Now, our use of ν introduces perturbations in the second line of first table and lines one and four of the second table in the proof of II in previous Lemma 5.2. Hence, by [START_REF] Scholtes | Convergence properties of a regularization scheme for mathematical programs with complementarity constraints[END_REF] We can now deduce the following result. Theorem 5.2 For any M-stationary point of (3) that satisfies MPCC -MFCQ, there exists a sequence of strong epsilon-stationary points of the regularization (10) that converges to that point. Proof Theorem 5.1 gives more relations between the parameters that are compatible with Corollary 5.1. Indeed for a chosen sequence of arbitrarily small parameters { k }, Corollary 5.1 requires that k = o( tk ) and Theorem 5.1 requires that tk > c and t k must be sufficiently small. Thus, a straightforward application of both of these results allows to conclude.

Conclusions

We studied in this paper a regularization method for the MPCC with approximate stationary points. Motivated by the recent discussion in [START_REF] Kanzow | The price of inexactness: Convergence properties of relaxation methods for mathematical programs with complementarity constraints revisited[END_REF], we introduced a new type of stationary point for the regularized subproblem, called strong epsilon-statonary point, whose sequence converges to an M-stationary point of the MPCC. Additionally, we proved that for any M-stationary point there exists a sequence of strong epsilon-statonary points in its neighborhood. These two theoretical results pave the way to a regularization method with strong properties.

The main perspective of this work is now to study an optimization algorithm able to compute a strong epsilon-stationary point for given t, t and . This algorithm must satisfy at each step exactly the feasibility of the regularized complementarity constraint as well as the complementarity condition of this constraint. We would like to mention that a possible direction is to adapt the active-set strategy studied in [START_REF] Kadrani | A globally convergent algorithm for mpcc[END_REF] since the active-set approach will be able to maintain the feasibility at each step.

Theorem 3 . 2 (

 32 Theorem 9 and 10 [19]

Example 5 . 2

 52 Consider the one-dimensional MPCC min x∈R -x s.t. 0 ≤ x ⊥ x ≥ 0.

  ∈ {1, . . . , m} and |h i (x)| ≤ t + (1 + c) ∀i ∈ {1, . . . , p}; Φ i (ŝ G , ŝH ; t) ≤ 0, νΦ i Φ i (ŝ G , ŝH ; t) = 0 ∀i ∈ {1, . . . , q}; VI. νG ≥ 0, νH ≥ 0, νΦ ≥ 0.

	II. ∇ s L R s t,t (x, ŝ, ν) III. |G i (x) -ŝG,i | ≤ t + (1 + c) , |H i (x) -ŝH,i | ≤ t + (1 + c) ∀i ∈ {1, . . . , q}; ≤ ; ∞
	IV. ŝG,i + t ≥ -, νG i (ŝ G,i + t) ≤ , ŝH,i + t ≥ -, νH i (ŝ H,i + t) ≤ ∀i ∈
	{1, . . . , q};
	V.

  with respect to s G (the case s H follows in the same way), we have∇ s G,i L R s t,t (x, s, ν) = ∇ s G,i L R s t,t (x, ŝ, ν)+ νΦ i α i . Therefore, ∇ s G,i L R s t,t (x, s, ν) < ∇ s G,i L R s t,t (x, ŝ, ν) + |ν Φ i |δ < . III. In Lemma 5.2, we proved that ∀i ∈ {1, . . . , q}, |G i (x) -ŝG,i | ≤ t + (1 + c) , |H i (x) -ŝH,i | ≤ t + (1 + c) . So, using that δ ≤ gives ∀i ∈ {1, . . . , q}, |G i (x) -sG,i | ≤ t + (2 + c) , |H i (x) -sH,i | ≤ t + (2 + c) , which proves III.IV. From Lemma 5.2, when ŝG,i + t ≥ 0 (resp. ŝH,i + t ≥ 0) we deduce sG,i + t ≥ -(resp. sH,i + t ≥ -) as δ < . Now, we verify the complementarity condition for i ∈ I 00 ++ ∪ I 0+

+ : νG i (s G,i + t) ≤ νG i (ŝ G,i + t) + νG i sG,i = δ λ G i ≤ 2 , for i ∈ I 0+ + , 2 + δ λ G i ≤ , for i ∈ I 00 ++ .

Acknowledgements This research was partially supported by a grant from "l'Ecole des Docteurs de l'UBL" and "le Conseil Régional de Bretagne". This research was partially supported by NSERC and FRQNT grants.