n-point Gaussian quadrature will produce exact results when
integrating a polynomial of degree (2n-1). It computes the
integral over the range [-1, 1] as a weighted sum of the
value of the integrand at points specified by the algorithm.
Thus, it is not suitable for integrating a function whose value
is known at (for example) evenly spaced points. If, however, the
integrand can be computed at any desired point, Gaussian quadrature
requires the integrand be evaluated a relatively small number of
times. As an example of its use, an electrostatic field can be
computed from a charge density by evaluation the integral
$$ \mathbf{E}(\mathbf{x}_2) = \int \rho(\mathbf{x}_1
\frac{\mathbf{x}_2 - \mathbf{x}_1}{|\mathbf{x}_2 - \mathbf{x}_1|^3}
d\mathbf{x}_1$$
While one could compute three integrals, one for each component,
the term |x2-x1|3|,
which is common to all three integrals, would be computed multiple
times. Using VectorValuedGLQ
instead of GLQuadrature
for such cases can lead to simpler and faster-executing code.
For the simplest case, one will create an anonymous inner
class implementing a method named mapping
as follows
(this example sets an array of length m):
To integrate, one then callsGLQuadrature glq = new GLQuadrature(n,m) { protected void mapping(double[] value, int m, double u) { double uu3 = u*u+3; for (int i = 0; i < m; i++) { value[i] = uu3 + i*u; } }
glq.integrate(vector, a, b)
where vector is an array of length m, and
a
and b
are the limits. Alternatively one
may call glq.integrate(vector, a, b, m1)
where m1 is
the number of subintervals on which to use Gauss-Legendre quadrature
(this is not necessary in the example above because Gauss-Legendre
quadrature is exact for polynomials of sufficiently low degree).
The type parameter P, if specified, can be used to provide a class
that will store the values of parameters used in computing the function
to be integrated. These values can be set by calling
setParameters(P)
and read by calling getParameters()
.
Parameters may also be provided explicitly when an "integrate" method
is called. The four-argument "mapping" should be implemented in this
case:
When used, this allows one to define parameterized functions - the parameters act as additional arguments that are typically constant during the integration. A parameter is stored, not copied, so it should not be modified while an integral is being computed (not an issue for single-threaded programs). If a series of integrals are computed for different values of the parameters, this allows one instance of GLQuadrature to be used, rather than a separate object for each value. An example of usage is:GLQuadrature<Data> glqp = new GLQuadrature<Data>(3, m) { protected void mapping(double[] value, int m, double u, Data data) { for (int i = 0; i < m; i++) { value[i] = u*u + data.value; } } }
or... Data data; double[] values = new double[m]; data.value = 3.0; glqp.setParameter(data); glqp.integrate(value, a, b); System.out.println(glqp.integrate(a, b));
... Data data; data.value = 3.0; glqp.integrateWithP(values, a, b, data));
In addition, the methods
getArguments(double,double)
and
integrate(double[], double[])
can be used to precompute function arguments in cases where those
will be used repetitively (e.g., for a series of integrals where
the parameters change but not the range or the number of points).
For example, consider the following code:
As an optimization, the code caches the arguments raised to the fifth power to reduce the number of multiplications needed, assumingGLQuadrature<Data> glq = new GLQuadrature<Data>(3, m) { double[] u5cache = null; public double[] getArguments(double a, double b) { double[] results = super.getArguments(a, b); int n = results.length - 1; u5cache = new double[n]; for (i = 0; i < n; i++) { double u = results[i]; double uu = u*u; u5cache[i] = uu*uu*u; result[i] = (double)i; } } protected void mapping(double[] values, int m, double u) { double u5 = u5cache[(int)(Math.round(u))]; Data data = getParameters(); for (int i = 0; i < m; i++) { values[i] = u5*data.value; } } }; double[] values = new double[m]; double[] args = glq.getArguments(a, b); Data data = new Data(); glq.setParameters(data); for (double value: list) { data.value = value; glq.integrate(values, args); ... }
list
contains multiple elements. This sort of optimization
would typically be used only in special cases where sufficiently
expensive terms in the function to be integrated can be precomputed.
Finally, the method newInstance(org.bzdev.math.RealToVectorMap, int, int)
will construct an instance of
GLQuadrature where the mapping is represented by a lambda expression.
For example,
shows how to useglq = VectorValuedGLQ.newInstance((result, m, u) -> { for (int i = 0; i < 10; i++) { result[i] = Math.sin(u*u + i); } }, 16, 10); glq.integrate(integrals, 0.0, Math.PI);
newInstance
with ECMAScript.
Note: some of the methods are named integrate
while
those that use parameters explicitly are named integrateWithP
.
While Java's type system treats primitive types and class types
differently, a scripting language may not.
-
Constructor Summary
Constructors -
Method Summary
Modifier and TypeMethodDescriptiondouble[]
getArguments
(double a, double b) Set up data for an integration.static double[]
getArguments
(double a, double b, int n) Set up data for an integration using an Gauss Legendre quadrature with a specified number of points.int
Get the number of points for which a function is evaluatedGet parameters.int
Get the minimum length of an array of double precision numbers that will be used hold the results of an integration.static double[]
getWeights
(double a, double b, int n) Get weights for explicit integration.void
integrate
(double[] results, double[] arguments) Integrate functions using precomputed arguments.void
integrate
(double[] results, double a, double b) Integrate the functions.void
integrate
(double[] results, double a, double b, int ms) Integrate the function using multiple subintervals.void
Integrate functions using precomputed arguments, accumulating the integral in an adder.void
Integrate the function, accumulating the integral in an adder.void
Integrate functions using multiple subintervals, accumulating the integral in adders.void
integrateWithP
(double[] results, double[] arguments, P p) Integrate a function using precomputed arguments with parameters.void
integrateWithP
(double[] results, double a, double b, int ms, P p) Integrate the function with explicit parameters using multiple subintervals.void
integrateWithP
(double[] results, double a, double b, P p) Integrate the function with explicit parameters.void
integrateWithP
(Adder[] adders, double[] arguments, P p) Integrate functions using precomputed arguments with parameters, accumulating the integral in adders.void
integrateWithP
(Adder[] adders, double a, double b, int ms, P p) Integrate functions with explicit parameters using multiple subintervals, accumulating the integrals in adders.void
integrateWithP
(Adder[] adders, double a, double b, P p) Integrate functions with explicit parameters, accumulating the integrals in adders.protected void
mapping
(double[] results, int m, double t) A mapping that specifies the functions to integrate.protected void
A mapping that specifies the functions (with parameters) to integrate.static VectorValuedGLQ
newInstance
(RealToVectorMap map, int n, int m) Create a new instance of VectorValuedGLQ that uses an instance of RealToVectorMap to map a value to the components of a vector.void
setParameters
(P parameters) Set parameters.
-
Constructor Details
-
VectorValuedGLQ
public VectorValuedGLQ(int n, int m) Constructor. Note: to determine weights used in the integration, the roots of Legendre polynomials or order n have to be calculated, and the accuracy decreases for large values of n: for n = 1000, for example, the value at a computed root may be around 10^-11 whereas for small values of n, it is around 10^-15, so the accuracy of the weights decreases for very large values of n. As a result, increasing n to excessively large values can be counterproductive beyond merely increasing computation time needlessly.- Parameters:
n
- the number of points to usem
- the number of values that will be provided as the result of an integration- Throws:
IllegalArgumentException
- n is less than 1
-
-
Method Details
-
setParameters
Set parameters. Parameters are used to provide additional values for computing the function being integrated. The argument is stored, not copied, and will typically not be changed while an integration is in progress.- Parameters:
parameters
- the parameters
-
getParameters
Get parameters. Parameters are used to provide additional values for computing the function being integrated.- Returns:
- an instance of the class representing a root finder's parameters (this will be the same instance passed to setParameters)
-
mapping
protected void mapping(double[] results, int m, double t) A mapping that specifies the functions to integrate. The default implementation calls the five-argument implementationmapping(double[],int,double,Object)
using the parameters supplied by a previous call tosetParameters(P)
. If parameters are* not used, this method may be overridden.- Parameters:
results
- an array holding the results, with the minimum width set by the constructor and read bygetResultLength()
m
- the number of elements in the results array that are usedt
- the argument for the integration
-
mapping
A mapping that specifies the functions (with parameters) to integrate. Typically, this mapping will be defined via an anonymous subclass. If the argument p is null, one will normally callgetParameters()
to obtain any parameters that were configured via thesetParameters(P)
method. The methods named integrateWithP will use this method.- Parameters:
results
- an array holding the results, with the minimum width set by the constructor and read bygetResultLength()
m
- the number of elements in the results array that are usedt
- the argument for the integrationp
- the parameters- Throws:
UnsupportedOperationException
- the method was needed but not implemented.
-
getNumberOfPoints
public int getNumberOfPoints()Get the number of points for which a function is evaluated- Returns:
- the number of points
-
getResultLength
public int getResultLength()Get the minimum length of an array of double precision numbers that will be used hold the results of an integration.- Returns:
- the length specified in the constructor
-
integrate
public void integrate(double[] results, double a, double b) throws IllegalArgumentException, UnsupportedOperationException Integrate the functions.- Parameters:
results
- an array holding the results, with the minimum width set by the constructor and read bygetResultLength()
a
- the lower limit of integrationb
- the upper limit of integration- Throws:
IllegalArgumentException
- the first argument was null or an array that was too shortUnsupportedOperationException
- a mapping was not implemented
-
integrate
public void integrate(Adder[] adders, double a, double b) throws IllegalArgumentException, UnsupportedOperationException Integrate the function, accumulating the integral in an adder.- Parameters:
adders
- the adders used to accumulate integralsa
- the lower limit of integrationb
- the upper limit of integration- Throws:
IllegalArgumentException
- the first argument was null or an array that was too shortUnsupportedOperationException
- a mapping was not implemented
-
integrateWithP
public void integrateWithP(double[] results, double a, double b, P p) throws IllegalArgumentException, UnsupportedOperationException Integrate the function with explicit parameters.- Parameters:
results
- an array holding the results, with the minimum width set by the constructor and read bygetResultLength()
a
- the lower limit of integrationb
- the upper limit of integrationp
- the parameters- Throws:
IllegalArgumentException
- the first argument was null or an array that was too shortUnsupportedOperationException
- a mapping was not implemented
-
integrateWithP
public void integrateWithP(Adder[] adders, double a, double b, P p) throws IllegalArgumentException, UnsupportedOperationException Integrate functions with explicit parameters, accumulating the integrals in adders.- Parameters:
adders
- the addersa
- the lower limit of integrationb
- the upper limit of integrationp
- the parameters- Throws:
IllegalArgumentException
- the first argument was null or an array that was too shortUnsupportedOperationException
- a mapping was not implemented
-
getArguments
public double[] getArguments(double a, double b) Set up data for an integration. This method is used in conjunction withintegrate(double[],double[])
for cases where a function may be integrated repeatedly over the same range, but with varying parameters.- Parameters:
a
- the lower limit of integrationb
- the upper limit of integration- Returns:
- an array of values at which the function should be evaluated, with the last element containing half the range of integration (the array length is one larger than the number of points provided to the constructor)
- See Also:
-
getArguments
public static double[] getArguments(double a, double b, int n) Set up data for an integration using an Gauss Legendre quadrature with a specified number of points. This method is used in conjunction withintegrate(double[],double[])
for cases where a function may be integrated repeatedly over the same range, but with varying parameters. The caller must ensure that the value passed as this method's final argument matches the value used in the constructor for the GLQuadrature instance used.- Parameters:
a
- the lower limit of integrationb
- the upper limit of integrationn
- the number of points at which the function will be evaluated- Returns:
- an array of values at which the function should be evaluated, with the last element containing half the range of integration (the array length is one larger than the number of points provided to the constructor)
- Throws:
IllegalArgumentException
- n is less than 1- See Also:
-
getWeights
public static double[] getWeights(double a, double b, int n) Get weights for explicit integration. To compute an integral from a to b explicitly (e.g., without creating an instance of GLQuadrature), one can use the following sequence of statements:
The weights returned by this method include a scaling factor for making the integral go from a to b instead of -1 to 1.double[] arguments = GLQuadrature.getArguments(a, b, n); double[] weights = GLQuadrature.getWeights(a, b, n); Adder adder = new Adder.Kahan(); for (int i = 0; i < n; i++) { double u = arguments[i]; ... y = [ an expression that depends on u ] ; adder.add(weights[i]*y])); } integral = adder.getSum();
- Parameters:
a
- the lower limit of integrationb
- the upper limit of integrationn
- the number of points at which a function will- Returns:
- an array of weights
- Throws:
IllegalArgumentException
- n is less than 1 be evaluated.
-
integrate
public void integrate(double[] results, double[] arguments) throws IllegalArgumentException, UnsupportedOperationException Integrate functions using precomputed arguments. This method is used in conjunction withgetArguments(double,double)
for cases where a function may be integrated repeatedly over the same range, but with varying parameters. It is provided as an optimization for a special case. Each integral's range is from a to b, where a and b are arguments togetArguments(double,double)
.- Parameters:
results
- an array holding the results, with the minimum width set by the constructor and read bygetResultLength()
arguments
- the array returned bygetArguments(a,b)
- Throws:
IllegalArgumentException
- the first argument was null or an array that was too shortUnsupportedOperationException
- a mapping was not implemented- See Also:
-
integrate
public void integrate(Adder[] adders, double[] arguments) throws IllegalArgumentException, UnsupportedOperationException Integrate functions using precomputed arguments, accumulating the integral in an adder. This method is used in conjunction withgetArguments(double,double)
for cases where a function may be integrated repeatedly over the same range, but with varying parameters. It is provided as an optimization for a special case.- Parameters:
adders
- the adders that will hold the resultsarguments
- the array returned bygetArguments(a,b)
- Throws:
IllegalArgumentException
- the first argument was null or an array that was too shortUnsupportedOperationException
- a mapping was not implemented- See Also:
-
integrateWithP
public void integrateWithP(double[] results, double[] arguments, P p) throws IllegalArgumentException, UnsupportedOperationException Integrate a function using precomputed arguments with parameters. This method is used in conjunction withgetArguments(double,double)
for cases where a function may be integrated repeatedly over the same range, but with varying parameters. It is provided as an optimization for a special case.- Parameters:
results
- an array holding the results, with the minimum width set by the constructor and read bygetResultLength()
arguments
- the array returned bygetArguments(a,b)
p
- the parameters- Throws:
IllegalArgumentException
- the first argument was null or an array that was too shortUnsupportedOperationException
- a mapping was not implemented- See Also:
-
integrateWithP
public void integrateWithP(Adder[] adders, double[] arguments, P p) throws IllegalArgumentException, UnsupportedOperationException Integrate functions using precomputed arguments with parameters, accumulating the integral in adders. This method is used in conjunction withgetArguments(double,double)
for cases where a function may be integrated repeatedly over the same range, but with varying parameters. It is provided as an optimization for a special case.- Parameters:
adders
- the addersarguments
- the array returned bygetArguments(a,b)
p
- the parameters- Throws:
IllegalArgumentException
- the first argument was null or an array that was too shortUnsupportedOperationException
- a mapping was not implemented- See Also:
-
integrate
public void integrate(double[] results, double a, double b, int ms) throws IllegalArgumentException, UnsupportedOperationException Integrate the function using multiple subintervals. The interval [a,b] is divided into m intervals, each of which is integrated separately and summed. This is useful when the argument n passed to the constructor is very large due to numerical accuracy issues in computing the values at a particular point for Legendre polynomials with a very high degree.- Parameters:
results
- an array holding the results, with the minimum width set by the constructor and read bygetResultLength()
a
- the lower limit of integrationb
- the upper limit of integrationms
- the number of subintervals- Throws:
IllegalArgumentException
- the first argument was null or an array that was too shortUnsupportedOperationException
- a mapping was not implemented
-
integrate
public void integrate(Adder[] adders, double a, double b, int ms) throws IllegalArgumentException, UnsupportedOperationException Integrate functions using multiple subintervals, accumulating the integral in adders. The interval [a,b] is divided into m intervals, each of which is integrated separately and summed. This is useful when the argument n passed to the constructor is very large due to numerical accuracy issues in computing the values at a particular point for Legendre polynomials with a very high degree.- Parameters:
adders
- the addersa
- the lower limit of integrationb
- the upper limit of integrationms
- the number of subintervals- Throws:
IllegalArgumentException
- the first argument was null or an array that was too shortUnsupportedOperationException
- a mapping was not implemented
-
integrateWithP
public void integrateWithP(double[] results, double a, double b, int ms, P p) throws IllegalArgumentException, UnsupportedOperationException Integrate the function with explicit parameters using multiple subintervals. The interval [a,b] is divided into m intervals, each of which is integrated separately and summed. This is useful when the argument n passed to the constructor is very large due to numerical accuracy issues in computing the values at a particular point for Legendre polynomials with a very high degree.- Parameters:
results
- an array holding the results, with the minimum width set by the constructor and read bygetResultLength()
a
- the lower limit of integrationb
- the upper limit of integrationms
- the number of subintervalsp
- the parameters- Throws:
IllegalArgumentException
- the first argument was null or an array that was too shortUnsupportedOperationException
- a mapping was not implemented
-
integrateWithP
public void integrateWithP(Adder[] adders, double a, double b, int ms, P p) throws IllegalArgumentException, UnsupportedOperationException Integrate functions with explicit parameters using multiple subintervals, accumulating the integrals in adders. The interval [a,b] is divided into m intervals, each of which is integrated separately and summed. This is useful when the argument n passed to the constructor is very large due to numerical accuracy issues in computing the values at a particular point for Legendre polynomials with a very high degree.- Parameters:
adders
- the addersa
- the lower limit of integrationb
- the upper limit of integrationms
- the number of subintervalsp
- the parameters- Throws:
IllegalArgumentException
- the first argument was null or an array that was too shortUnsupportedOperationException
- a mapping was not implemented
-
newInstance
Create a new instance of VectorValuedGLQ that uses an instance of RealToVectorMap to map a value to the components of a vector.- Parameters:
map
- the mapn
- the number of points to usem
- the number of elements in the results array that are used- Returns:
- the new instance of VectorValuedGLQ
- Throws:
IllegalArgumentException
- n is less than 1
-