Class VectorValuedGLQ<P>

java.lang.Object
org.bzdev.math.VectorValuedGLQ<P>

public abstract class VectorValuedGLQ<P> extends Object
Implementation of Gauss-Legendre quadrature that returns multiple integrals over the same range, with the results stored in an array.

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$$

E(x2) = ∫ ρ(x1) (x2-x1) /|x2-x1|3 dx1.
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):


      GLQuadrature 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;
        }
      }
 
To integrate, one then calls 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:


      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;
            }
        }
      }
 
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:

     ...
     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));
 
or

     ...
     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:


       GLQuadrature<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);
              ...
       }
 
As an optimization, the code caches the arguments raised to the fifth power to reduce the number of multiplications needed, assuming 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,


         glq = 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);
 
shows how to use 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
    Constructor
    Description
    VectorValuedGLQ(int n, int m)
    Constructor.
  • Method Summary

    Modifier and Type
    Method
    Description
    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 evaluated
    Get 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(Adder[] adders, double[] arguments)
    Integrate functions using precomputed arguments, accumulating the integral in an adder.
    void
    integrate(Adder[] adders, double a, double b)
    Integrate the function, accumulating the integral in an adder.
    void
    integrate(Adder[] adders, double a, double b, int ms)
    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
    mapping(double[] results, int m, double t, P p)
    A mapping that specifies the functions (with parameters) to integrate.
    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.

    Methods inherited from class java.lang.Object

    clone, equals, finalize, getClass, hashCode, notify, notifyAll, toString, wait, wait, wait
  • 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 use
      m - the number of values that will be provided as the result of an integration
      Throws:
      IllegalArgumentException - n is less than 1
  • Method Details

    • setParameters

      public void setParameters(P parameters)
      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

      public P 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 implementation mapping(double[],int,double,Object) using the parameters supplied by a previous call to setParameters(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 by getResultLength()
      m - the number of elements in the results array that are used
      t - the argument for the integration
    • mapping

      protected void mapping(double[] results, int m, double t, P p)
      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 call getParameters() to obtain any parameters that were configured via the setParameters(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 by getResultLength()
      m - the number of elements in the results array that are used
      t - the argument for the integration
      p - 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 by getResultLength()
      a - the lower limit of integration
      b - the upper limit of integration
      Throws:
      IllegalArgumentException - the first argument was null or an array that was too short
      UnsupportedOperationException - 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 integrals
      a - the lower limit of integration
      b - the upper limit of integration
      Throws:
      IllegalArgumentException - the first argument was null or an array that was too short
      UnsupportedOperationException - 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 by getResultLength()
      a - the lower limit of integration
      b - the upper limit of integration
      p - the parameters
      Throws:
      IllegalArgumentException - the first argument was null or an array that was too short
      UnsupportedOperationException - 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 adders
      a - the lower limit of integration
      b - the upper limit of integration
      p - the parameters
      Throws:
      IllegalArgumentException - the first argument was null or an array that was too short
      UnsupportedOperationException - 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 with integrate(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 integration
      b - 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 with integrate(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 integration
      b - the upper limit of integration
      n - 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:
           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();
       
      The weights returned by this method include a scaling factor for making the integral go from a to b instead of -1 to 1.
      Parameters:
      a - the lower limit of integration
      b - the upper limit of integration
      n - 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 with getArguments(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 to getArguments(double,double).
      Parameters:
      results - an array holding the results, with the minimum width set by the constructor and read by getResultLength()
      arguments - the array returned by getArguments(a,b)
      Throws:
      IllegalArgumentException - the first argument was null or an array that was too short
      UnsupportedOperationException - 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 with getArguments(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 results
      arguments - the array returned by getArguments(a,b)
      Throws:
      IllegalArgumentException - the first argument was null or an array that was too short
      UnsupportedOperationException - 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 with getArguments(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 by getResultLength()
      arguments - the array returned by getArguments(a,b)
      p - the parameters
      Throws:
      IllegalArgumentException - the first argument was null or an array that was too short
      UnsupportedOperationException - 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 with getArguments(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
      arguments - the array returned by getArguments(a,b)
      p - the parameters
      Throws:
      IllegalArgumentException - the first argument was null or an array that was too short
      UnsupportedOperationException - 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 by getResultLength()
      a - the lower limit of integration
      b - the upper limit of integration
      ms - the number of subintervals
      Throws:
      IllegalArgumentException - the first argument was null or an array that was too short
      UnsupportedOperationException - 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 adders
      a - the lower limit of integration
      b - the upper limit of integration
      ms - the number of subintervals
      Throws:
      IllegalArgumentException - the first argument was null or an array that was too short
      UnsupportedOperationException - 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 by getResultLength()
      a - the lower limit of integration
      b - the upper limit of integration
      ms - the number of subintervals
      p - the parameters
      Throws:
      IllegalArgumentException - the first argument was null or an array that was too short
      UnsupportedOperationException - 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 adders
      a - the lower limit of integration
      b - the upper limit of integration
      ms - the number of subintervals
      p - the parameters
      Throws:
      IllegalArgumentException - the first argument was null or an array that was too short
      UnsupportedOperationException - a mapping was not implemented
    • newInstance

      public 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.
      Parameters:
      map - the map
      n - the number of points to use
      m - the number of elements in the results array that are used
      Returns:
      the new instance of VectorValuedGLQ
      Throws:
      IllegalArgumentException - n is less than 1