From fe9178f9681af397c1bf70ce46a85d9719c1ff13 Mon Sep 17 00:00:00 2001 From: "Documenter.jl" Date: Thu, 22 Feb 2024 05:38:24 +0000 Subject: [PATCH] build based on 4593763 --- dev/.documenter-siteinfo.json | 2 +- dev/basics/FAQ/index.html | 2 +- dev/basics/IntegralFunction/index.html | 2 +- dev/basics/IntegralProblem/index.html | 2 +- dev/basics/SampledIntegralProblem/index.html | 8 ++++---- dev/basics/solve/index.html | 2 +- dev/index.html | 2 +- dev/search_index.js | 2 +- dev/solvers/IntegralSolvers/index.html | 2 +- dev/tutorials/caching_interface/index.html | 2 +- dev/tutorials/differentiating_integrals/index.html | 2 +- dev/tutorials/numerical_integrals/index.html | 2 +- 12 files changed, 15 insertions(+), 15 deletions(-) diff --git a/dev/.documenter-siteinfo.json b/dev/.documenter-siteinfo.json index 8467491..a2ca140 100644 --- a/dev/.documenter-siteinfo.json +++ b/dev/.documenter-siteinfo.json @@ -1 +1 @@ -{"documenter":{"julia_version":"1.10.1","generation_timestamp":"2024-02-22T04:57:02","documenter_version":"1.2.1"}} \ No newline at end of file +{"documenter":{"julia_version":"1.10.1","generation_timestamp":"2024-02-22T05:38:19","documenter_version":"1.2.1"}} \ No newline at end of file diff --git a/dev/basics/FAQ/index.html b/dev/basics/FAQ/index.html index 2727954..c367c1a 100644 --- a/dev/basics/FAQ/index.html +++ b/dev/basics/FAQ/index.html @@ -1,2 +1,2 @@ -Frequently Asked Questions · Integrals.jl
GitHub

Frequently Asked Questions

How should I use the in-place interface?

The in-place interface allows evaluating vector-valued integrands without allocating an output array. This can be beneficial for reducing allocations when integrating many functions simultaneously or to make use of existing in-place code. However, note that not all algorithms use in-place operations under the hood, i.e. HCubatureJL(), and may still allocate.

You can construct an IntegralFunction(f, prototype), where f is of the form f(y, u, p) where prototype is of the desired type and shape of y.

For small array outputs of a known size, consider using StaticArrays.jl for the return value of your integrand.

How should I use the batch interface?

The batch interface allows evaluating one (or more) integrals simultaneously at different points, which maximizes the parallelism for a given algorithm.

You can construct an out-of-place BatchIntegralFunction(bf) where bf is of the form bf(u, p) = stack(x -> f(x, p), eachslice(u; dims=ndims(u))), where f is the (unbatched) integrand.

You can construct an in-place BatchIntegralFunction(bf, prototype), where bf is of the form bf(y, u, p) = foreach((y,x) -> f(y,x,p), eachslice(y, dims=ndims(y)), eachslice(x, dims=ndims(x))).

Note that not all algorithms use in-place batched operations under the hood, i.e. QuadGKJL().

What should I do if my solution is not converged?

Certain algorithms, such as QuadratureRule used a fixed number of points to calculate an integral and cannot provide an error estimate. In this case, you have to increase the number of points and check the convergence yourself, which will depend on the accuracy of the rule you choose.

For badly-behaved integrands, such as (nearly) singular and highly oscillatory functions, most algorithms will fail to converge and either throw an error or silently return the incorrect result. In some cases Integrals.jl can provide an error code when things go wrong, but otherwise you can always check if the error estimate for the integral is less than the requested tolerance, e.g. sol.resid < max(abstol, reltol*norm(sol.u)). Sometimes using a larger tolerance or higher precision arithmetic may help.

How can I integrate arbitrarily-spaced data?

See SampledIntegralProblem.

How can I integrate on arbitrary geometries?

You can't, since Integrals.jl currently supports integration on hypercubes because that is what lower-level packages implement.

I don't see algorithm X or quadrature scheme Y ?

Fixed quadrature rules from other packages can be used with QuadratureRule. Otherwise, feel free to open an issue or pull request.

Can I take derivatives with respect to the limits of integration?

Currently this is not implemented.

+Frequently Asked Questions · Integrals.jl
GitHub

Frequently Asked Questions

How should I use the in-place interface?

The in-place interface allows evaluating vector-valued integrands without allocating an output array. This can be beneficial for reducing allocations when integrating many functions simultaneously or to make use of existing in-place code. However, note that not all algorithms use in-place operations under the hood, i.e. HCubatureJL, and may still allocate.

You can construct an IntegralFunction(f, prototype), where f is of the form f(y, u, p) where prototype is of the desired type and shape of y.

For small array outputs of a known size, consider using StaticArrays.jl for the return value of your integrand.

How should I use the batch interface?

The batch interface allows evaluating one (or more) integrals simultaneously at different points, which maximizes the parallelism for a given algorithm.

You can construct an out-of-place BatchIntegralFunction(bf) where bf is of the form bf(u, p) = stack(x -> f(x, p), eachslice(u; dims=ndims(u))), where f is the (unbatched) integrand. For interoperability with as many algorithms as possible, it is important that your out-of-place batch integrand accept an empty array of quadrature points and still return an output with a size and type consistent with the non-empty case.

You can construct an in-place BatchIntegralFunction(bf, prototype), where bf is of the form bf(y, u, p) = foreach((y,x) -> f(y,x,p), eachslice(y, dims=ndims(y)), eachslice(x, dims=ndims(x))).

Note that not all algorithms use in-place batched operations under the hood, i.e. QuadGKJL.

What should I do if my solution is not converged?

Certain algorithms, such as QuadratureRule used a fixed number of points to calculate an integral and cannot provide an error estimate. In this case, you have to increase the number of points and check the convergence yourself, which will depend on the accuracy of the rule you choose.

For badly-behaved integrands, such as (nearly) singular and highly oscillatory functions, most algorithms will fail to converge and either throw an error or silently return the incorrect result. In some cases Integrals.jl can provide an error code when things go wrong, but otherwise you can always check if the error estimate for the integral is less than the requested tolerance, e.g. sol.resid < max(abstol, reltol*norm(sol.u)). Sometimes using a larger tolerance or higher precision arithmetic may help.

How can I integrate arbitrarily-spaced data?

See SampledIntegralProblem.

How can I integrate on arbitrary geometries?

You can't, since Integrals.jl currently supports integration on hypercubes because that is what lower-level packages implement.

I don't see algorithm X or quadrature scheme Y ?

Fixed quadrature rules from other packages can be used with QuadratureRule. Otherwise, feel free to open an issue or pull request.

My integrand works with algorithm X but fails on algorithm Y

While bugs are not out of the question, certain algorithms, especially those implemented in C, are not compatible with arbitrary Julia types and have to return specific numeric types or arrays thereof. In some cases, such as ArblibJL, it is also expected that the integrand work with a custom quadrature point type. Moreover, some algorithms, such as VEGAS, only support scalar integrands. For more details see the solver page.

Can I take derivatives with respect to the limits of integration?

Currently this is not implemented.

diff --git a/dev/basics/IntegralFunction/index.html b/dev/basics/IntegralFunction/index.html index a914440..dacd60b 100644 --- a/dev/basics/IntegralFunction/index.html +++ b/dev/basics/IntegralFunction/index.html @@ -1,3 +1,3 @@ Integral Functions · Integrals.jl
GitHub

Integral Functions

SciMLBase.IntegralFunction — Type
IntegralFunction{iip,specialize,F,T} <: AbstractIntegralFunction{iip}

A representation of an integrand f defined by:

\[f(u, p)\]

For an in-place form of f see the iip section below for details on in-place or out-of-place handling.

IntegralFunction{iip,specialize}(f, [integrand_prototype])

Note that only f is required, and in the case of inplace integrands a mutable container integrand_prototype to store the result of the integrand. If integrand_prototype is present, f is interpreted as in-place, and otherwise f is assumed to be out-of-place.

iip: In-Place vs Out-Of-Place

Out-of-place functions must be of the form $y = f(u, p)$ and in-place functions of the form $f(y, u, p)$, where y is a number or array containing the output. Since f is allowed to return any type (e.g. real or complex numbers or arrays), in-place functions must provide a container integrand_prototype that is of the right type and size for the variable $y$, and the result is written to this container in-place. When in-place forms are used, in-place array operations, i.e. broadcasting, may be used by algorithms to reduce allocations. If integrand_prototype is not provided, f is assumed to be out-of-place.

specialize

This field is currently unused

Fields

The fields of the IntegralFunction type directly match the names of the inputs.

source
SciMLBase.BatchIntegralFunction — Type
BatchIntegralFunction{iip,specialize,F,T} <: AbstractIntegralFunction{iip}

A batched representation of an (non-batched) integrand f(u, p) that can be evaluated at multiple points simultaneously using threads, the gpu, or distributed memory defined by:

\[by = bf(bu, p)\]

Here we prefix variables with b to indicate they are batched variables, which implies that they are arrays whose last dimension is reserved for batching different evaluation points, or function values, and may be of a variable length. $bu$ is an array whose elements correspond to distinct evaluation points to f, and bf is a function to evaluate f 'point-wise' so that f(bu[..., i], p) == bf(bu, p)[..., i]. For example, a simple batching implementation of a scalar, univariate function is via broadcasting: bf(bu, p) = f.(bu, Ref(p)), although this interface exists in order to allow user parallelization. In general, the integration algorithm is allowed to vary the number of evaluation points between subsequent calls to bf.

For an in-place form of bf see the iip section below for details on in-place or out-of-place handling.

BatchIntegralFunction{iip,specialize}(bf, [integrand_prototype];
-                                     max_batch=typemax(Int))

Note that only bf is required, and in the case of inplace integrands a mutable container integrand_prototype to store a batch of integrand evaluations, with a last "batching" dimension.

The keyword max_batch is used to set a soft limit on the number of points to batch at the same time so that memory usage is controlled.

If integrand_prototype is present, bf is interpreted as in-place, and otherwise bf is assumed to be out-of-place.

iip: In-Place vs Out-Of-Place

Out-of-place functions must be of the form by = bf(bu, p) and in-place functions of the form bf(by, bu, p) where by is a batch array containing the output. Since the algorithm may vary the number of points to batch, the batching dimension can be of any length, including zero, and since bf is allowed to return arrays of any type (e.g. real or complex) or size, in-place functions must provide a container integrand_prototype of the desired type and size for by. If integrand_prototype is not provided, bf is assumed to be out-of-place.

In the out-of-place case, we require f(bu[..., i], p) == bf(bu, p)[..., i], and certain algorithms, such as those implemented in C, may infer the type or shape of by by calling bf with an empty array of input points, i.e. bu with size(bu)[end] == 0. Then it is expected for the resulting by to have the same type and size(by)[begin:end-1] for all subsequent calls.

When the in-place form is used, we require f(by[..., i], bu[..., i], p) == bf(by, bu, p)[..., i] and size(by)[begin:end-1] == size(integrand_prototype)[begin:end-1]. The algorithm should always pass the integrand by arrays that are similar to integrand_prototype, and may use views and in-place array operations to reduce allocations.

specialize

This field is currently unused

Fields

The fields of the BatchIntegralFunction type are f, corresponding to bf above, and integrand_prototype.

source
+ max_batch=typemax(Int))

Note that only bf is required, and in the case of inplace integrands a mutable container integrand_prototype to store a batch of integrand evaluations, with a last "batching" dimension.

The keyword max_batch is used to set a soft limit on the number of points to batch at the same time so that memory usage is controlled.

If integrand_prototype is present, bf is interpreted as in-place, and otherwise bf is assumed to be out-of-place.

iip: In-Place vs Out-Of-Place

Out-of-place functions must be of the form by = bf(bu, p) and in-place functions of the form bf(by, bu, p) where by is a batch array containing the output. Since the algorithm may vary the number of points to batch, the batching dimension can be of any length, including zero, and since bf is allowed to return arrays of any type (e.g. real or complex) or size, in-place functions must provide a container integrand_prototype of the desired type and size for by. If integrand_prototype is not provided, bf is assumed to be out-of-place.

In the out-of-place case, we require f(bu[..., i], p) == bf(bu, p)[..., i], and certain algorithms, such as those implemented in C, may infer the type or shape of by by calling bf with an empty array of input points, i.e. bu with size(bu)[end] == 0. Then it is expected for the resulting by to have the same type and size(by)[begin:end-1] for all subsequent calls.

When the in-place form is used, we require f(by[..., i], bu[..., i], p) == bf(by, bu, p)[..., i] and size(by)[begin:end-1] == size(integrand_prototype)[begin:end-1]. The algorithm should always pass the integrand by arrays that are similar to integrand_prototype, and may use views and in-place array operations to reduce allocations.

specialize

This field is currently unused

Fields

The fields of the BatchIntegralFunction type are f, corresponding to bf above, and integrand_prototype.

source diff --git a/dev/basics/IntegralProblem/index.html b/dev/basics/IntegralProblem/index.html index 2d1c0b4..4f10ce7 100644 --- a/dev/basics/IntegralProblem/index.html +++ b/dev/basics/IntegralProblem/index.html @@ -2,4 +2,4 @@ Integral Problems · Integrals.jl
GitHub

Integral Problems

SciMLBase.IntegralProblem — Type

Defines an integral problem. Documentation Page: https://docs.sciml.ai/Integrals/stable/

Mathematical Specification of an Integral Problem

Integral problems are multi-dimensional integrals defined as:

\[\int_{lb}^{ub} f(u,p) du\]

where p are parameters. u is a Number or AbstractVector whose geometry matches the space being integrated. This space is bounded by the lowerbound lb and upperbound ub, which are Numbers or AbstractVectors with the same geometry as u.

Problem Type

Constructors

IntegralProblem(f::AbstractIntegralFunction,domain,p=NullParameters(); kwargs...)
 IntegralProblem(f::AbstractIntegralFunction,lb,ub,p=NullParameters(); kwargs...)
 IntegralProblem(f,domain,p=NullParameters(); nout=nothing, batch=nothing, kwargs...)
-IntegralProblem(f,lb,ub,p=NullParameters(); nout=nothing, batch=nothing, kwargs...)
  • f: the integrand, callable function y = f(u,p) for out-of-place (default) or an IntegralFunction or BatchIntegralFunction for inplace and batching optimizations.
  • domain: an object representing an integration domain, i.e. the tuple (lb, ub).
  • lb: DEPRECATED: Either a number or vector of lower bounds.
  • ub: DEPRECATED: Either a number or vector of upper bounds.
  • p: The parameters associated with the problem.
  • nout: DEPRECATED (see IntegralFunction): length of the vector output of the integrand (by default the integrand is assumed to be scalar)
  • batch: DEPRECATED (see BatchIntegralFunction): number of points the integrand can evaluate simultaneously (by default there is no batching)
  • kwargs: Keyword arguments copied to the solvers.

Additionally, we can supply iip like IntegralProblem{iip}(...) as true or false to declare at compile time whether the integrator function is in-place.

Fields

The fields match the names of the constructor arguments.

source

The correct shape of the variables (inputs) u and the values (outputs) y of the integrand f depends on whether batching is used.

If batch == 0

single variable fmultiple variable f
scalar valued fu is a scalar, y is a scalaru is a vector, y is a scalar
vector valued fu is a scalar, y is a vectoru is a vector, y is a vector

If batch > 0

single variable fmultiple variable f
scalar valued fu is a vector, y is a vectoru is a matrix, y is a vector
vector valued fu is a vector, y is a matrixu is a matrix, y is a matrix

The last dimension is always used as the batching dimension, e.g., if u is a matrix, then u[:,i] is the ith point where the integrand will be evaluated.

+IntegralProblem(f,lb,ub,p=NullParameters(); nout=nothing, batch=nothing, kwargs...)

Additionally, we can supply iip like IntegralProblem{iip}(...) as true or false to declare at compile time whether the integrator function is in-place.

Fields

The fields match the names of the constructor arguments.

source

The correct shape of the variables (inputs) u and the values (outputs) y of the integrand f depends on whether batching is used.

If batch == 0

single variable fmultiple variable f
scalar valued fu is a scalar, y is a scalaru is a vector, y is a scalar
vector valued fu is a scalar, y is a vectoru is a vector, y is a vector

