Opened 7 years ago
Closed 7 years ago
#14496 closed enhancement (fixed)
unify the three implementations of gaussian qbinomial coefficients
Reported by:  chapoton  Owned by:  tbd 

Priority:  major  Milestone:  sage5.10 
Component:  combinatorics  Keywords:  gaussian binomial 
Cc:  fwclarke, andrew.mathas  Merged in:  sage5.10.beta2 
Authors:  Frédéric Chapoton  Reviewers:  Francis Clarke, Travis Scrimshaw 
Report Upstream:  N/A  Work issues:  
Branch:  Commit:  
Dependencies:  Stopgaps: 
Description
In sage 5.8, one can find the gaussian qbinomial coefficients in (at least) three places :
 sage.combinat.sf.kfpoly.q_bin
 sage.combinat.q_analogues.q_binomial
 sage.rings.arith.gaussian_binomial
The syntax is not quite the same for q_bin as for the two others.
Some timings:
sage: timeit("sage.combinat.sf.kfpoly.q_bin(7,5)") 625 loops, best of 3: 819 µs per loop sage: timeit("sage.combinat.q_analogues.q_binomial(12,5)") 625 loops, best of 3: 1.12 ms per loop sage: timeit("sage.rings.arith.gaussian_binomial(12,5)") 625 loops, best of 3: 294 µs per loop
The parents :
sage: sage.combinat.sf.kfpoly.q_bin(7,5).parent() Fraction Field of Univariate Polynomial Ring in t over Integer Ring sage: sage.combinat.q_analogues.q_binomial(12,5).parent() Univariate Polynomial Ring in q over Integer Ring sage: sage.rings.arith.gaussian_binomial(12,5).parent() Fraction Field of Univariate Polynomial Ring in q over Integer Ring
The behaviour with an added integer parameter
sage: sage.combinat.sf.kfpoly.q_bin(7,5,2) 114429029715 sage: sage.combinat.q_analogues.q_binomial(12,5,2) DOES NOT WORK sage: sage.rings.arith.gaussian_binomial(12,5,2) 114429029715
The behaviour with a polynomial parameter
sage: w=polygen(ZZ,'w') sage: sage.combinat.sf.kfpoly.q_bin(4,2,w) w^8 + w^7 + 2*w^6 + 2*w^5 + 3*w^4 + 2*w^3 + 2*w^2 + w + 1 sage: sage.combinat.q_analogues.q_binomial(6,2,w) w^8 + w^7 + 2*w^6 + 2*w^5 + 3*w^4 + 2*w^3 + 2*w^2 + w + 1 sage: sage.rings.arith.gaussian_binomial(6,2,w) w^8 + w^7 + 2*w^6 + 2*w^5 + 3*w^4 + 2*w^3 + 2*w^2 + w + 1 sage: sage.combinat.sf.kfpoly.q_bin(4,2,w).parent() Univariate Polynomial Ring in w over Integer Ring sage: sage.combinat.q_analogues.q_binomial(6,2,w).parent() Univariate Polynomial Ring in w over Integer Ring sage: sage.rings.arith.gaussian_binomial(6,2,w).parent() Fraction Field of Univariate Polynomial Ring in w over Integer Ring
Maybe it would be better to have a single function ?
Attachments (3)
Change History (35)
comment:1 Changed 7 years ago by
comment:2 Changed 7 years ago by
 Cc fwclarke added
 Status changed from new to needs_review
comment:3 Changed 7 years ago by
 Status changed from needs_review to needs_work
This looks good. It's certainly nice to have integer parameters giving integer results.
But n//d
can causes problems if q
takes a real value, and similarly in some other cases. Perhaps this should be a "try:
", with the exception giving n/d
or even gaussian_binomial(n, k)(q)
. [It is surely not necessary to introduce var
when this occurs before.]
I also don't like the way (inherited from the existing code) that n
plays two very different rôles; the line
n = prod([1  q**i for i in range(nk+1,n+1)])
is particularly bad style. It would make the code more readable to call them "numerator
" and "denominator
" (possibly abbreviated)
comment:4 Changed 7 years ago by
 Status changed from needs_work to needs_review
Here is a new version of the patch. I have modified the names d and n and introduced a try statement. I am not sure what kind of exceptions I need to catch.
Still it remains to choose if one keeps the gaussian_binomial or rather the q_binomial function.
comment:5 Changed 7 years ago by
 Component changed from PLEASE CHANGE to combinatorics
comment:6 Changed 7 years ago by
 Status changed from needs_review to needs_work
 There's a doctest failure at line 341 of
