3

The composite simpson's rule subdivides the interval into n equal subintervals, with n even.

Then $$\int_a^b f(x) dx \approx \frac{h}{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 $x_0 = a$ and $x_n=b$, $h=(b-a)/n$)

source

So what to do if given n equal intervals, with n odd?

One solution I implemented is to calculate using Simpson's rule for all but the last interval and use the trapezoid rule for the last interval. This seems to work well enough.

But I started thinking, how would I calculate the last interval using Simpson's rule as well?

Consider

$$\begin{aligned} I_1 &= \frac{h}{3}[f(x_{n-3}) + 4f(x_{n-2}) + f(x_{n-1})] \\ I_2 &= \frac{h}{3}[f(x_{n-2}) + 4f(x_{n-1}) + f(x_n)] \\ \end{aligned}$$

$x_n$ is the end of the integral interval. $I_1$ would be the sum corresponding to the last subinterval such that the total number of subintervals considered is even. $I_2$ is Simpson's rule applied to the final, odd subinterval, but half of the subinterval overlaps with $I_1$.

Is there a way to identify what's common to $I_1$ and $I_2$ and subtract it out? Just staring at it, I can't easily see what is the common part. Is it as simple as $(1/2)(h/3)[5f(x_{n-2}) + 5f(x_{n-1})]$?

bernie
  • 165
  • 1
  • 1
  • 6
  • $$\frac34I_1 + \frac34 I_2 = \int_{x_{n-3}}^{x_n} f(x),\mathrm{d}x - \frac14 f''h^3 + O(h^4),$$ which is obviously a worse remainder term. You can solve $\alpha I_1 + \beta I_2 = \int f + O(h^3)$ as a linear equation in two unknowns $\alpha, \beta$. – Kirill Nov 26 '16 at 16:06
  • I think if you take a linear combination of four two-interval integrals with coefficients $(\frac{15}{16},\frac5{16},\frac5{16},\frac{15}{16})$, then that gives an $\frac{-5}{144}f''''h^5$ remainder term. – Kirill Nov 26 '16 at 16:16
  • could you please explain where $(3/4)$ comes from? I tried to work out on paper, but could not figure it out. I probably need to look into the details of the derivation of simpson's rule eh? – bernie Nov 26 '16 at 22:02
  • I used this: Solve[ Thread[0 == Coefficient[ Series[x h/3 (f[0] + 4 f[h] + f[2 h]) + y h/3 (f[h] + 4 f[2 h] + f[3 h]), {h, 0, 6}] - Integrate[Series[f[h], {h, 0, 6}], {h, 0, 3 h}], h, {1, 2}]], {x, y}] – Kirill Nov 26 '16 at 22:47
  • Oh. I guess it's not so simple. Anyway, i am using simpson's rule + trapezoid in the case of odd intervals since it's not critical – bernie Nov 27 '16 at 17:27
  • That probably has the wrong order of error. [It's quite simple: expand everything in Taylor series and set equal to the true integral, expanded in Taylor series.] – Kirill Nov 27 '16 at 19:23

1 Answers1

5

A simple solution is to apply Simpson's (standard) rule to the first $n-3$ grid points, where $n-3$ is even for $n$ odd, and cover the remaining three gridpoints via the Simpson 3/8 formula:

$$I_{3/8} = \frac{3h}{8}[f(x_{n-3}) + 3f(x_{n-2}) + 3f(x_{n-1}) + f(x_n)].$$

Both have remainder terms of order $\mathcal O(h^5)$, so it keeps the order of the Simpson integration.

davidhigh
  • 3,127
  • 14
  • 15