Mathematical Problems in Engineering
Volume 2020 (2020), Article ID 4072651, 9 pages
https://doi.org/10.1155/2020/4072651

Research Article

A Cauchy Point Direction Trust Region Algorithm for Nonlinear Equations

Zhou Sheng¹ and Dan Luo²

¹Department of Mathematics, Nanjing University of Aeronautics and Astronautics, Nanjing 211106, Jiangsu, China
²School of Mathematics and Statistics, Baise University, Baise 533000, Guangxi, China

Correspondence should be addressed to Dan Luo; bsxyldan@163.com

Received 13 July 2020; Accepted 10 August 2020; Published 24 August 2020

Guest Editor: Maojun Zhang

Copyright © 2020 Zhou Sheng and Dan Luo. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

In this paper, a Cauchy point direction trust region algorithm is presented to solve nonlinear equations. The search direction is an optimal convex combination of the trust region direction and the Cauchy point direction with the sufficiently descent property and the automatic trust region property. The global convergence of the proposed algorithm is proven under some conditions. The preliminary numerical results demonstrate that the proposed algorithm is promising and has better convergence behaviors than the other two existing algorithms for solving large-scale nonlinear equations.

1. Introduction

The system of nonlinear equations is one of the most important mathematical models, with extensive applications in chemical equilibrium problems [1], problem of hydraulic circuits of power plants [2], L-shaped beam structures [3], and image restoration problems [4]. Some -norm regularized optimization problems in compressive sensing [5] are also obtained by reformulating the systems of nonlinear equations. In this paper, consider the following nonlinear system of equation:

where is continuously differentiable. Let

where . It is clear that the nonlinear system of equations (1) is equivalent to the unconstrained optimization problem as follows:

Many literature studies have studied a variety of numerical algorithms to solve the nonlinear system of equations, such as the conjugate gradient algorithm [6–9], Levenberg–Marquardt algorithm [10–13], and trust region algorithm [14–21]. The trust region method plays an important role in the area of nonlinear optimization. It was proposed by Powell in [22] to solve the nonlinear optimization problem by using an iterative structure. Recently, some trust region-based methods have been developed to solve nonlinear equations and have shown very promising results. For instance, Fan and Pan [16] discussed a regularized trust region subproblem with an adaptive trust region radius, which possesses quadratic convergence under the local error bound condition. Yuan et al. [17] proposed a trust region-based algorithm, in which the Jacobian is updated by the BFGS formula, for solving nonlinear equations. The algorithm showed global convergence at a quadratic convergence rate. This algorithm does not compute the Jacobian matrix in each interrogation, which significantly reduces the computation burden. Esmaeili and Kimiaei [19] studied an adaptive trust region algorithm with the moderate radius size by employing the nonmonotone technique, and the global and quadratic local convergence was given. However, the maximal dimension of the test problems was only 300. Qi et al. [23] proposed an asymptotic Newton search direction method which is an adequate combination of the projected gradient direction and the projected trust region direction for solving box-constrained nonsmooth equations and proved the global and quadratic local convergence. Motivated by the search direction by Qi et al., Yuan et al. [24] used the idea of the search direction to solve nonsmooth unconstrained optimization problems and also showed global convergence.

Inspired by the fact that the search direction by Qi et al. [23] speeds up the local convergence rate compared to the gradient direction, we propose an efficient Cauchy point direction trust region algorithm to solve (1) in this paper. The search direction is an optimal convex combination of the trust region direction and the Cauchy point direction with the automatic trust region property and the sufficiently descent property, which can improve the numerical performance. We also show the global convergence of the proposed algorithm. Numerical results indicate that this algorithm has better behavior convergence compared with that of the classical trust region algorithm and the adaptive regularized trust region algorithm. However, the local convergence is not given.

The remainder of this paper is organized as follows. In Section 2, we introduce our motivation and a new search direction. In Section 3, we give a Cauchy point direction trust region algorithm for solving problem (1). The global convergence of the proposed algorithm is proven in Section 4. Preliminary numerical results are reported in Section 5. We conclude our paper in Section 6.

Throughout this paper, unless otherwise specified, denotes the Euclidean norm of vectors or matrices.

2. Motivation and Search Direction

