Skip to main content

Numerical integration, quadrature for various domains

Project description

PyPi Version PyPI pyversions DOI GitHub stars PyPi downloads

Slack awesome

gh-actions codecov LGTM Code style: black

More than 1500 numerical integration schemes for line segments, circles, disks, triangles, quadrilaterals, spheres, balls, tetrahedra, hexahedra, wedges, pyramids, n-spheres, n-balls, n-cubes, n-simplices, the 1D half-space with weight functions exp(-r), the 2D space with weight functions exp(-r), the 3D space with weight functions exp(-r), the nD space with weight functions exp(-r), the 1D space with weight functions exp(-r2), the 2D space with weight functions exp(-r2), the 3D space with weight functions exp(-r2), and the nD space with weight functions exp(-r2), for fast integration of real-, complex-, and vector-valued functions.

For example, to numerically integrate any function over any given interval, install quadpy from the Python Package Index with

pip install quadpy

and do

import numpy
import quadpy


def f(x):
    return numpy.sin(x) - x


val, err = quadpy.quad(f, 0.0, 6.0)

This is like scipy with the addition that quadpy handles complex-, vector-, matrix-valued integrands, and "intervals" in spaces of arbitrary dimension.

To integrate over a triangle, do

import numpy
import quadpy


def f(x):
    return numpy.sin(x[0]) * numpy.sin(x[1])


triangle = numpy.array([[0.0, 0.0], [1.0, 0.0], [0.7, 0.5]])

# get a "good" scheme of degree 10
scheme = quadpy.t2.get_good_scheme(10)
val = scheme.integrate(f, triangle)

Most domains have get_good_scheme(degree). If you would like to use a particular scheme, you can pick one from the dictionary quadpy.t2.schemes.

All schemes have

scheme.points
scheme.weights
scheme.degree
scheme.source
scheme.test_tolerance

scheme.show()
scheme.integrate(
    # ...
)

and many have

scheme.points_symbolic
scheme.weights_symbolic

quadpy is fully vectorized, so if you like to compute the integral of a function on many domains at once, you can provide them all in one integrate() call, e.g.,

# shape (3, 5, 2), i.e., (corners, num_triangles, xy_coords)
triangles = numpy.stack(
    [
        [[0.0, 0.0], [1.0, 0.0], [0.0, 1.0]],
        [[1.2, 0.6], [1.3, 0.7], [1.4, 0.8]],
        [[26.0, 31.0], [24.0, 27.0], [33.0, 28]],
        [[0.1, 0.3], [0.4, 0.4], [0.7, 0.1]],
        [[8.6, 6.0], [9.4, 5.6], [7.5, 7.4]],
    ],
    axis=-2,
)

The same goes for functions with vectorized output, e.g.,

def f(x):
    return [numpy.sin(x[0]), numpy.sin(x[1])]

More examples under test/examples_test.py.

Read more about the dimensionality of the input/output arrays in the wiki.

Advanced topics:

Schemes

Line segment (C1)

See here for how to generate Gauss formulas for your own weight functions.

Example:

import numpy
import quadpy

scheme = quadpy.c1.gauss_patterson(5)
scheme.show()
val = scheme.integrate(lambda x: numpy.exp(x), [0.0, 1.0])

1D half-space with weight function exp(-r) (E1r)

Example:

import quadpy

scheme = quadpy.e1r.gauss_laguerre(5, alpha=0)
scheme.show()
val = scheme.integrate(lambda x: x ** 2)

1D space with weight function exp(-r2) (E1r2)

Example:

import quadpy

scheme = quadpy.e1r2.gauss_hermite(5)
scheme.show()
val = scheme.integrate(lambda x: x ** 2)

Circle (U2)

  • Krylov (1959, arbitrary degree)

Example:

import numpy
import quadpy

scheme = quadpy.u2.get_good_scheme(7)
scheme.show()
val = scheme.integrate(lambda x: numpy.exp(x[0]), [0.0, 0.0], 1.0)

Triangle (T2)

Apart from the classical centroid, vertex, and seven-point schemes we have

Example:

import numpy
import quadpy

scheme = quadpy.t2.get_good_scheme(12)
scheme.show()
val = scheme.integrate(lambda x: numpy.exp(x[0]), [[0.0, 0.0], [1.0, 0.0], [0.5, 0.7]])

Disk (S2)

Example:

import numpy
import quadpy

scheme = quadpy.s2.get_good_scheme(6)
scheme.show()
val = scheme.integrate(lambda x: numpy.exp(x[0]), [0.0, 0.0], 1.0)

Quadrilateral (C2)

Example:

import numpy
import quadpy

scheme = quadpy.c2.get_good_scheme(7)
val = scheme.integrate(
    lambda x: numpy.exp(x[0]),
    [[[0.0, 0.0], [1.0, 0.0]], [[0.0, 1.0], [1.0, 1.0]]],
)

The points are specified in an array of shape (2, 2, ...) such that arr[0][0] is the lower left corner, arr[1][1] the upper right. If your c2 has its sides aligned with the coordinate axes, you can use the convenience function

quadpy.c2.rectangle_points([x0, x1], [y0, y1])

to generate the array.

2D space with weight function exp(-r) (E2r)

Example:

import quadpy

scheme = quadpy.e2r.get_good_scheme(5)
scheme.show()
val = scheme.integrate(lambda x: x[0] ** 2)

2D space with weight function exp(-r2) (E2r2)

Example:

import quadpy

scheme = quadpy.e2r2.get_good_scheme(3)
scheme.show()
val = scheme.integrate(lambda x: x[0] ** 2)

Sphere (U3)

Example:

import numpy
import quadpy

