001 /*
002 * Licensed to the Apache Software Foundation (ASF) under one or more
003 * contributor license agreements. See the NOTICE file distributed with
004 * this work for additional information regarding copyright ownership.
005 * The ASF licenses this file to You under the Apache License, Version 2.0
006 * (the "License"); you may not use this file except in compliance with
007 * the License. You may obtain a copy of the License at
008 *
009 * http://www.apache.org/licenses/LICENSE-2.0
010 *
011 * Unless required by applicable law or agreed to in writing, software
012 * distributed under the License is distributed on an "AS IS" BASIS,
013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
014 * See the License for the specific language governing permissions and
015 * limitations under the License.
016 */
017 package org.apache.commons.math3.analysis.integration.gauss;
018
019 import org.apache.commons.math3.analysis.UnivariateFunction;
020 import org.apache.commons.math3.exception.DimensionMismatchException;
021 import org.apache.commons.math3.exception.NonMonotonicSequenceException;
022 import org.apache.commons.math3.util.MathArrays;
023 import org.apache.commons.math3.util.Pair;
024
025 /**
026 * Class that implements the Gaussian rule for
027 * {@link #integrate(UnivariateFunction) integrating} a weighted
028 * function.
029 *
030 * @since 3.1
031 * @version $Id: GaussIntegrator.java 1382197 2012-09-07 22:35:01Z erans $
032 */
033 public class GaussIntegrator {
034 /** Nodes. */
035 private final double[] points;
036 /** Nodes weights. */
037 private final double[] weights;
038
039 /**
040 * Creates an integrator from the given {@code points} and {@code weights}.
041 * The integration interval is defined by the first and last value of
042 * {@code points} which must be sorted in increasing order.
043 *
044 * @param points Integration points.
045 * @param weights Weights of the corresponding integration nodes.
046 * @throws NonMonotonicSequenceException if the {@code points} are not
047 * sorted in increasing order.
048 */
049 public GaussIntegrator(double[] points,
050 double[] weights)
051 throws NonMonotonicSequenceException {
052 if (points.length != weights.length) {
053 throw new DimensionMismatchException(points.length,
054 weights.length);
055 }
056
057 MathArrays.checkOrder(points, MathArrays.OrderDirection.INCREASING, true, true);
058
059 this.points = points.clone();
060 this.weights = weights.clone();
061 }
062
063 /**
064 * Creates an integrator from the given pair of points (first element of
065 * the pair) and weights (second element of the pair.
066 *
067 * @param pointsAndWeights Integration points and corresponding weights.
068 * @throws NonMonotonicSequenceException if the {@code points} are not
069 * sorted in increasing order.
070 *
071 * @see #GaussIntegrator(double[], double[])
072 */
073 public GaussIntegrator(Pair<double[], double[]> pointsAndWeights)
074 throws NonMonotonicSequenceException {
075 this(pointsAndWeights.getFirst(), pointsAndWeights.getSecond());
076 }
077
078 /**
079 * Returns an estimate of the integral of {@code f(x) * w(x)},
080 * where {@code w} is a weight function that depends on the actual
081 * flavor of the Gauss integration scheme.
082 * The algorithm uses the points and associated weights, as passed
083 * to the {@link #GaussIntegrator(double[],double[]) constructor}.
084 *
085 * @param f Function to integrate.
086 * @return the integral of the weighted function.
087 */
088 public double integrate(UnivariateFunction f) {
089 double s = 0;
090 double c = 0;
091 for (int i = 0; i < points.length; i++) {
092 final double x = points[i];
093 final double w = weights[i];
094 final double y = w * f.value(x) - c;
095 final double t = s + y;
096 c = (t - s) - y;
097 s = t;
098 }
099 return s;
100 }
101
102 /**
103 * @return the order of the integration rule (the number of integration
104 * points).
105 */
106 public int getNumberOfPoints() {
107 return points.length;
108 }
109 }