A Julia implementation of Algorithm NCL for constrained optimization
CCahier du GERAD G-2021-02
A Julia implementation of Algorithm NCLfor constrained optimization
January 27, 2021Ding Ma , Dominique Orban , and Michael A. Saunders Department of Management Science and Department of Marketing,College of Business, City University of Hong Kong, Hong Kong, [email protected] GERAD and Department of Mathematics and Industrial Engineering,Polytechnique Montréal, QC, Canada, [email protected] , https://dpo.github.io Systems Optimization Laboratory, Department of Management Science andEngineering, Stanford University, Stanford, CA, USA. [email protected] , http://stanford.edu/~saunders
Abstract.
Algorithm NCL is designed for general smooth optimizationproblems where first and second derivatives are available, including prob-lems whose constraints may not be linearly independent at a solution (i.e.,do not satisfy the LICQ). It is equivalent to the LANCELOT augmentedLagrangian method, reformulated as a short sequence of nonlinearlyconstrained subproblems that can be solved efficiently by IPOPT andKNITRO, with warm starts on each subproblem. We give numericalresults from a Julia implementation of Algorithm NCL on tax policymodels that do not satisfy the LICQ, and on nonlinear least-squaresproblems and general problems from the CUTEst test set.
Keywords:
Constrained optimization, second derivatives, Algorithm
NCL , Julia
Algorithm NCL (nonlinearly constrained augmented Lagrangian [15]) is designedfor smooth, constrained optimization problems for which first and second deriva-tives are available. Without loss of generality, we take the problem to beNCO minimize x ∈ R n φ ( x ) subject to c ( x ) = 0 , (cid:96) ≤ x ≤ u, where φ ( x ) is a scalar objective function and c ( x ) ∈ R m is a vector of linear ornonlinear constraints. Inequality constraints are accommodated by including slackvariables within x . We take the primal and dual solutions to be ( x ∗ , y ∗ , z ∗ ) . Wedenote the objective gradient by g ( x ) = ∇ φ ( x ) ∈ R n , and the constraint Jacobianby J ( x ) ∈ R m × n . The objective and constraint Hessians are H i ( x ) ∈ R n × n , i = 0 , , . . . , m . Cahier du GERAD G-2021-02 a r X i v : . [ m a t h . O C ] J a n D. Ma, D. Orban, and M. A. Saunders If J ( x ∗ ) has full row rank m , problem NCO satisfies the linear independenceconstraint qualification (LICQ) at x ∗ . Most constrained optimization solvershave difficulty if NCO does not satisfy the LICQ. An exception is LANCELOT [4,5,13]. Algorithm NCL inherits this desirable property by being equivalent tothe
LANCELOT algorithm. Assuming first and second derivatives are available,Algorithm NCL may be viewed as an efficient implementation of the
LANCELOT algorithm. Previously we have implemented Algorithm NCL in
AMPL [1,6,15] fortax policy problems [11,15] that could not otherwise be solved. Here we describeour implementation in Julia [3] and give results on the tax problems and on aset of nonlinear least-squares problems from the
CUTEst test set [9].
For problem NCO,
LANCELOT implements what we call a
BCL algorithm (bound-constrained augmented Lagrangian algorithm), which solves a sequence of about10 bound-constrained subproblemsBC k minimize x L ( x, y k , ρ k ) = φ ( x ) − y Tk c ( x ) + ρ k (cid:107) c ( x ) (cid:107) subject to (cid:96) ≤ x ≤ u for k = 0 , , , . . . , where y k is an estimate of the dual variable associated with c ( x ) = 0 , and ρ k > is a penalty parameter. Each subproblem BC k is solved(approximately) with a decreasing optimality tolerance ω k , giving an iterate ( x ∗ k , z ∗ k ) . If (cid:107) c ( x ∗ k ) (cid:107) is no larger than a decreasing feasibility tolerance η k , the dualvariable is updated to y k +1 = y k − ρ k c ( x ∗ k ) . Otherwise, the penalty parameter isincreased to ρ k +1 > ρ k .Optimality is declared if c ( x ∗ k ) ≤ η k and η k , ω k have already been decreasedto specified minimum values η ∗ , ω ∗ . Infeasibility is declared if c ( x ∗ k ) > η k and ρ k has already been increased to a specified maximum value ρ ∗ .If n is large and not many bounds are active at x ∗ , the BC k subproblemshave many degrees of freedom , and LANCELOT must optimize in high-dimensionsubspaces. The subproblems are therefore computationally expensive. The al-gorithm in
MINOS [16] (we call it an LCL algorithm) reduces this expense byincluding linearizations of the constraints within its subproblems:LC k minimize x L ( x, y k , ρ k ) = φ ( x ) − y Tk c ( x ) + ρ k (cid:107) c ( x ) (cid:107) subject to c ( x ∗ k − ) + J ( x ∗ k − )( x − x ∗ k − ) = 0 , (cid:96) ≤ x ≤ u. The SQP algorithm in
SNOPT [8] solves subproblems with the same linearizedconstraints and a quadratic approximation to the LC k objective. Complicationsarise for both MINOS and
SNOPT if the linearized constraints are infeasible. Available from https://github.com/optimizers/ncl
Cahier du GERAD G-2021-02
Julia Implementation of Algorithm NCL Algorithm NCL proceeds in the opposite way by introducing additionalvariable r ≡ − c ( x ) into subproblems LC k to obtain the NCL subproblemsNC k minimize x, r φ ( x ) + y Tk r + ρ k (cid:107) r (cid:107) subject to c ( x ) + r = 0 , (cid:96) ≤ x ≤ u. These subproblems have nonlinear constraints and far more degrees of freedomthan the original NCO! Indeed, the extra variables r make the subproblems moredifficult if they are solved by MINOS and
SNOPT . However, the subproblems satisfythe LICQ because of r . Also, interior solvers such as IPOPT [10] and
KNITRO [12] find r helpful because at each interior iteration p they update the currentprimal-dual point ( x p , r p , λ p ) by computing a search direction ( ∆x, ∆r, ∆λ ) froma linear system of the form ( H p + D p ) J Tp ρ k I IJ p I ∆x∆r∆λ = − g ( x p ) + J Tp λ p − z p + w p y k + ρ k r p + λ p c ( x p ) + r p , (1)where D p , z p and w p are an ill-conditioned positive-definite diagonal matrixand two vectors arising from the interior method, and each Lagrangian Hessian H p = H ( x p ) − (cid:80) i ( y k ) i H i ( x p ) may be altered to be more positive definite. Directmethods for solving each sparse system (1) are affected very little by the higherdimension caused by r , and they benefit significantly from (cid:0) J p I (cid:1) always havingfull row rank.If an optimal solution for NC k is ( x ∗ k , r ∗ k , y ∗ k , z ∗ k ) and the feasibility and opti-mality tolerances have decreased to their minimum values η ∗ and ω ∗ , a naturalstopping condition for Algorithm NCL is (cid:107) r ∗ k (cid:107) ∞ ≤ η ∗ , because the major itera-tions drive r toward zero and we see that if r = 0 , subproblem NC k is equivalentto the original problem NCO.We have found that Algorithm NCL is successful in practice because – there are only about 10 major iterations ( k = 1 , , . . . , ); – the search-direction computation (1) for interior solvers is more stable thanif the solvers are applied to NCO directly; – IPOPT and
KNITRO have run-time options that facilitate warm starts foreach subproblem NC k , k > . The above observations were confirmed by our
AMPL implementation of AlgorithmNCL in solving some large problems modeling taxation policy [11,15,17]. Theproblems have very many nonlinear inequality constraints c ( x ) ≥ in relativelyfew variables. They have the formTAX maximize c, y (cid:80) i λ i U i ( c i , y i ) subject to U i ( c i , y i ) − U i ( c j , y j ) ≥ for all i, jλ T ( y − c ) ≥ c, y ≥ , Cahier du GERAD G-2021-02
D. Ma, D. Orban, and M. A. Saunders where c i and y i are the consumption and income of taxpayer i , and λ is a vectorof positive weights. The utility functions U i ( c i , y i ) are each of the form U ( c, y ) = ( c − α ) − /γ − /γ − ψ ( y/w ) /η +1 /η + 1 , where w is the wage rate and α , γ , ψ and η are taxpayer heterogeneities. Moreprecisely, the utility functions are of the form U i,j,k,g,h ( c p,q,r,s,t , y p,q,r,s,t ) = ( c p,q,r,s,t − α k ) − /γ h − /γ h − ψ g ( y p,q,r,s,t /w i ) /η j +1 /η j + 1 , where ( i, j, k, g, h ) and ( p, q, r, s, t ) run over na wage types, nb elasticities of laborsupply, nc basic need types, nd levels of distaste for work, and ne elasticities ofdemand for consumption, with na , nb , nc , nd , ne determining the size of theproblem, namely m = T ( T − nonlinear constraints, n = 2 T variables, with T := na × nb × nc × nd × ne .To achieve reliability, we found it necessary to extend the AMPL model’sdefinition of U ( c, y ) to be a piecewise-continuous function that accommodatesnegative values of ( c − α ) .At a solution, a large proportion of the constraints are essentially active. Thefailure of LICQ causes numerical difficulties for MINOS , SNOPT , and
IPOPT . LANCELOT is more able to find a solution, except it is very slow on eachsubproblem NC k . For example, on the smallest problem of Table 2 with 32220constraints and 360 variables, LANCELOT running on NEOS [18] timed-out at anear-optimal point on the 11th major iteration after 8 hours of CPU.Note that when the constraints of NCO are inequalities c ( x ) ≥ as in problemTAX, the constraints of subproblem NC k become inequalities c ( x ) + r ≥ (andsimilarly for mixtures of equalities and inequalities). The inequalities mean “morefree variables” (more variables that are not on a bound). This increases theproblem difficulty for MINOS and
SNOPT , but has only a positive effect on theinterior solvers.In wishing to improve the efficiency of Algorithm NCL on larger tax problems,we found it possible to warm-start
IPOPT and
KNITRO on each NC k subproblem( k > ) by setting the run-time options shown in Table 1. These options were usedby NCL / IPOPT and
NCL / KNITRO to obtain the results in Table 2. We see that
NCL / IPOPT performed significantly better than
IPOPT itself, and similarly for
NCL / KNITRO compared to
KNITRO . The feasibility and optimality tolerances η k , ω k were fixed at η ∗ = ω ∗ = for all k . Our Julia implementation savescomputation by starting with larger η k , ω k and reducing them toward η ∗ , ω ∗ asin LANCELOT . In this section, ( c, y, λ ) refer to problem TAX, not the variables in Algorithm NCL. Cahier du GERAD G-2021-02
Julia Implementation of Algorithm NCL Table 1.
Run-time options for warm-starting IPOPT and KNITRO on subproblemNC k . IPOPT KNITRO k = 1 algorithm=1 k ≥ warm_start_init_point=yes bar_directinterval=0bar_initpt=2bar_murule=1 k = 2 , mu_init=1e-4 bar_initmu=1e-4bar_slackboundpush=1e-4 k = 4 , mu_init=1e-5 bar_initmu=1e-5bar_slackboundpush=1e-5 k = 6 , mu_init=1e-6 bar_initmu=1e-6bar_slackboundpush=1e-6 k = 8 , mu_init=1e-7 bar_initmu=1e-7bar_slackboundpush=1e-7 k ≥ mu_init=1e-8 bar_initmu=1e-8bar_slackboundpush=1e-8 Modeling languages such as
AMPL and
GAMS are domain-specific languages, asopposed to full-fledged, general-purpose programming languages like C or Java.In the terminology of Bentley [2], they are little languages . As such, they haveunderstandable, yet very real limitations, that make it difficult, impractical, andperhaps even impossible, to implement an algorithm such as Algorithm
NCL ina sufficiently generic manner so that it may be applied to arbitrary problems.Indeed, our
AMPL implementation of Algorithm
NCL is specific to the optimaltax policy problems, and it would be difficult to generalize it to other problems.One of the main motivations for implementing Algorithm
NCL in a languagesuch as Julia is to be able to solve a greater variety of optimization problems.We now describe the key features of our Julia implementation of Algorithm
NCL and show that it solves examples of the same tax problems more efficiently.We then give results on a set of nonlinear least-squares problems from the
CUTEst test set to indicate that Algorithm
NCL is a reliable solver for such problemswhere first and second derivatives are available for the interior solvers usedat each major iteration. To date, this means that Algorithm
NCL is effectivefor optimization problems modeled in
AMPL , GAMS , and
CUTEst . (We havenot made an implementation in
GAMS [7], but it would be possible to build amajor-iteration loop around calls to
IPOPT or KNITRO in the way that we didfor
AMPL [15].)
The main advantage of a Julia implementation over our original
AMPL imple-mentation is that we may take full advantage of our Julia software suite for
Cahier du GERAD G-2021-02
D. Ma, D. Orban, and M. A. Saunders
Table 2.
Solution of tax problems of increasing dimension using IPOPT and KNITROon the original problem (cold starts) and the AMPL implementation of Algorithm NCLwith IPOPT or KNITRO as subproblem solvers (warm-starting with the options inTable 1). The problem size increases with a problem parameter na . Other problemparameters are fixed at nb = nc = 3 , nd = ne = 2 . There are m nonlinear inequalityconstraints and n variables. For IPOPT, > indicates optimality was not achieved.IPOPT KNITRO NCL/IPOPT NCL/KNITRO na m n itns time itns time itns time itns time5 32220 360 449 217 168 53 322 146 339 639 104652 648 > >
928 825 655 1023 307 23911 156420 792 > > optimization, hosted under the JuliaSmoothOptimizers (JSO) organization [22].Our suite provides a general consistent API for solvers to interact with models byproviding flexible data types to represent the objective and constraint functions,to evaluate their derivatives, to examine bounds on the variables, to add slackvariables transparently, and to provide essentially any information that a solvermight request from a model. Thanks to interfaces to modeling languages such as
AMPL , CUTEst and
JuMP [14], solvers in JSO may be written without regardfor the language in which the model was written.The modules from our suite that are particularly useful in the context of ourimplementation of Algorithm NCL are the following. – NLPModels [24] is the main modeling package that defines the API on whichsolvers can rely to interact with models. Models are represented as instancesof a data type deriving from the base type
AbstractNLPModel , and solverscan evaluate the objective value by calling the obj() method, the gradientvector by calling the grad() method, and so forth. The main advantage of theconsistent API provided by NLPModels is that solvers need not worry aboutthe provenance of models. Other modules ensure communication betweenmodeling languages such as
AMPL , CUTEst or JuMP , and NLPModels. – AmplNLReader [20] is one such module, and, as the name indicates, allowsa solver written in Julia to interact with a model written in
AMPL . Thecommunication is made possible by the
AMPL
Solver Library (ASL) , whichrequires that the model be decoded as an nl file. – NLPModelsIpopt [23] is a thin translation layer between the low-level Juliainterface to
IPOPT provided by the IPOPT.jl package and NLPModels, andlets users solve any problem conforming to the NLPModels API with IPOPT . https://github.com/jump-dev/Ipopt.jl Cahier du GERAD G-2021-02
Julia Implementation of Algorithm NCL – NLPModelsKnitro [25] is similar to NLPModelsIpopt, but lets users solveproblems with
KNITRO via the low-level interface provided by KNITRO.jl .Julia is a convenient language built on top of state-of-the-art infrastructureunderlying modern compilers such as Clang. Julia may be used as an interactivelanguage for exploratory work in a read-eval-print loop similar to Matlab. However,Julia functions are transparently translated to low-level code and compiled thefirst time they are called. The net result is efficient compiled code whose efficiencyrivals that of binaries generated from standard compiled languages such as C andFortran. Though this last feature is not particularly important in the context ofAlgorithm NCL because the compiled solvers IPOPT and
KNITRO perform allthe work, it is paramount when implementing pure Julia optimization solvers.
The Julia implementation of Algorithm NCL, named NCL.jl [19], is in two parts.The first part defines a data type
NCLModel that derives from the basic datatype
AbstractNLPModel mentioned earlier and represents subproblem NC k . An NCLModel is a wrapper around the underlying problem NCO in which the currentvalues of ρ k and r can be updated efficiently. The second part is the solver itself,each iteration of which consists of a call to IPOPT or KNITRO , and parameterupdates. The solver takes an
NCLModel as input. If the input problem is not an
NCLModel , it is first converted into one. Parameters are initialized as η = 10 , ω = 10 , ρ = 100 , µ = 0 . , where µ is the initial barrier parameter for IPOPT or KNITRO . The initialvalues of x are those defined in the underlying model if any, or zero otherwise.We initialize r to zero and y to the vector of ones. When the subproblem solverreturns with NC k solution ( x ∗ k , r ∗ k , y ∗ k , z ∗ k ) , we check whether (cid:107) r ∗ k (cid:107) ≤ max( η k , η ∗ ) .If so, we decide that good progress has been made toward feasibility and update y k +1 = y k − ρ k r ∗ k , η k +1 = η k / , ω k +1 = ω k / , ρ k +1 = ρ k , where this definition of y k +1 is the first-order update of the multipliers. Otherwise,we keep most things the same but increase the penalty parameter: y k +1 = y k , η k +1 = η k , ω k +1 = ω k , ρ k +1 = min(10 ρ k , ρ ∗ ) , where ρ ∗ > is the threshold beyond which the user is alerted that the problemmay be infeasible. In our implementation, we use ρ ∗ = 10 .Note that updating the multipliers based on (cid:107) r ∗ k (cid:107) instead of (cid:107) c ( x ∗ k ) (cid:107) is adeparture from the classical augmented-Lagrangian update. From the optimalityconditions for NC k we can prove that the first-order update is equivalent tochoosing y k +1 = y ∗ k when NC k is solved accurately. We still have a choice betweenthe two updates because we use low accuracy for the early NC k . We could also https://github.com/jump-dev/KNITRO.jl Cahier du GERAD G-2021-02
D. Ma, D. Orban, and M. A. Saunders “trim” y k +1 (i.e., for inequality constraints c i ( x ) + r i ≥ or ≤ , set componentsof y k +1 with non-optimal sign to zero). These are topics for future research.With IPOPT as subproblem solver, we warm-start subproblem NC k +1 withthe options in Table 1 and ( y ∗ k , z ∗ k ) as initial values for the Lagrange multipliers.With KNITRO as subproblem solver, ( y ∗ k , z ∗ k ) as starting point did not help orharm KNITRO significantly. We allowed
KNITRO to determine its own initialmultipliers, and it proved to be significantly more reliable than
IPOPT in solvingthe NC k subproblems for the optimal tax policy problems. In the next sections,Algorithm NCL means our Julia implementation with KNITRO as subproblemsolver.
AMPL models of the optimal tax policy problems were input to the Julia im-plementation of Algorithm NCL. The notation 1D, 2D, 3D, 4D, 5D refers toproblem parameters na , nb , nc , nd , ne that define the utility function appearingin the objective and constraints. The subproblem solver was KNITRO
12 [12].Tables 3–7 illustrate that, as with our
AMPL implementation of AlgorithmNCL, about 10 major iterations are needed independent of the problem size. (Theproblems have increasing numbers of variables and greatly increasing numbers ofnonlinear inequality constraints.) In each iteration log, outer and inner refer to the NCL major iteration number k and thetotal number of KNITRO iterations for subproblems NC k ; NCL obj is the augmented Lagrangian objective value, which convergesto the objective value for the model; η and ω show the KNITRO feasibility and optimality tolerances η k and ω k decreasing from − to − ; (cid:107)∇ L (cid:107) is the size of the augmented Lagrangian gradient, namely (cid:107) g ( x ∗ k ) − J ( x ∗ k ) T y k +1 (cid:107) (a measure of the dual infeasibility at the end of majoriteration k ); ρ is the penalty parameter ρ k ; µ init is the initial value of KNITRO ’s barrier parameter; (cid:107) x (cid:107) is the size of the primal variable x ∗ k at the (approximate) solution ofNC k ; (cid:107) y (cid:107) is the size of the corresponding dual variable y ∗ k ; time is the number of seconds to solve NC k .We see from the decreasing inner iteration counts that KNITRO was able towarm-start each subproblem, and from the decreasing (cid:107) r (cid:107) and (cid:107)∇ L (cid:107) values thatit is sufficient to solve the early subproblems with low (but steadily increasing)accuracy. Cahier du GERAD G-2021-02
Julia Implementation of Algorithm NCL Table 3.
Tax1D problem with realistic data. NCL with KNITRO solving subproblems. (cid:7) (cid:4) julia> using NCLjulia> using AmplNLReaderjulia> tax1D = AmplModel("data/tax1D")Maximization problem data/tax1Dnvar = 24, ncon = 133 (1 linear)julia> NCLSolve(tax1D, outlev=0)outer inner NCL obj (cid:107) r (cid:107) η (cid:107)∇ L (cid:107) ω ρ µ init (cid:107) y (cid:107) (cid:107) x (cid:107) time1 5 -8.00e+02 9.7e-02 1.0e-02 7.6e-03 1.0e-02 1.0e+02 1.0e-01 1.0e+00 2.0e+02 0.132 12 -7.89e+02 4.2e-02 1.0e-02 4.3e-03 1.0e-02 1.0e+03 1.0e-03 1.0e+00 1.9e+02 0.003 7 -7.83e+02 5.7e-03 1.0e-02 1.0e-03 1.0e-02 1.0e+04 1.0e-03 1.0e+00 1.9e+02 0.004 3 -7.82e+02 1.3e-04 1.0e-03 1.0e-05 1.0e-03 1.0e+04 1.0e-05 5.8e+01 1.9e+02 0.005 2 -7.82e+02 2.3e-06 1.0e-04 1.0e-05 1.0e-04 1.0e+04 1.0e-05 5.9e+01 1.9e+02 0.006 2 -7.82e+02 9.3e-08 1.0e-05 1.0e-06 1.0e-05 1.0e+04 1.0e-06 5.9e+01 1.9e+02 0.007 2 -7.82e+02 7.7e-09 1.0e-06 1.0e-08 1.0e-06 1.0e+04 1.0e-06 5.9e+01 1.9e+02 0.00 (cid:6) (cid:5) Table 4.
Tax2D problem. NCL with KNITRO solving subproblems. (cid:7) (cid:4) julia> tax2D = AmplModel("data/tax2D")Maximization problem data/tax2Dnvar = 120, ncon = 3541 (1 linear)julia> NCLSolve(tax2D, outlev=0)outer inner NCL obj (cid:107) r (cid:107) η (cid:107)∇ L (cid:107) ω ρ µ init (cid:107) y (cid:107) (cid:107) x (cid:107) time1 16 -4.35e+03 6.1e-02 1.0e-02 4.2e-03 1.0e-02 1.0e+02 1.0e-01 1.0e+00 4.0e+02 0.152 15 -4.31e+03 2.5e-02 1.0e-02 2.7e-04 1.0e-02 1.0e+03 1.0e-03 1.0e+00 4.0e+02 0.133 16 -4.29e+03 7.8e-03 1.0e-02 3.5e-04 1.0e-02 1.0e+04 1.0e-03 1.0e+00 4.0e+02 0.164 15 -4.28e+03 5.1e-03 1.0e-03 1.0e-05 1.0e-03 1.0e+04 1.0e-05 7.9e+01 4.0e+02 0.145 32 -4.28e+03 1.2e-03 1.0e-03 1.0e-05 1.0e-03 1.0e+05 1.0e-05 7.9e+01 4.0e+02 0.326 12 -4.28e+03 1.5e-04 1.0e-03 1.5e-05 1.0e-03 1.0e+06 1.0e-06 7.9e+01 4.0e+02 0.157 4 -4.28e+03 1.8e-05 1.0e-04 2.7e-06 1.0e-04 1.0e+06 1.0e-06 2.0e+02 4.0e+02 0.068 4 -4.28e+03 1.2e-06 1.0e-05 1.3e-07 1.0e-05 1.0e+06 1.0e-07 2.0e+02 4.0e+02 0.059 3 -4.28e+03 3.5e-07 1.0e-06 1.0e-07 1.0e-06 1.0e+06 1.0e-07 2.0e+02 4.0e+02 0.05 (cid:6) (cid:5) Table 5.
Tax3D problem. NCL with KNITRO solving subproblems. (cid:7) (cid:4) julia> pTax3D = AmplModel("data/pTax3D")Maximization problem data/pTax3Dnvar = 216, ncon = 11557 (1 linear)julia> NCLSolve(pTax3D, outlev=0)outer inner NCL obj (cid:107) r (cid:107) η (cid:107)∇ L (cid:107) ω ρ µ init (cid:107) y (cid:107) (cid:107) x (cid:107) time1 9 -6.97e+03 4.5e-02 1.0e-02 9.1e-03 1.0e-02 1.0e+02 1.0e-01 1.0e+00 5.7e+02 0.542 18 -6.87e+03 1.7e-02 1.0e-02 2.4e-04 1.0e-02 1.0e+03 1.0e-03 1.0e+00 5.6e+02 0.993 16 -6.83e+03 7.8e-03 1.0e-02 1.7e-03 1.0e-02 1.0e+04 1.0e-03 1.0e+00 5.7e+02 1.014 17 -6.81e+03 5.2e-03 1.0e-03 1.5e-05 1.0e-03 1.0e+04 1.0e-05 7.9e+01 5.6e+02 0.995 54 -6.80e+03 2.6e-03 1.0e-03 1.2e-05 1.0e-03 1.0e+05 1.0e-05 7.9e+01 5.6e+02 3.156 22 -6.80e+03 4.5e-04 1.0e-03 8.0e-05 1.0e-03 1.0e+06 1.0e-06 7.9e+01 5.6e+02 1.307 9 -6.80e+03 1.1e-04 1.0e-04 1.0e-06 1.0e-04 1.0e+06 1.0e-06 5.2e+02 5.6e+02 0.568 8 -6.80e+03 1.1e-05 1.0e-04 1.1e-07 1.0e-04 1.0e+07 1.0e-07 5.2e+02 5.6e+02 0.499 5 -6.80e+03 1.1e-06 1.0e-05 1.0e-07 1.0e-05 1.0e+07 1.0e-07 5.3e+02 5.6e+02 0.3210 3 -6.80e+03 8.9e-08 1.0e-06 1.0e-08 1.0e-06 1.0e+07 1.0e-08 5.3e+02 5.6e+02 0.22 (cid:6) (cid:5) Cahier du GERAD G-2021-02
Table 6.
Tax4D problem. NCL with KNITRO solving subproblems. (cid:7) (cid:4) julia> pTax4D = AmplModel("data/pTax4D")Minimization problem data/pTax4Dnvar = 432, ncon = 46441 (1 linear)julia> NCLSolve(pTax4D, outlev=0)outer inner NCL obj (cid:107) r (cid:107) η (cid:107)∇ L (cid:107) ω ρ µ init (cid:107) y (cid:107) (cid:107) x (cid:107) time1 12 -1.34e+04 3.3e-02 1.0e-02 5.4e-03 1.0e-02 1.0e+02 1.0e-01 1.0e+00 7.2e+02 3.382 12 -1.31e+04 1.3e-02 1.0e-02 4.0e-03 1.0e-02 1.0e+03 1.0e-03 1.0e+00 7.2e+02 3.233 15 -1.30e+04 5.1e-03 1.0e-02 2.1e-04 1.0e-02 1.0e+04 1.0e-03 1.0e+00 7.1e+02 3.864 31 -1.30e+04 3.2e-03 1.0e-03 1.3e-05 1.0e-03 1.0e+04 1.0e-05 5.2e+01 7.0e+02 7.955 37 -1.30e+04 1.8e-03 1.0e-03 1.2e-05 1.0e-03 1.0e+05 1.0e-05 5.2e+01 7.0e+02 9.896 44 -1.29e+04 5.0e-04 1.0e-03 1.1e-06 1.0e-03 1.0e+06 1.0e-06 5.2e+01 7.0e+02 11.937 16 -1.29e+04 2.6e-04 1.0e-04 1.2e-05 1.0e-04 1.0e+06 1.0e-06 5.3e+02 7.0e+02 3.748 30 -1.29e+04 4.4e-05 1.0e-04 1.2e-07 1.0e-04 1.0e+07 1.0e-07 5.3e+02 7.0e+02 8.159 9 -1.29e+04 2.3e-05 1.0e-05 1.2e-07 1.0e-05 1.0e+07 1.0e-07 8.2e+02 7.0e+02 2.4910 11 -1.29e+04 3.8e-06 1.0e-05 1.0e-08 1.0e-05 1.0e+08 1.0e-08 8.2e+02 7.0e+02 3.0911 6 -1.29e+04 1.7e-07 1.0e-06 1.3e-08 1.0e-06 1.0e+08 1.0e-08 9.4e+02 7.0e+02 1.74 (cid:6) (cid:5) Table 7.
Tax5D problem. NCL with KNITRO solving subproblems. (cid:7) (cid:4) julia> pTax5D = AmplModel("data/pTax5D")Minimization problem data/pTax5Dnvar = 864, ncon = 186193 (1 linear)julia> NCLSolve(pTax5D, outlev=0)outer inner NCL obj (cid:107) r (cid:107) η (cid:107)∇ L (cid:107) ω ρ µ init (cid:107) y (cid:107) (cid:107) x (cid:107) time1 64 -1.76e+05 2.0e-01 1.0e-02 2.3e-03 1.0e-02 1.0e+02 1.0e-01 1.0e+00 1.1e+04 80.432 29 -1.74e+05 4.9e-02 1.0e-02 1.2e-03 1.0e-02 1.0e+03 1.0e-03 1.0e+00 1.1e+04 35.023 23 -1.74e+05 1.6e-02 1.0e-02 1.0e-03 1.0e-02 1.0e+04 1.0e-03 1.0e+00 1.1e+04 28.964 46 -1.74e+05 4.1e-03 1.0e-02 3.6e-05 1.0e-02 1.0e+05 1.0e-05 1.0e+00 1.1e+04 54.505 41 -1.74e+05 2.8e-03 1.0e-03 1.7e-05 1.0e-03 1.0e+05 1.0e-05 4.1e+02 1.1e+04 52.726 28 -1.74e+05 6.1e-04 1.0e-03 1.0e-06 1.0e-03 1.0e+06 1.0e-06 4.1e+02 1.1e+04 34.387 13 -1.74e+05 2.1e-04 1.0e-04 1.4e-06 1.0e-04 1.0e+06 1.0e-06 1.0e+03 1.1e+04 14.818 12 -1.74e+05 5.3e-05 1.0e-04 1.2e-07 1.0e-04 1.0e+07 1.0e-07 1.0e+03 1.1e+04 14.809 7 -1.74e+05 4.5e-06 1.0e-05 1.0e-07 1.0e-05 1.0e+07 1.0e-07 1.0e+03 1.1e+04 9.4910 5 -1.74e+05 8.0e-07 1.0e-06 1.2e-08 1.0e-06 1.0e+07 1.0e-08 1.0e+03 1.1e+04 7.02 (cid:6) (cid:5) Our Julia module CUTEst.jl [21] provides an interface with the
CUTEst [9]environment and problem collection. Its main feature is to let users instanti-ate problems from
CUTEst using the
CUTEstModel constructor so they can bemanipulated transparently or passed to a solver like any other
NLPModel .On a set of constrained problems with at least variables whoseconstraints are all nonlinear,
KNITRO solves and
NCL solves . Although oursimple implementation of
NCL is not competitive with plain
KNITRO in general,it does solve a few problems on which
KNITRO fails. Those are summarized inTables 8 and 9. The above results suggest that
NCL ’s strength might reside insolving difficult problems (rather than being the fastest), and that more researchis needed to improve its efficiency.
Cahier du GERAD G-2021-02
Julia Implementation of Algorithm NCL . Table 8.
KNITRO results on CUTEst constrained problems (a subset that failed). name nvar ncon f (cid:107)∇ L (cid:107) (cid:107) c (cid:107) t iter f ∇ f c ∇ c ∇ L statusCATENARY − . e +10 4 . e +00 2 . e +09 18 .
50 2000 7835 2002 7835 2002 2000 max_iterCOSHFUN − . e +17 5 . e −
01 0 . e +00 19 .
40 2000 2001 2001 2001 2001 2000 max_iterDRCAVTY1 . e +00 0 . e +00 2 . e −
03 456 .
00 2000 9191 2002 9191 2002 2000 max_iterEG3 . e +05 2 . e +03 3 . e −
01 7 .
18 51 55 52 55 52 52 infeasibleJUNKTURN . e −
03 1 . e −
02 6 . e −
07 123 .
00 1913 15051 1915 15051 1915 1914 unknownLUKVLE11 . e +04 5 . e +02 4 . e −
01 86 .
50 2000 6945 2001 6945 2001 2000 max_iterLUKVLE17 . e +04 1 . e −
02 9 . e −
07 47 .
00 2000 3190 2001 3190 2001 2000 max_iterLUKVLE18 . e +04 4 . e +01 2 . e −
08 83 .
90 2000 4190 2001 4190 2001 2000 max_iterORTHRDS2 . e +02 3 . e −
01 5 . e −
13 0 .
70 42 92 43 92 43 43 unknown
Table 9.
NCL results on the same problems (all successful). name nvar ncon f (cid:107)∇ L (cid:107) (cid:107) c (cid:107) t iter f ∇ f c ∇ c ∇ L statusCATENARY − . e +06 1 . e −
09 1 . e −
07 1 .
76 183 430 195 430 206 194 first_orderCOSHFUN − . e −
01 9 . e −
07 1 . e −
07 6 .
46 328 1712 337 1712 346 337 first_orderDRCAVTY1 . e +00 1 . e −
08 1 . e −
06 29 .
30 222 344 233 344 243 232 first_orderEG3 . e −
07 1 . e −
08 1 . e −
07 4 .
94 37 47 47 47 57 47 first_orderJUNKTURN . e −
06 6 . e −
07 1 . e −
07 5 .
29 108 131 119 131 129 118 first_orderLUKVLE11 . e +02 4 . e −
08 1 . e −
06 1 .
86 37 63 47 63 57 47 first_orderLUKVLE17 . e +04 1 . e −
08 1 . e −
07 2 .
45 60 94 78 94 96 78 first_orderLUKVLE18 . e +04 2 . e −
12 1 . e −
09 3 .
43 60 81 79 81 98 79 first_orderORTHRDS2 . e +02 9 . e −
08 1 . e −
07 1 .
23 47 63 62 63 77 62 first_order
An important class of problems worthy of special attention is nonlinear least-squares (NLS) problems of the form min x (cid:107) c ( x ) (cid:107) subject to (cid:96) ≤ x ≤ u, (2)where the Jacobian of c ( x ) is again J ( x ) , and the bounds are often empty. Suchproblems are not immediately meaningful to Algorithm NCL, but if they arepresented in the (probably infeasible) form min x subject to c ( x ) = 0 , (cid:96) ≤ x ≤ u, (3)the first NCL subproblem will beNC minimize x, r y T r + ρ (cid:107) r (cid:107) subject to c ( x ) + r = 0 , (cid:96) ≤ x ≤ u, which is well suited to KNITRO and is equivalent to (2) if y = 0 and ρ > .If we treat NLS problems as a special case, we can set y = 0 , ρ = 1 , η = η ∗ , ω = ω ∗ and obtain an optimal solution in one NCL iteration. In this sense,Algorithm
NCL is ideally suited to NLS problems (2).The
CUTEst collection features a number of NLS problems in both forms (2)and (3). While formulation (2) allows evaluation of the objective gradient J ( x ) T c ( x ) , it does not give access to J ( x ) itself. In contrast, a problem modeledas (3) allows solvers to access J ( x ) directly. Cahier du GERAD G-2021-02
The NLPModels modeling package allows us to formulate (2) from a problemgiven as (3) and fulfill requests for J ( x ) in (2) by returning the constraint Jacobianof (3). Alternatively, problem NC is easily created by the NCLModel constructor.The construction of both models is illustrated in Listing 1.1. Once a problemin the form (2) has been simulated in this way, it can be passed to
KNITRO ’snonlinear least-squares solver, which is a variant of the Levenberg-Marquardtmethod in which bound constraints are treated via an interior-point method.
Listing 1.1.
Formulating (2) from (3) (cid:7) (cid:4) julia> using CUTEstjulia> model = CUTEstModel("ARWHDNE") (3) julia> nls_model = FeasibilityResidual(model) (3) as representing (2) julia> knitro(nls_model) ρ =1.0) julia> knitro(ncl_model) (cid:6) (cid:5) We identified problems in the form (3) in
CUTEst . We solve each problemin two ways:Solver knitro_nls applies
KNITRO ’s nonlinear least-squares method to (2).Solver ncl_nls uses
KNITRO to perform a single
NCL iteration on NC .In both cases, KNITRO is given a maximum of iterations and minutes ofCPU time. Optimality and feasibility tolerances are set to − . knitro_nls solved problems to optimality, reached the iteration limit in cases and the time limit in cases, and failed for another reason in cases. ncl_nls solved problems to optimality, reached the iteration limit in casesand the time limit in cases, and failed for another reason in cases.Figure 1 shows Dolan-Moré performance profiles comparing the two solvers.The top and middle plots use the number of residual and residual Jacobianevaluations as metric, which, in the case of (3), corresponds to the number ofconstraint and constraint Jacobian evaluations. The bottom plot uses time asmetric. ncl_nls outperforms knitro_nls in all three measures and appearssubstantially more robust. It is important to keep in mind that a key differencebetween the two algorithms is that ncl_nls uses second-order information, andtherefore performs Hessian evaluations. Nevertheless, those evaluations are notso costly as to put NCL at a disadvantage in terms of run-time. For reference,Tables 10 and 11 in Appendix A give the detailed results.
Our
AMPL implementation of the tax policy models and Algorithm
NCL hasbeen the only way we could handle these particular problems reliably [15], with
KNITRO solving each subproblem accurately. Our Julia implementation of
NCL achieves greater efficiency on these
AMPL models by gradually tightening the
KNITRO feasibility and optimality tolerances. It also permits testing on a broadrange of problems, as illustrated on nonlinear least-squares problems and other
Cahier du GERAD G-2021-02
Julia Implementation of Algorithm NCL Fig. 1.
Performance profiles comparing knitro_nls (KNITRO’s NLS solver appliedto (2)) and ncl_nls (KNITRO solving NC ) on 127 nonlinear least squares problemsfrom CUTEst. ncl_nls is more efficient. Cahier du GERAD G-2021-02 problems from the
CUTEst test set. We believe Algorithm
NCL could become aneffective general-purpose optimization solver when first and second derivativesare available. It is especially useful when the LICQ is not satisfied at the solution.The current Julia implementation of
NCL (with
KNITRO as subproblem solver)is not quite competitive with
KNITRO itself on the general
CUTEst problems interms of run-time or number of evaluations, but it does solve some problems onwhich
KNITRO fails. An advantage is that the implementation is generic andmay be applied to problems from any collection adhering to the interface of theNLPModels.jl package [24].
Acknowledgments
We are deeply grateful to Professor Ken Judd and Dr Che-Lin Su for developingthe
AMPL tax policy model [11] that led to the development of Algorithm NCL[15], and to the developers of
AMPL , Julia,
IPOPT , and
KNITRO for makingthe implementation and evaluation of Algorithm NCL possible. In particular, wethank Dr Richard Waltz of Artelys for his help in finding runtime options forwarm-starting
KNITRO . We also give sincere thanks to Pierre-Élie Personnazfor obtaining the Julia/NCL results in Tables 3–7, and to Professor MehiddinAl-Baali and other organizers of the NAO-V conference
Numerical Analysis andOptimization at Sultan Qaboos University, Muscat, Oman, which brought theauthors together in January 2020. Finally, we are very grateful to the referees fortheir constructive questions and comments, and to Michael Friedlander for hishelpful discussions.
References
Commun. ACM , 29(8):711––721,August 1986.3. Jeff Bezanson, Alan Edelman, Stefan Karpinski, and Viral B Shah. Julia: A freshapproach to numerical computing.
SIAM Rev. , 59(1):65–98, 2017.4. A. R. Conn, N. I. M. Gould, and Ph. L. Toint. A globally convergent augmentedLagrangian algorithm for optimization with general constraints and simple bounds.
SIAM J. Numer. Anal. , 28:545–572, 1991.5. A. R. Conn, N. I. M. Gould, and Ph. L. Toint.
LANCELOT: A Fortran Packagefor Large-scale Nonlinear Optimization (Release A) . Lecture Notes in ComputationMathematics 17. Springer Verlag, Berlin, Heidelberg, New York, London, Paris andTokyo, 1992.6. R. Fourer, D. M. Gay, and B. W. Kernighan.
AMPL: A Modeling Language forMathematical Programming
SIAM Rev. , 47(1):99–131, 2005. SIGESTarticle.
Cahier du GERAD G-2021-02
Julia Implementation of Algorithm NCL
9. N. I. M. Gould, D. Orban, and Ph. L. Toint. CUTEst: a Constrained and Uncon-strained Testing Environment with safe threads.
Comput. Optim. Appl.
INFORMSJ. Comput. , 27(2), 2015.15. D. Ma, K. L. Judd, D. Orban, and M. A. Saunders. Stabilized optimization via anNCL algorithm. In M. Al-Baali et al., editor,
Numerical Analysis and Optimization,NAO-IV, Muscat, Oman, January 2017 , pages 173–191. Springer InternationalPublishing AG, 2018.16. B. A. Murtagh and M. A. Saunders. A projected Lagrangian algorithm and itsimplementation for sparse nonlinear constraints.
Math. Program. Study
A Detailed results for Julia/NCL on NLS problems
Table 10 reports the detailed results of KNITRO/Levenberg-Marquardt onproblems of the form (2) using the modeling mechanism of Section 5. In the tableheaders, “nvar” is the number of variables, “ncon” is the number of constraints(i.e., the number of least-squares residuals), f is the final objective value, (cid:107)∇ L (cid:107) is the final dual residual, t is the run-time in seconds, “iter” is the number ofiterations, “ c ” is the number of constraint (i.e, residual) evaluations, “ ∇ c ” is Cahier du GERAD G-2021-02 the number of constraint (i.e., residual) Jacobian evaluations, and “status” is thefinal solver status.Table 11 reports the results of Julia/NCL solving Problem NC for thesame models. In the interest of space, the second table does not repeat problemdimensions. The other columns are as follows: (cid:107) c (cid:107) is the final primal feasibility,and ∇ L is the number of Hessian evaluations. Table 10. knitro_nls results on
CUTEst nonlinear least-squares problems. problems were solved successfully. name nvar ncon f (cid:107)∇ L (cid:107) t iter c ∇ c statusARWHDNE
500 998 6 . e +01 8 . e −
06 0 .
46 22 167 23 first_orderBA-L1
57 12 1 . e −
25 3 . e −
10 0 .
00 5 6 6 first_orderBA-L16 . e +05 4 . e −
04 69 .
92 27 30 28 first_orderBA-L1SP
57 12 1 . e −
23 5 . e −
09 0 .
01 5 6 6 first_orderBA-L21 . e +05 7 . e −
05 129 .
28 119 160 120 first_orderBA-L49 . e +05 4 . e +05 1809 .
44 121 728 122 max_timeBA-L52 . e +06 3 . e +07 1808 .
91 37 200 38 max_timeBA-L73 . e +05 9 . e +06 1801 .
24 396 2312 397 max_timeBARDNE . e −
03 2 . e −
09 0 .
00 5 6 6 first_orderBDQRTICNE . e +04 1 . e −
01 91 .
67 500 3668 501 max_iterBEALENE . e −
25 2 . e −
12 0 .
00 6 8 7 first_orderBIGGS6NE . e −
17 5 . e −
09 0 .
01 30 124 31 first_orderBOX3NE . e −
19 3 . e −
10 0 .
00 5 6 6 first_orderBROWNBSNE . e +00 0 . e +00 0 .
00 12 53 13 first_orderBROWNDENE . e +04 2 . e +00 0 .
09 500 1697 501 max_iterBRYBNDNE . e −
21 8 . e −
11 1 .
16 6 7 7 first_orderCHAINWOONE . e +03 5 . e +01 63 .
96 500 1339 501 max_iterCHEBYQADNE
100 100 4 . e −
03 1 . e −
04 7 .
53 500 1994 501 max_iterCHNRSBNE
50 98 7 . e −
18 3 . e −
08 0 .
01 38 78 39 first_orderCHNRSNBMNE
50 98 1 . e −
20 1 . e −
09 0 .
02 54 132 55 first_orderCOATINGNE
134 252 2 . e −
01 5 . e −
07 0 .
01 9 11 10 first_orderCUBENE . e −
26 1 . e −
13 0 .
00 4 9 5 first_orderDECONVBNE
63 40 4 . e −
10 6 . e −
07 0 .
03 54 214 55 first_orderDECONVNE
63 40 2 . e −
16 8 . e −
10 0 .
00 2 3 3 first_orderDENSCHNBNE . e −
27 1 . e −
13 0 .
00 5 6 6 first_orderDENSCHNCNE . e −
22 8 . e −
11 0 .
00 7 8 8 first_orderDENSCHNDNE . e −
10 5 . e −
07 0 .
00 17 18 18 first_orderDENSCHNENE . e −
22 2 . e −
11 0 .
00 8 19 9 first_orderDENSCHNFNE . e −
26 2 . e −
12 0 .
00 5 6 6 first_orderDEVGLA1NE . e −
13 2 . e −
08 0 .
00 11 31 12 first_orderDEVGLA2NE . e −
15 5 . e −
07 0 .
00 9 16 10 first_orderEGGCRATENE . e +00 9 . e −
07 0 .
00 5 6 6 first_orderELATVIDUNE . e +01 3 . e −
06 0 .
00 15 30 16 first_orderENGVAL2NE . e −
32 2 . e −
16 0 .
00 10 14 11 first_orderERRINROSNE
50 98 2 . e +01 1 . e −
05 0 .
01 45 63 46 first_orderERRINRSMNE
50 98 1 . e +01 1 . e −
05 0 .
01 44 66 45 first_orderEXP2NE . e −
19 8 . e −
12 0 .
00 5 6 6 first_orderEXPFITNE . e −
01 4 . e −
07 0 .
00 10 12 11 first_orderEXTROSNBNE . e −
31 1 . e −
14 0 .
07 6 14 7 first_orderFBRAIN2NE . e −
01 6 . e −
07 0 .
08 8 11 9 first_orderFBRAINNE . e −
01 5 . e −
07 0 .
02 5 6 6 first_orderFREURONE . e +01 1 . e −
05 0 .
01 19 102 20 first_orderGENROSEBNE
500 998 7 . e +02 3 . e −
06 0 .
04 8 10 10 first_orderGENROSENE . e +02 6 . e +00 3 .
40 500 1567 501 max_iterGULFNE . e −
01 4 . e +06 0 .
27 500 2578 501 max_iterHATFLDANE . e −
17 3 . e −
09 0 .
00 8 10 9 first_orderHATFLDBNE . e −
03 1 . e −
07 0 .
00 7 8 8 first_orderHATFLDCNE
25 25 1 . e −
15 1 . e −
08 0 .
00 3 4 4 first_orderHATFLDDNE . e −
07 7 . e −
11 0 .
00 6 9 7 first_orderHATFLDENE . e −
06 8 . e −
12 0 .
00 5 6 6 first_orderHATFLDFLNE . e −
05 9 . e −
07 0 .
01 11 69 12 first_orderHELIXNE . e −
20 2 . e −
09 0 .
00 9 11 10 first_order
Continued on next page
Cahier du GERAD G-2021-02
Julia Implementation of Algorithm NCL Table 10. knitro_nls results on
CUTEst nonlinear least-squares problems. problems were solved successfully. name nvar ncon f (cid:107)∇ L (cid:107) t iter c ∇ c statusHIMMELBFNE . e +06 4 . e −
04 0 .
01 28 66 29 first_orderHS1NE . e −
17 4 . e −
09 0 .
00 9 20 10 first_orderHS25NE . e +01 9 . e −
09 0 .
00 0 1 1 first_orderHS2NE . e +00 2 . e −
08 0 .
00 6 8 8 first_orderINTEQNE
12 12 6 . e −
14 2 . e −
07 0 .
00 2 3 3 first_orderJENSMPNE . e +01 7 . e −
06 0 .
01 10 82 11 first_orderJUDGENE . e +00 1 . e −
06 0 .
00 9 10 10 first_orderKOEBHELBNE . e +01 7 . e −
07 0 .
10 194 757 195 first_orderKOWOSBNE . e −
04 6 . e −
07 0 .
00 16 30 17 first_orderLIARWHDNE . e −
27 4 . e −
12 1 .
07 5 6 6 first_orderLINVERSENE . e +02 1 . e −
02 13 .
28 500 1753 502 max_iterMANCINONE
100 100 7 . e −
22 2 . e −
08 0 .
10 5 6 6 first_orderMANNE . e +39 1 . e +19 95 .
71 167 168 168 unknownMARINE . e −
21 8 . e −
08 8 .
51 9 10 10 first_orderMEYER3NE . e +01 1 . e −
03 0 .
00 11 19 12 unknownMODBEALENE . e +00 8 . e −
07 109 .
80 38 73 39 first_orderMOREBVNE
10 10 1 . e −
14 1 . e −
08 0 .
00 2 3 3 first_orderMUONSINE . e +04 8 . e −
04 0 .
02 36 37 37 first_orderNGONE
200 5048 7 . e −
13 8 . e −
07 36 .
90 221 914 223 first_orderNONDIANE . e −
01 1 . e −
08 2 .
63 16 39 17 first_orderNONMSQRTNE . e +02 1 . e +01 154 .
67 500 2624 501 max_iterNONSCOMPNE . e −
06 5 . e −
07 11 .
31 80 217 81 first_orderOSCIGRNE . e −
24 8 . e −
09 3 .
62 7 8 8 first_orderOSCIPANE
10 10 5 . e −
01 1 . e −
01 0 .
11 500 2880 501 max_iterPALMER1ANE . e −
02 6 . e −
07 0 .
01 25 84 26 first_orderPALMER1BNE . e +00 2 . e −
08 0 .
00 6 7 7 first_orderPALMER1ENE . e −
01 1 . e +02 0 .
17 500 3471 501 max_iterPALMER1NE . e +03 8 . e +00 0 .
13 500 2561 501 max_iterPALMER2ANE . e −
03 5 . e −
07 0 .
02 58 248 59 first_orderPALMER2BNE . e −
01 5 . e −
08 0 .
00 8 10 9 first_orderPALMER2ENE . e −
02 8 . e −
07 0 .
01 9 68 10 first_orderPALMER2NE . e +03 4 . e −
06 0 .
01 26 92 27 unknownPALMER3ANE . e −
02 1 . e −
07 0 .
00 18 39 19 first_orderPALMER3BNE . e +00 1 . e −
07 0 .
00 11 14 12 first_orderPALMER3ENE . e −
02 1 . e +01 0 .
15 500 3348 501 max_iterPALMER3NE . e +03 5 . e −
02 0 .
15 500 3303 501 max_iterPALMER4ANE . e −
02 3 . e −
07 0 .
03 98 445 99 first_orderPALMER4BNE . e +00 6 . e −
07 0 .
00 14 18 15 first_orderPALMER4ENE . e −
02 9 . e −
07 0 .
04 69 667 70 first_orderPALMER4NE . e +03 2 . e −
05 0 .
01 35 163 36 unknownPALMER5ANE . e −
01 1 . e +00 0 .
13 500 3000 501 max_iterPALMER5BNE . e −
03 3 . e −
07 0 .
03 93 481 94 first_orderPALMER5ENE . e −
02 1 . e +00 0 .
16 500 3508 501 max_iterPALMER6ANE . e −
02 2 . e −
07 0 .
00 17 30 18 first_orderPALMER6ENE . e −
02 5 . e −
07 0 .
08 142 1654 143 first_orderPALMER7ANE . e +00 2 . e +00 0 .
13 500 2806 501 max_iterPALMER7ENE . e +00 1 . e +02 0 .
16 500 3243 500 max_iterPALMER8ANE . e −
02 5 . e −
07 0 .
02 72 303 73 first_orderPALMER8ENE . e −
01 2 . e −
04 0 .
22 500 5007 501 max_iterPENLT1NE
10 11 3 . e −
10 8 . e −
07 0 .
02 111 385 112 first_orderPENLT2NE . e −
11 9 . e −
07 0 .
03 163 613 164 first_orderPINENE . e −
17 1 . e −
07 1 .
72 2 10 3 first_orderPOWERSUMNE . e −
17 7 . e −
07 0 .
00 17 21 18 first_orderPRICE3NE . e −
22 3 . e −
10 0 .
00 7 8 8 first_orderPRICE4NE . e −
14 9 . e −
07 0 .
00 21 22 22 first_orderQINGNE
100 100 9 . e −
20 1 . e −
09 0 .
00 5 8 6 first_orderRSNBRNE . e −
30 3 . e −
15 0 .
00 11 34 12 first_orderS308NE . e −
01 9 . e −
07 0 .
00 34 36 35 first_orderSBRYBNDNE . e −
21 3 . e −
07 1 .
16 6 7 7 first_orderSINVALNE . e −
32 2 . e −
15 0 .
00 4 11 5 first_orderSPECANNE . e −
13 5 . e −
08 0 .
08 6 7 7 first_orderSROSENBRNE . e −
28 6 . e −
16 0 .
49 2 3 3 first_order
Continued on next page
Cahier du GERAD G-2021-02
Table 10. knitro_nls results on
CUTEst nonlinear least-squares problems. problems were solved successfully. name nvar ncon f (cid:107)∇ L (cid:107) t iter c ∇ c statusSSBRYBNDNE . e −
22 4 . e −
10 1 .
15 6 7 7 first_orderSTREGNE . e −
03 4 . e −
01 0 .
05 500 535 501 max_iterSTRTCHDVNE
10 9 3 . e −
10 4 . e −
07 0 .
00 8 9 9 first_orderTQUARTICNE . e −
26 2 . e −
13 0 .
39 1 2 2 first_orderTRIGON1NE
10 10 6 . e −
16 1 . e −
07 0 .
00 4 5 5 first_orderTRIGON2NE
10 31 1 . e +00 1 . e −
09 0 .
00 7 8 8 first_orderVARDIMNE
10 12 7 . e −
16 4 . e −
07 0 .
00 9 10 10 first_orderVIBRBEAMNE . e −
02 6 . e −
07 0 .
00 10 11 11 first_orderWATSONNE
12 31 1 . e −
15 1 . e −
13 0 .
00 4 5 5 first_orderWAYSEA1NE . e −
16 8 . e −
07 0 .
00 7 8 8 first_orderWAYSEA2NE . e −
18 1 . e −
08 0 .
00 11 17 12 first_orderWEEDSNE . e +00 3 . e −
07 0 .
00 15 25 16 first_orderWOODSNE . e −
01 0 . e +00 0 .
33 2 3 3 first_order
Table 11. ncl_nls results on
CUTEst nonlinear least-squaresproblems. 119 problems were solved successfully. name f (cid:107)∇ L (cid:107) (cid:107) c (cid:107) t iter c ∇ c ∇ L statusARWHDNE . e +01 1 . e −
14 1 . e −
16 2 .
65 30 105 32 31 first_orderBA-L1 . e −
31 9 . e −
16 2 . e −
10 0 .
01 5 6 7 6 first_orderBA-L16 . e +05 1 . e +03 4 . e −
07 937 .
15 88 454 90 89 unknownBA-L1SP . e −
23 7 . e −
12 8 . e −
05 0 .
01 4 5 6 5 first_orderBA-L21 . e +05 1 . e −
02 3 . e −
05 947 .
20 201 1093 203 203 max_timeBA-L49 . e +04 1 . e −
01 1 . e −
04 869 .
15 217 1189 219 219 max_timeBA-L52 . e +06 2 . e +02 8 . e −
01 1778 .
45 31 50 33 32 max_timeBA-L73 . e +05 1 . e +02 1 . e −
06 693 .
06 135 635 137 136 unknownBARDNE . e −
03 2 . e −
14 5 . e −
12 0 .
00 5 6 7 6 first_orderBDQRTICNE . e +04 3 . e −
10 1 . e −
11 0 .
47 17 25 19 18 first_orderBEALENE . e −
20 2 . e −
10 1 . e −
09 0 .
00 8 9 10 9 first_orderBIGGS6NE . e −
17 6 . e −
12 7 . e −
08 0 .
02 60 74 62 61 first_orderBOX3NE . e −
19 7 . e −
12 1 . e −
10 0 .
00 5 6 7 6 first_orderBROWNBSNE . e −
33 5 . e −
12 1 . e −
10 0 .
00 13 42 15 14 first_orderBROWNDENE . e +04 9 . e −
09 3 . e −
11 0 .
00 14 22 16 15 first_orderBRYBNDNE . e −
20 2 . e −
10 2 . e −
10 0 .
20 6 7 8 7 first_orderCHAINWOONE . e +03 1 . e +00 2 . e −
01 14 .
40 500 546 502 501 max_iterCHEBYQADNE . e −
03 5 . e −
07 2 . e −
07 0 .
60 33 51 35 34 first_orderCHNRSBNE . e −
21 1 . e −
10 7 . e −
11 0 .
02 34 51 36 35 first_orderCHNRSNBMNE . e −
23 4 . e −
11 3 . e −
10 0 .
02 38 57 40 39 first_orderCOATINGNE . e −
01 1 . e −
09 3 . e −
10 0 .
02 11 22 13 12 first_orderCUBENE . e −
33 9 . e −
17 6 . e −
15 0 .
00 2 7 4 3 first_orderDECONVBNE . e −
03 3 . e −
07 5 . e −
06 0 .
03 20 34 23 21 first_orderDECONVNE . e −
15 3 . e −
08 6 . e −
10 0 .
00 2 3 4 3 first_orderDENSCHNBNE . e −
17 2 . e −
09 5 . e −
09 0 .
00 6 8 8 7 first_orderDENSCHNCNE . e −
38 4 . e −
19 6 . e −
06 0 .
00 6 7 8 7 first_orderDENSCHNDNE . e −
28 4 . e −
17 3 . e −
04 0 .
00 15 16 17 16 first_orderDENSCHNENE . e −
21 9 . e −
11 1 . e −
06 0 .
00 15 20 17 16 first_orderDENSCHNFNE . e −
38 2 . e −
19 1 . e −
06 0 .
00 4 5 6 5 first_orderDEVGLA1NE . e −
13 2 . e −
08 4 . e −
11 0 .
01 16 45 18 17 first_orderDEVGLA2NE . e −
14 6 . e −
08 3 . e −
07 0 .
00 9 12 11 10 first_orderEGGCRATENE . e +00 3 . e −
08 2 . e −
09 0 .
00 4 5 6 5 first_orderELATVIDUNE . e +01 1 . e −
09 5 . e −
10 0 .
00 8 9 10 9 first_orderENGVAL2NE . e −
15 1 . e −
08 1 . e −
07 0 .
00 12 20 14 13 first_orderERRINROSNE . e +01 1 . e −
09 4 . e −
10 0 .
01 17 24 19 18 first_orderERRINRSMNE . e +01 1 . e −
08 1 . e −
09 0 .
02 25 38 27 26 first_orderEXP2NE . e −
19 7 . e −
13 1 . e −
13 0 .
00 5 6 7 6 first_orderEXPFITNE . e −
01 2 . e −
10 4 . e −
12 0 .
01 21 149 23 22 first_orderEXTROSNBNE − . e +00 1 . e −
06 1 . e −
09 1 .
16 250 1860 252 251 first_orderFBRAIN2NE . e −
01 6 . e −
10 3 . e −
11 0 .
11 5 6 7 6 first_orderFBRAINNE . e −
01 8 . e −
08 1 . e −
09 0 .
05 4 5 6 5 first_order
Continued on next page
Cahier du GERAD G-2021-02
Julia Implementation of Algorithm NCL Table 11. ncl_nls results on
CUTEst nonlinear least-squaresproblems. 119 problems were solved successfully. name f (cid:107)∇ L (cid:107) (cid:107) c (cid:107) t iter c ∇ c ∇ L statusFREURONE . e +01 1 . e −
10 5 . e −
12 0 .
00 7 15 9 8 first_orderGENROSEBNE . e +02 1 . e −
06 1 . e −
06 0 .
04 6 8 9 7 first_orderGENROSENE . e −
01 8 . e −
14 1 . e −
14 2 .
07 436 639 438 437 first_orderGULFNE . e −
20 3 . e −
11 1 . e −
10 0 .
02 20 28 22 21 first_orderHATFLDANE . e −
20 5 . e −
11 6 . e −
09 0 .
00 9 10 11 10 first_orderHATFLDBNE . e −
03 4 . e −
12 1 . e −
11 0 .
00 7 8 9 8 first_orderHATFLDCNE . e −
17 4 . e −
09 7 . e −
09 0 .
00 3 4 5 4 first_orderHATFLDDNE . e −
07 4 . e −
09 9 . e −
07 0 .
00 5 8 7 6 first_orderHATFLDENE . e −
06 5 . e −
11 1 . e −
09 0 .
00 5 6 7 6 first_orderHATFLDFLNE . e −
05 9 . e −
12 7 . e −
09 0 .
02 103 443 105 104 first_orderHELIXNE . e −
43 1 . e −
21 1 . e −
10 0 .
00 9 12 11 10 first_orderHIMMELBFNE . e +06 2 . e +02 2 . e −
04 0 .
09 500 720 502 501 max_iterHS1NE . e −
17 4 . e −
09 5 . e −
16 0 .
00 7 11 9 8 first_orderHS25NE . e +01 8 . e −
08 3 . e −
11 0 .
00 3 5 6 4 first_orderHS2NE . e +00 8 . e −
12 1 . e −
12 0 .
00 5 7 8 6 first_orderINTEQNE . e −
37 8 . e −
19 1 . e −
07 0 .
00 2 3 4 3 first_orderJENSMPNE . e +01 1 . e −
08 3 . e −
09 0 .
00 8 14 10 9 first_orderJUDGENE . e +00 4 . e −
07 2 . e −
07 0 .
00 5 6 7 6 first_orderKOEBHELBNE . e +01 1 . e −
08 4 . e −
09 0 .
13 119 164 121 121 first_orderKOWOSBNE . e −
04 2 . e −
10 2 . e −
09 0 .
00 11 14 13 12 first_orderLIARWHDNE . e −
21 4 . e −
11 4 . e −
11 0 .
17 6 7 8 7 first_orderLINVERSENE . e +02 3 . e −
07 1 . e −
07 0 .
31 19 22 22 20 first_orderMANCINONE . e −
32 1 . e −
16 2 . e −
09 0 .
12 4 5 6 5 first_orderMANNE − . e −
01 9 . e −
07 0 . e +00 10 .
27 344 346 347 345 first_orderMARINE . e +06 1 . e −
08 1 . e −
08 3 .
76 40 43 43 41 first_orderMEYER3NE . e +01 2 . e −
06 1 . e −
08 0 .
00 8 14 10 9 first_orderMODBEALENE . e −
18 1 . e −
09 2 . e −
09 1 .
80 11 12 13 12 first_orderMOREBVNE . e −
37 9 . e −
19 7 . e −
08 0 .
00 2 3 4 3 first_orderMUONSINE . e +04 2 . e −
11 6 . e −
14 0 .
02 10 15 12 11 first_orderNGONE − . e +00 4 . e −
09 3 . e −
11 2 .
46 106 108 109 107 first_orderNONDIANE . e −
01 1 . e −
16 1 . e −
07 0 .
16 6 7 8 7 first_orderNONMSQRTNE . e +02 1 . e −
07 2 . e −
06 4 .
01 19 30 21 20 first_orderNONSCOMPNE . e −
07 1 . e −
10 3 . e −
06 0 .
52 19 28 21 20 first_orderOSCIGRNE . e −
35 4 . e −
18 1 . e −
07 3 .
78 6 7 8 7 first_orderOSCIPANE . e −
01 9 . e −
06 5 . e −
04 0 .
16 500 2698 502 501 max_iterPALMER1ANE . e +01 2 . e −
06 6 . e −
08 0 .
01 12 16 14 13 first_orderPALMER1BNE . e +00 1 . e −
08 1 . e −
10 0 .
01 11 16 13 12 first_orderPALMER1ENE . e −
04 7 . e −
09 1 . e −
06 0 .
00 7 8 9 8 first_orderPALMER1NE . e +03 2 . e −
06 2 . e −
09 0 .
02 37 46 39 38 first_orderPALMER2ANE . e −
03 1 . e −
10 4 . e −
08 0 .
02 44 69 46 45 first_orderPALMER2BNE . e −
01 1 . e −
07 3 . e −
09 0 .
00 8 11 10 9 first_orderPALMER2ENE . e −
02 6 . e −
08 3 . e −
10 0 .
01 20 23 22 21 first_orderPALMER2NE . e +03 3 . e −
06 5 . e −
08 0 .
02 49 64 51 50 first_orderPALMER3ANE . e −
02 2 . e −
08 9 . e −
07 0 .
01 13 17 15 14 first_orderPALMER3BNE . e +00 6 . e −
12 1 . e −
14 0 .
01 27 38 29 28 first_orderPALMER3ENE . e −
05 4 . e −
11 1 . e −
09 0 .
01 14 19 16 15 first_orderPALMER3NE . e +03 1 . e −
05 2 . e −
08 0 .
00 6 10 8 7 first_orderPALMER4ANE . e −
02 2 . e −
08 1 . e −
07 0 .
01 11 20 13 12 first_orderPALMER4BNE . e +00 1 . e −
12 2 . e −
14 0 .
01 32 41 34 33 first_orderPALMER4ENE . e −
05 4 . e −
09 2 . e −
07 0 .
01 11 13 13 12 first_orderPALMER4NE . e +03 1 . e −
05 4 . e −
08 0 .
00 6 10 8 7 first_orderPALMER5ANE . e −
02 9 . e −
07 4 . e −
07 0 .
10 182 874 184 183 first_orderPALMER5BNE . e −
03 3 . e −
09 1 . e −
07 0 .
01 44 50 46 45 first_orderPALMER5ENE . e −
02 1 . e −
09 3 . e −
07 0 .
02 60 200 62 61 first_orderPALMER6ANE . e −
02 2 . e −
08 7 . e −
06 0 .
01 22 31 24 23 first_orderPALMER6ENE . e −
02 5 . e −
07 2 . e −
06 0 .
00 7 10 9 8 first_orderPALMER7ANE . e +00 1 . e −
07 2 . e −
05 0 .
16 355 1551 357 358 first_orderPALMER7ENE . e +00 1 . e −
06 4 . e −
10 0 .
03 54 274 56 55 first_orderPALMER8ANE . e −
02 9 . e −
08 3 . e −
07 0 .
01 19 31 21 20 first_orderPALMER8ENE . e −
01 6 . e −
08 3 . e −
09 0 .
01 25 35 27 26 first_orderPENLT1NE . e −
10 5 . e −
14 3 . e −
05 0 .
00 8 9 10 9 first_orderPENLT2NE . e −
11 1 . e −
18 5 . e −
08 0 .
00 5 10 7 6 first_order
Continued on next page
Cahier du GERAD G-2021-02
Table 11. ncl_nls results on
CUTEst nonlinear least-squaresproblems. 119 problems were solved successfully. name f (cid:107)∇ L (cid:107) (cid:107) c (cid:107) t iter c ∇ c ∇ L statusPINENE . e −
05 4 . e −
10 1 . e −
06 0 .
72 13 18 16 14 first_orderPOWERSUMNE . e −
22 1 . e −
11 6 . e −
10 0 .
01 45 48 47 46 first_orderPRICE3NE . e −
35 5 . e −
18 1 . e −
05 0 .
00 6 7 8 7 first_orderPRICE4NE . e −
32 3 . e −
19 1 . e −
04 0 .
00 13 14 15 14 first_orderQINGNE . e −
33 2 . e −
17 9 . e −
05 0 .
00 4 7 6 5 first_orderRSNBRNE . e −
32 3 . e −
16 3 . e −
14 0 .
00 1 3 3 2 first_orderS308NE . e −
01 1 . e −
09 9 . e −
10 0 .
00 10 16 12 11 first_orderSBRYBNDNE . e −
20 2 . e −
10 2 . e −
10 0 .
18 6 7 8 7 first_orderSINVALNE . e −
30 1 . e −
15 1 . e −
15 0 .
00 1 3 3 2 first_orderSPECANNE . e −
13 3 . e −
08 1 . e −
12 0 .
49 10 14 12 11 first_orderSROSENBRNE . e −
33 1 . e −
18 2 . e −
18 0 .
05 2 3 4 3 first_orderSSBRYBNDNE . e −
20 2 . e −
10 2 . e −
10 0 .
21 6 7 8 7 first_orderSTREGNE . e −
27 6 . e −
22 4 . e −
15 0 .
00 2 3 4 3 first_orderSTRTCHDVNE . e −
15 1 . e −
10 6 . e −
07 0 .
00 10 11 12 11 first_orderTQUARTICNE . e +00 0 . e +00 2 . e −
16 0 .
05 1 2 3 2 first_orderTRIGON1NE . e −
39 5 . e −
20 1 . e −
08 0 .
00 4 5 6 5 first_orderTRIGON2NE . e +00 2 . e −
09 1 . e −
09 0 .
00 11 12 13 12 first_orderVARDIMNE . e −
20 2 . e −
09 1 . e −
10 0 .
00 14 19 16 15 first_orderVIBRBEAMNE . e −
02 2 . e −
11 1 . e −
14 0 .
01 7 8 9 8 first_orderWATSONNE . e −
15 5 . e −
16 5 . e −
14 0 .
00 4 5 6 5 first_orderWAYSEA1NE . e +00 0 . e +00 2 . e −
08 0 .
00 7 8 9 8 first_orderWAYSEA2NE . e −
17 2 . e −
09 1 . e −
07 0 .
00 10 20 12 11 first_orderWEEDSNE . e +00 8 . e −
08 7 . e −
07 0 .
01 17 24 19 18 first_orderWOODSNE − . e +04 9 . e −
13 8 . e −
10 0 .
04 3 4 5 4 first_orderfirst_order