Opened 8 months ago

Last modified 3 days 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:

Status badges

Description (last modified by slelievre)

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 8 months ago by charpent

FWIW, Sympy doesn't seem to be able to pull this one:

sage: from sympy import Product                                                 
sage: from sympy.abc import j, p, q                                             
sage: foo=Product(1-q^j, [j,1,p]); foo                                          
Product(1 - q**j, (j, 1, p))
sage: type(foo)                                                                 
<class 'sympy.concrete.products.Product'>
sage: foo.doit()                                                                
Product(1 - q**j, (j, 1, p))

Mathematica isn't more helpful:

sage: bar=mathematica("Product[1-q^j, {j, 1, p}]"); bar                         
QPochhammer[q, q, p]

Giac's answer can't be believed :

age: var("j,q,p")                                                              
(j, q, p)
sage: gargs=list(map(giac, [1-q^j,j,1,p,1]))                                    
sage: gargs                                                                     
[-q^j+1, j, 1, p, 1]
sage: libgiac.product(*gargs)                                                   
-q+1

HTH, but doubting it...

comment:2 Changed 8 months ago by heluani

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.

EDIT: it may also be in the final conversion max_to_sr

Last edited 8 months ago by heluani (previous) (diff)

comment:3 Changed 8 months ago by chapoton

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 8 months ago by chapoton

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 7 months ago by mkoeppe

  • Milestone changed from sage-9.2 to sage-9.3

comment:6 Changed 7 weeks ago by slelievre

  • 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: Changed 7 weeks ago by gh-macrakis

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

Last edited 7 weeks ago by gh-macrakis (previous) (diff)

comment:8 Changed 7 weeks ago by slelievre

  • 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 7 weeks ago by kcrisman

  • What does max_simplify_product do?

Here:

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 7 weeks ago by gh-robert-dodier

"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 7 weeks ago by slelievre

  • Description modified (diff)

Thanks for the rephrasing suggestion and your thoughts on the problem.

comment:12 Changed 7 weeks ago by gh-DaveWitteMorris

@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: Changed 7 weeks ago by slelievre

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 7 weeks ago by gh-robert-dodier

  • 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: Changed 7 weeks ago by 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.

comment:16 in reply to: ↑ 15 Changed 7 weeks ago by slelievre

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 7 weeks ago by gh-robert-dodier

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 7 weeks ago by gh-DaveWitteMorris

Thanks!

comment:19 follow-up: Changed 3 weeks ago by slelievre

  • 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 3 weeks ago by gh-robert-dodier

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 3 days ago by mkoeppe

  • Milestone changed from sage-9.3 to sage-9.4

Moving to 9.4, as 9.3 has been released.

Note: See TracTickets for help on using tickets.