A Two-Step Newton Algorithm for the Weighted Complementarity Problem with Local Biquadratic Convergence

Abstract

We discuss the weighted complementarity problem, extending the nonlinear complementarity problem on $R^n$. In contrast to the NCP, many equilibrium problems in science, engineering, and economics can be transformed into WCPs for more efficient methods. Smoothing Newton algorithms, known for their at least locally superlinear convergence properties, have been widely applied to solve WCPs. We suggest a two-step Newton approach with a local biquadratic order convergence rate to solve the WCP. The new method needs to calculate two Newton equations at each iteration. We also insert a new term, which is of crucial importance for the local biquadratic convergence properties when solving the Newton equation. We demonstrate that the solution to the WCP is the accumulation point of the iterative sequence produced by the approach. We further demonstrate that the algorithm possesses local biquadratic convergence properties. Numerical results indicate the method to be practical and efficient.

1. Introduction

The weighted complementarity problem (WCP) was originally introduced by Potra [1], extending the standard nonlinear complementarity problem on $R^n$ (NCP). Because of its substantial applicability in various fields such as engineering, management, science, and market equilibrium, it has garnered a great deal of interest from researchers. In particular, problems such as linear programming, weighted centering problems [2], and Fisher market equilibrium problems [3] can be expressed using the model of the weighted linear complementarity problem (WLCP) as shown below:

$$x\in R_{+}^n,y\in R_{+}^n,xy=w,Ax+By+Cu=d,$$

(1)

and the WLCP model provides a more efficient approach to these problems than the NCP model does.

Here we consider a more general WCP model: to find a triple $(x,y,u)\in R^{2n}\times R^m$ satisfying

$$x\in R_{+}^n,y\in R_{+}^n,xy=w,G(x,y,u)=0,$$

(2)

in which $G(x,y,u):R^{2n+m}\to R^{n+m}$ is a nonlinear mapping and $w\in R_{+}^n$ is a weighted vector. When $w=0$, the WCP (2) becomes the NCP [4,5,6]. When the function $G(x,y,u)$ becomes linear, then the WCP (2) is reduced to (1) [7,8,9,10,11].

The WCP has been researched extensively, and numerous efficient algorithms have been proposed. Smoothing Newton algorithms have gained popularity among researchers because of their at least locally superlinear convergence properties [12,13,14,15,16]. With the aid of the complementarity function, the WCP problem can be transformed into an equivalent system of equations, which are then solved using Newton algorithms. This represents the core concept of employing Newton algorithms to solve the WCP. For WLCPs, Zhang [17] presents a smoothing Newton algorithm. Tang [18] offers a damped derivative-free Gauss–Newton method for a nonmonotone WLCP that is globally convergent and requires no problem assumptions.

When solving the nonlinear equations $G\left(z\right)=0$, it is widely known that the two-step Newton approach [19,20,21,22] typically achieves higher-order convergence, such as third-order or fourth-order convergence, than the classical Newton method. The two-step Newton approach for solving nonlinear equations has recently been used in an attempt to solve WCPs. Tang et al. [23] present a smoothing Newton approach to solve WCPs with local cubic convergence rates. The approach accelerates the convergence rate by adding an approximate Newton step, as the sequence of iterations is close to the solution, raising the local convergence rate from second-order to third-order. This approach can be viewed as an approximate two-step Newton algorithm. Liu et al. [24] propose a super-quadratic smoothing Newton algorithm for solving WCPs based on the following two-step Newton algorithm for solving nonlinear equations $G\left(z\right)=0$:

$$z^{k+1}=s^k-G{(z^k)}^{-1}G(s^k),\mathrm{where}{s}^k=z^k-G{(z^k)}^{-1}G(z^k).$$

When solving Newton equations with the right parameter selections, the algorithm in particular possesses the property of local cubic convergence. Argyros [25] and Magrenan Ruiz [26] discuss the two-step Newton algorithm, whose iteration sequence $\left\{z^k\right\}$ is generated by

$$z^{k+1}=s^k-G{(s^k)}^{-1}G(s^k),\mathrm{where}{s}^k=z^k-G{(z^k)}^{-1}G(z^k),$$

and they show that the two-step Newton algorithm has a fourth-order convergence. A natural question emerges whether we can apply this algorithm to solve the WCP to obtain an algorithm that is more efficient than the Newton method with local cubic convergence rate. With these considerations, we develop a two-step Newton iterative technique for solving the WCP (2) with a local biquadratic convergence property. The major contributions of the new algorithm are as follows.

• The algorithm computes two Newton equations directly to obtain the next iteration point, in contrast to the algorithm in [23]. If the value of the objective function meets a certain descent criterion, the algorithm takes the iteration point produced by the two Newton directions directly as the next iteration point; otherwise, the step size is confirmed via a derivative-free line search to find the next iteration point. By doing this, the computational efficiency of the algorithm is successfully improved without adding to the time investment.
• Compared with the algorithm in [24], we employ different Jacobian matrices for calculating Newton equations and add the new term ${\chi }_k=min\{1,{\xi }_k^4\}$ when doing so, in order to guarantee the local biquadratic convergence property. Due to this architecture, the new algorithm exhibits local biquadratic convergence under the right conditions.
• Because the nonlinear complementarity problem [4,5,6] and system of inequalities [27,28] can be transformed into an equivalent system of equations, the novel algorithm provides a fresh approach to solve these problems.

The paper proceeds as follows. Section 2 shows a smoothing function whose fundamental characteristics are also discussed. Section 3 proposes a new two-step smoothing Newton approach for WCPs and demonstrates its viability. Section 4 and Section 5 deal with the global and local strong convergence properties, respectively. Section 6 presents numerical experiments. Section 7 contains final remarks.

2. Preliminaries

In this paper, we deal with the WCP (2) by using smoothing Newton methods. To this end, we first introduce a smoothing function as

$${\phi }_c(\xi ,p,q)=\sqrt{p^2+q^2+2c+4{\xi }^2}-(p+q),$$

(3)

