Changeset 7485:c3c8fc4a0411

Timestamp:

11/27/07 06:43:17 (5 years ago)

Author:

William Stein <wstein@…>

Branch:

default

Message:

trac #1290 -- implement rencontres numbers in sage. (really by Dan Drake)

File:

• sage/combinat/combinat.py (modified) (2 diffs)

sage/combinat/combinat.py

@@ -193 +193 @@
 
 from sage.interfaces.all import gap, maxima
-from sage.rings.all import QQ, RR, ZZ
+from sage.rings.all import QQ, RR, ZZ, RealField
 from sage.rings.arith import binomial
 from sage.misc.sage_eval import sage_eval
@@ -2536 +2536 @@
     """
     return sage_eval(maxima.eval("bernpoly(x,%s)"%n), {'x':x})
+
+
+def rencontres_number(n, k):
+    r"""
+    Returns the Rencontres number D(n,k), the number of permutations of
+    {1, 2,..., n} with k fixed points.
+
+    EXAMPLES:
+    Because 312 and 231 are the two permutations of {1, 2, 3} with 0
+    fixed points, we have:
+    
+        sage: rencontres_number(3,0)
+        2
+
+    More examples:
+        sage: rencontres_number(6,1)
+        264
+          
+
+    REFERENCES:
+        http://en.wikipedia.org/wiki/Rencontres_number
+        
+        Sequence A008290 in the OEIS.
+
+    AUTHORS:
+        -- Dan Drake, 2007-11-27: post to sage-devel
+        -- William Stein, 2007-11-27: referee and clean up for inclusion in sage
+    """
+    F = factorial(n - k)
+    # Sufficient precision needed to guarantee the correct result (this is overkill somewhat)
+    R = RealField(10*F.bits() + 53)
+    e = R(1).exp()
+    if (n - k) % 2 == 0:
+        return binomial(n, k) * (F/e).ceil()
+    else:
+        return binomial(n, k) * (F/e).floor()
+