Ask Your Question
2

Possible bug in Sage-Giac integration interface needs confirmation

asked 2019-12-25 22:04:33 +0100

dsejas gravatar image

Hello, Sage Community! Trying to answer this related question, I think I have found a bug along with user @Nasser. Although my answer there explains what I think to be wrong, let me summarize it here.

When calling

integrate(x+1,x, algorithm="giac")

Sage works perfectly, but when calling

integrate((1-2*x^(1/3))^(3/4)/x,x, algorithm="giac")

we get the following traceback:

---------------------------------------------------------------------------
SyntaxError                               Traceback (most recent call last)
/Scientific/SageMath/local/lib/python3.7/site-packages/sage/interfaces/giac.py in _sage_(self, locals)
   1101                 return symbolic_expression_from_string(result, lsymbols,
-> 1102                     accept_sequence=True)
   1103 

/Scientific/SageMath/local/lib/python3.7/site-packages/sage/calculus/calculus.py in symbolic_expression_from_string(s, syms, accept_sequence)
   2384             _augmented_syms = syms
-> 2385             return parse_func(s)
   2386         finally:

/Scientific/SageMath/local/lib/python3.7/site-packages/sage/misc/parser.pyx in sage.misc.parser.Parser.parse_sequence (build/cythonized/sage/misc/parser.c:5479)()
    537 
--> 538     cpdef parse_sequence(self, s):
    539         """

/Scientific/SageMath/local/lib/python3.7/site-packages/sage/misc/parser.pyx in sage.misc.parser.Parser.parse_sequence (build/cythonized/sage/misc/parser.c:5369)()
    555         if tokens.next() != EOS:
--> 556             self.parse_error(tokens)
    557         if len(all) == 1 and isinstance(all, list):

/Scientific/SageMath/local/lib/python3.7/site-packages/sage/misc/parser.pyx in sage.misc.parser.Parser.parse_error (build/cythonized/sage/misc/parser.c:9742)()
   1006     cdef parse_error(self, Tokenizer tokens, msg="Malformed expression"):
-> 1007         raise SyntaxError(msg, tokens.s, tokens.pos)
   1008 

SyntaxError: Malformed expression

During handling of the above exception, another exception occurred:

NotImplementedError                       Traceback (most recent call last)
<ipython-input-92-987ddabbc645> in <module>()
----> 1 integrate((Integer(1)-Integer(2)*x**(Integer(1)/Integer(3)))**(Integer(3)/Integer(4))/x,x, algorithm="giac")

/Scientific/SageMath/local/lib/python3.7/site-packages/sage/misc/functional.py in integral(x, *args, **kwds)
    751     """
    752     if hasattr(x, 'integral'):
--> 753         return x.integral(*args, **kwds)
    754     else:
    755         from sage.symbolic.ring import SR

/Scientific/SageMath/local/lib/python3.7/site-packages/sage/symbolic/expression.pyx in sage.symbolic.expression.Expression.integral (build/cythonized/sage/symbolic/expression.cpp:64541)()
  12370                     R = ring.SR
  12371             return R(integral(f, v, a, b, **kwds))
> 12372         return integral(self, *args, **kwds)
  12373 
  12374     integrate = integral

/Scientific/SageMath/local/lib/python3.7/site-packages/sage/symbolic/integration/integral.py in integrate(expression, v, a, b, algorithm, hold)
    917         if not integrator:
    918             raise ValueError("Unknown algorithm: %s" % algorithm)
--> 919         return integrator(expression, v, a, b)
    920     if a is None:
    921         return indefinite_integral(expression, v, hold=hold)

/Scientific/SageMath/local/lib/python3.7/site-packages/sage/symbolic/integration/external.py in giac_integrator(expression, v, a, b)
    430         return expression.integrate(v, a, b, hold=True)
    431     else:
--> 432         return result._sage_()

/Scientific/SageMath/local/lib/python3.7/site-packages/sage/interfaces/giac.py in _sage_(self, locals)
   1103 
   1104             except Exception:
-> 1105                 raise NotImplementedError("Unable to parse Giac output: %s" % result)
   1106         else:
   1107             return [entry.sage() for entry in self]

NotImplementedError: Unable to parse Giac output: Evaluation time: 1.18
12*(1/4*ln(abs((-2*x^(1/3)+1)^(1/4)-1))-1/4*ln((-2*x^(1/3)+1)^(1/4)+1)+1/2*atan((-2*x^(1/3)+1)^(1/4))+1/3*((-2*x^(1/3)+1)^(1/4))^3)

Here is what I could deduce after reading a lot of source code and playing with examples: When we call something like this in Sage, the following process occurs at low level:

ex = (x+1)._giac_()
result = ex.integrate(x._giac_())
result._sage_()

The first line converts the Sage expression to Giac representation; the second calls the Giac integrate function, with respect to the variable x, which must also be converted to Giac representation (that's the x._giac_()); finally, the third line converts the result back to Sage representation (for example, replacing ln from Giac to log from Sage). This last step seems to be the problem. In this particular case, print(result) shows

x^2/2 + x

However, with a more complex expression, which takes more time to be integrated:

ex=((1-2*x^(1/3))^(3/4)/x)._giac_()
result =  ex.integrate(x._giac_())
print(result)

the same process gives a different result:

Evaluation time: 1.19
12*(1/4*ln(abs((-2*x^(1/3)+1)^(1/4)-1))-1/4*ln((-2*x^(1/3)+1)^(1/4)+1)+1/2*atan((-2*x^(1/3)+1)^(1/4))+1/3*((-2*x^(1/3)+1)^(1/4))^3)

Notice the "Evaluation time" bit that was absent in the previous example? This is Giac saying it took more time than usual to perform the integration. When that string makes presence is when Sage fails:

result._sage_()

will raise the exception in the traceback.

Can somebody confirm this is a bug?

edit retag flag offensive close merge delete

1 Answer

Sort by ยป oldest newest most voted
3

answered 2019-12-26 18:59:14 +0100

Emmanuel Charpentier gravatar image

updated 2021-01-03 19:53:59 +0100

I can repoduce the problem. However, it seems it is not (as I thought initially) in the Giac/Sage interface, since the following workaround works:

sage: Ex=(1-2*x^(1/3))^(3/4)/x
sage: from giacpy_sage import *
sage: Ig=libgiac.integrate(Ex,x).sage();Ig
4*(-2*x^(1/3) + 1)^(3/4) + 6*arctan((-2*x^(1/3) + 1)^(1/4)) - 3*log((-2*x^(1/3) + 1)^(1/4) + 1) + 3*log(abs((-2*x^(1/3) + 1)^(1/4) - 1))

The presence of the timing information given by Giac in the expression to convert is highly unwelcome...

This is now Trac#28913.

EDIT : Trac#28913 is now fixed.

edit flag offensive delete link more

Comments

Thank you very much, @Emmanuel Charpentier! I am wondering whether my hypothesis about the "Evaluation time" message causing this bug was in fact correct.

dsejas gravatar imagedsejas ( 2021-01-04 06:10:13 +0100 )edit

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

1 follower

Stats

Asked: 2019-12-25 22:04:33 +0100

Seen: 702 times

Last updated: Jan 03 '21