>
> Sage punts numerical integrals to QUADPACK or a translation of it,
>
What does GSL use? I forgot that Scipy also has quadrature, in addition to
Maxima... wealth of riches.
> right? QUADPACK is based on Gauss-Kronrod rules which are essentially
> Gaussian integration + an efficient refinement scheme. To determine
> whether to refine the estimate, the integrand is evaluated at some
> points. The difficulty with 1/sqrt(x^4 + 2) is that it is nearly zero
> for much of the interval (1, 10^11). If it is close enough to zero
> at the evaluation points, the refinement heuristic can be fooled.
> Probably other quadrature methods can be fooled in a similar way.
>
> I guess this is a bug in the sense that the result is far from correct
> even though the code is a correct implementation of the method.
> I wonder what other heuristic could be invented -- maybe since we
> can see that integrals on shorter intervals are much larger and
> the integrand is nonnegative, the integral on the whole interval
> must be at least that large -- perhaps that can be formalized.
>
> I am trying the examples with QUADPACK functions in Maxima --
> I find that quad_qags(f2, x, 1, 10^11) fails (with error=5, "integral
> is probably divergent or slowly convergent") but
> quad_qag(f2, x, 1, 10^11, 4) succeeds, likewise quad_qagi(f2, x, 1, inf)
> succeeds. If Sage is indeed calling QUADPACK, perhaps at least the
> error number can be reported?
>
That would be the "nintegrate()" method on symbolic expressions, and indeed
sage: f2.nintegrate(x,1,1000000)
(0.8815680604421181, 1.6586910832421555e-12, 819, 0)
sage: f2.nintegrate(x,1,10000000)
(-1.0000000002652395e-07, 9.9265675869388e-15, 861, 5)
so the "5" in the latter is the error code, which we even document as
- ``5`` - integral is probably divergent or slowly
convergent
Also, gp seems to handle it.
sage: gp.eval('intnum(x=1,100000000000,1/sqrt(x^4+2))')
'0.88156906043147435374520375967552406680'
sage: gp.eval('intnum(x=1,1000000000000,1/sqrt(x^4+2))')
'0.88156906044791138558085421922579474969'
Harald, thanks. I think now we have six or seven ways to do quadrature in
Sage!
According
to
https://www.gnu.org/software/gsl/manual/html_node/Numerical-Integration.html#Numerical-Integration
GSL does indeed reimplement QUADPACK, so I'm not sure what we should do
here.
--
You received this message because you are subscribed to the Google Groups
"sage-support" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to sage-support+unsubscr...@googlegroups.com.
To post to this group, send email to sage-support@googlegroups.com.
Visit this group at http://groups.google.com/group/sage-support.
For more options, visit https://groups.google.com/d/optout.
