For the simplest case, one will create an anonymous inner
class implementing a method named function
as follows:
To integrate, one then callsGLQuadrature glq = new GLQuadrature(3) { protected double function(double u) { return u*u + 3; } }
glq.integrate(a, b)
where
a
and b
are the limits. Alternatively one
may call glq.integrate(a, b, m)
where m 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 two-argument "function" 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) { protected double function(double u, Data data) { return u*u + data.value; } }
or... Data data; data.value = 3.0; glqp.setParameter(data); System.out.println(glqp.integrate(a, b));
... Data data; data.value = 3.0; System.out.println(glqp.integrateWithP(a, b, data));
In addition, the methods
getArguments(double,double)
and
integrate(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) { 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 double function(double u) { double u5 = u5cache[(int)(Math.round(u))]; Data data = getParameters(); return u5 * data.value; } }; double[] args = glq.getArguments(a, b); Data data = new Data(); glq.setParameters(data); for (double value: list) { data.value = value; System.our.println(glq.integrate(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.RealValuedFunctOps, int)
will construct an instance of
GLQuadrature where the function is represented as an instance of
RealValuedFunction
. While the implementation is trivial in Java,
newInstance
simplifies the use of this class
from a scripting language. For example,
shows how to usefs = {valueAt: function(u) {return Math.sin(u);}} rvf = new RealValuedFunction(scripting, fs); glq = GLQuadrature.newInstance(rvf, 16); integral = glq.integrate(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 TypeMethodDescriptionprotected double
function
(double t) The function to integrate.protected double
The function (with parameters) to integrate.double[]
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.static double[]
getWeights
(double a, double b, int n) Get weights for explicit integration.double
integrate
(double[] arguments) Integrate a function using precomputed arguments.double
integrate
(double a, double b) Integrate the function.double
integrate
(double a, double b, int m) Integrate the function using multiple subintervals.void
Integrate a function using precomputed arguments, accumulating the integral in an adder.void
Integrate the function, accumulating the integral in an adder.void
Integrate the function using multiple subintervals, accumulating the integral in an adder.double
integrateWithP
(double[] arguments, P p) Integrate a function using precomputed arguments with parameters.double
integrateWithP
(double a, double b, int m, P p) Integrate the function with explicit parameters using multiple subintervals.double
integrateWithP
(double a, double b, P p) Integrate the function with explicit parameters.void
integrateWithP
(Adder adder, double[] arguments, P p) Integrate a function using precomputed arguments with parameters, accumulating the integral in an adder.void
integrateWithP
(Adder adder, double a, double b, int m, P p) Integrate the function with explicit parameters using multiple subintervals, accumulating the integral in an adder.void
integrateWithP
(Adder adder, double a, double b, P p) Integrate the function with explicit parameters, accumulating the integral in an adder.static GLQuadrature
newInstance
(RealValuedFunctOps f, int n) Create a new instance of GLQuadrature that uses an instance of RealValuedFunction or RealValueFunctOps as its function.void
setParameters
(P parameters) Set parameters.
-
Constructor Details
-
GLQuadrature
public GLQuadrature(int n) 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 use- 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)
-
function
protected double function(double t) The function to integrate. The default implementation calls the two-argument implementation offunction(double,Object)
using the parameters supplied by a previous call tosetParameters(P)
. If parameters are not used, this method may be overridden.- Parameters:
t
- the function's argument- Returns:
- the value of the function given its argument
-
function
The function (with parameters) to integrate. Typically, this function 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:
t
- the function's argumentp
- the parameters- Returns:
- the value of the function given its arguments
- 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
-
integrate
public double integrate(double a, double b) Integrate the function.- Parameters:
a
- the lower limit of integrationb
- the upper limit of integration- Returns:
- the definite integral from a to b
-
integrate
Integrate the function, accumulating the integral in an adder.- Parameters:
adder
- the adder used to accumulate integralsa
- the lower limit of integrationb
- the upper limit of integration
-
integrateWithP
Integrate the function with explicit parameters.- Parameters:
a
- the lower limit of integrationb
- the upper limit of integrationp
- the parameters- Returns:
- the definite integral from a to b
-
integrateWithP
Integrate the function with explicit parameters, accumulating the integral in an adder.- Parameters:
adder
- the addera
- the lower limit of integrationb
- the upper limit of integrationp
- the parameters
-
getArguments
public double[] getArguments(double a, double b) Set up data for an integration. This method is used in conjunction withintegrate(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[])
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 double integrate(double[] arguments) Integrate a function 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.- Parameters:
arguments
- the array returned bygetArguments(a,b)
- Returns:
- the definite integral from a to b
- See Also:
-
integrate
Integrate a function 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:
adder
- the adderarguments
- the array returned bygetArguments(a,b)
- See Also:
-
integrateWithP
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:
arguments
- the array returned bygetArguments(a,b)
p
- the parameters- Returns:
- the definite integral from a to b
- See Also:
-
integrateWithP
Integrate a function using precomputed arguments with parameters, 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:
adder
- the adderarguments
- the array returned bygetArguments(a,b)
p
- the parameters- See Also:
-
integrate
public double integrate(double a, double b, int m) 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 not very large due to numerical accuracy issues in computing the values at a particular point for Legendre polynomials with a very high degree.- Parameters:
a
- the lower limit of integrationb
- the upper limit of integrationm
- the number of subintervals- Returns:
- the definite integral from a to b
-
integrate
Integrate the function using multiple subintervals, accumulating the integral in an adder. 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:
adder
- the addera
- the lower limit of integrationb
- the upper limit of integrationm
- the number of subintervals
-
integrateWithP
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 not very large due to numerical accuracy issues in computing the values at a particular point for Legendre polynomials with a very high degree.- Parameters:
a
- the lower limit of integrationb
- the upper limit of integrationm
- the number of subintervalsp
- the parameters- Returns:
- the definite integral from a to b
-
integrateWithP
Integrate the function with explicit parameters using multiple subintervals, accumulating the integral in an adder. 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:
adder
- the addera
- the lower limit of integrationb
- the upper limit of integrationm
- the number of subintervalsp
- the parameters
-
newInstance
Create a new instance of GLQuadrature that uses an instance of RealValuedFunction or RealValueFunctOps as its function.- Parameters:
f
- the functionn
- the number of points to use- Returns:
- a new instance of GLQuadrature
- Throws:
IllegalArgumentException
- n is less than 1
-