Opened 10 months ago
Last modified 6 weeks ago
#30520 new defect
Error in the sign of a product
Reported by: | tmonteil | Owned by: | |
---|---|---|---|
Priority: | major | Milestone: | sage-9.4 |
Component: | symbolics | Keywords: | maxima |
Cc: | gh-DaveWitteMorris, gh-macrakis, gh-robert-dodier, kcrisman, nbruin, slelievre | Merged in: | |
Authors: | Reviewers: | ||
Report Upstream: | N/A | Work issues: | |
Branch: | Commit: | ||
Dependencies: | Stopgaps: |
Description (last modified by )
Reported at Ask Sage question 53345.
The output for this product has an incorrect leading minus sign:
sage: N, q, k = SR.var('N, q, k') sage: product(1 - q^k, k, 1, N) -(-1)^N*product(q^k - 1, k, 1, N)
Expected:
sage: N, q, k = SR.var('N, q, k') sage: product(1 - q^k, k, 1, N) (-1)^N*product(q^k - 1, k, 1, N)
Fixed in Maxima by this commit:
Change History (21)
comment:1 Changed 10 months ago by
comment:2 Changed 9 months ago by
for what it's worth, it's either in the interface to maxima in maxima_eval
or in maxima itself, what ultimately gets called is:
sage: var('n,q,i') (n, q, i) sage: expr = 1 - q**i sage: from sage.interfaces import maxima_lib sage: sage.interfaces.maxima_lib.maxima.sr_prod(expr, i, 1, n) -(-1)^n*product(q^i - 1, i, 1, n)
I couldn't reproduce this on pure maxima, but I have never used it before, so it may well be there.
comment:3 Changed 9 months ago by
A simple case
sage: from sage.calculus.calculus import symbolic_product sage: var('x,n') (x, n) sage: symbolic_product(-x**2,x,1,n) -(-1)^n*factorial(n)^2
comment:4 Changed 9 months ago by
This seems to happen in max_simplify_prod
. Before:
ipdb> p maxima_eval(([max_prod],[sr_to_max(SR(a)) for a in args])) <ECL: ((%PRODUCT SIMP) ((MTIMES SIMP) -1 ((MEXPT SIMP) |$_SAGE_VAR_x| 2)) |$_SAGE_VAR_x| 1 |$_SAGE_VAR_n|)>
and after:
ipdb> p maxima_eval([max_simplify_prod, ([max_prod],[sr_to_max(SR(a)) for a in a rgs])]) <ECL: ((MTIMES SIMP) ((MEXPT SIMP) -1 ((MPLUS SIMP) -1 |$_SAGE_VAR_n|)) ((MEXPT SIMP) ((MFACTORIAL SIMP) |$_SAGE_VAR_n|) 2))>
comment:5 Changed 8 months ago by
- Milestone changed from sage-9.2 to sage-9.3
comment:6 Changed 3 months ago by
- Cc gh-DaveWitteMorris gh-macrakis gh-robert-dodier slelievre added
- Description modified (diff)
Any insights from Maxima or Sage symbolics experts?
comment:7 follow-up: ↓ 9 Changed 3 months ago by
Hi, thanks for the Maxima problem report.
I was unable to reproduce this in Maxima 5.44.0/SBCL 2.0.0/Windows 10.
Could you please report the details necessary to reproduce the problem:
- The exact version of Maxima/SBCL/OS you're using.
- The exact commands that were passed to Maxima.
- Any non-default global settings that Sage uses when calling Maxima.
- What does max_simplify_product do?
- What the var(...) declaration/statement does in Maxima terms. Does it declare anything special about those variables?
- What does "before this ticket" and "after this ticket" mean? Was it run on different versions of Maxima? If so, which ones?
Thanks,
-s
comment:8 Changed 3 months ago by
- Cc kcrisman nbruin added
- Keywords maxima added
Apologies for cc-ing you. Rather than in Maxima, the bug might well be in Sage's interface to it.
In Sage, the bug can lead to wrong results as in
sage: from sage.calculus.calculus import symbolic_product sage: x, n = SR.var('x,n') sage: symbolic_product(-x**2, x, 1, n) -(-1)^n*factorial(n)^2
or cause an error to be raised as in
sage: symbolic_product(-x, x, 1, n) Traceback (most recent call last) ... RuntimeError: ECL says: Error executing code in Maxima: factorial: factorial of negative integer -1 not defined.
Not sure what exactly gets sent to Maxima in those cases.
I have Maxima 5.44.0, ECL 21.2.1, macOS 10.14.6.
sage: pv = installed_packages() # package versions sage: print(', '.join(f'{p} {pv[p]}' for p in ['maxima', 'ecl'])) maxima 5.44.0, ecl 21.2.1
Sage calls Maxima with the following global settings:
besselexpand : true ; display2d : false ; domain : complex ; keepfloat : true ; load(to_poly_solve) ; load(simplify_sum) ; load(diag) ; nolabels : true ;
"Before this ticket" is the observed bug. "After this ticket" is the goal (so far not reached).
comment:9 in reply to: ↑ 7 Changed 3 months ago by
- What does max_simplify_product do?
max_simplify_prod=EclObject("$SIMPLIFY_PRODUCT")
which is presumably in the contributed simplify_sum
package, more precisely here. I had trouble finding references in the official documentation, though.
- What the var(...) declaration/statement does in Maxima terms. Does it declare anything special about those variables?
No.
comment:10 Changed 3 months ago by
"Before this ticket" is the observed bug. "After this ticket" is the goal (so far not reached).
I know this is just a convention, but, um, this is pretty confusing; I don't see how anyone who was not already among the cognoscenti would know what it means. "Observed" and "expected behavior" seems adequate.
I can't be sure about what's going on here, but I speculate the error is coming out of simplify_product
in the simplify_sum
code.
simplify_product ('product(1 - f(k), k, 1, N)); => (-1)^(N-1)*'product(f(k)-1,k,1,N)
The leading term is off by 1 (right?) -- I think this is the origin of the stray -1 which was observed.
comment:11 Changed 3 months ago by
- Description modified (diff)
Thanks for the rephrasing suggestion and your thoughts on the problem.
comment:12 Changed 3 months ago by
@gh-robert-dodier: Yes, I think you are right. In line 1112 of the file suggested by kcrisman, the code for simplify_product
says:
(-1)^(hi-lo)*simplify_product(apply(nounify('product), [-term, %n, lo, hi]))
If I understand correctly, then this is an off-by-one error: the number of terms in the product is hi - lo + 1
because hi
and lo
are both included.
comment:13 follow-up: ↓ 15 Changed 3 months ago by
I located the copy of that file in my Sage installation:
$ cd $(sage -c 'print(SAGE_ROOT)') $ find . -name "solve_rec.mac" ./local/share/maxima/5.44.0/share/solve_rec/solve_rec.mac
and changed line 1112 as suggested by Dave.
Now I get the correct product from the ticket description:
sage: N, q, k = SR.var('N, q, k') sage: product(1 - q^k, k, 1, N) (-1)^N*product(q^k - 1, k, 1, N)
and for the sum of -x^2
:
sage: from sage.calculus.calculus import symbolic_product sage: x, n = SR.var('x, n') sage: symbolic_product(-x**2, x, 1, n) (-1)^n*factorial(n)^2
but this still gives trouble:
sage: from sage.calculus.calculus import symbolic_product sage: x, n = SR.var('x, n') sage: symbolic_product(-x, x, 1, n) #0: simplify_product(product='product(-_SAGE_VAR_x,_SAGE_VAR_x,1,_SAGE_VAR_n)) ... Traceback (most recent call last) ... RuntimeError: ECL says: Error executing code in Maxima: factorial: factorial of negative integer -1 not defined.
It would be nice to fix that too.
comment:14 Changed 3 months ago by
- Description modified (diff)
I've come to the same conclusion. I think this patch fixes it:
$ git diff diff --git a/share/solve_rec/solve_rec.mac b/share/solve_rec/solve_rec.mac index 2727db808..8979e3f1d 100644 --- a/share/solve_rec/solve_rec.mac +++ b/share/solve_rec/solve_rec.mac @@ -1109,7 +1109,7 @@ simplify_product(prod) := block( prod ) else if part(term, 0)="-" and length(args(term))=1 then - (-1)^(hi-lo)*simplify_product(apply(nounify('product), [-term, %n, lo, hi])) + (-1)^(hi-lo+1)*simplify_product(apply(nounify('product), [-term, %n, lo, hi])) /* Give up! */ else p%*apply(nounify('product), [factor(term), %n, lo, hi])
Although the existing, incorrect code is exercised by the unit tests (rtest_simplify_sum and rtest_solve_rec), fixing it does not change any results! Not sure what's going on there; I guess it is a happy accident. Or not.
Anyway I will add a few unit tests and commit the patch.
comment:15 in reply to: ↑ 13 ; follow-up: ↓ 16 Changed 3 months ago by
Replying to slelievre:
but this still gives trouble:
sage: from sage.calculus.calculus import symbolic_product sage: x, n = SR.var('x, n') sage: symbolic_product(-x, x, 1, n) #0: simplify_product(product='product(-_SAGE_VAR_x,_SAGE_VAR_x,1,_SAGE_VAR_n)) ... Traceback (most recent call last) ... RuntimeError: ECL says: Error executing code in Maxima: factorial: factorial of negative integer -1 not defined.It would be nice to fix that too.
Well, that's a different bug, although it's in the same file; looks like it is in the stanza which begins if freeof(%n,aa) and freeof(%n,bb)
around line 1096. In the interest of clarity, I would suggest opening a separate bug report about it.
comment:16 in reply to: ↑ 15 Changed 3 months ago by
Replying to gh-robert-dodier:
Replying to slelievre:
but this still gives trouble:
sage: from sage.calculus.calculus import symbolic_product sage: x, n = SR.var('x, n') sage: symbolic_product(-x, x, 1, n) #0: simplify_product(product='product(-_SAGE_VAR_x,_SAGE_VAR_x,1,_SAGE_VAR_n)) ... Traceback (most recent call last) ... RuntimeError: ECL says: Error executing code in Maxima: factorial: factorial of negative integer -1 not defined.It would be nice to fix that too.
Well, that's a different bug, although it's in the same file; looks like it is in the stanza which begins
if freeof(%n,aa) and freeof(%n,bb)
around line 1096. In the interest of clarity, I would suggest opening a separate bug report about it.
You are right. This is now tracked at #31557.
comment:17 Changed 3 months ago by
Fixed (to the best of my knowledge; no doubt you'll want to verify for Sage) by commit 762fd85 on maxima-code/master which applies the patch for simplify_product shown above.
comment:18 Changed 3 months ago by
Thanks!
comment:19 follow-up: ↓ 20 Changed 2 months ago by
- Description modified (diff)
Should we apply the upstream commit fixing this as a patch or wait for the next Maxima release?
comment:20 in reply to: ↑ 19 Changed 2 months ago by
Replying to slelievre:
Should we apply the upstream commit fixing this as a patch or wait for the next Maxima release?
For what it's worth, a new Maxima version is on the horizon; it is being discussed on the mailing list, although no work on it has been announced. So my guess is that a new version will appear in, let's say, one to three months. I don't know whether that implies you should wait or not.
comment:21 Changed 6 weeks ago by
- Milestone changed from sage-9.3 to sage-9.4
Moving to 9.4, as 9.3 has been released.
FWIW, Sympy doesn't seem to be able to pull this one:
Mathematica isn't more helpful:
Giac's answer can't be believed :
HTH, but doubting it...