sagemath / sage

Main repository of SageMath
https://www.sagemath.org
Other
1.28k stars 439 forks source link

Unevaluated integrals to infinity have nonsense numeric value #22567

Open pelegm opened 7 years ago

pelegm commented 7 years ago

Running

integrate(floor(x), x, 0, infinity, algorithm='sympy')

returns

integrate(floor(x), x, 0, +Infinity)

and trying

integrate(floor(x), x, 0, infinity, algorithm='sympy').n()

returns (!!!) -679.7441466712775. This is awfully wrong. It shows with any unevaluated integral.

CC: @kcrisman @rwst

Component: calculus

Keywords: integrate

Issue created by migration from https://trac.sagemath.org/ticket/22567

pelegm commented 7 years ago

Description changed:

--- 
+++ 
@@ -2,6 +2,11 @@

 ```python
 integrate(floor(x), x, 0, infinity, algorithm='sympy')
+```
+returns
+
+```python
+integrate(floor(x), x, 0, +Infinity)

and trying

rwst commented 7 years ago

Description changed:

--- 
+++ 
@@ -13,4 +13,4 @@
 ```python
 integrate(floor(x), x, 0, infinity, algorithm='sympy').n()

-returns (!!!) -679.7441466712775. This is awfully wrong. +returns (!!!) -679.7441466712775. This is awfully wrong. It shows with any unevaluated integral.

pelegm commented 7 years ago
comment:4

Note that it doesn't happen without algorithm='sympy'.

rwst commented 7 years ago
comment:5

I don't think this is correct:

sage: integrate((1+x+x^2)^(1/3), x, 0, infinity)
integrate((x^2 + x + 1)^(1/3), x, 0, +Infinity)
sage: _.n()
8.835093500042741e+20
kcrisman commented 7 years ago

Changed keywords from floor, integrate, sympy to integrate

kcrisman commented 7 years ago
comment:6

I want to point out that all these integrals involve infinity. I can't remember any more exactly where to find it (did poke about a bit) but there should be code, and not original code, for how to deal with this situation - I think it just sends it to the GSL algorithm numerical_integral which would likely give nonsense for infinity. And indeed that is the case - compare this Sage cell:

sage: print numerical_integral(floor(x), 0, +Infinity)
sage: print numerical_integral((x^2 + x + 1)^(1/3), 0, +Infinity)
(-679.7441466712775, 632.307547415802)
(8.835093500042741e+20, 1.6991730463958232e+21)