742 Stimmen

Geschwindigkeitsvergleich mit Projekt Euler: C vs. Python vs. Erlang vs. Haskell

Ich habe die Problem Nr. 12 de Projekt Euler als Programmierübung und zum Vergleich meiner (sicherlich nicht optimalen) Implementierungen in C, Python, Erlang und Haskell. Um eine höhere Ausführungszeit zu erreichen, suche ich nach der ersten Dreieckszahl mit mehr als 1000 Teilern statt 500, wie im ursprünglichen Problem angegeben.

Das Ergebnis ist das folgende:

C:

lorenzo@enzo:~/erlang$ gcc -lm -o euler12.bin euler12.c
lorenzo@enzo:~/erlang$ time ./euler12.bin
842161320

real    0m11.074s
user    0m11.070s
sys 0m0.000s

Python:

lorenzo@enzo:~/erlang$ time ./euler12.py 
842161320

real    1m16.632s
user    1m16.370s
sys 0m0.250s

Python mit PyPy:

lorenzo@enzo:~/Downloads/pypy-c-jit-43780-b590cf6de419-linux64/bin$ time ./pypy /home/lorenzo/erlang/euler12.py 
842161320

real    0m13.082s
user    0m13.050s
sys 0m0.020s

Erlang:

lorenzo@enzo:~/erlang$ erlc euler12.erl 
lorenzo@enzo:~/erlang$ time erl -s euler12 solve
Erlang R13B03 (erts-5.7.4) [source] [64-bit] [smp:4:4] [rq:4] [async-threads:0] [hipe] [kernel-poll:false]

Eshell V5.7.4  (abort with ^G)
1> 842161320

real    0m48.259s
user    0m48.070s
sys 0m0.020s

Haskell:

lorenzo@enzo:~/erlang$ ghc euler12.hs -o euler12.hsx
[1 of 1] Compiling Main             ( euler12.hs, euler12.o )
Linking euler12.hsx ...
lorenzo@enzo:~/erlang$ time ./euler12.hsx 
842161320

real    2m37.326s
user    2m37.240s
sys 0m0.080s

Zusammenfassung:

  • C: 100%
  • Python: 692% (118% mit PyPy)
  • Erlang: 436% (135% Dank an RichardC)
  • Haskell: 1421%

Ich nehme an, dass C einen großen Vorteil hat, da es für die Berechnungen lange und nicht beliebig lange ganze Zahlen verwendet, wie die anderen drei. Außerdem muss es nicht erst eine Laufzeitumgebung laden (wie die anderen?).

Frage 1: Verlieren Erlang, Python und Haskell durch die Verwendung von Ganzzahlen beliebiger Länge an Geschwindigkeit, oder verlieren sie nicht, solange die Werte kleiner als MAXINT ?

Frage 2: Warum ist Haskell so langsam? Gibt es ein Compiler-Flag, das die Bremsen ausschaltet, oder liegt es an meiner Implementierung? (Letzteres ist sehr wahrscheinlich, da Haskell für mich ein Buch mit sieben Siegeln ist).

Frage 3: Können Sie mir Tipps geben, wie ich diese Implementierungen optimieren kann, ohne die Art und Weise zu ändern, wie ich die Faktoren bestimme? Optimierung in jeder Hinsicht: schöner, schneller, "nativer" in der Sprache.

EDITAR:

Frage 4: Erlauben meine funktionalen Implementierungen LCO (Last Call Optimization, auch bekannt als Tail Recursion Elimination) und vermeiden damit das Hinzufügen unnötiger Frames auf dem Call Stack?

Ich habe wirklich versucht, denselben Algorithmus so ähnlich wie möglich in den vier Sprachen zu implementieren, obwohl ich zugeben muss, dass meine Kenntnisse in Haskell und Erlang sehr begrenzt sind.


Verwendete Quellcodes:

#include <stdio.h>
#include <math.h>

int factorCount (long n)
{
    double square = sqrt (n);
    int isquare = (int) square;
    int count = isquare == square ? -1 : 0;
    long candidate;
    for (candidate = 1; candidate <= isquare; candidate ++)
        if (0 == n % candidate) count += 2;
    return count;
}