In this section, motivated by the definition of the Cauchy point which is parallel to the steepest descent direction, we present a search direction with the Cauchy point. Let be the trust region radius and .

2.1. Cauchy Point Direction

Noting that is a semidefinite matrix, we define the Cauchy point at the current point :

where

, and is the Jacobian of or its approximation.

The following lemma gives a nice property of the Cauchy point direction .

Proposition 1 (automatic trust region property). Suppose that is defined by (4), then .

Proof. It is easy to get the result of this proposition from . This completes the proof.
From the definition of and the above proposition, we can determine the descent property and the automatic trust region property of the Cauchy point direction .

2.2. Trust Region Direction

Let the trial step be a solution of the following trust region subproblem:

where is the same as the above definition.

From the famous result of Powell [22], we give the following lemma and omit its proof.

Lemma 1. Let be a solution of the trust region subproblem (6), then

As mentioned by Qi et al. [23] and Yuan et al. [24], when is far from the optimization solution, the trust region direction may not possess the descent property. However, it follows from the definition of the Cauchy point direction that always possesses the descent property. On the contrary, if the search direction is defined by using a combination of and , it will be a descent direction. Hence, we give a new search direction by using a convex combination form in the next section.

2.3. Search Direction

Let be a solution of the one-dimensional minimization problem as follows:

and the search direction is defined as

For the above one-dimensional minimization problem, we can use the gold section method, Fibonacci method, parabolic method [25], etc. The idea of the search direction follows from Qi et al. [23]. However, we use the Cauchy point direction which possesses the automatic trust region property to get optimal by using the gold section method for problem (6) in this paper.

3. New Algorithm

In this section, let be the optimal solution of (3). Then, the actual reduction is defined as

and the predicted reduction is

Thus, is defined by the ratio between and :

We now list the detailed steps of the Cauchy point direction trust region algorithm in Algorithm 1.

Algorithm 1. Cauchy point direction trust region algorithm for solving (3).

 Step 1. Given initial point and parameters , , , , and . Let . Step 2. If the termination condition is satisfied at the iteration point , then stop. Otherwise, go to Step 3. Step 3. Solve the trust region subproblem (6) to obtain and compute and from (4) and (5), respectively. Step 4. To solve the one-dimensional minimization problem (8) giving by using Algorithm 2, calculate the search direction from (6). Step 5. Compute , if

 are satisfied; let , and update by

 Step 6. If , return to Step 3, and let . Otherwise, return to Step 2. Let .

In Step 4 of Algorithm 1, we use the gold section algorithm as given in Algorithm 2.

Algorithm 2. Gold section algorithm [25] for solving problem (8).

 Step 1. For the given initial parameters , compute

 and . Let . Step 2. If , go to Step 3. Otherwise, go to Step 4. Step 3. If the termination condition is satisfied, stop. Otherwise, let

 and compute and go to Step 5. Step 4. If the termination condition is satisfied, stop. Otherwise, let

 and compute and go to Step 5. Step 5. Let ; return to Step 2.

4. Global Convergence

In this section, the global convergence of Algorithm 1 is proven. Firstly, it is easy to have the following lemma from Lemma 1.

Lemma 2. Let be defined by (9), where is a solution of problem (8). Then,

Proof. From the definition of , we have . It follows from the definition of and Lemma 1 that holds. Thus, (16) holds, and therefore, the proof is completed.
The next lemma indicates an important property of the Cauchy point direction , which is the sufficiently descent property.

Lemma 3. holds, for any .

Proof. For any , from (9), we have

and note that ; thus,

This completes the proof.
To derive the global convergence of Algorithm 1, the following assumption is considered.

Assumption 1. The Jacobian of is bounded; i.e., there exists a positive constant such that , for any .
In the sequel, we show that “Steps 3–6” of Algorithm 1 do not cycle in the inner loop infinitely.

Lemma 4. Suppose that Assumption 1 holds, let be generated by Algorithm 1. If is not a stationary point of problem (3), then Algorithm 1 is well defined; i.e., the inner cycle of Algorithm 1 can be left in a finite number of internal iterates.

Proof. It follows from which is not a stationary point that there exists a positive constant such that

Since and Assumption 1 holds, there exists a positive constant such that

Letting , then from (4) and (22), we obtain

and noting (21), we have , which implies

