Simpson s Rule


Simpson's Rule

The trapezoidal algorithm is a second-order function integration technique. There are many other integration schemes with the same general form as Eq. (21.4). A commonly used example is the fourth-order technique known as Simpson's rule.

Equation 21.12

graphics/21equ12.gif


In Eq. (21.12), the 4/3 and 2/3 factors alternate for the interior points and D x = ( b a )/( N “ 1). The number of divisions, N , must be odd since Simpson's rule uses parabolas rather than lines to approximate the intervals.

Simpson's rule shares a common characteristic with the trapezoidal algorithm in that one approximation can be used as a building block toward a more accurate approximation. In the case of Simpson's rule, this feature is accomplished by adding three times the points that were added in the previous approximation. Another interesting feature of Simpson's rule is that the expression for any given level of approximation is related to the corresponding trapezoidal rule expressions. For example, consider an integration range between x = 0 and x = 1. If you multiply 4/3 times Eq. (21.6) and subtract 1/3 times Eq. (21.5) you obtain the following expression.

Equation 21.13

graphics/21equ13.gif


Eq. (21.13) is the five-point version of Simpson's rule. Indeed, the Simpson's rule value at any given level of refinement can always be expressed in terms of the corresponding trapezoidal expressions.

Equation 21.14

graphics/21equ14.gif


The S term in Eq. (21.14) is Simpson's rule's evaluation of the integral and the y terms are the trapezoidal approximations.

Using Eq. (21.14) it is very easy to implement Simpson's rule with only slight modifications to the trapezoidal() method. We will call our Simpson's rule method simpsonsRule() and place it inside the Integrator class. The method computes successive approximations to the integral function with the trapezoidal rule and uses the trapezoidal results to compute Simpson's rule approximations according to Eq. (21.14).

The method takes the same input arguments as the trapezoidal() method. After declaring some variables , the method computes initial values of the trapezoidal and Simpson's rule approximations. The simpsonsRule() method then goes through the same iteration sequence as the trapezoidal() method. Points are added to the integration domain, each iteration adding twice as many points as the one before it. The Simpson's rule approximation is computed from the current and previous trapezoidal values. When the Simpson's rule values converge to within a certain tolerance, the method returns.

The simpsonsRule() method source code is ”

 public static double simpsonsRule(Function function,          double start, double end, double tolerance) {   double dx = end - start;   double value, oldValue, sum, scale;   double valueSimpson, oldValueSimpson;   int maxIter = 20;   //  Set initial values for the trapezoidal   //  method (value) and Simpson's rule (valueSimpson)   value = 0.5*dx*( function.getValue(start) +                    function.getValue(end) );   valueSimpson = 1.0e+6;   //  Subdivide the integration range by adding more   //  points, first one point, then two more, then   //  four more, and so on. Recompute the extended   //  trapezoidal value and the Simpson's rule value.   //  If the Simpson's rule value converges, return   //  the value.   for(int i=0; i<maxIter; ++i) {     sum = 0.0;     scale = Math.pow(0.5,i+1);     for(int j=0; j<Math.pow(2,i); ++j) {       sum += function.getValue(                   start + scale*(2*j + 1.0)*dx);     }     oldValue = value;     value = 0.5*value + scale*dx*sum;     oldValueSimpson = valueSimpson;     valueSimpson = (4.0*value - oldValue)/3.0;     if ( Math.abs(valueSimpson-oldValueSimpson) <=          tolerance*Math.abs(oldValueSimpson) ) {       return valueSimpson;     }   }   //  Solution did not converge   return 12345.6; } 

Let's apply the simpsonRule() method to the same two problems as were solved by the trapezoidal() method. The driver to do this computation is called TestSimpson.java and its source code is shown here.

 import TechJava.MathLib.*; public class TestSimpson {   public static void main(String args[]) {     double value;     //  Integrate the function y = x*ln(x) from     //  x=2 to x=4 with an error tolerance of 1.0e-4     LogFunction function = new LogFunction();     value = Integrator.simpsonsRule(                            function,2.0,4.0,1.0e-4);     System.out.println("value = " + value);     //  Integrate the function y = sqrt(x^2 + 2^2) from     //  x=0 to x=3 with an error tolerance of 1.0e-4     SqrtFunction function2 = new SqrtFunction(2.0);     value = Integrator.simpsonsRule(                           function2,0.0,3.0,1.0e-4);     System.out.println("value = " + value);   } } 

Output ”

 value = 6.704064541919408 value = 7.7978469366358025 

The simpsonsRule() method provides better accuracy than the trapezoidal() method with fewer iterations and function evaluations required to achieve a converged solution. The simpsonsRule() method results match the exact results to eight and five decimal places respectively. The simpsonsRule() method converged in three iterations requiring nine function evaluations. The trapezoidal() method, on the other hand, required five iterations and 65 function evaluations to achieve convergence.



Technical Java. Applications for Science and Engineering
Technical Java: Applications for Science and Engineering
ISBN: 0131018159
EAN: 2147483647
Year: 2003
Pages: 281
Authors: Grant Palmer

flylib.com © 2008-2017.
If you may any questions please contact us: flylib@qtcs.net