[Avoid cfunction closure] GSL root finding problem with function passed as argument

General Usage
1

I have a usage scenario in which the functions that represent the equations for the multidimentional root-finding problem are generated in run time. I use GSL.jl to solve the root-finding problem, which requires wrapping the functions for the C library.

What I used to do is to make use of the @cfunction runtime enclosure. But, now I realize that this is not supported on ARM platforms. The way how I do that is illustrated below. So, I need to somehow circumvent this issue with @cfunction. How would that be done?

What I used to do is something like below (what matters the most is just the beginning part). The use of @cfunction is giving me error on Apple Silicon.

ERROR: cfunction: closures are not supported on this platform

using GSL
# For Mac
using LinearAlgebra, OpenBLAS32_jll
LinearAlgebra.BLAS.lbt_forward(OpenBLAS32_jll.libopenblas_path)

function solve!(f!, x0::AbstractArray, maxiter)
    # Do what GSL.@gsl_multiroot_function does to wrap the function
    function _f(x_vec, p, y_vec)
        x = GSL.wrap_gsl_vector(x_vec)
        y = GSL.wrap_gsl_vector(y_vec)
        f!(y, x)
        return Cint(GSL.GSL_SUCCESS)
    end
    # Does not work on ARM !
    cfunc = @cfunction($_f, Cint, (Ptr{GSL.gsl_vector}, Ptr{Cvoid}, Ptr{GSL.gsl_vector}))
    N = length(x0)
    gsl_func = GSL.gsl_multiroot_function(Base.unsafe_convert(Ptr{Cvoid}, cfunc), N, 0)

    # The remaining code solves the root-finding problem
    vinit = GSL.vector_alloc(N)
    solver = GSL.multiroot_fsolver_alloc(GSL.gsl_multiroot_fsolver_hybrids, N)

    for (i, v) in enumerate(x0)
        v = convert(Cdouble, v)
        GSL.vector_set(vinit, i-1, v)
    end

    GSL.multiroot_fsolver_set(solver, gsl_func, vinit)

    iter = 0
    while iter < maxiter
        iter += 1
        status = GSL.multiroot_fsolver_iterate(solver)
    end
    x = GSL.multiroot_fsolver_root(solver)
    # Copy the solution
    sol = Vector{Float64}(unsafe_wrap(Array{Cdouble}, GSL.vector_ptr(x, 0), N))
    GSL.multiroot_fsolver_free(solver)
    GSL.vector_free(vinit)
    return sol
end

# The function for the actual scenario is only known in run time
function f!(x::Array{Float64}, out::Array{Float64})
    out[1:5] =  (x .- collect(1.0:5.0)).^2
end

x0 = [1.0, 5.0, 2.0, 1.5, -1.0]

solve!(f!, x0, 5)

Thank you for any advice!

2

OK. I figured it out by myself after studying this amazing blog: Passing Julia Callback Functions to C (julialang.org)

The trick is to make use of the parameter argument available for GSL functions, which avoids the need for using any cfunction closure.