Integrating the sgn() function can produce incorrect results #11590

orlitzky · 2011-07-12T20:05:19Z

Actual result:

sage: integrate(x * sgn(x^2 - 1/4), x, -1, 0)
1/2

Since the argument to sgn() has only one root, -1/2, on (-1, 0), there are only two subintervals on which sgn() can have different values. In particular,

sage: sgn(x^2 - 1/4)(x = -0.75)
1

sage: sgn(x^2 - 1/4)(x = -0.25)
-1

Now, the original, actual result should be equivalent to the sum of the following:

sage: integrate(x, x, -1, -0.5)
-0.375

sage: integrate(-x, x, -0.5, 0)
0.125

So, something went wrong during the initial integration.

Upstream: Fixed upstream, in a later stable release.

CC: @nbruin @fchapoton

Component: calculus

Stopgaps: #12731

Author: Michael Orlitzky

Branch/Commit: bf1b50b

Reviewer: Frédéric Chapoton

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

The text was updated successfully, but these errors were encountered:

kcrisman · 2011-07-29T18:30:49Z

comment:2

This appears to be in Maxima, and is reported at their bug tracker.


(%i11) integrate(x*signum(x^2-1/4),x,-1,0);
                                       1
(%o11)                                 -
                                       2
(%i15) integrate(x*signum(x^2-1/4),x);
                                   ! 2   1!
                                   !x  - -!
                                   !     4!
(%o15)                             --------
                                      2

kcrisman · 2011-07-29T18:30:49Z

Upstream: Reported upstream. Little or no feedback.

kcrisman · 2012-01-28T02:50:59Z

comment:3

Even with #11483, this isn't working right. It should, though - see Barton's post at the Maxima ticket:

By the way:

(%i4) load(abs_integrate)$

Correct antiderivative:

(%i5) 'integrate(x*signum(x^2-1/4),x);
(%o5) abs(x^2-1/4)/2

Correct definite integral

(%i6) 'integrate(x*signum(x^2-1/4),x,-1,0);
(%o6) -1/4

so I'm not sure why this is still returning the "wrong" thing. Probably something about the integration code...

kcrisman · 2012-01-28T02:50:59Z

Changed upstream from Reported upstream. Little or no feedback. to Reported upstream. Developers acknowledge bug.

orlitzky · 2012-01-30T20:37:37Z

comment:4

I'm stumped on this one. We get the correct antiderivative:

sage: integrate(x*sgn(x^2-1/4),x)     
1/2*abs(x^2 - 1/4)

And ECL gives us the right answer, so it's not an environment setting:

sage: from sage.interfaces.maxima_lib import ecl_eval
sage: ecl_eval("#$'integrate(x*signum(x^2-1/4),x,-1,0);$")                           
<ECL: ((RAT SIMP) -1 4)>

But going through maxima_eval is still trouble:

sage: integrate(x*sgn(x^2-1/4),x,-1,0)                                            
1/2
sage: from sage.interfaces.maxima_lib import maxima_eval
sage: a = '($INTEGRATE ((MTIMES SIMP) $X ((%SIGNUM SIMP) ((MPLUS SIMP) ((RAT SIMP) (- 1) 4) ((MEXPT SIMP) $X 2))) ) $X -1 0)'
sage: maxima_eval(a)
<ECL: ((RAT SIMP) 1 2)>

kcrisman · 2012-01-30T20:41:59Z

comment:5

Well, I'm glad it wasn't just me!

Could it be that because we get an answer without abs_integrate it just returns the 'regular' answer? But I don't recall the integrate code doing that, I think once we turn on abs_integrate it should just 'be on'...

orlitzky · 2012-01-30T21:23:15Z

comment:6

Replying to @kcrisman:

Well, I'm glad it wasn't just me!

Could it be that because we get an answer without abs_integrate it just returns the 'regular' answer? But I don't recall the integrate code doing that, I think once we turn on abs_integrate it should just 'be on'...

It's possible. The abs_integrate code looks like it defines extra definite integration methods, and then loops through them until 'integrate is gone from the expression. Within maxima, the extra methods obviously get tried first, because we get the right answer. But I suppose something in the ECL integration could be trying the default definite integration first.

