NLsolve with StaticArrays

Specific Domains Numerics
1

I am wondering if there is a way to make NLsolve work with StaticArrays, because my function is written to take SVector.

I thought if I use the keyword argument inplace=false as described in here, I would be able to use a function and Jacobian that do not require in-place overwriting of intermediate states, and therefore I could define them as F = f(x::SVector) and J = j(x::SVector). I tried this idea by modifying the example shown in NLsolve’s README as follows:

using NLsolve
using StaticArrays

function f(x::SVector{2})
    return @SVector [(x[1]+3)*(x[2]^3-7)+18, sin(x[2]*exp(x[1])-1)]
end

function j(x::SVector{2})
    u = exp(x[1])*cos(x[2]*exp(x[1])-1)

    return @SMatrix [x[2]^3-7 3*x[2]^2*(x[1]+3); x[2]*u u]
end

nlsolve(f, j, @SVector([0.1,1.2]), inplace=false)

However, NLsolve still attempts to overwrite the solution:

julia> nlsolve(f, j, @SVector([0.1,1.2]), inplace=false)
ERROR: setindex!(::SVector{2, Float64}, value, ::Int) is not defined.
Stacktrace:
  [1] error(s::String)
    @ Base ./error.jl:33
  [2] setindex!(a::SVector{2, Float64}, value::Float64, i::Int64)
    @ StaticArrays ~/.julia/packages/StaticArrays/OWJK7/src/indexing.jl:3
[...]

So, I tried to change SVector and SMatrix in the above example to MVector and MMatrix, so that the vectors and matrices are modifiable:

using NLsolve
using StaticArrays

function f(x::MVector{2})
    return @MVector [(x[1]+3)*(x[2]^3-7)+18, sin(x[2]*exp(x[1])-1)]
end

function j(x::MVector{2})
    u = exp(x[1])*cos(x[2]*exp(x[1])-1)

    return @MMatrix [x[2]^3-7 3*x[2]^2*(x[1]+3); x[2]*u u]
end

nlsolve(f, j, @MVector([0.1,1.2]), inplace=false)

This change eliminates the errors, but generates a wrong solution: the solution doesn’t change from the initial guess [0.1, 1.2]:

julia> nlsolve(f, j, @MVector([0.1,1.2]), inplace=false)
Results of Nonlinear Solver Algorithm
 * Algorithm: Trust-region with dogleg and autoscaling
 * Starting Point: [0.1, 1.2]
 * Zero: [0.1, 1.2]
 * Inf-norm of residuals: 1.656800
 * Iterations: 1000
 * Convergence: false
   * |x - x'| < 0.0e+00: false
   * |f(x)| < 1.0e-08: false
 * Function Calls (f): 2
 * Jacobian Calls (df/dx): 1

Any suggestions?

2

Here you are creating a new vector, so you lost the track of the original one. Probably this should be:

function f(x::MVector{2})
    x .= @MVector [(x[1]+3)*(x[2]^3-7)+18, sin(x[2]*exp(x[1])-1)]
    return x # redundant here
end
3

Thanks for the suggestion. It makes my second example to find the correct solution. But isn’t this effectively making f in-place? Such a change doesn’t seem to serve the purpose of using the keyword argument inplace=false.

4

Yes, it is. No, it does not solve your original problem. I would bet (without really knowing) that inplace=false, if working, is only creating a copy of the original array to start with. You can do that yourself in a simple interface that copies the SVector to an MVector if the goal is to preserve the initial point.