int main ()
{
    long triangle = 1;
    int index = 1;
    while (factorCount (triangle) < 1001)
    {
        index ++;
        triangle += index;
    }
    printf ("%ld\n", triangle);
}

#! /usr/bin/env python3.2

import math

def factorCount (n):
    square = math.sqrt (n)
    isquare = int (square)
    count = -1 if isquare == square else 0
    for candidate in range (1, isquare + 1):
        if not n % candidate: count += 2
    return count

triangle = 1
index = 1
while factorCount (triangle) < 1001:
    index += 1
    triangle += index

print (triangle)

-module (euler12).
-compile (export_all).

factorCount (Number) -> factorCount (Number, math:sqrt (Number), 1, 0).

factorCount (_, Sqrt, Candidate, Count) when Candidate > Sqrt -> Count;

factorCount (_, Sqrt, Candidate, Count) when Candidate == Sqrt -> Count + 1;

factorCount (Number, Sqrt, Candidate, Count) ->
    case Number rem Candidate of
        0 -> factorCount (Number, Sqrt, Candidate + 1, Count + 2);
        _ -> factorCount (Number, Sqrt, Candidate + 1, Count)
    end.

nextTriangle (Index, Triangle) ->
    Count = factorCount (Triangle),
    if
        Count > 1000 -> Triangle;
        true -> nextTriangle (Index + 1, Triangle + Index + 1)  
    end.

solve () ->
    io:format ("~p~n", [nextTriangle (1, 1) ] ),
    halt (0).

factorCount number = factorCount' number isquare 1 0 - (fromEnum $ square == fromIntegral isquare)
    where square = sqrt $ fromIntegral number
          isquare = floor square

factorCount' number sqrt candidate count
    | fromIntegral candidate > sqrt = count
    | number `mod` candidate == 0 = factorCount' number sqrt (candidate + 1) (count + 2)
    | otherwise = factorCount' number sqrt (candidate + 1) count

nextTriangle index triangle
    | factorCount triangle > 1000 = triangle
    | otherwise = nextTriangle (index + 1) (triangle + index + 1)

main = print $ nextTriangle 1 1

9voto

user5994461 Punkte 3652

Einige weitere Zahlen und Erklärungen für die C-Version. Anscheinend hat das in all den Jahren niemand gemacht. Vergessen Sie nicht, diese Antwort hochzuvoten, damit sie für alle sichtbar und lernbar wird.

Schritt Eins: Benchmarking der Programme des Autors

Laptop-Spezifikationen:

  • CPU i3 M380 (931 MHz - maximaler Batteriesparmodus)
  • 4 GB Speicher
  • Win7 64 Bit
  • Microsoft Visual Studio 2012 Ultimate
  • Cygwin mit gcc 4.9.3
  • Python 2.7.10

Befehle:

compiling on VS x64 command prompt > `for /f %f in ('dir /b *.c') do cl /O2 /Ot /Ox %f -o %f_x64_vs2012.exe`
compiling on cygwin with gcc x64   > `for f in ./*.c; do gcc -m64 -O3 $f -o ${f}_x64_gcc.exe ; done`
time (unix tools) using cygwin > `for f in ./*.exe; do  echo "----------"; echo $f ; time $f ; done`

.

----------
$ time python ./original.py

real    2m17.748s
user    2m15.783s
sys     0m0.093s
----------
$ time ./original_x86_vs2012.exe

real    0m8.377s
user    0m0.015s
sys     0m0.000s
----------
$ time ./original_x64_vs2012.exe

real    0m8.408s
user    0m0.000s
sys     0m0.015s
----------
$ time ./original_x64_gcc.exe

real    0m20.951s
user    0m20.732s
sys     0m0.030s

Die Dateinamen lauten: integertype_architecture_compiler.exe

  • Integer-Typ ist vorerst das gleiche wie das Originalprogramm (mehr dazu später)
  • Architektur ist x86 oder x64, abhängig von den Compiler-Einstellungen
  • Compiler ist gcc oder vs2012

Schritt zwei: Untersuchen, verbessern und erneut bewerten

VS ist 250% schneller als gcc. Die beiden Compiler sollten eine ähnliche Geschwindigkeit liefern. Offensichtlich stimmt etwas mit dem Code oder den Compiler-Optionen nicht. Lasst uns das untersuchen!

