Factoriel et log factoriel

Factoriel et log factoriel

Factoriel:

Quand vous avez besoin de calculer n!, il existe plusieurs solutions: 

  1. La solution « rapide »: qui utilise une boucle ou une fonction récursive: 
    def factorial_for(n):
      r = 1
      for i in range(2, n + 1):
        r *= i
      return(r)
    
    def factorial_rec(n):
      if n > 1:
        return(n * factorial_rec(n - 1))
      else:
        return(1)

    Ici, la multiplication séquentielle de chaque nombre va générer un nombre très grand très rapidement, ce qui n’est pas optimum. En effet, les ordinateurs sont plus rapides quand ils multiplient des petits nombres (120×30240) par rapport à la même multiplication avec un grand nombre (362880×10).

  2. La solution « intelligente »: Ici, l’idée est de séparer une séquence de nombres à multiplier en deux par le milieu, de façon récursive.
    10! = 1x2x3x4x5x6x7x8x9 x 10 = 362880 x 10 (rush solution)
    = 1x2x3x4x5 x 6x7x8x9x10 = 120 x 30240 (1st split)
    = 1x2x3 x 4×5 x 6x7x8 x 9×10 = 6 x 20 x 336 x 90 (2nd split)
    def range_prod(lo, hi):
      if lo + 1 < hi:
        mid = (hi + lo) // 2
        return(range_prod(lo, mid) * range_prod(mid + 1, hi))
      if lo == hi:
        return(lo)
      return(lo * hi)
    
    def factorial_tree(n):
      if n < 2:
        return(1)
      return(range_prod(1, n))
  3. La solution « paresseuse » ou « très intelligente » (tout dépend du point de vue): consiste à calculer une estimation du factoriel avec l’approximation de Stirling. 
    import math
    def stirling(n):
      return(math.sqrt(2 * math.pi * n) * (n / math.e) ** n)

    Cette implémentation est limitée à n<170. Toutefois, comme vous pouvez le constater plus bas, avec n=170, le résultat est énorme.
    Avec ce genre de gros nombre, nous préférons de toute façon travailler avec les log factoriels.

    In [01]: stirling(170)
    Out[01]: 7.253858934542908e+306
    In [02]: stirling(171)
    Out[02]: inf
    

     

Log Factorial:

Pour les deux premières solutions, vous n’avez qu’à ajouter math.log() (en faisant import math d’abord) :

def log_fact_rec(n):
  return(math.log(factorial_rec(n)))
def log_fact_tree(n):
  return(mat.log(factorial_tree(n)))
def log_fact_for(n):
  return(math.log(factorial_for(n)))

Le log factoriel avec l’approximation de Stirling est implémentée comme ceci :

def log_stirling(n):
 fn = float(n)
 fn = n * math.log(fn) + math.log(2 * math.pi * fn) / 2 - fn + \
      (fn ** -1) / 12 - (fn ** -3) / 360 + (fn ** -5) / 1260 - \
      (fn ** -7) / 1680 + (fn ** -9) / 1188
 return(fn)

Note:

Dans le cas où vous avez de multiples log factoriels à calculer, il peut être utile de sauvegarder et réutiliser des résultats précédents et ainsi gagner du temps (par exemple, 10! = 6!x7x8x9x10) comme ceci :

def log_fact_exact_saved(n, value_saved):
  select = max((i for i in value_saved.keys() if i <= n))
  if select != n:
    value_saved[n] = range_prod(select + 1, n) * value_saved[select]
  return(math.log(value_saved[n]))

def log_fact_approx_saved(n, value_saved):
  if n not in value_saved:
    value_saved[n] = log_stirling(n)
  return(value_saved[n])

Test:

Nous pouvons maintenant tester et comparer chacune de ces solutions en utilisant les librairies math et scipy pour voir laquelle est la plus rapide :