orlitzky · 2012-01-31T01:14:11Z

comment:7

Something else is weird here. In standalone maxima-5.26, we don't get back the wrong answer. But through sage, without abs_integrate, we get 1/2:

sage: integrate(x * sgn(x^2 - 1/4), x, -1, 0) 
1/2
sage: maxima.console()
;;; Loading #P"/home/mjo/src/sage-5.0.beta1/local/lib/ecl/sb-bsd-sockets.fas"
;;; Loading #P"/home/mjo/src/sage-5.0.beta1/local/lib/ecl/sockets.fas"
;;; Loading #P"/home/mjo/src/sage-5.0.beta1/local/lib/ecl/defsystem.fas"
;;; Loading #P"/home/mjo/src/sage-5.0.beta1/local/lib/ecl/cmp.fas"
Maxima 5.26.0 http://maxima.sourceforge.net
using Lisp ECL 11.1.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) display2d:false;

(%o1) false
(%i2) 'integrate(x*signum(x^2-1/4),x);

(%o2) 'integrate(x*signum(x^2-1/4),x)
(%i3) 'integrate(x*signum(x^2-1/4),x,-1,0);

(%o3) 'integrate(x*signum(x^2-1/4),x,-1,0)

Intriguing!

nbruin · 2012-01-31T21:29:55Z

comment:8

Some random observations that may or may not be relevant. First, let's separate the "Reader" (parser) and the "Evaluator".

sage: from sage.libs.ecl import *
sage: I1=EclObject("#$integrate(x*signum(x^2-1/4),x,-1,0);$").cadr().cdr().car() 
sage: I2=EclObject("#$'integrate(x*signum(x^2-1/4),x,-1,0);$").cadr().cdr().car()
sage: I1
<ECL: (($INTEGRATE)
 ((MTIMES) $X
  ((%SIGNUM) ((MPLUS) ((MEXPT) $X 2) ((MMINUS) ((MQUOTIENT) 1 4)))))
 $X ((MMINUS) 1) 0)>
sage: I2
<ECL: ((%INTEGRATE)
 ((MTIMES) $X
  ((%SIGNUM) ((MPLUS) ((MEXPT) $X 2) ((MMINUS) ((MQUOTIENT) 1 4)))))
 $X ((MMINUS) 1) 0)>

Note the subtle difference between $INTEGRATE (function) %INTEGRATE (inert integral form). Sage's sr_integral produces the $INTEGRATE expression, so the two alternatives tested are not equivalent. I suspect those trigger different code-paths. You could try changing

interfaces/maxima_lib.py:205:max_integrate=EclObject("$INTEGRATE")

to max_integrate=EclObject("%INTEGRATE") instead and see if that solves more than it messes up.

Slightly less likely to make a difference: If you look at the full expression resulting for EclObject("#$1+2$"), you'll notice an MEVAL*, whereas maxima_eval is defined with MEVAL. I've never worked out which one is the proper one to use. So, you should test all combinations (these skip setting up the proper maxima-error-catching environment):

sage: EclObject(["meval*",["quote", I1]]).eval()
sage: EclObject(["meval*",["quote", I2]]).eval()
sage: EclObject(["meval",["quote", I1]]).eval()
sage: EclObject(["meval",["quote", I2]]).eval()

orlitzky · 2012-01-31T23:48:07Z

comment:9

