NLsolve stops quickly

Hey everyone,

I have a multivariate non-linear root-finding problem for which I would like to use NLsolve.jl. It does not seem to work, however:

using NLsolve
function f!(F, β)
    F = mbagg(β, par) - mc.(β, par = par)
nlsolve(f!, rand(par.n))

The output that I receive is

julia> nlsolve(f!, rand(par.n))
Results of Nonlinear Solver Algorithm
 * Algorithm: Trust-region with dogleg and autoscaling
 * Starting Point: [0.6188824549955863, 0.5222000432685836, 0.6009845689936288, 0.5690479861423647, 0.9369342312698823, 0.2516124762445371, 0.362404246259066]
 * Zero: [NaN, NaN, NaN, NaN, NaN, NaN, NaN]
 * Inf-norm of residuals: 0.936934
 * Iterations: 1
 * Convergence: false
   * |x - x'| < 0.0e+00: false
   * |f(x)| < 1.0e-08: false
 * Function Calls (f): 2
 * Jacobian Calls (df/dx): 1

Which is not really informative about what is happening. Can someone help me and point me to the right direction?

Sorry for not being able to provide a MWE, but there is quite a bit in mbagg.

Maybe it has something to do that my function is only defined over positive values for beta. What can I do to solve such a constrained problem?

I am now trying the following:

f(β,par) = sum((mbagg(β, par) - mc.(β, par = par)).^2)
g(β) = f(β, par)

# Let's see if this works instead
lower = [1e-6, 1e-6, 1e-6, 1e-6, 1e-6, 1e-6, 1e-6]
upper = [Inf, Inf, Inf, Inf, Inf, Inf, Inf]
optimize(g, lower, upper, β_n)

But I receive the following error:

julia> optimize(g, lower, upper, β_n)
ERROR: MethodError: objects of type Array{Float64,1} are not callable
Use square brackets [] for indexing an Array.
 [1] (::NLSolversBase.var"#fg!#8"{typeof(g),Array{Float64,1}})(::Array{Float64,1}, ::Array{Float64,1}) at C:\Users\Ilja\.julia\packages\NLSolversBase\5oIMo\src\objective_types\abstract.jl:13.jl:13
 [2] value_gradient!!(::OnceDifferentiable{Float64,Array{Float64,1},Array{Float64,1}}, ::Array{Float64,1}) at C:\Users\Ilja\.julia\packages\NLSolversBase\5oIMo\src\interface.jl:82     at64,1},Array{Float64,1}}, ::Array{Float64,1}
 [3] initial_state(::BFGS{InitialStatic{Float64},HagerZhang{Float64,Base.RefValue{Bool}},Nothing,Nothing,Flat}, ::Optim.Options{Float64,Nothing}, ::OnceDifferentiable{Float64,Array{Float64,1},Array{Float64,1}}, ::Array{Float64,1}) at C:\Users\Ilja\.julia\packages\Optim\TNmSw\src\multivariate\solvers\first_order\bfgs.jl:66
 [4] optimize at C:\Users\Ilja\.julia\packages\Optim\TNmSw\src\multivariate\optimize\optimize.jl:33 [inlined]
 [5] #optimize#92 at C:\Users\Ilja\.julia\packages\Optim\TNmSw\src\multivariate\optimize\interface.jl:136 [inlined]                                                                     im.Options{Float64,Nothing}) at C:\Users\Ilja
 [6] optimize(::Function, ::Array{Float64,1}, ::Array{Float64,1}, ::Array{Float64,1}, ::BFGS{InitialStatic{Float64},HagerZhang{Float64,Base.RefValue{Bool}},Nothing,Nothing,Flat}, ::Optim.Options{Float64,Nothing}) at C:\Users\Ilja\.julia\packages\Optim\TNmSw\src\multivariate\optimize\interface.jl:134 (repeats 2 times)
 [7] top-level scope at REPL[55]:1
 [8] include_string(::Function, ::Module, ::String, ::String) at .\loading.jl:1088

I do not really understand how I can fix the error.

As far as I know, NLsolve has no option to bound parameters. So, you would have to transform your equation to avoid NaN results. For example, substitute beta by c.^2, solve for c, and calculate beta=sqrt.(abs.(c))

Optim takes a scalar function, you could try g(β) = sum( f(β, par).^2 )

1 Like

So I tried Optim.jl, it takes first long to run and then gives me the error seen above. Which is strange, since g(β) works just fine.

The part with beta = c.^2 and NLsolve sounds interesting, I’ll give it a try.

Its hard to say without MWE: I get this error when the function returns an array, instead of a scalar. Maybe try restarting the REPL?
You could insert some print statements to track down which parameters cause the problem


I guess optimize thinks I am passing a Jacobian and Hessian and, therefore, tries to do something like upper(beta). In that case I might be calling optimize wrongly. Let’s see if I can find the mistake. Especially since optimize(g, β_n) seems to work now. In the end I am going with

f(β,par) = sum((mbagg(β, par) - mc.(β, par = par)).^2)
g(β)     = f(β, par)
res      = optimize(g, β_n)

Anyways, does anybody know how to solve this problem with bounds on beta, i.e., call optimize correctly?

Final update, this here seems to work, but it is really slow compared to the above. The above, however, errors out at the moment when it tries to use negative values. One problem was that bounds are not allowed to be integers, it seems

lower = 1e-4.*ones(par.n)
upper = 10.0.*ones(par.n)

optimize(g, lower, upper, rand(par.n), Fminbox(NelderMead()))

In the end I’ll go with the substitution method you describe, but you have to square the solution in the end instead of taking the root :slight_smile:

f(c,par) = sum((mbagg(c.^2, par) - mc.(c.^2, par = par)).^2)
g(c)     = f(c, par)
(optimize(g, rand(par.n)).minimizer).^2
1 Like