View Javadoc

1   /*
2    * Licensed to the Apache Software Foundation (ASF) under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * The ASF licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    *
9    *      http://www.apache.org/licenses/LICENSE-2.0
10   *
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
16   */
17  package org.apache.commons.math.analysis.integration;
18  
19  import org.apache.commons.math.FunctionEvaluationException;
20  import org.apache.commons.math.MathRuntimeException;
21  import org.apache.commons.math.MaxIterationsExceededException;
22  import org.apache.commons.math.analysis.UnivariateRealFunction;
23  
24  /**
25   * Implements the <a href="http://mathworld.wolfram.com/TrapezoidalRule.html">
26   * Trapezoidal Rule</a> for integration of real univariate functions. For
27   * reference, see <b>Introduction to Numerical Analysis</b>, ISBN 038795452X,
28   * chapter 3.
29   * <p>
30   * The function should be integrable.</p>
31   *  
32   * @version $Revision: 799857 $ $Date: 2009-08-01 09:07:12 -0400 (Sat, 01 Aug 2009) $
33   * @since 1.2
34   */
35  public class TrapezoidIntegrator extends UnivariateRealIntegratorImpl {
36  
37      /** Intermediate result. */
38      private double s;
39  
40      /**
41       * Construct an integrator for the given function.
42       * 
43       * @param f function to integrate
44       * @deprecated as of 2.0 the integrand function is passed as an argument
45       * to the {@link #integrate(UnivariateRealFunction, double, double)}method.
46       */
47      @Deprecated
48      public TrapezoidIntegrator(UnivariateRealFunction f) {
49          super(f, 64);
50      }
51  
52      /**
53       * Construct an integrator.
54       */
55      public TrapezoidIntegrator() {
56          super(64);
57      }
58  
59      /**
60       * Compute the n-th stage integral of trapezoid rule. This function
61       * should only be called by API <code>integrate()</code> in the package.
62       * To save time it does not verify arguments - caller does.
63       * <p>
64       * The interval is divided equally into 2^n sections rather than an
65       * arbitrary m sections because this configuration can best utilize the
66       * alrealy computed values.</p>
67       *
68       * @param f the integrand function
69       * @param min the lower bound for the interval
70       * @param max the upper bound for the interval
71       * @param n the stage of 1/2 refinement, n = 0 is no refinement
72       * @return the value of n-th stage integral
73       * @throws FunctionEvaluationException if an error occurs evaluating the
74       * function
75       */
76      double stage(final UnivariateRealFunction f,
77                   final double min, final double max, final int n)
78          throws FunctionEvaluationException {
79          
80          long i, np;
81          double x, spacing, sum = 0;
82          
83          if (n == 0) {
84              s = 0.5 * (max - min) * (f.value(min) + f.value(max));
85              return s;
86          } else {
87              np = 1L << (n-1);           // number of new points in this stage
88              spacing = (max - min) / np; // spacing between adjacent new points
89              x = min + 0.5 * spacing;    // the first new point
90              for (i = 0; i < np; i++) {
91                  sum += f.value(x);
92                  x += spacing;
93              }
94              // add the new sum to previously calculated result
95              s = 0.5 * (s + sum * spacing);
96              return s;
97          }
98      }
99  
100     /** {@inheritDoc} */
101     @Deprecated
102     public double integrate(final double min, final double max)
103         throws MaxIterationsExceededException, FunctionEvaluationException, IllegalArgumentException {
104         return integrate(f, min, max);
105     }
106 
107     /** {@inheritDoc} */
108     public double integrate(final UnivariateRealFunction f,
109                             final double min, final double max)
110         throws MaxIterationsExceededException, FunctionEvaluationException, IllegalArgumentException {
111         
112         int i = 1;
113         double t, oldt;
114         
115         clearResult();
116         verifyInterval(min, max);
117         verifyIterationCount();
118 
119         oldt = stage(f, min, max, 0);
120         while (i <= maximalIterationCount) {
121             t = stage(f, min, max, i);
122             if (i >= minimalIterationCount) {
123                 final double delta = Math.abs(t - oldt);
124                 final double rLimit =
125                     relativeAccuracy * (Math.abs(oldt) + Math.abs(t)) * 0.5; 
126                 if ((delta <= rLimit) || (delta <= absoluteAccuracy)) {
127                     setResult(t, i);
128                     return result;
129                 }
130             }
131             oldt = t;
132             i++;
133         }
134         throw new MaxIterationsExceededException(maximalIterationCount);
135     }
136 
137     /** {@inheritDoc} */
138     @Override
139     protected void verifyIterationCount() throws IllegalArgumentException {
140         super.verifyIterationCount();
141         // at most 64 bisection refinements
142         if (maximalIterationCount > 64) {
143             throw MathRuntimeException.createIllegalArgumentException(
144                     "invalid iteration limits: min={0}, max={1}",
145                     0, 64);
146         }
147     }
148 }