**Department of Mathematics**

Abstract

In this paper, we present a neural network approach for solving nonlinear complementarity problems. The neural network model is derived from an unconstrained minimization reformulation of the complementarity problem. The existence and the convergence of the trajectory of the neural network are addressed in detail. In addition, we also explore the stability properties, such as the stability in the sense of Lyapunov, the asymptotic stability and the exponential stability, for the neural network model. The theory developed here is also valid for neural network models derived from a number of reformulation methods for nonlinear complementarity problems. Simulation results are also reported. c 2001 Elsevier Science B.V. All rights reserved. MSC: 90C33; 90C30; 65H10 Keywords: Neural network; Nonlinear complementarity problem; Stability; Reformulation

1. Introduction

We are interested in Ã¿nding a solution of the nonlinear complementarity problem (NCP(F)) Fi (x)¿0; xi ¿0 and xi Fi (x) = 0 i = 1; : : : ; n; where F : Rn Rn is assumed to be continuously di erentiable. When F(x) = Mx + q for some M Rn×n and q Rn , NCP(F) is reduced to the linear complementarity problem (LCP(M; q)). There are rich sources for complementarity problems, see the book [9] and a recent review paper [14].

This work is supported in part by grant FRG=97-98=II-42 of Hong Kong Baptist University and the Australian Research Council. Corresponding author. E-mail addresses: liliao@hkbu.edu.hk (L.-Z. Liao), hdqi@maths.unsw.edu.au (H. Qi), l.qi@unsw.edu.au (L. Qi). 0377-0427/01/$ - see front matter c 2001 Elsevier Science B.V. All rights reserved. PII: S 0 3 7 7 - 0 4 2 7 ( 0 0 ) 0 0 2 6 2 - 4

344

In the past decade, there have been growing interests in solving NCP(F) via various unconstrained minimization methods, among which are so called reformulation methods see [17]. In this paper, we will pay a particular attention to a reformulation approach based on the FischerBurmeister NCP function. The introduction of neural networks (or artiÃ¿cial neural networks) in optimization started in 1980s (see [5,19]). Since then, signiÃ¿cant research results have been achieved for various optimization problems, such as linear programming [35], quadratic programming [1], linear complementarity problems [23], and nonlinear programming [30]. The essence of neural network approach for optimization is to establish an energy function (nonnegative) and a dynamic system which is a representation of an (artiÃ¿cial) neural network. The dynamic system is normally in the form of Ã¿rst-order ordinary differential equations. It is expected that for an initial state, the dynamic system will approach its static state (or equilibrium point) which corresponds the solution of the underlying optimization problem. An important requirement is that the energy function decreases monotonically as the dynamic system approaches an equilibrium point. In this paper, following a reformulation of NCP(F), an unconstrained optimization problem is formulated. As a result, we are able to establish an energy function and a neural network. However, the function deÃ¿ning the dynamic system is continuous only, not necessarily di erentiable. Hence classical results of stabilities on di erentiable case are not applicable to our dynamic system. Despite of these di culties, we could establish the existence, the convergence of the trajectory, as well as some stability results such as the stability in the sense of Lyapunov and the asymptotic stability. Especially, we establish the relationship between the exponential stability and the regularity condition of NCP(F). In the previous study [23] on LCP(M; q), we restrict M to be both P0 matrix and R0 matrix so that the trajectory is guaranteed to be convergent to a static state. We do not extend this result to its nonlinear correspondence here because the argument is similar to the linear case [23]. The rest of the paper is organized as follows. In Section 2, some preliminary results are provided. Section 3 is devoted to the neural network architecture and implementation of our new method. Convergence and stability results are discussed in Section 4. Simulation results of the new method are reported in Section 5. Finally, some concluding remarks are drawn in Section 6. 2. Preliminaries In this section, we recall some basic classes of functions and matrices, the properties of some special NCP functions, as well as some stability concepts in di erential equations. We also present some related results which will be used later. Some of these results are only discovered recently. 2.1. Functions and matrices Throughout the paper, we assume that F : Rn Rn is a continuously di erentiable function. Let M be a matrix in Rn×n . DeÃ¿nition 2.1. A matrix M Rn×n is said to be · a P0 matrix if all of its principal minors are nonnegative; · a P matrix if all of its principal minors are positive; positive semideÃ¿nite if x; Mx ¿0 for all x Rn , and positive deÃ¿nite if x; Mx ¿ 0 for all x (= 0) Rn . Obviously, a positive-deÃ¿nite matrix is a P matrix, and a positive-semideÃ¿nite matrix is a P0 matrix. Let N = {1; : : : ; n}. DeÃ¿nition 2.2. The function F : Rn Rn is said to be a · P0 -function if for all x; y Rn with x = y max (xi - yi )[Fi (x) - Fi (y)]¿0;

