#3921 Expanded subtracted from unexpanded with e^ix does not integrate to zero

[kcrisman]

The following expressions to be integrated is just zero - an expanded subtracted from expanded.

Maxima 5.44.0 http://maxima.sourceforge.net
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:%e^-(3*%i*x)*(sin(3*x)+x)-%e^-(3*%i*x)*sin(3*x)-x*%e^-(3*%i*x);
          - 3 %i x                    - 3 %i x                - 3 %i x
(%o2)   %e         (sin(3 x) + x) - %e         sin(3 x) - x %e
(%i3) integrate(f,x,0,2*%pi);
(%o3)                               %i %pi

Presumably some kind of branch cut issue? Note even more bizarrely:

(%i3) integrate(f,x,0,2);        
(%o3)                                  0
(%i4) integrate(f,x,0,1);
(%o4)                                  0
(%i5) integrate(f,x,0,2*%pi);
(%o5)                               %i %pi
(%i6) integrate(f,x,0,-%pi);
                                     %i %pi
(%o6)                              - ------
                                       3
(%i7) integrate(f,x,0,3*%pi);
                                   4 %i %pi
(%o7)                              --------
                                      3
(%i8) integrate(f,x,-%pi,%pi);
(%o8)                                  0
(%i9) integrate(f,x,-%pi,2*%pi);
                                   4 %i %pi
(%o9)                              --------
(%i10) integrate(f,x,-1,1);
(%o10)                                 0

This is relatively recent:

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;
(%o1)                                real
(%i2) domain:complex;
(%o2)                               complex
(%i3) f:%e^-(3*%i*x)*(sin(3*x)+x)-%e^-(3*%i*x)*sin(3*x)-x*%e^-(3*%i*x);
          - 3 %i x                    - 3 %i x                - 3 %i x
(%o3)   %e         (sin(3 x) + x) - %e         sin(3 x) - x %e
(%i5) integrate(f,x,0,2*%pi);
(%o5)                                  0

See downstream at Sage ticket 33034.

[Robert Dodier]

Hmm, it's not clear to me that domain:complex makes a difference, maybe I'm mistaken.

I see the bug in 5.43.2 + Clisp, but not in 5.37.2 + Clisp. Not sure what's going on there.

Anyway it seems feasible to do git bisect given that the bug seems to have appeared in the last few years.

[kcrisman]

Yeah, the domain probably doesn't - I just included it because Sage automatically invokes it. Since it doesn't happen in 5.42.2 but does in 5.43.2, that would sound like a fairly narrow window!

[Robert Dodier]

After some investigation, so far it looks to me like the bug is in function SCRAT in src/defint.lisp. SCRAT tries to construct a rational function and then computes residues. I don't know yet what part of the process is going off the tracks here.

SCRAT is only called when the limits of integration are 0 to pi or 2 pi, so if the limits of integration are variables, some other function is called, and the integral turns out correctly.

(%i2) trace(?scrat);
(%o2)                        [scrat]
(%i3) defint(%e^(%i*x)*sin(x), x, a, b);
                         2 %i a     2 %i b
                       %e         %e
(%o3)                  -------- - --------
                          4          4
(%i4) %, a=0, b=2*%pi;
(%o4)                           0

But if the limits are 0 to 2 pi, then SCRAT is called, and it returns an incorrect result.

(%i5) defint(%e^(%i*x)*sin(x), x, 0, 2*%pi);
                 %i x
1 Enter scrat [%e     sin(x), 2 %pi]
1 Exit  scrat %i %pi
(%o5)                        %i %pi

I don't know enough about the residue calculus to figure this out, so I'll ask for help on the mailing list.

SCRAT hasn't been touched in many years. If the behavior varies by Maxima versions, I suspect it's because stuff farther up the call chain did or did not call SCRAT.

[Robert Dodier]

False alarm -- the SCRAT case mentioned isn't the problem, it's a correct result. So the error must be elsewhere. I'll keep looking.

[Robert Dodier]

I think I have a bug fix for SUPEREXPT in src/sin.lisp. I'll test this out.

diff --git a/src/sin.lisp b/src/sin.lisp
index e6d78daea..33ca59dc6 100644
--- a/src/sin.lisp
+++ b/src/sin.lisp
@@ -766,7 +766,7 @@
       (exptflag nil)) ; When T, the substitution is not possible.

   (defun superexpt (exp var bas1 pow1)
-    (prog (y (new-var (gensym "NEW-VAR-")))
+    (prog (y ($logabs nil) (new-var (gensym "NEW-VAR-")))
       (putprop new-var t 'internal)
       (setq base bas1
             pow pow1

[Robert Dodier]

I've applied the patch shown above as commit d243bce. Now the various cases shown in the original report all evaluate to 0 as expected (these have been added to the test suite). Closing this report as fixed.

[kcrisman]

Wow, fast work, thanks! I'll let downstream know - let us know when this is in a stable Maxima release.