and therefore, (13) holds. Next, we prove that (14) is also satisfied. By contradiction, we assume that “Step 3–Step 6” of Algorithm 1 do cycle in the inner loop infinitely; then, .
From the definition of , we get

It follows from (13) and Lemma 3 that there exists a positive constant such that

By the definition , we have , and noting Proposition 1 and , further, we get

Thus, from (26) and (27), we have

and therefore , which means that when is sufficiently small. However, the definition of in Step 5 of Algorithm 1 contradicts the fact that . Hence, the cycling between Step 3 and Step 6 terminates finitely. This completes the proof.
We can get the global convergence of Algorithm 1 based on the above lemmas.

Theorem 1. Suppose that Assumption 1 holds and be generated by Algorithm 1, then

Proof. Suppose that the conclusion is not true, then there exists a positive constant such that for any ,

From Lemma 2, Assumption 1, and the above inequation, we have

Noting the updating rule of , , (31), and Lemma 4, we deduce that

which yields a contradiction that this means (30) is not true; i.e., the result (29) of the theorem holds, and therefore, the proof is completed.

5. Numerical Results

In this section, the numerical results of the proposed algorithm for solving problem (1) are reported. We denote the proposed Algorithm 1 by CTR and compare it with some existing algorithms. The adaptive trust region algorithm by Fan and Pan [16] and the classical trust region algorithm [26] are denoted by ATR and TTR, respectively. In numerical experiments, all parameters of Algorithm 1 are chosen as follows: , , , , , and . We choose in Algorithm 2. Steihaug method [27] is employed to obtain the approximative trial step for solving the trust region subproblem (6). We will terminate algorithms if the number of iterations is larger than 2000.

All codes were implemented in MATLAB R2010b. The numerical experiments were performed on a PC with an Intel Pentium (R) Dual-Core CPU at 3.20 GHz and 2.00 GB of RAM and using the Windows 7 operating system.

We list the test problems in Table 1 that can also be found in [28]. In Tables 1–3, “No.” denotes the number of the test problem, “Problem name” denotes the name of the test problem, “” denotes the initial point, “” denotes the dimension of the problem, “NI” is the total number of iterations, “NF” is the number of the function evaluations, “” refers to the optimum function values, and “” denotes the mean value of . It is worth mentioning that “” represents the stopping of the algorithms in situations that the number of iterations exceeds the maximal iterations, but the terminal conditions cannot be satisfied yet.

Table 1: Test problems.

No.
Problem name

1
Exponential function 1

2
Exponential function 2

3
Extended Rosenbrock function

4
Trigonometric function

5
Singular function

6
Logarithmic function

7
Broyden tridiagonal function

8
Trigexp function

9
Function 15

10
Strictly convex function 1

11
Strictly convex function 2

12
Zero Jacobian function

13
Linear function-full rank

14
Penalty function

15
Brown almost-linear function

16
Extended Powell singular function

17
Tridiagonal system

18
Extended Freudenstein and Roth function

19
Discrete boundary value problem

20
Troesch problem

Table 2:  Numerical results.

No.
n
CTR
ATR
TTR

NI/NF

NI/NF

NI/NF

1
500
6/7
3.18E − 06
3.32E − 07
7/8
6.73E − 06
6/7
3.18E − 06

1000
5/6
8.95E − 06
3.32E − 07
7/8
4.74E − 06
5/6
8.95E − 06

2000
5/6
6.32E − 06
3.32E − 07
7/8
3.34E − 06
5/6
6.32E − 06

2
500
2/3
3.50E − 06
1.57E − 06
2/3
3.38E − 06
2/3
3.38E − 06

1000
2/3
6.67E − 06
3.32E − 07
2/3
6.69E − 06
2/3
6.69E − 06

2000
3/4
9.13E − 06
3.32E − 07
3/8
8.39E − 06
3/5
8.44E − 06

3
500
12/16
2.16E − 07
3.31E − 01
12/20
3.89E − 06
16/21
9.44E − 10

1000
11/15
5.92E − 07
2.96E − 01
12/20
5.51E − 06
16/21
1.39E − 09

2000
11/15
1.78E − 07
3.50E − 01
12/20
7.82E − 06
16/21
1.84E − 09

4
100
9/11
7.33E − 06
5.28E − 01
2001/2293

2001/2270

