- Email: [email protected]
A filter method for solving nonlinear complementarity problems
Applied Mathematics and Computation 167 (2005) 677–694 www.elsevier.com/locate/amc
A filter method for solving nonlinear complementarity problems Pu-yan Nie Department of Mathematics, College of Information Science and Technology, Jinan University, Guangzhou 510632, PR China
Abstract Filter methods are extensively studied to handle nonlinear programming problems recently. Because of good numerical results, filter techniques are attached importance to. In this paper, filter approaches are employed to tackle nonlinear complementarity problems (NCPs). Firstly, NCP conditions are transformed into a nonlinear programming problem. Then, to obtain a trial step, the corresponding nonlinear programming problems are solved by some existing strategies. Moreover, filter criterion is utilized to evaluate a trial iterate. The purpose in this paper is to employ filter approaches to attack NCPs. In essence, multi-objective view is utilized to attack NCPs because the idea of filter methods stems from multi-objective problems. Furthermore, a new filter method, based on the special two objects which differs from others, is brought forward. Moreover, Maratos effects are overcome in our new filter approach by weakening acceptable conditions. 2004 Elsevier Inc. All rights reserved. Keywords: Filter methods; Nonlinear complementarity problems; Nonlinear programming; Global convergence
E-mail address: [email protected] 0096-3003/$ - see front matter 2004 Elsevier Inc. All rights reserved. doi:10.1016/j.amc.2004.06.125
678
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
1. Introduction In this work, the system of nonlinear complementarity problem (NCP), with the following form, is considered F ðxÞ P 0 x P 0
and
xT F ðxÞ ¼ 0;
ð1:1Þ
where x 2 Rn, F : Rn ! Rn is a given function whose properties will be specified later and n is the dimension of the variables. The NCP has been utilized as a general framework for quadratic programming, linear complementarity, and the other mathematical programming. Different concepts have been developed to treat this problem. In the last few years, growing attention has been paid to approaches which employ a reformulation of NCP as a system of nonlinear equations or a minimization problem [4,9,11,12,15,16]. Theoretical results have been established to present conditions under which a stationary point of the minimization problem is a solution to the complementarity problem. Many researchers reformulate (1.1) as a system of nonsmooth equations. Then, many kinds of NCP functions appear to lead to a system of equations. We have the NCP functions with the form U(x, F(x)) and aim to find x satisfying Uðxi ; F i ðxÞÞ ¼ 0;
i ¼ 1; 2; . . . ; n;
ð1:2Þ
where U : Rn ! Rn and (1.2) is equivalent to (1.1) to a certain degree. Many algorithms based on (1.2) have been proposed. For example, NewtonÕs methods and generalized NewtonÕs approaches (see in [4,10,12]). To handle (1.1) by a minimization problem, the constraints are C ¼ fx P 0; F ðxÞ P 0g; and the objective function is kXF(x)k2 or something similar, where X = diag{x1,x2, . . ., xn}. Some optimization techniques can be employed. This type of method is re-visited in this paper. Moreover, the filter technique is employed to attack (1.1). In this work, filter methods are utilized to deal with NCPs. In the filter method, at each step, after a subproblem is solved, the filter criterion is employed to determine whether to accept the trial point or not. To give the subproblem, (1.1) is described as another equivalent form Xn minimize F 2 ðxÞx2i i¼1 i subject to
F j ðxÞ P 0; x P 0:
The algorithm is proposed, which are based on (1.3) to tackle (1.1).
ð1:3Þ
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
679
Filter approaches are introduced as follows. Filter method is initially put forward by Fletcher and Leyffer in 2002 [5]. Filter strategies do very well to balance the objective function and constraint conditions. Because of good numerical results, filter methods have been combined with trust region approaches, SQP techniques, pattern search method, interior point strategy and composite-like step methods [2,5–8,13,17]. In [18], super-linear local convergence is achieved for filter-SQP methods. In an optimization approach, we hope to find a satisfying point, which relates not only to objective function but also to constraint conditions. Two functions related close to constraint conditions and objective function are thus given as follows hðxÞ ¼
n X
maxf0; F j ðxÞg;
ð1:4Þ
j¼1
pðxÞ ¼
n X
2
x2i F i ðxÞ þ rhðxÞ;
ð1:5Þ
j¼1
where r is a constant. Moreover, x P 0 is always satisfied. At the kth point, if a new point xk + sk is obtained, a filter rule to judge k x + sk is presented. Certainly, we hope that p(xk + sk) < p(xi) or h(xk + sk) < h(xi) where i 2 Fk and Fk is a set consisted of ‘‘good’’ points. To obtain the convergent results, the strong conditions are required pðxk Þ 6 chðxk Þ þ pðxi Þ
ð1:6Þ
hðxk Þ 6 bhðxi Þ;
ð1:7Þ
or
where i 2 Fk and 1 > b > c > 0 are two constants. Certainly, whether a point is accepted by the filter is related to when it is generated (see [5] in detail). The filter is based on (1.6) and (1.7). xk + sk is accepted by the filter if (1.6) and (1.7) are satisfied for all i 2 Fk . Then, let xk+1 = xk + sk and update F. For convenience, denote Dk ðxkþ1 Þ ¼ fiji 2 Fk ; pðxi Þ P pðxkþ1 Þ and hðxi Þ P hðxkþ1 Þg: kþ1
Let Dk ðx Then,
Þ be the set of the points which is dominated by x
Fkþ1 ¼ Fk
[
fk þ 1g=Dk ðxkþ1 Þ:
k+1
ð1:8Þ k
:= x + sk. ð1:9Þ
We point out that the function of p(x) in this paper is different from those in the other papers. The choice of r is flexible. With suitable r < 0, Maratos effect can be overcome because the acceptable criterion is relaxed.
680
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
The paper is organized as follows: In Section 2, an algorithm is put forward. The algorithm is analyzed in Section 3. Some numerical results and remarks are listed in Section 4.
2. A filter algorithm for NCPs To handle (1.1) based on (1.3), we hope to find a better point at the beginning of the kth iteration. A step is obtained by solving the following subproblem " #T n X 1 2 T minimize WH ðxk F ðxk ÞÞi d k ðdÞ ¼ d Bk d þ r 2 i¼1 subject to
rF i ðxk Þd þ F i ðxk Þ P 0;
ð2:1Þ
kdk 6 D; d þ xk P 0; where F(xk) is a vector and Bk is a Hessian of Lagrangian function to (1.3) or an approximate Hessian matrix. The new point is defined as xk ðDÞ ¼ xk þ d k ; where dk is generated by (2.1). To determine whether to accept the trial step or not, some criterion is utilized. The following merit function is thus employed. mk ðxÞ ¼
n X 2 ððxF ðxÞÞi Þ:
ð2:2Þ
i¼1
The predicted reduction and actual reduction, to judge the new trial point, are presented as follows: Aredðsk Þ ¼ ðmk ðxk þ sk Þ mk ðxk ÞÞ; ð2:3Þ 1 Predðsk Þ ¼ ðsk ÞT Bk sk þ gTk sk ; ð2:4Þ 2 P where gk ¼ r ni¼1 ðxk F ðxk ÞÞ2i . When (2.1) is inconsistent, a restoration algorithm is employed to decrease the constraint violation degree and to obtain a better point by solving ($Fi(xk)sk + Fi(xk)) = 0 where ($Fi(xk)sk + Fi(xk)) = min{0, $Fi(xk)sk + Fi(xk)}. For convenience, the following notation is applied throughout this paper.
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
681
Definition 1 hIk ¼ minfhi jhi > 0; i ¼ 1; 2; . . . ; kg; pFk ¼ minfpi jhi ¼ 0; 1 6 i 6 kg; and pIk is the corresponding value to hIk . Namely, some way, to let (h(x 0 ), p(x 0 )) lying near (0, pF) or (hI, pI), is utilized. In the algorithm of this paper, the restoration algorithm is employed. For the update of Dk, the way of YuanÕs [19,14] is employed. The algorithm is formally stated as follows. Algorithm 1 (Trust Region Filter Algorithm For NCPs) Step 0. Choose D0, x0, r, a1, a2, b, g, g2, c, where D0 > 0, 0 < a1, a2, b, g, g2 < 1 and c 2 ð0; 12Þ; b 2 ð12 ; 1Þ. Set k :¼ 0; F ¼ fx0 g. Step 1. Compute hIk ; pIk ; pFk . Step 2. Solve (2.1) to get sk. Step 3. If (2.1) has no solution or sk = 0 but hk > 0, utilize Restoration Algorithm to get skr . Set xk :¼ xk þ skr , update hIk and goto Step 2. If sk = 0 and hk = 0 , then stop. Step 4. Compute rk ¼
Aredðsk Þ : Predðsk Þ
If rk 6 g, then, let xkþ1 :¼ xk ; Dkþ1 :¼ 12 Dk and k := k + 1 and go to Step 2. Compute ^ pk ¼ pðxk þ sk Þ; ^ hk ¼ hðxk þ sk Þ. If xk + sk is not acceptable to the filter, then Dkþ1 :¼ 12 Dk and goto Step 2. If hðxk Þ 6 Dk minfg2 ; a1 Dak 2 g, go to Step 5. If hðxk Þ P Dk minfg2 ; a1 Dak 2 g, call Restoration Algorithm to produce a point xkr :¼ xk þ skr such that: (A) xkr is acceptable to the filter; (B) hðxkr Þ 6 g2 minfhIk ; a1 D2k g. Let xk :¼ xkr and goto Step 2. Step 5. xk+1 := xk + sk is accepted by the filter. pkþ1 ¼ ^pk ; hkþ1 ¼ ^hk and remove the points dominated by (pk+1, hk+1) from the filter. Set Dk+1 := 2Dk and generate Bk+1. 2 Step 6. If hðxkþ1 Þ 6 Dkþ1 minfg2 ; a1 Dakþ1 g, then, set k := k + 1 and goto Step 1. Otherwise, k := k + 1 and call Restoration Algorithm to produce a point xkr ¼ xk þ skr satisfying (A) and (B). Let xk :¼ xkr and goto Step 2. If hðxk Þ > Dk minfg2 ; a1 Dak 2 g, a restoration algorithm is utilized to compute x such that (A) and (B) are all met. kr
682
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
In a restoration algorithm, we therefore aim to decrease the value of h(x). The direct way is employed Newton method or similar ways to solve F(x + sr) = 0. Let M jk ðdÞ ¼ ðhðxkj Þ hðxkj þ dÞÞ: A restoration algorithm is then given, where
Wjk ðdÞ ¼ kðrF ðxkj Þd þ F ðxkj Þ k1 hðxkj Þ;
ð2:5Þ
and rjk ¼
M jk ðskj Þ : Wjk ðskj Þ
ð2:6Þ
Algorithm 2 (Restoration Algorithm) Step 0. Let xk0 :¼ xk ; D0k :¼ Dk ; j :¼ 0; g; g2 2 ð0; 1Þ. Step 1. If (A) and (B) are met, then, let xkr :¼ xkj and stop. Step 2. Compute minimize
Wjk ðdÞ
subject to
kdk1 6 Dk
ð2:7Þ
to obtain skj . Calculate rjk . Step 3. If rjk 6 g, then, let xkjþ1 :¼ xkj ; Djþ1 :¼ 12 Djk ; j :¼ j þ 1 and goto Step 2. k jþ1 Step 4. Otherwise, let xkjþ1 :¼ xkj þ skj ; Dk :¼ 2Djk . Update Ajþ1 k ; j :¼ j þ 1 and goto Step 1. We point out that xkj is the jth restoration point following xk and xkr is the final restoration point following xk. The above Restoration Algorithm is a NewtonÕs method for F(x) = 0. This method is frequently utilized. Of course, there are other Restoration Algorithms, for example, the interior point restoration algorithm in [17] and the SLP restoration algorithm in [2]. We can also solve a quadratic programming to obtain the result. In fact, Algorithm 2 is a linear programming approach to reduce the constraint violation degree. Actually, Pred(sk) and Ared(sk) are permitted to be less than 0 at the same time. Therefore, this is a nonmonotonic algorithm. We hope the restoration algorithm to terminate finitely. In Section 3, this property will be shown. About r, it will be discussed in Section 4.
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
683
3. The convergence properties Just as that in [1,2,7,8], the analysis to the algorithm are based on the following standard assumptions. Further, to obtain the convergence, the sufficient reduction plays a crucial role throughout. Assumption 1 (1) The set {xk} 2 X is nonempty and bounded. (2) The function F(x) is twice continuously differentiable on an open set containing X. (3) When solving (2.7), we have
Wjk ðdÞ ¼ hðxkj Þ k½rF ðxkj Þd þ F ðxkj Þ k1 P b2 minfhðxkj Þ; Djk g; where b2 > 0 is a constant. (4) The matrix sequence {Bk} is bounded. (5) The Restoration Algorithm has a solution satisfying kskr k 6 s0 hk . (1) and (2) are the standard assumptions. (3) is the sufficient reduction condition, which is a very weak condition because Cauchy step satisfies this condition. It is regarded as a condition in this paper. In a trust region method, (3) guarantees the global convergence. (4) plays an important role to obtain the convergent result. But it has minor effects to the local convergent rate. The following results are based on Assumption 1. Analyzing the restoration algorithm, we obtain Lemma 1. The Restoration Algorithm 2 terminates finitely under Assumption 1. Proof. Assume the restoration algorithm does not terminate finitely. The termination criterion will be satisfied if limk!1 hðxkj Þ ¼ 0. It is thus reasonable that there exists > 0 with hðxkj Þ > for all j. We show it by contradiction. Let K ¼ fj : rjk > gg: According to Assumption 1 and (2.6), we have 1>
1 1 X X X k j1 ðhðxkj1 Þ hðxkj ÞÞ P g Wj1 Þ P gb2 minf; Djk g: k ðs j¼1
j¼1
k2K
Therefore, Djk ! 0; j 2 K. By the algorithm, we have Djk ! 0, for all j. Moreover,
684
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
Wjk ðskj Þ ¼ M jk ðskj Þ þ oðDjk Þ; P Djk from Step 4 of Restoration Algorithm. Namely, when Djk ! 0. Thus, Djþ1 k j fDk g is increased when Djk is very small, which contradicts Djk ! 0. Consequently, the result is obtained and the proof is complete. h The optimal properties of Algorithm 1 considered, the following result is obtained immediately. Lemma 2. Every new iteration xk+1 5 xk is acceptable to the filter set F. Proof. From Algorithm 1, a new iteration is produced in Steps 5 and 6. In both cases, xk+1 is accepted by the filter. The result therefore holds and the proof is complete. h Theorem 1. Under Assumption 1, suppose there are infinitely many points added to the filter. Then lim hðxk Þ ¼ 0:
k!1
Proof. See [7] or [13].
h
If there are finitely many points added to the filter, then, the following conclusions hold. Theorem 2. Under Assumption 1, suppose there are finitely many points added to the filter. Then hðxk Þ ¼ 0; for certain k1 and k > k1. Proof. The result is obvious from the algorithm.
h
Similar to [2,7], the global convergence theorem concerns Kuhn–Tucker (K–T) necessary conditions under a Mangasarian–Fromowitz constraint qualification (MFCQ). This is an extended form of Fritz–John conditions for a problem that includes equality constraints. We analyze the properties of solution to (2.1). Without loss of generality, assume that k$Fik 6 M, k$2Fik 6 M2 for all x 2 X and i = 1, 2, . . ., n. Then, we consider Fi(xk + sk). From Taylor expansion, if xk + sk is the solution to (2.1) there is T
T
F i ðxk þ sk Þ ¼ F i ðxk Þ þ rF i ðxk Þ sk þ ðsk Þ r2 F i ðxk þ hðxk þ sk xk ÞÞsk ;
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
685
where h 2 (0, 1). When Dk is small enough and Fi(xk) = min{0, Fi(xk)}, we have Fi(xk) + $Fi(xk)Tsk P 0 for i = 1, 2, . . ., n. We therefore obtain n X kðsk ÞT r2 F i ðxk þ hðxk þ sk xk ÞÞsk k hðxk þ sk Þ 6 i¼1
6 nkr2 F i ðxk þ hðxk þ sk xk ÞÞkD2k 6 2n2 M 2 D2k :
ð3:1Þ
Lemma 3. If the (k + j)th step relates to Restoration Algorithm, j is large enough and Dk ! 0, then, 1þa2 2 a1 D1þa kþj 6 hkþj 6 4a1 Dkþj :
ð3:2Þ
Proof. The first part of the inequality is obvious by the algorithm. Consider the second part of the inequality. If xk+j1 is refused by the filter, then Dk+j1 6 2Dk+j, xk+j = xk+j1 and there is no Restoration Algorithm at the (k + j 1)th 1 step. Accordingly, the result is true. If xk+j1 is accepted and Dkþj 6 ð4a21 Þ1a2 , according to Step 5 of Algorithm 1, we have Dk+j1 6 Dk+j and 2 hkþj 6 2D2kþj1 6 2D2kþj 6 4a1 D1þa kþj :
The result is established.
h
About (2.1), Fletcher et al. do very excellent analysis in [7], which is listed as follows Theorem 3. Let Assumption 1 hold and xw 2 X be a feasible point of problem (1.3) at which MFCQ holds, but which is not a K–T point. Then there exists a neighborhood N0 of xw and positive constants n1, n2, n3 such that for all x 2 N0 ˙ X and all Dk for which n2 hðcðxÞÞ 6 Dk 6 n3 ;
ð3:3Þ
it follows that (3.2) has a feasible solution d at which the predicted reduction satisfies 1 Wk ðdÞ P n1 Dk : 3 If Dk 6 (1 g3)n1/3nM2, then n X 2 2 ððxk Þi F 2i ðxk Þ ðxk þ dÞi F 2i ðxk þ dÞÞ P g3 Wk ðdÞ; i¼1
where g < g3. qffiffiffiffiffiffiffiffi k k If hk > 0 and Dk 6 n2bh 2 M , then, h(x + d) 6 bhk. 2
ð3:4Þ
ð3:5Þ
686
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
Proof. See Lemmas 4 and 5 in [7]. h Lemma 4. Assume that Algorithm 1 generates a sequence {xk} and its accumulation xw satisfying the conditions of Theorem 3, if there are infinitely many points enter into the filter, then, there are only finitely many points related with Restoration Algorithm. Proof. Under conditions of Theorem 3, if there are infinitely many points related to Restoration Algorithm, Dk ! 0. (See Theorem 7 in this paper.) According to Theorem 3, consider the conditions when a point is accepted. If r 5 0, we show the result. When ( 1 sffiffiffiffiffiffiffiffiffiffiffi) ð1 g3 Þn1 g 3 n 1 a2 2bhk Dk 6 min n3 ; ; ; ; ð3:6Þ 3nM 2 3a1 jrj n2 M 2 Dk P
1 1þa1 2 h ; a1 k
ð3:7Þ
and hk > 0, then, (3.4) and (3.5) are true. Moreover, according to Theorem 3, we obtain h(xk + sk) 6 bhk. xk + sk will be therefore accepted by the filter because n n X X ðxk Þ2i F 2i ðxk Þ ðxk þ sk Þ2i F 2i ðxk þ sk Þ pðxk Þ pðxk þ sk Þ ¼ i¼1
i¼1
þ rðhk hðxk þ sk ÞÞ 1 2 >0 P g3 Wk ðsk Þ jrjhk P g3 n1 Dk jrja1 D1þa k 3 and h(xk + sk) 6 bhk. By the way, the interval between (3.6) and (3.7) 1is not empty if D is small 1þa
enough. Further, under (3.6) and (3.7), Dkþj P a11 hk 2 is satisfied for j P 1 if k is large enough. Namely, every trial step satisfies (3.7) for j P 1. If there is a Restoration Algorithm, from (3.4), there exists a constant 13 n1 P 0 > 0 such that kþj WH Þ P 0 kskþj k: kþj ðs
From (3.2) and the above analysis, we have kþj pðxkþj Þ pðxkþj þ skþj Þ P g3 WH Þ þ rhðxkþj Þ rhðxkþj þ skþj Þ kþj ðs kþj P g3 WH Þ jrjhðxkþj Þ kþj ðs 1þa2 2 2 P g3 0 Dkþj jrja1 D1þa kþj P 2cn M 2 Dkþj
P 2cn2 M 2 D2kþj P chðxkþj ðDkþj ÞÞ;
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
687
g3 0 k+j r 2 if Dakþj 6 jrja1 þ2cn )k 2 M . The above inequality holds because limk!1k(s 2
6 limk!1s0hk+j = 0 and Dk ! 0. The above inequality means that p is monotonic 1 1þa decreasing and there exists a trial point accepted by the filter for Dkþj P a11 hk 2 and 1 sffiffiffiffiffiffiffiffiffiffiffi) ð1 g3 Þn1 g3 n1 a 2 2bhk 6 min n3 ; ; ; : 3nM 2 3a1 jrj n2 M 2 (
Dkþj
1 1þa
Therefore, Dkþjþ1 P a11 hk 2 . We now show that (3.7) is reasonable. If (3.7) were always false, there 1 1þa2
should be Dk < a11 hk
and limk!1Dk = 0 according to limk!1hk = 0. Further1
more, from (3.2) when 2Dk+1 P Dk and Dkþ1 6 ð8na2 1M 2 Þ1þa2 for sufficient large k, hkþ1 6 2n2 M 2 D2k 6 8n2 M 2 D2kþ1 : 1 1þa
1 1þa2
Thus, a11 hkþ12 6 Dkþ1 , which contradicts Dk < a11 hk If r = 0, (3.6) is substituted for (
ð1 g3 Þn1 ; Dk 6 min n3 ; 3nM 2
because of limk!1Dk = 0.
sffiffiffiffiffiffiffiffiffiffiffi) 2bhk : n2 M 2
Similar result is obtained. The result is accordingly obtained and the proof is complete. h We point out that the result in Lemma 4 holds for any bounded constant r. The global convergence can be obtained combining the results of Theorems 1–3 and Lemma 4. Theorem 4. Let Assumption 1 hold and let the sequence {xk}, where {xk} 2 X is generated by Algorithm 1, converges to a feasible point of problem (1.3) at which MFCQ holds. Then, {xk} has an accumulation which is a K–T point of (1.3). Proof. If Algorithm 1 terminates finitely, the result holds apparently. Assume Algorithm 1 terminates infinitely. If the result were false, from Theorem 1 and the above analysis, there would exist an integer k0 such that hk < n2 is small H enough, mk ðxk Þ mk ðxk þ dÞ P g½WH k ð0Þ Wk ðdÞ for all k > k0. Furthermore, the Restoration Algorithm does not appear when k > k0 from Lemma 4. For convenience, we denote H K 2 ¼ fkjk > k 0 ; and mk ðxk Þ mk ðxk þ dÞ P g½WH k ð0Þ Wk ðdÞg:
ð3:8Þ
688
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
Moreover, from Assumption 1 there exists X X 1> ½pðxk Þ pðxkþ1 Þ P ½g3 Wk ðdÞ þ rðhk hkþ1 Þ K2
P
X 1 K2
3
K2
g3 n1 Dk þ rðhk hkþ1 Þ :
Consequently, limk!1Dk = 0 for k 2 K2 because hk ! 0. Taking into account the choice of Dk, limk!1Dk = 0 for all k. There thus exists Dk such that 2 Da > k
6 ðjrj þ maxfð1 bÞ; cgÞ g3 n1 a1
ð3:9Þ
for sufficiently large k and 2 P hðxk þ sk Þ: a1 D1þa k
ð3:10Þ
(3.10) is issued from Step 4 of Algorithm 1. We consider the next trial step xk + sk. By the above analysis, p(xk) is monotonic decreasing and 1 pðxk Þ pðxk þ sk Þ P g3 n1 Dk þ rðhk hðxk þ sk ÞÞ 3 1 1þa2 2 P g3 n1 Da þ rðhk hðxk þ sk ÞÞ k Dk 3 2 P ðjrj þ maxfð1 bÞ; cgÞ2a1 D1þa þ rðhk hðxk þ sk ÞÞ k P maxfð1 bÞ; cghðxk þ sk Þ; k
ð3:11Þ
k
which means that x + s will be accepted by the filter. Dk will not be reduced when Dk is small enough and k > k0, which contradicts limk!1Dk = 0. The result thus holds and the proof is complete. h The K–T points of (1.3) but not the solutions to (1.1) are defined as local infeasibility points as follows Definition 2. If x* is a K–T point but not the solution to (1.1), we call x* a local infeasibility point. In fact, there may be no solution or we may fail to find the solution to (1.1) in some region. For global case, it is very difficult to describe it just as the other global problems. Local infeasibility points are therefore defined. The following result is obviously obtained from Theorem 4. Theorem 5. Under assumptions of Theorem 4, the algorithm termites finitely or has an accumulation point which is the solution of (1.1) or a local infeasibility point.
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
689
Under certain conditions, the stationary points of (1.3) are the solutions of (1.1). For example, if the solution to (1.3) is zero, it is the solution of (1.1). The global minimization is just the solution of (1.1) if there exists. We hope that we can obtain the global minimum to (1.3). We give a condition to guarantee this result. Definition 3. The mapping F is a P0 function, i.e., for every x1, x2 2 Rn with x1 5 x2, there exists an index i 2 {1, 2, . . ., n} such that x1i 6¼ x2i
and ðx1i x2i ÞðF i ðx1 Þ F i ðx2 ÞÞ P 0:
Correspondingly, we define the P0-matrix. The gradient of P0 function is a P0matrix. M 2 Rn·n is said to be a P0-matrix if all the principal minors are nonnegative. P0-matrix has the following properties. Lemma 5. M 2 Rn·n, the following two statements are equivalent (1) M is a P0-matrix. (2) For any x 2 Rn 5 0, there exists an index xi such that xi(Mx)i P 0.
Proof. See Theorem 3.4.2 in [3].
h
Now, the result is presented as follows. Theorem 6. If F is a P0 function, then xw is a global minimum if and only if it is a K–T point. Proof. The ‘‘only if’’ part is obvious. We show the ‘‘if’’ part. If xw is not the global minimum to (1.3), without loss of generality, assume that xH S > 0
and
F S ðxH Þ > 0;
where S is a subset of {1, 2, . . ., n}. According to K–T conditions, we thus have kS = 0, lS = 0, where kS and lS are the corresponding Lagrangian to the constraints. Further, from K–T conditions there exists X 2 H H H 2 H ðxH ð3:12Þ i Þ F i ðx Þrxj F i ðx Þ þ xj F j ðx Þ ¼ 0 i2S 2 H H 2 H for any j 2 S. Denoting ui ¼ ðxH k i Þ F k i ðx Þ > 0 and vi ¼ xk i F k i ðx Þ where ki 2 S, from (3.12), we have
rs Fu þ v ¼ 0:
690
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
By v > 0, we have $sFu < 0. This contradicts the fact that $sF is a P0-matrix issuing from Lemma 5. (The sub-matrix of a P0-matrix is also a P0-matrix.) The result is consequently established. h But when there is no solution to (1.1), some other algorithms may run slowly near a local infeasibility point while the algorithm may avoid it because its second information is made full use of. Now, the following result, which is employed in Lemma 4, is proved. Theorem 7. If there are infinitely many points related to Restoration Algorithm and they are far away from K–T points, then, limk!1Dk = 0. Proof. Apparently there is a subsequence of {Dk} satisfying limj!1 Dtj ¼ 0 because hk ! 0 and there are infinitely many points related to Restoration Algorithm. If xk is far away from K–T points, from the result in [19], there exist k two positive constants s and such that WH s minfD0 ; g. Further, k ðs Þ P when (
Dkþj
ð1 g3 Þn1 6 min n3 ; ; 3nM 2
sffiffiffiffiffiffiffiffiffiffiffi) 2bhk ; n2 M 2
this trial step will be accepted by filter in terms of Theorem 3. It is therefore impossible to have any subsequence of {Dk} satisfying limj!1 Dtj ¼ 0, which contracts the above fact. Therefore, limk!1Dk = 0 and the proof is complete. h The global convergence properties have been obtained. Actually, n0 is a very small real. Algorithm 1 is an SQP trust-region method combined with a new filter technique for NCPs. The approach in this paper has the following advantages: • This method supports a new approach to get the solution of (1.1). In early filter techniques, objective function and constraint violation degree function are reduced. While in the algorithm in this paper, we aim to reduce the objective and some other function, which is a combination of constraint violation and objective function. When (1.1) has no solution, a technique to terminate is proposed. • It is a method without multipliers and it is very difficult to deal with multipliers. • This approach is a nonmonotonic one, which may be beneficial to obtain global solutions. • The subproblem employed in this paper is based on a smoothing problem.
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
691
• The method in this paper has no requirement about feasible initial point, which is convenient.
4. Conclusion remarks and numerical results The researchers concentrate on filter techniques because filter approaches can efficiently balance the objective function and constraint violation function. The combination of various objective function with filter technique brings out a new method. Of course, diverse kinds of p can be used in various methods. As for r, it can be positive or negative because the filter method can guarantee hk ! 0. Therefore, it has no effect on the convergence of the algorithm. In the other papers about filter method, to attack (1.3), authors assume that r = 0. When r > 0, the accepting condition becomes strict. But, when r < 0 the accepting condition becomes moderate. Our filter approach is therefore regarded as a general version about the choices of r in filter method family. When a suitable r is employed, some bad cases may be overcome. Bk is an approximate Hessian. Some approaches, such as BFGS, can be applied to update. It has minor effect on the convergent properties. We also do not discuss the initial trust region radius, although it plays an important role to obtain good numerical results. In Assumption 1, it is required that {Bk} be bounded. Actually, when kBkk 6 qk where q is a positive constant, the convergent results also hold. We hope that we can moderate the conditions. But MFCQ seems necessary. Some examples illustrate the algorithm. (In the following examples, the toleration is 1.0e7 and b = 0.05, c = 0.98, a1 = a2 = 0.5, g = 0.25, g2 = 0.75). Our algorithm compares with that in [10]. The toleration is 1.0e6 in [10]. Example 1. Let F(x) = Mx + q, where 2
4
1
6 6 1 4 6 6 6 6 0 1 6 M ¼6 6 6 6 6 6 0 0 4 0
0
0
1 4
0
0
and q = (1, 1, . . . , 1)T.
0
3
7 0 7 7 7 7 0 7 7 7 7 7 7 7 1 7 5 4
692
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
The results of Example 1 lie in the following table, where the initial point is x0 = (0,0, . . ., 0)T.
NF NG IT(JQÕs)
n=4
n=8
n = 16
n = 32
n = 64
n = 128
n = 256
2 2 4
2 2 4
2 2 4
2 2 4
2 2 4
2 2 4
2 2 4
where NF and NG mean the number of calculations of f and gradient, respectively. IT (JQÕs) means the number of iterates in [10]. This problem has no relation with r. But it works very well. Example 2. Let F 1 ðxÞ ¼ 3x21 þ 2x1 x2 þ 2x22 þ x3 þ 3x4 6; F 2 ðxÞ ¼ 2x21 þ x1 þ x22 þ 10x3 þ 2x4 2; F 3 ðxÞ ¼ 3x21 þ x1 x2 þ 2x22 þ 2x3 þ 9x4 9; F 4 ðxÞ ¼ x21 þ 3x22 þ 2x3 þ 3x4 3:
pffiffi This example has one degenerate solution ð 26 ; 0; 0; 12Þ and one nondegenerate solution (1, 0, 3, 0). The results of Example 2 is listed in the following table if we let r = 0.
Intial point
NF
Ng
IT (JQÕs)
(0, 0, 0) (1, 1, 1)
12 16
8 10
13 8
Example 2 appears in [12]. The following example and Example 1 come from [9]. Example 3. Let F(x) = Mx + q, where M ¼ diag
1 2 ; ;...;1 n n
and q = (1, 1, . . ., 1)T.
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
693
The initial point is x0 = (0, 0, . . ., 0)T. We give the results.
n = 10 n = 20 n = 40 n = 80 n = 160 n = 320 No trust region restrict in RA
NF
3
3
3
3
3
3
NG Initial trust region NF radius is 1-norm of gradient NG IT(JQÕs)
2 2
2 2
2 2
2 2
2 2
2 2
2 4
2 4
2 4
2 4
2 4
2 4
From Examples 1–3, the algorithm do very well compared with other methods to deal with NCP function. We let the constant r be flexible. Certainly, we can give some criterions to choose r or to update r. But, the monotonic property of r is not necessary if r is allowed for change. The value of r is not required to be very large. In a restoration algorithm, the instinctive thought is to minimize h(x). To achieve the convergence properties of Algorithm 1, the sufficient reduction condition is necessary because both in filter method and in trust region approach the sufficient reduction condition is very important. the algorithm, the merit PIn n k 2 2 k function is necessary to enable the information of P ðx Þi F i ðx Þ to keep in fit i¼1 n k with the quadratic model. Moreover, we can use ðx Þi F i ðxk Þ instead of i¼1 Pn k 2 2 k i¼1 ðx Þi F i ðx Þ. The similar results can be obtained. About the parameters b and c, where 1 > b > c > 0. b is required to close enough to 1, while c sufficiently close to 0. There are several advantages: (1) The acceptable criterion is not very strong. (2) The relation of p and h becomes more implicit. Certainly, some conditions in the algorithm can be modified to get the global convergence. In the algorithm, we hope to avoid the linkage between p and h. But it is by no means yet. Because when the relation is lost, the convergence properties are not guaranteed. In this way will a nondescent direction be produced. As for the norm of the constraint violation function, the other forms can be chosen, such as l1-norm and l2-norm. For example, l1-norm is used in [1]. About the filter method, there is much work to do. For example, filter methods with three dimensions are interesting. Meanwhile, the choice of alternative for the optimal measure is still an interesting topic. In the numerical results, when this technique is taken steps, it plays better to a certain degree. It is well known, filter method has the advantage of high-efficiency.
694
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
Acknowledgment The Project Sponsored by the Scientific Research Foundation for the Returned Overseas Chinese Scholars, State Education Ministry.
References [1] C. Audet, J.E. Dennis Jr., Combining patter search and filter algorithms for derivative free optimization, SIAM Journal on Optimization 14 (2004) 980–1010. [2] C.M. Chin, R. Fletcher, On the global convergence of an SLP-filter algorithm that takes EQP steps, Mathematical Programming 96 (2003) 161–177. [3] R.W. Cottle, J.S. Pang, R.E. Stone, The Linear Complementarity Problem, Academic Press, New York, NY, 1992. [4] A. Fischer, Solution of monotone complementarity problems with locally Lipschitzian functions, Mathematical Programming 76 (1997) 513–532. [5] R. Fletcher, S. Leyffer, Nonlinear programming without a penalty function, Mathematical Programming 91 (2002) 239–269. [6] R. Fletcher, N.I.M. Gould, S. Leyffer et al., Global convergence of a trust-region SQP-filter algorithm for general nonlinear, SIAM Journal on Optimization 13 (2003) 635–659. [7] R. Fletcher, S. Leyffer, P.L. Toint, On the global convergence of a filter SQP algorithm, SIAM Journal on Optimization 13 (2002) 44–59. [8] R. Fletcher, S. Leyffer, A bundle filter method for nonsmooth nonlinear optimization, Tech. Report NA/195, Department of Mathematics, University of Dundee, 1999. [9] C. Geiger, C. Kanzow, On the solution of monotone complementarity problems, Computational Optimization and Application 5 (1996) 155–173. [10] H.Y. Jiang, L.Q. Qi, A new nonsmooth equations approach to nonlinear complementarity problems, SIAM Journal on Control and Optimization 35 (1997) 178–193. [11] C. Kanzow, Some equation-based methods for the nonlinear complementarity problem, Optimization Methods Software 3 (1994) 327–340. [12] J.J. More´, Global methods for nonlinear complementarity problems, Mathematics of Operations Research 21 (1996) 589–614. [13] P.Y. Nie, Composite-step like filter methods for equality constraint problems, Journal of Computational Mathematics 21 (2003) 613–624. [14] P.Y. Nie, A null space method for system of nonlinear equations, Applied Mathematics and Computation 149 (2004) 215–226. [15] J.S. Pang, S.A. Gabriel, NE/SQP: A robust algorithm for nonlinear complementarity problem, Mathematical Programming 60 (1993) 295–337. [16] J.S. Pang, A B-differentiable equations based, globally and locally quadratically convergent algorithm for nonlinear programming, complementarity, and variational inequality problems, Mathematical Programming 51 (1991) 101–131. [17] M. Ulbrich, S. Ulbrich, L.N. Vicente, A globally convergent primal-dual interior-point filter method for nonconvex nonlinear programming, Mathematical Programming 100 (2004) 379– 410. [18] S. Ulbrich, On the superlinear local convergence of a filter-SQP methods, Mathematical Programming 100 (2004) 217–245. [19] Y.X. Yuan, On the convergence of a new trust region algorithm, Numerische Mathematik 70 (1995) 515–539.
A filter method for solving nonlinear complementarity problems Pu-yan Nie Department of Mathematics, College of Information Science and Technology, Jinan University, Guangzhou 510632, PR China
Abstract Filter methods are extensively studied to handle nonlinear programming problems recently. Because of good numerical results, filter techniques are attached importance to. In this paper, filter approaches are employed to tackle nonlinear complementarity problems (NCPs). Firstly, NCP conditions are transformed into a nonlinear programming problem. Then, to obtain a trial step, the corresponding nonlinear programming problems are solved by some existing strategies. Moreover, filter criterion is utilized to evaluate a trial iterate. The purpose in this paper is to employ filter approaches to attack NCPs. In essence, multi-objective view is utilized to attack NCPs because the idea of filter methods stems from multi-objective problems. Furthermore, a new filter method, based on the special two objects which differs from others, is brought forward. Moreover, Maratos effects are overcome in our new filter approach by weakening acceptable conditions. 2004 Elsevier Inc. All rights reserved. Keywords: Filter methods; Nonlinear complementarity problems; Nonlinear programming; Global convergence
E-mail address: [email protected] 0096-3003/$ - see front matter 2004 Elsevier Inc. All rights reserved. doi:10.1016/j.amc.2004.06.125
678
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
1. Introduction In this work, the system of nonlinear complementarity problem (NCP), with the following form, is considered F ðxÞ P 0 x P 0
and
xT F ðxÞ ¼ 0;
ð1:1Þ
where x 2 Rn, F : Rn ! Rn is a given function whose properties will be specified later and n is the dimension of the variables. The NCP has been utilized as a general framework for quadratic programming, linear complementarity, and the other mathematical programming. Different concepts have been developed to treat this problem. In the last few years, growing attention has been paid to approaches which employ a reformulation of NCP as a system of nonlinear equations or a minimization problem [4,9,11,12,15,16]. Theoretical results have been established to present conditions under which a stationary point of the minimization problem is a solution to the complementarity problem. Many researchers reformulate (1.1) as a system of nonsmooth equations. Then, many kinds of NCP functions appear to lead to a system of equations. We have the NCP functions with the form U(x, F(x)) and aim to find x satisfying Uðxi ; F i ðxÞÞ ¼ 0;
i ¼ 1; 2; . . . ; n;
ð1:2Þ
where U : Rn ! Rn and (1.2) is equivalent to (1.1) to a certain degree. Many algorithms based on (1.2) have been proposed. For example, NewtonÕs methods and generalized NewtonÕs approaches (see in [4,10,12]). To handle (1.1) by a minimization problem, the constraints are C ¼ fx P 0; F ðxÞ P 0g; and the objective function is kXF(x)k2 or something similar, where X = diag{x1,x2, . . ., xn}. Some optimization techniques can be employed. This type of method is re-visited in this paper. Moreover, the filter technique is employed to attack (1.1). In this work, filter methods are utilized to deal with NCPs. In the filter method, at each step, after a subproblem is solved, the filter criterion is employed to determine whether to accept the trial point or not. To give the subproblem, (1.1) is described as another equivalent form Xn minimize F 2 ðxÞx2i i¼1 i subject to
F j ðxÞ P 0; x P 0:
The algorithm is proposed, which are based on (1.3) to tackle (1.1).
ð1:3Þ
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
679
Filter approaches are introduced as follows. Filter method is initially put forward by Fletcher and Leyffer in 2002 [5]. Filter strategies do very well to balance the objective function and constraint conditions. Because of good numerical results, filter methods have been combined with trust region approaches, SQP techniques, pattern search method, interior point strategy and composite-like step methods [2,5–8,13,17]. In [18], super-linear local convergence is achieved for filter-SQP methods. In an optimization approach, we hope to find a satisfying point, which relates not only to objective function but also to constraint conditions. Two functions related close to constraint conditions and objective function are thus given as follows hðxÞ ¼
n X
maxf0; F j ðxÞg;
ð1:4Þ
j¼1
pðxÞ ¼
n X
2
x2i F i ðxÞ þ rhðxÞ;
ð1:5Þ
j¼1
where r is a constant. Moreover, x P 0 is always satisfied. At the kth point, if a new point xk + sk is obtained, a filter rule to judge k x + sk is presented. Certainly, we hope that p(xk + sk) < p(xi) or h(xk + sk) < h(xi) where i 2 Fk and Fk is a set consisted of ‘‘good’’ points. To obtain the convergent results, the strong conditions are required pðxk Þ 6 chðxk Þ þ pðxi Þ
ð1:6Þ
hðxk Þ 6 bhðxi Þ;
ð1:7Þ
or
where i 2 Fk and 1 > b > c > 0 are two constants. Certainly, whether a point is accepted by the filter is related to when it is generated (see [5] in detail). The filter is based on (1.6) and (1.7). xk + sk is accepted by the filter if (1.6) and (1.7) are satisfied for all i 2 Fk . Then, let xk+1 = xk + sk and update F. For convenience, denote Dk ðxkþ1 Þ ¼ fiji 2 Fk ; pðxi Þ P pðxkþ1 Þ and hðxi Þ P hðxkþ1 Þg: kþ1
Let Dk ðx Then,
Þ be the set of the points which is dominated by x
Fkþ1 ¼ Fk
[
fk þ 1g=Dk ðxkþ1 Þ:
k+1
ð1:8Þ k
:= x + sk. ð1:9Þ
We point out that the function of p(x) in this paper is different from those in the other papers. The choice of r is flexible. With suitable r < 0, Maratos effect can be overcome because the acceptable criterion is relaxed.
680
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
The paper is organized as follows: In Section 2, an algorithm is put forward. The algorithm is analyzed in Section 3. Some numerical results and remarks are listed in Section 4.
2. A filter algorithm for NCPs To handle (1.1) based on (1.3), we hope to find a better point at the beginning of the kth iteration. A step is obtained by solving the following subproblem " #T n X 1 2 T minimize WH ðxk F ðxk ÞÞi d k ðdÞ ¼ d Bk d þ r 2 i¼1 subject to
rF i ðxk Þd þ F i ðxk Þ P 0;
ð2:1Þ
kdk 6 D; d þ xk P 0; where F(xk) is a vector and Bk is a Hessian of Lagrangian function to (1.3) or an approximate Hessian matrix. The new point is defined as xk ðDÞ ¼ xk þ d k ; where dk is generated by (2.1). To determine whether to accept the trial step or not, some criterion is utilized. The following merit function is thus employed. mk ðxÞ ¼
n X 2 ððxF ðxÞÞi Þ:
ð2:2Þ
i¼1
The predicted reduction and actual reduction, to judge the new trial point, are presented as follows: Aredðsk Þ ¼ ðmk ðxk þ sk Þ mk ðxk ÞÞ; ð2:3Þ 1 Predðsk Þ ¼ ðsk ÞT Bk sk þ gTk sk ; ð2:4Þ 2 P where gk ¼ r ni¼1 ðxk F ðxk ÞÞ2i . When (2.1) is inconsistent, a restoration algorithm is employed to decrease the constraint violation degree and to obtain a better point by solving ($Fi(xk)sk + Fi(xk)) = 0 where ($Fi(xk)sk + Fi(xk)) = min{0, $Fi(xk)sk + Fi(xk)}. For convenience, the following notation is applied throughout this paper.
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
681
Definition 1 hIk ¼ minfhi jhi > 0; i ¼ 1; 2; . . . ; kg; pFk ¼ minfpi jhi ¼ 0; 1 6 i 6 kg; and pIk is the corresponding value to hIk . Namely, some way, to let (h(x 0 ), p(x 0 )) lying near (0, pF) or (hI, pI), is utilized. In the algorithm of this paper, the restoration algorithm is employed. For the update of Dk, the way of YuanÕs [19,14] is employed. The algorithm is formally stated as follows. Algorithm 1 (Trust Region Filter Algorithm For NCPs) Step 0. Choose D0, x0, r, a1, a2, b, g, g2, c, where D0 > 0, 0 < a1, a2, b, g, g2 < 1 and c 2 ð0; 12Þ; b 2 ð12 ; 1Þ. Set k :¼ 0; F ¼ fx0 g. Step 1. Compute hIk ; pIk ; pFk . Step 2. Solve (2.1) to get sk. Step 3. If (2.1) has no solution or sk = 0 but hk > 0, utilize Restoration Algorithm to get skr . Set xk :¼ xk þ skr , update hIk and goto Step 2. If sk = 0 and hk = 0 , then stop. Step 4. Compute rk ¼
Aredðsk Þ : Predðsk Þ
If rk 6 g, then, let xkþ1 :¼ xk ; Dkþ1 :¼ 12 Dk and k := k + 1 and go to Step 2. Compute ^ pk ¼ pðxk þ sk Þ; ^ hk ¼ hðxk þ sk Þ. If xk + sk is not acceptable to the filter, then Dkþ1 :¼ 12 Dk and goto Step 2. If hðxk Þ 6 Dk minfg2 ; a1 Dak 2 g, go to Step 5. If hðxk Þ P Dk minfg2 ; a1 Dak 2 g, call Restoration Algorithm to produce a point xkr :¼ xk þ skr such that: (A) xkr is acceptable to the filter; (B) hðxkr Þ 6 g2 minfhIk ; a1 D2k g. Let xk :¼ xkr and goto Step 2. Step 5. xk+1 := xk + sk is accepted by the filter. pkþ1 ¼ ^pk ; hkþ1 ¼ ^hk and remove the points dominated by (pk+1, hk+1) from the filter. Set Dk+1 := 2Dk and generate Bk+1. 2 Step 6. If hðxkþ1 Þ 6 Dkþ1 minfg2 ; a1 Dakþ1 g, then, set k := k + 1 and goto Step 1. Otherwise, k := k + 1 and call Restoration Algorithm to produce a point xkr ¼ xk þ skr satisfying (A) and (B). Let xk :¼ xkr and goto Step 2. If hðxk Þ > Dk minfg2 ; a1 Dak 2 g, a restoration algorithm is utilized to compute x such that (A) and (B) are all met. kr
682
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
In a restoration algorithm, we therefore aim to decrease the value of h(x). The direct way is employed Newton method or similar ways to solve F(x + sr) = 0. Let M jk ðdÞ ¼ ðhðxkj Þ hðxkj þ dÞÞ: A restoration algorithm is then given, where
Wjk ðdÞ ¼ kðrF ðxkj Þd þ F ðxkj Þ k1 hðxkj Þ;
ð2:5Þ
and rjk ¼
M jk ðskj Þ : Wjk ðskj Þ
ð2:6Þ
Algorithm 2 (Restoration Algorithm) Step 0. Let xk0 :¼ xk ; D0k :¼ Dk ; j :¼ 0; g; g2 2 ð0; 1Þ. Step 1. If (A) and (B) are met, then, let xkr :¼ xkj and stop. Step 2. Compute minimize
Wjk ðdÞ
subject to
kdk1 6 Dk
ð2:7Þ
to obtain skj . Calculate rjk . Step 3. If rjk 6 g, then, let xkjþ1 :¼ xkj ; Djþ1 :¼ 12 Djk ; j :¼ j þ 1 and goto Step 2. k jþ1 Step 4. Otherwise, let xkjþ1 :¼ xkj þ skj ; Dk :¼ 2Djk . Update Ajþ1 k ; j :¼ j þ 1 and goto Step 1. We point out that xkj is the jth restoration point following xk and xkr is the final restoration point following xk. The above Restoration Algorithm is a NewtonÕs method for F(x) = 0. This method is frequently utilized. Of course, there are other Restoration Algorithms, for example, the interior point restoration algorithm in [17] and the SLP restoration algorithm in [2]. We can also solve a quadratic programming to obtain the result. In fact, Algorithm 2 is a linear programming approach to reduce the constraint violation degree. Actually, Pred(sk) and Ared(sk) are permitted to be less than 0 at the same time. Therefore, this is a nonmonotonic algorithm. We hope the restoration algorithm to terminate finitely. In Section 3, this property will be shown. About r, it will be discussed in Section 4.
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
683
3. The convergence properties Just as that in [1,2,7,8], the analysis to the algorithm are based on the following standard assumptions. Further, to obtain the convergence, the sufficient reduction plays a crucial role throughout. Assumption 1 (1) The set {xk} 2 X is nonempty and bounded. (2) The function F(x) is twice continuously differentiable on an open set containing X. (3) When solving (2.7), we have
Wjk ðdÞ ¼ hðxkj Þ k½rF ðxkj Þd þ F ðxkj Þ k1 P b2 minfhðxkj Þ; Djk g; where b2 > 0 is a constant. (4) The matrix sequence {Bk} is bounded. (5) The Restoration Algorithm has a solution satisfying kskr k 6 s0 hk . (1) and (2) are the standard assumptions. (3) is the sufficient reduction condition, which is a very weak condition because Cauchy step satisfies this condition. It is regarded as a condition in this paper. In a trust region method, (3) guarantees the global convergence. (4) plays an important role to obtain the convergent result. But it has minor effects to the local convergent rate. The following results are based on Assumption 1. Analyzing the restoration algorithm, we obtain Lemma 1. The Restoration Algorithm 2 terminates finitely under Assumption 1. Proof. Assume the restoration algorithm does not terminate finitely. The termination criterion will be satisfied if limk!1 hðxkj Þ ¼ 0. It is thus reasonable that there exists > 0 with hðxkj Þ > for all j. We show it by contradiction. Let K ¼ fj : rjk > gg: According to Assumption 1 and (2.6), we have 1>
1 1 X X X k j1 ðhðxkj1 Þ hðxkj ÞÞ P g Wj1 Þ P gb2 minf; Djk g: k ðs j¼1
j¼1
k2K
Therefore, Djk ! 0; j 2 K. By the algorithm, we have Djk ! 0, for all j. Moreover,
684
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
Wjk ðskj Þ ¼ M jk ðskj Þ þ oðDjk Þ; P Djk from Step 4 of Restoration Algorithm. Namely, when Djk ! 0. Thus, Djþ1 k j fDk g is increased when Djk is very small, which contradicts Djk ! 0. Consequently, the result is obtained and the proof is complete. h The optimal properties of Algorithm 1 considered, the following result is obtained immediately. Lemma 2. Every new iteration xk+1 5 xk is acceptable to the filter set F. Proof. From Algorithm 1, a new iteration is produced in Steps 5 and 6. In both cases, xk+1 is accepted by the filter. The result therefore holds and the proof is complete. h Theorem 1. Under Assumption 1, suppose there are infinitely many points added to the filter. Then lim hðxk Þ ¼ 0:
k!1
Proof. See [7] or [13].
h
If there are finitely many points added to the filter, then, the following conclusions hold. Theorem 2. Under Assumption 1, suppose there are finitely many points added to the filter. Then hðxk Þ ¼ 0; for certain k1 and k > k1. Proof. The result is obvious from the algorithm.
h
Similar to [2,7], the global convergence theorem concerns Kuhn–Tucker (K–T) necessary conditions under a Mangasarian–Fromowitz constraint qualification (MFCQ). This is an extended form of Fritz–John conditions for a problem that includes equality constraints. We analyze the properties of solution to (2.1). Without loss of generality, assume that k$Fik 6 M, k$2Fik 6 M2 for all x 2 X and i = 1, 2, . . ., n. Then, we consider Fi(xk + sk). From Taylor expansion, if xk + sk is the solution to (2.1) there is T
T
F i ðxk þ sk Þ ¼ F i ðxk Þ þ rF i ðxk Þ sk þ ðsk Þ r2 F i ðxk þ hðxk þ sk xk ÞÞsk ;
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
685
where h 2 (0, 1). When Dk is small enough and Fi(xk) = min{0, Fi(xk)}, we have Fi(xk) + $Fi(xk)Tsk P 0 for i = 1, 2, . . ., n. We therefore obtain n X kðsk ÞT r2 F i ðxk þ hðxk þ sk xk ÞÞsk k hðxk þ sk Þ 6 i¼1
6 nkr2 F i ðxk þ hðxk þ sk xk ÞÞkD2k 6 2n2 M 2 D2k :
ð3:1Þ
Lemma 3. If the (k + j)th step relates to Restoration Algorithm, j is large enough and Dk ! 0, then, 1þa2 2 a1 D1þa kþj 6 hkþj 6 4a1 Dkþj :
ð3:2Þ
Proof. The first part of the inequality is obvious by the algorithm. Consider the second part of the inequality. If xk+j1 is refused by the filter, then Dk+j1 6 2Dk+j, xk+j = xk+j1 and there is no Restoration Algorithm at the (k + j 1)th 1 step. Accordingly, the result is true. If xk+j1 is accepted and Dkþj 6 ð4a21 Þ1a2 , according to Step 5 of Algorithm 1, we have Dk+j1 6 Dk+j and 2 hkþj 6 2D2kþj1 6 2D2kþj 6 4a1 D1þa kþj :
The result is established.
h
About (2.1), Fletcher et al. do very excellent analysis in [7], which is listed as follows Theorem 3. Let Assumption 1 hold and xw 2 X be a feasible point of problem (1.3) at which MFCQ holds, but which is not a K–T point. Then there exists a neighborhood N0 of xw and positive constants n1, n2, n3 such that for all x 2 N0 ˙ X and all Dk for which n2 hðcðxÞÞ 6 Dk 6 n3 ;
ð3:3Þ
it follows that (3.2) has a feasible solution d at which the predicted reduction satisfies 1 Wk ðdÞ P n1 Dk : 3 If Dk 6 (1 g3)n1/3nM2, then n X 2 2 ððxk Þi F 2i ðxk Þ ðxk þ dÞi F 2i ðxk þ dÞÞ P g3 Wk ðdÞ; i¼1
where g < g3. qffiffiffiffiffiffiffiffi k k If hk > 0 and Dk 6 n2bh 2 M , then, h(x + d) 6 bhk. 2
ð3:4Þ
ð3:5Þ
686
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
Proof. See Lemmas 4 and 5 in [7]. h Lemma 4. Assume that Algorithm 1 generates a sequence {xk} and its accumulation xw satisfying the conditions of Theorem 3, if there are infinitely many points enter into the filter, then, there are only finitely many points related with Restoration Algorithm. Proof. Under conditions of Theorem 3, if there are infinitely many points related to Restoration Algorithm, Dk ! 0. (See Theorem 7 in this paper.) According to Theorem 3, consider the conditions when a point is accepted. If r 5 0, we show the result. When ( 1 sffiffiffiffiffiffiffiffiffiffiffi) ð1 g3 Þn1 g 3 n 1 a2 2bhk Dk 6 min n3 ; ; ; ; ð3:6Þ 3nM 2 3a1 jrj n2 M 2 Dk P
1 1þa1 2 h ; a1 k
ð3:7Þ
and hk > 0, then, (3.4) and (3.5) are true. Moreover, according to Theorem 3, we obtain h(xk + sk) 6 bhk. xk + sk will be therefore accepted by the filter because n n X X ðxk Þ2i F 2i ðxk Þ ðxk þ sk Þ2i F 2i ðxk þ sk Þ pðxk Þ pðxk þ sk Þ ¼ i¼1
i¼1
þ rðhk hðxk þ sk ÞÞ 1 2 >0 P g3 Wk ðsk Þ jrjhk P g3 n1 Dk jrja1 D1þa k 3 and h(xk + sk) 6 bhk. By the way, the interval between (3.6) and (3.7) 1is not empty if D is small 1þa
enough. Further, under (3.6) and (3.7), Dkþj P a11 hk 2 is satisfied for j P 1 if k is large enough. Namely, every trial step satisfies (3.7) for j P 1. If there is a Restoration Algorithm, from (3.4), there exists a constant 13 n1 P 0 > 0 such that kþj WH Þ P 0 kskþj k: kþj ðs
From (3.2) and the above analysis, we have kþj pðxkþj Þ pðxkþj þ skþj Þ P g3 WH Þ þ rhðxkþj Þ rhðxkþj þ skþj Þ kþj ðs kþj P g3 WH Þ jrjhðxkþj Þ kþj ðs 1þa2 2 2 P g3 0 Dkþj jrja1 D1þa kþj P 2cn M 2 Dkþj
P 2cn2 M 2 D2kþj P chðxkþj ðDkþj ÞÞ;
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
687
g3 0 k+j r 2 if Dakþj 6 jrja1 þ2cn )k 2 M . The above inequality holds because limk!1k(s 2
6 limk!1s0hk+j = 0 and Dk ! 0. The above inequality means that p is monotonic 1 1þa decreasing and there exists a trial point accepted by the filter for Dkþj P a11 hk 2 and 1 sffiffiffiffiffiffiffiffiffiffiffi) ð1 g3 Þn1 g3 n1 a 2 2bhk 6 min n3 ; ; ; : 3nM 2 3a1 jrj n2 M 2 (
Dkþj
1 1þa
Therefore, Dkþjþ1 P a11 hk 2 . We now show that (3.7) is reasonable. If (3.7) were always false, there 1 1þa2
should be Dk < a11 hk
and limk!1Dk = 0 according to limk!1hk = 0. Further1
more, from (3.2) when 2Dk+1 P Dk and Dkþ1 6 ð8na2 1M 2 Þ1þa2 for sufficient large k, hkþ1 6 2n2 M 2 D2k 6 8n2 M 2 D2kþ1 : 1 1þa
1 1þa2
Thus, a11 hkþ12 6 Dkþ1 , which contradicts Dk < a11 hk If r = 0, (3.6) is substituted for (
ð1 g3 Þn1 ; Dk 6 min n3 ; 3nM 2
because of limk!1Dk = 0.
sffiffiffiffiffiffiffiffiffiffiffi) 2bhk : n2 M 2
Similar result is obtained. The result is accordingly obtained and the proof is complete. h We point out that the result in Lemma 4 holds for any bounded constant r. The global convergence can be obtained combining the results of Theorems 1–3 and Lemma 4. Theorem 4. Let Assumption 1 hold and let the sequence {xk}, where {xk} 2 X is generated by Algorithm 1, converges to a feasible point of problem (1.3) at which MFCQ holds. Then, {xk} has an accumulation which is a K–T point of (1.3). Proof. If Algorithm 1 terminates finitely, the result holds apparently. Assume Algorithm 1 terminates infinitely. If the result were false, from Theorem 1 and the above analysis, there would exist an integer k0 such that hk < n2 is small H enough, mk ðxk Þ mk ðxk þ dÞ P g½WH k ð0Þ Wk ðdÞ for all k > k0. Furthermore, the Restoration Algorithm does not appear when k > k0 from Lemma 4. For convenience, we denote H K 2 ¼ fkjk > k 0 ; and mk ðxk Þ mk ðxk þ dÞ P g½WH k ð0Þ Wk ðdÞg:
ð3:8Þ
688
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
Moreover, from Assumption 1 there exists X X 1> ½pðxk Þ pðxkþ1 Þ P ½g3 Wk ðdÞ þ rðhk hkþ1 Þ K2
P
X 1 K2
3
K2
g3 n1 Dk þ rðhk hkþ1 Þ :
Consequently, limk!1Dk = 0 for k 2 K2 because hk ! 0. Taking into account the choice of Dk, limk!1Dk = 0 for all k. There thus exists Dk such that 2 Da > k
6 ðjrj þ maxfð1 bÞ; cgÞ g3 n1 a1
ð3:9Þ
for sufficiently large k and 2 P hðxk þ sk Þ: a1 D1þa k
ð3:10Þ
(3.10) is issued from Step 4 of Algorithm 1. We consider the next trial step xk + sk. By the above analysis, p(xk) is monotonic decreasing and 1 pðxk Þ pðxk þ sk Þ P g3 n1 Dk þ rðhk hðxk þ sk ÞÞ 3 1 1þa2 2 P g3 n1 Da þ rðhk hðxk þ sk ÞÞ k Dk 3 2 P ðjrj þ maxfð1 bÞ; cgÞ2a1 D1þa þ rðhk hðxk þ sk ÞÞ k P maxfð1 bÞ; cghðxk þ sk Þ; k
ð3:11Þ
k
which means that x + s will be accepted by the filter. Dk will not be reduced when Dk is small enough and k > k0, which contradicts limk!1Dk = 0. The result thus holds and the proof is complete. h The K–T points of (1.3) but not the solutions to (1.1) are defined as local infeasibility points as follows Definition 2. If x* is a K–T point but not the solution to (1.1), we call x* a local infeasibility point. In fact, there may be no solution or we may fail to find the solution to (1.1) in some region. For global case, it is very difficult to describe it just as the other global problems. Local infeasibility points are therefore defined. The following result is obviously obtained from Theorem 4. Theorem 5. Under assumptions of Theorem 4, the algorithm termites finitely or has an accumulation point which is the solution of (1.1) or a local infeasibility point.
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
689
Under certain conditions, the stationary points of (1.3) are the solutions of (1.1). For example, if the solution to (1.3) is zero, it is the solution of (1.1). The global minimization is just the solution of (1.1) if there exists. We hope that we can obtain the global minimum to (1.3). We give a condition to guarantee this result. Definition 3. The mapping F is a P0 function, i.e., for every x1, x2 2 Rn with x1 5 x2, there exists an index i 2 {1, 2, . . ., n} such that x1i 6¼ x2i
and ðx1i x2i ÞðF i ðx1 Þ F i ðx2 ÞÞ P 0:
Correspondingly, we define the P0-matrix. The gradient of P0 function is a P0matrix. M 2 Rn·n is said to be a P0-matrix if all the principal minors are nonnegative. P0-matrix has the following properties. Lemma 5. M 2 Rn·n, the following two statements are equivalent (1) M is a P0-matrix. (2) For any x 2 Rn 5 0, there exists an index xi such that xi(Mx)i P 0.
Proof. See Theorem 3.4.2 in [3].
h
Now, the result is presented as follows. Theorem 6. If F is a P0 function, then xw is a global minimum if and only if it is a K–T point. Proof. The ‘‘only if’’ part is obvious. We show the ‘‘if’’ part. If xw is not the global minimum to (1.3), without loss of generality, assume that xH S > 0
and
F S ðxH Þ > 0;
where S is a subset of {1, 2, . . ., n}. According to K–T conditions, we thus have kS = 0, lS = 0, where kS and lS are the corresponding Lagrangian to the constraints. Further, from K–T conditions there exists X 2 H H H 2 H ðxH ð3:12Þ i Þ F i ðx Þrxj F i ðx Þ þ xj F j ðx Þ ¼ 0 i2S 2 H H 2 H for any j 2 S. Denoting ui ¼ ðxH k i Þ F k i ðx Þ > 0 and vi ¼ xk i F k i ðx Þ where ki 2 S, from (3.12), we have
rs Fu þ v ¼ 0:
690
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
By v > 0, we have $sFu < 0. This contradicts the fact that $sF is a P0-matrix issuing from Lemma 5. (The sub-matrix of a P0-matrix is also a P0-matrix.) The result is consequently established. h But when there is no solution to (1.1), some other algorithms may run slowly near a local infeasibility point while the algorithm may avoid it because its second information is made full use of. Now, the following result, which is employed in Lemma 4, is proved. Theorem 7. If there are infinitely many points related to Restoration Algorithm and they are far away from K–T points, then, limk!1Dk = 0. Proof. Apparently there is a subsequence of {Dk} satisfying limj!1 Dtj ¼ 0 because hk ! 0 and there are infinitely many points related to Restoration Algorithm. If xk is far away from K–T points, from the result in [19], there exist k two positive constants s and such that WH s minfD0 ; g. Further, k ðs Þ P when (
Dkþj
ð1 g3 Þn1 6 min n3 ; ; 3nM 2
sffiffiffiffiffiffiffiffiffiffiffi) 2bhk ; n2 M 2
this trial step will be accepted by filter in terms of Theorem 3. It is therefore impossible to have any subsequence of {Dk} satisfying limj!1 Dtj ¼ 0, which contracts the above fact. Therefore, limk!1Dk = 0 and the proof is complete. h The global convergence properties have been obtained. Actually, n0 is a very small real. Algorithm 1 is an SQP trust-region method combined with a new filter technique for NCPs. The approach in this paper has the following advantages: • This method supports a new approach to get the solution of (1.1). In early filter techniques, objective function and constraint violation degree function are reduced. While in the algorithm in this paper, we aim to reduce the objective and some other function, which is a combination of constraint violation and objective function. When (1.1) has no solution, a technique to terminate is proposed. • It is a method without multipliers and it is very difficult to deal with multipliers. • This approach is a nonmonotonic one, which may be beneficial to obtain global solutions. • The subproblem employed in this paper is based on a smoothing problem.
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
691
• The method in this paper has no requirement about feasible initial point, which is convenient.
4. Conclusion remarks and numerical results The researchers concentrate on filter techniques because filter approaches can efficiently balance the objective function and constraint violation function. The combination of various objective function with filter technique brings out a new method. Of course, diverse kinds of p can be used in various methods. As for r, it can be positive or negative because the filter method can guarantee hk ! 0. Therefore, it has no effect on the convergence of the algorithm. In the other papers about filter method, to attack (1.3), authors assume that r = 0. When r > 0, the accepting condition becomes strict. But, when r < 0 the accepting condition becomes moderate. Our filter approach is therefore regarded as a general version about the choices of r in filter method family. When a suitable r is employed, some bad cases may be overcome. Bk is an approximate Hessian. Some approaches, such as BFGS, can be applied to update. It has minor effect on the convergent properties. We also do not discuss the initial trust region radius, although it plays an important role to obtain good numerical results. In Assumption 1, it is required that {Bk} be bounded. Actually, when kBkk 6 qk where q is a positive constant, the convergent results also hold. We hope that we can moderate the conditions. But MFCQ seems necessary. Some examples illustrate the algorithm. (In the following examples, the toleration is 1.0e7 and b = 0.05, c = 0.98, a1 = a2 = 0.5, g = 0.25, g2 = 0.75). Our algorithm compares with that in [10]. The toleration is 1.0e6 in [10]. Example 1. Let F(x) = Mx + q, where 2
4
1
6 6 1 4 6 6 6 6 0 1 6 M ¼6 6 6 6 6 6 0 0 4 0
0
0
1 4
0
0
and q = (1, 1, . . . , 1)T.
0
3
7 0 7 7 7 7 0 7 7 7 7 7 7 7 1 7 5 4
692
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
The results of Example 1 lie in the following table, where the initial point is x0 = (0,0, . . ., 0)T.
NF NG IT(JQÕs)
n=4
n=8
n = 16
n = 32
n = 64
n = 128
n = 256
2 2 4
2 2 4
2 2 4
2 2 4
2 2 4
2 2 4
2 2 4
where NF and NG mean the number of calculations of f and gradient, respectively. IT (JQÕs) means the number of iterates in [10]. This problem has no relation with r. But it works very well. Example 2. Let F 1 ðxÞ ¼ 3x21 þ 2x1 x2 þ 2x22 þ x3 þ 3x4 6; F 2 ðxÞ ¼ 2x21 þ x1 þ x22 þ 10x3 þ 2x4 2; F 3 ðxÞ ¼ 3x21 þ x1 x2 þ 2x22 þ 2x3 þ 9x4 9; F 4 ðxÞ ¼ x21 þ 3x22 þ 2x3 þ 3x4 3:
pffiffi This example has one degenerate solution ð 26 ; 0; 0; 12Þ and one nondegenerate solution (1, 0, 3, 0). The results of Example 2 is listed in the following table if we let r = 0.
Intial point
NF
Ng
IT (JQÕs)
(0, 0, 0) (1, 1, 1)
12 16
8 10
13 8
Example 2 appears in [12]. The following example and Example 1 come from [9]. Example 3. Let F(x) = Mx + q, where M ¼ diag
1 2 ; ;...;1 n n
and q = (1, 1, . . ., 1)T.
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
693
The initial point is x0 = (0, 0, . . ., 0)T. We give the results.
n = 10 n = 20 n = 40 n = 80 n = 160 n = 320 No trust region restrict in RA
NF
3
3
3
3
3
3
NG Initial trust region NF radius is 1-norm of gradient NG IT(JQÕs)
2 2
2 2
2 2
2 2
2 2
2 2
2 4
2 4
2 4
2 4
2 4
2 4
From Examples 1–3, the algorithm do very well compared with other methods to deal with NCP function. We let the constant r be flexible. Certainly, we can give some criterions to choose r or to update r. But, the monotonic property of r is not necessary if r is allowed for change. The value of r is not required to be very large. In a restoration algorithm, the instinctive thought is to minimize h(x). To achieve the convergence properties of Algorithm 1, the sufficient reduction condition is necessary because both in filter method and in trust region approach the sufficient reduction condition is very important. the algorithm, the merit PIn n k 2 2 k function is necessary to enable the information of P ðx Þi F i ðx Þ to keep in fit i¼1 n k with the quadratic model. Moreover, we can use ðx Þi F i ðxk Þ instead of i¼1 Pn k 2 2 k i¼1 ðx Þi F i ðx Þ. The similar results can be obtained. About the parameters b and c, where 1 > b > c > 0. b is required to close enough to 1, while c sufficiently close to 0. There are several advantages: (1) The acceptable criterion is not very strong. (2) The relation of p and h becomes more implicit. Certainly, some conditions in the algorithm can be modified to get the global convergence. In the algorithm, we hope to avoid the linkage between p and h. But it is by no means yet. Because when the relation is lost, the convergence properties are not guaranteed. In this way will a nondescent direction be produced. As for the norm of the constraint violation function, the other forms can be chosen, such as l1-norm and l2-norm. For example, l1-norm is used in [1]. About the filter method, there is much work to do. For example, filter methods with three dimensions are interesting. Meanwhile, the choice of alternative for the optimal measure is still an interesting topic. In the numerical results, when this technique is taken steps, it plays better to a certain degree. It is well known, filter method has the advantage of high-efficiency.
694
P.-y. Nie / Appl. Math. Comput. 167 (2005) 677–694
Acknowledgment The Project Sponsored by the Scientific Research Foundation for the Returned Overseas Chinese Scholars, State Education Ministry.
References [1] C. Audet, J.E. Dennis Jr., Combining patter search and filter algorithms for derivative free optimization, SIAM Journal on Optimization 14 (2004) 980–1010. [2] C.M. Chin, R. Fletcher, On the global convergence of an SLP-filter algorithm that takes EQP steps, Mathematical Programming 96 (2003) 161–177. [3] R.W. Cottle, J.S. Pang, R.E. Stone, The Linear Complementarity Problem, Academic Press, New York, NY, 1992. [4] A. Fischer, Solution of monotone complementarity problems with locally Lipschitzian functions, Mathematical Programming 76 (1997) 513–532. [5] R. Fletcher, S. Leyffer, Nonlinear programming without a penalty function, Mathematical Programming 91 (2002) 239–269. [6] R. Fletcher, N.I.M. Gould, S. Leyffer et al., Global convergence of a trust-region SQP-filter algorithm for general nonlinear, SIAM Journal on Optimization 13 (2003) 635–659. [7] R. Fletcher, S. Leyffer, P.L. Toint, On the global convergence of a filter SQP algorithm, SIAM Journal on Optimization 13 (2002) 44–59. [8] R. Fletcher, S. Leyffer, A bundle filter method for nonsmooth nonlinear optimization, Tech. Report NA/195, Department of Mathematics, University of Dundee, 1999. [9] C. Geiger, C. Kanzow, On the solution of monotone complementarity problems, Computational Optimization and Application 5 (1996) 155–173. [10] H.Y. Jiang, L.Q. Qi, A new nonsmooth equations approach to nonlinear complementarity problems, SIAM Journal on Control and Optimization 35 (1997) 178–193. [11] C. Kanzow, Some equation-based methods for the nonlinear complementarity problem, Optimization Methods Software 3 (1994) 327–340. [12] J.J. More´, Global methods for nonlinear complementarity problems, Mathematics of Operations Research 21 (1996) 589–614. [13] P.Y. Nie, Composite-step like filter methods for equality constraint problems, Journal of Computational Mathematics 21 (2003) 613–624. [14] P.Y. Nie, A null space method for system of nonlinear equations, Applied Mathematics and Computation 149 (2004) 215–226. [15] J.S. Pang, S.A. Gabriel, NE/SQP: A robust algorithm for nonlinear complementarity problem, Mathematical Programming 60 (1993) 295–337. [16] J.S. Pang, A B-differentiable equations based, globally and locally quadratically convergent algorithm for nonlinear programming, complementarity, and variational inequality problems, Mathematical Programming 51 (1991) 101–131. [17] M. Ulbrich, S. Ulbrich, L.N. Vicente, A globally convergent primal-dual interior-point filter method for nonconvex nonlinear programming, Mathematical Programming 100 (2004) 379– 410. [18] S. Ulbrich, On the superlinear local convergence of a filter-SQP methods, Mathematical Programming 100 (2004) 217–245. [19] Y.X. Yuan, On the convergence of a new trust region algorithm, Numerische Mathematik 70 (1995) 515–539.