scheme = quadpy.u3.get_good_scheme(19)
# scheme.show()
val = scheme.integrate(lambda x: numpy.exp(x[0]), [0.0, 0.0, 0.0], 1.0)

Integration on the sphere can also be done for functions defined in spherical coordinates:

import numpy
import quadpy


def f(theta_phi):
    theta, phi = theta_phi
    return numpy.sin(phi) ** 2 * numpy.sin(theta)


scheme = quadpy.u3.get_good_scheme(19)
val = scheme.integrate_spherical(f)

Ball (S3)

Example:

import numpy
import quadpy

scheme = quadpy.s3.get_good_scheme(4)
# scheme.show()
val = scheme.integrate(lambda x: numpy.exp(x[0]), [0.0, 0.0, 0.0], 1.0)

Tetrahedron (T3)

Example:

import numpy
import quadpy

scheme = quadpy.t3.get_good_scheme(5)
# scheme.show()
val = scheme.integrate(
    lambda x: numpy.exp(x[0]),
    [[0.0, 0.0, 0.0], [1.0, 0.0, 0.0], [0.5, 0.7, 0.0], [0.3, 0.9, 1.0]],
)

Hexahedron (C3)

Example:

import numpy
import quadpy

scheme = quadpy.c3.product(quadpy.c1.newton_cotes_closed(3))
# scheme.show()
val = scheme.integrate(
    lambda x: numpy.exp(x[0]),
    quadpy.c3.cube_points([0.0, 1.0], [-0.3, 0.4], [1.0, 2.1]),
)

Pyramid (P3)

  • Felippa (2004, 9 schemes up to degree 5)

Example:

import numpy
import quadpy

scheme = quadpy.p3.felippa_5()

val = scheme.integrate(
    lambda x: numpy.exp(x[0]),
    [
        [0.0, 0.0, 0.0],
        [1.0, 0.0, 0.0],
        [0.5, 0.7, 0.0],
        [0.3, 0.9, 0.0],
        [0.0, 0.1, 1.0],
    ],
)

Wedge (W3)

Example:

import numpy
import quadpy

scheme = quadpy.w3.felippa_3()
val = scheme.integrate(
    lambda x: numpy.exp(x[0]),
    [
        [[0.0, 0.0, 0.0], [1.0, 0.0, 0.0], [0.5, 0.7, 0.0]],
        [[0.0, 0.0, 1.0], [1.0, 0.0, 1.0], [0.5, 0.7, 1.0]],
    ],
)

3D space with weight function exp(-r) (E3r)

Example:

import quadpy

scheme = quadpy.e3r.get_good_scheme(5)
# scheme.show()
val = scheme.integrate(lambda x: x[0] ** 2)

3D space with weight function exp(-r2) (E3r2)

Example:

import quadpy

scheme = quadpy.e3r2.get_good_scheme(6)
# scheme.show()
val = scheme.integrate(lambda x: x[0] ** 2)

n-Simplex (Tn)

Example:

import numpy
import quadpy

dim = 4
scheme = quadpy.tn.grundmann_moeller(dim, 3)
val = scheme.integrate(
    lambda x: numpy.exp(x[0]),
    numpy.array(
        [
            [0.0, 0.0, 0.0, 0.0],
            [1.0, 2.0, 0.0, 0.0],
            [0.0, 1.0, 0.0, 0.0],
            [0.0, 3.0, 1.0, 0.0],
            [0.0, 0.0, 4.0, 1.0],
        ]
    ),
)

n-Sphere (Un)

Example:

import numpy
import quadpy

dim = 4
scheme = quadpy.un.dobrodeev_1978(dim)
val = scheme.integrate(lambda x: numpy.exp(x[0]), numpy.zeros(dim), 1.0)

n-Ball (Sn)

Example:

import numpy
import quadpy

dim = 4
scheme = quadpy.sn.dobrodeev_1970(dim)
val = scheme.integrate(lambda x: numpy.exp(x[0]), numpy.zeros(dim), 1.0)

n-Cube (Cn)

Example:

import numpy
import quadpy

dim = 4
scheme = quadpy.cn.stroud_cn_3_3(dim)
val = scheme.integrate(
    lambda x: numpy.exp(x[0]),
    quadpy.cn.ncube_points([0.0, 1.0], [0.1, 0.9], [-1.0, 1.0], [-1.0, -0.5]),
)

nD space with weight function exp(-r) (Enr)

Example:

import quadpy

dim = 4
scheme = quadpy.enr.stroud_enr_5_4(dim)
val = scheme.integrate(lambda x: x[0] ** 2)

nD space with weight function exp(-r2) (Enr2)

Example:

import quadpy

dim = 4
scheme = quadpy.enr2.stroud_enr2_5_2(dim)
val = scheme.integrate(lambda x: x[0] ** 2)

Installation

quadpy is available from the Python Package Index, so with

pip install quadpy

you can install.

Testing

To run the tests, check out this repository and type

MPLBACKEND=Agg pytest

License

This software is published under the GPLv3 license.

Project details


Download files

Download the file for your platform. If you're not sure which to choose, learn more about installing packages.

Files for quadpy, version 0.16.4
Filename, size File type Python version Upload date Hashes
Filename, size quadpy-0.16.4-py3-none-any.whl (1.3 MB) File type Wheel Python version py3 Upload date Hashes View
Filename, size quadpy-0.16.4.tar.gz (1.1 MB) File type Source Python version None Upload date Hashes View

Supported by

Pingdom Pingdom Monitoring Google Google Object Storage and Download Analytics Sentry Sentry Error logging AWS AWS Cloud computing DataDog DataDog Monitoring Fastly Fastly CDN DigiCert DigiCert EV certificate StatusPage StatusPage Status page