5
500
22/25
4.39E − 06
2.96E − 01
18/19
5.84E − 06
19/20
8.32E − 06

1000
25/27
6.37E − 06
4.35E − 01
23/58
5.68E − 06
23/24
4.33E − 06

2000
24/28
5.24E − 06
4.50E − 01
19/20
5.53E − 06
21/23
5.23E − 06

6
500
5/6
2.12E − 06
9.81E − 04
5/6
4.91E − 08
7/8
2.74E − 08

1000
5/6
1.44E − 06
5.05E − 01
5/6
1.91E − 07
7/8
1.58E − 07

2000
6/7
7.69E − 07
6.08E − 01
5/6
7.56E − 08
7/8
5.08E − 08

7
500
6/7
2.57E − 09
1.59E − 01
6/7
1.51E − 09
6/7
4.28E − 11

1000
6/7
2.63E − 06
1.63E − 01
6/7
1.03E − 08
6/7
3.83E − 06

2000
6/7
5.50E − 06
2.26E − 01
6/7
2.60E − 08
6/7
5.43E − 06

8
100
14/16
4.31E − 08
2.92E − 01
2001/2281

2001/2268

9
1000
14/18
5.42E − 06
1.46E − 02
14/20
2.76E − 07
16/22
1.36E − 09

10
500
6/7
5.54E − 08
4.14E − 01
6/7
1.33E − 07
6/7
1.09E − 07

1000
6/7
6.04E − 08
3.71E − 01
6/7
2.49E − 07
6/7
9.91E − 07

2000
6/7
1.16E − 07
4.47E − 01
6/7
4.33E − 07
7/8
2.32E − 07

11
100
8/9
3.22E − 07
2.21E − 01
8/9
3.38E − 08
8/9
1.15E − 06

12
500
5/6
3.62E − 12
8.36E − 01
16/17
5.81E − 06
18/19
5.96E − 06

1000
5/6
4.15E − 11
8.82E − 01
16/17
7.38E − 06
18/19
9.31E − 06

2000
6/8
3.42E − 07
8.29E − 01
16/17
8.23E − 06
19/20
2.82E − 06

13
500
9/10
2.41E − 07
7.19E − 01
4/5
7.67E − 06
9/10
1.04E − 08

1000
9/10
3.93E − 07
7.62E − 01
4/5
8.32E − 06
9/10
3.37E − 08

2000
10/11
7.01E − 06
7.30E − 01
5/6
1.00E − 10
10/11
1.35E − 07

14
500
4/5
6.96E − 06
1.09E − 02
7/10
5.22E − 06
5/6
5.34E − 06

1000
12/13
1.79E − 07
1.77E − 01
8/10
2.39E − 06
6/8
2.65E − 06

2000
10/11
5.39E − 06
1.15E − 01
11/16
1.70E − 06
10/11
5.32E − 07

15
500
3/4
2.85E − 09
3.61E − 01
3/4
2.48E − 06
4/5
2.49E − 06

1000
4/5
1.32E − 07
5.56E − 01
3/4
6.34E − 07
6/7
8.80E − 08

2000
4/5
1.40E − 08
5.04E − 01
4/5
3.74E − 07
6/7
3.63E − 07

16
500
1/2
1.10E − 06
3.32E − 07
1/2
1.10E − 06
1/2
1.10E − 06

1000
1/2
1.55E − 06
3.32E − 07
1/2
1.55E − 06
1/2
1.55E − 06

2000
1/2
2.20E − 06
3.32E − 07
1/2
2.20E − 06
1/2
2.20E − 06

17
1000
34/44
9.22E − 09
3.58E − 01
90/219
6.85E − 09
63/75
4.89E − 10

Table 3:  Numerical results (continued).

No.
n
CTR
ATR
TTR

NI/NF

NI/NF

NI/NF

18
500
6/8
3.57E − 07
4.13E − 01
7/8
5.84E − 08
7/8
6.12E − 09

1000
8/10
1.14E − 06
3.17E − 01
7/8
7.87E − 08
8/9
7.10E − 08

2000
8/10
3.08E − 06
3.17E − 01
7/8
1.07E − 07
8/9
2.29E − 08

19
500
2/3
9.77E − 09
3.32E − 07
3/4
6.77E − 09
2/3
9.77E − 09

