quadpack fcns need abs err argument, was: quad_qag error

Logged In: YES
user_id=28849
Originator: NO

Maxima says integrate(legendre_p(2,x),x,0,1) is 0. quad_qag says [2.0e-17, 4.28e-15, 31, 2]. That is, the integral is 2e-17.

The abnormal return code is 2, which is excessive round-off. I think this is correct because the default relative error is 1e-8 but the absolute error is 0.

I think quad_qag (and friends) probably should allow another arg to let the user specify the desired absolute error. And quad_qag probably shouldn't set the absolute error to 0; perhaps something more like double-float-epsilon (about 1e-16) should be used instead.

Logged In: YES
user_id=501686
Originator: NO

Agreed w/ Ray T on both counts: current behavior is correct, and quad_qag and friends should accept an absolute error argument. Since the original Fortran QUADPACK functions accept an absolute error argument, all that Maxima quadpack need to do is to pass the argument. Changing title of this report accordingly.

Logged In: YES
user_id=28849
Originator: NO

Further investigation indicates that the problem is not caused by epsabs being 0. The first call to dqk31 (because key = 3) is telling quadpack that the integral suffers from roundoff.

Using the adaptive quad_qags routine doesn't help either because the first call to dqk21 indicates roundoff error already and no further progress is made.

Nevertheless, we should allow an epsabs argument.

Logged In: YES
user_id=28849
Originator: NO

Support for abs err argument added. See quadpack.lisp, rev 1.10 and related checkins.

Now we get

quad_qag(legendre_p(2,x),x,0,1,3,'epsabs = 1d-12) -> [2.0e-17, 4.28e-15, 31, 0]

Closing report.