I took a much stupider approach, but I did essentially try replacing $INTEGRATE with %INTEGRATE. Once I discovered the (#$$) syntax, I tried replacing the call to maxima_eval in sr_integral with a simple string substitution. It causes more problems than it solves, I think.

If I replace meval with meval* in a few places, I'm able to break some additional things but not fix any =)

I guess the Right Thing to do would be to go read the maxima source...

nbruin · 2012-02-01T01:06:36Z

comment:10

Replying to @orlitzky:

I tried replacing the call to maxima_eval in sr_integral with a simple string substitution. It causes more problems than it solves, I think.

Yes, that would make things much worse (not to mention a major step back in other ways). The sr_integral interface tries to avoid string-based conversion, but maxima_lib in general is still capable of it in a much better way than crafting your own string-mangling based on #$...$. For instance:

sage.calculus.calculus.maxima(sage.calculus.calculus.dummy_integrate(x * sgn(x^2 - 1/4), x, -1, 0))

should work and compares nicely to the pexpect-based

maxima(sage.calculus.calculus.dummy_integrate(x * sgn(x^2 - 1/4), x, -1, 0))

You should really try the commands I suggested earlier to see where the problem is, though (no changes to the source required!). I don't have ready access to a sage build with the right patches, but you obviously do. Even on a clearly insufficiently patched sage I observe different behaviour:

sage: maxima_eval(I1)
<ECL: ((RAT SIMP) 1 2)>
sage: maxima_eval(I2)
RuntimeError: ECL says: Error executing code in Maxima: first: empty argument.

orlitzky · 2012-02-01T01:29:28Z

comment:11

Sorry, I should have been more clear: I did try everything you suggested. Replacing $INTEGRATE with %INTEGRATE causes the same test failures in maxima_lib.py that I got when I swapped out the whole thing for string substitution. It does fix this particular answer, though.

Here's a session with just the abs_integrate patch applied, maxima-5.26. I left the error at the top intact, for whatever reason it doesn't work at all until after I've integrated something.

sage: from sage.libs.ecl import *                                               
sage: I1=EclObject("#$integrate(x*signum(x^2-1/4),x,-1,0);$").cadr().cdr().car()
---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)

/home/mjo/src/sage-5.0.beta1/devel/sage-devel/<ipython console> in <module>()

/home/mjo/src/sage-5.0.beta1/local/lib/python2.7/site-packages/sage/libs/ecl.so in sage.libs.ecl.EclObject.cadr (sage/libs/ecl.c:5528)()

TypeError: cadr can only be applied to a cons
sage: integrate(x*sgn(x^2-1/4),x,-1,0)                                          
1/2
sage: I1=EclObject("#$integrate(x*signum(x^2-1/4),x,-1,0);$").cadr().cdr().car()
sage: I2=EclObject("#$'integrate(x*signum(x^2-1/4),x,-1,0);$").cadr().cdr().car()
sage: I1
<ECL: (($INTEGRATE)
 ((MTIMES) $X
  ((%SIGNUM) ((MPLUS) ((MEXPT) $X 2) ((MMINUS) ((MQUOTIENT) 1 4)))))
 $X ((MMINUS) 1) 0)>
sage: I2
<ECL: ((%INTEGRATE)
 ((MTIMES) $X
  ((%SIGNUM) ((MPLUS) ((MEXPT) $X 2) ((MMINUS) ((MQUOTIENT) 1 4)))))
 $X ((MMINUS) 1) 0)>
sage: EclObject(["meval*",["quote", I1]]).eval()
<ECL: ((RAT SIMP) 1 2)>
sage: EclObject(["meval*",["quote", I2]]).eval()
<ECL: ((RAT SIMP) -1 4)>
sage: EclObject(["meval",["quote", I1]]).eval()
<ECL: ((RAT SIMP) 1 2)>
sage: EclObject(["meval",["quote", I2]]).eval()
<ECL: ((RAT SIMP) -1 4)>

And here's a session with $INTEGRATE switched to %INTEGRATE:

sage: from sage.libs.ecl import *                                                
sage: I1=EclObject("#$integrate(x*signum(x^2-1/4),x,-1,0);$").cadr().cdr().car()
---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)

/home/mjo/src/sage-5.0.beta1/devel/sage-devel/<ipython console> in <module>()

/home/mjo/src/sage-5.0.beta1/local/lib/python2.7/site-packages/sage/libs/ecl.so in sage.libs.ecl.EclObject.cadr (sage/libs/ecl.c:5528)()

TypeError: cadr can only be applied to a cons
sage: integrate(x*sgn(x^2-1/4),x,-1,0)                                           
-1/4
sage: I1=EclObject("#$integrate(x*signum(x^2-1/4),x,-1,0);$").cadr().cdr().car() 
sage: I2=EclObject("#$'integrate(x*signum(x^2-1/4),x,-1,0);$").cadr().cdr().car()
sage: I1
<ECL: (($INTEGRATE)
 ((MTIMES) $X
  ((%SIGNUM) ((MPLUS) ((MEXPT) $X 2) ((MMINUS) ((MQUOTIENT) 1 4)))))
 $X ((MMINUS) 1) 0)>
