division by zero from integrate to giac. Where does long giac output come from?
I found number of integrals where sagemath says ValueError: power::eval(): division by zero:
when using giac
algorithm.
But when I try the same integrate command in giac directly, it actually gives no output. Not able to integrate it. giac just says Done
.
So I wonder why then sagemath gives division by zero? It also gives very large output below that, which I do not know where it came from, since giac does not show any output.
At first I thought giac may be does not display on the screen an output which is too large? I could not find googling around anything about it. But sagemath
claims that giac
has returned this vey very long result from integrate. How is that possible?
Using giac 1.7
and sagemath 9.3
on Linux
Here is the command inside giac
10>> integrate((e*x+d)^(9/2)/(a*d*e+(a*e^2+c*d^2)*x+c*d*e*x^2),x)
Warning, need to choose a branch for the root of a polynomial with parameters. This might be wrong.
The choice was done assuming [a,c,d,exp(1),exp(2)]=[60,-94,37,-33,49]
Warning, need to choose a branch for the root of a polynomial with parameters. This might be wrong.
The choice was done assuming [a,c,d,exp(1),exp(2)]=[-3,-98,-46,47,-64]
Warning, need to choose a branch for the root of a polynomial with parameters. This might be wrong.
The choice was done assuming [a,c,d,exp(1),exp(2)]=[-22,22,-65,22,-81]
Warning, need to choose a branch for the root of a polynomial with parameters. This might be wrong.
The choice was done assuming [a,c,d,exp(1),exp(2)]=[1,6,75,-85,-12]
Evaluation time: 10.44
Done
// Time 10.44
And inside sagemath
var('e x d a c')
sage: integrate((e*x+d)^(9/2)/(a*d*e+(a*e^2+c*d^2)*x+c*d*e*x^2),x, algorithm="giac")
---------------------------------------------------------------------------
ValueError Traceback (most recent call last)
/usr/lib/python3.9/site-packages/sage/symbolic/expression.pyx in sage.symbolic.expression.Expression._div_ (build/cythonized/sage/symbolic/expression.cpp:25839)()
3894 else:
-> 3895 x = left._gobj / _right._gobj
3896 return new_Expression_from_GEx(left._parent, x)
ValueError: power::eval(): division by zero
During handling of the above exception, another exception occurred:
ZeroDivisionError Traceback (most recent call last)
/usr/lib/python3.9/site-packages/sage/interfaces/giac.py in _sage_(self, locals)
1131 try:
-> 1132 return symbolic_expression_from_string(result, lsymbols,
1133 accept_sequence=True)
/usr/lib/python3.9/site-packages/sage/calculus/calculus.py in symbolic_expression_from_string(s, syms, accept_sequence)
2407 if isinstance(v,Function)})
-> 2408 return parse_func(s)
2409
/usr/lib/python3.9/site-packages/sage/misc/parser.pyx in sage.misc.parser.Parser.parse_sequence (build/cythonized/sage/misc/parser.c:5837)()
549
--> 550 cpdef parse_sequence(self, s):
551 """
/usr/lib/python3.9/site-packages/sage/misc/parser.pyx in sage.misc.parser.Parser.parse_sequence (build/cythonized/sage/misc/parser.c:5702)()
565 cdef Tokenizer tokens = Tokenizer(s)
--> 566 all = self.p_sequence(tokens)
567 if tokens.next() != EOS:
/usr/lib/python3.9/site-packages/sage/misc/parser.pyx in sage.misc.parser.Parser.p_sequence (build/cythonized/sage/misc/parser.c:6450)()
631 elif token == '(':
--> 632 obj = self.p_tuple(tokens)
633 elif token == EOS:
/usr/lib/python3.9/site-packages/sage/misc/parser.pyx in sage.misc.parser.Parser.p_tuple (build/cythonized/sage/misc/parser.c:7237)()
699 tokens.reset(start)
--> 700 return self.p_eqn(tokens)
701
/usr/lib/python3.9/site-packages/sage/misc/parser.pyx in sage.misc.parser.Parser.p_eqn (build/cythonized/sage/misc/parser.c:7394)()
728 """
--> 729 lhs = self.p_expr(tokens)
730 cdef int op = tokens.next()
/usr/lib/python3.9/site-packages/sage/misc/parser.pyx in sage.misc.parser.Parser.p_expr (build/cythonized/sage/misc/parser.c:7746)()
768 cdef int op
--> 769 operand1 = self.p_term(tokens)
770 op = tokens.next()
/usr/lib/python3.9/site-packages/sage/misc/parser.pyx in sage.misc.parser.Parser.p_term (build/cythonized/sage/misc/parser.c:8150)()
812 else:
--> 813 operand1 = operand1 / operand2
814 op = tokens.next()
/usr/lib/python3.9/site-packages/sage/structure/element.pyx in sage.structure.element.Element.__truediv__ (build/cythonized/sage/structure/element.c:13162)()
1734 if HAVE_SAME_PARENT(cl):
-> 1735 return (<Element>left)._div_(right)
1736 if BOTH_ARE_ELEMENT(cl):
/usr/lib/python3.9/site-packages/sage/symbolic/expression.pyx in sage.symbolic.expression.Expression._div_ (build/cythonized/sage/symbolic/expression.cpp:25919)()
3900 if 'division by zero' in str(msg):
-> 3901 raise ZeroDivisionError("symbolic division by zero")
3902 else:
ZeroDivisionError: symbolic division by zero
During handling of the above exception, another exception occurred:
NotImplementedError Traceback (most recent call last)
<ipython-input-6-754567723580> in <module>
----> 1 integrate((e*x+d)**(Integer(9)/Integer(2))/(a*d*e+(a*e**Integer(2)+c*d**Integer(2))*x+c*d*e*x**Integer(2)),x, algorithm="giac")
/usr/lib/python3.9/site-packages/sage/misc/functional.py in integral(x, *args, **kwds)
757 """
758 if hasattr(x, 'integral'):
--> 759 return x.integral(*args, **kwds)
760 else:
761 from sage.symbolic.ring import SR
/usr/lib/python3.9/site-packages/sage/symbolic/expression.pyx in sage.symbolic.expression.Expression.integral (build/cythonized/sage/symbolic/expression.cpp:66867)()
12645 R = SR
12646 return R(integral(f, v, a, b, **kwds))
> 12647 return integral(self, *args, **kwds)
12648
12649 integrate = integral
/usr/lib/python3.9/site-packages/sage/symbolic/integration/integral.py in integrate(expression, v, a, b, algorithm, hold)
988 if not integrator:
989 raise ValueError("Unknown algorithm: %s" % algorithm)
--> 990 return integrator(expression, v, a, b)
991 if a is None:
992 return indefinite_integral(expression, v, hold=hold)
/usr/lib/python3.9/site-packages/sage/symbolic/integration/external.py in giac_integrator(expression, v, a, b)
446 return expression.integrate(v, a, b, hold=True)
447 else:
--> 448 return result._sage_()
/usr/lib/python3.9/site-packages/sage/interfaces/giac.py in _sage_(self, locals)
1134
1135 except Exception:
-> 1136 raise NotImplementedError("Unable to parse Giac output: %s" % result)
1137 else:
1138 return [entry.sage() for entry in self]
NotImplementedError: Unable to parse Giac output: ((-4*a^8*c^2*d^2*exp(2)^8+2*a^8*sqrt(-c^2*d^3-c*d*sqrt(c^2*d^4-4*a*c*d^2*exp(1)^2+a^2*exp(2)^2+2*a*c*d^2*exp(2))+a*c*d*exp(2))*sqrt(2)*sqrt(c^2*d^4-4*a*c*d^2*exp(1)^2+a^2*exp(2)^2+2*a*c*d^2*exp(2))*exp(2)^8+44*a^7*c^3*d^4*exp(1)^2*exp(2)^6-12*a^7*c^3*d^4*exp(2)^7-22*a^7*c*d^2*sqrt(-c^2*d^3-c*d*sqrt(c^2*d^4-4*a*c*d^2*exp(1)^2+a^2*exp(2)^2+2*a*c*d^2*exp(2))+a*c*d*exp(2))*sqrt(2)*sqrt(c^2*d^4-
... pages and pages of output like the above follows....
p(1))*d^9*c^6-4*sqrt(d+x*exp(1))*d^7*exp(1)^2*c^5*a-2*sqrt(d+x*exp(1))*d^7*c^5*a*exp(2)+4*sqrt(d+x*exp(1))*d^5*exp(1)^2*c^4*a^2*exp(2)+2*sqrt(d+x*exp(1))*d^5*c^4*a^2*exp(2)^2-2*sqrt(d+x*exp(1))*d^3*c^3*a^3*exp(2)^3)/d^7/c^7
sage:
My question is: Where did all this output come from? The ones starting after Unable to parse Giac output:
.
I know because giac has warning messages, this could have confused sagemath in parsing the return value, and this might explain the division by zero. There is a question on this already, when giac returns warning messages.
But my question here is different. Since giac
did not give any output from the integrate command, I was wondering then where did this very long output came from?
I cannot answer your main question, but have a few remarks. First, this works fine with libgiac. Secondly, there is a problem that the Giac interface cannot distinguish between the number $e$ and the variable $e$, so you may want to avoid using $e$ as a variable; see https://trac.sagemath.org/ticket/30133. Finally, this integral also works fine with Giac 1.6.0.47 (which comes with the Sage 9.3 distribution), apparently because $e$ seems to be treated differently there than with Giac 1.7.0.
this works fine with libgiac
hummm. I am using 9.3 sagemath with giac 1.7 and it does not work for me. It throws an exceptionAnd when doing
print(anti)
nothing is displayed.As for changing the letter
e
, I am running the Rubi integration test suite as is. This is the same input used for all CAS systems and not practical for me to edit these files and change them differently for each CAS and change letters. Also, insidegiac
itself, it works usinge
. But as I mentioned, I see no output in giac.Indeed, I get the same error with libgiac now. I am not sure what I have done differently before. Though, if I replace
e
by something else, then it works.Regarding the letter
e
, I agree this is a problem. Hopefully it will be fixed by the ticket I linked to above, if someone has enough interest in fixing this bug. Feel free to take a look at the source code.