Sehr schneller Primzahl-Finder

[GTA-Place]

Auch euch will ich meinen Source nicht vorbehalten:

Delphi-Quellcode:

function Prim(Zahl: Cardinal): Boolean;  
var  
  Teiler: PCardinal;  
  Wurzel: Cardinal;  
begin  
  Result := True;    // Result = True  
  if not odd(Zahl) OR (Zahl <= 5) then    // Ist die Zahl nich ungerade oder kleiner als 5, dann  
  begin  
    if (Zahl <> 2) AND (Zahl <> 3) AND (Zahl <> 5) then    // Ist die Zahl nicht 2 und nicht 3 und nicht 5, dann  
      Result := False;    // Result = False  
    Exit;                 // Exit  
  end;  
  Teiler := @PrimS[0];    // Teiler = @PrimS[0]  
  Wurzel := Trunc(sqrt(Zahl));    // Wurzel aus der Zahl  
  while Teiler^ <= Wurzel do    // Solange Teiler^ <= Wurzel ist, mache  
  begin      
    if Zahl mod Teiler^ = 0 then    // Ist Zahl / Teiler^ = Rest 0, dann  
    begin  
      Result := False;    // Result = False  
      Break;    // Break  
    end;  
    Inc(Teiler);    // Teiler erhöhen um 1  
  end;  
end;

Delphi-Quellcode:

procedure TMainForm.StartButClick(Sender: TObject);  
var  
  Start, Ende: Real;  
  Z: PCardinal;  
  X, LPrim: Cardinal;  
  PrimF: TStringList;  
begin  
  try  
    Von := StrToInt(VonEdit.Text);    // Start  
    Bis := StrToInt(BisEdit.Text);    // Endwert  

 
    if Bis > 10 then  
      SetLength(PrimS, Trunc(0.4*Bis-(Bis/4)))    // Größe des Array: 0.4*Bis-(Bis/4)  
    else  
      SetLength(PrimS, 4);  

 
    LPrim := 0;    // Letze Prim = 0  
    Z := @PrimS[0];     // Gefundene Prims = 0  

 
    if (Von >= 0) AND (Bis >= 0) AND (Von < Bis) then    // Von >= 0; Bis >= 0; Von < Bis;  
    begin  
      Start := GetTickCount;    // Start-Zeit  

 
      for X := Von to Bis do    // Schleife: Startwert -> Endwert  
      begin  
        if Prim(X) then    // Funktion ausführen, wenn Prim dann  
        begin  
          Z^ := X;    // Prim in Array schreiben  
          Inc(Z);    // Z erhöhen um 1  
          LPrim := X;    // Letze Prim = X  
        end;  
        if X mod 20000 = 0 then    // Wenn X : 20.000 = Rest 0, dann  
        begin  
          Application.ProcessMessages;    // Anwendung aktualisieren  
          PrimLab.Caption := 'Aktuelle Primzahl: ' + IntToStr(LPrim);    // Akt. Primzahl anzeigen  
        end;  
      end;  

 
      Ende := GetTickCount;    // End-Zeit  
      DauerLab.Caption := 'Diese Überprüfung hat ' + FloatToStr((Ende - Start) / 1000) + ' Sekunden gedauert.';    // Dauer anzeigen  

 
      PrimLab.Caption := 'Speichern...';    // "Speichern..." anzeigen  

 
      Z := @PrimS[0];    // Z auf 0 stellen  
      PrimF := TStringList.Create;    // Stringlist erzeugen  

 
      for X := 0 to Length(PrimS)-1 do    // Von 0 bis Größe des Array  
      begin  
        if Z^ = 0 then    // Wenn Z^ = 0, dann  
          Break;    // Break  
        PrimF.Add(IntToStr(Z^));    // Prim in Stringlist schreiben  
        Inc(Z);    // Z erhöhen um 1  
      end;  

 
      PrimF.SaveToFile('Prim.txt');    // Stringlist speichern  
      PrimF.Free;    // Stringlist freigeben  

 
      PrimLab.Caption := 'Aktuelle Primzahl: ' + IntToStr(LPrim);    // Akt. Primzahl anzeigen  
    end  
    else  
      ShowMessage('Ungültige Eingabe(n)!');    // Bei falschen Eingaben, Nachricht anzeigen  
  except  
    ShowMessage('Es ist ein Fehler aufgetreten!');    // Wenn Fehler auftritt, Nachricht anzeigen  
  end;  