Der erste interessante Punkt sind die Integer-Typen. Konvertierungen können teuer sein und Konsistenz ist wichtig für eine bessere Codegenerierung und Optimierung. Alle Ganzzahlen sollten vom gleichen Typ sein.

Es ist ein Mischmasch aus int y long gerade jetzt. Das werden wir noch verbessern. Welcher Typ soll verwendet werden? Den schnellsten. Wir müssen sie alle miteinander vergleichen!

----------
$ time ./int_x86_vs2012.exe

real    0m8.440s
user    0m0.016s
sys     0m0.015s
----------
$ time ./int_x64_vs2012.exe

real    0m8.408s
user    0m0.016s
sys     0m0.015s
----------
$ time ./int32_x86_vs2012.exe

real    0m8.408s
user    0m0.000s
sys     0m0.015s
----------
$ time ./int32_x64_vs2012.exe

real    0m8.362s
user    0m0.000s
sys     0m0.015s
----------
$ time ./int64_x86_vs2012.exe

real    0m18.112s
user    0m0.000s
sys     0m0.015s
----------
$ time ./int64_x64_vs2012.exe

real    0m18.611s
user    0m0.000s
sys     0m0.015s
----------
$ time ./long_x86_vs2012.exe

real    0m8.393s
user    0m0.015s
sys     0m0.000s
----------
$ time ./long_x64_vs2012.exe

real    0m8.440s
user    0m0.000s
sys     0m0.015s
----------
$ time ./uint32_x86_vs2012.exe

real    0m8.362s
user    0m0.000s
sys     0m0.015s
----------
$ time ./uint32_x64_vs2012.exe

real    0m8.393s
user    0m0.015s
sys     0m0.015s
----------
$ time ./uint64_x86_vs2012.exe

real    0m15.428s
user    0m0.000s
sys     0m0.015s
----------
$ time ./uint64_x64_vs2012.exe

real    0m15.725s
user    0m0.015s
sys     0m0.015s
----------
$ time ./int_x64_gcc.exe

real    0m8.531s
user    0m8.329s
sys     0m0.015s
----------
$ time ./int32_x64_gcc.exe

real    0m8.471s
user    0m8.345s
sys     0m0.000s
----------
$ time ./int64_x64_gcc.exe

real    0m20.264s
user    0m20.186s
sys     0m0.015s
----------
$ time ./long_x64_gcc.exe

real    0m20.935s
user    0m20.809s
sys     0m0.015s
----------
$ time ./uint32_x64_gcc.exe

real    0m8.393s
user    0m8.346s
sys     0m0.015s
----------
$ time ./uint64_x64_gcc.exe

real    0m16.973s
user    0m16.879s
sys     0m0.030s

Ganzzahlige Typen sind int long int32_t uint32_t int64_t y uint64_t de #include <stdint.h>

Es gibt viele Integer-Typen in C, plus einige vorzeichenbehaftete/unvorzeichenbehaftete, mit denen man spielen kann, plus die Möglichkeit, als x86 oder x64 zu kompilieren (nicht zu verwechseln mit der tatsächlichen Integer-Größe). Das sind eine Menge Versionen zum Kompilieren und Ausführen ^^

Dritter Schritt: Die Zahlen verstehen

Endgültige Schlussfolgerungen:

  • 32-Bit-Ganzzahlen sind ~200% schneller als 64-Bit-Äquivalente
  • ohne Vorzeichen 64 Bits Ganzzahlen sind 25 % schneller als vorzeichenbehaftet 64 Bit (Leider habe ich keine Erklärung dafür)

Fangfrage: "Was sind die Größen von int und long in C?"
Die richtige Antwort lautet: Die Größe von int und long in C ist nicht genau definiert!

Aus der C-Spezifikation:

int besteht aus mindestens 32 Bits
long ist mindestens ein int

Aus der gcc-Manualseite (-m32- und -m64-Flags):

Die 32-Bit-Umgebung setzt int, long und pointer auf 32 Bit und erzeugt Code, der auf jedem i386-System läuft.
Die 64-Bit-Umgebung setzt int auf 32 Bit und long und pointer auf 64 Bit und generiert Code für die x86-64-Architektur von AMD.

