F # parallélisation question lors du calcul de nombres parfaits?

Question

Je suis en train d'optimiser un petit programme qui calcule les nombres parfaits d'un exposant donné.

Les pistes de programme (presque) parfaitement, mais quand j'ouvre le gestionnaire de tâches, il fonctionne toujours sur un seul thread. Cela signifie que je dois faire quelque chose de mal, mais ma connaissance de F # est encore dans une phase « commencement ».

Je vais essayer de poser cette question aussi claire que je le peux, mais si je ne faisant, s'il vous plaît laissez-moi savoir.

Un nombre parfait est un nombre dont la somme de tous les diviseurs de (sauf pour le nombre lui-même) est égal au nombre lui-même (par exemple 6 est parfaite, puisque la somme des diviseurs 1, 2 et 3 de celui-ci sont 6).

J'utilise les nombres premiers pour accélérer le calcul, c'est que je ne suis pas intéressé par les listes (énormes) où tous les diviseurs sont stockés. Pour ce faire, j'utilise la formule qu'Euclide avéré être correct: (2 * (puissance de num - 1)) * (2 * (puissance de num - 1)) lorsque celle-ci est un Mersenne premier. J'ai utilisé un algorithme très rapide de stackoverflow (par @Juliet) pour déterminer si un nombre donné est premier.

Comme je l'ai lu dans plusieurs articles (je ne l'ai pas encore acheté un bon livre, si la honte sur moi) sur Internet, j'ai découvert que des séquences de meilleures performances que les listes. Voilà pourquoi j'ai commencé à créer une fonction qui génère une séquence de nombres parfaits:

   let perfectNumbersTwo (n : int) =  
    seq { for i in 1..n do 
           if (PowShift i) - 1I |> isPrime 
           then yield PowShift (i-1) * ((PowShift i)-1I)
        } 

Le PowShift de helperfunction est mis en œuvre comme suit:

    let inline PowShift (exp:int32) = 1I <<< exp ;;

J'utilise un opérateur de décalage de bits, depuis la base de tous les calculs de puissance est de 2, donc cela pourrait être un moyen facile. Bien sûr, je suis toujours reconnaissant pour les contributions sur la question que je posais à ce sujet sur: F # problèmes d'alimentation qui accepte les arguments à bigints> F # problèmes d'alimentation qui accepte les arguments à bigints

La fonction Juliet créé ( ici empruntèrent ) est comme suit:

let isPrime ( n : bigint) = 
    let maxFactor = bigint(sqrt(float n))
    let rec loop testPrime tog =
        if testPrime > maxFactor then true
        elif n % testPrime = 0I then false
        else loop (testPrime + tog) (6I - tog)
    if n = 2I || n = 3I || n = 5I then true
    elif n <= 1I || n % 2I = 0I || n % 3I = 0I || n % 5I = 0I then false
    else loop 7I 4I;;

En utilisant ce code, sans parallèle, il faut environ 9 minutes sur mon ordinateur portable pour trouver le 9 nombre parfait (qui se compose de 37 chiffres, et peut être trouvé avec la valeur 31 pour l'exposant). Depuis mon ordinateur portable possède un processeur à deux noyaux, et une seule est en cours d'exécution à 50 pour cent (pleine charge pour une base) I mais que je puisse accélérer les calculs en calculant les résultats parallèles.

Alors j'ai changé ma fonction perfectnumber comme suit:

//Now the function again, but async for parallel computing
let perfectNumbersAsync ( n : int) =
    async {
        try
            for x in 1.. n do
                if PowShift x - 1I |> isPrime then
                    let result = PowShift (x-1) * ((PowShift x)-1I)
                    printfn "Found %A as a perfect number" result
        with
            | ex -> printfn "Error%s" (ex.Message);
    }

Pour appeler cette fonction, j'utilise une petite fonction d'aide pour l'exécuter:

 let runPerfects n =
    [n]
        |> Seq.map perfectNumbersAsync
        |> Async.Parallel
        |> Async.RunSynchronously
        |> ignore

Lorsque le résultat du calcul est ignoré async, depuis que je suis dans le l'afficher fonction perfectNumbersAsync.

Le code ci-dessus compile et il fonctionne, mais il utilise toujours qu'un seul noyau (bien qu'il fonctionne 10 secondes plus vite lors du calcul du 9 nombre parfait). Je crains qu'il doit faire quelque chose avec les fonctions d'aide PowShift et isPrime, mais je ne suis pas certain. Est-ce que je dois mettre le code de ces fonctions d'aide dans le bloc de async perfectNumbersAsync? Il n'améliore pas ...