end;

Das Programm überprüft 10.000.000 Zahlen in erstaunlichen 6 Sekunden.
Und das Speichern geht so schnell, dass man "Speichern..." gar nicht sieht.

[GTA-Place]

Das DF (besonders Phantom1) hat mal weiter probiert und Phantom1 hat diese Funktion geschrieben:

Delphi-Quellcode:

procedure SavePrimes(MaxPrime: Cardinal; const FileName: string = '');  
const  
  CACHE = 64*1024;  
  STEMPEL: array[0..7] of Byte = (1, 7, 11, 13, 17, 19, 23, 29);  
  MODS: array[0..29] of Byte = (0, 1, 0, 0, 0, 0, 0, 2, 0, 0, 0, 4, 0, 8, 0, 0,  
                                0, 16, 0, 32, 0, 0, 0, 64, 0, 0, 0, 0, 0, 128);  
var  
  Primes, PrimesLUT: array of Byte;  
  i, j, k, PrimeLen, PrimeBits, Num, Num2, m, mbit, s: Cardinal;  
  f: TextFile;  
begin  
  if FileName<>'' then begin  
    AssignFile(f, FileName);  
    ReWrite(f);  
    WriteLn(f, '2'+#10#13+'3'+#10#13+'5');  
  end;  
 
  SetLength(PrimesLUT, Trunc(Sqrt(MaxPrime)/30)); // max 2184 Byte für 2^32 ;-)  
  PrimesLUT[0]:=$01;  
  PrimeLen:=Length(PrimesLUT);  
  PrimeBits:=PrimeLen*30;  
  for i:=0 to Trunc(Sqrt(PrimeBits)/30) do  
    for j:=0 to 7 do  
      if PrimesLUT[i] and (1 shl j)=0 then begin  
        s:=STEMPEL[j];  
        Num:=i*30+s;  
        Num2:=Num*Num;  
        mbit:=Num2 mod 30;  
        m:=(Num2-mbit) div 30;  
        while m<PrimeLen do begin  
          PrimesLUT[m]:=PrimesLUT[m] or MODS[mbit];  
          Inc(m, i);  
          Inc(mbit, s);  
          if mbit>29 then begin  
            Dec(mbit, 30);  
            Inc(m);  
          end;  
        end;  
      end;  

  SetLength(Primes, CACHE);  
  PrimeLen:=Length(Primes);  
  PrimeBits:=PrimeLen*30;  
  for k:=0 to MaxPrime div PrimeBits do begin  
    FillChar(Primes[0], PrimeLen, 0);  
    for i:=0 to Trunc(Sqrt((k+1)*PrimeBits)/30) do  
      for j:=0 to 7 do  
        if PrimesLUT[i] and (1 shl j)=0 then begin  
          s:=STEMPEL[j];  
          Num:=i*30+s;  
          if k=0 then  
            Num2:=Num*Num  
          else  
            Num2:=Trunc(k*PrimeBits/Num)*Num+Num;  
          mbit:=Num2 mod 30;  
          m:=(Num2-mbit) div 30-k*PrimeLen;  
          while m<PrimeLen do begin  
            primes[m]:=Primes[m] or MODS[mbit];  
            Inc(m, i);  
            Inc(mbit, s);  
            if mbit>29 then begin  
              Dec(mbit, 30);  
              Inc(m);  
            end;  
          end;  
        end;  
    if FileName<>'' then  
      for i:=0 To PrimeLen-1 do  
        for j:=0 to 7 do begin  
          if k*PrimeBits+i*30+STEMPEL[j]>MaxPrime then  
            Break;  
          if not ((i=0) and (j=0) and (k=0)) and (Primes[i] and (1 shl j)=0) then  
            WriteLn(f, IntToStr(k*PrimeBits+i*30+STEMPEL[j]));  
        end;  
  end;  
  if FileName<>'' then  
    CloseFile(f);  
end;

Um 500.000.000 Zahlen zu überprüfen, braucht diese Funktion nur ca. 4 Sekunden.


EDIT: Code wurde von Phantom1 verbessert: 500 mio Zahlen in ca. 2.65 Sekunden.

[negaH]

Toll, das ist in fact ein 8/30 Comb Sieve, so wie es in meinem Posting im Attachment zu finden ist, sehr wahrscheinlich sogar von meinen Source abgekupfert.

Zitat:

 SetLength(PrimesLUT, Trunc(Sqrt(MaxPrime)/30)); // max 2184 Byte für 2^32 ;-)

Kannst du mir mathematisch erklären werum im Remark "max 2184 Bytes für 2^32" steht ? bzw. wie sich diese Schranke auf Grund des Verfahrens ergibt ?


Der Unterschied ist halt nur in der Performance und Flexibilität zu suchen.

Zitat:

EDIT: Code wurde von Phantom1 verbessert: 500 mio Zahlen in ca. 2.65 Sekunden.

Auf welchem Rechner ermittelt ?

Mein obiges Sieb schafft die 500Mio in 2.7 Sekunden auf einem P4 1.5GHz.
Falls die 2.65 Sekunden auch auf einem P4 1.5GHz ermittelt wurden so wundert mich dies schon, da mein Sieb mit hochoptimierten Assembler Routinen arbeitet die im Grunde wesentlich schneller als reiner PASCAL Source sein sollten.

Gruß Hagen

[Phantom1]

Zitat von negaH:

Toll, das ist in fact ein 8/30 Comb Sieve, so wie es in meinem Posting im Attachment zu finden ist, sehr wahrscheinlich sogar von meinen Source abgekupfert.

Ich hab mir dein Code bis gerade eben noch nicht angeschaut. Sieht ziemlich komplex aus, schade das ich kein Assembler kann.

Zitat von negaH:

Kannst du mir mathematisch erklären werum im Remark "max 2184 Bytes für 2^32" steht ? bzw. wie sich diese Schranke auf Grund des Verfahrens ergibt ?

Mein Programm sucht nur primzahlen bis 2^32. Da man nur die Primzahlen bis zur Wurzel(2^32) braucht um die restlichen wegzustreichen zu können, wären das 65536. Da ich ein 30er Stempel verwende ->65536/30 = 2184 Bytes nötig um alle Primzahlen bis 65536 (als bits) zu speichern. Dieses Array (PrimeLUT) benutze ich dann für die spätere Berechnung der einzelnen Array-Abschnitte und spare mir dadurch viele unnötige Streichungen.

Zitat von negaH:

Auf welchem Rechner ermittelt ?
Mein obiges Sieb schafft die 500Mio in 2.7 Sekunden auf einem P4 1.5GHz.
Falls die 2.65 Sekunden auch auf einem P4 1.5GHz ermittelt wurden so wundert mich dies schon, da mein Sieb mit hochoptimierten Assembler Routinen arbeitet die im Grunde wesentlich schneller als reiner PASCAL Source sein sollten.

Mein Rechner ist ein Athlon XP-M@2000MHz. Gemessen wurde mit RealTime-Priority.

mfg
Phantom1

[PierreB]

Zitat:

Gemessen wurde mit RealTime-Priority.

Aähm nur mal kurz ne ganz kleine Frage, ich steh grad aufm Schlauch: Wie kann ich denn die Zeit die ein Algorithmus braucht messen ??

[negaH]

Zitat:

Ich hab mir dein Code bis gerade eben noch nicht angeschaut. Sieht ziemlich komplex aus, schade das ich kein Assembler kann.

Alle relevanten Assemblerparts sind auch auskommentiert in reinem PASCAL in diesem Source enthalten.
Verzichtet man in meiner Implementation auf alle zusätzlichen Features, wie:
- Index Berechnungen von Primzahlen
- speichern, laden, erzeugen einer Cache Datei
- Anpassungen des Sources für D3 etc.
- direkte Berechnung der Residue Tabelle zu einem beliebigen Startwert

und packt dies alles in eine einzigste Funktion, so dürfte ein Source rauskommen der zu deinem Source fast 1 zu 1 identisch ist.

Es verwundert mich eben schon weil ich damals definitiv im Netz keinen einzigsten Source finden konnte der dieses 8/30Comb Sieb implemntierte. Meine Referenz war ein PostScript von D.J.Bernstein, eine mathematiche Abhandlung über dieses Sieb.

Gruß Hagen

[Kara]

<OT>Wisst ihr, was mir zu diesem Thread einfällt?
"Meiner ist läää-schneller" </OT>

[negaH]

@Kara:

die Schnelligkeit einer Implementierung hängt von verschiedenen Faktoren ab:

1.) der mathematische Algorithmus bestimmt die grundsätzliche Laufzeitkomplexität der Implementation
2.) die programmiertechnische Umsetzung der Datenstrukturen und Algorithmen bestimmt die Gesamtlaufzeit

