2
$\begingroup$

For a given function $f(x)$, I have tried to find its numerical integral using Simpson's 1/3 and Simpson's 3/8 rules.

I then compare the solution from the numerical quadratures to the analytical integral to calculate the errors. I find that the error for the 3/8 rule is double that of the 1/3 rule, while the opposite is expected. Where did I go wrong?

The Python3 function for the 1/3 rule:

def simpsons13(a, b, N):
    """
    Calculates the numerical integral of a function f(x) using the Simpson's 1/3rd rule:

    F(x) = Σ(0 to (N-2)/2) Δx/3 * (f(x(2i)) + 4f(x(2i + 1)) + f(x(2i + 2)))

    Parameters:

    a:      The lower limit of the definite integral (real)
    b:      The upper limit of the definite integral (real)
    N:      A positive, even integer to denote the number of intervals of the function for the integral

    Returns I, the numerical integral calculated through Simpson's 1/3rd rule
    """

    if N%2 != 0:
        return "N must be even"     

    dx = (b - a)/N
    narr = np.linspace(a, b, N+1)   # N intervals corresponds to N + 1 points
    I = 0

    for i in range(int((N - 2) / 2) + 1):

        I = I + dx/3*(f(narr[2*i]) + 4*f(narr[2*i + 1]) + f(narr[2*i + 2]))

    return I

The Python3 function for the 3/8 rule:

def simpsons38(a, b, N):
    """
    Calculates the numerical integral of a function f(x) using the Simpson's 3/8 rule:

    F(x) = Σ(0 to (N-3)/3) 3Δx/8 * (f(x(3i)) + 3f(x(3i + 1)) + 3f(x(3i + 2)) + f(x(3i + 3)))

    Parameters:

    a:      The lower limit of the definite integral (real)
    b:      The upper limit of the definite integral (real)
    N:      A positive, even integer to denote the number of intervals of the function for the integral

    Returns I, the numerical integral calculated through Simpson's 3/8 rule
    """

    if N%3 != 0:
        return "N must be divisible by 3"

    dx = (b - a)/N
    narr = np.linspace(a, b, N+1)
    I = 0

    for i in range(int((N-3)/3) + 1):

        I = I + 3*dx/8 * (f(narr[3*i]) + 3*f(narr[3*i + 1]) + 3*f(narr[3*i + 2]) + f(narr[3*i + 3]))

    return I 

I consider the simple function $f(x) = xe^x$ and find its numerical integral from $[-\pi, \pi]$ by considering $N = 24$ intervals. I calculate the errors between the numerical and analytical solutions as:

error = np.abs(ans - approx)

I get the error by the 1/3 method to be $0.003669436$ and by the 3/8 method to be $0.00816864$. The 3/8 rule has more than double the error of the 1/3 rule. Why is this happening?

$\endgroup$
1
  • 1
    $\begingroup$ Could you check if the number of sub-intervals resp. step size in the error formulas refers to the actual sample points or to segments with internal points. That is, in the latter view you would have applied 1/3 to 12 segments and 3/8 to 8 segments. /// Comparing for equal numbers of function evaluations is of course a fair measure. $\endgroup$ Aug 16 at 8:33
5
$\begingroup$

The simple segment-$[a,b]$-with-midpoints error formulas are for the 1/3 rule $$ E=-\frac{(b-a)^5}{90·2^5}f^{(4)}(\zeta) $$ and for the 3/8 rule $$ E=-\frac{(b-a)^5}{90·72}f^{(4)}(\zeta) $$ Now the first rule has 2 sub-intervals while the second has 3. To get comparable values for equal numbers of function evaluations apply this to an interval with 6 sub-intervals, giving 3 segments for the 1/3 rule and 2 segments for the 3/8 rule. Then the composite 1/3 rule has an error of $$ E=-3·\frac{((b-a)/3)^5}{90·2^5}f^{(4)}(\zeta) =-\frac{(b-a)^5}{90·2^5·3^4}f^{(4)}(\zeta) =-\frac{((b-a)/6)^5}{30}f^{(4)}(\zeta) $$ and the composite 3/8 rule $$ E=-2·\frac{((b-a)/2)^5}{90·72}f^{(4)}(\zeta) =-\frac{(b-a)^5}{40·2^5·3^4} =-\frac{((b-a)/6)^5}{40/3}f^{(4)}(\zeta) $$ The midpoints $\zeta$ can and will be different for each formula.

In total the error of the 3/8 rule is about $\frac94$ of the error of the 1/3 rule for the same number of sampling points. This is also exactly what you observed, as $36·\frac94=81$.


Note that if this gets expanded to Runge-Kutta methods, where the "classical" method is based on the 1/3 rule and the 3/8 method one the 3/8 rule, both methods have 4 stages, that is, 4 function evaluations. This means that both methods can be compared based on the one-segment error formulas, where thus the 3/8 method has the smaller step error. (One would have to confirm that with the actual one-step Runge-Kutta error, which was indeed done already by Kutta in the original paper.)

$\endgroup$
1
  • $\begingroup$ Great point! I didn't think about the unequal sub-intervals. This solves the issue. Thanks! $\endgroup$
    – justauser
    Aug 16 at 12:57

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.