Issue with symbolic multiple integration (upstream?)

[kedlaya]

This symbolic double integral works fine:

sage: var('x,y')
(x, y)
sage: integral(integral(1, (y, 0, sqrt((1-x^2)/2))), (x,0,1))
1/8*sqrt(2)*pi

but the analogous triple integral breaks:

sage: var('x,y,z')
(x, y, z)
sage: assume(x<1,x>-1) # Complains otherwise
sage: integral(integral(integral(1, (z, 0, sqrt((1-x^2-2*y^2)/3))), (y,0,sqrt((1-x^2)/2))), (x, 0, 1))
RuntimeError                              Traceback (most recent call last)
...
RuntimeError: ECL says: Error executing code in Maxima: expt: undefined: 0 to a negative exponent.

On the other hand, even for the same region of integration this does not arise for other choices of integrand:

sage: var('x,y,z')
(x, y, z)
sage: integral(integral(integral(x*y, (z, 0, sqrt((1-x^2-2*y^2)/3))), (y,0,sqrt((1-x^2)/2))), (x, 0, 1))
1/90*sqrt(3)

Is this an upstream issue?

Upstream: Reported upstream. Developers acknowledge bug.

Component: calculus

Keywords: integral

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

[kcrisman]

comment:1

It seems that this might be a Maxima antiderivative issue.

Z = integral(integral(1, (z, 0, sqrt((1-x^2-2*y^2)/3))), (y,0,sqrt((1-x^2)/2)))
integral(Z,x)  
<same explosion>

If we imitate the Sage Maxima,

domain:complex;
X:integrate(1,z,0,sqrt((1-x^2-2*y^2)/3));
assume(x<1); 
assume(x>-1);
Y:integrate(X,y,0,sqrt((1-x^2)/2));
display2d:false;
integrate(Y,x,0,1);

This reveals

%pi/(3*2^(3/2))

Is this what you would expect?

I think we got rid of the abs_integrate package so I would not expect that to be part of the problem (try load(abs_integrate) to test), and I did get some related errors in some intermediate steps, but I am not finding this particular error in Maxima - yet. Hopefully this initial debug helps somebody track it down further.

[mkoeppe]

comment:2

Moving tickets to milestone sage-9.2 based on a review of last modification date, branch status, and severity.

[mkoeppe]

comment:4

Moving to 9.4, as 9.3 has been released.

[fchapoton]

Changed keywords from none to integral

[fchapoton]

comment:8

isolation of the issue in maxima integration:

sage: f=-1/12*(sqrt(3)*sqrt(2)*x^2 - sqrt(3)*sqrt(2))*arcsin(4/3*sqrt(-1/2*x^2 +1/2)/sqrt(-8/9*x^2 + 8/9))
sage: f.integrate(x,0,1)
1/36*sqrt(6)*pi
sage: f.integrate(x,0,1,algorithm='maxima')
BOUM!!

[fchapoton]

comment:9

simplified bug

sage: f = x^2*arcsin(1/3*sqrt(-x^2 + 1)/sqrt(-1/9*x^2 + 1/9))
sage: f.integrate(x,0,1,algorithm='maxima')

where in fact the square roots inside do cancel each other.

[kcrisman]

comment:10

If you can try this inside of the most recent Sage Maxima (I think mine might be a little old), then we can report that upstream.

Maxima 5.42.2 http://maxima.sourceforge.net
using Lisp ECL 16.1.2
Distributed under the GNU Public License. See the file COPYING.
Dedicated to the memory of William Schelter.
The function bug_report() provides bug reporting information.
(%i1) domain:complex;
(%o1)                               complex
(%i2) f:x^2*asin(1/3*sqrt(-x^2 + 1)/sqrt(-1/9*x^2 + 1/9));
                                               2
                             2       sqrt(1 - x )
(%o2)                       x  asin(--------------)
                                                2
                                           1   x
                                    3 sqrt(- - --)
                                           9   9
(%i3) integrate(f,x,0,1);

expt: undefined: 0 to a negative exponent.
 -- an error. To debug this try: debugmode(true);

Nice tracking this down.

[fchapoton]

comment:11

here it goes

Maxima 5.45.0 https://maxima.sourceforge.io
using Lisp ECL 21.2.1
Distributed under the GNU Public License. See the file COPYING.
Dedicated to the memory of William Schelter.
The function bug_report() provides bug reporting information.
(%i1) domain:complex;
(%o1)                               complex
(%i2) f:x^2*asin(1/3*sqrt(-x^2 + 1)/sqrt(-1/9*x^2 + 1/9));
                                               2
                             2       sqrt(1 - x )
(%o2)                       x  asin(--------------)
                                                2
                                           1   x
                                    3 sqrt(- - --)
                                           9   9
(%i3) integrate(f,x,0,1);

expt: undefined: 0 to a negative exponent.
 -- an error. To debug this try: debugmode(true);

[kcrisman]

Upstream: Reported upstream. No feedback yet.

[kcrisman]

comment:12

Thanks! This is now Maxima bug 3906.

[fchapoton]

comment:13

one can even remove the x**2 factor, the bug remains the same.

[fchapoton]

comment:14

on the other hand,

maxima: f:asin(2*sqrt(-x + 1)/sqrt(-4*x + 4));
asin((2*sqrt(1-x))/sqrt(4-4*x))
maxima: integrate(f,x,0,1);

gives another error:

CODE:
    integrate(f,x,0,1);
Maxima ERROR:

`quotient' by `zero'

[kcrisman]

Changed upstream from Reported upstream. No feedback yet. to Reported upstream. Developers acknowledge bug.