Aus der MSDN-Dokumentation (Data Type Ranges) https://msdn.microsoft.com/en-us/library/s3f49ktz%28v=vs.110%29.aspx :

int, 4 Bytes, auch bekannt als signed
long, 4 Bytes, auch bekannt als long int und signed long int

Um dies abzuschließen: Gelernte Lektionen

  • 32-Bit-Ganzzahlen sind schneller als 64-Bit-Ganzzahlen.

  • Standard-Ganzzahlentypen sind weder in C noch in C++ gut definiert, sie variieren je nach Compiler und Architektur. Wenn Sie Konsistenz und Vorhersagbarkeit benötigen, verwenden Sie die uint32_t Ganzzahlige Familie aus #include <stdint.h> .

  • Geschwindigkeitsprobleme gelöst. Alle anderen Sprachen sind wieder Hunderte von Prozent im Rückstand, C & C++ gewinnen wieder! Das tun sie immer. Die nächste Verbesserung wird Multithreading mit OpenMP sein :D

8voto

Mark Washeim Punkte 101

Frage 1: Verlieren Erlang, Python und Haskell an Geschwindigkeit durch die Verwendung von durch die Verwendung von Ganzzahlen beliebiger Länge oder nicht, solange die Werte kleiner als MAXINT sind?

Frage eins kann für Erlang mit Nein beantwortet werden. Die letzte Frage wird durch eine entsprechende Verwendung von Erlang beantwortet, wie in:

http://bredsaal.dk/learning-erlang-using-projecteuler-net

Da es schneller ist als Ihr ursprüngliches C-Beispiel, würde ich vermuten, dass es zahlreiche Probleme gibt, wie andere bereits ausführlich beschrieben haben.

Dieses Erlang-Modul wird auf einem billigen Netbook in etwa 5 Sekunden ausgeführt ... Es verwendet das Netzwerk-Thread-Modell in Erlang und demonstriert so, wie man die Vorteile des Ereignismodells nutzen kann. Es könnte über viele Knoten verteilt werden. Und es ist schnell. Nicht mein Code.

-module(p12dist).  
-author("Jannich Brendle, jannich@bredsaal.dk, http://blog.bredsaal.dk").  
-compile(export_all).

server() ->  
  server(1).

server(Number) ->  
  receive {getwork, Worker_PID} -> Worker_PID ! {work,Number,Number+100},  
  server(Number+101);  
  {result,T} -> io:format("The result is: \~w.\~n", [T]);  
  _ -> server(Number)  
  end.

worker(Server_PID) ->  
  Server_PID ! {getwork, self()},  
  receive {work,Start,End} -> solve(Start,End,Server_PID)  
  end,  
  worker(Server_PID).

start() ->  
  Server_PID = spawn(p12dist, server, []),  
  spawn(p12dist, worker, [Server_PID]),  
  spawn(p12dist, worker, [Server_PID]),  
  spawn(p12dist, worker, [Server_PID]),  
  spawn(p12dist, worker, [Server_PID]).

solve(N,End,_) when N =:= End -> no_solution;

solve(N,End,Server_PID) ->  
  T=round(N*(N+1)/2),
  case (divisor(T,round(math:sqrt(T))) > 500) of  
    true ->  
      Server_PID ! {result,T};  
    false ->  
      solve(N+1,End,Server_PID)  
  end.

divisors(N) ->  
  divisor(N,round(math:sqrt(N))).

divisor(_,0) -> 1;  
divisor(N,I) ->  
  case (N rem I) =:= 0 of  
  true ->  
    2+divisor(N,I-1);  
  false ->  
    divisor(N,I-1)  
  end.

Der nachstehende Test fand an einem: Intel(R) Atom(TM) CPU N270 @ 1.60GHz

~$ time erl -noshell -s p12dist start

The result is: 76576500.

^C

BREAK: (a)bort (c)ontinue (p)roc info (i)nfo (l)oaded
       (v)ersion (k)ill (D)b-tables (d)istribution
a

real    0m5.510s
user    0m5.836s
sys 0m0.152s

6voto

thanos Punkte 713

Versuchen Sie GO:

package main

import "fmt"
import "math"