import sys
sys.setrecursionlimit(1000000) // increase the depth 
import math
import scipy.special
def test_saved(saved, approx):
    value = [10, 50, 100, 1000, 1000, 2000, 5000, 5000, 5000,
             5000, 10000, 11000, 11000, 11000, 20000, 20000]
    if approx:
        if saved:
            value_saved = {1: 1}
            for num in value:
                log_fact_approx_saved(num, value_saved)
        else:
            for num in value:
                log_stirling(num)
    else:
        if saved:
            value_saved = {1: 1}
            for num in value:
                log_fact_exact_saved(num, value_saved)
        else:
            for num in value:
                log_fact_tree(num)

def log_fact_math(n):
    return(math.log(math.factorial(n)))

def log_fact_scipy(n, exact=False):
    return(math.log(scipy.special.factorial(n, exact)))

import timeit
num = 10000
rep = 200
print("rec:" + str(timeit.timeit("log_fact_rec(" + str(num) + ")",
      setup="from __main__ import log_fact_rec", number=rep)))
print("for:" + str(timeit.timeit("log_fact_for(" + str(num) + ")",
      setup="from __main__ import log_fact_for", number=rep)))
print("scipy exact:" + str(timeit.timeit("log_fact_scipy(" + str(num) + ", True)",
      setup="from __main__ import log_fact_scipy", number=rep)))
print("math:" + str(timeit.timeit("log_fact_math(" + str(num) + ")",
      setup="from __main__ import log_fact_math", number=rep)))
print("tree:" + str(timeit.timeit("log_fact_tree(" + str(num) + ")",
      setup="from __main__ import log_fact_tree", number=rep)))

print("scipy approx:" + str(timeit.timeit("log_fact_scipy(" + str(num) + ")",
      setup="from __main__ import log_fact_scipy", number=rep)))
print("log_stirling:" + str(timeit.timeit("log_stirling(" + str(num) + ")",
      setup="from __main__ import log_stirling", number=rep)))

print("value log(100!) tree:" + str(log_fact_tree(100)))
print("value log(100!) scipy approx:" + str(log_fact_scipy(100)))
print("value log(100!) log_stirling:" + str(log_stirling(100)))

print("severals factorial:" + str(timeit.timeit("test_saved(False, False)",
      setup="from __main__ import test_saved", number=rep)))
print("severals factorial saved:" + str(timeit.timeit("test_saved(True, False)",
      setup="from __main__ import test_saved", number=rep)))
print("severals factorial approx:" + str(timeit.timeit("test_saved(False, True)",
      setup="from __main__ import test_saved", number=rep)))
print("severals factorial approx saved:" + str(timeit.timeit("test_saved(True, True)",
      setup="from __main__ import test_saved", number=rep)))

En utilisant python 2.7.6 et scipy 0.17.0, nous avons :

rec:4.52898597717
for:4.09089803696
scipy exact:4.08975815773
math:4.02955198288
tree:1.2138531208
Nous pouvons voir que la solution « intelligente » est la plus rapide,

même lorsque comparée à scipy

scipy approx:0.00166201591492
log_stirling:0.000301837921143
 Encore ici, scipy est le plus lent
value log(100!) tree:363.739375556
value log(100!) scipy approx:363.739375556
value log(100!) log_stirling:363.739375556
Dans ce cas-ci, toutes les méthodes sont équivalentes pour n>=100
several factorials:14.3231370449
several factorials saved:3.51164102554
several factorials approx:0.00481510162354
several factorials approx saved:0.00360107421875
Sauvegarder des résultats et les réutiliser peut sauver beaucoup de temps,

mais prendra plus de mémoire.

J’espère que cet article vous aidera à choisir la meilleure solution (adaptée à vos besoins!) pour calculer des nombres factoriels et log factoriels.

By | 2017-04-29T17:06:29+00:00 19 février 2016|Categories: Performance, Python|Tags: |0 Commentaires

About the Author:

Informaticien de formation, j’ai vite compris que la bioinformatique regorge d’égnimes à résoudre. Comme dans le « Sommet des dieux »(Jiro Taniguchi), il y a toujours un nouveau sommet à gravir ou un itinéraire plus direct à tenter.

Laisser un commentaire