HackerRank : piège à code
Je me suis fait avoir une fois de plus. Après avoir résolu une septantaine* de problèmes de Project Euler de plus en plus mathématiquement ardus, suis tombé sur HackerRank.com , qui m’avait l’air plus orienté programmation, avec des petits problèmes solubles en quelques minutes, voire une petite soirée. Du moins en apparence…
Par exemple, j’ai attaqué le problème “nCr ” très confiant : étant donnés deux entiers n et r, combien y’a-t-il de manières de choisir r objets parmi n ?
Fâââcile : c’est le coefficient binomial :
$${n \choose k} = C_n^k\\, = \frac{n!}{k!(n-k)!}$$\[mathjax\]
Avec l’habitude on sait qu’une formule avec des factorielles peut faire exploser n’importe quel ordinateur si on n’y prend pas garde donc on se méfie, on trouve sur RosettaCode ou directement sur la Wikipédia un bout de code Python de ce genre :
\[python\]def binomial_coefficient(n,k): if k < 0 or k > n: return 0 if k == 0 or k == n: return 1 k = min(k, n - k) # take advantage of symmetry c = 1 for i in range(k): c = c * (n - i) // (i + 1) # // dénote la division entière en Python return c
\[/python\]On copie-colle ça directement dans l’éditeur de code intégré à HackerRank, on ajoute deux lignes pour lire les données du “test case” et “printer” le résultat que HackerRank va comparer aux bonnes réponses. Yapluka cliquer le bouton “Run Code”… Petite bulle : j’ai oublié le %142857 pour fournir les résultats modulo 142857 comme demandé. Correction, re-click sur “Run Code” : ça passe. Après il faut cliquer sur “Submit Code” pour exécuter le programme avec d’autres “test cases” que ceux fournis dans la donnée…
Et là, ça ne va plus bien du tout, parce que la donnée précise que “1 <= n <= 1000000000 , 0 <= r <= n “, ce qui signifie en clair que le programme doit pouvoir dans le pire cas calculer en un temps raisonnable le produit des entiers de 500'000'001 à 1'000'000'000, sans faire exploser la mémoire.
Avec le code ci-dessus, binomial_coefficient(10000,5000) ne prend que 0.029 secondes sur mon PC, mais 4.16 secondes pour binomial_coefficient(10000,5000). Ca peut surprendre car la fonction semble O(k)*** , mais les multiplications et divisons sur des nombres de millions de chiffres, ça prend aussi pas mal de temps. Donc binomial_coefficient(1000000000,500000000), on oublie.
Une tentation serait d’utiliser une approximation basée sur la formule de Stirling
$$\log {n\choose k} \approx (n+\frac{1}{2})\log n - (k+\frac{1}{2})\log k - (n-k+\frac{1}{2})\log (n-k) - \frac{1}{2}\log 2\pi $$extrêmement rapide à calculer: log_binomial_coefficient(1000000000,500000000) ne met que 10 microsecondes à donner le résultat, 693147169.9725189 . Donc
$${1000000000\choose 500000000} \approx e^{693147169.9725189}$$mais en raison de l’approximation, il n’y a aucune chance que les derniers chiffres de ce nombre soient corrects, or c’est ceux qui m’intéressent, en raison du modulo 142857…
Mais bien sur ! La ligne la plus importante de la donnée c’est “Output all answers modulo 142857.” ! J’ai été piégé par l’anodin “output all answers …”. En fait c’est tout le programme qu’il faut penser en arithmétique modulaire !
L’arithmétique modulaire
Pour ceux qui sont arrivés jusqu’ici en ne se rappelant plus trop ce que “modulo” veut dire, un nombre x “modulo” n est égal au reste de la division x/n et se note x mod n. Par exemple 45 mod 7 = 3 parce que 45 = 6x7 + 3
L’arithmétique modulaire, c’est le calcul avec des nombres “modulo quelque chose”. C’est très utile dans certains domaines et facile car on peut facilement définir certaines opérations
- (a + b) mod n= ((a mod n) + (b mod n)) mod n
- (a - b) mod n= ((a mod n) - (b mod n)) mod n
- (a x b) mod n= ((a mod n) x (b mod n)) mod n
mais pour la division, dont on a besoin pour résoudre le problème, c’est plus compliqué. Tellement que je traduis ci-dessous [1] , la page la plus claire que j’ai pu trouver sur le sujet:
La division en arithmétique modulaire
Il y a trois manières de calculer x= a / b (mod n):
- Essai et Erreur : Est-ce que 1.b = a ?, ou 2.b = a ?, ou 3.b = a ?, etc ….
- L’algorithme d’Euclide inversé : résoudre bx + ny = a pour trouver x.
- Multiplier par 1/b (mod n), si vous savez ce qu’est 1/b (mod n)
La méthode à utiliser dépend de l’information que vous avez.
1. Essai/Erreur:
Par exemple, si on vous donne une table de multiplication, alors la méthode par essai et erreur est clairement la plus simple. Toutes les valeurs 1.b, 2.b, 3.b, 4.b, 5.b, … sont déjà dans la colonne b de votre table, donc tout ce que vous avez à faire est de chercher dans cette colonne jusqu’à ce que vous trouviez a, et x est alors le numéro de la ligne correspondante.
La méthode marche aussi sans table de multiplication SI n est petit. Par exemple si n=5, alors il faut utiliser la méthode essai/erreur, car il n’y a que 4 nombres à vérifier.
Une autre version souvent plus rapide de l’essai/erreur pour la division consiste à examiner les valeurs successives de a (mod n) en ajoutant continuellement n à a jusqu’à trouver un a divisible par b, ou qui permette de réduire la fraction a/b. Par exemple pour 3/16 (mod 53) on a:
$$\frac{3}{16} = \frac{56 (=3+53)}{16} = \frac{7}{2} = \frac{60 (=7+53)}{2} = 30 (mod 53)$$2) Algorithme d’Euclide étendu:
Si n est plus grand, vous devriez considérer l’algorithme d’euclide étendu . Voici pourquoi: a/b = x (mod n) 11/14 = x (mod 31) a = bx (mod n) 11 = 14x (mod 31) a = bx + ny for some y. 11 = 14x + 31y for some y. This is one linear equation in 2 variables, so it is precisely the sort of equation we learned how to solve using the Backwards Euclidean Algorithm. Doing this gives us a solution for x and y, and all we care about is x. Let’s try this out with the above example with 11/14 (mod 31).
To solve 11 = 14x + 31y, we first do the Euc. Alg for 14, 31.
a) 31 = 2 x 14 + 3 —> 3 = 31 - 2 x 14 – – - - – –
b) 14 = 4 x 3 + 2 —> 2 = 14 - 4 x 3 – - - - – -
c) 3 = 1 x 2 + 1 —> 1 = 3 - 1 x 2 - - - - - -
so gcd(31,14) = 1.
Backwards:
Start with (c). 1 = 3 - 1 x 2 - - -
Sub (b) for 2. = 3 - 1 x (14 - 4 x 3) - - – -
Collect terms. = -1 x 14 + 5 x 3 – -
Sub (a) for 3. = -1 x 14 + 5 x (31 - 2 x 14) - – – –
Collect terms. 1 = 5 x 31 + (-11) x 14. - – –
Recall that we wanted to solve 11 = 14x + 31y, so we multiply the last equation by 11 to obtain 11 = 55 x 31 + (-121) x 14. Thus we get x = -121.
Finally, we must represent x mod 31. -121 = -4 x 31 + 3, so x = 121 = 3 (mod 31).
Let us check this. We claimed that 11/14 = 3 (mod 31). Then we must have 11 = 14 x 3 (mod 31). This works correctly, because 14 x 3 = 42 = 31 + 11.
- Multiply by 1/b (mod n). ————————— It is rare that you will already be told what 1/b is, but this is precisely what you are told in Homework 16, problem 3, so make the most of it.
Recall that a/b = a x (1/b), so, for example, if we know that 1/5 = 3 (mod 7), then 2/5 = 2 x 1/5 = 2 x 3 = 6 (mod 7). Check that this works: 5 x 6 = 30 = 4 x 7 + 2. Good.
This is sort of like multiplying by .5 instead of dividing by 2, except that in this case, multiplying is a billion times easier than dividing, because that’s just the way that modular arithmetic works.
When might you start out knowing 1/b if you weren’t already told the answer like in Homework 16, problem 3? Well, if you want to divide by b a lot in some fixed base, then you should really only do a hard division problem _once_. Figure out 1/b once, and then for all the other a/b examples you want to calculate, just multiply a x (1/b).
For example, in part (2) of this discussion, when we calculated that 11/14 = 3 (mod 31), we first obtained the equation
1 = 5 x 31 + (-11) x 14. - – –
Well, this is precisely the equation we need for finding the reciprocal of 14. This equation says that (-11) x 14 = 1 (mod 31), so 1/14 = -11 = 20 (mod 31).
Given that 1/14 = 20 (mod 31), it should be easy to calculate anything/14 (mod 31). For example,
11/14 = 11 x 20 = 220 = 7 x 31 + 3 = 3 (mod 31) (as shown before)
28/14 = 28 x 20 = 560 = 18 x 31 + 2 = 2 (mod 31).
Okay, that summarizes the 3 techniques for division. Deciding which one to use is up to you, but usually it’s best to use the method that will take the least amount of time.
Back to the code
Notes:
- * : septantaine , c’est tout de même plus logique que soixante-dizaine, non ?
- ** : devoir pour demain : comment j’ai fait pour le savoir en moins d’une minute ?
- ***: comme il n’y a qu’une boucle for, le temps de calcul devrait être proportionnel à k
Références
- Sarah Dean Rasmussen, “Division in Modular Arithmetic ”, 2004, sur math.harward.edu
- http://www.doc.ic.ac.uk/~mrh/330tutor/ch03.html
- https://fr.wikibooks.org/wiki/Approfondissements _de_lyc%C3%A9e/Arithm%C3%A9tique_modulaire
- http://math.stackexchange.com/questions/95491/n-choose-k-bmod-m-using-chinese-remainder-theorem
- A. Granville “Arithmetic properties of binomial coefficients. I. Binomial coefficients modulo prime powers .” 1997, 253-276. CMS Conf. Proc. 20, Amer. Math. Soc., Providence