Plus je joue avec F #, plus j'apprends à apprécier cette langue, mais comme ce cas, je suis dans le besoin de certains experts parfois.)

Merci d'avance pour la lecture, j'espère seulement que je me suis fait un peu clai ...

Robert.

Answer 1

Commentaire de @Jeffrey Sax est certainement intéressant, donc je pris un certain temps pour faire une petite expérience. Le test Lucas-Lehmer est écrit comme suit:

let lucasLehmer p =
    let m = (PowShift p) - 1I
    let rec loop i acc =
        if i = p-2 then acc
        else loop (i+1) ((acc*acc - 2I)%m)
    (loop 0 4I) = 0I

Avec le test Lucas-Lehmer, je peux obtenir quelques premiers nombres parfaits très rapide:

let mersenne (i: int) =     
    if i = 2 || (isPrime (bigint i) && lucasLehmer i) then
        let p = PowShift i
        Some ((p/2I) * (p-1I))
    else None

let runPerfects n =
    seq [1..n]
        |> Seq.choose mersenne
        |> Seq.toArray

let m1 = runPerfects 2048;; // Real: 00:00:07.839, CPU: 00:00:07.878, GC gen0: 112, gen1: 2, gen2: 1

Le test Lucas-Lehmer aide à réduire le temps de vérifier les nombres premiers. Au lieu de tester la divisibilité de 2 ^ p-1 qui prend O(sqrt(2^p-1)), nous utilisons le test de primalité qui est au plus O(p^3). Avec n = 2048, je suis en mesure de trouver les premiers 15 numéros Mersenne en 7,83 secondes. Le 15 nombre de Mersenne est celui avec i = 1279 et il se compose de 770 chiffres.

J'ai essayé de paralléliser runPerfects en utilisant le module PSeq F # PowerPack . PSeq ne conserve pas l'ordre de la séquence d'origine, afin d'être honnête, j'ai trié la séquence de sortie. Étant donné que le test de primalité est un équilibre tout à fait entre les indices, le résultat est tout à fait encourageant:

#r "FSharp.Powerpack.Parallel.Seq.dll"    
open Microsoft.FSharp.Collections

let runPerfectsPar n =
    seq [1..n]
        |> PSeq.choose mersenne
        |> PSeq.sort (* align with sequential version *)
        |> PSeq.toArray 

let m2 = runPerfectsPar 2048;; // Real: 00:00:02.288, CPU: 00:00:07.987, GC gen0: 115, gen1: 1, gen2: 0

Avec la même entrée, la version parallèle a pris 2,28 secondes, ce qui est équivalent à 3,4x sur ma machine speedup quad-core. Je crois que le résultat pourrait être améliorée si vous utilisez construction Parallel.For et la partition de la plage d'entrée sensiblement.

Answer 2

Un commentaire rapide sur la vitesse et parallelisability,

Votre isPrime est O (sqrt (n)), et chaque successive n est d'environ 2 x aussi grand que le dernier, donc prendra environ 1,5 x aussi longtemps à calculer, ce qui signifie que le calcul des derniers numéros prendra beaucoup plus

J'ai fait quelques bidouillages avec le test de primalité et certaines choses que j'ai trouvé qui sont utiles sont:

1. Pour grand N, (vous les numéros des tests avec 20 chiffres), la densité première est en fait assez faible, de sorte que vous allez faire beaucoup de divisions par des nombres composés. Une meilleure approche consiste à précalculer une table de nombres premiers (en utilisant un tamis) jusqu'à une limite maximale (probablement déterminée par la quantité de mémoire). Notez que vous êtes le plus susceptible de trouver des facteurs avec un petit nombre. Une fois que vous manquez de mémoire pour votre table, vous pouvez tester le reste des numéros avec votre fonction existante, avec un point de départ plus grande.

2. Une autre approche consiste à utiliser plusieurs threads dans la vérification. Par exemple, vous vérifiez actuellement x,x+4,x+6... facteurs. En étant un peu plus intelligent, vous pouvez faire le congruent numérique à 1 mod 3 en 1 fil et le nombre congru à 2 mod 3 dans un autre thread.

Non

. 2 est le plus simple, mais n ° 1 est plus efficace et offre un potentiel pour faire le flux de contrôle avec OutOfMemoryExceptions qui peut toujours être intéressant

EDIT: Donc, je mis en œuvre ces deux idées, il trouve 2305843008139952128 presque instantanément, trouver 2658455991569831744654692615953842176 prend 7 minutes sur mon ordinateur (quad core AMD 3200). La plupart du temps est consacré à la vérification 2 ^ 61 est premier, donc un meilleur algorithme serait probablement préférable de vérifier les nombres premiers: Code ici

let swatch = new System.Diagnostics.Stopwatch()
swatch.Start()
let inline PowShift (exp:int32) = 1I <<< exp ;;
let limit = 10000000 //go to a limit, makes table gen slow, but should pay off
printfn "making table"
//returns an array of all the primes up to limit
let table =
    let table = Array.create limit true //use bools in the table to save on memory
    let tlimit = int (sqrt (float limit)) //max test no for table, ints should be fine
    table.[1] <- false //special case
    [2..tlimit] 
    |> List.iter (fun t -> 
        if table.[t]  then //simple optimisation
            let mutable v = t*2
            while v < limit do
                table.[v] <- false
                v <- v + t)
    let out = Array.create (50847534) 0I //wolfram alpha provides pi(1 billion) - want to minimize memory
    let mutable idx = 0
    for x in [1..(limit-1)] do
        if table.[x] then
            out.[idx] <- bigint x
            idx <- idx + 1
    out |> Array.filter (fun t -> t <> 0I) //wolfram no is for 1 billion as limit, we use a smaller number
printfn "table made"

let rec isploop testprime incr max n=
    if testprime > max then true
    else if n % testprime = 0I then false
    else isploop (testprime + incr) incr max n

let isPrime ( n : bigint) = 
    //first test the table
    let maxFactor = bigint(sqrt(float n))
    match table |> Array.tryFind (fun t -> n % t = 0I && t <= maxFactor) with
    |Some(t) -> 
        false
    |None -> //now slow test
        //I have 4 cores so
        let bases = [|limit;limit+1;limit+3;limit+4|] //uses the fact that 10^x congruent to 1 mod 3
        //for 2 cores, drop last 2 terms above and change 6I to 3I
        match bases |> Array.map (fun t -> async {return isploop (bigint t) 6I maxFactor n}) |> Async.Parallel |> Async.RunSynchronously |> Array.tryFind (fun t -> t = false) with
        |Some(t) -> false
        |None -> true


let pcount = ref 0
let perfectNumbersTwo (n : int) =  
    seq { for i in 2..n do 
           if (isPrime (bigint i)) then
                if (PowShift i) - 1I |> isPrime then
                    pcount := !pcount + 1
                    if !pcount = 9 then
                        swatch.Stop()
                        printfn "total time %f seconds, %i:%i m:s"  (swatch.Elapsed.TotalSeconds) (swatch.Elapsed.Minutes) (swatch.Elapsed.Seconds)
                    yield PowShift (i-1) * ((PowShift i)-1I)
        } 


perfectNumbersTwo 62 |> Seq.iter (printfn "PERFECT: %A") //62 gives 9th number

printfn "done"
System.Console.Read() |> ignore