All content

This is the list of all content in torchquad. The type backend tensor in the documentation is a placeholder for the tensor type of the current numerical backend, for example numpy.array or torch.Tensor.

We are continuously implementing new content in our library. For the code, please see the code page or check out our full code and latest news at https://github.com/esa/torchquad.

class torchquad.BaseIntegrator[source]

Bases: object

The (abstract) integrator that all other integrators inherit from. Provides no explicit definitions for methods.

static evaluate_integrand(fn, points, weights=None, args=None)[source]

Evaluate the integrand function at the passed points

Parameters:
  • fn (function) – Integrand function

  • points (backend tensor) – Integration points

  • weights (backend tensor, optional) – Integration weights. Defaults to None.

  • args (list or tuple, optional) – Any arguments required by the function. Defaults to None.

Returns:

Integrand function output int: Number of evaluated points

Return type:

backend tensor

integrate()[source]
class torchquad.Boole[source]

Bases: NewtonCotes

Boole’s rule. See https://en.wikipedia.org/wiki/Newton%E2%80%93Cotes_formulas#Closed_Newton%E2%80%93Cotes_formulas .

integrate(fn, dim, N=None, integration_domain=None, backend=None)[source]

Integrates the passed function on the passed domain using Boole’s rule.

Parameters:
  • fn (func) – The function to integrate over.

  • dim (int) – Dimensionality of the integration domain.

  • N (int, optional) – Total number of sample points to use for the integration. N has to be such that N^(1/dim) - 1 % 4 == 0. Defaults to 5 points per dimension if None is given.

  • integration_domain (list or backend tensor, optional) – Integration domain, e.g. [[-1,1],[0,1]]. Defaults to [-1,1]^dim. It can also determine the numerical backend.

  • backend (string, optional) – Numerical backend. Defaults to integration_domain’s backend if it is a tensor and otherwise to the backend from the latest call to set_up_backend or “torch” for backwards compatibility.

Returns:

Integral value

Return type:

backend-specific number

class torchquad.GaussLegendre[source]

Bases: Gaussian

Gauss Legendre quadrature rule in torch for any domain [a,b]. See https://en.wikipedia.org/wiki/Gaussian_quadrature#Gauss%E2%80%93Legendre_quadrature.

Examples

>>> gl=torchquad.GaussLegendre()
>>> integral = gl.integrate(lambda x:np.sin(x), dim=1, N=101, integration_domain=[[0,5]]) #integral from 0 to 5 of np.sin(x)
|TQ-INFO| Computed integral was 0.7163378000259399 #analytic result = 1-np.cos(5)
class torchquad.Gaussian[source]

Bases: GridIntegrator

Base method for Gaussian Quadrature. Different Gaussian methods should inherit from this class, and override as necessary methods. Default behaviour is Gauss-Legendre quadrature on [-1,1] (i.e., this “parent” class should __not__ be used directly with other integration domains, and for this parent class integration_domain as an argument to integrate is ignored internally).

For an example of how to properly override the behavior to acheive different Gaussian Integration methods, please see the Custom Integrators section of the Tutorial or the implementation of GaussLegendre.

The primary methods/attributes of interest to override are _root_fn (for different polynomials, like numpy.polynomial.legendre.leggauss), _apply_composite_rule (as in other integration methods), and _resize_roots (for handling different integration domains).

name

A human-readable name for the integral.

Type:

str

_root_fn

A function that returns roots and weights like numpy.polynomial.legendre.leggauss.

Type:

function

_root_args

a way of adding information to be passed into _root_fn as needed. This is then used when caching roots/weights to potentially distinguish different calls to _root_fn based on arguments.

Type:

tuple

_cache

a cache for roots and weights, used internally.

Type:

dict

integrate(fn, dim, N=8, integration_domain=None, backend=None)[source]

Integrates the passed function on the passed domain using a Gaussian rule (Gauss-Legendre on [-1,1] as a default).

Parameters:
  • fn (func) – The function to integrate over.

  • dim (int) – Dimensionality of the integration domain.

  • N (int, optional) – Total number of sample points to use for the integration. Should be odd. Defaults to 3 points per dimension if None is given.

  • integration_domain (list or backend tensor, optional) – Integration domain, e.g. [[-1,1],[0,1]]. Defaults to [-1,1]^dim. It also determines the numerical backend if possible.

  • backend (string, optional) – Numerical backend. This argument is ignored if the backend can be inferred from integration_domain. Defaults to the backend from the latest call to set_up_backend or “torch” for backwards compatibility.

Returns:

Integral value

Return type:

backend-specific number

class torchquad.GridIntegrator[source]

Bases: BaseIntegrator

The abstract integrator that grid-like integrators (Newton-Cotes and Gaussian) integrators inherit from