func main() {
    var n, m, c int
    for i := 1; ; i++ {
        n, m, c = i * (i + 1) / 2, int(math.Sqrt(float64(n))), 0
        for f := 1; f < m; f++ {
            if n % f == 0 { c++ }
    }
    c *= 2
    if m * m == n { c ++ }
    if c > 1001 {
        fmt.Println(n)
        break
        }
    }
}

Ich verstehe:

ursprüngliche C-Version: 9.1690 100%
gehen: 8.2520 111%

Aber mit:

package main

import (
    "math"
    "fmt"
 )

// Sieve of Eratosthenes
func PrimesBelow(limit int) []int {
    switch {
        case limit < 2:
            return []int{}
        case limit == 2:
            return []int{2}
    }
    sievebound := (limit - 1) / 2
    sieve := make([]bool, sievebound+1)
    crosslimit := int(math.Sqrt(float64(limit))-1) / 2
    for i := 1; i <= crosslimit; i++ {
        if !sieve[i] {
            for j := 2 * i * (i + 1); j <= sievebound; j += 2*i + 1 {
                sieve[j] = true
            }
        }
    }
    plimit := int(1.3*float64(limit)) / int(math.Log(float64(limit)))
    primes := make([]int, plimit)
    p := 1
    primes[0] = 2
    for i := 1; i <= sievebound; i++ {
        if !sieve[i] {
            primes[p] = 2*i + 1
            p++
            if p >= plimit {
                break
            }
        }
    }
    last := len(primes) - 1
    for i := last; i > 0; i-- {
        if primes[i] != 0 {
            break
        }
        last = i
    }
    return primes[0:last]
}

func main() {
    fmt.Println(p12())
}
// Requires PrimesBelow from utils.go
func p12() int {
    n, dn, cnt := 3, 2, 0
    primearray := PrimesBelow(1000000)
    for cnt <= 1001 {
        n++
        n1 := n
        if n1%2 == 0 {
            n1 /= 2
        }
        dn1 := 1
        for i := 0; i < len(primearray); i++ {
            if primearray[i]*primearray[i] > n1 {
                dn1 *= 2
                break
            }
            exponent := 1
            for n1%primearray[i] == 0 {
                exponent++
                n1 /= primearray[i]
            }
            if exponent > 1 {
                dn1 *= exponent
            }
            if n1 == 1 {
                break
            }
        }
        cnt = dn * dn1
        dn = dn1
    }
    return n * (n - 1) / 2
}

Ich verstehe:

ursprüngliche C-Version: 9.1690 100%
thaumkid's c Version: 0.1060 8650%
erste Version: 8.2520 111%
zweite Go-Version: 0.0230 39865%

Ich habe auch Python3.6 und pypy3.3-5.5-alpha ausprobiert:

ursprüngliche C-Version: 8.629 100%
Version c von thaumkid: 0.109 7916%
Python3.6: 54.795 16%
pypy3.3-5.5-alpha: 13.291 65%

und dann mit folgendem Code erhielt ich:

ursprüngliche C-Version: 8.629 100%
Version c von thaumkid: 0.109 8650%
Python3.6: 1.489 580%
pypy3.3-5.5-alpha: 0.582 1483%

def D(N):
    if N == 1: return 1
    sqrtN = int(N ** 0.5)
    nf = 1
    for d in range(2, sqrtN + 1):
        if N % d == 0:
            nf = nf + 1
    return 2 * nf - (1 if sqrtN**2 == N else 0)

L = 1000
Dt, n = 0, 0

while Dt <= L:
    t = n * (n + 1) // 2
    Dt = D(n/2)*D(n+1) if n%2 == 0 else D(n)*D((n+1)/2)
    n = n + 1

print (t)

5voto

user3125280 Punkte 2759

C++11, < 20ms für mich - Führen Sie es hier aus

Ich verstehe, dass Sie Tipps zur Verbesserung Ihrer sprachspezifischen Kenntnisse haben möchten, aber da das hier bereits gut behandelt wurde, dachte ich, ich würde etwas Kontext für Leute hinzufügen, die sich den Mathematica-Kommentar zu Ihrer Frage usw. angesehen und sich gefragt haben, warum dieser Code so viel langsamer war.

Diese Antwort dient hauptsächlich dazu, Kontext zu liefern, um hoffentlich den Leuten zu helfen, den Code in Ihrer Frage / anderen Antworten leichter zu bewerten.

