Using OpenQuad#

This guide describes the basic steps for using the key functionality of OpenQuad, illustrated with a concrete example. If you are not sure how to use this package, this guide is a good starting point.

Initialization#

As a first step, import the class for the desired geometry. For example, openquad.S2 for the integral over the surface of the 2d unit sphere, \(\mathrm{S}^2\):

>>> import numpy as np
>>> from openquad import S2

We will also use NumPy later on.

Initiate the class with the method specification (see S2) in the form of a list of tuples, containing the a string with the name of the method and a dictionary with required parameters.

Either use a method tailored to this geometry

>>> quad = S2([('S2-Design-Womersley', degree=7)])

or combine lower-dimensional methods with compatible geometries

>>> quad = S2([
...     ('GaussLegendre', dict(degree=7))]),
...     ('Trapezoid', dict(size=8)),
... ], polar_sampling='cos')

Some classes take additional parameters, like polar_sampling for openquad.S2.

Common functionality#

All geometry classes share a set of common attributes and functions, supplemented by functionality specific to certain geometries.

Sample points and weights#

Sample points and weights are stored as NumPy arrays in points and weights, respecitvely.

>>> quad.weights.shape
(32,)
>>> quad.points.shape
(2, 32)

Weights are normalized such that their sum equals the volume of the integration domain. For \(\mathrm{S}^2\), this is the area \(4\pi\):

>>> np.sum(quad.weights) / (4 * np.pi)
np.float(1.0)

Sample points in points are given in the default coordinates of the selected integration domain. For S2, these are spherical polar angles. Other coordinates might be available, e.g. angles or xyz.

>>> np.array_equal(quad.angles, quad.points)
True
>>> quad.xyz.shape
(3, 32)

Exporting quadratures#

You can save quadrature points and weights as a textfile with savetxt().

>>> quad.savetxt('points_and_weights.txt')

Integration#

Each class is equipped with the integrate() function, which can handle arrays and Python callables.

Suppose the integrand \(f(x)\) is a Python function, e.g.

>>> def f(theta, phi):
...     """Spherical harmonic |Y_{2, 1}|^2."""
...     return ((np.sin(2*theta) * np.cos(phi))*np.sqrt(15/(16*np.pi)))**2

To perform the integral of this function directly

>>> quad.integrate(f)
np.float64(1.0000000000000002)

In some situations it may be desirable or necessary to access the function values available on the quadrature grid.

>>> f_values = f(*quad.angles)
>>> f_values.shape
(32,)

You can perform the integration on the array data at a later point with

>>> quad.integrate(f_values)
np.float64(1.0000000000000002)

Other parameters#

Other attributes that are available for all top-level classes include:

  • dim: the dimension of the domain \(\mathcal{D}\).

  • size: the number of sample points.

  • shape: the shape of points.

See the API reference for details.