Few years ago I have published some ideas on how to improve numerical stability of the Newton-Cotes formulas of closed type. Today I want to apply the same ideas to so-called “open” NC-formulas when boundary points are not used for integral approximation.

Just to remind the basic facts and ideas. Newton-Cotes rules use polynomial interpolation to approximate the function under integral sign. Then integral of the function is assumed to be equal to the integral of the approximating polynomial (which can be computed explicitly).

Idea is simple enough, but it has one serious drawback.

Although polynomials are easy to work with, they are very unpredictable/problematic beasts. Enough to say that even such simple task as root finding of polynomials is inherently ill-conditioned (e.g. Wilkinson’s polynomial).

Similar issues with stability arise when we use polynomial interpolation over equidistant nodes (used by Newton-Cotes). The classic example of why this is a bad idea is Runge’s phenomenon.

Because of the issues, Gauss-family quadrature are usually used instead of Newton-Cotes. However in some applications we stuck with uniform sampling and have no other choice but to improve the Newton-Cotes.

One of the most obvious ways of doing that is to use **least-squares approximation** instead of interpolation. Another useful tool is studying properties of the formulas in spectral domain.

There is additional information in my older posts (1, 2, etc.) – here I just provide resulting formulas for “stable” Newton-Cotes rules of open type.

## Task

We consider the problem of numerical approximation of integral of a function $f(x)$ on $x\in[a,b]$:

having only discrete values ${f_i=f(x_i)}$ sampled at equidistant points $x_i$ with step $h$:

Our goal here is to use only interior nodes for approximation (ignoring $f_0$ and $f_{N-1}$).

Closed rules (when all nodes are used) were derived here.

## Least-squares Newton-Cotes formulas (open type)

$O(h^5)$:

$O(h^7)$:

Spectral plots show how sensitive every formula to noise (e.g. rounding errors, etc.) in the data. Overall as closer curve to axis X as better the stability of the corresponding rule. In combination with $l_1$ norm it gives us full picture for stability assessment.

Classic rules based on interpolation show strong amplification of rounding errors (or any noise in the data in general) – the red curves.

Not surprisingly we see the substantial improvement when we use least-squares. However to reach the stability we have to use much longer filters. It looks like the “open” case is ill-conditioned in its nature (interpolatory formulas have very high $l_1$ norm right from the onset).

We advise using rule with at least $N=11$ if you consider $O(h^7)$-family.