Dieser Code verwendet nur ein paar (hässliche) Optimierungen, die nichts mit der verwendeten Sprache zu tun haben und auf der Grundlage von:

  1. jede Traingle-Zahl ist von der Form n(n+1)/2
  2. n und n+1 sind koprim
  3. die Anzahl der Teiler ist eine multiplikative Funktion

    include <iostream>

    include <cmath>

    include <tuple>

    include <chrono>

    using namespace std;

    // Calculates the divisors of an integer by determining its prime factorisation.

    int get_divisors(long long n) { int divisors_count = 1;

    for(long long i = 2;
        i <= sqrt(n);
        /* empty */)
    {
        int divisions = 0;
        while(n % i == 0)
        {
            n /= i;
            divisions++;
        }
    
        divisors_count *= (divisions + 1);
    
        //here, we try to iterate more efficiently by skipping
        //obvious non-primes like 4, 6, etc
        if(i == 2)
            i++;
        else
            i += 2;
    }
    
    if(n != 1) //n is a prime
        return divisors_count * 2;
    else
        return divisors_count;

    }

    long long euler12() { //n and n + 1 long long n, n_p_1;

    n = 1; n_p_1 = 2;
    
    // divisors_x will store either the divisors of x or x/2
    // (the later iff x is divisible by two)
    long long divisors_n = 1;
    long long divisors_n_p_1 = 2;
    
    for(;;)
    {
        /* This loop has been unwound, so two iterations are completed at a time
         * n and n + 1 have no prime factors in common and therefore we can
         * calculate their divisors separately
         */
    
        long long total_divisors;                 //the divisors of the triangle number
                                                  // n(n+1)/2
    
        //the first (unwound) iteration
    
        divisors_n_p_1 = get_divisors(n_p_1 / 2); //here n+1 is even and we
    
        total_divisors =
                  divisors_n
                * divisors_n_p_1;
    
        if(total_divisors > 1000)
            break;
    
        //move n and n+1 forward
        n = n_p_1;
        n_p_1 = n + 1;
    
        //fix the divisors
        divisors_n = divisors_n_p_1;
        divisors_n_p_1 = get_divisors(n_p_1);   //n_p_1 is now odd!
    
        //now the second (unwound) iteration
    
        total_divisors =
                  divisors_n
                * divisors_n_p_1;
    
        if(total_divisors > 1000)
            break;
    
        //move n and n+1 forward
        n = n_p_1;
        n_p_1 = n + 1;
    
        //fix the divisors
        divisors_n = divisors_n_p_1;
        divisors_n_p_1 = get_divisors(n_p_1 / 2);   //n_p_1 is now even!
    }
    
    return (n * n_p_1) / 2;

    }

    int main() { for(int i = 0; i < 1000; i++) { using namespace std::chrono; auto start = high_resolution_clock::now(); auto result = euler12(); auto end = high_resolution_clock::now();

        double time_elapsed = duration_cast<milliseconds>(end - start).count();
    
        cout << result << " " << time_elapsed << '\n';
    }
    return 0;

    }

Das dauert durchschnittlich 19 ms auf meinem Desktop und 80 ms auf meinem Laptop, was weit von den meisten anderen Codes entfernt ist, die ich hier gesehen habe. Und es sind zweifellos noch viele Optimierungen möglich.

1voto

gnasher729 Punkte 49452

Ich bin davon ausgegangen, dass die Anzahl der Faktoren nur dann groß ist, wenn die beteiligten Zahlen viele kleine Faktoren haben. Daher habe ich den hervorragenden Algorithmus von thaumkid verwendet, aber zunächst eine Annäherung an die Faktorenzahl vorgenommen, die niemals zu klein ist. Das ist ganz einfach: Man sucht nach Primfaktoren bis 29, dann prüft man die verbleibende Zahl und berechnet eine Obergrenze für die Anzahl der Faktoren. Daraus errechnet man eine Obergrenze für die Anzahl der Faktoren, und wenn diese Zahl hoch genug ist, berechnet man die genaue Anzahl der Faktoren.

Der folgende Code braucht diese Annahme nicht, um korrekt zu sein, sondern um schnell zu sein. Es scheint zu funktionieren; nur etwa eine von 100.000 Zahlen ergibt eine Schätzung, die hoch genug ist, um eine vollständige Prüfung zu erfordern.

