Проблема с распараллеливанием в F # при вычислении точных чисел?

Комментарий @Jeffrey Sax определенно интересен, поэтому мне потребовалось время, чтобы провести небольшой эксперимент. Тест Лукаса-Лемера записывается следующим образом:

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

С помощью теста Лукаса-Лемера я могу очень быстро получить несколько первых идеальных чисел:

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

Тест Лукаса-Лемера помогает сократить время проверки простых чисел. Вместо проверки делимости 2 ^ p-1, для которой требуется O(sqrt(2^p-1)), мы используем тест на простоту, который составляет не более O(p^3). С n = 2048 я могу найти первые 15 чисел Мерсенна за 7,83 секунды. Пятнадцатое число Мерсенна - это число с i = 1279 и состоит из 770 цифр.

Я попытался распараллелить runPerfects с помощью модуля PSeq в F # Powerpack. PSeq не сохраняет порядок исходной последовательности, поэтому, честно говоря, я отсортировал выходную последовательность. Поскольку тест на простоту представляет собой довольно сбалансированное соотношение показателей, результат весьма обнадеживает:

#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

При том же вводе параллельная версия заняла 2,28 секунды, что эквивалентно увеличению скорости в 3,4 раза на моем четырехъядерном компьютере. Я считаю, что результат можно улучшить, если использовать конструкцию Parallel.For и разумно разделить диапазон ввода.

Один быстрый комментарий по скорости и параллелизму,

Ваш isPrime равен O (sqrt (n)), и каждое последующее n примерно в 2 раза больше, чем последнее, поэтому для расчета потребуется примерно 1,5 раза больше времени, что означает, что вычисление последних чисел займет гораздо больше времени.

Я проделал несколько хакерских работ с тестированием на простоту и обнаружил некоторые полезные вещи:

1. Для большого N (вы тестируете числа с 20 цифрами) плотность простых чисел на самом деле довольно низкая, поэтому вы будете делать много делений на составные числа. Лучшим подходом является предварительное вычисление таблицы простых чисел (с использованием сита) до некоторого максимального предела (вероятно, определяемого объемом памяти). Обратите внимание, что вы, скорее всего, найдете множители с небольшими числами. Как только у вас закончится память для таблицы, вы можете протестировать остальные числа с помощью существующей функции с большей начальной точкой.

2. Другой подход - использовать при проверке несколько потоков. Например, в настоящее время вы проверяете x,x+4,x+6... как факторы. Будучи немного умнее, вы можете сделать число, совпадающее с 1 модулем 3 в 1 потоке, и числа, совпадающее с 2 модулем 3, в другом потоке.

№2 - самый простой, но №1 - более эффективный и предоставляет возможность для выполнения потока управления с OutOfMemoryExceptions, что всегда может быть интересно.

РЕДАКТИРОВАТЬ: Итак, я реализовал обе эти идеи, он находит 2305843008139952128 почти мгновенно, а поиск 2658455991569831744654692615953842176 занимает 7 минут на моем компьютере (четырехъядерный AMD 3200). Большую часть времени тратится на проверку простых чисел 2 ^ 61, поэтому для проверки простых чисел, вероятно, будет лучше лучший алгоритм: Код здесь

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