Man kann sich immer hinstellen und die Arbeit anderer diskreditieren ohne es real besser machen zu können oder eben nicht im entferntesten verstehen zu können welche Arbeit hinter einer solchen Umsetzung tatsächlich steckt.

Der schlußendlich Vergleich der Effizienz einer jeweiligen Implementierung hat also rein garnichts "Penis-Längen" zu tuen, sondern mit dem Ehrgeiz eines Programmierers seine Arbeit so gut wie nur möglich zu machen. Das sich daraus direkt auch der Stolz auf die eigene Arbeit ergibt, die man natürlicherweise auch mit der Arbeit anderer vergleichen will, ist doch wohl ein rein menchliches Bedürfnis.

Wenn man nun bedenkt das die Primzahlfindung als Operation wiederum die Laufzeiten der darauf basierenden Algorithmen beeinflusst so ist es ein wichtiges Bedürfnis diese Primzahlberechnungen so effizient wie nur möglich hinzubekommen.

Sich also, wenn auch offtopic, sich über solche Vergleiche lustig zu machen, deutet nur daraufhin das man das grundsätzliche Problem überhaupt nicht verstehen konnte. (schade eigentlich da du anscheinend nicht bereit bist an deinen Schwächen zu arbeiten).

Denn ich behaupte mal, das jeder Poster in diesem Thread, mit seiner eigenen Implementierung sich viel Mühe und Arbeit gemacht hat sie
1.) mathematisch zu begreifen
2.) sich im WEB nach relevanten Informationen auf Suche zu begeben
3.) lange an der Implementierungen rumgefeilt hat um sie so effizient wie möglich hinzubekommen

