[ANN] Cuba.jl: library for multidimensional numerical integration

[Mosè Giordano]

Dear all,

I am proud to announce Cuba.jl a library for multidimensional numerical integration with four independent algorithms: Vegas, Suave, Divonne, and Cuhre (this algorithm is the same used in Cubature.jl).  This package is just a wrapper around Cuba Library, written in C by Thomas Hahn.

Cuba.jl is a registered Julia package, so you can install it with the package manager:

Pkg.add("Cuba")

The package is usable, but I must admit user interface is not optimal.  One has to define a function of this type:

function integrand(ndim::Cint, xx::Ptr{Cdouble}, ncomp::Cint, ff::Ptr{Cdouble},
                   userdata::Ptr{Void})
    # Take arrays from "xx" and "ff" pointers.
    x = pointer_to_array(xx, (ndim,))
    f = pointer_to_array(ff, (ncomp,))
    # Do calculations on "f" here
    #   ...
    # Store back the results to "ff"
    ff = pointer_from_objref(f)
return Cint(0)::Cint
end

and then call one of the four integrator functions available with this syntax:

Vegas(integrand, ndim, ncomp[, keywords...])
Suave(integrand, ndim, ncomp[, keywords...])
Divonne(integrand, ndim, ncomp[, keywords...])
Cuhre(integrand, ndim, ncomp[, keywords...])

Issue #3 tracks this problem, if someone wants to help on this is warmly welcome.

Documentation of the package is available at https://cubajl.readthedocs.org/ and you can also download the PDF version of the manual from https://media.readthedocs.org/pdf/cubajl/latest/cubajl.pdf  In particular, there is a section with some useful examples: https://cubajl.readthedocs.org/en/latest/#examples

Even though Cuba.jl does not support parallelization (see issue #1), its performance is comparable with those of equivalent codes written in C or Fortran relying on Cuba Library: https://github.com/giordano/Cuba.jl#performance

Cuba.jl is released under the terms of LGPLv3 and is available for GNU/Linux and OS X (Windows support is currently missing, Cubature.jl is a better alternative for that platform).

Feel free to share your comments and suggestions on this package!

Cheers,
Mosè

[Mosè Giordano]

A quick update: with new version 0.0.5 the syntax of the integrand function has been changed and now the integrand function should be passed as in Cubature,jl, for the case of vector-valued integrands.  So, for example, to integrate log(x)/sqrt(x) for x in [0, 1], one can use one of the following commands:

Vegas( (x,f) -> f[1] = log(x[1])/sqrt(x[1]), 1, 1)
Suave( (x,f) -> f[1] = log(x[1])/sqrt(x[1]), 1, 1)
Divonne( (x,f) -> f[1] = log(x[1])/sqrt(x[1]), 1, 1)
Cuhre( (x,f) -> f[1] = log(x[1])/sqrt(x[1]), 1, 1)

More details and examples are available in the manual: http://cubajl.readthedocs.org/

Let me thank Steven G. Johnson for his help.

Bye,
Mosè

[Chris Rackauckas]

Nice! I was looking for a Vegas function awhile ago and the GSL.jl one isn't bound correctly yet. This will be a nice addition to the Julia package list. Good job!

[Oliver Schulz]

Very nice, Mose! I thought I might have to write a Julia wrapper for Cuba at some point, but you beat me to it (which I'm grateful for!). :-)

I noticed that you're maintaining a modified version of Cuba for this - I guess because the upstream Cuba doesn't support "make shared"? Are you in contact with Thomas? He may be willing to add that - if you like, I can ask, his office is two floors up from mine. ;-)

[Kaj Wiik]

Very nice work indeed!

Is it possible to pass 'userdata' to integrand in Julia? It is mentioned in sources but is not listed in keywords.

Thanks,
Kaj

[Steven G. Johnson]

On Tuesday, April 12, 2016 at 4:51:07 PM UTC-4, Kaj Wiik wrote:

Very nice work indeed!

Is it possible to pass 'userdata' to integrand in Julia? It is mentioned in sources but is not listed in keywords.

You don't need an explicit userdata argument -- this is just a crude simulation of closures in C, but Julia has true closures and lexical scoping.

For example

myvar = 7

Vegas( (x,f) -> f[1] = log(x[1] + myvar)/sqrt(x[1]), 1, 1)

will "capture" the value of myvar and pass it through Vegas to the integrand, via the closure.

[Mosè Giordano]

Hi Chris,

thanks for your comments :-)  Regarding integration algorithm, unless you have specific reasons to use Vegas I warmly recommend to use Cuhre, which in my experience is the best in terms of speed and precision.  There must be a reason if Cubature,jl, the only widely used numerical integration packages so far, implements the same algorithm ;-)  Only once I found Divonne more useful because I had a wild function with peaks in known positions and you can tell Divonne to increase sampling around those points.

Bye,

Mosè

[Mosè Giordano]

Hi Oliver,

2016-04-12 21:38 GMT+02:00 Oliver Schulz <oliver...@tu-dortmund.de>:
> Very nice, Mose! I thought I might have to write a Julia wrapper for Cuba at
> some point, but you beat me to it (which I'm grateful for!). :-)

Thanks!

> I noticed that you're maintaining a modified version of Cuba for this - I
> guess because the upstream Cuba doesn't support "make shared"?

Yes, exactly, these are the differences with vanilla Cuba Library:
https://github.com/giordano/cuba/compare/julia Only configure and
makefile are involved.

> Are you in
> contact with Thomas? He may be willing to add that - if you like, I can ask,
> his office is two floors up from mine. ;-)

I had been in contact with him in the past, and I think I'll contact
him very soon again.

Ciao,
Mosè

[Kaj Wiik]

On Wednesday, April 13, 2016 at 6:14:29 AM UTC+3, Steven G. Johnson wrote:

 

You don't need an explicit userdata argument -- this is just a crude simulation of closures in C, but Julia has true closures and lexical scoping.

For example

myvar = 7

Vegas( (x,f) -> f[1] = log(x[1] + myvar)/sqrt(x[1]), 1, 1)

will "capture" the value of myvar and pass it through Vegas to the integrand, via the closure.

Ah, of course. Some Julia magic :-).

Thanks!

Kaj
 

[Mosè Giordano]

Hi Kaj,

thanks for your question, it served me as inspiration to add another example to the manual: https://cubajl.readthedocs.org/en/latest/#passing-data-to-the-integrand-function

Bye,
Mosè

[Kaj Wiik]

Perfect!

The documentation looks very good, thanks!

Cheers,
Kaj

[Mosè Giordano]

Another relevant update for this package: Cuba.jl can now be used on Windows (GNU/Linux and Mac OS were already supported).  I managed to cross compile the shared library for Windows on GNU/Linux, so during building of the package a prebuilt dll is downloaded for Windows users.

If you want to check out Cuba.jl, remember to issue Pkg.update() before installing it, in order to get the latest version (there have been two releases today, v0.1.3 being the latest one).

The only open issue is support for parallelization (https://github.com/giordano/Cuba.jl/issues/1), a feature provided by the native Cuba library in C, which however doesn't seem to play nicely with Julia because of the use of fork() function.

Cheers,
Mosè