1000
2/3
9.10E − 10
3.32E − 07
3/4
6.17E − 10
2/3
9.09E − 10

2000
2/3
4.47E − 10
3.32E − 07
3/4
3.02E − 10
2/3
4.47E − 10

20
1000
4/5
1.90E − 07
3.32E − 07
8/9
8.94E − 08
4/5
1.87E − 07

It is shown from Tables 2 and 3 that the proposed algorithm is efficient to solve these test problems and is competitive with the other two algorithms in NI and NF for most problems. For problems 1, 2, 16, 19, and 20, is close to zero, which implies that the trust region direction was used many times. Moreover, the number of iteration is the same as that of algorithm TTR, which means that our algorithm speeds up the local convergence of the iterates. It can be seen that is not too large for other test problems. However, for problems 4 and 8, algorithms ATR and TTR are not efficient to get the optimal solution when . We used the performance profiles by Dolan and Moré [29] to compare these algorithms with respect to NI and NF again, as shown in Figures 1 and 2. Therefore, it is noted from Tables 2 and 3 and Figures 1 and 2 that the proposed algorithm is well promising and has better performance in test problems in comparison with two existing algorithms. In order to investigate the behavior of the convergence rate of three algorithms, we plot the curve on problem 10 with and in Figure 3, where is the optimal point. Figure 3 indicates that Algorithm 1 converges faster than others.

Figure 1: Performance profiles of the algorithms for the test problems (NI).

Figure 2: Performance profiles of the algorithms for the test problems (NF).

(a)

(b)

Figure 3: The behavior of the convergence: (a) n = 1000 and (b) n = 2000.

6. Conclusions

In this paper, we presented a Cauchy point direction trust region algorithm for solving the nonlinear system of equations. This algorithm combines the trust region direction and the Cauchy point direction with the descent property and the automatic trust region. The optimal convex combination parameters were determined by using the gold section algorithm which is an exact one-dimensional linear search algorithm. The new algorithm was proven to have the descent property and the trust region property. We also showed the global convergence of the proposed algorithm under suitable condition. More importantly, the numerical results demonstrated the faster convergence of the proposed algorithm over two existing methods, which is very promising to solve large-scale nonlinear equations. However, the convergence rate of the proposed algorithm remains unclear. Therefore, the demonstration of the theoretical convergence rate of the proposed algorithm, application in image restoration problems, and compressive sensing will be the topics among others in our future study.

Data Availability

The data used to support the findings of this study are included within the article.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (11261006 and 11661009) and the Guangxi Natural Science Foundation (2020GXNSFAA159069).

References

K. Meintjes and A. P. Morgan, “Chemical equilibrium systems as numerical test problems,” ACM Transactions on Mathematical Software, vol. 16, no. 2, pp. 143–151, 1990.

A. A. Levin, V. Chistyakov, E. A. Tairov, and V. F. Chistyakov, “On application of the structure of the nonlinear equations system, describing hydraulic circuits of power plants, in computations,” Bulletin of the South Ural State University. Series “Mathematical Modelling, Programming and Computer Software, vol. 9, no. 4, pp. 53–62, 2016.

F. Georgiades, “Nonlinear equations of motion of L-shaped beam structures,” European Journal of Mechanics—A/Solids, vol. 65, pp. 91–122, 2017.

G. Yuan, T. Li, and W. Hu, “A conjugate gradient algorithm for large-scale nonlinear equations and image restoration problems,” Applied Numerical Mathematics, vol. 147, pp. 129–141, 2020.

Y. Xiao and H. Zhu, “A conjugate gradient method to solve convex constrained monotone equations with applications in compressive sensing,” Journal of Mathematical Analysis and Applications, vol. 405, no. 1, pp. 310–319, 2013.

Z. Papp and S. Rapajić, “FR type methods for systems of large-scale nonlinear monotone equations,” Applied Mathematics and Computation, vol. 269, pp. 816–823, 2015.

W. Zhou and F. Wang, “A PRP-based residual method for large-scale monotone nonlinear equations,” Applied Mathematics and Computation, vol. 261, pp. 1–7, 2015.

G. Yuan and M. Zhang, “A three-terms Polak-Ribière-Polyak conjugate gradient algorithm for large-scale nonlinear equations,” Journal of Computational and Applied Mathematics, vol. 286, pp. 186–195, 2015.