If batch > 0

single variable fmultiple variable f
scalar valued fu is a vector, y is a vectoru is a matrix, y is a vector
vector valued fu is a vector, y is a matrixu is a matrix, y is a matrix

The last dimension is always used as the batching dimension, e.g., if u is a matrix, then u[:,i] is the ith point where the integrand will be evaluated.

diff --git a/dev/basics/SampledIntegralProblem/index.html b/dev/basics/SampledIntegralProblem/index.html index 475b3ff..b751096 100644 --- a/dev/basics/SampledIntegralProblem/index.html +++ b/dev/basics/SampledIntegralProblem/index.html @@ -24,13 +24,13 @@ 1.0

Now, we can integrate this data set as follows:

problem = SampledIntegralProblem(y, x)
 method = TrapezoidalRule()
 solve(problem, method)
retcode: Success
-u: 0.33379501385041543

The exact answer is of course $ 1/3 $.

Details

SciMLBase.SampledIntegralProblem — Type

Defines a integral problem over pre-sampled data. Documentation Page: https://docs.sciml.ai/Integrals/stable/

Mathematical Specification of a data Integral Problem

Sampled integral problems are defined as:

\[\sum_i w_i y_i\]

where y_i are sampled values of the integrand, and w_i are weights assigned by a quadrature rule, which depend on sampling points x.

Problem Type

Constructors

SampledIntegralProblem(y::AbstractArray, x::AbstractVector; dim=ndims(y), kwargs...)
  • y: The sampled integrand, must be a subtype of AbstractArray. It is assumed that the values of y along dimension dim correspond to the integrand evaluated at sampling points x
  • x: Sampling points, must be a subtype of AbstractVector.
  • dim: Dimension along which to integrate. Defaults to the last dimension of y.
  • kwargs: Keyword arguments copied to the solvers.

Fields

The fields match the names of the constructor arguments.

source
CommonSolve.solve — Method
solve(prob::SampledIntegralProblem, alg::SciMLBase.AbstractIntegralAlgorithm; kwargs...)

Keyword Arguments

There are no keyword arguments used to solve SampledIntegralProblems

source

Non-equidistant grids

If the sampling points x are provided as an AbstractRange (constructed with the range function for example), faster methods are used that take advantage of the fact that the points are equidistantly spaced. Otherwise, general methods are used for non-uniform grids.

Example:

f = x -> x^7
+u: 0.33379501385041543

The exact answer is of course $ 1/3 $.

Details

SciMLBase.SampledIntegralProblem — Type

Defines a integral problem over pre-sampled data. Documentation Page: https://docs.sciml.ai/Integrals/stable/

Mathematical Specification of a data Integral Problem

Sampled integral problems are defined as:

\[\sum_i w_i y_i\]

where y_i are sampled values of the integrand, and w_i are weights assigned by a quadrature rule, which depend on sampling points x.

Problem Type

Constructors

SampledIntegralProblem(y::AbstractArray, x::AbstractVector; dim=ndims(y), kwargs...)
  • y: The sampled integrand, must be a subtype of AbstractArray. It is assumed that the values of y along dimension dim correspond to the integrand evaluated at sampling points x
  • x: Sampling points, must be a subtype of AbstractVector.
  • dim: Dimension along which to integrate. Defaults to the last dimension of y.
  • kwargs: Keyword arguments copied to the solvers.

Fields

The fields match the names of the constructor arguments.

source
CommonSolve.solve — Method
solve(prob::SampledIntegralProblem, alg::SciMLBase.AbstractIntegralAlgorithm; kwargs...)

Keyword Arguments

There are no keyword arguments used to solve SampledIntegralProblems

source

Non-equidistant grids

If the sampling points x are provided as an AbstractRange (constructed with the range function for example), faster methods are used that take advantage of the fact that the points are equidistantly spaced. Otherwise, general methods are used for non-uniform grids.

Example:

f = x -> x^7
 x = [0.0; sort(rand(1000)); 1.0]
 y = f.(x)
 problem = SampledIntegralProblem(y, x)
 method = TrapezoidalRule()
 solve(problem, method)
retcode: Success
-u: 0.1250030260094972

Evaluating multiple integrals at once

If the provided data set y is a multidimensional array, the integrals are evaluated across only one of its axes. For performance reasons, the last axis of the array y is chosen by default, but this can be modified with the dim keyword argument to the problem definition.

f1 = x -> x^2
+u: 0.1250040468299159

Evaluating multiple integrals at once

If the provided data set y is a multidimensional array, the integrals are evaluated across only one of its axes. For performance reasons, the last axis of the array y is chosen by default, but this can be modified with the dim keyword argument to the problem definition.

f1 = x -> x^2
 f2 = x -> x^3
 f3 = x -> x^4
 x = range(0, 1, length = 20)
@@ -47,10 +47,10 @@
 y = f.(x)
 problem = SampledIntegralProblem(y, x)
 method = TrapezoidalRule()
-solve(problem, method)
source
Integrals.SimpsonsRule — Type
SimpsonsRule

Struct for evaluating an integral via the Simpson's composite 1/3-3/8 rule over AbstractRanges (evenly spaced points) and Simpson's composite 1/3 rule for non-equidistant grids.

Example with equidistant data:

using Integrals
+solve(problem, method)
source
Integrals.SimpsonsRule — Type
SimpsonsRule

Struct for evaluating an integral via the Simpson's composite 1/3-3/8 rule over AbstractRanges (evenly spaced points) and Simpson's composite 1/3 rule for non-equidistant grids.

Example with equidistant data:

using Integrals
 f = x -> x^2
 x = range(0, 1, length=20)
 y = f.(x)
 problem = SampledIntegralProblem(y, x)
 method = SimpsonsRule()
-solve(problem, method)
source
+solve(problem, method)source diff --git a/dev/basics/solve/index.html b/dev/basics/solve/index.html index 70c3865..eddfc25 100644 --- a/dev/basics/solve/index.html +++ b/dev/basics/solve/index.html @@ -1,2 +1,2 @@ -Common Solver Options (Solve Keyword Arguments) · Integrals.jl
GitHub

Common Solver Options (Solve Keyword Arguments)

CommonSolve.solve — Method
solve(prob::IntegralProblem, alg::SciMLBase.AbstractIntegralAlgorithm; kwargs...)

Keyword Arguments

The arguments to solve are common across all of the quadrature methods. These common arguments are:

  • maxiters (the maximum number of iterations)
  • abstol (absolute tolerance in changes of the objective value)
  • reltol (relative tolerance in changes of the objective value)
source

Additionally, the extra keyword arguments are splatted to the library calls, so see the documentation of the integrator library for all the extra details. These extra keyword arguments are not guaranteed to act uniformly.

+Common Solver Options (Solve Keyword Arguments) · Integrals.jl
GitHub

Common Solver Options (Solve Keyword Arguments)

CommonSolve.solve — Method
solve(prob::IntegralProblem, alg::SciMLBase.AbstractIntegralAlgorithm; kwargs...)

Keyword Arguments

The arguments to solve are common across all of the quadrature methods. These common arguments are:

  • maxiters (the maximum number of iterations)
  • abstol (absolute tolerance in changes of the objective value)
  • reltol (relative tolerance in changes of the objective value)
source

Additionally, the extra keyword arguments are splatted to the library calls, so see the documentation of the integrator library for all the extra details. These extra keyword arguments are not guaranteed to act uniformly.

diff --git a/dev/index.html b/dev/index.html index dedb030..a0836ae 100644 --- a/dev/index.html +++ b/dev/index.html @@ -171,4 +171,4 @@ [83775a58] Zlib_jll v1.2.13+1 [8e850b90] libblastrampoline_jll v5.8.0+1 [8e850ede] nghttp2_jll v1.52.0+1 - [3f19e933] p7zip_jll v17.4.0+2

You can also download the manifest file and the project file.

+ [3f19e933] p7zip_jll v17.4.0+2

You can also download the manifest file and the project file.