calculate_grid(N, integration_domain, disable_integration_domain_check=False)[source]

Calculate grid points, widths and N per dim

Parameters:
  • N (int) – Number of points

  • integration_domain (backend tensor) – Integration domain

  • disable_integration_domain_check (bool) – Disbaling integration domain checks (default False)

Returns:

Grid points backend tensor: Grid widths int: Number of grid slices per dimension

Return type:

backend tensor

calculate_result(**kwargs)
get_jit_compiled_integrate(dim, N=None, integration_domain=None, backend=None)[source]

Create an integrate function where the performance-relevant steps except the integrand evaluation are JIT compiled. Use this method only if the integrand cannot be compiled. The compilation happens when the function is executed the first time. With PyTorch, return values of different integrands passed to the compiled function must all have the same format, e.g. precision.

Parameters:
  • dim (int) – Dimensionality of the integration domain.

  • N (int, optional) – Total number of sample points to use for the integration. See the integrate method documentation for more details.

  • integration_domain (list or backend tensor, optional) – Integration domain, e.g. [[-1,1],[0,1]]. Defaults to [-1,1]^dim. It can also determine the numerical backend.

  • backend (string, optional) – Numerical backend. Defaults to integration_domain’s backend if it is a tensor and otherwise to the backend from the latest call to set_up_backend or “torch” for backwards compatibility.

Returns:

JIT compiled integrate function where all parameters except the integrand and domain are fixed

Return type:

function(fn, integration_domain)

integrate(fn, dim, N, integration_domain, backend)[source]

Integrate the passed function on the passed domain using a Composite Newton Cotes rule. The argument meanings are explained in the sub-classes.

Returns:

integral value

Return type:

float

class torchquad.IntegrationGrid(N, integration_domain, grid_func=<function grid_func>, disable_integration_domain_check=False)[source]

Bases: object

This class is used to store the integration grid for methods like Trapezoid or Simpsons, which require a grid.

h = None
points = None
class torchquad.MonteCarlo[source]

Bases: BaseIntegrator

Monte Carlo integration

calculate_result(**kwargs)
calculate_sample_points(N, integration_domain, seed=None, rng=None)[source]

Calculate random points for the integrand evaluation

Parameters:
  • N (int) – Number of points

  • integration_domain (backend tensor) – Integration domain

  • seed (int, optional) – Random number generation seed for the sampling point creation, only set if provided. Defaults to None.

  • rng (RNG, optional) – An initialised RNG; this can be used when compiling the function for Tensorflow

Returns:

Sample points

Return type:

backend tensor

get_jit_compiled_integrate(dim, N=1000, integration_domain=None, seed=None, backend=None)[source]

Create an integrate function where the performance-relevant steps except the integrand evaluation are JIT compiled. Use this method only if the integrand cannot be compiled. The compilation happens when the function is executed the first time. With PyTorch, return values of different integrands passed to the compiled function must all have the same format, e.g. precision.

Parameters:
  • dim (int) – Dimensionality of the integration domain.

  • N (int, optional) – Number of sample points to use for the integration. Defaults to 1000.

  • integration_domain (list or backend tensor, optional) – Integration domain, e.g. [[-1,1],[0,1]]. Defaults to [-1,1]^dim. It can also determine the numerical backend.

  • seed (int, optional) – Random number generation seed for the sequence of sampling point calculations, only set if provided. The returned integrate function calculates different points in each invocation with and without specified seed. Defaults to None.

  • backend (string, optional) – Numerical backend. Defaults to integration_domain’s backend if it is a tensor and otherwise to the backend from the latest call to set_up_backend or “torch” for backwards compatibility.

Returns:

JIT compiled integrate function where all parameters except the integrand and domain are fixed

Return type:

function(fn, integration_domain)

integrate(fn, dim, N=1000, integration_domain=None, seed=None, rng=None, backend=None)[source]

Integrates the passed function on the passed domain using vanilla Monte Carlo Integration.

Parameters:
  • fn (func) – The function to integrate over.

  • dim (int) – Dimensionality of the function’s domain over which to integrate.

  • N (int, optional) – Number of sample points to use for the integration. Defaults to 1000.

  • integration_domain (list or backend tensor, optional) – Integration domain, e.g. [[-1,1],[0,1]]. Defaults to [-1,1]^dim. It can also determine the numerical backend.

  • seed (int, optional) – Random number generation seed to the sampling point creation, only set if provided. Defaults to None.

  • rng (RNG, optional) – An initialised RNG; this can be used when compiling the function for Tensorflow

  • backend (string, optional) – Numerical backend. Defaults to integration_domain’s backend if it is a tensor and otherwise to the backend from the latest call to set_up_backend or “torch” for backwards compatibility.

Raises:

ValueError – If len(integration_domain) != dim