where $\xi \in (0,1)$ and $c\in R_{+}$. The following lemma shows some basic properties of ${\phi }_c(\xi ,p,q)$, whose proof is obvious by simple calculations.

Lemma 1. 

Let ${\phi }_c(\xi ,p,q)$ be defined by (3). Then

1.

$p\ge 0,q\ge 0,pq=c\iff {\phi }_c(0,p,q)=0.$

2.

If $\xi >0$, then ${\phi }_c(\xi ,p,q)$ is continuously differentiable for any $(p,q)\in R^2$.

With the smoothing function ${\phi }_c(\xi ,p,q)$ defined by (3), for the WCP (2), we define a function $F(\xi ,x,y,u):R\times R^{2n+m}\to R\times R^{2n+m}$ by

$$F(\xi ,x,y,u)=\left(\begin{array}{c}\xi \\ {\phi }_w(\xi ,x,y)\\ G(x,y,u)\end{array}\right),$$

(4)

where

$${\phi }_w(\xi ,x,y)=\left(\begin{array}{c}{\phi }_{w_1}(\xi ,x_1,y_1)\\ ⋮\\ {\phi }_{w_n}(\xi ,x_n,y_n)\end{array}\right).$$

(5)

Let $z=(\xi ,x,y,u)$ to simplify the notation, then,

$$F\left(z\right)=F(\xi ,x,y,u)=0\iff \xi =0,G(x,y,u)=0,{\phi }_w(\xi ,x,y)=0.$$

It is simple to demonstrate that $F\left(z\right)$ is continuously differentiable on $R^{2n+m}$ for any $\xi >0$ by using Lemma 1. By simple calculations, we have the Jacobian matrix for $F\left(z\right)$ as below:

$$F^{\prime }\left(z\right)=\left(\begin{array}{cccc}1& 0& 0& 0\\ J_1& J_2& J_3& 0\\ 0& G_x^{\prime }& G_y^{\prime }& G_u^{\prime }\end{array}\right),$$

(6)

where

$$\begin{array}{cc}\hfill J_1& =\mathrm{vec}\left\{{\displaystyle \frac{4\xi }{\sqrt{x_i^2+y_i^2+2w_i+4{\xi }^2}}}\right\},i=1,2,\dots ,n,\hfill \end{array}$$

(7)

$$\begin{array}{cc}\hfill J_2& =\mathrm{diag}\left\{{\displaystyle \frac{x_i}{\sqrt{x_i^2+y_i^2+2w_i+4{\xi }^2}}}-1\right\},i=1,2,\dots ,n,\hfill \end{array}$$

(8)

$$\begin{array}{cc}\hfill J_3& =\mathrm{diag}\left\{{\displaystyle \frac{s_i}{\sqrt{x_i^2+y_i^2+2w_i+4{\xi }^2}}}-1\right\},i=1,2,\dots ,n.\hfill \end{array}$$

(9)

We next discuss the nonsingularity of Jacobian matrix $H^{\prime }\left(z\right)$. For this purpose, we need an assumption.

Assumption 1. 

Suppose that $Rank\left(G_y^{\prime }\right)=m$, for any $(\Delta x,\Delta y,\Delta u)\in R^{2n+m}$, if 

$$G_x^{\prime }\Delta x+G_y^{\prime }\Delta y+G_u^{\prime }\Delta u=0,$$

(10)

then $\langle \Delta x,\Delta y\rangle \ge 0$.

If $G(x,y,u)$ is linear, i.e., $G(x,y,u)=Ax+By+Cu-d$, then (10) reduces to

$$A\Delta x+B\Delta y+C\Delta u=0,$$

which means that the associated WLCP is monotone [9,17,29]. Similar to Lemma 1 in [17], we can draw a conclusion as follows.

Theorem 1. 

If Assumption 1 is true, then $F^{\prime }\left(z\right)$ is invertible for any $\xi >0$.

3. Description of the Method

This section presents and illustrates the feasibility of a new two-step smoothing Newton approach. We start with the formal explanation of the new approach.

Remark 1. 

1.

It is worth noting that the Newton direction $d_1^k+d_2^k$ obtained by computing Newton equations twice is not necessarily the descent direction of the objective function. Therefore, to guarantee the global convergence properties, we introduce a derivation-free line search. When the objective function satisfies a certain descent quantity, we can use $d_1^k+d_2^k$ directly as a descent direction. Otherwise, we utilize (15) to generate a step length for obtaining the next iteration point.

2.

In Step 2, the additional term ${\chi }_k=min\{1,{\xi }_k^4\}$ is added to Newton Equations (11) and (12), in contrast to the current smoothing Newton methods [17,23,30]. The property of local biquadratic convergence of Algorithm 1 depends on this particular perturbation term. Algorithm 1 has a similar computational cost to the traditional Newton approach even though it computes the Newton direction twice.

3.

The main distinction between Algorithm 1 and the accelerated algorithm in [23] as two-step Newton algorithms is that Algorithm 1 employs two Newton directions from the beginning, whereas the accelerated algorithm in [23] begins with one Newton direction and adds a second Newton direction when certain conditions are met. We also note that Algorithm 1 is able to solve the WCP better than the accelerated method in the following section of numerical experiments.

4.

In addition to the difference between Algorithm 1 and the algorithm in [24] in solving Newton directions, another difference is that when the Newton equation is used as the descent direction for the line search, there are different choices for the descent direction. As described in Step 4, the descent direction is chosen to be the sum of two Newton directions under certain conditions. In the subsequent discussion, it will be shown that this choice is made to ensure global convergence of the algorithm.

Algorithm 1 A Two-Step Newton Method

Input parameters: Required stopping criterion $\delta >0$, ${\xi }_0>0$, $c,l,\rho \in (0,1)$, $\gamma ,{\sigma }_1,{\sigma }_2\in (0,1)$, $\kappa \ge 0$ such that $\gamma +\gamma \kappa +\kappa c<1$, and $h={(\gamma ,0,0,0)}^T\in R\times R^{2n+m}$. $\left\{{\eta }_k\right\}\subseteq R_{+}$ satisfies that ${\sum }_{k=0}^{\infty }{\eta }_k\le \eta <\infty$ and $\underset{k\to \infty }{lim}{\eta }_k=0$, and starting point $(x^0,y^0,u^0)\in R^{2n+m}$.