diff --git a/dev/search_index.js b/dev/search_index.js index e4d88d8..d16b86e 100644 --- a/dev/search_index.js +++ b/dev/search_index.js @@ -1,3 +1,3 @@ var documenterSearchIndex = {"docs": -[{"location":"basics/SampledIntegralProblem/#Integrating-pre-sampled-data","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"","category":"section"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"In some cases, instead of a function that acts as integrand, one only possesses a list of data points y at a set of sampling locations x, that must be integrated. This package contains functionality for doing that.","category":"page"},{"location":"basics/SampledIntegralProblem/#Example","page":"Integrating pre-sampled data","title":"Example","text":"","category":"section"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"Say, by some means we have generated a dataset x and y:","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"using Integrals # hide\nf = x -> x^2\nx = range(0, 1, length = 20)\ny = f.(x)","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"Now, we can integrate this data set as follows:","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"problem = SampledIntegralProblem(y, x)\nmethod = TrapezoidalRule()\nsolve(problem, method)","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"The exact answer is of course $ 1/3 $.","category":"page"},{"location":"basics/SampledIntegralProblem/#Details","page":"Integrating pre-sampled data","title":"Details","text":"","category":"section"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"SciMLBase.SampledIntegralProblem\nsolve(::SampledIntegralProblem, ::SciMLBase.AbstractIntegralAlgorithm)","category":"page"},{"location":"basics/SampledIntegralProblem/#SciMLBase.SampledIntegralProblem","page":"Integrating pre-sampled data","title":"SciMLBase.SampledIntegralProblem","text":"Defines a integral problem over pre-sampled data. Documentation Page: https://docs.sciml.ai/Integrals/stable/\n\nMathematical Specification of a data Integral Problem\n\nSampled integral problems are defined as:\n\nsum_i w_i y_i\n\nwhere y_i are sampled values of the integrand, and w_i are weights assigned by a quadrature rule, which depend on sampling points x.\n\nProblem Type\n\nConstructors\n\nSampledIntegralProblem(y::AbstractArray, x::AbstractVector; dim=ndims(y), kwargs...)\n\ny: The sampled integrand, must be a subtype of AbstractArray. It is assumed that the values of y along dimension dim correspond to the integrand evaluated at sampling points x\nx: Sampling points, must be a subtype of AbstractVector.\ndim: Dimension along which to integrate. Defaults to the last dimension of y.\nkwargs: Keyword arguments copied to the solvers.\n\nFields\n\nThe fields match the names of the constructor arguments.\n\n\n\n\n\n","category":"type"},{"location":"basics/SampledIntegralProblem/#CommonSolve.solve-Tuple{SampledIntegralProblem, SciMLBase.AbstractIntegralAlgorithm}","page":"Integrating pre-sampled data","title":"CommonSolve.solve","text":"solve(prob::SampledIntegralProblem, alg::SciMLBase.AbstractIntegralAlgorithm; kwargs...)\n\nKeyword Arguments\n\nThere are no keyword arguments used to solve SampledIntegralProblems\n\n\n\n\n\n","category":"method"},{"location":"basics/SampledIntegralProblem/#Non-equidistant-grids","page":"Integrating pre-sampled data","title":"Non-equidistant grids","text":"","category":"section"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"If the sampling points x are provided as an AbstractRange (constructed with the range function for example), faster methods are used that take advantage of the fact that the points are equidistantly spaced. Otherwise, general methods are used for non-uniform grids.","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"Example:","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"using Integrals # hide\nf = x -> x^7\nx = [0.0; sort(rand(1000)); 1.0]\ny = f.(x)\nproblem = SampledIntegralProblem(y, x)\nmethod = TrapezoidalRule()\nsolve(problem, method)","category":"page"},{"location":"basics/SampledIntegralProblem/#Evaluating-multiple-integrals-at-once","page":"Integrating pre-sampled data","title":"Evaluating multiple integrals at once","text":"","category":"section"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"If the provided data set y is a multidimensional array, the integrals are evaluated across only one of its axes. For performance reasons, the last axis of the array y is chosen by default, but this can be modified with the dim keyword argument to the problem definition.","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"using Integrals # hide\nf1 = x -> x^2\nf2 = x -> x^3\nf3 = x -> x^4\nx = range(0, 1, length = 20)\ny = [f1.(x) f2.(x) f3.(x)]\nproblem = SampledIntegralProblem(y, x; dim = 1)\nmethod = SimpsonsRule()\nsolve(problem, method)","category":"page"},{"location":"basics/SampledIntegralProblem/#Supported-methods","page":"Integrating pre-sampled data","title":"Supported methods","text":"","category":"section"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"Right now, only the TrapezoidalRule and SimpsonsRule are supported.","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"TrapezoidalRule\nSimpsonsRule","category":"page"},{"location":"basics/SampledIntegralProblem/#Integrals.TrapezoidalRule","page":"Integrating pre-sampled data","title":"Integrals.TrapezoidalRule","text":"TrapezoidalRule\n\nStruct for evaluating an integral via the trapezoidal rule.\n\nExample with sampled data:\n\nusing Integrals\nf = x -> x^2\nx = range(0, 1, length=20)\ny = f.(x)\nproblem = SampledIntegralProblem(y, x)\nmethod = TrapezoidalRule()\nsolve(problem, method)\n\n\n\n\n\n","category":"type"},{"location":"basics/SampledIntegralProblem/#Integrals.SimpsonsRule","page":"Integrating pre-sampled data","title":"Integrals.SimpsonsRule","text":"SimpsonsRule\n\nStruct for evaluating an integral via the Simpson's composite 1/3-3/8 rule over AbstractRanges (evenly spaced points) and Simpson's composite 1/3 rule for non-equidistant grids.\n\nExample with equidistant data:\n\nusing Integrals\nf = x -> x^2\nx = range(0, 1, length=20)\ny = f.(x)\nproblem = SampledIntegralProblem(y, x)\nmethod = SimpsonsRule()\nsolve(problem, method)\n\n\n\n\n\n","category":"type"},{"location":"tutorials/caching_interface/#Integrals-with-Caching-Interface","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"","category":"section"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"Often, integral solvers allocate memory or reuse quadrature rules for solving different problems. For example, if one is going to solve the same integral for several parameters","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"using Integrals\n\nprob = IntegralProblem((x, p) -> sin(x * p), (0, 1), 14.0)\nalg = QuadGKJL()\n\nsolve(prob, alg)\n\nprob = remake(prob, p = 15.0)\nsolve(prob, alg)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"then it would be more efficient to allocate the heap used by quadgk across several calls, shown below by directly calling the library","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"using QuadGK\nsegbuf = QuadGK.alloc_segbuf()\nquadgk(x -> sin(14x), 0, 1, segbuf = segbuf)\nquadgk(x -> sin(15x), 0, 1, segbuf = segbuf)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"Integrals.jl's caching interface automates this process to reuse resources if an algorithm supports it and if the necessary types to build the cache can be inferred from prob. To do this with Integrals.jl, you simply init a cache, solve!, replace p, and solve again. This uses the SciML init interface","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"using Integrals\n\nprob = IntegralProblem((x, p) -> sin(x * p), (0, 1), 14.0)\nalg = QuadGKJL()\n\ncache = init(prob, alg)\nsol1 = solve!(cache)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"cache.p = 15.0\nsol2 = solve!(cache)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"The caching interface is intended for updating p and domain. Note that the types of these variables is not allowed to change. If it is necessary to change the integrand f instead of defining a new IntegralProblem, consider using FunctionWrappers.jl.","category":"page"},{"location":"tutorials/caching_interface/#Caching-for-sampled-integral-problems","page":"Integrals with Caching Interface","title":"Caching for sampled integral problems","text":"","category":"section"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"For sampled integral problems, it is possible to cache the weights and reuse them for multiple data sets.","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"using Integrals\n\nx = 0.0:0.1:1.0\ny = sin.(x)\n\nprob = SampledIntegralProblem(y, x)\nalg = TrapezoidalRule()\n\ncache = init(prob, alg)\nsol1 = solve!(cache)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"cache.y = cos.(x) # use .= to update in-place\nsol2 = solve!(cache)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"If the grid is modified, the weights are recomputed.","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"cache.x = 0.0:0.2:2.0\ncache.y = sin.(cache.x)\nsol3 = solve!(cache)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"For multi-dimensional datasets, the integration dimension can also be changed","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"using Integrals\n\nx = 0.0:0.1:1.0\ny = sin.(x) .* cos.(x')\n\nprob = SampledIntegralProblem(y, x)\nalg = TrapezoidalRule()\n\ncache = init(prob, alg)\nsol1 = solve!(cache)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"cache.dim = 1\nsol2 = solve!(cache)","category":"page"},{"location":"tutorials/numerical_integrals/#Numerically-Solving-Integrals","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"","category":"section"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"For basic multidimensional quadrature, we can construct and solve a IntegralProblem. The integral we want to evaluate is:","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"int_1^3int_1^3int_1^3 sum_1^3 sin(u_i) du_1du_2du_3","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"We can numerically approximate this integral:","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"using Integrals\nf(u, p) = sum(sin.(u))\ndomain = (ones(3), 3ones(3)) # (lb, ub)\nprob = IntegralProblem(f, domain)\nsol = solve(prob, HCubatureJL(); reltol = 1e-3, abstol = 1e-3)\nsol.u","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"where the first argument of IntegralProblem is the integrand, the second argument is the lower bound, and the third argument is the upper bound. p are the parameters of the integrand. In this case, there are no parameters, but still f must be defined as f(x,p) and not f(x). For an example with parameters, see the next tutorial. The first argument of solve is the problem we are solving, the second is an algorithm to solve the problem with. Then there are keywords which provides details how the algorithm should work, in this case tolerances how precise the numerical approximation should be.","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"We can also evaluate multiple integrals at once. We could create two IntegralProblems for this, but that is wasteful if the integrands share a lot of computation. For example, we also want to evaluate:","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"int_1^3int_1^3int_1^3 sum_1^3 cos(u_i) du_1du_2du_3","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"using Integrals\nf(u, p) = [sum(sin.(u)), sum(cos.(u))]\ndomain = (ones(3), 3ones(3)) # (lb, ub)\nprob = IntegralProblem(f, domain)\nsol = solve(prob, HCubatureJL(); reltol = 1e-3, abstol = 1e-3)\nsol.u","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"Another way to think about this is that the integrand is now a vector valued function. In general, we should be able to integrate any type that is in a vector space and supports addition and scalar multiplication, although Integrals.jl allows scalars and arrays. In the above example, the integrand was defined out-of-position. This means that a new output vector is created every time the function f is called. If we do not want these allocations, we can also define f in-position.","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"using Integrals, Cubature\nfunction f(y, u, p)\n y[1] = sum(sin.(u))\n y[2] = sum(cos.(u))\nend\nprototype = zeros(2)\ndomain = (ones(3), 3ones(3)) # (lb, ub)\nprob = IntegralProblem(IntegralFunction(f, prototype), domain)\nsol = solve(prob, CubatureJLh(); reltol = 1e-3, abstol = 1e-3)\nsol.u","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"where y is a cache to store the evaluation of the integrand and prototype is an instance of y with the desired type and shape. We needed to change the algorithm to CubatureJLh() because HCubatureJL() does not support in-position under the hood. f evaluates the integrand at a certain point, but most adaptive quadrature algorithms need to evaluate the integrand at multiple points in each step of the algorithm. We would thus often like to parallelize the computation. The batch interface allows us to compute multiple points at once. For example, here we do allocation-free multithreading with Cubature.jl:","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"using Integrals, Cubature, Base.Threads\nfunction f(y, u, p)\n Threads.@threads for i in 1:size(u, 2)\n y[1, i] = sum(sin, @view(u[:, i]))\n y[2, i] = sum(cos, @view(u[:, i]))\n end\nend\nprototype = zeros(2, 0)\ndomain = (ones(3), 3ones(3)) # (lb, ub)\nprob = IntegralProblem(BatchIntegralFunction(f, prototype), domain)\nsol = solve(prob, CubatureJLh(); reltol = 1e-3, abstol = 1e-3)\nsol.u","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"Both u and y changed from vectors to matrices, where each column is respectively a point the integrand is evaluated at or the evaluation of the integrand at the corresponding point. The prototype now has an extra dimension for batching that can be of size zero. Try to create yourself an out-of-position version of the above problem. For the full details of the batching interface, see the problem page.","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"If we would like to compare the results against Cuba.jl's Cuhre method, then the change is a one-argument change:","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"using Integrals\nusing Cuba\nf(u, p) = sum(sin.(u))\ndomain = (ones(3), 3ones(3)) # (lb, ub)\nprob = IntegralProblem(f, domain)\nsol = solve(prob, CubaCuhre(); reltol = 1e-3, abstol = 1e-3)\nsol.u","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"However, Cuhre does not support vector valued integrands. The solvers page gives an overview of which arguments each algorithm can handle.","category":"page"},{"location":"tutorials/numerical_integrals/#One-dimensional-integrals","page":"Numerically Solving Integrals","title":"One-dimensional integrals","text":"","category":"section"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"Integrals.jl also has specific solvers for integrals in a single dimension, such as QuadGKJL. For example, we can create our own sine function by integrating the cosine function from 0 to x.","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"using Integrals\nmy_sin(x) = solve(IntegralProblem((x, p) -> cos(x), (0.0, x)), QuadGKJL()).u\nx = 0:0.1:(2 * pi)\n@. my_sin(x) ≈ sin(x)","category":"page"},{"location":"tutorials/numerical_integrals/#Infinity-handling","page":"Numerically Solving Integrals","title":"Infinity handling","text":"","category":"section"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"Integrals.jl can also handle infinite integration bounds. For infinite upper bounds u is substituted with a+fract1-t, and the integral is thus transformed to:","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"int_a^infty f(u)du = int_0^1 fleft(a+fract1-tright)frac1(1-t)^2dt","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"Integrals with an infinite lower bound are handled in the same way. If both upper and lower bound are infinite, u is substituted with fract1-t^2,","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"int_-infty^infty f(u)du = int_-1^1 fleft(fract1-t^2right)frac1+t^2(1-t^2)^2dt","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"For multidimensional integrals, each variable with infinite bounds is substituted the same way. The details of the math behind these transforms can be found here.","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"As an example, let us integrate the standard bivariate normal probability distribution over the area above the horizontal axis, which should be equal to 05.","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"using Distributions\nusing Integrals\ndist = MvNormal(ones(2))\nf = (x, p) -> pdf(dist, x)\ndomain = ([-Inf, 0.0], [Inf, Inf]) # (lb, ub)\nprob = IntegralProblem(f, domain)\nsolve(prob, HCubatureJL(), reltol = 1e-3, abstol = 1e-3)","category":"page"},{"location":"basics/IntegralProblem/#prob","page":"Integral Problems","title":"Integral Problems","text":"","category":"section"},{"location":"basics/IntegralProblem/","page":"Integral Problems","title":"Integral Problems","text":"SciMLBase.IntegralProblem","category":"page"},{"location":"basics/IntegralProblem/#SciMLBase.IntegralProblem","page":"Integral Problems","title":"SciMLBase.IntegralProblem","text":"Defines an integral problem. Documentation Page: https://docs.sciml.ai/Integrals/stable/\n\nMathematical Specification of an Integral Problem\n\nIntegral problems are multi-dimensional integrals defined as:\n\nint_lb^ub f(up) du\n\nwhere p are parameters. u is a Number or AbstractVector whose geometry matches the space being integrated. This space is bounded by the lowerbound lb and upperbound ub, which are Numbers or AbstractVectors with the same geometry as u.\n\nProblem Type\n\nConstructors\n\nIntegralProblem(f::AbstractIntegralFunction,domain,p=NullParameters(); kwargs...)\nIntegralProblem(f::AbstractIntegralFunction,lb,ub,p=NullParameters(); kwargs...)\nIntegralProblem(f,domain,p=NullParameters(); nout=nothing, batch=nothing, kwargs...)\nIntegralProblem(f,lb,ub,p=NullParameters(); nout=nothing, batch=nothing, kwargs...)\n\nf: the integrand, callable function y = f(u,p) for out-of-place (default) or an IntegralFunction or BatchIntegralFunction for inplace and batching optimizations.\ndomain: an object representing an integration domain, i.e. the tuple (lb, ub).\nlb: DEPRECATED: Either a number or vector of lower bounds.\nub: DEPRECATED: Either a number or vector of upper bounds.\np: The parameters associated with the problem.\nnout: DEPRECATED (see IntegralFunction): length of the vector output of the integrand (by default the integrand is assumed to be scalar)\nbatch: DEPRECATED (see BatchIntegralFunction): number of points the integrand can evaluate simultaneously (by default there is no batching)\nkwargs: Keyword arguments copied to the solvers.\n\nAdditionally, we can supply iip like IntegralProblem{iip}(...) as true or false to declare at compile time whether the integrator function is in-place.\n\nFields\n\nThe fields match the names of the constructor arguments.\n\n\n\n\n\n","category":"type"},{"location":"basics/IntegralProblem/","page":"Integral Problems","title":"Integral Problems","text":"The correct shape of the variables (inputs) u and the values (outputs) y of the integrand f depends on whether batching is used.","category":"page"},{"location":"basics/IntegralProblem/","page":"Integral Problems","title":"Integral Problems","text":"If batch == 0","category":"page"},{"location":"basics/IntegralProblem/","page":"Integral Problems","title":"Integral Problems","text":" single variable f multiple variable f\nscalar valued f u is a scalar, y is a scalar u is a vector, y is a scalar\nvector valued f u is a scalar, y is a vector u is a vector, y is a vector","category":"page"},{"location":"basics/IntegralProblem/","page":"Integral Problems","title":"Integral Problems","text":"If batch > 0","category":"page"},{"location":"basics/IntegralProblem/","page":"Integral Problems","title":"Integral Problems","text":" single variable f multiple variable f\nscalar valued f u is a vector, y is a vector u is a matrix, y is a vector\nvector valued f u is a vector, y is a matrix u is a matrix, y is a matrix","category":"page"},{"location":"basics/IntegralProblem/","page":"Integral Problems","title":"Integral Problems","text":"The last dimension is always used as the batching dimension, e.g., if u is a matrix, then u[:,i] is the ith point where the integrand will be evaluated.","category":"page"},{"location":"basics/solve/#Common-Solver-Options-(Solve-Keyword-Arguments)","page":"Common Solver Options (Solve Keyword Arguments)","title":"Common Solver Options (Solve Keyword Arguments)","text":"","category":"section"},{"location":"basics/solve/","page":"Common Solver Options (Solve Keyword Arguments)","title":"Common Solver Options (Solve Keyword Arguments)","text":"solve(prob::IntegralProblem, alg::SciMLBase.AbstractIntegralAlgorithm)","category":"page"},{"location":"basics/solve/#CommonSolve.solve-Tuple{IntegralProblem, SciMLBase.AbstractIntegralAlgorithm}","page":"Common Solver Options (Solve Keyword Arguments)","title":"CommonSolve.solve","text":"solve(prob::IntegralProblem, alg::SciMLBase.AbstractIntegralAlgorithm; kwargs...)\n\nKeyword Arguments\n\nThe arguments to solve are common across all of the quadrature methods. These common arguments are:\n\nmaxiters (the maximum number of iterations)\nabstol (absolute tolerance in changes of the objective value)\nreltol (relative tolerance in changes of the objective value)\n\n\n\n\n\n","category":"method"},{"location":"basics/solve/","page":"Common Solver Options (Solve Keyword Arguments)","title":"Common Solver Options (Solve Keyword Arguments)","text":"Additionally, the extra keyword arguments are splatted to the library calls, so see the documentation of the integrator library for all the extra details. These extra keyword arguments are not guaranteed to act uniformly.","category":"page"},{"location":"tutorials/differentiating_integrals/#Differentiating-Integrals","page":"Differentiating Integrals","title":"Differentiating Integrals","text":"","category":"section"},{"location":"tutorials/differentiating_integrals/","page":"Differentiating Integrals","title":"Differentiating Integrals","text":"Integrals.jl is a fully differentiable quadrature library. Thus, it adds the ability to perform automatic differentiation over any of the libraries that it calls. It integrates with ForwardDiff.jl for forward-mode automatic differentiation and Zygote.jl for reverse-mode automatic differentiation. For example:","category":"page"},{"location":"tutorials/differentiating_integrals/","page":"Differentiating Integrals","title":"Differentiating Integrals","text":"using Integrals, ForwardDiff, FiniteDiff, Zygote, Cuba\nf(x, p) = sum(sin.(x .* p))\ndomain = (ones(2), 3ones(2)) # (lb, ub)\np = ones(2)\n\nfunction testf(p)\n prob = IntegralProblem(f, domain, p)\n sin(solve(prob, CubaCuhre(), reltol = 1e-6, abstol = 1e-6)[1])\nend\ntestf(p)\ndp1 = Zygote.gradient(testf, p)\ndp2 = FiniteDiff.finite_difference_gradient(testf, p)\ndp3 = ForwardDiff.gradient(testf, p)\ndp1[1] ≈ dp2 ≈ dp3","category":"page"},{"location":"basics/FAQ/#Frequently-Asked-Questions","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"","category":"section"},{"location":"basics/FAQ/#How-should-I-use-the-in-place-interface?","page":"Frequently Asked Questions","title":"How should I use the in-place interface?","text":"","category":"section"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"The in-place interface allows evaluating vector-valued integrands without allocating an output array. This can be beneficial for reducing allocations when integrating many functions simultaneously or to make use of existing in-place code. However, note that not all algorithms use in-place operations under the hood, i.e. HCubatureJL(), and may still allocate.","category":"page"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"You can construct an IntegralFunction(f, prototype), where f is of the form f(y, u, p) where prototype is of the desired type and shape of y.","category":"page"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"For small array outputs of a known size, consider using StaticArrays.jl for the return value of your integrand.","category":"page"},{"location":"basics/FAQ/#How-should-I-use-the-batch-interface?","page":"Frequently Asked Questions","title":"How should I use the batch interface?","text":"","category":"section"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"The batch interface allows evaluating one (or more) integrals simultaneously at different points, which maximizes the parallelism for a given algorithm.","category":"page"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"You can construct an out-of-place BatchIntegralFunction(bf) where bf is of the form bf(u, p) = stack(x -> f(x, p), eachslice(u; dims=ndims(u))), where f is the (unbatched) integrand.","category":"page"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"You can construct an in-place BatchIntegralFunction(bf, prototype), where bf is of the form bf(y, u, p) = foreach((y,x) -> f(y,x,p), eachslice(y, dims=ndims(y)), eachslice(x, dims=ndims(x))).","category":"page"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"Note that not all algorithms use in-place batched operations under the hood, i.e. QuadGKJL().","category":"page"},{"location":"basics/FAQ/#What-should-I-do-if-my-solution-is-not-converged?","page":"Frequently Asked Questions","title":"What should I do if my solution is not converged?","text":"","category":"section"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"Certain algorithms, such as QuadratureRule used a fixed number of points to calculate an integral and cannot provide an error estimate. In this case, you have to increase the number of points and check the convergence yourself, which will depend on the accuracy of the rule you choose.","category":"page"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"For badly-behaved integrands, such as (nearly) singular and highly oscillatory functions, most algorithms will fail to converge and either throw an error or silently return the incorrect result. In some cases Integrals.jl can provide an error code when things go wrong, but otherwise you can always check if the error estimate for the integral is less than the requested tolerance, e.g. sol.resid < max(abstol, reltol*norm(sol.u)). Sometimes using a larger tolerance or higher precision arithmetic may help.","category":"page"},{"location":"basics/FAQ/#How-can-I-integrate-arbitrarily-spaced-data?","page":"Frequently Asked Questions","title":"How can I integrate arbitrarily-spaced data?","text":"","category":"section"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"See SampledIntegralProblem.","category":"page"},{"location":"basics/FAQ/#How-can-I-integrate-on-arbitrary-geometries?","page":"Frequently Asked Questions","title":"How can I integrate on arbitrary geometries?","text":"","category":"section"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"You can't, since Integrals.jl currently supports integration on hypercubes because that is what lower-level packages implement.","category":"page"},{"location":"basics/FAQ/#I-don't-see-algorithm-X-or-quadrature-scheme-Y-?","page":"Frequently Asked Questions","title":"I don't see algorithm X or quadrature scheme Y ?","text":"","category":"section"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"Fixed quadrature rules from other packages can be used with QuadratureRule. Otherwise, feel free to open an issue or pull request.","category":"page"},{"location":"basics/FAQ/#Can-I-take-derivatives-with-respect-to-the-limits-of-integration?","page":"Frequently Asked Questions","title":"Can I take derivatives with respect to the limits of integration?","text":"","category":"section"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"Currently this is not implemented.","category":"page"},{"location":"#Integrals.jl:-Unified-Integral-Approximation-Interface","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"","category":"section"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"Integrals.jl is a unified interface for the numerical approximation of integrals (quadrature) in Julia. It interfaces with other packages of the Julia ecosystem to make it easy to test alternative solver packages and pass small types to control algorithm swapping.","category":"page"},{"location":"#Installation","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Installation","text":"","category":"section"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"To install Integrals.jl, use the Julia package manager:","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"using Pkg\nPkg.add(\"Integrals\")","category":"page"},{"location":"#Contributing","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Contributing","text":"","category":"section"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"Please refer to the SciML ColPrac: Contributor's Guide on Collaborative Practices for Community Packages for guidance on PRs, issues, and other matters relating to contributing to SciML.\nSee the SciML Style Guide for common coding practices and other style decisions.\nThere are a few community forums:\nThe #diffeq-bridged and #sciml-bridged channels in the Julia Slack\nThe #diffeq-bridged and #sciml-bridged channels in the Julia Zulip\nOn the Julia Discourse forums\nSee also SciML Community page","category":"page"},{"location":"#Reproducibility","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Reproducibility","text":"","category":"section"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"
The documentation of this SciML package was built using these direct dependencies,","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"using Pkg # hide\nPkg.status() # hide","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"
","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"
and using this machine and Julia version.","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"using InteractiveUtils # hide\nversioninfo() # hide","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"
","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"
A more complete overview of all dependencies and their versions is also provided.","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"using Pkg # hide\nPkg.status(; mode = PKGMODE_MANIFEST) # hide","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"
","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"using TOML\nusing Markdown\nversion = TOML.parse(read(\"../../Project.toml\", String))[\"version\"]\nname = TOML.parse(read(\"../../Project.toml\", String))[\"name\"]\nlink_manifest = \"https://github.com/SciML/\" * name * \".jl/tree/gh-pages/v\" * version *\n \"/assets/Manifest.toml\"\nlink_project = \"https://github.com/SciML/\" * name * \".jl/tree/gh-pages/v\" * version *\n \"/assets/Project.toml\"\nMarkdown.parse(\"\"\"You can also download the\n[manifest]($link_manifest)\nfile and the\n[project]($link_project)\nfile.\n\"\"\")","category":"page"},{"location":"basics/IntegralFunction/#func","page":"Integral Functions","title":"Integral Functions","text":"","category":"section"},{"location":"basics/IntegralFunction/","page":"Integral Functions","title":"Integral Functions","text":"SciMLBase.IntegralFunction\nSciMLBase.BatchIntegralFunction","category":"page"},{"location":"basics/IntegralFunction/#SciMLBase.IntegralFunction","page":"Integral Functions","title":"SciMLBase.IntegralFunction","text":"IntegralFunction{iip,specialize,F,T} <: AbstractIntegralFunction{iip}\n\nA representation of an integrand f defined by:\n\nf(u p)\n\nFor an in-place form of f see the iip section below for details on in-place or out-of-place handling.\n\nIntegralFunction{iip,specialize}(f, [integrand_prototype])\n\nNote that only f is required, and in the case of inplace integrands a mutable container integrand_prototype to store the result of the integrand. If integrand_prototype is present, f is interpreted as in-place, and otherwise f is assumed to be out-of-place.\n\niip: In-Place vs Out-Of-Place\n\nOut-of-place functions must be of the form y = f(u p) and in-place functions of the form f(y u p), where y is a number or array containing the output. Since f is allowed to return any type (e.g. real or complex numbers or arrays), in-place functions must provide a container integrand_prototype that is of the right type and size for the variable y, and the result is written to this container in-place. When in-place forms are used, in-place array operations, i.e. broadcasting, may be used by algorithms to reduce allocations. If integrand_prototype is not provided, f is assumed to be out-of-place.\n\nspecialize\n\nThis field is currently unused\n\nFields\n\nThe fields of the IntegralFunction type directly match the names of the inputs.\n\n\n\n\n\n","category":"type"},{"location":"basics/IntegralFunction/#SciMLBase.BatchIntegralFunction","page":"Integral Functions","title":"SciMLBase.BatchIntegralFunction","text":"BatchIntegralFunction{iip,specialize,F,T} <: AbstractIntegralFunction{iip}\n\nA batched representation of an (non-batched) integrand f(u, p) that can be evaluated at multiple points simultaneously using threads, the gpu, or distributed memory defined by:\n\nby = bf(bu p)\n\nHere we prefix variables with b to indicate they are batched variables, which implies that they are arrays whose last dimension is reserved for batching different evaluation points, or function values, and may be of a variable length. bu is an array whose elements correspond to distinct evaluation points to f, and bf is a function to evaluate f 'point-wise' so that f(bu[..., i], p) == bf(bu, p)[..., i]. For example, a simple batching implementation of a scalar, univariate function is via broadcasting: bf(bu, p) = f.(bu, Ref(p)), although this interface exists in order to allow user parallelization. In general, the integration algorithm is allowed to vary the number of evaluation points between subsequent calls to bf.\n\nFor an in-place form of bf see the iip section below for details on in-place or out-of-place handling.\n\nBatchIntegralFunction{iip,specialize}(bf, [integrand_prototype];\n max_batch=typemax(Int))\n\nNote that only bf is required, and in the case of inplace integrands a mutable container integrand_prototype to store a batch of integrand evaluations, with a last \"batching\" dimension.\n\nThe keyword max_batch is used to set a soft limit on the number of points to batch at the same time so that memory usage is controlled.\n\nIf integrand_prototype is present, bf is interpreted as in-place, and otherwise bf is assumed to be out-of-place.\n\niip: In-Place vs Out-Of-Place\n\nOut-of-place functions must be of the form by = bf(bu, p) and in-place functions of the form bf(by, bu, p) where by is a batch array containing the output. Since the algorithm may vary the number of points to batch, the batching dimension can be of any length, including zero, and since bf is allowed to return arrays of any type (e.g. real or complex) or size, in-place functions must provide a container integrand_prototype of the desired type and size for by. If integrand_prototype is not provided, bf is assumed to be out-of-place.\n\nIn the out-of-place case, we require f(bu[..., i], p) == bf(bu, p)[..., i], and certain algorithms, such as those implemented in C, may infer the type or shape of by by calling bf with an empty array of input points, i.e. bu with size(bu)[end] == 0. Then it is expected for the resulting by to have the same type and size(by)[begin:end-1] for all subsequent calls.\n\nWhen the in-place form is used, we require f(by[..., i], bu[..., i], p) == bf(by, bu, p)[..., i] and size(by)[begin:end-1] == size(integrand_prototype)[begin:end-1]. The algorithm should always pass the integrand by arrays that are similar to integrand_prototype, and may use views and in-place array operations to reduce allocations.\n\nspecialize\n\nThis field is currently unused\n\nFields\n\nThe fields of the BatchIntegralFunction type are f, corresponding to bf above, and integrand_prototype.\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#solvers","page":"Integral Solver Algorithms","title":"Integral Solver Algorithms","text":"","category":"section"},{"location":"solvers/IntegralSolvers/","page":"Integral Solver Algorithms","title":"Integral Solver Algorithms","text":"The following algorithms are available:","category":"page"},{"location":"solvers/IntegralSolvers/","page":"Integral Solver Algorithms","title":"Integral Solver Algorithms","text":"QuadGKJL: Uses QuadGK.jl, which supports one-dimensional integration of scalar and array-valued integrands with in-place or batched forms. Integrands that are both in-place and batched are implemented in the wrapper but are not supported under the hood.\nHCubatureJL: Uses HCubature.jl, which supports scalar and array-valued integrands and works best in low dimensions, e.g. ≤ 8. In-place integrands are implemented in the wrapper but are not supported under the hood. Batching is not supported.\nVEGAS: Uses MonteCarloIntegration.jl, which requires scalar, Float64-valued integrands and works in any number of dimensions.\nVEGASMC: Uses MCIntegration.jl. Requires using MCIntegration. Doesn't support batching.\nCubatureJLh: h-Cubature from Cubature.jl. Requires using Cubature.\nCubatureJLp: p-Cubature from Cubature.jl. Requires using Cubature.\nCubaVegas: Vegas from Cuba.jl. Requires using Cuba.\nCubaSUAVE: SUAVE from Cuba.jl. Requires using Cuba.\nCubaDivonne: Divonne from Cuba.jl. Requires using Cuba. Works only for >1-dimensional integrations.\nCubaCuhre: Cuhre from Cuba.jl. Requires using Cuba. Works only for >1-dimensional integrations.\nGaussLegendre: Performs fixed-order Gauss-Legendre quadrature. Requires using FastGaussQuadrature.\nQuadratureRule: Accepts a user-defined function that returns nodes and weights.\nArblibJL: real- and complex-valued univariate integration of holomorphic and meromorphic functions from Arblib.jl. Requires using Arblib.","category":"page"},{"location":"solvers/IntegralSolvers/","page":"Integral Solver Algorithms","title":"Integral Solver Algorithms","text":"QuadGKJL\nHCubatureJL\nCubatureJLp\nCubatureJLh\nVEGAS\nVEGASMC\nCubaVegas\nCubaSUAVE\nCubaDivonne\nCubaCuhre\nGaussLegendre\nQuadratureRule\nArblibJL","category":"page"},{"location":"solvers/IntegralSolvers/#Integrals.QuadGKJL","page":"Integral Solver Algorithms","title":"Integrals.QuadGKJL","text":"QuadGKJL(; order = 7, norm=norm)\n\nOne-dimensional Gauss-Kronrod integration from QuadGK.jl. This method also takes the optional arguments order and norm. Which are the order of the integration rule and the norm for calculating the error, respectively\n\nReferences\n\n@article{laurie1997calculation, title={Calculation of Gauss-Kronrod quadrature rules}, author={Laurie, Dirk}, journal={Mathematics of Computation}, volume={66}, number={219}, pages={1133–1145}, year={1997} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.HCubatureJL","page":"Integral Solver Algorithms","title":"Integrals.HCubatureJL","text":"HCubatureJL(; norm=norm, initdiv=1)\n\nMultidimensional \"h-adaptive\" integration from HCubature.jl. This method also takes the optional arguments initdiv and norm. Which are the initial number of segments each dimension of the integration domain is divided into, and the norm for calculating the error, respectively.\n\nReferences\n\n@article{genz1980remarks, title={Remarks on algorithm 006: An adaptive algorithm for numerical integration over an N-dimensional rectangular region}, author={Genz, Alan C and Malik, Aftab Ahmad}, journal={Journal of Computational and Applied mathematics}, volume={6}, number={4}, pages={295–302}, year={1980}, publisher={Elsevier} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.CubatureJLp","page":"Integral Solver Algorithms","title":"Integrals.CubatureJLp","text":"CubatureJLp(; error_norm=Cubature.INDIVIDUAL)\n\nMultidimensional p-adaptive integration from Cubature.jl. This method is based on repeatedly doubling the degree of the cubature rules, until convergence is achieved. The used cubature rule is a tensor product of Clenshaw–Curtis quadrature rules. error_norm specifies the convergence criterion for vector valued integrands. Defaults to Cubature.INDIVIDUAL, other options are Cubature.PAIRED, Cubature.L1, Cubature.L2, or Cubature.LINF.\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.CubatureJLh","page":"Integral Solver Algorithms","title":"Integrals.CubatureJLh","text":"CubatureJLh(; error_norm=Cubature.INDIVIDUAL)\n\nMultidimensional h-adaptive integration from Cubature.jl. error_norm specifies the convergence criterion for vector valued integrands. Defaults to Cubature.INDIVIDUAL, other options are Cubature.PAIRED, Cubature.L1, Cubature.L2, or Cubature.LINF.\n\nReferences\n\n@article{genz1980remarks, title={Remarks on algorithm 006: An adaptive algorithm for numerical integration over an N-dimensional rectangular region}, author={Genz, Alan C and Malik, Aftab Ahmad}, journal={Journal of Computational and Applied mathematics}, volume={6}, number={4}, pages={295–302}, year={1980}, publisher={Elsevier} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.VEGAS","page":"Integral Solver Algorithms","title":"Integrals.VEGAS","text":"VEGAS(; nbins = 100, ncalls = 1000, debug=false)\n\nMultidimensional adaptive Monte Carlo integration from MonteCarloIntegration.jl. Importance sampling is used to reduce variance. This method also takes three optional arguments nbins, ncalls and debug which are the initial number of bins each dimension of the integration domain is divided into, the number of function calls per iteration of the algorithm, and whether debug info should be printed, respectively.\n\nLimitations\n\nThis algorithm can only integrate Float64-valued functions\n\nReferences\n\n@article{lepage1978new, title={A new algorithm for adaptive multidimensional integration}, author={Lepage, G Peter}, journal={Journal of Computational Physics}, volume={27}, number={2}, pages={192–203}, year={1978}, publisher={Elsevier} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.VEGASMC","page":"Integral Solver Algorithms","title":"Integrals.VEGASMC","text":"VEGASMC(; kws...)\n\nMarkov-chain based Vegas algorithm from MCIntegration.jl\n\nRefer to MCIntegration.integrate for documentation on the keywords, which are passed directly to the solver with a set of defaults that works for conforming integrands.\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.CubaVegas","page":"Integral Solver Algorithms","title":"Integrals.CubaVegas","text":"CubaVegas()\n\nMultidimensional adaptive Monte Carlo integration from Cuba.jl. Importance sampling is used to reduce variance.\n\nReferences\n\n@article{lepage1978new, title={A new algorithm for adaptive multidimensional integration}, author={Lepage, G Peter}, journal={Journal of Computational Physics}, volume={27}, number={2}, pages={192–203}, year={1978}, publisher={Elsevier} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.CubaSUAVE","page":"Integral Solver Algorithms","title":"Integrals.CubaSUAVE","text":"CubaSUAVE()\n\nMultidimensional adaptive Monte Carlo integration from Cuba.jl. Suave stands for subregion-adaptive VEGAS. Importance sampling and subdivision are thus used to reduce variance.\n\nReferences\n\n@article{hahn2005cuba, title={Cuba—a library for multidimensional numerical integration}, author={Hahn, Thomas}, journal={Computer Physics Communications}, volume={168}, number={2}, pages={78–95}, year={2005}, publisher={Elsevier} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.CubaDivonne","page":"Integral Solver Algorithms","title":"Integrals.CubaDivonne","text":"CubaDivonne()\n\nMultidimensional adaptive Monte Carlo integration from Cuba.jl. Stratified sampling is used to reduce variance.\n\nReferences\n\n@article{friedman1981nested, title={A nested partitioning procedure for numerical multiple integration}, author={Friedman, Jerome H and Wright, Margaret H}, journal={ACM Transactions on Mathematical Software (TOMS)}, volume={7}, number={1}, pages={76–92}, year={1981}, publisher={ACM New York, NY, USA} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.CubaCuhre","page":"Integral Solver Algorithms","title":"Integrals.CubaCuhre","text":"CubaCuhre()\n\nMultidimensional h-adaptive integration from Cuba.jl.\n\nReferences\n\n@article{berntsen1991adaptive, title={An adaptive algorithm for the approximate calculation of multiple integrals}, author={Berntsen, Jarle and Espelid, Terje O and Genz, Alan}, journal={ACM Transactions on Mathematical Software (TOMS)}, volume={17}, number={4}, pages={437–451}, year={1991}, publisher={ACM New York, NY, USA} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.GaussLegendre","page":"Integral Solver Algorithms","title":"Integrals.GaussLegendre","text":"GaussLegendre{C, N, W}\n\nStruct for evaluating an integral via (composite) Gauss-Legendre quadrature. The field C will be true if subintervals > 1, and false otherwise.\n\nThe fields nodes::N and weights::W are defined by nodes, weights = gausslegendre(n) for a given number of nodes n.\n\nThe field subintervals::Int64 = 1 (with default value 1) defines the number of intervals to partition the original interval of integration [a, b] into, splitting it into [xⱼ, xⱼ₊₁] for j = 1,…,subintervals, where xⱼ = a + (j-1)h and h = (b-a)/subintervals. Gauss-Legendre quadrature is then applied on each subinterval. For example, if [a, b] = [-1, 1] and subintervals = 2, then Gauss-Legendre quadrature will be applied separately on [-1, 0] and [0, 1], summing the two results.\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.QuadratureRule","page":"Integral Solver Algorithms","title":"Integrals.QuadratureRule","text":"QuadratureRule(q; n=250)\n\nAlgorithm to construct and evaluate a quadrature rule q of n points computed from the inputs as x, w = q(n). It assumes the nodes and weights are for the standard interval [-1, 1]^d in d dimensions, and rescales the nodes to the specific hypercube being solved. The nodes x may be scalars in 1d or vectors in arbitrary dimensions, and the weights w must be scalar. The algorithm computes the quadrature rule sum(w .* f.(x)) and the caller must check that the result is converged with respect to n.\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.ArblibJL","page":"Integral Solver Algorithms","title":"Integrals.ArblibJL","text":"ArblibJL(; check_analytic=false, take_prec=false, warn_on_no_convergence=false, opts=C_NULL)\n\nOne-dimensional adaptive Gauss-Legendre integration using rigorous error bounds and precision ball arithmetic. Generally this assumes the integrand is holomorphic or meromorphic, which is the user's responsibility to verify. The result of the integral is not guaranteed to satisfy the requested tolerances, however the result is guaranteed to be within the error estimate.\n\nArblib.jl only supports integration of univariate real- and complex-valued functions with both inplace and out-of-place forms. See their documentation for additional details the algorithm arguments and on implementing high-precision integrands. Additionally, the error estimate is included in the return value of the integral, representing a ball.\n\n\n\n\n\n","category":"type"}] +[{"location":"basics/SampledIntegralProblem/#Integrating-pre-sampled-data","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"","category":"section"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"In some cases, instead of a function that acts as integrand, one only possesses a list of data points y at a set of sampling locations x, that must be integrated. This package contains functionality for doing that.","category":"page"},{"location":"basics/SampledIntegralProblem/#Example","page":"Integrating pre-sampled data","title":"Example","text":"","category":"section"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"Say, by some means we have generated a dataset x and y:","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"using Integrals # hide\nf = x -> x^2\nx = range(0, 1, length = 20)\ny = f.(x)","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"Now, we can integrate this data set as follows:","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"problem = SampledIntegralProblem(y, x)\nmethod = TrapezoidalRule()\nsolve(problem, method)","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"The exact answer is of course $ 1/3 $.","category":"page"},{"location":"basics/SampledIntegralProblem/#Details","page":"Integrating pre-sampled data","title":"Details","text":"","category":"section"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"SciMLBase.SampledIntegralProblem\nsolve(::SampledIntegralProblem, ::SciMLBase.AbstractIntegralAlgorithm)","category":"page"},{"location":"basics/SampledIntegralProblem/#SciMLBase.SampledIntegralProblem","page":"Integrating pre-sampled data","title":"SciMLBase.SampledIntegralProblem","text":"Defines a integral problem over pre-sampled data. Documentation Page: https://docs.sciml.ai/Integrals/stable/\n\nMathematical Specification of a data Integral Problem\n\nSampled integral problems are defined as:\n\nsum_i w_i y_i\n\nwhere y_i are sampled values of the integrand, and w_i are weights assigned by a quadrature rule, which depend on sampling points x.\n\nProblem Type\n\nConstructors\n\nSampledIntegralProblem(y::AbstractArray, x::AbstractVector; dim=ndims(y), kwargs...)\n\ny: The sampled integrand, must be a subtype of AbstractArray. It is assumed that the values of y along dimension dim correspond to the integrand evaluated at sampling points x\nx: Sampling points, must be a subtype of AbstractVector.\ndim: Dimension along which to integrate. Defaults to the last dimension of y.\nkwargs: Keyword arguments copied to the solvers.\n\nFields\n\nThe fields match the names of the constructor arguments.\n\n\n\n\n\n","category":"type"},{"location":"basics/SampledIntegralProblem/#CommonSolve.solve-Tuple{SampledIntegralProblem, SciMLBase.AbstractIntegralAlgorithm}","page":"Integrating pre-sampled data","title":"CommonSolve.solve","text":"solve(prob::SampledIntegralProblem, alg::SciMLBase.AbstractIntegralAlgorithm; kwargs...)\n\nKeyword Arguments\n\nThere are no keyword arguments used to solve SampledIntegralProblems\n\n\n\n\n\n","category":"method"},{"location":"basics/SampledIntegralProblem/#Non-equidistant-grids","page":"Integrating pre-sampled data","title":"Non-equidistant grids","text":"","category":"section"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"If the sampling points x are provided as an AbstractRange (constructed with the range function for example), faster methods are used that take advantage of the fact that the points are equidistantly spaced. Otherwise, general methods are used for non-uniform grids.","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"Example:","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"using Integrals # hide\nf = x -> x^7\nx = [0.0; sort(rand(1000)); 1.0]\ny = f.(x)\nproblem = SampledIntegralProblem(y, x)\nmethod = TrapezoidalRule()\nsolve(problem, method)","category":"page"},{"location":"basics/SampledIntegralProblem/#Evaluating-multiple-integrals-at-once","page":"Integrating pre-sampled data","title":"Evaluating multiple integrals at once","text":"","category":"section"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"If the provided data set y is a multidimensional array, the integrals are evaluated across only one of its axes. For performance reasons, the last axis of the array y is chosen by default, but this can be modified with the dim keyword argument to the problem definition.","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"using Integrals # hide\nf1 = x -> x^2\nf2 = x -> x^3\nf3 = x -> x^4\nx = range(0, 1, length = 20)\ny = [f1.(x) f2.(x) f3.(x)]\nproblem = SampledIntegralProblem(y, x; dim = 1)\nmethod = SimpsonsRule()\nsolve(problem, method)","category":"page"},{"location":"basics/SampledIntegralProblem/#Supported-methods","page":"Integrating pre-sampled data","title":"Supported methods","text":"","category":"section"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"Right now, only the TrapezoidalRule and SimpsonsRule are supported.","category":"page"},{"location":"basics/SampledIntegralProblem/","page":"Integrating pre-sampled data","title":"Integrating pre-sampled data","text":"TrapezoidalRule\nSimpsonsRule","category":"page"},{"location":"basics/SampledIntegralProblem/#Integrals.TrapezoidalRule","page":"Integrating pre-sampled data","title":"Integrals.TrapezoidalRule","text":"TrapezoidalRule\n\nStruct for evaluating an integral via the trapezoidal rule.\n\nExample with sampled data:\n\nusing Integrals\nf = x -> x^2\nx = range(0, 1, length=20)\ny = f.(x)\nproblem = SampledIntegralProblem(y, x)\nmethod = TrapezoidalRule()\nsolve(problem, method)\n\n\n\n\n\n","category":"type"},{"location":"basics/SampledIntegralProblem/#Integrals.SimpsonsRule","page":"Integrating pre-sampled data","title":"Integrals.SimpsonsRule","text":"SimpsonsRule\n\nStruct for evaluating an integral via the Simpson's composite 1/3-3/8 rule over AbstractRanges (evenly spaced points) and Simpson's composite 1/3 rule for non-equidistant grids.\n\nExample with equidistant data:\n\nusing Integrals\nf = x -> x^2\nx = range(0, 1, length=20)\ny = f.(x)\nproblem = SampledIntegralProblem(y, x)\nmethod = SimpsonsRule()\nsolve(problem, method)\n\n\n\n\n\n","category":"type"},{"location":"tutorials/caching_interface/#Integrals-with-Caching-Interface","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"","category":"section"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"Often, integral solvers allocate memory or reuse quadrature rules for solving different problems. For example, if one is going to solve the same integral for several parameters","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"using Integrals\n\nprob = IntegralProblem((x, p) -> sin(x * p), (0, 1), 14.0)\nalg = QuadGKJL()\n\nsolve(prob, alg)\n\nprob = remake(prob, p = 15.0)\nsolve(prob, alg)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"then it would be more efficient to allocate the heap used by quadgk across several calls, shown below by directly calling the library","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"using QuadGK\nsegbuf = QuadGK.alloc_segbuf()\nquadgk(x -> sin(14x), 0, 1, segbuf = segbuf)\nquadgk(x -> sin(15x), 0, 1, segbuf = segbuf)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"Integrals.jl's caching interface automates this process to reuse resources if an algorithm supports it and if the necessary types to build the cache can be inferred from prob. To do this with Integrals.jl, you simply init a cache, solve!, replace p, and solve again. This uses the SciML init interface","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"using Integrals\n\nprob = IntegralProblem((x, p) -> sin(x * p), (0, 1), 14.0)\nalg = QuadGKJL()\n\ncache = init(prob, alg)\nsol1 = solve!(cache)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"cache.p = 15.0\nsol2 = solve!(cache)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"The caching interface is intended for updating p and domain. Note that the types of these variables is not allowed to change. If it is necessary to change the integrand f instead of defining a new IntegralProblem, consider using FunctionWrappers.jl.","category":"page"},{"location":"tutorials/caching_interface/#Caching-for-sampled-integral-problems","page":"Integrals with Caching Interface","title":"Caching for sampled integral problems","text":"","category":"section"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"For sampled integral problems, it is possible to cache the weights and reuse them for multiple data sets.","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"using Integrals\n\nx = 0.0:0.1:1.0\ny = sin.(x)\n\nprob = SampledIntegralProblem(y, x)\nalg = TrapezoidalRule()\n\ncache = init(prob, alg)\nsol1 = solve!(cache)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"cache.y = cos.(x) # use .= to update in-place\nsol2 = solve!(cache)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"If the grid is modified, the weights are recomputed.","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"cache.x = 0.0:0.2:2.0\ncache.y = sin.(cache.x)\nsol3 = solve!(cache)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"For multi-dimensional datasets, the integration dimension can also be changed","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"using Integrals\n\nx = 0.0:0.1:1.0\ny = sin.(x) .* cos.(x')\n\nprob = SampledIntegralProblem(y, x)\nalg = TrapezoidalRule()\n\ncache = init(prob, alg)\nsol1 = solve!(cache)","category":"page"},{"location":"tutorials/caching_interface/","page":"Integrals with Caching Interface","title":"Integrals with Caching Interface","text":"cache.dim = 1\nsol2 = solve!(cache)","category":"page"},{"location":"tutorials/numerical_integrals/#Numerically-Solving-Integrals","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"","category":"section"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"For basic multidimensional quadrature, we can construct and solve a IntegralProblem. The integral we want to evaluate is:","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"int_1^3int_1^3int_1^3 sum_1^3 sin(u_i) du_1du_2du_3","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"We can numerically approximate this integral:","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"using Integrals\nf(u, p) = sum(sin.(u))\ndomain = (ones(3), 3ones(3)) # (lb, ub)\nprob = IntegralProblem(f, domain)\nsol = solve(prob, HCubatureJL(); reltol = 1e-3, abstol = 1e-3)\nsol.u","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"where the first argument of IntegralProblem is the integrand, the second argument is the lower bound, and the third argument is the upper bound. p are the parameters of the integrand. In this case, there are no parameters, but still f must be defined as f(x,p) and not f(x). For an example with parameters, see the next tutorial. The first argument of solve is the problem we are solving, the second is an algorithm to solve the problem with. Then there are keywords which provides details how the algorithm should work, in this case tolerances how precise the numerical approximation should be.","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"We can also evaluate multiple integrals at once. We could create two IntegralProblems for this, but that is wasteful if the integrands share a lot of computation. For example, we also want to evaluate:","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"int_1^3int_1^3int_1^3 sum_1^3 cos(u_i) du_1du_2du_3","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"using Integrals\nf(u, p) = [sum(sin.(u)), sum(cos.(u))]\ndomain = (ones(3), 3ones(3)) # (lb, ub)\nprob = IntegralProblem(f, domain)\nsol = solve(prob, HCubatureJL(); reltol = 1e-3, abstol = 1e-3)\nsol.u","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"Another way to think about this is that the integrand is now a vector valued function. In general, we should be able to integrate any type that is in a vector space and supports addition and scalar multiplication, although Integrals.jl allows scalars and arrays. In the above example, the integrand was defined out-of-position. This means that a new output vector is created every time the function f is called. If we do not want these allocations, we can also define f in-position.","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"using Integrals, Cubature\nfunction f(y, u, p)\n y[1] = sum(sin.(u))\n y[2] = sum(cos.(u))\nend\nprototype = zeros(2)\ndomain = (ones(3), 3ones(3)) # (lb, ub)\nprob = IntegralProblem(IntegralFunction(f, prototype), domain)\nsol = solve(prob, CubatureJLh(); reltol = 1e-3, abstol = 1e-3)\nsol.u","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"where y is a cache to store the evaluation of the integrand and prototype is an instance of y with the desired type and shape. We needed to change the algorithm to CubatureJLh() because HCubatureJL() does not support in-position under the hood. f evaluates the integrand at a certain point, but most adaptive quadrature algorithms need to evaluate the integrand at multiple points in each step of the algorithm. We would thus often like to parallelize the computation. The batch interface allows us to compute multiple points at once. For example, here we do allocation-free multithreading with Cubature.jl:","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"using Integrals, Cubature, Base.Threads\nfunction f(y, u, p)\n Threads.@threads for i in 1:size(u, 2)\n y[1, i] = sum(sin, @view(u[:, i]))\n y[2, i] = sum(cos, @view(u[:, i]))\n end\nend\nprototype = zeros(2, 0)\ndomain = (ones(3), 3ones(3)) # (lb, ub)\nprob = IntegralProblem(BatchIntegralFunction(f, prototype), domain)\nsol = solve(prob, CubatureJLh(); reltol = 1e-3, abstol = 1e-3)\nsol.u","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"Both u and y changed from vectors to matrices, where each column is respectively a point the integrand is evaluated at or the evaluation of the integrand at the corresponding point. The prototype now has an extra dimension for batching that can be of size zero. Try to create yourself an out-of-position version of the above problem. For the full details of the batching interface, see the problem page.","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"If we would like to compare the results against Cuba.jl's Cuhre method, then the change is a one-argument change:","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"using Integrals\nusing Cuba\nf(u, p) = sum(sin.(u))\ndomain = (ones(3), 3ones(3)) # (lb, ub)\nprob = IntegralProblem(f, domain)\nsol = solve(prob, CubaCuhre(); reltol = 1e-3, abstol = 1e-3)\nsol.u","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"However, Cuhre does not support vector valued integrands. The solvers page gives an overview of which arguments each algorithm can handle.","category":"page"},{"location":"tutorials/numerical_integrals/#One-dimensional-integrals","page":"Numerically Solving Integrals","title":"One-dimensional integrals","text":"","category":"section"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"Integrals.jl also has specific solvers for integrals in a single dimension, such as QuadGKJL. For example, we can create our own sine function by integrating the cosine function from 0 to x.","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"using Integrals\nmy_sin(x) = solve(IntegralProblem((x, p) -> cos(x), (0.0, x)), QuadGKJL()).u\nx = 0:0.1:(2 * pi)\n@. my_sin(x) ≈ sin(x)","category":"page"},{"location":"tutorials/numerical_integrals/#Infinity-handling","page":"Numerically Solving Integrals","title":"Infinity handling","text":"","category":"section"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"Integrals.jl can also handle infinite integration bounds. For infinite upper bounds u is substituted with a+fract1-t, and the integral is thus transformed to:","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"int_a^infty f(u)du = int_0^1 fleft(a+fract1-tright)frac1(1-t)^2dt","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"Integrals with an infinite lower bound are handled in the same way. If both upper and lower bound are infinite, u is substituted with fract1-t^2,","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"int_-infty^infty f(u)du = int_-1^1 fleft(fract1-t^2right)frac1+t^2(1-t^2)^2dt","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"For multidimensional integrals, each variable with infinite bounds is substituted the same way. The details of the math behind these transforms can be found here.","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"As an example, let us integrate the standard bivariate normal probability distribution over the area above the horizontal axis, which should be equal to 05.","category":"page"},{"location":"tutorials/numerical_integrals/","page":"Numerically Solving Integrals","title":"Numerically Solving Integrals","text":"using Distributions\nusing Integrals\ndist = MvNormal(ones(2))\nf = (x, p) -> pdf(dist, x)\ndomain = ([-Inf, 0.0], [Inf, Inf]) # (lb, ub)\nprob = IntegralProblem(f, domain)\nsolve(prob, HCubatureJL(), reltol = 1e-3, abstol = 1e-3)","category":"page"},{"location":"basics/IntegralProblem/#prob","page":"Integral Problems","title":"Integral Problems","text":"","category":"section"},{"location":"basics/IntegralProblem/","page":"Integral Problems","title":"Integral Problems","text":"SciMLBase.IntegralProblem","category":"page"},{"location":"basics/IntegralProblem/#SciMLBase.IntegralProblem","page":"Integral Problems","title":"SciMLBase.IntegralProblem","text":"Defines an integral problem. Documentation Page: https://docs.sciml.ai/Integrals/stable/\n\nMathematical Specification of an Integral Problem\n\nIntegral problems are multi-dimensional integrals defined as:\n\nint_lb^ub f(up) du\n\nwhere p are parameters. u is a Number or AbstractVector whose geometry matches the space being integrated. This space is bounded by the lowerbound lb and upperbound ub, which are Numbers or AbstractVectors with the same geometry as u.\n\nProblem Type\n\nConstructors\n\nIntegralProblem(f::AbstractIntegralFunction,domain,p=NullParameters(); kwargs...)\nIntegralProblem(f::AbstractIntegralFunction,lb,ub,p=NullParameters(); kwargs...)\nIntegralProblem(f,domain,p=NullParameters(); nout=nothing, batch=nothing, kwargs...)\nIntegralProblem(f,lb,ub,p=NullParameters(); nout=nothing, batch=nothing, kwargs...)\n\nf: the integrand, callable function y = f(u,p) for out-of-place (default) or an IntegralFunction or BatchIntegralFunction for inplace and batching optimizations.\ndomain: an object representing an integration domain, i.e. the tuple (lb, ub).\nlb: DEPRECATED: Either a number or vector of lower bounds.\nub: DEPRECATED: Either a number or vector of upper bounds.\np: The parameters associated with the problem.\nnout: DEPRECATED (see IntegralFunction): length of the vector output of the integrand (by default the integrand is assumed to be scalar)\nbatch: DEPRECATED (see BatchIntegralFunction): number of points the integrand can evaluate simultaneously (by default there is no batching)\nkwargs: Keyword arguments copied to the solvers.\n\nAdditionally, we can supply iip like IntegralProblem{iip}(...) as true or false to declare at compile time whether the integrator function is in-place.\n\nFields\n\nThe fields match the names of the constructor arguments.\n\n\n\n\n\n","category":"type"},{"location":"basics/IntegralProblem/","page":"Integral Problems","title":"Integral Problems","text":"The correct shape of the variables (inputs) u and the values (outputs) y of the integrand f depends on whether batching is used.","category":"page"},{"location":"basics/IntegralProblem/","page":"Integral Problems","title":"Integral Problems","text":"If batch == 0","category":"page"},{"location":"basics/IntegralProblem/","page":"Integral Problems","title":"Integral Problems","text":" single variable f multiple variable f\nscalar valued f u is a scalar, y is a scalar u is a vector, y is a scalar\nvector valued f u is a scalar, y is a vector u is a vector, y is a vector","category":"page"},{"location":"basics/IntegralProblem/","page":"Integral Problems","title":"Integral Problems","text":"If batch > 0","category":"page"},{"location":"basics/IntegralProblem/","page":"Integral Problems","title":"Integral Problems","text":" single variable f multiple variable f\nscalar valued f u is a vector, y is a vector u is a matrix, y is a vector\nvector valued f u is a vector, y is a matrix u is a matrix, y is a matrix","category":"page"},{"location":"basics/IntegralProblem/","page":"Integral Problems","title":"Integral Problems","text":"The last dimension is always used as the batching dimension, e.g., if u is a matrix, then u[:,i] is the ith point where the integrand will be evaluated.","category":"page"},{"location":"basics/solve/#Common-Solver-Options-(Solve-Keyword-Arguments)","page":"Common Solver Options (Solve Keyword Arguments)","title":"Common Solver Options (Solve Keyword Arguments)","text":"","category":"section"},{"location":"basics/solve/","page":"Common Solver Options (Solve Keyword Arguments)","title":"Common Solver Options (Solve Keyword Arguments)","text":"solve(prob::IntegralProblem, alg::SciMLBase.AbstractIntegralAlgorithm)","category":"page"},{"location":"basics/solve/#CommonSolve.solve-Tuple{IntegralProblem, SciMLBase.AbstractIntegralAlgorithm}","page":"Common Solver Options (Solve Keyword Arguments)","title":"CommonSolve.solve","text":"solve(prob::IntegralProblem, alg::SciMLBase.AbstractIntegralAlgorithm; kwargs...)\n\nKeyword Arguments\n\nThe arguments to solve are common across all of the quadrature methods. These common arguments are:\n\nmaxiters (the maximum number of iterations)\nabstol (absolute tolerance in changes of the objective value)\nreltol (relative tolerance in changes of the objective value)\n\n\n\n\n\n","category":"method"},{"location":"basics/solve/","page":"Common Solver Options (Solve Keyword Arguments)","title":"Common Solver Options (Solve Keyword Arguments)","text":"Additionally, the extra keyword arguments are splatted to the library calls, so see the documentation of the integrator library for all the extra details. These extra keyword arguments are not guaranteed to act uniformly.","category":"page"},{"location":"tutorials/differentiating_integrals/#Differentiating-Integrals","page":"Differentiating Integrals","title":"Differentiating Integrals","text":"","category":"section"},{"location":"tutorials/differentiating_integrals/","page":"Differentiating Integrals","title":"Differentiating Integrals","text":"Integrals.jl is a fully differentiable quadrature library. Thus, it adds the ability to perform automatic differentiation over any of the libraries that it calls. It integrates with ForwardDiff.jl for forward-mode automatic differentiation and Zygote.jl for reverse-mode automatic differentiation. For example:","category":"page"},{"location":"tutorials/differentiating_integrals/","page":"Differentiating Integrals","title":"Differentiating Integrals","text":"using Integrals, ForwardDiff, FiniteDiff, Zygote, Cuba\nf(x, p) = sum(sin.(x .* p))\ndomain = (ones(2), 3ones(2)) # (lb, ub)\np = ones(2)\n\nfunction testf(p)\n prob = IntegralProblem(f, domain, p)\n sin(solve(prob, CubaCuhre(), reltol = 1e-6, abstol = 1e-6)[1])\nend\ntestf(p)\ndp1 = Zygote.gradient(testf, p)\ndp2 = FiniteDiff.finite_difference_gradient(testf, p)\ndp3 = ForwardDiff.gradient(testf, p)\ndp1[1] ≈ dp2 ≈ dp3","category":"page"},{"location":"basics/FAQ/#Frequently-Asked-Questions","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"","category":"section"},{"location":"basics/FAQ/#How-should-I-use-the-in-place-interface?","page":"Frequently Asked Questions","title":"How should I use the in-place interface?","text":"","category":"section"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"The in-place interface allows evaluating vector-valued integrands without allocating an output array. This can be beneficial for reducing allocations when integrating many functions simultaneously or to make use of existing in-place code. However, note that not all algorithms use in-place operations under the hood, i.e. HCubatureJL, and may still allocate.","category":"page"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"You can construct an IntegralFunction(f, prototype), where f is of the form f(y, u, p) where prototype is of the desired type and shape of y.","category":"page"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"For small array outputs of a known size, consider using StaticArrays.jl for the return value of your integrand.","category":"page"},{"location":"basics/FAQ/#How-should-I-use-the-batch-interface?","page":"Frequently Asked Questions","title":"How should I use the batch interface?","text":"","category":"section"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"The batch interface allows evaluating one (or more) integrals simultaneously at different points, which maximizes the parallelism for a given algorithm.","category":"page"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"You can construct an out-of-place BatchIntegralFunction(bf) where bf is of the form bf(u, p) = stack(x -> f(x, p), eachslice(u; dims=ndims(u))), where f is the (unbatched) integrand. For interoperability with as many algorithms as possible, it is important that your out-of-place batch integrand accept an empty array of quadrature points and still return an output with a size and type consistent with the non-empty case.","category":"page"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"You can construct an in-place BatchIntegralFunction(bf, prototype), where bf is of the form bf(y, u, p) = foreach((y,x) -> f(y,x,p), eachslice(y, dims=ndims(y)), eachslice(x, dims=ndims(x))).","category":"page"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"Note that not all algorithms use in-place batched operations under the hood, i.e. QuadGKJL.","category":"page"},{"location":"basics/FAQ/#What-should-I-do-if-my-solution-is-not-converged?","page":"Frequently Asked Questions","title":"What should I do if my solution is not converged?","text":"","category":"section"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"Certain algorithms, such as QuadratureRule used a fixed number of points to calculate an integral and cannot provide an error estimate. In this case, you have to increase the number of points and check the convergence yourself, which will depend on the accuracy of the rule you choose.","category":"page"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"For badly-behaved integrands, such as (nearly) singular and highly oscillatory functions, most algorithms will fail to converge and either throw an error or silently return the incorrect result. In some cases Integrals.jl can provide an error code when things go wrong, but otherwise you can always check if the error estimate for the integral is less than the requested tolerance, e.g. sol.resid < max(abstol, reltol*norm(sol.u)). Sometimes using a larger tolerance or higher precision arithmetic may help.","category":"page"},{"location":"basics/FAQ/#How-can-I-integrate-arbitrarily-spaced-data?","page":"Frequently Asked Questions","title":"How can I integrate arbitrarily-spaced data?","text":"","category":"section"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"See SampledIntegralProblem.","category":"page"},{"location":"basics/FAQ/#How-can-I-integrate-on-arbitrary-geometries?","page":"Frequently Asked Questions","title":"How can I integrate on arbitrary geometries?","text":"","category":"section"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"You can't, since Integrals.jl currently supports integration on hypercubes because that is what lower-level packages implement.","category":"page"},{"location":"basics/FAQ/#I-don't-see-algorithm-X-or-quadrature-scheme-Y-?","page":"Frequently Asked Questions","title":"I don't see algorithm X or quadrature scheme Y ?","text":"","category":"section"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"Fixed quadrature rules from other packages can be used with QuadratureRule. Otherwise, feel free to open an issue or pull request.","category":"page"},{"location":"basics/FAQ/#My-integrand-works-with-algorithm-X-but-fails-on-algorithm-Y","page":"Frequently Asked Questions","title":"My integrand works with algorithm X but fails on algorithm Y","text":"","category":"section"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"While bugs are not out of the question, certain algorithms, especially those implemented in C, are not compatible with arbitrary Julia types and have to return specific numeric types or arrays thereof. In some cases, such as ArblibJL, it is also expected that the integrand work with a custom quadrature point type. Moreover, some algorithms, such as VEGAS, only support scalar integrands. For more details see the solver page.","category":"page"},{"location":"basics/FAQ/#Can-I-take-derivatives-with-respect-to-the-limits-of-integration?","page":"Frequently Asked Questions","title":"Can I take derivatives with respect to the limits of integration?","text":"","category":"section"},{"location":"basics/FAQ/","page":"Frequently Asked Questions","title":"Frequently Asked Questions","text":"Currently this is not implemented.","category":"page"},{"location":"#Integrals.jl:-Unified-Integral-Approximation-Interface","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"","category":"section"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"Integrals.jl is a unified interface for the numerical approximation of integrals (quadrature) in Julia. It interfaces with other packages of the Julia ecosystem to make it easy to test alternative solver packages and pass small types to control algorithm swapping.","category":"page"},{"location":"#Installation","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Installation","text":"","category":"section"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"To install Integrals.jl, use the Julia package manager:","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"using Pkg\nPkg.add(\"Integrals\")","category":"page"},{"location":"#Contributing","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Contributing","text":"","category":"section"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"Please refer to the SciML ColPrac: Contributor's Guide on Collaborative Practices for Community Packages for guidance on PRs, issues, and other matters relating to contributing to SciML.\nSee the SciML Style Guide for common coding practices and other style decisions.\nThere are a few community forums:\nThe #diffeq-bridged and #sciml-bridged channels in the Julia Slack\nThe #diffeq-bridged and #sciml-bridged channels in the Julia Zulip\nOn the Julia Discourse forums\nSee also SciML Community page","category":"page"},{"location":"#Reproducibility","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Reproducibility","text":"","category":"section"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"
The documentation of this SciML package was built using these direct dependencies,","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"using Pkg # hide\nPkg.status() # hide","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"
","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"
and using this machine and Julia version.","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"using InteractiveUtils # hide\nversioninfo() # hide","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"
","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"
A more complete overview of all dependencies and their versions is also provided.","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"using Pkg # hide\nPkg.status(; mode = PKGMODE_MANIFEST) # hide","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"
","category":"page"},{"location":"","page":"Integrals.jl: Unified Integral Approximation Interface","title":"Integrals.jl: Unified Integral Approximation Interface","text":"using TOML\nusing Markdown\nversion = TOML.parse(read(\"../../Project.toml\", String))[\"version\"]\nname = TOML.parse(read(\"../../Project.toml\", String))[\"name\"]\nlink_manifest = \"https://github.com/SciML/\" * name * \".jl/tree/gh-pages/v\" * version *\n \"/assets/Manifest.toml\"\nlink_project = \"https://github.com/SciML/\" * name * \".jl/tree/gh-pages/v\" * version *\n \"/assets/Project.toml\"\nMarkdown.parse(\"\"\"You can also download the\n[manifest]($link_manifest)\nfile and the\n[project]($link_project)\nfile.\n\"\"\")","category":"page"},{"location":"basics/IntegralFunction/#func","page":"Integral Functions","title":"Integral Functions","text":"","category":"section"},{"location":"basics/IntegralFunction/","page":"Integral Functions","title":"Integral Functions","text":"SciMLBase.IntegralFunction\nSciMLBase.BatchIntegralFunction","category":"page"},{"location":"basics/IntegralFunction/#SciMLBase.IntegralFunction","page":"Integral Functions","title":"SciMLBase.IntegralFunction","text":"IntegralFunction{iip,specialize,F,T} <: AbstractIntegralFunction{iip}\n\nA representation of an integrand f defined by:\n\nf(u p)\n\nFor an in-place form of f see the iip section below for details on in-place or out-of-place handling.\n\nIntegralFunction{iip,specialize}(f, [integrand_prototype])\n\nNote that only f is required, and in the case of inplace integrands a mutable container integrand_prototype to store the result of the integrand. If integrand_prototype is present, f is interpreted as in-place, and otherwise f is assumed to be out-of-place.\n\niip: In-Place vs Out-Of-Place\n\nOut-of-place functions must be of the form y = f(u p) and in-place functions of the form f(y u p), where y is a number or array containing the output. Since f is allowed to return any type (e.g. real or complex numbers or arrays), in-place functions must provide a container integrand_prototype that is of the right type and size for the variable y, and the result is written to this container in-place. When in-place forms are used, in-place array operations, i.e. broadcasting, may be used by algorithms to reduce allocations. If integrand_prototype is not provided, f is assumed to be out-of-place.\n\nspecialize\n\nThis field is currently unused\n\nFields\n\nThe fields of the IntegralFunction type directly match the names of the inputs.\n\n\n\n\n\n","category":"type"},{"location":"basics/IntegralFunction/#SciMLBase.BatchIntegralFunction","page":"Integral Functions","title":"SciMLBase.BatchIntegralFunction","text":"BatchIntegralFunction{iip,specialize,F,T} <: AbstractIntegralFunction{iip}\n\nA batched representation of an (non-batched) integrand f(u, p) that can be evaluated at multiple points simultaneously using threads, the gpu, or distributed memory defined by:\n\nby = bf(bu p)\n\nHere we prefix variables with b to indicate they are batched variables, which implies that they are arrays whose last dimension is reserved for batching different evaluation points, or function values, and may be of a variable length. bu is an array whose elements correspond to distinct evaluation points to f, and bf is a function to evaluate f 'point-wise' so that f(bu[..., i], p) == bf(bu, p)[..., i]. For example, a simple batching implementation of a scalar, univariate function is via broadcasting: bf(bu, p) = f.(bu, Ref(p)), although this interface exists in order to allow user parallelization. In general, the integration algorithm is allowed to vary the number of evaluation points between subsequent calls to bf.\n\nFor an in-place form of bf see the iip section below for details on in-place or out-of-place handling.\n\nBatchIntegralFunction{iip,specialize}(bf, [integrand_prototype];\n max_batch=typemax(Int))\n\nNote that only bf is required, and in the case of inplace integrands a mutable container integrand_prototype to store a batch of integrand evaluations, with a last \"batching\" dimension.\n\nThe keyword max_batch is used to set a soft limit on the number of points to batch at the same time so that memory usage is controlled.\n\nIf integrand_prototype is present, bf is interpreted as in-place, and otherwise bf is assumed to be out-of-place.\n\niip: In-Place vs Out-Of-Place\n\nOut-of-place functions must be of the form by = bf(bu, p) and in-place functions of the form bf(by, bu, p) where by is a batch array containing the output. Since the algorithm may vary the number of points to batch, the batching dimension can be of any length, including zero, and since bf is allowed to return arrays of any type (e.g. real or complex) or size, in-place functions must provide a container integrand_prototype of the desired type and size for by. If integrand_prototype is not provided, bf is assumed to be out-of-place.\n\nIn the out-of-place case, we require f(bu[..., i], p) == bf(bu, p)[..., i], and certain algorithms, such as those implemented in C, may infer the type or shape of by by calling bf with an empty array of input points, i.e. bu with size(bu)[end] == 0. Then it is expected for the resulting by to have the same type and size(by)[begin:end-1] for all subsequent calls.\n\nWhen the in-place form is used, we require f(by[..., i], bu[..., i], p) == bf(by, bu, p)[..., i] and size(by)[begin:end-1] == size(integrand_prototype)[begin:end-1]. The algorithm should always pass the integrand by arrays that are similar to integrand_prototype, and may use views and in-place array operations to reduce allocations.\n\nspecialize\n\nThis field is currently unused\n\nFields\n\nThe fields of the BatchIntegralFunction type are f, corresponding to bf above, and integrand_prototype.\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#solvers","page":"Integral Solver Algorithms","title":"Integral Solver Algorithms","text":"","category":"section"},{"location":"solvers/IntegralSolvers/","page":"Integral Solver Algorithms","title":"Integral Solver Algorithms","text":"The following algorithms are available:","category":"page"},{"location":"solvers/IntegralSolvers/","page":"Integral Solver Algorithms","title":"Integral Solver Algorithms","text":"QuadGKJL: Uses QuadGK.jl, which supports one-dimensional integration of scalar and array-valued integrands with in-place or batched forms. Integrands that are both in-place and batched are implemented in the wrapper but are not supported under the hood.\nHCubatureJL: Uses HCubature.jl, which supports scalar and array-valued integrands and works best in low dimensions, e.g. ≤ 8. In-place integrands are implemented in the wrapper but are not supported under the hood. Batching is not supported.\nVEGAS: Uses MonteCarloIntegration.jl, which requires scalar, Float64-valued integrands and works in any number of dimensions.\nVEGASMC: Uses MCIntegration.jl. Requires using MCIntegration. Doesn't support batching.\nCubatureJLh: h-Cubature from Cubature.jl. Requires using Cubature.\nCubatureJLp: p-Cubature from Cubature.jl. Requires using Cubature.\nCubaVegas: Vegas from Cuba.jl. Requires using Cuba.\nCubaSUAVE: SUAVE from Cuba.jl. Requires using Cuba.\nCubaDivonne: Divonne from Cuba.jl. Requires using Cuba. Works only for >1-dimensional integrations.\nCubaCuhre: Cuhre from Cuba.jl. Requires using Cuba. Works only for >1-dimensional integrations.\nGaussLegendre: Performs fixed-order Gauss-Legendre quadrature. Requires using FastGaussQuadrature.\nQuadratureRule: Accepts a user-defined function that returns nodes and weights.\nArblibJL: real- and complex-valued univariate integration of holomorphic and meromorphic functions from Arblib.jl. Requires using Arblib.","category":"page"},{"location":"solvers/IntegralSolvers/","page":"Integral Solver Algorithms","title":"Integral Solver Algorithms","text":"QuadGKJL\nHCubatureJL\nCubatureJLp\nCubatureJLh\nVEGAS\nVEGASMC\nCubaVegas\nCubaSUAVE\nCubaDivonne\nCubaCuhre\nGaussLegendre\nQuadratureRule\nArblibJL","category":"page"},{"location":"solvers/IntegralSolvers/#Integrals.QuadGKJL","page":"Integral Solver Algorithms","title":"Integrals.QuadGKJL","text":"QuadGKJL(; order = 7, norm=norm)\n\nOne-dimensional Gauss-Kronrod integration from QuadGK.jl. This method also takes the optional arguments order and norm. Which are the order of the integration rule and the norm for calculating the error, respectively\n\nReferences\n\n@article{laurie1997calculation, title={Calculation of Gauss-Kronrod quadrature rules}, author={Laurie, Dirk}, journal={Mathematics of Computation}, volume={66}, number={219}, pages={1133–1145}, year={1997} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.HCubatureJL","page":"Integral Solver Algorithms","title":"Integrals.HCubatureJL","text":"HCubatureJL(; norm=norm, initdiv=1)\n\nMultidimensional \"h-adaptive\" integration from HCubature.jl. This method also takes the optional arguments initdiv and norm. Which are the initial number of segments each dimension of the integration domain is divided into, and the norm for calculating the error, respectively.\n\nReferences\n\n@article{genz1980remarks, title={Remarks on algorithm 006: An adaptive algorithm for numerical integration over an N-dimensional rectangular region}, author={Genz, Alan C and Malik, Aftab Ahmad}, journal={Journal of Computational and Applied mathematics}, volume={6}, number={4}, pages={295–302}, year={1980}, publisher={Elsevier} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.CubatureJLp","page":"Integral Solver Algorithms","title":"Integrals.CubatureJLp","text":"CubatureJLp(; error_norm=Cubature.INDIVIDUAL)\n\nMultidimensional p-adaptive integration from Cubature.jl. This method is based on repeatedly doubling the degree of the cubature rules, until convergence is achieved. The used cubature rule is a tensor product of Clenshaw–Curtis quadrature rules. error_norm specifies the convergence criterion for vector valued integrands. Defaults to Cubature.INDIVIDUAL, other options are Cubature.PAIRED, Cubature.L1, Cubature.L2, or Cubature.LINF.\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.CubatureJLh","page":"Integral Solver Algorithms","title":"Integrals.CubatureJLh","text":"CubatureJLh(; error_norm=Cubature.INDIVIDUAL)\n\nMultidimensional h-adaptive integration from Cubature.jl. error_norm specifies the convergence criterion for vector valued integrands. Defaults to Cubature.INDIVIDUAL, other options are Cubature.PAIRED, Cubature.L1, Cubature.L2, or Cubature.LINF.\n\nReferences\n\n@article{genz1980remarks, title={Remarks on algorithm 006: An adaptive algorithm for numerical integration over an N-dimensional rectangular region}, author={Genz, Alan C and Malik, Aftab Ahmad}, journal={Journal of Computational and Applied mathematics}, volume={6}, number={4}, pages={295–302}, year={1980}, publisher={Elsevier} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.VEGAS","page":"Integral Solver Algorithms","title":"Integrals.VEGAS","text":"VEGAS(; nbins = 100, ncalls = 1000, debug=false)\n\nMultidimensional adaptive Monte Carlo integration from MonteCarloIntegration.jl. Importance sampling is used to reduce variance. This method also takes three optional arguments nbins, ncalls and debug which are the initial number of bins each dimension of the integration domain is divided into, the number of function calls per iteration of the algorithm, and whether debug info should be printed, respectively.\n\nLimitations\n\nThis algorithm can only integrate Float64-valued functions\n\nReferences\n\n@article{lepage1978new, title={A new algorithm for adaptive multidimensional integration}, author={Lepage, G Peter}, journal={Journal of Computational Physics}, volume={27}, number={2}, pages={192–203}, year={1978}, publisher={Elsevier} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.VEGASMC","page":"Integral Solver Algorithms","title":"Integrals.VEGASMC","text":"VEGASMC(; kws...)\n\nMarkov-chain based Vegas algorithm from MCIntegration.jl\n\nRefer to MCIntegration.integrate for documentation on the keywords, which are passed directly to the solver with a set of defaults that works for conforming integrands.\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.CubaVegas","page":"Integral Solver Algorithms","title":"Integrals.CubaVegas","text":"CubaVegas()\n\nMultidimensional adaptive Monte Carlo integration from Cuba.jl. Importance sampling is used to reduce variance.\n\nReferences\n\n@article{lepage1978new, title={A new algorithm for adaptive multidimensional integration}, author={Lepage, G Peter}, journal={Journal of Computational Physics}, volume={27}, number={2}, pages={192–203}, year={1978}, publisher={Elsevier} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.CubaSUAVE","page":"Integral Solver Algorithms","title":"Integrals.CubaSUAVE","text":"CubaSUAVE()\n\nMultidimensional adaptive Monte Carlo integration from Cuba.jl. Suave stands for subregion-adaptive VEGAS. Importance sampling and subdivision are thus used to reduce variance.\n\nReferences\n\n@article{hahn2005cuba, title={Cuba—a library for multidimensional numerical integration}, author={Hahn, Thomas}, journal={Computer Physics Communications}, volume={168}, number={2}, pages={78–95}, year={2005}, publisher={Elsevier} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.CubaDivonne","page":"Integral Solver Algorithms","title":"Integrals.CubaDivonne","text":"CubaDivonne()\n\nMultidimensional adaptive Monte Carlo integration from Cuba.jl. Stratified sampling is used to reduce variance.\n\nReferences\n\n@article{friedman1981nested, title={A nested partitioning procedure for numerical multiple integration}, author={Friedman, Jerome H and Wright, Margaret H}, journal={ACM Transactions on Mathematical Software (TOMS)}, volume={7}, number={1}, pages={76–92}, year={1981}, publisher={ACM New York, NY, USA} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.CubaCuhre","page":"Integral Solver Algorithms","title":"Integrals.CubaCuhre","text":"CubaCuhre()\n\nMultidimensional h-adaptive integration from Cuba.jl.\n\nReferences\n\n@article{berntsen1991adaptive, title={An adaptive algorithm for the approximate calculation of multiple integrals}, author={Berntsen, Jarle and Espelid, Terje O and Genz, Alan}, journal={ACM Transactions on Mathematical Software (TOMS)}, volume={17}, number={4}, pages={437–451}, year={1991}, publisher={ACM New York, NY, USA} }\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.GaussLegendre","page":"Integral Solver Algorithms","title":"Integrals.GaussLegendre","text":"GaussLegendre{C, N, W}\n\nStruct for evaluating an integral via (composite) Gauss-Legendre quadrature. The field C will be true if subintervals > 1, and false otherwise.\n\nThe fields nodes::N and weights::W are defined by nodes, weights = gausslegendre(n) for a given number of nodes n.\n\nThe field subintervals::Int64 = 1 (with default value 1) defines the number of intervals to partition the original interval of integration [a, b] into, splitting it into [xⱼ, xⱼ₊₁] for j = 1,…,subintervals, where xⱼ = a + (j-1)h and h = (b-a)/subintervals. Gauss-Legendre quadrature is then applied on each subinterval. For example, if [a, b] = [-1, 1] and subintervals = 2, then Gauss-Legendre quadrature will be applied separately on [-1, 0] and [0, 1], summing the two results.\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.QuadratureRule","page":"Integral Solver Algorithms","title":"Integrals.QuadratureRule","text":"QuadratureRule(q; n=250)\n\nAlgorithm to construct and evaluate a quadrature rule q of n points computed from the inputs as x, w = q(n). It assumes the nodes and weights are for the standard interval [-1, 1]^d in d dimensions, and rescales the nodes to the specific hypercube being solved. The nodes x may be scalars in 1d or vectors in arbitrary dimensions, and the weights w must be scalar. The algorithm computes the quadrature rule sum(w .* f.(x)) and the caller must check that the result is converged with respect to n.\n\n\n\n\n\n","category":"type"},{"location":"solvers/IntegralSolvers/#Integrals.ArblibJL","page":"Integral Solver Algorithms","title":"Integrals.ArblibJL","text":"ArblibJL(; check_analytic=false, take_prec=false, warn_on_no_convergence=false, opts=C_NULL)\n\nOne-dimensional adaptive Gauss-Legendre integration using rigorous error bounds and precision ball arithmetic. Generally this assumes the integrand is holomorphic or meromorphic, which is the user's responsibility to verify. The result of the integral is not guaranteed to satisfy the requested tolerances, however the result is guaranteed to be within the error estimate.\n\nArblib.jl only supports integration of univariate real- and complex-valued functions with both inplace and out-of-place forms. See their documentation for additional details the algorithm arguments and on implementing high-precision integrands. Additionally, the error estimate is included in the return value of the integral, representing a ball.\n\n\n\n\n\n","category":"type"}] } diff --git a/dev/solvers/IntegralSolvers/index.html b/dev/solvers/IntegralSolvers/index.html index 58955ea..871b66f 100644 --- a/dev/solvers/IntegralSolvers/index.html +++ b/dev/solvers/IntegralSolvers/index.html @@ -1,2 +1,2 @@ -Integral Solver Algorithms · Integrals.jl
GitHub

Integral Solver Algorithms

The following algorithms are available:

  • QuadGKJL: Uses QuadGK.jl, which supports one-dimensional integration of scalar and array-valued integrands with in-place or batched forms. Integrands that are both in-place and batched are implemented in the wrapper but are not supported under the hood.
  • HCubatureJL: Uses HCubature.jl, which supports scalar and array-valued integrands and works best in low dimensions, e.g. ≤ 8. In-place integrands are implemented in the wrapper but are not supported under the hood. Batching is not supported.
  • VEGAS: Uses MonteCarloIntegration.jl, which requires scalar, Float64-valued integrands and works in any number of dimensions.
  • VEGASMC: Uses MCIntegration.jl. Requires using MCIntegration. Doesn't support batching.
  • CubatureJLh: h-Cubature from Cubature.jl. Requires using Cubature.
  • CubatureJLp: p-Cubature from Cubature.jl. Requires using Cubature.
  • CubaVegas: Vegas from Cuba.jl. Requires using Cuba.
  • CubaSUAVE: SUAVE from Cuba.jl. Requires using Cuba.
  • CubaDivonne: Divonne from Cuba.jl. Requires using Cuba. Works only for >1-dimensional integrations.
  • CubaCuhre: Cuhre from Cuba.jl. Requires using Cuba. Works only for >1-dimensional integrations.
  • GaussLegendre: Performs fixed-order Gauss-Legendre quadrature. Requires using FastGaussQuadrature.
  • QuadratureRule: Accepts a user-defined function that returns nodes and weights.
  • ArblibJL: real- and complex-valued univariate integration of holomorphic and meromorphic functions from Arblib.jl. Requires using Arblib.
Integrals.QuadGKJL — Type
QuadGKJL(; order = 7, norm=norm)

One-dimensional Gauss-Kronrod integration from QuadGK.jl. This method also takes the optional arguments order and norm. Which are the order of the integration rule and the norm for calculating the error, respectively

References

@article{laurie1997calculation, title={Calculation of Gauss-Kronrod quadrature rules}, author={Laurie, Dirk}, journal={Mathematics of Computation}, volume={66}, number={219}, pages={1133–1145}, year={1997} }

source
Integrals.HCubatureJL — Type
HCubatureJL(; norm=norm, initdiv=1)

Multidimensional "h-adaptive" integration from HCubature.jl. This method also takes the optional arguments initdiv and norm. Which are the initial number of segments each dimension of the integration domain is divided into, and the norm for calculating the error, respectively.

References

@article{genz1980remarks, title={Remarks on algorithm 006: An adaptive algorithm for numerical integration over an N-dimensional rectangular region}, author={Genz, Alan C and Malik, Aftab Ahmad}, journal={Journal of Computational and Applied mathematics}, volume={6}, number={4}, pages={295–302}, year={1980}, publisher={Elsevier} }

source
Integrals.CubatureJLp — Type
CubatureJLp(; error_norm=Cubature.INDIVIDUAL)

Multidimensional p-adaptive integration from Cubature.jl. This method is based on repeatedly doubling the degree of the cubature rules, until convergence is achieved. The used cubature rule is a tensor product of Clenshaw–Curtis quadrature rules. error_norm specifies the convergence criterion for vector valued integrands. Defaults to Cubature.INDIVIDUAL, other options are Cubature.PAIRED, Cubature.L1, Cubature.L2, or Cubature.LINF.

source
Integrals.CubatureJLh — Type
CubatureJLh(; error_norm=Cubature.INDIVIDUAL)

Multidimensional h-adaptive integration from Cubature.jl. error_norm specifies the convergence criterion for vector valued integrands. Defaults to Cubature.INDIVIDUAL, other options are Cubature.PAIRED, Cubature.L1, Cubature.L2, or Cubature.LINF.

References

@article{genz1980remarks, title={Remarks on algorithm 006: An adaptive algorithm for numerical integration over an N-dimensional rectangular region}, author={Genz, Alan C and Malik, Aftab Ahmad}, journal={Journal of Computational and Applied mathematics}, volume={6}, number={4}, pages={295–302}, year={1980}, publisher={Elsevier} }

source
Integrals.VEGAS — Type
VEGAS(; nbins = 100, ncalls = 1000, debug=false)

Multidimensional adaptive Monte Carlo integration from MonteCarloIntegration.jl. Importance sampling is used to reduce variance. This method also takes three optional arguments nbins, ncalls and debug which are the initial number of bins each dimension of the integration domain is divided into, the number of function calls per iteration of the algorithm, and whether debug info should be printed, respectively.

Limitations

This algorithm can only integrate Float64-valued functions

References

@article{lepage1978new, title={A new algorithm for adaptive multidimensional integration}, author={Lepage, G Peter}, journal={Journal of Computational Physics}, volume={27}, number={2}, pages={192–203}, year={1978}, publisher={Elsevier} }

source
Integrals.VEGASMC — Type
VEGASMC(; kws...)

Markov-chain based Vegas algorithm from MCIntegration.jl

Refer to MCIntegration.integrate for documentation on the keywords, which are passed directly to the solver with a set of defaults that works for conforming integrands.

source
Integrals.CubaVegas — Type
CubaVegas()

Multidimensional adaptive Monte Carlo integration from Cuba.jl. Importance sampling is used to reduce variance.

References

@article{lepage1978new, title={A new algorithm for adaptive multidimensional integration}, author={Lepage, G Peter}, journal={Journal of Computational Physics}, volume={27}, number={2}, pages={192–203}, year={1978}, publisher={Elsevier} }

source
Integrals.CubaSUAVE — Type
CubaSUAVE()

Multidimensional adaptive Monte Carlo integration from Cuba.jl. Suave stands for subregion-adaptive VEGAS. Importance sampling and subdivision are thus used to reduce variance.

References

@article{hahn2005cuba, title={Cuba—a library for multidimensional numerical integration}, author={Hahn, Thomas}, journal={Computer Physics Communications}, volume={168}, number={2}, pages={78–95}, year={2005}, publisher={Elsevier} }

source
Integrals.CubaDivonne — Type
CubaDivonne()

Multidimensional adaptive Monte Carlo integration from Cuba.jl. Stratified sampling is used to reduce variance.

References

@article{friedman1981nested, title={A nested partitioning procedure for numerical multiple integration}, author={Friedman, Jerome H and Wright, Margaret H}, journal={ACM Transactions on Mathematical Software (TOMS)}, volume={7}, number={1}, pages={76–92}, year={1981}, publisher={ACM New York, NY, USA} }

source
Integrals.CubaCuhre — Type
CubaCuhre()

Multidimensional h-adaptive integration from Cuba.jl.

References

@article{berntsen1991adaptive, title={An adaptive algorithm for the approximate calculation of multiple integrals}, author={Berntsen, Jarle and Espelid, Terje O and Genz, Alan}, journal={ACM Transactions on Mathematical Software (TOMS)}, volume={17}, number={4}, pages={437–451}, year={1991}, publisher={ACM New York, NY, USA} }

source
Integrals.GaussLegendre — Type
GaussLegendre{C, N, W}

Struct for evaluating an integral via (composite) Gauss-Legendre quadrature. The field C will be true if subintervals > 1, and false otherwise.

The fields nodes::N and weights::W are defined by nodes, weights = gausslegendre(n) for a given number of nodes n.

The field subintervals::Int64 = 1 (with default value 1) defines the number of intervals to partition the original interval of integration [a, b] into, splitting it into [xⱼ, xⱼ₊₁] for j = 1,…,subintervals, where xⱼ = a + (j-1)h and h = (b-a)/subintervals. Gauss-Legendre quadrature is then applied on each subinterval. For example, if [a, b] = [-1, 1] and subintervals = 2, then Gauss-Legendre quadrature will be applied separately on [-1, 0] and [0, 1], summing the two results.

source
Integrals.QuadratureRule — Type

QuadratureRule(q; n=250)

Algorithm to construct and evaluate a quadrature rule q of n points computed from the inputs as x, w = q(n). It assumes the nodes and weights are for the standard interval [-1, 1]^d in d dimensions, and rescales the nodes to the specific hypercube being solved. The nodes x may be scalars in 1d or vectors in arbitrary dimensions, and the weights w must be scalar. The algorithm computes the quadrature rule sum(w .* f.(x)) and the caller must check that the result is converged with respect to n.

source
Integrals.ArblibJL — Type
ArblibJL(; check_analytic=false, take_prec=false, warn_on_no_convergence=false, opts=C_NULL)

One-dimensional adaptive Gauss-Legendre integration using rigorous error bounds and precision ball arithmetic. Generally this assumes the integrand is holomorphic or meromorphic, which is the user's responsibility to verify. The result of the integral is not guaranteed to satisfy the requested tolerances, however the result is guaranteed to be within the error estimate.

Arblib.jl only supports integration of univariate real- and complex-valued functions with both inplace and out-of-place forms. See their documentation for additional details the algorithm arguments and on implementing high-precision integrands. Additionally, the error estimate is included in the return value of the integral, representing a ball.

source
+Integral Solver Algorithms · Integrals.jl
GitHub

Integral Solver Algorithms

The following algorithms are available:

  • QuadGKJL: Uses QuadGK.jl, which supports one-dimensional integration of scalar and array-valued integrands with in-place or batched forms. Integrands that are both in-place and batched are implemented in the wrapper but are not supported under the hood.
  • HCubatureJL: Uses HCubature.jl, which supports scalar and array-valued integrands and works best in low dimensions, e.g. ≤ 8. In-place integrands are implemented in the wrapper but are not supported under the hood. Batching is not supported.
  • VEGAS: Uses MonteCarloIntegration.jl, which requires scalar, Float64-valued integrands and works in any number of dimensions.
  • VEGASMC: Uses MCIntegration.jl. Requires using MCIntegration. Doesn't support batching.
  • CubatureJLh: h-Cubature from Cubature.jl. Requires using Cubature.
  • CubatureJLp: p-Cubature from Cubature.jl. Requires using Cubature.
  • CubaVegas: Vegas from Cuba.jl. Requires using Cuba.
  • CubaSUAVE: SUAVE from Cuba.jl. Requires using Cuba.
  • CubaDivonne: Divonne from Cuba.jl. Requires using Cuba. Works only for >1-dimensional integrations.
  • CubaCuhre: Cuhre from Cuba.jl. Requires using Cuba. Works only for >1-dimensional integrations.
  • GaussLegendre: Performs fixed-order Gauss-Legendre quadrature. Requires using FastGaussQuadrature.
  • QuadratureRule: Accepts a user-defined function that returns nodes and weights.
  • ArblibJL: real- and complex-valued univariate integration of holomorphic and meromorphic functions from Arblib.jl. Requires using Arblib.
Integrals.QuadGKJL — Type
QuadGKJL(; order = 7, norm=norm)

One-dimensional Gauss-Kronrod integration from QuadGK.jl. This method also takes the optional arguments order and norm. Which are the order of the integration rule and the norm for calculating the error, respectively

References

@article{laurie1997calculation, title={Calculation of Gauss-Kronrod quadrature rules}, author={Laurie, Dirk}, journal={Mathematics of Computation}, volume={66}, number={219}, pages={1133–1145}, year={1997} }

source
Integrals.HCubatureJL — Type
HCubatureJL(; norm=norm, initdiv=1)

Multidimensional "h-adaptive" integration from HCubature.jl. This method also takes the optional arguments initdiv and norm. Which are the initial number of segments each dimension of the integration domain is divided into, and the norm for calculating the error, respectively.

References

@article{genz1980remarks, title={Remarks on algorithm 006: An adaptive algorithm for numerical integration over an N-dimensional rectangular region}, author={Genz, Alan C and Malik, Aftab Ahmad}, journal={Journal of Computational and Applied mathematics}, volume={6}, number={4}, pages={295–302}, year={1980}, publisher={Elsevier} }

source
Integrals.CubatureJLp — Type
CubatureJLp(; error_norm=Cubature.INDIVIDUAL)

Multidimensional p-adaptive integration from Cubature.jl. This method is based on repeatedly doubling the degree of the cubature rules, until convergence is achieved. The used cubature rule is a tensor product of Clenshaw–Curtis quadrature rules. error_norm specifies the convergence criterion for vector valued integrands. Defaults to Cubature.INDIVIDUAL, other options are Cubature.PAIRED, Cubature.L1, Cubature.L2, or Cubature.LINF.

source
Integrals.CubatureJLh — Type
CubatureJLh(; error_norm=Cubature.INDIVIDUAL)

Multidimensional h-adaptive integration from Cubature.jl. error_norm specifies the convergence criterion for vector valued integrands. Defaults to Cubature.INDIVIDUAL, other options are Cubature.PAIRED, Cubature.L1, Cubature.L2, or Cubature.LINF.

References

@article{genz1980remarks, title={Remarks on algorithm 006: An adaptive algorithm for numerical integration over an N-dimensional rectangular region}, author={Genz, Alan C and Malik, Aftab Ahmad}, journal={Journal of Computational and Applied mathematics}, volume={6}, number={4}, pages={295–302}, year={1980}, publisher={Elsevier} }

source
Integrals.VEGAS — Type
VEGAS(; nbins = 100, ncalls = 1000, debug=false)

Multidimensional adaptive Monte Carlo integration from MonteCarloIntegration.jl. Importance sampling is used to reduce variance. This method also takes three optional arguments nbins, ncalls and debug which are the initial number of bins each dimension of the integration domain is divided into, the number of function calls per iteration of the algorithm, and whether debug info should be printed, respectively.

Limitations

This algorithm can only integrate Float64-valued functions

References

@article{lepage1978new, title={A new algorithm for adaptive multidimensional integration}, author={Lepage, G Peter}, journal={Journal of Computational Physics}, volume={27}, number={2}, pages={192–203}, year={1978}, publisher={Elsevier} }

source
Integrals.VEGASMC — Type
VEGASMC(; kws...)

Markov-chain based Vegas algorithm from MCIntegration.jl

Refer to MCIntegration.integrate for documentation on the keywords, which are passed directly to the solver with a set of defaults that works for conforming integrands.

source
Integrals.CubaVegas — Type
CubaVegas()

Multidimensional adaptive Monte Carlo integration from Cuba.jl. Importance sampling is used to reduce variance.

References

@article{lepage1978new, title={A new algorithm for adaptive multidimensional integration}, author={Lepage, G Peter}, journal={Journal of Computational Physics}, volume={27}, number={2}, pages={192–203}, year={1978}, publisher={Elsevier} }

source
Integrals.CubaSUAVE — Type
CubaSUAVE()

Multidimensional adaptive Monte Carlo integration from Cuba.jl. Suave stands for subregion-adaptive VEGAS. Importance sampling and subdivision are thus used to reduce variance.

References

@article{hahn2005cuba, title={Cuba—a library for multidimensional numerical integration}, author={Hahn, Thomas}, journal={Computer Physics Communications}, volume={168}, number={2}, pages={78–95}, year={2005}, publisher={Elsevier} }

source
Integrals.CubaDivonne — Type
CubaDivonne()

Multidimensional adaptive Monte Carlo integration from Cuba.jl. Stratified sampling is used to reduce variance.

References

@article{friedman1981nested, title={A nested partitioning procedure for numerical multiple integration}, author={Friedman, Jerome H and Wright, Margaret H}, journal={ACM Transactions on Mathematical Software (TOMS)}, volume={7}, number={1}, pages={76–92}, year={1981}, publisher={ACM New York, NY, USA} }

source
Integrals.CubaCuhre — Type
CubaCuhre()

Multidimensional h-adaptive integration from Cuba.jl.

References

@article{berntsen1991adaptive, title={An adaptive algorithm for the approximate calculation of multiple integrals}, author={Berntsen, Jarle and Espelid, Terje O and Genz, Alan}, journal={ACM Transactions on Mathematical Software (TOMS)}, volume={17}, number={4}, pages={437–451}, year={1991}, publisher={ACM New York, NY, USA} }

source
Integrals.GaussLegendre — Type
GaussLegendre{C, N, W}

Struct for evaluating an integral via (composite) Gauss-Legendre quadrature. The field C will be true if subintervals > 1, and false otherwise.

The fields nodes::N and weights::W are defined by nodes, weights = gausslegendre(n) for a given number of nodes n.

The field subintervals::Int64 = 1 (with default value 1) defines the number of intervals to partition the original interval of integration [a, b] into, splitting it into [xⱼ, xⱼ₊₁] for j = 1,…,subintervals, where xⱼ = a + (j-1)h and h = (b-a)/subintervals. Gauss-Legendre quadrature is then applied on each subinterval. For example, if [a, b] = [-1, 1] and subintervals = 2, then Gauss-Legendre quadrature will be applied separately on [-1, 0] and [0, 1], summing the two results.

source
Integrals.QuadratureRule — Type

QuadratureRule(q; n=250)

Algorithm to construct and evaluate a quadrature rule q of n points computed from the inputs as x, w = q(n). It assumes the nodes and weights are for the standard interval [-1, 1]^d in d dimensions, and rescales the nodes to the specific hypercube being solved. The nodes x may be scalars in 1d or vectors in arbitrary dimensions, and the weights w must be scalar. The algorithm computes the quadrature rule sum(w .* f.(x)) and the caller must check that the result is converged with respect to n.

source
Integrals.ArblibJL — Type
ArblibJL(; check_analytic=false, take_prec=false, warn_on_no_convergence=false, opts=C_NULL)

One-dimensional adaptive Gauss-Legendre integration using rigorous error bounds and precision ball arithmetic. Generally this assumes the integrand is holomorphic or meromorphic, which is the user's responsibility to verify. The result of the integral is not guaranteed to satisfy the requested tolerances, however the result is guaranteed to be within the error estimate.

Arblib.jl only supports integration of univariate real- and complex-valued functions with both inplace and out-of-place forms. See their documentation for additional details the algorithm arguments and on implementing high-precision integrands. Additionally, the error estimate is included in the return value of the integral, representing a ball.

source
diff --git a/dev/tutorials/caching_interface/index.html b/dev/tutorials/caching_interface/index.html index 6094902..2c2d011 100644 --- a/dev/tutorials/caching_interface/index.html +++ b/dev/tutorials/caching_interface/index.html @@ -68,4 +68,4 @@ 0.35130314420012265 0.32000752789011855 0.28551450214186125 - 0.24816870986674897 + 0.24816870986674897 diff --git a/dev/tutorials/differentiating_integrals/index.html b/dev/tutorials/differentiating_integrals/index.html index 49077f4..bb5c009 100644 --- a/dev/tutorials/differentiating_integrals/index.html +++ b/dev/tutorials/differentiating_integrals/index.html @@ -12,4 +12,4 @@ dp1 = Zygote.gradient(testf, p) dp2 = FiniteDiff.finite_difference_gradient(testf, p) dp3 = ForwardDiff.gradient(testf, p) -dp1[1] ≈ dp2 ≈ dp3
true
+dp1[1] ≈ dp2 ≈ dp3
true
diff --git a/dev/tutorials/numerical_integrals/index.html b/dev/tutorials/numerical_integrals/index.html index 444b7fe..60c50f7 100644 --- a/dev/tutorials/numerical_integrals/index.html +++ b/dev/tutorials/numerical_integrals/index.html @@ -71,4 +71,4 @@ domain = ([-Inf, 0.0], [Inf, Inf]) # (lb, ub) prob = IntegralProblem(f, domain) solve(prob, HCubatureJL(), reltol = 1e-3, abstol = 1e-3)
retcode: Success
-u: 0.49995156348205083
+u: 0.49995156348205083