Returns:

Integral value

Return type:

backend-specific number

class torchquad.RNG(backend, seed=None, torch_save_state=False)[source]

Bases: object

A random number generator helper class for multiple numerical backends

Notes

jax_get_key()[source]

Get the current PRNGKey. This function is needed for non-determinism when JIT-compiling with JAX.

jax_set_key(key)[source]

Set the PRNGKey. This function is needed for non-determinism when JIT-compiling with JAX.

uniform(size, dtype)[source]

Generate uniform random numbers in [0, 1) for the given numerical backend. This function is backend-specific; its definitions are in the constructor.

Parameters:
  • size (list) – The shape of the generated numbers tensor

  • dtype (backend dtype) – The dtype for the numbers, e.g. torch.float32

Returns:

A tensor with random values for the given numerical backend

Return type:

backend tensor

class torchquad.Simpson[source]

Bases: NewtonCotes

Simpson’s rule. See https://en.wikipedia.org/wiki/Newton%E2%80%93Cotes_formulas#Closed_Newton%E2%80%93Cotes_formulas .

integrate(fn, dim, N=None, integration_domain=None, backend=None)[source]

Integrates the passed function on the passed domain using Simpson’s rule.

Parameters:
  • fn (func) – The function to integrate over.

  • dim (int) – Dimensionality of the integration domain.

  • N (int, optional) – Total number of sample points to use for the integration. Should be odd. Defaults to 3 points per dimension if None is given.

  • integration_domain (list or backend tensor, optional) – Integration domain, e.g. [[-1,1],[0,1]]. Defaults to [-1,1]^dim. It can also determine the numerical backend.

  • backend (string, optional) – Numerical backend. Defaults to integration_domain’s backend if it is a tensor and otherwise to the backend from the latest call to set_up_backend or “torch” for backwards compatibility.

Returns:

Integral value

Return type:

backend-specific number

class torchquad.Trapezoid[source]

Bases: NewtonCotes

Trapezoidal rule. See https://en.wikipedia.org/wiki/Newton%E2%80%93Cotes_formulas#Closed_Newton%E2%80%93Cotes_formulas .

integrate(fn, dim, N=1000, integration_domain=None, backend=None)[source]

Integrates the passed function on the passed domain using the trapezoid rule.

Parameters:
  • fn (func) – The function to integrate over.

  • dim (int) – Dimensionality of the function to integrate.

  • N (int, optional) – Total number of sample points to use for the integration. Defaults to 1000.

  • integration_domain (list or backend tensor, optional) – Integration domain, e.g. [[-1,1],[0,1]]. Defaults to [-1,1]^dim. It can also determine the numerical backend.

  • backend (string, optional) – Numerical backend. Defaults to integration_domain’s backend if it is a tensor and otherwise to the backend from the latest call to set_up_backend or “torch” for backwards compatibility.

Returns:

Integral value

Return type:

backend-specific number

class torchquad.VEGAS[source]

Bases: BaseIntegrator

VEGAS Enhanced. Refer to https://arxiv.org/abs/2009.05112 . Implementation inspired by https://github.com/ycwu1030/CIGAR/ . EQ <n> refers to equation <n> in the above paper. JAX and Tensorflow are unsupported. For Tensorflow there exists a VEGAS+ implementation called VegasFlow: https://github.com/N3PDF/vegasflow

integrate(fn, dim, N=10000, integration_domain=None, seed=None, rng=None, use_grid_improve=True, eps_rel=0, eps_abs=0, max_iterations=20, use_warmup=True, backend=None)[source]

Integrates the passed function on the passed domain using VEGAS.

If the integrand output is far away from zero, i.e. lies within [b, b+c] for a constant b with large absolute value and small constant c, VEGAS does not adapt well to the integrand. Shifting the integrand so that it is close to zero may improve the accuracy of the calculated integral in this case. This method does not support multi-dimensional/vectorized integrands (i.e., integrating an integrand repeatedly over a grid of points).

Parameters:
  • fn (func) – The function to integrate over.

  • dim (int) – Dimensionality of the function’s domain over which to integrate.

  • N (int, optional) – Approximate maximum number of function evaluations to use for the integration. This value can be exceeded if the vegas stratification distributes evaluations per hypercube very unevenly. Defaults to 10000.

  • integration_domain (list, optional) – Integration domain, e.g. [[-1,1],[0,1]]. Defaults to [-1,1]^dim.

  • seed (int, optional) – Random number generation seed for the sampling point creation; only set if provided. Defaults to None.

  • rng (RNG, optional) – An initialised RNG; this can be used as alternative to the seed argument and to avoid problems with integrand functions which reset PyTorch’s RNG seed.

  • use_grid_improve (bool, optional) – If True, improve the vegas map after each iteration. Defaults to True.

  • eps_rel (float, optional) – Relative error to abort at. Defaults to 0.

  • eps_abs (float, optional) – Absolute error to abort at. Defaults to 0.

  • max_iterations (int, optional) – Maximum number of vegas iterations to perform. The number of performed iterations is usually lower than this value because the number of sample points per iteration increases every fifth iteration. Defaults to 20.

  • use_warmup (bool, optional) – If True, execute a warmup to initialize the vegas map. Defaults to True.

  • backend (string, optional) – Numerical backend. “jax” and “tensorflow” are unsupported. Defaults to integration_domain’s backend if it is a tensor and otherwise to the backend from the latest call to set_up_backend or “torch” for backwards compatibility.

