Moyen rapide de calculer n! mod m où m est premier?

Question

J'étais curieux de savoir si il y avait une bonne façon de le faire. Mon code actuel est quelque chose comme:

def factorialMod(n, modulus):
    ans=1
    for i in range(1,n+1):
        ans = ans * i % modulus    
    return ans % modulus

Mais il semble assez lent!

Je ne peux pas calculer n! et puis appliquer le premier module, car parfois n est si grand que n! n'est tout simplement pas possible de calculer explicitement.

J'ai aussi vu http://en.wikipedia.org/wiki/Stirling%27s_approximation et je me demande si cela peut être utilisé à tous ici, d'une certaine façon?

Ou, comment pourrais-je créer une récursif, memoized fonction en C++?

Answer 1

n peut être arbitrairement grand

Eh bien, n ne peut pas être arbitrairement grand - si n >= m, alors n! ≡ 0 (mod m) (parce qu' m est l'un des facteurs, par la définition de la factorielle).

---

En supposant n << m et vous avez besoin d'une exacte valeur, votre algorithme ne peut aller plus vite, à ma connaissance. Toutefois, si n > m/2, vous pouvez utiliser l'identité suivante (le théorème de Wilson - Merci @Daniel Fischer!)

pour limiter le nombre de multiplications à propos de m-n

(m-1)! ≡ -1 (mod m)
1 * 2 * 3 * ... * (n-1) * n * (n+1) * ... * (m-2) * (m-1) ≡ -1 (mod m)
n! * (n+1) * ... * (m-2) * (m-1) ≡ -1 (mod m)
n! ≡ -[(n+1) * ... * (m-2) * (m-1)]-1 (mod m)

Cela nous donne une façon simple de calculer n! (mod m) en m-n-1 multiplications, plus modulaire inverse:

def factorialMod(n, module):
ans=1
 si n <= module//2:
 #calculer la factorielle normalement (à droite argument de gamme() est exclusive)
 for i in range(1,n+1):
 ans = (ans * i) % module 
autre chose:
 #Fancypants méthode pour n grand
 for i in range(n+1,le module):
 ans = (ans * i) % module
 ans = modinv(sna, le module)
 ans = -1*ans + module
 retour sna % module

---

Nous pouvons reformuler l'équation ci-dessus, d'une autre façon, qui peut ou peut-effectuez pas un peu plus rapide. À l'aide de l'identité suivante:

nous pouvons reformuler l'équation

