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