HAdaptiveIntegration

Overview

HAdaptiveIntegration.jl is a Julia package for approximating integrals of functions over various predefined AbstractDomains. It uses embedded cubature rules to build error estimates, and refines the integration domain by splitting its mesh elements until a certain tolerance is reached. Features include:

  • Adaptive integration over simplices and orthotope of any dimension
  • Use of efficient (tabulated) cubatures for low-dimensional simplices (triangle and tetrahedron) and orthotope (rectangle and cuboid)
  • Support for custom cubature rules
  • Arbitrary precision arithmetic

Installation

The package can be installed with the Julia package manager. From the Julia REPL, type ] to enter the Pkg REPL mode and run

pkg> add HAdaptiveIntegration

Or, equivalently, via the Julia Pkg API

julia> import Pkg; Pkg.add("HAdaptiveIntegration")

Then use

using HAdaptiveIntegration

Basic usage

The main function exported by this package is integrate(f, Ω), which is used to approximate

\[I = \int_{\Omega} f(x) \, \mathrm{d}x\]

where $\Omega \subset \mathbb{R}^d$ is a AbstractDomain object, and $f \colon \mathbb{R}^d \to \mathbb{F}$ is a function. Here is a simple example: first, we define a function,

fct = x -> cis(sum(x)) / (sum(abs2, x) + 1e-2)
Function signature

The function f must accept a single argument x which is an abstract vector of length d, the dimension of the domain (concretely, f is called through f(::SVector)). The return type T of f can be any type that supports the operations +(T, T), norm(T), and multiplication by a scalar (e.g. vectors or matrices).

Domains are constructed using the following functions (see their respective docstrings for more details):

Simplices

To integrate the above function over a $d$-dimensional simplices (triangle, tetrahedron, ...), defined by their vertices, we can use

  • Triangle

    I, E = integrate(fct, Triangle((0, 0), (1, 0), (0, 1)))
    (2.9073045602673253 + 1.2009254989232099im, 4.646511607369907e-8)

    The result I is the integral of f over a triangle with vertices (0,0), (1,0), and (0,1), and E is an error estimate.

  • Tetrahedron

    I, E = integrate(fct, Tetrahedron((0, 0, 0), (1, 0, 0), (0, 1, 0), (0, 0, 1)))
    (0.6767480949657807 + 0.4437444630481043im, 1.2058253541977129e-8)
  • $4$-simplex

    I, E = integrate(
             fct,
             Simplex((0, 0, 0, 0), (1, 0, 0, 0), (0, 1, 0, 0), (0, 0, 1, 0), (0, 0, 0, 1));
             rtol=1e-4
           )
    (0.14934972051376968 + 0.1243949385231301im, 1.773289010714962e-5)

    The keyword arguments atol and rtol can be used to control the desired absolute and relative error tolerances, respectively.

Orthotopes (a.k.a. hyperrectangle)

To integrate the same function over a $d$-dimensional axis-aligned orthotope (rectangle, cuboid, and hyperrectangle), defined by their low and high corners, we can use

  • Rectangle

    I, E = integrate(fct, Rectangle((0, 0), (1, 1)))
    (3.061858492661162 + 1.7034087551176278im, 4.8525654841517435e-8)
  • Cuboid

    I, E = integrate(fct, Cuboid((0, 0, 0), (1, 1, 1)))
    (0.7332718951640882 + 1.2361267744925115im, 2.083830752837712e-8)
  • $4$-orthotope (hyperrectangle)

    I, E = integrate(fct, Orthotope((0, 0, 0, 0), (1, 1, 1, 1)); rtol=1e-4)
    (-0.02397615741770491 + 0.8405904960226667im, 8.287878623984825e-5)

Going further

In the previous examples we covered the basic usage of the integrate function. There are, however, other options that can be passed to integrate in order to customize various aspects of the underlying algorithm (e.g. passing a different cubature rule, using a buffer to avoid memory allocations, etc.). For more details, see the docstring of the integrate function, as well as the next section on advanced usage.