sage: I2
<ECL: ((%INTEGRATE)
 ((MTIMES) $X
  ((%SIGNUM) ((MPLUS) ((MEXPT) $X 2) ((MMINUS) ((MQUOTIENT) 1 4)))))
 $X ((MMINUS) 1) 0)>
sage: EclObject(["meval*",["quote", I1]]).eval()                                 
<ECL: ((RAT SIMP) 1 2)>
sage: EclObject(["meval*",["quote", I2]]).eval()                                 
<ECL: ((RAT SIMP) -1 4)>
sage: EclObject(["meval",["quote", I1]]).eval()                                  
<ECL: ((RAT SIMP) 1 2)>
sage: EclObject(["meval",["quote", I2]]).eval()                                  
<ECL: ((RAT SIMP) -1 4)>

The same thing, except the original integration actually works. Problem is, other tests start to fail with %INTEGRATE. Here's an example:

$ sage -t sage/interfaces/maxima_lib.py
sage -t  "devel/sage-devel/sage/interfaces/maxima_lib.py"   
**********************************************************************
File "/home/mjo/src/sage-5.0.beta1/devel/sage-devel/sage/interfaces/maxima_lib.py", line 661:
    sage: integral(x^n,x)
Expected:
    Traceback (most recent call last):
    ...
    ValueError: Computation failed since Maxima requested additional
    constraints; using the 'assume' command before integral evaluation
    *may* help (example of legal syntax is 'assume(n+1>0)',
    see `assume?` for more details)
    Is  n+1  zero or nonzero?
Got:
    x^(n + 1)/(n + 1)
**********************************************************************
File "/home/mjo/src/sage-5.0.beta1/devel/sage-devel/sage/interfaces/maxima_lib.py", line 685:
    sage: integrate(sgn(x) - sgn(1-x), x)
Expected:
    abs(x - 1) + abs(x)
Got:
    (x - 1)*sgn(x - 1) + x*sgn(x)
**********************************************************************
File "/home/mjo/src/sage-5.0.beta1/devel/sage-devel/sage/interfaces/maxima_lib.py", line 700:
    sage: integrate(cos(x + abs(x)), x)
Expected:
    1/4*(sgn(x) + 1)*sin(2*x) - 1/2*x*sgn(x) + 1/2*x
Got:
    -1/4*(2*x - sin(2*x))*sgn(x) + 1/2*x + 1/4*sin(2*x)
**********************************************************************
File "/home/mjo/src/sage-5.0.beta1/devel/sage-devel/sage/interfaces/maxima_lib.py", line 711:
    sage: integral(abs(cos(x))*sin(x),(x,pi/2,pi))
Expected:
    1/2
Got:
    -1/2

nbruin · 2012-02-01T01:59:44Z

comment:12

Replying to @orlitzky:

I left the error at the top intact, for whatever reason it doesn't work at all until after I've integrated something.

sage: from sage.libs.ecl import *

This doesn't load maxima inside ecl yet, but once you evaluate an integral, maxima is loaded and the #$...$ macro works. Sorry about that.

nbruin · 2012-02-01T02:06:40Z

comment:13

Replying to @kcrisman:

Even with #11483, this isn't working right. It should, though - see Barton's post at the Maxima ticket:
By the way:

(%i4) load(abs_integrate)$

Correct antiderivative:

(%i5) 'integrate(x*signum(x^2-1/4),x);
(%o5) abs(x^2-1/4)/2

Correct definite integral

(%i6) 'integrate(x*signum(x^2-1/4),x,-1,0);
(%o6) -1/4
so I'm not sure why this is still returning the "wrong" thing.

Did you try it with integrate rather than 'integrate? Given what we've seen elsewhere in this ticket, I suspect it still gives the wrong answer. In sage, we are interfacing with integrate. If that is still broken, then the bug is not fixed as far as sage is concerned.

If maxima's position is that we should use a different integration routine (i.e., 'integrate) then we need a heuristic on when to use what routine ... isn't it maxima's job to figure this out?

