[SciPy-Dev] Incorrect results from scipy.integrate.romberg

[Thomas Robitaille]

Hi,

I think there is some kind of bug with the romberg integration function. In the following code, romberg integration give an incorrect answer that is too large by a factor of two:

from scipy.integrate import quad, romberg

def powerlaw(x):
    return x**3
    
print "%12.3e" % quad(powerlaw, 1e1, 1e5)[0] # gives 2.500.e+19
print "%12.3e" % romberg(powerlaw, 1e1, 1e5) # gives 5.000.e+19

The correct answer is 2.5e19, but romberg gives twice that.

Also, as in my previous email, it seems the tol= argument in romberg() is the absolute tolerance, not relative. Would it not be a good idea to be able to specify a relative tolerance?

Cheers,

Thomas