Unlike the Riemann sum methods, Stewart does not give a sigma formula for Simpson's Rule. Instead, the textbook defines it as a series:
$$ \int_{a}^b f(x) dx \approx S_n = \frac{\Delta x}{3} [f(x_0) + 4f(x_1) + 2f(x_2) + 4f(x_3) + ... + 2f(x_{n-2}) + 4f(x_{n-1}) + f(x_n)] $$Where $n$ is even and $\Delta x = \frac{b - a}{n} $.
The authors say "note the pattern of the coefficients" (pg. 520), which feels like a cop-out. What does $x_0$ mean? How do I find the actual values I'm supposed to plug into the function being integrated? Attempting to convert the formula into sigma form made these blindspots obvious.
$$ \frac{\Delta x}{3} \sum_{i = 0}^{n} \Delta x \ S(i) $$Where $n$ is even and $\Delta x = \frac{b - a}{n} $.
While it is possible to start the summation at 1 (to match the formulas for the left and right Riemann sums), I think it is simpler to start at 0.
What is $S(i)$? This is a special piecewise function that determines the coefficient (if any) to multiply $f(x_{i})$ by.
Note that precedence is important here. I don't know how to represent this in math notation, but programatically the logic is:
# Will not actually run because n, x_i, and f(x) are not defined
def S(i):
if (i == 0) or (i == n):
return f(x_i)
elif (i % 2 == 0):
return 2 * f(x_i)
else:
# odd
return 4 * f(x_i)
Finally, $x_{i}$ is a function that returns the value to plug into $f(x)$ given the current index $i$.
$$x_{i} = a + i \ \Delta x $$The notation, $x_{i}$, is somewhat confusing because it looks like a variable but is actually a function.
This is from an old uOttawa Calculus I exam.
Approximate the integral $ \int_{0}^4 \frac{1}{x+1} dx $ numerically to 4 decimal places with n = 4 subdivisions with Simpson's Rule.
# Lower bound
a = 0
# Upper bound
b = 4
# Number of intervals
n = 4
# Difference in x between intervals
deltaX = (b - a) / n
# Function being integrated
def f(x):
return 1.0 / (x + 1.0)
# Determines what value to plug into f(x) for each i
def currentX(i):
return a + (i * deltaX)
# Function that applies Simpson's coefficients (based on i)
def S(i):
# if zeroth and last values
if (i == 0) or (i == n):
# no coefficient
return f(currentX(i))
# if even
elif (i % 2 == 0):
# coefficient of 2
return 2 * f(currentX(i))
# if odd
else:
# coefficient of 4
return 4 * f(currentX(i))
Note that the return statements are different from the previously defined $S(i)$ function, here $currentX()$ is used to determine the value to submit to f(x).
# Computation
simpsonSum = 0.0
# Loop from 0 to n (inclusive)
for i in range(0, n+1):
simpsonSum += S(i)
# Don't forget to complete the formula by multiplying the sum by (deltaX / 3.0)!
(deltaX / 3.0)*simpsonSum
1.622222222222222
Which is roughly equivalent to Wolfram Alpha's answer.