kcrisman · 2012-02-01T02:28:03Z

comment:14

Did you try it with integrate rather than 'integrate?

(%i5) display2d:false;

(%o5) false
(%i6) integrate(x*signum(x^2-1/4),x,-1,0);

(%o6) 1/2
(%i7) 'integrate(x*signum(x^2-1/4),x,-1,0);

(%o7) -1/4

Note that this was just Barton Willis (Maxima dev) example.

Given what we've seen elsewhere in this ticket, I suspect it still gives the wrong answer. In sage, we are interfacing with integrate. If that is still broken, then the bug is not fixed as far as sage is concerned.

I don't think they claimed it was fixed.

If maxima's position is that we should use a different integration routine (i.e., 'integrate) then we need a heuristic on when to use what routine ... isn't it maxima's job to figure this out?

Yes, I am a little stumped on this. I thought that the apostrophe just meant "don't evaluate nounform", but apparently I was mistaken. I just reported this example on the Maxima ticket, but I wouldn't hold my breath waiting for it, since there is a way to get it to work correctly there.

roed314 · 2012-03-22T19:01:56Z

Stopgaps: #12731

kcrisman · 2014-12-08T15:00:45Z

comment:16

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Integrating the sgn() function can produce incorrect results #11590

Integrating the sgn() function can produce incorrect results #11590

orlitzky commented Jul 12, 2011

kcrisman commented Jul 29, 2011

kcrisman commented Jul 29, 2011

kcrisman commented Jan 28, 2012

kcrisman commented Jan 28, 2012

orlitzky commented Jan 30, 2012

kcrisman commented Jan 30, 2012

orlitzky commented Jan 30, 2012

orlitzky commented Jan 31, 2012

nbruin commented Jan 31, 2012

orlitzky commented Jan 31, 2012

nbruin commented Feb 1, 2012

orlitzky commented Feb 1, 2012

nbruin commented Feb 1, 2012

nbruin commented Feb 1, 2012

kcrisman commented Feb 1, 2012

roed314 commented Mar 22, 2012

kcrisman commented Dec 8, 2014

kcrisman commented Dec 8, 2014

sagetrac-jakobkroeker mannequin commented Feb 10, 2015

sagetrac-robert-dodier mannequin commented Apr 14, 2016

kcrisman commented Apr 14, 2016

fchapoton commented Aug 26, 2019

orlitzky commented Dec 18, 2020

orlitzky commented Dec 18, 2020

orlitzky commented Dec 18, 2020

orlitzky commented Dec 18, 2020

fchapoton commented Dec 19, 2020

fchapoton commented Dec 19, 2020

vbraun commented Dec 27, 2020

Integrating the sgn() function can produce incorrect results #11590

Integrating the sgn() function can produce incorrect results #11590

Comments

orlitzky commented Jul 12, 2011

kcrisman commented Jul 29, 2011

kcrisman commented Jul 29, 2011

kcrisman commented Jan 28, 2012

kcrisman commented Jan 28, 2012

orlitzky commented Jan 30, 2012

kcrisman commented Jan 30, 2012

orlitzky commented Jan 30, 2012

orlitzky commented Jan 31, 2012

nbruin commented Jan 31, 2012

orlitzky commented Jan 31, 2012

nbruin commented Feb 1, 2012

orlitzky commented Feb 1, 2012

nbruin commented Feb 1, 2012

nbruin commented Feb 1, 2012

kcrisman commented Feb 1, 2012

roed314 commented Mar 22, 2012

kcrisman commented Dec 8, 2014

kcrisman commented Dec 8, 2014

sagetrac-jakobkroeker mannequin commented Feb 10, 2015

sagetrac-robert-dodier mannequin commented Apr 14, 2016

kcrisman commented Apr 14, 2016

fchapoton commented Aug 26, 2019

orlitzky commented Dec 18, 2020

orlitzky commented Dec 18, 2020

orlitzky commented Dec 18, 2020

orlitzky commented Dec 18, 2020

fchapoton commented Dec 19, 2020

fchapoton commented Dec 19, 2020

vbraun commented Dec 27, 2020