Définition

Un entier positif nest un nombre pratique (séquence OEIS A005153 ) si tous les entiers positifs plus petits peuvent être représentés comme des sommes de diviseurs distincts de n.

Par exemple, 18est un nombre pratique: ses diviseurs sont 1, 2, 3, 6, 9 et 18, et les autres entiers positifs inférieurs à 18 peuvent être formés comme suit:

 4 = 1 + 3          5 = 2 + 3           7 = 1 + 6
 8 = 2 + 6          10 = 1 + 9         11 = 2 + 9
12 = 3 + 9 = 1 + 2 + 9 = 1 + 2 + 3 + 6
13 = 1 + 3 + 9      14 = 2 + 3 + 9      15 = 6 + 9
16 = 1 + 6 + 9      17 = 2 + 6 + 9

Mais ce 14n'est pas un nombre pratique: ses diviseurs sont 1, 2, 7 et 14, et il n'y a aucun sous-ensemble de ceux-ci qui s'ajoute à 4, 5, 6, 11, 12 ou 13.

Défi

Écrivez un programme, une fonction ou un verbe qui prend en entrée un entier positif xet retourne ou imprime le x ^ème nombre pratique, indexé de 1 pour cohérence avec OEIS. Votre code doit être suffisamment efficace pour pouvoir gérer jusqu'à 250000 entrées en moins de deux minutes sur un ordinateur de bureau raisonnable. (NB mon implémentation de référence en Java gère 250000 en moins de 0,5 seconde, et mon implémentation de référence en Python le gère en 12 secondes).

Cas de test

Input        Expected output
1            1
8            18
1000         6500
250000       2764000
1000000      12214770
3000000      39258256

J (99 caractères)

f=:3 :0
'n c'=.0 1
while.c<y do.
'p e'=.__ q:n=.n+2
c=.c+*/(}.p)<:}:1+*/\(<:p^e+1)%<:p
end.
n+n=0
)

Étant donné que l'énoncé du problème demande un "programme, fonction ou verbe ", quelqu'un a dû faire une soumission J. J les gens remarqueront que je n'ai pas vraiment joué au golf (!) Ou optimisé cela. Comme les autres entrées, j'ai utilisé le théorème de Stewart, mentionné sur le lien OEIS, pour tester si chaque nombre pair est pratique ou non.

Je n'ai pas facilement accès à un "ordinateur de bureau raisonnable" avec J installé. Sur mon netbook de six ans, il f 250000calcule en 120,6 secondes, ce qui n'est pas tout à fait en moins de deux minutes, mais probablement sur n'importe quel ordinateur un peu plus raisonnable, cela se termine à temps.

Mathematica, 126 121 caractères

Merci à Bélisaire.

Utilisation de la formule sur wikipedia.

f=(i=j=1;While[j<#,If[And@@Thread[#[[;;,1]]<2+Most@DivisorSum[FoldList[#Power@@#2&,1,#],#&]&@FactorInteger@++i],j++]];i)&

Exemples:

f[1]

1

f[8]

18

f[250000]

2764000

Il a fallu des années 70 pour calculer f[250000]sur mon ordinateur.

Haskell - 329

s 1=[]
s n=p:(s$div n p)where d=dropWhile((/=0).mod n)[2..ceiling$sqrt$fromIntegral n];p=if null d then n else head d
u=foldr(\v l@((n,c):q)->if v==n then(n,c+1):q else(v,1):l)[(0,1)]
i z=(z<2)||(head w==2)&&(and$zipWith(\(n,_)p->n-1<=p)(tail n)$scanl1(*)$map(\(n,c)->(n*n^c-1)`div`(n-1))n)where w=s z;n=u w
f=((filter i[0..])!!)

Exemples:

> f 1
1
> f 13
32
> f 1000
6500

Voici une petite suite de tests (avant celle ci-dessus):

import Data.Time.Clock
import System.IO

test x = do
    start <- getCurrentTime
    putStr $ (show x) ++ " -> " ++ (show $ f x)
    finish <- getCurrentTime
    putStrLn $ " [" ++ (show $ diffUTCTime finish start) ++ "]"

main = do
    hSetBuffering stdout NoBuffering
    mapM_ test [1, 8, 1000, 250000, 1000000, 3000000]

Résultats des tests après avoir été compilés avec ghc -O3:

1 -> 1 [0.000071s]
8 -> 18 [0.000047s]
1000 -> 6500 [0.010045s]
250000 -> 2764000 [29.084049s]
1000000 -> 12214770 [201.374324s]
3000000 -> 39258256 [986.885397s]

Javascript, 306 307 282B

function y(r){for(n=r-1,k=1;n;k++)if(p=[],e=[],c=0,P=s=1,!((x=k)%2|1==x)){while(x>1){for(f=x,j=2;j<=Math.sqrt(f);j++)if(f%j==0){f=j;break}f!=p[c-1]?(p.push(f),e.push(2),c++):e[c-1]++,x/=f}for(i=0;c>i;i++){if(p[i]>P+1){s=0;break}P*=(Math.pow(p[i],e[i])-1)/(p[i]-1)}s&&n--}return k-1}

250k env. 6s sur mon ordinateur portable.

Code non golfé commenté: http://jsfiddle.net/82xb9/3/ maintenant avec un meilleur test sigma et une meilleure condition si (merci commentaires)

Versions de pré-édition: http://jsfiddle.net/82xb9/ http://jsfiddle.net/82xb9/1/