Hier ist der Code:

// Return at least the number of factors of n.
static uint64_t approxfactorcount (uint64_t n)
{
    uint64_t count = 1, add;

#define CHECK(d)                            \
    do {                                    \
        if (n % d == 0) {                   \
            add = count;                    \
            do { n /= d; count += add; }    \
            while (n % d == 0);             \
        }                                   \
    } while (0)

    CHECK ( 2); CHECK ( 3); CHECK ( 5); CHECK ( 7); CHECK (11); CHECK (13);
    CHECK (17); CHECK (19); CHECK (23); CHECK (29);
    if (n == 1) return count;
    if (n < 1ull * 31 * 31) return count * 2;
    if (n < 1ull * 31 * 31 * 37) return count * 4;
    if (n < 1ull * 31 * 31 * 37 * 37) return count * 8;
    if (n < 1ull * 31 * 31 * 37 * 37 * 41) return count * 16;
    if (n < 1ull * 31 * 31 * 37 * 37 * 41 * 43) return count * 32;
    if (n < 1ull * 31 * 31 * 37 * 37 * 41 * 43 * 47) return count * 64;
    if (n < 1ull * 31 * 31 * 37 * 37 * 41 * 43 * 47 * 53) return count * 128;
    if (n < 1ull * 31 * 31 * 37 * 37 * 41 * 43 * 47 * 53 * 59) return count * 256;
    if (n < 1ull * 31 * 31 * 37 * 37 * 41 * 43 * 47 * 53 * 59 * 61) return count * 512;
    if (n < 1ull * 31 * 31 * 37 * 37 * 41 * 43 * 47 * 53 * 59 * 61 * 67) return count * 1024;
    if (n < 1ull * 31 * 31 * 37 * 37 * 41 * 43 * 47 * 53 * 59 * 61 * 67 * 71) return count * 2048;
    if (n < 1ull * 31 * 31 * 37 * 37 * 41 * 43 * 47 * 53 * 59 * 61 * 67 * 71 * 73) return count * 4096;
    return count * 1000000;
}

// Return the number of factors of n.
static uint64_t factorcount (uint64_t n)
{
    uint64_t count = 1, add;

    CHECK (2); CHECK (3);

    uint64_t d = 5, inc = 2;
    for (; d*d <= n; d += inc, inc = (6 - inc))
        CHECK (d);

    if (n > 1) count *= 2; // n must be a prime number
    return count;
}

// Prints triangular numbers with record numbers of factors.
static void printrecordnumbers (uint64_t limit)
{
    uint64_t record = 30000;

    uint64_t count1, factor1;
    uint64_t count2 = 1, factor2 = 1;

    for (uint64_t n = 1; n <= limit; ++n)
    {
        factor1 = factor2;
        count1 = count2;

        factor2 = n + 1; if (factor2 % 2 == 0) factor2 /= 2;
        count2 = approxfactorcount (factor2);

        if (count1 * count2 > record)
        {
            uint64_t factors = factorcount (factor1) * factorcount (factor2);
            if (factors > record)
            {
                printf ("%lluth triangular number = %llu has %llu factors\n", n, factor1 * factor2, factors);
                record = factors;
            }
        }
    }
}

Dreieckszahl mit 13824 Faktoren in etwa 0,7 Sekunden, die 879.207.615. Dreieckszahl mit 61.440 Faktoren in 34 Sekunden, die 12.524.486.975. Dreieckszahl mit 138.240 Faktoren in 10 Minuten und 5 Sekunden und die 26.467.792.064. Dreieckszahl mit 172.032 Faktoren in 21 Minuten und 25 Sekunden (2,4GHz Core2 Duo), so dass dieser Code im Durchschnitt nur 116 Prozessorzyklen pro Zahl benötigt. Die letzte Dreieckszahl selbst ist größer als 2^68, also

CodeJaeger.com

CodeJaeger ist eine Gemeinschaft für Programmierer, die täglich Hilfe erhalten..
Wir haben viele Inhalte, und Sie können auch Ihre eigenen Fragen stellen oder die Fragen anderer Leute lösen.

Powered by:

X