Output: an approximate solution $(x^k,s^k,u^k)$ to the WCP (2);

Step 0. Let $z^0=({\xi }_0,x^0,y^0,u^0)$ and $k=0$.

Step 1. If $\parallel F\left(z^k\right)\parallel \le \delta$, stop.

Step 2. a. Calculate $d_1^k$ by $$F^{\prime }\left(z^k\right)d_1^k=-F\left(z^k\right)+{\chi }_{k}h$$ (11) where ${\chi }_k=min\{1,{\xi }_k^4\}$. Let $s^k=z^k+d_1^k$. b. Calculate $d_2^k$ by $$F^{\prime }\left(s^k\right)d_2^k=-F\left(s^k\right)+{\chi }_{k}h.$$ (12)

Step 3. If $$\parallel F(z^k+d_1^k+d_2^k)\parallel >l·\parallel F\left(z^k\right)\parallel ,$$ (13) go to Step 4. Else, set $d^k=d_1^k+d_2^k$ and ${\beta }_k=1$, go to Step 5.

Step 4. Let $$d^k=\left\{\begin{array}{cc}{d}_1^k+d_2^k,\hfill & \mathrm{if}\parallel F\left(s^k\right)\parallel \le c\parallel F\left(z^k\right)\parallel \mathrm{and}\parallel F^{\prime }\left(z^k\right)F^{\prime }{\left(s^k\right)}^{-1}\parallel \le \kappa ,\hfill \\ d_1^k+{\beta }_k{d}_2^k,\hfill & \mathrm{otherwise},\hfill \end{array}\right.$$ (14) and let ${\beta }_k$ be the maximum of $\{{\rho }^0,\rho ,{\rho }^2,\cdots \}$, satisfying: $$\parallel F(z^k+{\rho }^{m\left(k\right)}{d}^k){\parallel }^2\le (1+{\eta }_k)\parallel F\left(z^k\right){\parallel }^2-{\sigma }_1{\left({\rho }^{m\left(k\right)}\right)}^2\parallel d^k{\parallel }^2-{\sigma }_2{\left({\rho }^{m\left(k\right)}\right)}^2{\parallel F\left(z^k\right)\parallel }^2.$$ (15)

Step 5. Set $z^{k+1}=z^k+{\beta }_k{d}^k$ and $k=k+1$. Return to Step 1.

To investigate the convergence of Algorithm 1, we first demonstrate that it is clearly defined.

Theorem 2. 

Supposing Assumption 1 is true, then Algorithm 1 is well-defined.

Proof of Theorem 2. 

As $F^{\prime }\left(z\right)$ is a nonsingular duo to Theorem 1, Step 2 is feasible. By the definition of ${\chi }_k$, we have that

$${\chi }_k=min\{1,{\xi }_k^4\}=1\le {\xi }_k\le \parallel F\left(z^k\right)\parallel ,\mathrm{if}{\xi }_k\ge 1$$

or

$${\chi }_k=min\{1,{\xi }_k^4\}={\xi }_k^4\le {\xi }_k\le \parallel F\left(z^k\right)\parallel ,\mathrm{if}{\xi }_k<1.$$

Thus,

$${\chi }_k\le \parallel F\left(z^k\right)\parallel .$$

(16)

We next discuss the following two cases.

I.

$d^k=d_1^k+d_2^k$

We obtain from (11), (12), (14), and (16) that

$$\begin{array}{cc}\hfill & F{\left(z^k\right)}^T{F}^{\prime }\left(z^k\right)d^k\hfill \\ \hfill =& F{\left(z^k\right)}^T{F}^{\prime }\left(z^k\right)(d_1^k+d_2^k)\hfill \\ \hfill =& F{\left(z^k\right)}^T\{-F\left(z^k\right)+{\chi }_{k}h+F^{\prime }\left(z^k\right)\left[F^{\prime }{\left(s^k\right)}^{-1}(-F\left(s^k\right)+{\chi }_{k}h)\right]\}\hfill \\ \hfill \le & -(1-\gamma )\parallel F\left(z^k\right){\parallel }^2-F{\left(z^k\right)}^T{F}^{\prime }\left(z^k\right)F^{\prime }{\left(s^k\right)}^{-1}F\left(s^k\right)\hfill \\ \hfill & +{\chi }_{k}F{\left(z^k\right)}^T{F}^{\prime }\left(z^k\right)F^{\prime }{\left(s^k\right)}^{-1}h\hfill \\ \hfill \le & -(1-\gamma -\kappa c-\kappa \gamma )\parallel F\left(z^k\right){\parallel }^2,\hfill \end{array}$$

(17)

which together with $\gamma +\gamma \kappa +\kappa c<1$ yields that

$$F{\left(z^k\right)}^T{F}^{\prime }\left(z^k\right)d^k<0.$$

(18)

Note that, for any $k\ge 0$, if the line search (15) is not satisfied, then

$$\begin{array}{cc}\hfill & \parallel F(z^k+{\beta }_k{d}^k){\parallel }^2\hfill \\ \hfill >& (1+{\eta }_k)\parallel F\left(z^k\right){\parallel }^2-{\sigma }_1{\beta }_k^2\parallel d^k{\parallel }^2-{\sigma }_2{\beta }_k^2{\parallel F\left(z^k\right)\parallel }^2\hfill \\ \hfill \ge & \parallel F\left(z^k\right){\parallel }^2-{\sigma }_1{\beta }_k^2\parallel d^k{\parallel }^2-{\sigma }_2{\beta }_k^2{\parallel F\left(z^k\right)\parallel }^2.\hfill \end{array}$$

(19)

Then,

$${\displaystyle \frac{\parallel F(z^k+{\beta }_k{d}^k){\parallel }^2-{\parallel F\left(z^k\right)\parallel }^2}{{\beta }_k}}>-{\sigma }_1{\beta }_k\parallel d^k{\parallel }^2-{\sigma }_2{\beta }_k{\parallel F\left(z^k\right)\parallel }^2.$$

(20)

By letting $k\to \infty$ on (20), we obtain

$$F{\left(z^k\right)}^T{F}^{\prime }\left(z^k\right)d^k\ge 0,$$

which is in contradiction with (18). We can thus derive a step size ${\beta }_k$ that satisfies (15).

II.

$d^k=d_1^k+{\beta }_k{d}_2^k$ We obtain from (11) and (16) that

$$\begin{array}{cc}\hfill F{\left(z^k\right)}^T{F}^{\prime }\left(z^k\right)d_1^k=& F{\left(z^k\right)}^T(-F\left(z^k\right)+{\chi }_{k}h)\hfill \\ \hfill \le & -(1-\gamma )\parallel F\left(z^k\right){\parallel }^2\hfill \\ \hfill <& 0.\hfill \end{array}$$

(21)

On the other hand, if (15) is not satisfied, then

$$\begin{array}{cc}\hfill & \parallel F(z^k+{\beta }_k(d_1^k+{\beta }_k{d}_2^k)){\parallel }^2\hfill \\ \hfill >& (1+{\eta }_k)\parallel F\left(z^k\right){\parallel }^2-{\sigma }_1{\beta }_k^2\parallel d_1^k+{\beta }_k{d}_2^k{\parallel }^2-{\sigma }_2{\beta }_k^2{\parallel F\left(z^k\right)\parallel }^2\hfill \\ \hfill \ge & \parallel F\left(z^k\right){\parallel }^2-{\sigma }_1{\beta }_k^2\parallel d_1^k+{\beta }_k{d}_2^k{\parallel }^2-{\sigma }_2{\beta }_k^2{\parallel F\left(z^k\right)\parallel }^2.\hfill \end{array}$$

(22)

Then,

$${\displaystyle \frac{\parallel F(z^k+{\beta }_k(d_1^k+{\beta }_k{d}_2^k)){\parallel }^2-{\parallel F\left(z^k\right)\parallel }^2}{{\beta }_k}}>-{\sigma }_1{\beta }_k\parallel d_1^k+{\beta }_k{d}_2^k{\parallel }^2-{\sigma }_2{\beta }_k{\parallel F\left(z^k\right)\parallel }^2.$$

(23)

By letting $k\to \infty$ on (23), we obtain

$$F{\left(z^k\right)}^T{F}^{\prime }\left(z^k\right)d_1^k\ge 0,$$

which is in contradiction with (21). We can thus derive a step size ${\beta }_k$ that satisfies (15).

In conclusion, Algorithm 1 is well-defined. □

4. Global Convergence

We start defining the set $\mathrm{\Omega }\left(z\right)$ as

$$\mathrm{\Omega }\left(z\right)=\{z\in R_{+}\times R^{2n+m}|\parallel F\left(z\right)\parallel \le e^{\frac{\xi }{2}}\parallel F\left(z^0\right)\parallel \}.$$

Theorem 3. 

If Assumption 1 is true, then Algorithm 1 generates a sequence $z^k=({\xi }_k,x^k,s^k,u^k)$ satisfying $0\le {\xi }_{k+1}\le {\xi }_k\le \cdots \le {\xi }_0$ and $z^k\in \mathrm{\Omega }\left(z\right)$.

Proof of Theorem 3. 

We show ${\xi }_k\ge 0$ by induction. Supposing that ${\xi }_k\ge 0$, it follows from (11) and (12) that

$$\Delta {\xi }_k^1=-{\xi }_k+{\chi }_k\gamma$$

(24)

and

$$\Delta {\xi }_k^2=-({\xi }_k+\Delta {\xi }_k^1)+{\chi }_k\gamma .$$

(25)

Then, by Step 5, we have

$$\begin{array}{cc}\hfill {\xi }_{k+1}& ={\xi }_k+{\beta }_k(\Delta {\xi }_k^1+{\beta }_k\Delta {\xi }_k^2)\hfill \\ \hfill & ={\xi }_k+{\beta }_k\{-{\xi }_k+{\chi }_k\gamma +{\beta }_k[-({\xi }_k+\Delta {\xi }_k^1)+{\chi }_k\gamma ]\}\hfill \\ \hfill & =(1-{\beta }_k){\xi }_k+{\beta }_k{\chi }_k\gamma ,\hfill \end{array}$$

(26)

indicating that ${\xi }_{k+1}\ge 0$ due to the fact that $0\le {\beta }_k\le 1$. Moreover, combining with (16) yields

$${\xi }_{k+1}-{\xi }_k=(1-{\beta }_k){\xi }_k+{\beta }_k{\chi }_k\gamma -{\xi }_k={\beta }_k({\chi }_k\gamma -{\xi }_k)\le 0,$$

i.e., ${\xi }_{k+1}\le {\xi }_k$ for any $k\ge 0$.

It follows from (15) that

$$\begin{array}{cc}\hfill \parallel F\left(z^{k+1}\right){\parallel }^2& \le (1+{\eta }_k){\parallel F\left(z^k\right)\parallel }^2\hfill \\ \hfill & \le (1+{\eta }_k)·(1+{\eta }_{k-1})\cdots (1+{\eta }_0){\parallel F\left(z^0\right)\parallel }^2\hfill \\ \hfill & \le {\left[{\displaystyle \sum _{j=0}^k}{\displaystyle \frac{1}{k+1}}(1+{\eta }_j)\right]}^{k+1}{\parallel F\left(z^0\right)\parallel }^2\hfill \\ \hfill & \le {\left(1+{\displaystyle \frac{\eta }{k+1}}\right)}^{k+1}{\parallel F\left(z^0\right)\parallel }^2\hfill \\ \hfill & \le e^{\eta }{\parallel F\left(z^0\right)\parallel }^2,\hfill \end{array}$$

(27)

which indicates that $z^k\in \mathrm{\Omega }\left(z\right)$. □

Theorem 4. 

If Assumption 1 is satisfied and $\left\{z^k\right\}$ is bounded, then $\underset{k\to \infty }{lim}{\xi }_k=0$.

Proof of Theorem 4. 

We have that $\left\{{\xi }_k\right\}$ is monotonically nonincreasing and bounded from Theorem 3, and is therefore convergent. Let $\underset{k\to \infty }{lim}{\xi }_k={\xi }_{*}\ge 0$. If ${\xi }_{*}=0$, the conclusion is clearly valid. Supposing that ${\xi }_{*}>0$, we next show a contradiction.

Letting $\underset{k\to \infty }{lim}{z}^k=z^{*}=({\xi }_{*},x^{*},y^{*},u^{*})$, then $\underset{k\to \infty }{lim}\parallel F\left(z^k\right)\parallel =\parallel F\left(z^{*}\right)\parallel \ge {\xi }_{*}>0.$ From (15), we have

$$\parallel F(z^k+{\beta }_k{d}^k){\parallel }^2\le (1+{\eta }_k)\parallel F\left(z^k\right){\parallel }^2-{\sigma }_1{\beta }_k^2\parallel d^k{\parallel }^2-{\sigma }_2{\beta }_k^2{\parallel F\left(z^k\right)\parallel }^2.$$

Since $\underset{k\to \infty }{lim}{\eta }_k=0$, by letting $k\to \infty$, we obtain that

$$\parallel F\left(z^{*}\right){\parallel }^2\le \parallel F\left(z^{*}\right){\parallel }^2-{\sigma }_1{\beta }_{*}^2\parallel d^{*}{\parallel }^2-{\sigma }_2{\beta }_{*}^2{\parallel F\left(z^{*}\right)\parallel }^2,$$

i.e.,

$${\beta }_{*}^2({\sigma }_1\parallel d^{*}{\parallel }^2+{\sigma }_2\parallel F\left(z^{*}\right){\parallel }^2)\le 0,$$

which indicates that ${\beta }_{*}=0$ due to $F\left(z^{*}\right)>0$ and ${\sigma }_1,{\sigma }_2>0$. We proceed to the discussion of two cases.

Case 1: $\parallel F\left(s^k\right)\parallel \le c\parallel F\left(z^k\right)\parallel$ and $\parallel F^{\prime }\left(z^k\right)F^{\prime }{\left(s^k\right)}^{-1}\parallel \le \kappa$.

Letting $\widehat{\beta }={\displaystyle \frac{{\beta }_k}{\rho }}$, it holds that

$$\begin{array}{cc}\hfill & \parallel F(z^k+\widehat{\beta }{d}^k{\parallel }^2\hfill \\ \hfill >& (1+{\eta }_k)\parallel F\left(z^k\right){\parallel }^2-{\sigma }_1{\widehat{\beta }}^2\parallel d^k{\parallel }^2-{\sigma }_2{\widehat{\beta }}^2{\parallel F\left(z^k\right)\parallel }^2\hfill \\ \hfill \ge & \parallel F\left(z^k\right){\parallel }^2-{\sigma }_1{\widehat{\beta }}^2\parallel d^k{\parallel }^2-{\sigma }_2{\widehat{\beta }}^2{\parallel F\left(z^k\right)\parallel }^2\hfill \end{array}$$

(28)

for sufficiently large k. Since

$$\begin{array}{cc}\hfill & \parallel F(z^k+\widehat{\beta }{d}^k){\parallel }^2\hfill \\ \hfill =& \parallel F\left(z^k\right)+\widehat{\beta }{F}^{\prime }\left(z^k\right)d^k{\parallel }^2+o\left(\widehat{\beta }\right)\hfill \\ \hfill =& \parallel F\left(z^k\right){\parallel }^2+2\widehat{\beta }F{\left(z^k\right)}^T{F}^{\prime }\left(z^k\right)d^k+o\left(\widehat{\beta }\right),\hfill \end{array}$$

(29)

combining (28) with (29), we obtain that

$$2F{\left(z^k\right)}^T{F}^{\prime }\left(z^k\right)d^k+o\left(\widehat{\beta }\right)>-\widehat{\beta }({\sigma }_1\parallel d^k{\parallel }^2+{\sigma }_2\parallel F\left(z^k\right){\parallel }^2).$$

(30)

By (11), (12), and (14), we have

$$\begin{array}{cc}\hfill & F{\left(z^k\right)}^T{F}^{\prime }\left(z^k\right)d^k\hfill \\ \hfill =& F{\left(z^k\right)}^T{F}^{\prime }\left(z^k\right)(d_1^k+d_2^k)\hfill \\ \hfill =& F{\left(z^k\right)}^T\{-F\left(z^k\right)+{\chi }_{k}h+F^{\prime }\left(z^k\right)\left[F^{\prime }{\left(s^k\right)}^{-1}(-F\left(s^k\right)+{\chi }_{k}h)\right]\}\hfill \\ \hfill \le & -\parallel F\left(z^k\right){\parallel }^2+{\chi }_k\gamma \parallel F\left(z^k\right)\parallel +\kappa c\parallel F\left(z^k\right){\parallel }^2+{\chi }_k\kappa \gamma \parallel F\left(z^k\right)\parallel .\hfill \end{array}$$

(31)

Then, from (30) and (31), we obtain

$$\begin{array}{cc}\hfill & 2[(-1+\kappa c)\parallel F\left(z^k\right){\parallel }^2+{\chi }_k(\gamma +\kappa \gamma )\parallel F\left(z^k\right)\parallel ]\hfill \\ \hfill >& -\widehat{\beta }({\sigma }_1\parallel d^k{\parallel }^2+{\sigma }_2\parallel F\left(z^k\right){\parallel }^2).\hfill \end{array}$$

(32)

By letting $k\to \infty$ in (32), it holds that

$$2[(-1+\kappa c)\parallel F\left(z^{*}\right)\parallel +(\gamma +\kappa \gamma ){\chi }_{*}]\ge 0,$$

where ${\chi }_{*}=min\{1,{\xi }_{*}^4\}$. Then, we have

$${\chi }_{*}\ge {\displaystyle \frac{1-\kappa c}{\gamma +\gamma \kappa }}\parallel F\left(z^{*}\right)\parallel >{\displaystyle \frac{\gamma +\gamma \kappa }{\gamma +\gamma \kappa }}\parallel F\left(z^{*}\right)\parallel =\parallel F\left(z^{*}\right)\parallel ,$$

which contradicts (16). Thus, ${\xi }_{*}=0$.

Case 2: If the condition that $\parallel F\left(s^k\right)\parallel \le c\parallel F\left(z^k\right)\parallel$ and $\parallel F^{\prime }\left(z^k\right)F^{\prime }{\left(s^k\right)}^{-1}\parallel \le \kappa$ is not satisfied, we obtain from Step 4 that $d^k=d_1^k+{\beta }_k{d}_2^k$. Similarly to Case 1, it can be deduced that

$$\gamma {\chi }_{*}\ge \parallel F\left(z^{*}\right)\parallel$$

and then

$${\chi }_{*}\ge {\displaystyle \frac{\parallel F\left(z^{*}\right)\parallel }{\gamma }}>\parallel F\left(z^{*}\right)\parallel ,$$

which is a contradiction. Therefore, we have ${\xi }_{*}=0$. □

Theorem 5. 

If Assumption 1 is satisfied, then the sequence of iterations $\left\{z^k\right\}$ produced by Algorithm 1 converges to a solution to the WCP (2).

Proof of Theorem 5. 

From (15), we obtain

$$\parallel F\left(z^{k+1}\right){\parallel }^2\le (1+{\eta }_k){\parallel F\left(z^k\right)\parallel }^2.$$

(33)

Since ${\sum }_{k=0}^{\infty }{\eta }_k\le \xi <\infty$, $\{\parallel F\left(z^k\right){\parallel }^2\}$ is convergent according to Lemma 3.3 in [31], and $\{\parallel F\left(z^k\right)\parallel \}$ is also convergent as a result.

Suppose that $\underset{k\to \infty }{lim}{z}^k=z^{*}=({\xi }_{*},x^{*},y^{*},u^{*})$ without loss of generality. We only need to verify that $\parallel F\left(z^{*}\right)\parallel =0$. If not, we can obtain by a similar proof to that in Theorem 4,

$$1-\gamma -\kappa c-\gamma \kappa \le 0,$$

which is a contradiction, or

$$(1-\gamma )\parallel F\left(z^{*}\right){\parallel }^2\le 1,$$

which is a contradiction. Hence, $\parallel F\left(z^{*}\right)\parallel =0$. □

5. Local Convergence

We deal with the local biquadratic convergence property in this section.

Theorem 6. 

If Assumption 1 is true, all $D\in \partial F\left(z^{*}\right)$ are nonsingular. If the conditions hold that both $G^{\prime }(x,y,u)$ and $F^{\prime }(x,y,u)$ are of Lipschitz continuity near $z^{*}$, then $d^k=d_1^k+d_2^k$ for any sufficiently large k, and $\left\{z^k\right\}$ converges locally biquadratically to $z^{*}$.

Proof of Theorem 6. 

Since $z^{*}$ is the solution to the WCP (2), the Jacobian matrix $F^{\prime }\left(z^k\right)$ is invertible for any $z^k$ sufficiently close to $z^{*}$ according to Theorem 1. We obtain, for any sufficiently large k,

$$\parallel F^{\prime }{\left(z^k\right)}^{-1}\parallel =O\left(1\right)$$

(34)

from the condition that all $D\in \partial F\left(z^{*}\right)$ are nonsingular and Proposition 4.1 in [32]. Additionally, since $F\left(z\right)$ is strongly semismooth and locally Lipschitz continuous,

$$\parallel F(z^k)-F(z^{*})-F^{\prime }(z^k)(z^k-z^{*})\parallel =O(\parallel z^k-z^{*}{\parallel }^2)$$

(35)

and

$$\parallel F(z^k)\parallel = \parallel F(z^k)-F\left(z^{*}\right)\parallel =O(\parallel z^k-z^{*}\parallel ).$$

(36)

Since

$$\parallel {\chi }_{k}h\parallel \le \gamma {\xi }_k^4\le {\parallel F(z^k)\parallel }^4,$$

(37)

by (11) and (34)–(37), we have

$$\begin{array}{cc}\hfill \parallel s^k-z^{*}\parallel & = \parallel z^k+d_1^k-z^{*}\parallel \hfill \\ \hfill & = \parallel z^k+F^{\prime }{\left(z^k\right)}^{-1}(-F\left(z^k\right)+{\chi }_{k}h)-z^{*}\parallel \hfill \\ \hfill & =O\left(\parallel {\chi }_{k}h+F\left(z^k\right)-F\left(z^{*}\right)-F^{\prime }\left(z^k\right)(z^k-z^{*})\parallel \right)\hfill \\ \hfill & \le O(\parallel F\left(z^k\right){\parallel }^4)+O(\parallel z^k-z^{*}{\parallel }^2)\hfill \\ \hfill & =O(\parallel z^k-z^{*}{\parallel }^2),\hfill \end{array}$$

(38)

implying that $s^k$ is sufficiently close to $z^{*}$ as is $z^k$. By (36) and (38), we obtain

$$\parallel F\left(s^k\right)\parallel =O(\parallel s^k-z^{*}\parallel )=O(\parallel F\left(z^k\right){\parallel }^2).$$

(39)

Hence, combining (12), (34), (37), and (39) yields

$$\begin{array}{cc}\hfill \parallel d_2^k\parallel & =∥F^{\prime }{\left(y\right)}^{-1}\left(-F\left(s^k\right)+{\chi }_{k}h\right)∥\hfill \\ \hfill & \le O\left(\parallel F\left(s^k\right)\parallel +\parallel {\chi }_{k}h\parallel \right)\hfill \\ \hfill & =O(\parallel F\left(z^k\right){\parallel }^2),\hfill \end{array}$$

(40)

and combining with (36) and (38) yields

$$\begin{array}{cc}\hfill \parallel z^k+d_1^k+d_2^k-\Delta z^{*}\parallel & \le \parallel s^k-z^{*}\parallel +\parallel d_2^k\parallel \hfill \\ \hfill & =O(\parallel z^k-z^{*}{\parallel }^2),\hfill \end{array}$$

(41)

indicating that $z^k+d_1^k+d_2^k$ is sufficiently close to $z^{*}$ as is $z^k$.

Therefore, by combining (36) and (41), we have

$$\begin{array}{cc}\hfill & \parallel F(z^k+d_1^k+d_2^k)\parallel \hfill \\ \hfill =& O(\parallel z^k+d_1^k+d_2^k-z^{*}\parallel )\hfill \\ \hfill =& O(\parallel z^k-z^{*}{\parallel }^2)\hfill \\ \hfill =& o(\parallel F\left(z^k\right)\parallel )\hfill \\ \hfill =& l_k\parallel F\left(z^k\right)\parallel ,\hfill \end{array}$$

(42)

where $l_k\to 0$. This means that (13) holds, indicating that ${\beta }_k\equiv 1$ for any sufficiently large k, i.e.,

$$z^{k+1}=z^k+d_1^k+d_2^k.$$

(43)

Upon the condition that $F^{\prime }\left(z\right)$ is Lipschitz continuous, we have from (12), (34), (37), and (43) that

$$\begin{array}{cc}\hfill & \parallel z^{k+1}-z^{*}\parallel \hfill \\ \hfill =& \parallel s^k+d_2^k-z^{*}\parallel \hfill \\ \hfill =& \parallel s^k-z^{*}+F^{\prime }{\left(s^k\right)}^{-1}(-F\left(s^k\right)+{\chi }_{k}h)\parallel \hfill \\ \hfill \le & O\left(\parallel F\left(s^k\right)-F\left(z^{*}\right)-F^{\prime }\left(s^k\right)(s^k-z^{*})\parallel +\parallel {\chi }_{k}h\parallel \right)\hfill \\ \hfill =& O(\parallel s^k-z^{*}{\parallel }^2)+O(\parallel F\left(z^k\right){\parallel }^4)\hfill \\ \hfill =& O(\parallel z^k-z^{*}{\parallel }^4).\hfill \end{array}$$

(44)

Additionally, we have from (36) that

$$\parallel F\left(z^{k+1}\right)\parallel =O(\parallel z^{k+1}-z^{*}\parallel )=O(\parallel z^k-z^{*}{\parallel }^4)=O(\parallel F\left(z^k\right){\parallel }^4),$$

implying that $\left\{z^k\right\}$ converges locally biquadratically to $z^{*}$. □

6. Numerical Experiments

This section reports the computational efficiency of Algorithm 1, denoted as SNQ_L, for WLCPs and WCPs. The numerical experiments are performed on a PC with 16 GB RAM running MATLAB R2018b. In our tests, we set $n=2m$, $\delta ={10}^{-6}, \rho =0.6,$$\gamma =0.01, l=0.1, c=0.01, {\sigma }_1={\sigma }_2=0.001, {\xi }_0=0.1, \mathrm{and} {\eta }_k=1/2^{k+2}.$

Furthermore, we use the algorithm in [23], denoted as ANM_Tang, and that in [24], denoted as TSN_Liu, and compare them to SNQ_L. We utilize the same parameters for ANM_Tang and TSN_Liu as those in [23,24], respectively.

6.1. Numerical Tests for a WLCP

Consider the following WLCP,

$$x\in R_{+}^n,y\in R_{+}^n,xy=w,Ax+By+Cu=d,$$

where

$$A=\left(\begin{array}{c}M\\ N\end{array}\right),B=\left(\begin{array}{c}0\\ -I\end{array}\right),C=\left(\begin{array}{c}0\\ -M^T\end{array}\right),d=\left(\begin{array}{c}Mf\\ g\end{array}\right).$$

We set $M\in R^{m\times n}$ to be generated from a standard normal distribution and $N=P^{T}P/\parallel P^{T}P\parallel$, where $P\in R^{n\times n}$ is uniformly distributed over the interval $(0,1)$. The elements of $f\in R^n$ and $g\in R^n$, respectively, follow uniform distributions over the intervals $(0,1)$ and $(-1,0)$. Then, w is generated by $w=\widehat{x}\widehat{y}$ with $\widehat{y}=N\widehat{x}-g$ and $\widehat{x}=\mathtt{rand}(\mathtt{n},\mathtt{1})$. The starting points $(x^0,y^0,u^0)={(0,0,\dots ,0)}^T$. In the following, Aveero reflects the average value of $\parallel F\left(z^k\right)\parallel$ at the end of the iteration, Avek stands for the average number of iterations, and AveCPU stands for the average runtime of the associated algorithm in seconds.

To verify the local convergence rates of TSN_Liu, ANM_Tang, and SNQ_L, we first test them on a specific case with $n=3000$.

Table 1. Display of the variation in $\parallel F\left(z^k\right)\parallel$ as the number of iterations k increases.

k	SNQ_L	ANM_Tang	TSN_Liu
1	$3.5540\times {10}^1$	$1.5520\times {10}^1$	$6.6272\times {10}^0$
2	$3.7347\times {10}^{-2}$	$1.6755\times {10}^0$	$2.6368\times {10}^{-1}$
3	$4.2037\times {10}^{-8}$	$5.1659\times {10}^{-3}$	$3.6580\times {10}^{-4}$
4	∖	$3.2440\times {10}^{-7}$	$6.6823\times {10}^{-12}$

illustrates the variations in $\parallel F\left(z^k\right)\parallel$ with an increase in the number of iterations. From Table 1, it is evident that our algorithm SNQ_L exhibits local fourth-order convergence and, compared to ANM_Tang and TSN_Liu, it is able to converge to the solution more rapidly.

Next, we randomly generate 10 instances of varying size for testing on each problem. The test results are displayed in

Table 2. Numerical comparison results of the three algorithms for solving the WLCP.

n	SNQ_L			ANM_Tang			TSN_Liu
	Avek	AveCPU	Aveero	Avek	AveCPU	Aveero	Avek	AveCPU	Aveero
1000	3.0	1.4839	$2.2325\times {10}^{-7}$	3.5	1.5091	$9.4296\times {10}^{-8}$	4.0	1.8323	$6.1994\times {10}^{-11}$
2000	3.0	9.8266	$4.3843\times {10}^{-7}$	3.6	9.4651	$3.9083\times {10}^{-7}$	4.0	15.3614	$1.1060\times {10}^{-10}$
3000	3.0	36.4712	$2.9185\times {10}^{-7}$	4.0	37.7936	$3.2919\times {10}^{-7}$	4.0	51.4373	$1.6420\times {10}^{-10}$
4000	3.1	89.5122	$5.9255\times {10}^{-7}$	4.0	101.6152	$6.0811\times {10}^{-7}$	4.0	110.9661	$4.0618\times {10}^{-10}$
5000	3.2	165.3057	$4.9597\times {10}^{-7}$	4.0	160.4525	$7.7175\times {10}^{-7}$	4.0	216.3353	$1.4777\times {10}^{-10}$
6000	3.3	289.1295	$4.1117\times {10}^{-7}$	4.0	296.4240	$3.8819\times {10}^{-7}$	4.0	368.0030	$7.5609\times {10}^{-9}$
7000	3.4	478.9888	$7.6231\times {10}^{-11}$	4.0	758.7281	$1.1165\times {10}^{-11}$	4.0	584.5569	$5.9611\times {10}^{-11}$
8000	3.4	731.9228	$4.7637\times {10}^{-7}$	4.0	1174.9030	$1.3886\times {10}^{-11}$	4.0	871.0305	$3.6116\times {10}^{-10}$

, showing that SNQ_L consistently requires fewer iterations and typically utilizes less CPU time to achieve the stopping tolerance compared to ANM_Tang and TSN_Liu. Furthermore, we can observe that SNQ_L requires even fewer iterations and saves more CPU time compared to ANM_Tang and TSN_Liu as the problem size increases. This is due to the local fourth-order convergence rate exhibited by SNQ_L.

6.2. Numerical Tests for a WCP

Consider the following WCP,

$$x\in R_{+}^n,y\in R_{+}^n,xy=w,G(x,y,u)=0,$$

with

$$G(x,y,u)=\left(\begin{array}{c}Rx-P^{T}u-y+d\\ P(x-f)\end{array}\right),$$

where $R=D^{T}D/\parallel D^{T}D\parallel$, $D\in R^{n\times n}$ is uniformly distributed over the interval $(0,1)$, and $P\in R^{m\times n}$ is obtained from a uniform distribution on the interval $(0,1)$. The vectors $d,f,w\in R^n$ are all uniformly distributed over the interval $(0,1)$.

We perform 10 random experiments for each dimension. The average test results are displayed in

Table 3. Numerical comparison results of the three algorithms for solving the WCP.

n	SNQ_L			ANM_Tang			TSN_Liu
	Avek	AveCPU	Aveero	Avek	AveCPU	Aveero	Avek	AveCPU	Aveero
500	3.3	0.5099	$1.5306\times {10}^{-7}$	4.0	0.4196	$1.4293\times {10}^{-7}$	4.2	0.3881	$8.6678\times {10}^{-8}$
1000	3.7	2.0400	$1.2644\times {10}^{-7}$	4.0	2.1728	$1.2036\times {10}^{-7}$	4.7	1.9977	$3.1280\times {10}^{-8}$
1500	3.9	5.4378	$9.5575\times {10}^{-8}$	4.0	5.3286	$1.6225\times {10}^{-8}$	4.9	6.5876	$8.3102\times {10}^{-8}$
2000	3.8	12.2531	$1.0501\times {10}^{-7}$	4.3	17.5133	$1.2893\times {10}^{-7}$	4.9	17.1392	$2.4996\times {10}^{-7}$
2500	4.0	25.5304	$1.1410\times {10}^{-11}$	4.3	26.1118	$5.1913\times {10}^{-7}$	4.9	35.1112	$3.1818\times {10}^{-8}$
3000	4.0	44.2121	$6.2264\times {10}^{-9}$	5.0	42.3906	$1.0581\times {10}^{-7}$	5.1	56.0453	$5.7393\times {10}^{-8}$
3500	4.0	74.8664	$1.4966\times {10}^{-7}$	5.0	70.4032	$1.4966\times {10}^{-7}$	5.3	91.4523	$2.8774\times {10}^{-9}$
4000	4.0	109.1118	$5.8961\times {10}^{-11}$	6.0	149.8035	$1.6964\times {10}^{-8}$	5.4	139.3034	$4.1105\times {10}^{-8}$
4500	4.0	144.2821	$1.3936\times {10}^{-10}$	6.0	164.8998	$1.4673\times {10}^{-8}$	5.4	198.5842	$2.9454\times {10}^{-7}$
5000	4.0	210.1606	$6.9677\times {10}^{-9}$	7.0	259.7001	$7.6293\times {10}^{-9}$	5.2	278.2668	$4.1158\times {10}^{-8}$

. These numerical results also demonstrate that SNQ_L is more stable and effective compared to ANM_Tang and TSN_Liu. Furthermore, as the dimension of the problem increases, SNQ_L requires less time and fewer iterations.

7. Conclusions

We suggested a two-step Newton algorithm for the WCP by combining a two-step Newton approach for nonlinear systems of equations with a classical Newton approach for WCPs. The algorithm equivalently describes the WCP as a nonlinear set of equations. The new algorithm computes an additional Newton direction to obtain the next iteration point at each iteration. When the objective function value satisfies a specific descent requirement, the algorithm uses the iteration point created by the two Newton directions directly as the next iteration point; otherwise, the step size is decided by a derivative-free line search to find the next iteration point. Under certain assumptions, the global and local biquadratic convergence properties are verified. In numerical tests, we also compared the algorithm with a two-step Newton algorithm and an accelerated Newton algorithm, and the test results demonstrate that the computational efficiency of the new algorithm is effectively improved without increasing the time cost and that the algorithm has the property of local biquadratic convergence when the sequence of iterations is close to the solution to the WCP.