iN xi =yi· P-function if for all x; y Rn with x = y max(xi - yi )[Fi (x) - Fi (y)] ¿ 0;

iN · uniform P-function if there exists a constant max(xi - yi )[Fi (x) - Fi (y)]¿

iN 0 such that for all x; y Rn x - y 2; · monotone function if for all x; y Rn (x - y)T (F(x) - F(y))¿0: It is known that F is a P0 function if and only if the Jacobian matrix F (x) is a P0 matrix for all x Rn , and if F (x) is a P matrix for all x Rn then F must be a P function. 2.2. The FischerBurmeister NCP function The FischerBurmeister function (a; b) = a2 + b2 - a - b: : R2 R is given by

The function was Ã¿rst used by Fischer [15] to construct a Newton-type method for constrained optimization problem, and later was extensively used to solve the nonlinear complementarity problem, see [10,12,13,16,20,28,33]. Here are some properties of this function. These results can be found in the papers just mentioned. For the deÃ¿nition of semismoothness and strong semismoothness, one can refer to the paper [29] by Qi and Sun. Proposition 2.3. (i) (a; b) = 0 if and only if a¿0; b¿0; ab = 0. (ii) The square of (a; b) is continuously di erentiable. (iii) is twice continuously di erentiable everywhere except at the origin; but it is strongly semismooth at the origin.