G. Yuan, Z. Meng, and Y. Li, “A modified Hestenes and Stiefel conjugate gradient algorithm for large-scale nonsmooth minimizations and nonlinear equations,” Journal of Optimization Theory and Applications, vol. 168, no. 1, pp. 129–152, 2016.

J. Fan and J. Pan, “Convergence properties of a self-adaptive Levenberg-Marquardt algorithm under local error bound condition,” Computational Optimization and Applications, vol. 34, no. 1, pp. 47–62, 2006.

J. Fan, “Accelerating the modified Levenberg-Marquardt method for nonlinear equations,” Mathematics of Computation, vol. 83, no. 287, pp. 1173–1187, 2013.

X. Yang, “A higher-order Levenberg-Marquardt method for nonlinear equations,” Applied Mathematics and Computation, vol. 219, no. 22, pp. 10682–10694, 2013.

K. Amini and F. Rostami, “A modified two steps Levenberg-Marquardt method for nonlinear equations,” Journal of Computational and Applied Mathematics, vol. 288, pp. 341–350, 2015.

J.-l. Zhang and Y. Wang, “A new trust region method for nonlinear equations,” Mathematical Methods of Operations Research (ZOR), vol. 58, no. 2, pp. 283–298, 2003.

J. Fan, “Convergence rate of the trust region method for nonlinear equations under local error bound condition,” Computational Optimization and Applications, vol. 34, no. 2, pp. 215–227, 2006.

J. Fan and J. Pan, “An improved trust region algorithm for nonlinear equations,” Computational Optimization and Applications, vol. 48, no. 1, pp. 59–70, 2011.

G. Yuan, Z. Wei, and X. Lu, “A BFGS trust-region method for nonlinear equations,” Computing, vol. 92, no. 4, pp. 317–333, 2011.

H. Esmaeili and M. Kimiaei, “A new adaptive trust-region method for system of nonlinear equations,” Applied Mathematical Modelling, vol. 38, no. 11-12, pp. 3003–3015, 2014.

H. Esmaeili and M. Kimiaei, “An efficient adaptive trust-region method for systems of nonlinear equations,” International Journal of Computer Mathematics, vol. 92, no. 1, pp. 151–166, 2015.

M. Kimiaei and H. Esmaeili, “A trust-region approach with novel filter adaptive radius for system of nonlinear equations,” Numerical Algorithms, vol. 73, no. 4, pp. 999–1016, 2016.

H. Esmaeili and M. Kimiaei, “A trust-region method with improved adaptive radius for systems of nonlinear equations,” Mathematical Methods of Operations Research, vol. 83, no. 1, pp. 109–125, 2016.

M. J. D. Powell, “Convergence properties of a class of minimization algorithms,” in Nonlinear Programming, pp. 1–27, Elsevier, Amsterdam, Netherlands, 1975.

L. Qi, X. J. Tong, and D. H. Li, “Active-set projected trust-region algorithm for box-constrained nonsmooth equations,” Journal of Optimization Theory and Applications, vol. 120, no. 3, pp. 601–625, 2004.

G. Yuan, Z. Wei, and Z. Wang, “Gradient trust region algorithm with limited memory BFGS update for nonsmooth convex minimization,” Computational Optimization and Applications, vol. 54, no. 1, pp. 45–64, 2013.

Y. X. Yuan and W. Sun, Optimization Theory and Methods, Science Press, Beijing, China, 1997.

A. R. Conn, N. I. M. Gould, and P. L. Toint, Trust-region Methods, Society for Industrial and Applied Mathematics, Philadelphia, PA, USA, 2000.

T. Steihaug, “The conjugate gradient method and trust regions in large scale optimization,” SIAM Journal on Numerical Analysis, vol. 20, no. 3, pp. 626–637, 1983.

W. La Cruz, J. M. Martínez, and M. Raydan, “Spectral residual method without gradient information for solving large-scale nonlinear systems of equations,” Mathematics of Computation, vol. 75, no. 255, pp. 1429–1449, 2006.

E. D. Dolan and J. J. Moré, “Benchmarking optimization software with performance profiles,” Mathematical Programming, vol. 91, no. 2, pp. 201–213, 2002.