Es steckt also in jedem Posting hier viel Zeit drinnen.
Statt davon zu lernen kann man sich auch darüber lustig machen, bravo

EDIT:
4.) sie dann hier in der DP der Allgmeinheit zur Verfügung stellt damit Leute wie DU daraus lernen können.

Gruß Hagen

[Kara]

Ist ja gut. Du brauchst dich durch meinen kleinen Einwurf nicht angegriffen fühlen.

[negaH]

Warum ist dann überhaupt ein solcher Einwurf nötig gewesen ?

Einfach nur um was zu posten, oder einfach nur um überhaupt was zu sagen, oder einfach aus Frust weil man es nicht besser kann ?

Alles nur Anzeichen dafür das Du unreif bist.

Gruß Hagen

PS: ich hoffe das meine "Gegenkritik" an deinem Verhalten dich nun ermuntert mal 90 Grad den Kopf nach unten zu kippen und zu schauen wie lang er wirklich ist.

[Kara]

Ich hoffe, dass ich nichts sehe...

Und ich frage mich grad, ob du einen schlechten Tag hast, oder warum meine zwei Zeilen dich zu solchen Antworten treiben.
Aber bevor wir noch länger OT bleiben, ziehe ich mich geschlagen aus diesem Thread zurück und überlasse ihn wieder euren Expertengesprächen.