Subject: how to "fool" Gauss-Kronrod 21-point rule
From: Raymond Toy
Date: Fri, 13 Apr 2012 10:30:12 -0700
On Thu, Apr 12, 2012 at 11:22 AM, Edwin Woollett <woollett at charter.net>wrote:
>
> Robert,
>
> In your paper
> http://riso.sourceforge.net/**docs/heterogeneous-polytree.**pdf<http://riso.sourceforge.net/docs/heterogeneous-polytree.pdf>
>
> which I found via http://en.scientificcommons.**org/robert_dodier<http://en.scientificcommons.org/robert_dodier>
>
> on page 5, in section 5: Numerical Subtleties,
> Integration Algorithm,
>
> you say:
>
> "Integration algorithm. Numerical integrations in
> more than one dimensions are difficult. In the current
> implementation, multidimensional integrations are reduced
> to repeated one-dimensional integrations. The
> one-dimensional integrations are carried out by an
> adaptive region-splitting algorithm based on a Gauss-
> Kronrod 21-point rule. (The code is a translation
> of the QAGS algorithm from quadpack, a collection
> of quadrature algorithms available from www.-
> netlib.org.) The adaptive quadrature algorithm can
> be "fooled" if the integrand varies on a scale much
> smaller than I/42, where I is the length of the interval
> of integration. For this reason, riso tries to find
> the smallest effective support of the integrand, as described
> under the preceding heading."
>
> I would be interested in some examples which illustrate
> this ability to fool quad_qags.
>
Not exactly sure what you're looking for but here is one example where
quad_qags gives no error or warning messages, but the result is totally
wrong:
quad_qags(4^(-16)/((x - %pi/4)^2 + 16^(-16)),x,0,1);
[1.4602196559157534e-6, 1.8681728113217153e-6, 693, 5]
Maxima can actually compute the integral:
ratsimp(integrate(4^(-16)/((x - %pi/4)^2 + 16^(-16)),x,0,1));
atan(1073741824 %pi) - atan(1073741824 %pi - 4294967296)
float(%)
3.141592652208402
(This example from the quadpack book, also listed in
numerical/slatec/quadpack.lisp.
Ray