n! ≡ -[(n+1) * ... * (m-2) * (m-1)]-1 (mod m)
n! ≡ -[(n+1-m) * ... * (m-2-m) * (m-1-m)]-1 (mod m)
 (inverser l'ordre des termes)
n! ≡ -[(-1) * (-2) * ... * -(m-n-2) * -(m-n-1)]-1 (mod m)
n! ≡ -[(1) * (2) * ... * (m-n-2) * (m-n-1) * (-1)(m-n-1)]-1 (mod m)
n! ≡ [(m-n-1)!]-1 * (-1)(m-n) (mod m)

Cela peut être écrit en Python comme suit:

def factorialMod(n, module):
ans=1
 si n <= module//2:
 #calculer la factorielle normalement (à droite argument de gamme() est exclusive)
 for i in range(1,n+1):
 ans = (ans * i) % module 
autre chose:
 #Fancypants méthode pour n grand
 for i in range(1,module-n):
 ans = (ans * i) % module
 ans = modinv(sna, le module)

 #Puisque m est impair-le premier, (-1)^(m-n) = -1 si n est pair, +1 si n est impair
 si n % 2 == 0:
 ans = -1*ans + module
 retour sna % module

---

Si vous n'avez pas besoin d'une exacte valeur, la vie devient un peu plus facile - vous pouvez utiliser Stirling approximation pour calculer une valeur approchée de la valeur en O(log n) du temps (en utilisant l'exponentiation par la quadrature).

---

Enfin, je dois mentionner que si c'est le temps critique et vous êtes à l'aide de Python, essayez de passer à C++. Par expérience personnelle, vous devriez vous attendre environ un ordre de grandeur augmentation de la vitesse ou plus, simplement parce que c'est exactement le genre de CPU serré en boucle que natif du code compilé excelle à (aussi, pour quelque raison que ce soit, GMP, semble beaucoup plus fine que Python Bignum).

Answer 2

L'expansion de mon commentaire à la réponse:

Oui, il y a des moyens plus efficaces pour ce faire. Mais ils sont extrêmement pénible.

Donc, à moins que vous vraiment besoin de plus de performances, je vous suggère de ne pas essayer de les appliquer.

---

La clé est à noter que le module (qui est essentiellement une division) va être le goulot d'étranglement de l'opération. Heureusement, il y a des algorithmes rapides qui vous permettront d'effectuer le module sur le même nombre de nombreuses fois.

• La Division par l'Invariant de nombres Entiers à l'aide de Multiplication
• Montgomery Réduction

Ces méthodes sont rapides, car ils éliminent essentiellement le module.

---

Ces seules les méthodes devrait vous donner une accélération modérée. Pour être vraiment efficace, vous pouvez avoir besoin de dérouler la boucle afin de permettre une meilleure CIB:

Quelque chose comme ceci:

ans0 = 1
ans1 = 1
for i in range(1,(n+1) / 2):
    ans0 = ans0 * (2*i + 0) % modulus    
    ans1 = ans1 * (2*i + 1) % modulus    

return ans0 * ans1 % modulus

mais en tenant compte d'une étrange nbre d'itérations et de le combiner avec l'une des méthodes que j'ai lié ci-dessus.

Certains pourraient faire valoir que la boucle de dérouler devrait être laissé au compilateur. Je vais à contre-argumenter que les compilateurs ne sont actuellement pas assez intelligent pour dérouler cette boucle. Regarder de plus près et vous verrez pourquoi.

---

Notez que bien que ma réponse est indépendant de la langue, il est destiné principalement à C ou C++.

Answer 3

n! mod m peut être calculé en opérations O (n ^{1/2 + ε} ) à la place du naïf O (n). Cela nécessite l'utilisation de la multiplication polynomiale FFT, et ne vaut que pour les très grands n, par exemple n> 10 ⁴ .

Un aperçu de l'algorithme et quelques timings peuvent être vus ici: http://fredrikj.net/blog/2012/03/factorials-mod-n-and-wilsons-theorem/

Answer 4

Si nous voulons calculer M = a*(a+1) * ... * (b-1) * b (mod p), nous pouvons utiliser l'approche suivante, si nous supposons que nous pouvons ajouter, soustraire et multiplier rapidement (mod p), et d'obtenir un temps de fonctionnement de la complexité de l' O( sqrt(b-a) * polylog(b-a) ).

Pour des raisons de simplicité, supposons (b-a+1) = k^2, est un carré. Maintenant, nous pouvons diviser notre produit en k parties, c'est à dire M = [a*..*(a+k-1)] *...* [(b-k+1)*..*b]. Chacun des facteurs de ce produit est de la forme p(x)=x*..*(x+k-1), pour x.

À l'aide d'un algorithme de multiplication rapide de polynômes, comme Schönhage–Strassen algorithme, dans un divide & conquer manière, on peut trouver les coefficients du polynôme p(x) in O( k * polylog(k) ). Maintenant, apparemment, il existe un algorithme de substitution k points au même degré-k polynôme en O( k * polylog(k) ), ce qui signifie, nous pouvons calculer le p(a), p(a+k), ..., p(b-k+1) rapide.

Cet algorithme de substitution nombre de points en un seul polynôme est décrit dans le livre "nombres premiers", par C. Pomerance et R. Crandall. Finalement, lorsque vous avez ces k valeurs, vous pouvez les multiplier en O(k) et obtenir la valeur souhaitée.

Notez que l'ensemble de nos activités et prises (mod p). Le fonctionnement exact de temps est - O(sqrt(b-a) * log(b-a)^2 * log(log(b-a))).

Answer 5

En développant mon commentaire, cela prend environ 50% du temps pour tout n dans [100, 100007] où m = (117 | 1117):

 Function facmod(n As Integer, m As Integer) As Integer
    Dim f As Integer = 1
    For i As Integer = 2 To n
        f = f * i
        If f > m Then
            f = f Mod m
        End If
    Next
    Return f
End Function