# C callback for GSL

**URL:** <https://discourse.julialang.org/t/c-callback-for-gsl/7058>\
**Category:** Numerics\
**Created:** [November 14, 2017, 5:41am UTC](https://discourse.julialang.org/t/c-callback-for-gsl/7058 "2017-11-14T05:41:06Z")\
**Posts on this page:** 2\
**Page:** 1

<div class="post-metadata">

**Author:** ![mzaffalon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mzaffalon/32/214168_2.png) [@mzaffalon](https://discourse.julialang.org/u/mzaffalon)\
**Post date:** [November 14, 2017, 5:41am UTC](https://discourse.julialang.org/t/c-callback-for-gsl/7058/1 "2017-11-14T05:41:06Z")

</div>

I would like to add [fixed points quadrature](https://www.gnu.org/software/gsl/doc/html/integration.html#fixed-point-quadratures) to `GSL.jl` and to understand how callbacks for GSL work, I read the [blog post](https://julialang.org/blog/2013/05/callback#passing-closures-in-data-structures) as linked [in the documentation](https://docs.julialang.org/en/stable/manual/calling-c-and-fortran-code/#More-About-Callbacks-1).

Here is the code for Gauss-Legendre adaptive quadrature routine that includes the changes to reflect the recent PR by @yuyichao:

```julia
function gsl_function_wrap(x::Cdouble, gsl_func::Ptr{Void})
    jfunc = unsafe_pointer_to_objref(gsl_func)::Function
    convert(Cdouble, jfunc(x))::Cdouble
end

struct GSL_Function
    cfunc::Ptr{Void}
    jfunc::Any
    function GSL_Function(jfunc::Function)
        cfunc = cfunction(gsl_function_wrap, Cdouble, (Cdouble, Ptr{Void}))
        new(cfunc, jfunc)
    end
end

function integration_qag(julia_func::Function, a::Real, b::Real)
    w = ccall((:gsl_integration_workspace_alloc, "libgsl"), Ptr{Void}, (Csize_t,), 10^7)
    result = Ref{Cdouble}(0)
    abserr = Ref{Cdouble}(0)
    errno = ccall((:gsl_integration_qag, "libgsl"), Cint,
                  (Ref{GSL_Function}, Cdouble, Cdouble, Cdouble, Cdouble, Csize_t,
                   Cint, Ptr{Void}, Ref{Cdouble}, Ref{Cdouble}),
                  GSL_Function(julia_func), a, b, 1e-12, 1e-12, 10^3,
                  1, w, result, abserr)
    ccall((:gsl_integration_workspace_free, "libgsl"), Void, (Ptr{Void},), w)
    return result[], abserr[]
end

```

Note that `ccall((:gsl_integration_qag,...), ...)` in the blog post is missing the `epsabs` argument.

Now for the question: the code above seems to work, for instance try

```julia
f(x)=x^2
integration_qag(f, -2.0, 1.0)

```

but I would like to make sure this is correct. In particular, I understand that

- `GSL_Function(julia_func)` should be passed as `Ref{GSL_Function}`;
- the workspace `w` allocated by the C function `gsl_integration_workspace_alloc` can be passed as `Ptr{Void}`,
- the type assertion in `gsl_function_wrap` is required `jfunc = unsafe_pointer_to_objref(gsl_func)::Function` and
- putting `const cfunc = cfunction(gsl_function_wrap, Cdouble, (Cdouble, Ptr{Void}))` in the global scope returns a NULL pointer if the code above is in a module.

The blog post should be updated (I find it too useful to be removed from the documentation) since it uses the `&` construct that is deprecated and passes Julia objects as `Ptr` when `Ref` instead should be used. I can submit a PR, if somebody can confirm the code above is correct (I will add the other missing parameters to the Julia wrapper) and [if the syntax has settled](https://github.com/stevengj/Cubature.jl/pull/23#issuecomment-270642596).

---

<div class="post-metadata">

**Author:** ![mzaffalon](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mzaffalon/32/214168_2.png) [@mzaffalon](https://discourse.julialang.org/u/mzaffalon)\
**Post date:** [November 14, 2017, 6:41am UTC](https://discourse.julialang.org/t/c-callback-for-gsl/7058/2 "2017-11-14T06:41:03Z")

</div>

I should add that before writing the post, I was too missing the extra parameter `abserr` and getting segfaults, which disappeared once I used the correct C function signature.
