Counterexample of Fubini - Strange Result with Maxima

[Michael Jung]

Dear fellow developers,

I've encountered a really strange result in Sage while using Maxima.

sage: f(x,y) = (x^2-y^2)/(x^2+y^2)^2
sage: integrate(integrate(abs(f(x,y)), x, 0, 1), y, 0, 1)        

-1/4*pi   

This is really weird. At least, the result should be positive! SymPy however yields the correct result:

sage: integrate(integrate(abs(f(x,y)), x, 0, 1, algorithm='sympy'), y, 0, 1, algorithm='sympy')

+Infinity

What is the issue here? Did I do something wrong, or is this a bug with Maxima?

Best
Michael

[Michael Orlitzky]

On 3/18/20 11:48 AM, Michael Jung wrote:
> Dear fellow developers,
> I've encountered a really strange result in Sage while using Maxima.
>
> |
> sage:f(x,y)=(x^2-y^2)/(x^2+y^2)^2
> sage:integrate(integrate(abs(f(x,y)),x,0,1),y,0,1)
>
> -1/4*pi
>
> |||
>
> This is really weird. At least, the result should be positive! SymPy
> however yields the correct result:

> ...


For posterity:

Maxima is unable to integrate pretty much anything with an absolute
value in it. I tried to fix this a decade ago in

https://trac.sagemath.org/ticket/11483

but the cure was worse than the disease, and it had to be rolled back:

https://trac.sagemath.org/ticket/12731

There is still a lot of room for improvement. SymPy could be tried first
when integrating expressions containing an absolute value, for one. We
already _fall back_ to giac/sympy if maxima throws an error; but when it
simply returns garbage, the problem goes unnoticed.

[Sébastien Labbé]

There is still a lot of room for improvement. SymPy could be tried first
when integrating expressions containing an absolute value, for one. We
already _fall back_ to giac/sympy if maxima throws an error; but when it
simply returns garbage, the problem goes unnoticed.

Why do we use maxima first as opposed to giac/sympy? Is it because it is faster than giac/sympy? Is it because it returns answers that are correct but for which giac/sympy returns incorrect results?

It seems to me that a library returning wrong answer as opposed to no answer should not be the default.

Sébastien

 

[Michael Orlitzky]

On 12/18/20 2:55 AM, Sébastien Labbé wrote:
>
> Why do we use maxima first as opposed to giac/sympy? Is it because it is
> faster than giac/sympy? Is it because it returns answers that are
> correct but for which giac/sympy returns incorrect results?
>

I personally have never tried it because I'm afraid to check the
mathematical equivalence of the 289075206 doctests that would change.

[Dima Pasechnik]

isn't pynac/ginac used a lot?

--
You received this message because you are subscribed to the Google Groups "sage-devel" group.
To unsubscribe from this group and stop receiving emails from it, send an email to sage-devel+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/sage-devel/4a61a913-8e06-1247-63d4-4113ebaa870e%40orlitzky.com.

[Michael Orlitzky]

On 12/18/20 9:53 AM, Dima Pasechnik wrote:
> isn't pynac/ginac used a lot?

At least for symbolic integration, the current list in
symbolic/integration/integral.py is,

self.integrators = [external.maxima_integrator,
external.giac_integrator,
external.sympy_integrator]

[Emmanuel Charpentier]

It is difficult to know before the fact which library will returnthe expected result. Consider :

```

sage: var("x, a, b")

(x, a,

b)
sage: dgamma(x, a, b)=x^(a-1)*(1-x)^(b-1)/beta(a, b)
sage: with assuming(a>0, b>0): dgamma(x, a, b).integrate(x,0,1)
1

This is quasi-immediate. And of course so is :

sage: with assuming(a>0, b>0): dgamma(x, a+1, b+1).integrate(x,0,1)
1

But

sage: with assuming(a>0, b>0): dgamma(x, a+1, b+1).integrate(x,0,1, algorithm=”sympy”)
gamma(a + 1)hypergeometric((-b, a + 1), (a + 2,), 1)/(beta(a + 1, b + 1)gamma(a + 2))

This is slow ; furthermore, simplifying this expression to 1 isn't trivialin `sage`, the best way being `sympy.simplify`... Worse :

sage: with assuming(a>0, b>0): dgamma(x, a, b).integrate(x,0,1, algorithm=”sympy”)
```
never returns…

I wouldn’t (currently) switch default-to-maxima for default-to-sympy…

HTH,

​