Usually a function satisfying (i) in Proposition 2.3 is called an NCP function. By this property, NCP(F) can be equivalently reformulated as Ã¿nding a solution of the following equation: (x1 ; F1 (x)) . . (x) = = 0: . (x n ; Fn (x) (1)

We note that is locally Lipschitz continuous everywhere, so that Clarke's [8] generalized Jacobian @ (x) is well deÃ¿ned at any point. To ease our later presentation, we deÃ¿ne E(x) =

1 2 (x) 2 :

We note that E(x)¿0 for all x Rn and x solves NCP(F) if and only if E(x) = 0. Hence solving NCP(F) is equivalent to Ã¿nding the global minimizer of the following unconstrained minimization problem min E(x): n xR

(2) , see [26,29,25]; The Ã¿rst result in the next proposition follows directly from the semismoothness of and the proofs of the last two can be found in [10]. Proposition 2.4. (i) For any x Rn ; we have (x + d) - (x) - Vd = o( d ) for d 0 and V @ (x + d):

(ii) The function E is continuously di erentiable with E(x) = V T (x) for an arbitrary element V @ (x). (iii) If F is a P0 function; then any stationary point x of optimization problem (2); i.e.; E(x)=0; is a solution to NCP(F). In the statement of Proposition 2.4, we make use of the Landau symbol o(·): If { k } and {Ã¿k } are two sequences of positive numbers, then k = o(Ã¿k ) if limk k =Ã¿k = 0 for Ã¿k 0. It is well known that linear programming and convex quadratic programming problems can be equivalently reformulated as monotone linear complementarity problems, and the general convex programming can be reformulated as a monotone nonlinear complementarity problem. Hence the class of P0 NCP (where the function F is a P0 function) covers a large range of problems, as well as of the problems arising from economy and engineering [14,25]. 2.3. Stability in di erential equations Now we recall some stability results from [34] on the following di erential equation: x(t) = f(x(t)); x(t0 ) = x0 Rn : (3) The following classical results on the existence and uniqueness of the solution to (3) hold. Theorem 2.5. Assume that f is a continuous mapping from Rn to Rn . Then for arbitrary t0 ¿0 and x0 Rn there exists a local solution x(t); t [t0 ; ) to (3) for some t0 . If furthermore f is locally Lipschitzian continuous at x0 then the solution is unique; and if f is Lipschitzian continuous in Rn then can be extended to +. If a local solution deÃ¿ned on [t0 ; ) cannot be extended to a local solution on a larger interval [t0 ; 1 ), 1 ¿ , then it is called a maximal solution, and the interval [t0 ; ) is the maximal interval of existence. An arbitrary local solution has an extension to a maximal one. The maximal interval of existence associated with x0 is often denoted by [t0 ; (x0 )). Theorem 2.6. Assume that f is a continuous mapping from Rn to Rn . If x(t); t [t0 ; (x0 )); is a maximal solution and (x0 ) ¡ + then

t (x0 ) lim x(t) = +:

The results stated in Theorems 2.5 and 2.6 can be found in [34, p. 74, Theorems 1:1 and 1:2]. A point x Rn is called an equilibrium point of (3) if f(x ) = 0. The following are some common deÃ¿nitions on stability. DeÃ¿nition 2.7 (Stability in the sense of Lyapunov). Let x(t) be a solution of (3). An isolated equilibrium point x is Lyapunov stable if for any x0 = x(t0 ) and any scalar ¿ 0 there exists a ¿ 0 so that if x(t0 ) - x ¡ then x(t) - x ¡ for t¿t0 . DeÃ¿nition 2.8 (Asymptotic stability). An isolated equilibrium point x is said to be asymptotic stable if in addition to being a Lyapunov stable it has the property that x(t) x as t +, if x(t0 ) - x ¡ : DeÃ¿nition 2.9 (Lyapunov function). Let Rn be an open neighborhood of x. A continuously di erentiable function E : Rn R is said to be a Lyapunov function at the state x (over the set ) for Eq. (3) if E(x) = 0; E(x) ¿ 0 for x ; x = x; dE(x(t)) = [x(t) E(x(t))]T f(x(t))60; x : dt (4)

A Lyapunov function is often called an energy function for (3). The next result addresses the relationship between stabilities and a Lyapunov function, see [34,24,4]. Theorem 2.10. (i) An isolated equilibrium point x is Lyapunov stable if there exists a Lyapunov function over some neighborhood of x . (ii) An isolated equilibrium point x is asymptotically stable if there exists a Lyapunov function over some neighborhood of of x satisfying dE(x(t)) ¡ 0 x(t) ; x(t) = x : dt A strong notion than the Lyapunov stability is the so called exponential stability.

DeÃ¿nition 2.11. An isolated equilibrium point x is exponentially stable for (3) if there exist ! ¡ 0; Ã„ ¿ 0; ¿ 0 such that arbitrary solution x(t) of (3), with the initial condition x(t0 ) = x0 , x(t0 ) - x ¡ , is deÃ¿ned on [0; ) and satisÃ¿es x(t) - x 6Ã„e!t x(t0 ) - x ; t¿t0 :

It is clear that exponentially stable equilibria are asymptotic stable. In the next section we will propose a dynamic system for NCP(F) and study its stabilities. 3. Neural network model As we discussed in Section 2, NCP(F) can be reformulated as an unconstrained minimization problem (2). The objective function E is continuously di erentiable for all x Rn . Hence it is natural to use the steepest descent-based neural network model for problem (2) d x(t) (5) = - E(x(t)); x(0) = x0 ; ¿ 0; dt where is a scaling factor. ¿ 1 indicates that a longer step could be taken. For simplicity of our analysis, we let = 1. From Proposition 2.4(ii), E(x) can be easily obtained. Since we are considering the general NCP(F), the computation of E(x(t)) will determine how the neural network (5) can be implemented on hardware. The discussion in [6], especially Appendix B in [6], indicates that for certain very nonlinear problems, computer assisted neural network implementation of (5) might be necessary. Fig. 1 provides an indication of how neural network (5) would be implemented on hardware. A procedure to evaluate an element V @ (x) is provided in the following, for a proof see [10]. Algorithm 3.1 (Procedure to evaluate an element V @ (x)). (S:0) Let x Rn be given; and let Vi denote the ith row of a matrix V Rn×n .

(S:1) Set S1 = {i | xi = 0 = Fi (x)}. (S:2) Set z Rn such that zi = 0 for i S1 and zi = 1 for i S1 . (S:3) For i S1 ; set Vi = zi - 1 eT + i (zi ; Fi (x)T z) xi - 1 eT + i (xi ; Fi (x)) Fi (x)T z - 1 Fi (x)T : (zi ; Fi (x)T z)

(S:4) For i S1 ; set Vi = Fi (x) - 1 Fi (x)T : (xi ; Fi (x))

There are other ways to evaluate elements in @ (x), for example see [28]. We stress that any element in @ (x) enjoys the following structure: V = Da (x) + Db (x)F (x); where Da (x) and Db (x) are diagonal matrices, which can be calculated explicitly, see [10] for details.

4. Stability analysis In this section, we address the stability issues on the neural network (5) in two aspects. First we investigate the general behavior of the solution trajectory of (5). This study includes the uniqueness, the convergence, and other various properties of the trajectory. Then we focus on a particular case where the equilibrium point is isolated. 4.1. Existence of the trajectory Let S denote the solution set of NCP(F) and let x S. Then it is easy to see consequently E(x) = 0 (Proposition 2.4(ii)). Hence we have (x) = 0,

Proposition 4.1. Every solution to NCP(F) is an equilibrium point of the neural network (5). Conversely; if x Rn is an equilibrium of (5) and F(x) is a P0 function; then x S. Proof. We only need to address the second part of the proposition. In optimization, an equilibrium point x is called a stationary point of (2). It is well known [10] that F(x) being a P0 function is a su cient condition for a stationary point to be a solution to NCP(F). Let L(x0 ) denote the level set associated with the initial state x0 and be given by L(x0 ) = {x Rn | E(x)6E(x0 )}: The next result addresses the existence of the solution trajectory of (5).

Theorem 4.2. (i) For an arbitrary initial state x0 ; there exists exactly one maximal solution x(t); t [t0 ; (x0 )) of the Eq. (5). (ii) If the level set L(x0 ) is bounded or F(x) is Lipschitzian continuous; then (x0 ) = +. Proof. (i) It is known that E(x) is continuous (not necessarily di erentiable). Indeed it is locally Lipschitzian continuous by following, for example, a similar argument in [18]. Hence Theorem 2.5 implies the maximal solution is unique. (ii) If F is Lipschitzian continuous in Rn , then E(x) is Lipschitzian continuous in Rn [18]. It follows from Theorem 2.5 that (x0 ) = +. Now we assume that the level set L(x0 ) is bounded. If (x0 ) ¡ + , then it follows from Theorem 2.6 that

t (x0 ) lim x(t) = +: = inf {s¿0 | s ¡ (x0 ); x(s) Lc (x0 )} ¡ +

Let 0 where Lc (x0 ) is the complement of the set L(x0 ) in Rn . Since L(x0 ) is closed due to the continuity of E(x) and bounded by the assumption, we have x( 0 ) L(x0 ) and 0 ¡ (x0 ). Hence for some s ( 0 ; (x0 )), E(x(s)) ¿ E(x( 0 )): However, it follows from (2) and (5) that dE(x(t)) = E(x(t))T (x(t)) = - E(x(t)) 2 60; dt i.e., E(x(·)) is nonincreasing on [t0 ; (x0 )), a contradiction to (6). This completes (ii). (7) (6)

Although Theorem 4.2 states the existence of the solution trajectory of (5), it does not state its convergence. The following corollary provides such a result. Corollary 4.3. (i) Let x(t): [t0 ; (x0 )) be the unique maximal solution to (5). If (x0 ) = + and {x(t)} is bounded; then

t+

lim E(x(t)) = 0:

(ii) If F(x) is a uniform P function; then L(x0 ) is bounded and every accumulation point of the trajectory x(t) is a solution to NCP(F). Proof. (i) It is proved in (7) that E(x(t)) is a monotonically decreasing function in t. We also note that E(x(t)) is a nonnegative function over Rn , i.e., E(x(t)) is bounded from below. The dynamic system (5) corresponds the steepest descent dynamic model for the unconstrained minimization problem (2). Hence the analysis carried on this model in [21] implies that the trajectory of (5) would reach a steady state. (ii) It is known that if F is a uniform P function in Rn , then the level set L(x0 ) is bounded [10,13]. It follows from Theorem 4.2 that (x0 ) = + and the trajectory {x(t)} L(x0 ) by the proof of Theorem 4.2(ii). Let x be any accumulation point of the trajectory, we have from (i) of this theorem, E(x ) = 0; i.e., x is an equilibrium point of (5). Since F is a P function, it is also a P0 function, Proposition 4.1 implies that x is a solution of NCP(F). As an example shown in [24], the conclusion in Corollary 4.3(i) may not hold if {x(t)} is unbounded. The boundedness of the level set L(x0 ) guarantees the boundedness of {x(t)}. The class of uniform P functions are su cient for these boundedness. There are modiÃ¿cations of the FischerBurmeister function aiming at guaranteeing the boundedness of corresponding level sets for general class of functions. For example, a penalized FischerBurmeister NCP function is introduced in [2], which is deÃ¿ned by (a; b) = - (a; b) + (1 - )a+ b+ ; where (0; 1) is an arbitrary but Ã¿xed parameter, a+ = max{a; 0} and b+ = max{b; 0}. Correspondingly, we can deÃ¿ne (x), the objective function E (x), and the level set L (x0 ). It is known that if F is monotone and there is a strictly feasible point, i.e., x ¿ 0; F(x) ¿ 0, then the level set L (x0 ) is bounded [2, Proposition 3:10]. Hence by Corollary 4.3 the corresponding neural network model based on the penalized FischerBurmeister NCP function can solve monotone problems with a strictly feasible point. Although the penalized FischerBurmeister function enjoys some advantages over the FischerBurmeister function, we restrict our discussion on the FischerBurmeister function only, because the results for the penalized FischerBurmeister function can be easily obtained from the results for the FischerBurmeister function. 4.2. Stability of an isolated equilibrium In this subsection, we will discuss various stabilities introduced for di erential equations in Section 2 to our di erential equation (5). Now let x be a solution to NCP(F). Obviously x is an equilibrium point of (5) by Proposition 4.1. To discuss the stability at x we assume that x is an isolated equilibrium point of (5), that is there is a neighborhood Rn of x such that E(x ) = 0; E(x) = 0 x and x = x :

Now we state our Ã¿rst result in this subsection. Theorem 4.4. If x is an isolated equilibrium point of (5); x is asymptotically stable for (5). Proof. First we show that E(x) is a Lyapunov function over the set for equation (5). We note that E(x) is nonnegative function over Rn . Since x is a solution to NCP(F), obviously E(x ) = 0. For any x \{x }, we claim that E(x) ¿ 0. Otherwise if there is an x \{x } satisfying E(x) = 0 (this implies (x) = 0), then for any V @ (x), we have from Proposition 2.4 E(x) = V T (x) = 0: Hence x is an equilibrium point of (5), contradicting with the isolatedness of x in check the second condition in (4). dE(x(t)) = [x(t) E(x(t))]T x(t) = - x(t) E(x(t)) 2 60: dt

Now we (8)

Hence the function E(x) is a Lyapunov function for (5) over the set of x , we have from (8) that Because the isolatedness dE(x(t)) ¡ 0; x(t) and x(t) = x : dt It follows from Theorem 2.10(ii) that x is asymptotically stable for (5). We recall that x is said to be a regular solution to NCP(F) if every element V @ (x ) is nonsingular. It is also well known that the mapping of the generalized Jacobian (·) is upper semicontinuous [8]. Let B(x ; ) denote a ball centered at x with radius , i.e., B(x ; ) = {x Rn | x - x 6 }: From the upper semicontinuity of @ (·), we have the following property for a regular solution. It can be stated similarly as [29, Proposition 3:1]. Proposition 4.5. If x is a regular solution to NCP(F); then there is a neighborhood B(x ; ) of x and a constant C such that for all x B(x ; ) and any V (x); V is nonsingular and max{ V ; V -1 }6C: Theorem 4.6. If x is a regular solution to NCP(F); then it is exponentially stable for (5). Proof. It is known that if x is a regular solution, then it is an isolated solution for the equation (x) = 0 [29,26]. From Proposition 2.4(ii), x is an isolated equilibrium of (5). Hence x is asymptotically stable by Theorem 4.4. Let ¿ 0 be su ciently small such that for any x(t0 ) B(x ; ); x(t) x as t +, and the results in Proposition 4.5 hold. Hence there exists Ã„1 ¿ 0 and Ã„2 ¿ 0 such that Ã„1 v 2 6vT V T Vv6Ã„2 v 2 ; Since the function Reducing x B(x ; ): V @ (x): (9) (10) (11) (x) is semismooth, we have from Proposition 2.4 that (x) = (x ) + V (x - x ) + o( x - x ); |o( x - x )|6 x - x for some 0 ¡ ¡ Ã„1 and x B(x ; ). Now let (t) = x(t) - x 2 ; Then t [t0 ; ): if necessary, we assume that the last term in (10) satisÃ¿es d (t) = 2[x(t) - x ]T x(t) = -2[x(t) - x ]T E(x(t)) = -2[x(t) - x ]T V T (x(t)) dt for all V @ (x(t)). Let = inf {t [t0 ; ) | x(t) - x ¿ } (12) be the Ã¿rst exit time of the solution from the ball B(x ; ); (inf = +). Substituting (10) into (12), and noticing (9), and (x ) = 0, we have for all t I = [t0 ; ) d (t) 6 -2[x(t) - x ]T V T V (x(t) - x ) + x(t) - x dt 6 (-2Ã„1 + ) x(t) - x 2 = (-2Ã„1 + ) (t): By [34, p. 95, Corollary 2:1] (t)6e(-2Ã„1 + )t (t0 ); or equivalently x(t) - x 6e!t x(t0 ) - x ; where ! = -Ã„1 + =2 ¡ 0. If 2

t I; t I; (13) + , then

6lim x(t) - x 6e! x(t0 ) - x ¡ ; a contradiction. Thus exponential stability. = + and relation (13) completes the proof by noticing the deÃ¿nition of To better understand the importance of Theorem 4.6, we recall the di erential equation (3). We have the classical result [34, p. 101, Theorem 2:3]. Theorem 4.7. Assume that a continuous function f is di erentiable at an equilibrium point x. Then x is exponentially stable for (3) if and only if the Jacobi matrix f (x) is stable. A matrix A Rn×n is stable if and only if !(A) = sup{Re | (A)} ¡ 0; where (A) is the set of all eigenvalues of A, and Re is the real part of (A). This result cannot be directly used for our problem (5) because the gradient mapping of E(x) may not be di erentiable at a solution. Hence Theorem 4.6 can be viewed as an extension of Theorem 4.7 from the di erentiable case to a semismooth case. Such extension deserves a further study to see if the semismoothness requirement is essential. There are various conditions guaranteeing the regularity at a solution x . We do not intend to review these conditions here. Interested readers may refer to [10,11,28,13,33]. But we do want to point out that the regularity often leads to the superlinear convergence of some optimization methods. Here we relate the regularity to the exponential stability, which we believe is quite interesting. 5. Numerical examples In this section, the following four nonlinear complementarity problems were tested on neural network (5).

Example 5.1. This is the fourth example of Watson [32] deÃ¿ned by 5 i=1

F(x) = 2 exp (xi - i + 2)2 x1 + 1 x2 x3 - 1 : x4 - 2 x5 - 3

F(x) is monotone on the positive orthant but not P0 on Rn . The solution for this problem is (0; 0; 1; 2; 3). Example 5.2 (KojimaShindo nondegenerate NCP test problem [22]).

2 2 3x1 + 2x1 x2 + 2x2 + x3 + 3x4 - 6 2x2 + x + x2 + 3x + 2x - 2 1 3 4 2 F(x) = 2 1 : 2 3x1 + x1 x2 + 2x2 + 2x3 + 3x4 - 1 2 2 x1 + 3x2 + 2x3 + 3x4 - 3

This problem has the solution x = ( 6=2; 0; 0; 0:5); F(x ) = (0; 2 + 6=2; 5; 0): Example 5.3 (KojimaShindo degenerate NCP test problem [22]).

2 2 3x1 + 2x1 x2 + 2x2 + x3 + 3x4 - 6 2x2 + x + x2 + 10x + 2x - 2 1 3 4 2 F(x) = 2 1 : 2 3x1 + x1 x2 + 2x2 + 2x3 + 9x4 - 9 2 2 x1 + 3x2 + 2x3 + 3x4 - 3

This problem has the following two solutions: x = ( 6=2; 0; 0; 0:5); F(x ) = (0; 2 + 6=2; 0; 0) and x = (1; 0; 3; 0); F(x ) = (0; 31; 0; 4):

Example 5.4. This problem is NCP reformulation of a quadratic problem taken from [11, p. 25]. F(x) = 2x1 + 4x2 2x2 + 4x1 :

It is easy to see that x = (0; 0) is the only solution of this problem. However it is not a regular solution since @ (x ) contains the zero matrix. Our simulation is conducted on Matlab version 5.2. The ordinary di erential equation solver engaged is ode23. The initial state used in all implementations is x0 = (0; : : : ; 0)T for Examples 5:15:3, and x0 = (1; 1)T for Example 5.4. The scaling factor in our tests is chosen at 103 and

L.-Z. Liao et al. / Journal of Computational and Applied Mathematics 131 (2001) 343359 Table 1 Simulation results for Examples 5.15.4 Example 5.1 5.2 5.3 5.4 103 106 103 106 103 106 103 106 tf (seconds) 0.011 1.1e-5 0.022 2.2e-5 0.022 2.2e-5 0.014 1.4e-5 x(tf ) - x 9.6e-6 9.6e-6 1.8e-5 1.8e-5 1.7e-5 1.7e-5 1.4e-5 1.4e-5 E(x(tf )) 9.6e-6 9.6e-6 8.1e-6 8.1e-6 8.2e-6 8.2e-6 1.2e-6 1.2e-6 E(x(tf )) 4.6e-11 4.6e-11 7.3e-11 7.3e-11 7.0e-11 7.0e-11 8.4e-11 8.4e-11

Fig. 2. Evolution of x(t) for Example 5.1 with

= 103. 106 . The stopping criterion in all our tests is either E(x(t)) 610-5 or E(x(t))610-10 . Table 1 summarizes the simulation results. The tf in Table 1 represents the Ã¿nal time when the stopping criterion is met. For Example 5.3, the trajectory of the neural network model (5) converges to ( 6=2; 0; 0; 0:5). The evolutions of the solution trajectories for these examples are illustrated in Figs. 25. From Figs. 25 and Table 1, we can make the following observations: (i) All trajectories converge to their corresponding static states, respectively. This convergence is faster when a larger scaling factor is applied. (ii) The trajectories of some variables may go outside of the feasible region, i.e., {x | x¿0}. For example, the trajectory of x1 in Example 5.1 goes down as low as -0:74, while the trajectories of x3 in Examples 5.25.3 stay negative in most of the times. So it is important for our neural network (5) that the function F is deÃ¿ned everywhere. However, in some real-life

Fig. 3. Evolution of x(t) for Example 5.2 with = 103 .

Fig. 4. Evolution of x(t) for Example 5.3 with = 103 .

applications, the function F is often deÃ¿ned in a speciÃ¿ed region, not in the whole space. We leave the research along this direction to a possible future study. (iii) Although Example 5.3 is degenerate while Example 5.2 is not, the convergence rate for both problems are almost the same. This observation raises the following question: Does the degeneracy slow down the convergence of the neural network proposed here? Once again we leave this problem to the future research. (iv) Example 5.4 is tested for the case where the solution is not regular, and yet isolated. From Fig. 5, we can see that the convergence is really fast. These limited numerical experiments veriÃ¿ed the stability of the proposed neural network. 6. Conclusion remarks In this paper, we have studied a neural network approach for solving nonlinear complementarity problems. The neural network is obtained from a well-known unconstrained minimization reformulation of the complementarity problem via the FischerBurmeister NCP function. Besides the study on the existence and the convergence of the trajectory of the neural network, we pay much attention to the study of the stability aspects, among which are the stability in the sense of Lyapunov, the asymptotic stability as well as the exponential stability. Especially, it is shown that the regularity condition, which often implies the superlinear convergence of the optimization methods, implies the exponential stability of the neural network. Numerical experiments are provided to demonstrate the e ciency of our neural network model. We should point out that most of the results developed here are also valid for neural networks derived from a number of NCP functions except the result on exponential stability. This is because that for di erent NCP functions, we need di erent regularity conditions to guarantee the nonsingularity of the generalized Jacobian, see [31]. Now we discuss some directions for future research. In some real applications, the function F is only deÃ¿ned on some speciÃ¿ed region (feasible region), and is undeÃ¿ned outside. However, our neural network (5) is proposed on the whole space Rn . So it is ideal to modify the neural network to cover such class of problems. In other words, the modiÃ¿ed neural network should keep the trajectory within the speciÃ¿ed region so that the function F is well-deÃ¿ned during the implementation. Another possible research is to extend the neural network (5) for NCP(F) to a more general case, the box constrained variational inequality problem. This could be studied by using some kind of BVIP functions from [27]. It should be pointed out that the box constraints can be fulÃ¿lled by employing limiting integrators with nonlinear (hardware) limiters at their output [7]. In such an approach, all box constraints are "hard", i.e., the constraints must not be violated either at the Ã¿nal solution or during the optimization process. Unlike the above considerations, Chen and Fang [4] transform the box constrained region into an unbounded region by some nonlinear transformations. All these progresses of treating the box constraints should be taken into account while proposing neural networks for the box constrained variational inequality problem. Finally, we like to point out that there are ways such as the smoothing [3] to be used to reduce the computation complexity of the element in (x). This would result in a class of smoothing neural networks, which would have a di erent structure from the one proposed in this paper.

References

[1] A. Bouzerdoum, T.R. Pattison, Neural network for quadratic optimization with bound constraints, IEEE Trans. Neural Networks 4 (1993) 293304. [2] B. Chen, X. Chen, C. Kanzow, A penalized FischerBurmeister NCP function: theoretical investigation and numerical results, Applied Mathematics Report AMR 97=28, School of Mathematics, University of New South Wales (Revised 1999). [3] X. Chen, L. Qi, D. Sun, Global and superlinear convergence of the smoothing Newton method and its application to general box constrained variational inequalities, Math. Comput. 67 (1998) 519540. [4] Y.-H. Chen, S.-C. Fang, Solving convex programming problems with equality constraints by neural networks, Comput. Math. Appl. 36 (1998) 4168. [5] L.O. Chua, G.-N. Lin, Nonlinear programming without computation, IEEE Trans. Circuits Systems CAS-31 (1984) 182188. [6] A. Cichocki, R. Unbehauen, Neural Networks for Optimization and Signal Processing, Wiley, New York, 1993. [7] A. Cichocki, R. Unbehauen, K. Weinzierl, R. Holzel, A new neural network for solving linear programming problems, Eur. J. Oper. Res. 93 (1996) 244256. [8] F.H. Clarke, Optimization and Nonsmooth Analysis, Wiley, New York, 1983. [9] R.W. Cottle, J.-S. Pang, R.E. Stone, The Linear Complementarity Problem, Academic Press, Boston, 1992. [10] T. De Luca, F. Facchinei, C. Kanzow, A semismooth equation approach to the solution of nonlinear complementarity problems, Math. Programm. 75 (1996) 407439. [11] F. Facchinei, A. Fischer, C. Kanzow, On the accurate identiÃ¿cation of active constraints, SIAM J. Optim. 9 (1999) 1432. [12] F. Facchinei, C. Kanzow, A nonsmooth inexact Newton method for the solution of large-scale nonlinear complementarity problems, Math. Programm. 76 (1997) 493512. [13] F. Facchinei, J. Soares, A new merit function for nonlinear complementarity problems and a related algorithm, SIAM J. Optim. 7 (1997) 225247. [14] M.C. Ferris, J.-S. Pang, Engineering and economic applications of complementarity problems, SIAM Rev. 39 (1997) 669713. [15] A. Fischer, A special Newton-type optimization method, Optimization 24 (1992) 269284. [16] A. Fischer, An NCP-function and its use for the solution of complementarity problems, in: D. Du, L. Qi, R. Womersley (Eds.), Recent Advance in Nonsmooth Optimization, World ScientiÃ¿c Publishers, New Jersey, 1995, pp. 88105. [17] M. Fukushima, L. Qi (Eds.), Reformulation -- Nonsmooth, Piecewise Smooth, Semismooth and Smoothing Methods, Kluwer Academic Publisher, Nowell, Maryland, 1999.

[18] C. Geiger, C. Kanzow, On the resolution of monotone complementarity problems, Comput. Optim. Appl. 5 (1996) 155173. [19] J.J. HopÃ¿eld, D.W. Tank, Neural computation of decisions in optimization problems, Biol. Cybernet. 52 (1985) 141152. [20] H. Jiang, L. Qi, A new nonsmooth equation approach to nonlinear complementarity problems, SIAM J. Control Optim. 35 (1997) 178193. [21] M.P. Kennedy, L.O. Chua, Neural network for nonlinear programming, IEEE Trans. Circuits Systems 35 (1988) 554562. [22] M. Kojima, S. Shindo, Extensions of Newton and quasi-Newton methods to systems of PC 1 equations, J. Oper. Res. Soc. Jpn. 29 (1986) 352374. [23] L.-Z. Liao, H.-D. Qi, A neural network for the linear complementarity problem, Math. Comput. Modelling 29 (1999) 918. [24] W.E. Lillo, M.-H. Loh, S. Hui, S.H. Zak, On solving constrained optimization problems with neural networks: a penalty method approach, IEEE Trans. Neural Networks 4 (1993) 931940. [25] J.-S. Pang, L. Qi, Nonsmooth equations: motivation and algorithms, SIAM J. Optim. 3 (1993) 443465. [26] L. Qi, Convergence analysis of some algorithms for solving nonsmooth equations, Math. Oper. Res. 18 (1993) 227244. [27] L. Qi, Regular pseudo-smooth NCP and BVIP functions and globally and quadratically convergent generalized Newton methods for complementarity and variational inequality problems, Math. Oper. Res. 24 (1999) 440471. [28] L. Qi, H. Jiang, Semismooth Karush-Kuhn-Tucker equations and convergence analysis of Newton and quasi-Newton methods for solving these equations, Math. Oper. Res. 22 (1997) 301325. [29] L. Qi, J. Sun, A nonsmooth version of Newton's method, Math. Programm. 58 (1993) 353368. [30] A. RodrÃƒ iguez-VÃƒ zquez, R. DomÃƒ a inguez-Castro, A. Rueda, J.L. Huertas, E. SÃƒ nchez-Sinencio, Nonlinear a switched-capacitor `neural' networks for optimization problems, IEEE Trans. Circuits Systems 37 (1990) 384398. [31] D. Sun, L. Qi, On NCP functions, Comput. Optim. Appl. 13 (1999) 201220. [32] L.T. Watson, Solving the nonlinear complementarity problem by a homotopy method, SIAM J. Control Optim. 17 (1979) 3646. [33] N. Yamashita, M. Fukushima, ModiÃ¿ed Newton methods for solving semismooth reformulations of monotone complementarity problems, Math. Programm. 76 (1997) 469491. [34] J. Zabczyk, Mathematical Control Theory: An Introduction, Birkhauser, Boston, 1992. [35] S.H. Zak, V. Upatising, S. Hui, Solving linear programming problems with neural networks: a comparative study, IEEE Trans. Neural Networks 6 (1995) 94104.