Adaptive Simpson and Romberg integration methods

[robve]

(03-24-2021 11:11 PM)Albert Chan Wrote:  

(03-24-2021 08:54 PM)robve Wrote:  Improper integrals require integration over an open interval to avoid evaluating the endpoints, for example using Romberg with midpoint quadratures, see e.g. Numerical Recipes 2nd Ed. Ch.4.4 qromo and as implemented in the HP-71B Math Pac.

No. closed-intervals is fine, if we can assume end-points really is zero.
Normally, u-transformation will achieve it, but not always ...

A counter-example, from one of my old thread. Expected: \(\displaystyle \int _{-1}^1 {dx \over \sqrt{1-x^2}} = asin(x)\bigg| _{-1}^1 = \pi\)

(02-06-2020 11:16 PM)Albert Chan Wrote:  A simple test showed that HP71B INTEGRAL do extrapolations from trapezoids, not mid-points rule.

>10 DISP INTEGRAL(-1, 1, 1E-5, 1/SQRT(1-IVAR^2)), IBOUND
>RUN
3.14156045534              -3.1415397954E-5

this failure to converge (65535 sample points !) is due to missing end points evaluation.

\(\large \int _{-1} ^1 f(x) dx
= \int _{-1} ^ 1 {3 \over 2}(1-u^2) f \left({ u (3-u^2) \over 2} \right) du
= \int _{-1} ^ 1 g(u) du \)

\(f(x) = {1\over\sqrt{1-x^2}}\quad → \quad
g(u) = {3 (1-u^2) \over \sqrt{(4-u^2)(1-u^2)^2}}\)

\(\displaystyle{\lim_{u^2 \to 1^-} g(u)} = \displaystyle{\lim_{u^2 \to 1^-}{3\over\sqrt{4-u^2}}} = \sqrt3 ≠ 0\)

Makes sense. Let me add to this that if the value at the endpoint is zero then it can certainly be "ignored". A closed interval integration method can be used by ignoring such floating point exceptions. So yes, that works in those cases.

However, this is a specific case for zero. We're considering any finite limiting value at one or both of the endpoints. A integral is improper and may require open intervals if its integrand goes to a finite limiting value at its finite upper and lower limits (integral endpoint), but cannot be evaluated right on one of those limits, for example \( \frac{\sin x}{x} \) at \( x=0 \) since \( \lim_{x\rightarrow 0} \frac{\sin x}{x} = 1 \).

Weird that HP-71B Math Pac uses trapezoidal quadratures when trying to exclude endpoints. The recommended Romberg modification is to use an open quadrature as the basis to replace (the popular) trapezoidal, while thereby preserving the nice properties of Romberg. Numerical Recipes 2nd ed thus gives a generalization of the Romberg algorithm with open quadratures, but there is no specific implementation given, which with some effort can be derived, e.g. as qromo in the post (my qromo is not the same as the generic qromo in NR 2nd ed, and differ significantly structurally).

The midpoint version qromo essentially performs high-degree polynomial extrapolation to obtain the finite values at the (missing) endpoints of the integration interval. The integral's precision is about as high as the trapezoidal version, despite the fact that the integrand is never evaluated at the endpoints, when considering the same number of points. Note the ratio of n subdivisions is about 2:3 for midpoint versus trapezoidal to get very roughly the same number of points in each (see post).

- Rob

Edit: fixed typo and last sentence added.