Raises:

ValueError – If the integration_domain or backend argument is invalid

Returns:

Integral value

Return type:

backend-specific float

torchquad.enable_cuda(data_type='float32')[source]

This function sets torch’s default device to CUDA if possible. Call before declaring any variables! The default precision can be set here initially, or using set_precision later.

Parameters:

data_type ("float32", "float64" or None, optional) – Data type to use. If None, skip the call to set_precision. Defaults to “float32”.

torchquad.plot_convergence(evals, fvals, ground_truth, labels, dpi=150)[source]

Plots errors vs. function evaluations (fevals) and shows the convergence rate.

Parameters:
  • evals (list of np.array) – Number of evaluations, for each method a np.array of ints.

  • fvals (list of np.array) – Function values for evals.

  • ground_truth (np.array) – Ground truth values.

  • labels (list) – Method names.

  • dpi (int, optional) – Plot dpi. Defaults to 150.

torchquad.plot_runtime(evals, runtime, labels, dpi=150, y_axis_name='Runtime [s]')[source]

Plots the runtime vs. function evaluations (fevals).

Parameters:
  • evals (list of np.array) – Number of evaluations, for each method a np.array of fevals.

  • runtime (list of np.array) – Runtime for evals.

  • labels (list) – Method names.

  • dpi (int, optional) – Plot dpi. Defaults to 150.

  • y_axis_name (str, optional) – Name for y axis. Deafults to “Runtime [s]”.

torchquad.set_log_level(log_level: str)[source]

Set the log level for the logger. The preset log level when initialising Torchquad is the value of the TORCHQUAD_LOG_LEVEL environment variable, or ‘WARNING’ if the environment variable is unset.

Parameters:

log_level (str) – The log level to set. Options are ‘TRACE’,’DEBUG’, ‘INFO’, ‘SUCCESS’, ‘WARNING’, ‘ERROR’, ‘CRITICAL’

torchquad.set_precision(data_type='float32', backend='torch')[source]

This function allows the user to set the default precision for floating point numbers for the given numerical backend. Call before declaring your variables. NumPy and Tensorflow don’t have global dtypes: https://github.com/numpy/numpy/issues/6860 https://github.com/tensorflow/tensorflow/issues/26033 Therefore, torchquad sets the dtype argument for these two when initialising the integration domain.

Parameters:
  • data_type (string, optional) – Data type to use, either “float32” or “float64”. Defaults to “float32”.

  • backend (string, optional) – Numerical backend for which the data type is changed. Defaults to “torch”.

torchquad.set_up_backend(backend, data_type=None, torch_enable_cuda=True)[source]

Configure a numerical backend for torchquad.

With the torch backend, this function calls torchquad.enable_cuda unless torch_enable_cuda is False. With the tensorflow backend, this function enables tensorflow’s numpy behaviour, which is a requirement for torchquad. If a data type is passed, set the default floating point precision with torchquad.set_precision.

Parameters:
  • backend (string) – Numerical backend, e.g. “torch”

  • data_type ("float32", "float64" or None, optional) – Data type which is passed to set_precision. If None, do not call set_precision except if CUDA is enabled for torch. Defaults to None.

  • torch_enable_cuda (Bool, optional) – If True and backend is “torch”, call enable_cuda. Defaults to True.

class torchquad.integration.newton_cotes.NewtonCotes[source]

Bases: GridIntegrator

The abstract integrator that Composite Newton Cotes integrators inherit from

class torchquad.integration.base_integrator.BaseIntegrator[source]

Bases: object

The (abstract) integrator that all other integrators inherit from. Provides no explicit definitions for methods.

static evaluate_integrand(fn, points, weights=None, args=None)[source]

Evaluate the integrand function at the passed points

Parameters:
  • fn (function) – Integrand function

  • points (backend tensor) – Integration points

  • weights (backend tensor, optional) – Integration weights. Defaults to None.

  • args (list or tuple, optional) – Any arguments required by the function. Defaults to None.

Returns:

Integrand function output int: Number of evaluated points

Return type:

backend tensor

integrate()[source]