sage/combinat/sf/kfpoly.py
. I don't know what's being computed here, so I don't know what's going wrong.
 I thought the standard spacing was as
def gaussian_binomial(n, k, q=None):
 I think
TypeError
is enough; I never saw anything else when I was testing your code.
 Is it intentional not to alter
sage/combinat/q_analogues.py
?
So it's nearly ready for a positive review.
comment:7 Changed 7 years ago by
 Status changed from needs_work to needs_info
Here is a new patch.
Point 1 is indeed a problem. I do not understand either what is computed in the weight function. The procedure q_bin was returning 1 for any negative values of n and any k. This is not a very natural behaviour.
Point 4 : I have chosen that the name to be kept would be q_binomial.
There is another problem remaining. The existing algorithm in sage.combinat.q_analogues uses cyclotomic polynomials and become faster than the naive one when the arguments are big (and of similar size ?), but is much slower in some cases (when n is small or k is small ?).
sage: timeit("sage.combinat.q_analogues.q_binomial(85,3)") 625 loops, best of 3: 246 µs per loop sage: timeit("sage.combinat.q_analogues.q_binomial_via_cyclotomic(85,3)") 125 loops, best of 3: 2.78 ms per loop sage: timeit("sage.combinat.q_analogues.q_binomial(85,43)") 25 loops, best of 3: 9.21 ms per loop sage: timeit("sage.combinat.q_analogues.q_binomial_via_cyclotomic(85,43)") 25 loops, best of 3: 8.18 ms per loop
So one needs to keep both and find a strategy of choice..
comment:8 Changed 7 years ago by
I'll take a look at why the weight function in SF fails; but the best way to fix it is in the weight, have a check for negative n
or k
and just circumvent the call to q_binomial()
and return 1.
As for the speed of the naive vs. cyclotomic polys, I would pick a cutoff value which you think is about the breakeven point(s) for the two, and call which ever function is faster based on that.
Also, have you tried implementing it by computing the qbinomial, then substituting in a user defined q
value to the result (in particular, have you compared the speed)? I'm thinking this might be slower than just computing it outright if q
is in GF(3)
for example... If you haven't, I can write it up and run the tests.
Thanks for noticing this and working on it,
Travis
comment:9 Changed 7 years ago by
Good idea, I have done just that in the new patch, so there is no doctest failure now.
There remains to choose a strategy..
comment:10 followup: ↓ 11 Changed 7 years ago by
Some timings were given in the ticket #13166
According to my own timings, it seems that for n <= ~ 78 the naive algorithm is faster than the cyclotomic algorithm for computing the gaussian polynomial.
comment:11 in reply to: ↑ 10 Changed 7 years ago by
Replying to chapoton:
According to my own timings, it seems that for n <= ~ 78 the naive algorithm is faster than the cyclotomic algorithm for computing the gaussian polynomial.
I get much the same. For gaussian_binomial(2*k, k)
I found the naive method faster for k < 30
, and the cyclotomic method faster for `k > 30'.
One point about the naive method:
sage: %timeit gaussian_binomial(50, 4) 1000 loops, best of 3: 287 µs per loop sage: %timeit gaussian_binomial(50, 46) 100 loops, best of 3: 3.41 ms per loop
Since the code contains two loops of length k
, it's obviously better to evaluate gaussian_binomial(n, k)
as gaussian_binomial(n, nk)
if k > n/2
.
comment:12 Changed 7 years ago by
Here is a new patch
 using the minimum of k and nk in the naive algo
 implementing the following strategy :
if n <= 70 or k <= N/4 or q is not a polynomial, use the naive algo
otherwise use the cyclotomic algo
I have found this strategy using some timings (I only timed the polynomial case).
comment:13 Changed 7 years ago by
 Reviewers set to Francis Clarke, Travis Scrimshaw
 Status changed from needs_info to needs_review
comment:14 Changed 7 years ago by
 Type changed from task to enhancement
comment:15 Changed 7 years ago by
 Status changed from needs_review to needs_work
Two points remain, I think:
 On speed, I see that you are not using the cyclotomic method when
q
is not a polynomial. But by usingcyclotomic_value
the faster method can be efficient for "scalar" values too.
 I don't see any reason why
gaussian_binomial
should be deprecated. This terminology is in common use, and a mathematician who knows only this name should not be told that it is "wrong". She should be able to use it without necessarily ever knowing that some people use a different name. Doesn't deprecation mean that eventually the function will disappear from Sage? This would only cause confusion. Surely the two should be synonyms, as are, for example, the numberfield methodsring_of_integers
andmaximal_order
.
Changed 7 years ago by
comment:16 Changed 7 years ago by
Here is a new patch.
For point 2, I have now defined gaussian_binomial as being an alias to q_binomial. The global namespace only contains gaussian_binomial, as before the patch.
For point 1, I have made some changes in the choice of algo, mainly about symbolic inputs. But I do not think that it is really worth spending more time on trying harder to enhance the speed for arbitrary inputs. This function is supposed to be used mainly without 3rd argument.
comment:17 Changed 7 years ago by
 Status changed from needs_work to needs_review
comment:18 Changed 7 years ago by
Hey Frederic and Francis,
I've uploaded a review patch which does two things I wanted, an input for the algorithm and handling objects which may not have division (by computing it using a formal variable and then substituting in q
). If you're happy with my changes, feel free to set this to a positive review.
Best,
Travis
comment:19 Changed 7 years ago by
hello,
it seems to me that at the end you can write simply
return q_binomial(n, k)(q)
because then the polygen import will be done inside the function
Changed 7 years ago by
comment:20 Changed 7 years ago by
Yep, that works. I've uploaded the updated patch.
comment:21 Changed 7 years ago by
 Status changed from needs_review to needs_work
I was about to say positive review, when I noticed two things:

gaussian_polynomial(5, 2, 7r)
fails. But this caused by a feature ofcyclotomic_value
, so should probably be ignored.

gaussian_binomial(4, 2, Zmod(6)(2), algorithm='naive')
fails. This can be very simply rectified by excepting aZeroDivisionError
as well as aTypeError
just before the "last attempt".
comment:22 followup: ↓ 23 Changed 7 years ago by
 Status changed from needs_work to needs_review
Point 2 is done. I vote to ignore point 1 (anyway, who cares about that ?)
comment:23 in reply to: ↑ 22 Changed 7 years ago by
 Status changed from needs_review to positive_review
Replying to chapoton:
Point 2 is done. I vote to ignore point 1 (anyway, who cares about that ?)
Agreed. Positive review.
comment:24 Changed 7 years ago by
 Status changed from positive_review to needs_work
dochtml.log:[combinat ] /mazur/release/merger/sage5.10.beta2/local/lib/python2.7/sitepackages/sage/combinat/q_analogues.py:docstring of sage.combinat.q_analogues.q_binomial:141: WARNING: Duplicate explicit target name: "ch2006". dochtml.log:[combinat ] /mazur/release/merger/sage5.10.beta2/local/lib/python2.7/sitepackages/sage/combinat/q_analogues.py:docstring of sage.combinat.q_analogues.q_binomial:141: WARNING: duplicate citation CH2006, other instance in /mazur/release/merger/sage5.10.beta2/devel/sage/doc/en/reference/combinat/sage/combinat/q_analogues.rst
comment:25 Changed 7 years ago by
ah, well, sorry.
This is because I defined an alias in an incorrect way, using
gaussian_binomial = q_binomial
But what is the correct way ?
comment:26 Changed 7 years ago by
Hum, it seems that I have used the correct way to create an alias. But as the doc is duplicated (instead of saying that one function is an alias for the other), it will necessarily create a conflict any time it contains a reference. I am surprised that this problem has never been met before !
One can compare the situation to the definition of charpoly and characteristic_polynomial in sage.misc.functional, which does not contain a reference, hence dos not raise this problem.
comment:27 Changed 7 years ago by
There are two workarounds I can see:
 Put the reference in a place that isn't duplicated such as at the module level.
 Write in
gaussian_binomial()
which callsq_binomial()
and it's doc just says "seeq_binomial()
"
Option 2 is more work and less clean IMO. Option 1 isn't perfect to me either but is the lesser of 2 evils. Up to you; either way would be okay with me.
comment:28 Changed 7 years ago by
 Cc andrew.mathas added
comment:29 Changed 7 years ago by
 Status changed from needs_work to positive_review
ok, this should be good now. Back to positive review
comment:30 Changed 7 years ago by
Hey Frederic,
Two minor things on the doc for gaussian_binomial()
:
 There should be periods '.' at the end of the sentences.
 The 'this function' is slightly ambiguous to me, I would explicitly put
:func:`q_binomial()`
to avoid this.
Otherwise I'm also happy with the patch. Thanks.
Changed 7 years ago by
comment:31 Changed 7 years ago by
Done, and thanks for the review
comment:32 Changed 7 years ago by
 Merged in set to sage5.10.beta2
 Resolution set to fixed
 Status changed from positive_review to closed
Here is already a seemingly faster function. There remains to deprecate the two others. The most annoying one is q_bin, with a different syntax and with a global variable inside.