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 */ 017package org.apache.commons.math3.analysis.integration.gauss; 018 019import org.apache.commons.math3.analysis.UnivariateFunction; 020import org.apache.commons.math3.exception.DimensionMismatchException; 021import org.apache.commons.math3.exception.NonMonotonicSequenceException; 022import org.apache.commons.math3.util.Pair; 023 024/** 025 * This class's implements {@link #integrate(UnivariateFunction) integrate} 026 * method assuming that the integral is symmetric about 0. 027 * This allows to reduce numerical errors. 028 * 029 * @since 3.3 030 */ 031public class SymmetricGaussIntegrator extends GaussIntegrator { 032 /** 033 * Creates an integrator from the given {@code points} and {@code weights}. 034 * The integration interval is defined by the first and last value of 035 * {@code points} which must be sorted in increasing order. 036 * 037 * @param points Integration points. 038 * @param weights Weights of the corresponding integration nodes. 039 * @throws NonMonotonicSequenceException if the {@code points} are not 040 * sorted in increasing order. 041 * @throws DimensionMismatchException if points and weights don't have the same length 042 */ 043 public SymmetricGaussIntegrator(double[] points, 044 double[] weights) 045 throws NonMonotonicSequenceException, DimensionMismatchException { 046 super(points, weights); 047 } 048 049 /** 050 * Creates an integrator from the given pair of points (first element of 051 * the pair) and weights (second element of the pair. 052 * 053 * @param pointsAndWeights Integration points and corresponding weights. 054 * @throws NonMonotonicSequenceException if the {@code points} are not 055 * sorted in increasing order. 056 * 057 * @see #SymmetricGaussIntegrator(double[], double[]) 058 */ 059 public SymmetricGaussIntegrator(Pair<double[], double[]> pointsAndWeights) 060 throws NonMonotonicSequenceException { 061 this(pointsAndWeights.getFirst(), pointsAndWeights.getSecond()); 062 } 063 064 /** 065 * {@inheritDoc} 066 */ 067 @Override 068 public double integrate(UnivariateFunction f) { 069 final int ruleLength = getNumberOfPoints(); 070 071 if (ruleLength == 1) { 072 return getWeight(0) * f.value(0d); 073 } 074 075 final int iMax = ruleLength / 2; 076 double s = 0; 077 double c = 0; 078 for (int i = 0; i < iMax; i++) { 079 final double p = getPoint(i); 080 final double w = getWeight(i); 081 082 final double f1 = f.value(p); 083 final double f2 = f.value(-p); 084 085 final double y = w * (f1 + f2) - c; 086 final double t = s + y; 087 088 c = (t - s) - y; 089 s = t; 090 } 091 092 if (ruleLength % 2 != 0) { 093 final double w = getWeight(iMax); 094 095 final double y = w * f.value(0d) - c; 096 final double t = s + y; 097 098 s = t; 099 } 100 101 return s; 102 } 103}