Autor |
Beitrag |
Mathematiker
Beiträge: 2622
Erhaltene Danke: 1447
Win 7, 8.1, 10
Delphi 5, 7, 10.1
|
Verfasst: Mi 23.08.17 17:10
Hallo,
im Ergebnis der Verbesserung des Programms durch Horst_H läuft jetzt ein Projekt zur Suche nach den kleinsten n-stelligen Primzahlvierlingen, d.h. wir suchen für jede Stellenzahl n = 1, ..., 1000 einen Summanden a, so dass
10^{n-1}+a, 10^{n-1}+a+2, 10^{n-1}+a+6, 10^{n-1}+a+8 jeweils Primzahlen sind.
Im Moment sind außer Horst_H auch zwei vom Matroids Matheplaneten am "rechnen", natürlich ich auch.
Sollte jemand von euch mitmachen wollen, so gibt's unter mathematikalpha.de/primzahlvierlinge ein kleines Programm und eine genauere Erklärung der Suche.
Bevor man mitmacht ist es aber ratsam, sich die aktuelle Liste der Vierlinge anzusehen, da regelmäßig neue Ergebnisse eingehen.
Außer einer lobenden Erwähnung gibt es aber leider nichts zu gewinnen.
Großes Ziel ist es bis zu 1000 Stellen zu kommen. Mehr wäre schön, aber leider steigt der Rechenaufwand immer mehr.
So eine Berechnung wurde bisher noch nicht durchgeführt. D.h., jeder, der einen neuen Wert findet, geht in die "Mathematikgeschichte" ein.
Beste Grüße und (vielleicht) erfolgreiches Suchen
Steffen
_________________ Töten im Krieg ist nach meiner Auffassung um nichts besser als gewöhnlicher Mord. Albert Einstein
|
|
Gammatester
Beiträge: 328
Erhaltene Danke: 101
|
Verfasst: Fr 25.08.17 13:01
Hier ein Vierling bei 10^(620-1)+49370415361 etc. Noch mal separat mit mp_pprime
getestet.
Gruß Gammatester
Edit nach Horst's Beitrag: Noch ein paar Anregungen zum Programm. Wenn noch Zahlen offen, sind sollten diese sofort angezeigt werden, nicht erst nach Suche. Und wichtig: Bei großen Zahlen ist es mM besser, nicht ganze 20-er Blöcke zu bearbeiten, da man muß schon großes Glück haben. Ich habe gestern ohne Erfolg 620-639 bearbeitet, bis mir klar wurde, daß ich ca 20-mal schneller bin, wenn ich nur 620 beackere, habe dann Dein Programm ausgetrickst
Edit 2 Horst_H: ich habe nur Anhalts-Punkte für die Zeit: Ich bin heute um ca 9:00 mit 620,7228500000 gestartet und vor ca 45 min fertig geworden. Vielleicht eine weitere Anregung für die nächste Programmversion.
Zuletzt bearbeitet von Gammatester am Fr 25.08.17 13:26, insgesamt 4-mal bearbeitet
Für diesen Beitrag haben gedankt: Horst_H, Mathematiker
|
|
Horst_H
Beiträge: 1653
Erhaltene Danke: 243
WIN10,PuppyLinux
FreePascal,Lazarus
|
Verfasst: Fr 25.08.17 13:12
Hallo,
da würde mich die Laufzeit interessieren.
Gerade eben auf i3-4330 3.5 Ghz in fast 12h fertig:
Quelltext 1:
| 553 1447073137021 11:53:33.600 |
Ich schätze um die 30 min == 12h *( 620/553)^2.2 ) * 49.37/1447
Gruß Horst
Für diesen Beitrag haben gedankt: Mathematiker
|
|
Mathematiker
Beiträge: 2622
Erhaltene Danke: 1447
Win 7, 8.1, 10
Delphi 5, 7, 10.1
|
Verfasst: Fr 25.08.17 14:02
Gammatester hat folgendes geschrieben : | Hier ein Vierling bei 10^(620-1)+49370415361 etc. Noch mal separat mit mp_pprime
getestet. |
Danke. Ist schon in die Liste eingetragen.
Steffen
_________________ Töten im Krieg ist nach meiner Auffassung um nichts besser als gewöhnlicher Mord. Albert Einstein
|
|
Horst_H
Beiträge: 1653
Erhaltene Danke: 243
WIN10,PuppyLinux
FreePascal,Lazarus
|
Verfasst: So 27.08.17 14:41
Hallo,
ein kleines Testprogramm, um GMP und mpArith zu vergleichen:
1: 2: 3: 4: 5: 6: 7: 8: 9: 10: 11: 12: 13: 14: 15: 16: 17: 18: 19: 20: 21: 22: 23: 24: 25: 26: 27: 28: 29: 30: 31: 32: 33: 34: 35: 36: 37: 38: 39: 40: 41: 42: 43: 44: 45: 46: 47: 48: 49: 50: 51: 52: 53: 54: 55: 56: 57: 58: 59: 60: 61: 62: 63: 64: 65: 66: 67: 68: 69: 70: 71: 72: 73: 74: 75: 76: 77: 78: 79: 80: 81: 82: 83: 84: 85: 86: 87: 88: 89: 90: 91: 92: 93: 94: 95: 96: 97: 98: 99: 100: 101: 102: 103: 104: 105: 106: 107: 108: 109: 110: 111:
| program GMP_Test; {$IFNDEF FPC} {$Apptype CONSOLE} {$ELSE} {$MODE DELPHI} {$ENDIF}
uses sysutils,gmp,mp_base, mp_numth, mp_prng, mp_prime, mp_types; type tTestWerte = record twSZ, twVs : NativeUint; end; const cTestwerte : array[1..10] of tTestWerte = ((twSz: 100;twVs:2191), (twSz: 200;twVs:5491), (twSz: 300;twVs:4681), (twSz: 400;twVs:7081), (twSz: 500;twVs:15721), (twSz: 600;twVs:2161), (twSz: 700;twVs:20251), (twSz: 800;twVs:9511), (twSz: 900;twVs:6991), (twSz:1000;twVs:19411)); {$IFDEF CPU32} runden = 25; {$ENDIF} {$IFDEF CPU64} runden = 100; {$ENDIF} var x :mpz_t; y: mp_int;
procedure TestlaufMpArith(var TestZahl:mp_int); var T1,T0: TDateTime; i: integer; Begin t0:= time; For i := 1 to runden do s_mp_is_psp2(TestZahl); t1:= time; Writeln(86400*1000*(T1-T0)/runden:8:6,' ms'); end;
procedure TestlaufGMP(var TestZahl:mpz_t); var T1,T0: TDateTime; i: integer; Begin t0:= time; i := 1; For i := 1 to runden do mpz_probab_prime_p(TestZahl, 0); t1:= time; Writeln(86400*1000*(T1-T0)/runden:8:6,' ms'); end;
var i : integer; begin
mpz_init(x); writeln('GMP'); For i := low(cTestwerte) to High(cTestwerte) do Begin mpz_set_ui(x,10); with cTestwerte[i] do Begin mpz_pow_ui(x,x,twSz-1); mpz_add_ui(x,x,twVs); write(twSz:5,' '); end; TestlaufGMP(x); end; mpz_clear(x);
writeln; writeln('MpArith'); mp_init(y); For i := low(cTestwerte) to High(cTestwerte) do Begin mp_set(y, 10); with cTestwerte[i] do Begin mp_expt_int(y,twSz-1,y); mp_add_d(y,twVs,y); write(twSz:5,' '); end; TestlaufMpArith(y); end; mp_clear(y);
end. |
Laufzeit bei mir, während die Suche auf einem Thread weiterlief.
1: 2: 3: 4: 5: 6: 7: 8: 9: 10: 11: 12: 13: 14: 15: 16: 17: 18: 19: 20: 21: 22: 23: 24: 25: 26: 27: 28: 29: 30: 31: 32: 33: 34: 35: 36: 37: 38: 39: 40: 41: 42: 43: 44: 45: 46: 47: 48: 49: 50: 51:
| CPU64 # ./GMP_Test GMP 100 0.030000 ms 200 0.170000 ms 300 0.500000 ms 400 1.080000 ms 500 2.050000 ms 600 3.590000 ms 700 5.340000 ms 800 7.630000 ms 900 10.660000 ms 1000 14.320000 ms
MpArith 100 0.150000 ms 200 0.880000 ms 300 3.400000 ms 400 6.990000 ms 500 12.950000 ms 600 29.610000 ms 700 39.780000 ms 800 52.050000 ms 900 77.090000 ms 1000 104.640000 ms
CPU32 # ./GMP_Test GMP 100 0.120000 ms 200 0.880000 ms 300 2.760000 ms 400 5.920000 ms 500 10.720000 ms 600 18.560000 ms 700 28.560000 ms 800 33.480000 ms 900 50.280000 ms 1000 82.480000 ms
MpArith 100 0.440000 ms 200 3.000000 ms 300 10.720000 ms 400 22.400000 ms 500 40.600000 ms 600 184.160000 ms 700 215.800000 ms 800 278.840000 ms 900 554.600000 ms 1000 674.480000 ms |
Was bleibt gmp ist fast 7x schneller.
Statt 13 Stunden dann 2 ist doch ein Wort. 560 Stellen rödelt schon 16h bei 1.922e12.
Gruß Horst
EDIT: Programm geändert für 32-Bit
Zuletzt bearbeitet von Horst_H am So 27.08.17 15:48, insgesamt 1-mal bearbeitet
Für diesen Beitrag haben gedankt: Mathematiker
|
|
Mathematiker
Beiträge: 2622
Erhaltene Danke: 1447
Win 7, 8.1, 10
Delphi 5, 7, 10.1
|
Verfasst: So 27.08.17 14:49
Hallo Horst,
Horst_H hat folgendes geschrieben : | Hallo,
Was bleibt gmp ist fast 7x schneller.
Statt 13 Stunden dann 2 ist doch ein Wort. 560 Stellen rödelt schon 16h bei 1.922e12.
|
Das ist hervorragend.
Mein Problem nun, wie komme ich an eine Exe zum Ausprobierenheran (64 Bit?).
Liebe Grüße
Steffen
_________________ Töten im Krieg ist nach meiner Auffassung um nichts besser als gewöhnlicher Mord. Albert Einstein
|
|
Delphi-Laie
Beiträge: 1600
Erhaltene Danke: 232
Delphi 2 - RAD-Studio 10.1 Berlin
|
Verfasst: So 27.08.17 15:18
|
|
Horst_H
Beiträge: 1653
Erhaltene Danke: 243
WIN10,PuppyLinux
FreePascal,Lazarus
|
Verfasst: So 27.08.17 16:00
Hallo,
das Programm "program GMP_Test;" 2 Antworten höher habe ich jetzt geändert, einfach mögliche Primzahlen suchen lassen.
Unter 32-Bit sind die Laufzeiten aber erheblich "schlechter", da ist ja nochmal ein Faktor 6 gegenüber 64-Bit drin.
Besonders auffällig ist der Sprung bei mparith bei CPU32 bei Stellenzahl 500-> 600 40->184 ms.
Die Laufzeiten sind ja fast proportional der Stellenzahl n^3 (1000/500 -> 82/10,7 = 2^2,94/ oder 104/12,95 = 2^ 3,03 )
Gruß Horst
Für diesen Beitrag haben gedankt: Mathematiker
|
|
Mathematiker
Beiträge: 2622
Erhaltene Danke: 1447
Win 7, 8.1, 10
Delphi 5, 7, 10.1
|
Verfasst: So 27.08.17 16:26
Hallo,
ich stelle mich wohl wieder zu doof an.
Ich bin auf der Seite gmplib.org/ und lade die Datei gmp-6.1.2.tar.lz herunter. Ob das für Delphi ist, weiß ich nicht, ist auch egal, da ich das Ding nicht entpacken kann, auch nicht mit dem 7-Zip-Dateimanager.
Und ohne gmp geht es nicht.
Beste Grüße
Steffen
Hast sich erledigt, ich bin eben doch zu doof.
_________________ Töten im Krieg ist nach meiner Auffassung um nichts besser als gewöhnlicher Mord. Albert Einstein
|
|
Delphi-Laie
Beiträge: 1600
Erhaltene Danke: 232
Delphi 2 - RAD-Studio 10.1 Berlin
|
Verfasst: So 27.08.17 17:56
gmp-6.1.2.tar.bz2 läßt sich mit Winrar öffnen.
Allerdings ist das eine C-Bibliothek und mithin für mich uninteressant.
Oder wie handhabt Ihr diese Quelltexte?
|
|
Lelf
Beiträge: 42
Erhaltene Danke: 21
|
Verfasst: So 27.08.17 20:49
Hier gibt es einen Gmp-wrapper-for-delphi:
code.google.com/arch...loads-host/downloads
Gruß
Lelf
Für diesen Beitrag haben gedankt: Delphi-Laie, Mathematiker
|
|
Gammatester
Beiträge: 328
Erhaltene Danke: 101
|
Verfasst: Mo 28.08.17 11:18
Für diesen Beitrag haben gedankt: Horst_H, Mathematiker
|
|
Gammatester
Beiträge: 328
Erhaltene Danke: 101
|
Verfasst: Mo 28.08.17 16:01
Gammatester hat folgendes geschrieben : | Der Sprung könnte beim Wechsel von Karatsuba auf Toom-3 verursacht werden.
|
Horst_H, ich kann Deinen Sprung hier nicht nachvollziehen. Hier meine Werte fur den Test mit expt_mod auf i3/2.3 Ghz (Delphi 3 rechnet mit mp_digit bits 15/16, alle anderen mit mp_digit bits 31/32, FPC32/T=4K: FPC32 mit Toom3 = 4*Karatsuba).
Quelltext 1: 2: 3: 4: 5: 6: 7: 8: 9: 10: 11:
| D3/16 FPC32 D18/32 FPC64 FPC32/T=4K 100 1.342 1.938 2.459 1.093 1.950 200 6.541 6.042 7.057 4.074 5.551 300 17.410 14.102 20.856 7.424 14.079 400 39.333 28.126 45.369 14.153 28.101 500 77.455 52.632 86.937 24.796 52.558 600 128.781 95.760 150.708 48.929 98.861 700 199.945 142.172 227.475 69.307 147.441 800 293.470 205.109 336.141 98.148 210.664 900 424.481 327.736 476.608 170.889 287.117 1000 568.588 456.911 642.765 207.466 385.704 |
Für FPC32 scheint das Anheben der Toom-3-Grenze etwas zu bringen.
Gruß Gammatester
Edit: Ich habe die Suche für 621 aufgegeben, nachdem nach tagelanger Suche die Suchlistenanzeige auf "1,-1" gesprungen ist (Programmversion vom 23.08. 07:20). Ist das ein bekanntes Problem?
|
|
Mathematiker
Beiträge: 2622
Erhaltene Danke: 1447
Win 7, 8.1, 10
Delphi 5, 7, 10.1
|
Verfasst: Mo 28.08.17 16:31
Hallo Gammatester,
Gammatester hat folgendes geschrieben : |
Ich habe die Suche für 621 aufgegeben, nachdem nach tagelanger Suche die Suchlistenanzeige auf "1,-1" gesprungen ist (Programmversion vom 23.08. 07:20). Ist das ein bekanntes Problem? |
Das ist neu. Im Moment kann ich mir den Grund noch nicht erklären, werde aber den Fehler suchen.
Unter mathematikalpha.de/primzahlvierlinge habe ich eine neuere Variante, die aller drei Minuten den Stand unter suchliste_kopie.txt sicherheitshalber speichert.
Beste Grüße
Steffen
_________________ Töten im Krieg ist nach meiner Auffassung um nichts besser als gewöhnlicher Mord. Albert Einstein
|
|
Horst_H
Beiträge: 1653
Erhaltene Danke: 243
WIN10,PuppyLinux
FreePascal,Lazarus
|
Verfasst: Di 29.08.17 11:22
Hallo,
ich versuche ja, das sieben zu beschleunigen, da ich den PrimzahlTest selbst nicht beschleunigen kann.
Und, wer hätte das gedacht, das geht auch
4e9 durchsiebt er bei mir in 1,5 Sekunden. ( Edit4 Freepascal Linux 64-Bit oder 32-Bit gleich )
Berechnung der nächsten Position stark vereinfacht, aber nicht unbedingt verständlicher.
1: 2: 3: 4: 5: 6: 7: 8: 9: 10: 11: 12: 13: 14: 15: 16: 17: 18: 19: 20: 21: 22: 23: 24: 25: 26: 27: 28: 29: 30: 31: 32: 33: 34: 35: 36: 37: 38: 39: 40: 41: 42: 43: 44: 45: 46: 47: 48: 49: 50: 51: 52: 53: 54: 55: 56: 57: 58: 59: 60: 61: 62: 63: 64: 65: 66: 67: 68: 69: 70: 71: 72: 73: 74: 75: 76: 77: 78: 79: 80: 81: 82: 83: 84: 85: 86: 87: 88: 89: 90: 91: 92: 93: 94: 95: 96: 97: 98: 99: 100: 101: 102: 103: 104: 105: 106: 107: 108: 109: 110: 111: 112: 113: 114: 115: 116: 117: 118: 119: 120: 121: 122: 123: 124: 125: 126: 127: 128: 129: 130: 131: 132: 133: 134: 135: 136: 137: 138: 139: 140: 141: 142: 143: 144: 145: 146: 147: 148: 149: 150: 151: 152: 153: 154: 155: 156: 157: 158: 159: 160: 161: 162: 163: 164: 165: 166: 167: 168: 169: 170: 171: 172: 173: 174: 175: 176: 177: 178: 179: 180: 181: 182: 183: 184: 185: 186: 187: 188: 189: 190: 191: 192: 193: 194: 195: 196: 197: 198: 199: 200: 201: 202: 203: 204: 205: 206: 207: 208: 209: 210: 211: 212: 213: 214: 215: 216: 217: 218: 219: 220: 221: 222: 223: 224: 225: 226: 227: 228: 229: 230: 231: 232: 233: 234: 235: 236: 237: 238: 239: 240: 241: 242: 243: 244: 245: 246: 247: 248: 249: 250: 251: 252: 253: 254: 255: 256: 257: 258: 259: 260: 261: 262: 263: 264: 265: 266: 267: 268: 269: 270: 271: 272: 273: 274: 275: 276: 277: 278: 279: 280: 281: 282: 283: 284: 285: 286: 287: 288: 289: 290: 291: 292: 293: 294: 295: 296: 297: 298: 299: 300: 301: 302: 303: 304: 305: 306: 307: 308: 309: 310: 311: 312:
| program FindeQuadrupel;
{$IFDEF FPC} {$MODE DELPHI} {$OPTIMIZATION ON,REGVAR,CSE,PEEPHOLE} {$ELSE} {$APPTYPE CONSOLE} {$ENDIF} uses sysutils; type tSieveElem = NativeUint;
tDelta30 = array[0..3] of NativeUint; tMulOfs = record moMul,moOfs: byte; end; tMulPosArr = array[0..3] of tMulOfs;
tPrimMod30 =packed record pm30Offset : LongWord; pmNextSivbNr : word; pm30AktIdx, pm30Delta : byte; end; tpPrimMod30 = ^tPrimMod30;
const cMulOfs : array [0..7]of tMulPosArr = (((moMul: 2;moOfs:0),(moMul: 4;moOfs: 0),(moMul: 2;moOfs:0),(moMul:22;moOfs: 1)), ((moMul: 4;moOfs:1),(moMul: 8;moOfs: 2),(moMul: 4;moOfs:1),(moMul:14;moOfs: 3)), ((moMul: 6;moOfs:2),(moMul:16;moOfs: 6),(moMul: 6;moOfs:2),(moMul: 2;moOfs: 1)), ((moMul:12;moOfs:5),(moMul: 4;moOfs: 2),(moMul:12;moOfs:5),(moMul: 2;moOfs: 1)), ((moMul:12;moOfs:7),(moMul: 4;moOfs: 2),(moMul:12;moOfs:7),(moMul: 2;moOfs: 1)), ((moMul: 6;moOfs:4),(moMul:16;moOfs:10),(moMul: 6;moOfs:4),(moMul: 2;moOfs: 1)), ((moMul: 4;moOfs:3),(moMul: 8;moOfs: 6),(moMul: 4;moOfs:3),(moMul:14;moOfs:11)), ((moMul: 2;moOfs:2),(moMul: 4;moOfs: 4),(moMul: 2;moOfs:2),(moMul:22;moOfs:21))); cMod30ToIdx : array[0..29] of byte = (111,0,111,111,111,111,111,1,111,111,111,2,111,3,111, 111,111,4,111,5,111,111,111,6,111,111,111,111,111,7);
cIdxToMod30 : array[0..7] of byte = (1,7,11,13,17,19,23,29);
cBitSize = SizeOf(tSieveElem)*8; cAndMask = cBitSize-1;
cSiebGroesse = 64*30*130208; cPrimeMax =63245;var SiebMod30 : array[0..(cSiebGroesse-1) DIV 30 DIV cBitSize] of tSieveElem; Prim30List :array[0..11895200] of tPrimMod30; PrimSieb : array of boolean; T0,T1: TDAteTime; gblP30, Prim30idx, LastPrime : longWord;
tmpDelta30:tDelta30;
function Auszaehlen:Longint; var i:NativeUInt; elem:tSieveElem; begin result := 0; Begin For i := Low(SiebMod30) to High(SiebMod30)-1 do Begin elem := NOT SiebMod30[i]; repeat inc(result,elem AND 1); elem := elem shr 1; until elem = 0; end; elem := NOT SiebMod30[High(SiebMod30)]; i := High(SiebMod30)*30*cBitsize; while i < cSiebGroesse do Begin inc(result,elem AND 1); elem := elem shr 1; inc(i,30); end; end; end;
procedure MulListErstellen(p30:LongWord); var pMulPos : ^tMulOfs; begin pMulPos := @cMulOfs[p30 AND 7][0]; p30 := p30 shr 3; with pMulPos^ do tmpDelta30[0] := moOfs+p30*moMul; inc(pMulPos); with pMulPos^ do tmpDelta30[1] := moOfs+p30*moMul; inc(pMulPos); with pMulPos^ do tmpDelta30[2] := moOfs+p30*moMul; inc(pMulPos); with pMulPos^ do tmpDelta30[3] := moOfs+p30*moMul; end; procedure PrimSieben; var i,p: NativeInt; Begin p := 2; repeat i := p*p; while i <= cPrimeMax do Begin PrimSieb[i] := true; inc(i,p); end; repeat inc(p); until NOT PrimSieb[p] until (p*p) >cPrimeMax; end;
procedure PrimListeErstellen; var i: NativeInt; p30,p30Idx,DIV30,Aktidx,Last30Prime: NativeUInt; Begin T0 := now; setlength(PrimSieb,cPrimeMax+1); PrimSieben;
Prim30idx := 0; LastPrime := 1; Last30Prime := 0; For i := 7 to cPrimeMax do IF Not PrimSieb[i] then Begin DIV30 := i DIV 30; p30Idx := cMod30ToIdx[i-30*DIV30]; p30 := DIV30 shl 3+p30Idx; MulListErstellen(p30); IF cMulOfs[p30Idx][3].moMul = 2 then Begin inc(DIV30,tmpDelta30[0]); Aktidx:= 1 end else Begin DIV30:= tmpDelta30[3] shr 1; Aktidx:= 0; end; with Prim30List[Prim30idx] do Begin pm30Offset:= Div30; pm30Delta := p30-Last30Prime; pm30AktIdx:= Aktidx; end; inc(Prim30idx); if Prim30idx>High(Prim30List) then Begin BREAK; end; Last30Prime:= p30; LastPrime := i; end; dec(Prim30idx); setlength(PrimSieb,0); T1 := now; writeln('Anzahl Primzahlen ',Prim30idx,' Letzte PrimZahl ', lastPrime); Writeln('Laufzeit sieben',FormatDateTime('S.ZZZ',T1-T0)); end;
procedure EinmalSieben(var Prim30:tPrimMod30); var p30,mytmp,AktIdx: NativeUInt; begin with Prim30 do Begin MulListErstellen(gblP30); p30 := pm30Offset; AktIdx := pm30AktIdx; end; while p30 < cSiebGroesse DIV 30 do Begin mytmp := p30 DIV cBitSize; SiebMod30[mytmp]:= SiebMod30[mytmp] OR (tSieveElem(1) shl (p30 AND cAndMask)); inc(p30,tmpDelta30[AktIdx]); AktIdx := (AktIdx + 1) AND 3; end; with Prim30 do Begin pm30Offset:= p30-cSiebGroesse DIV 30; pm30AktIdx := AktIdx; end;
end;
function AlleSieben:longint; var pPL : tpPrimMod30; i,k,cnt : NativeInt;
begin k := 0; cnt := 0; repeat fillchar(SiebMod30,SizeOf(SiebMod30),0); i := 0; pPL := @Prim30List[0]; gblP30 := 0; repeat inc(gblP30,pPL^.pm30Delta); If pPL^.pm30Offset < cSiebGroesse DIV 30 then EinmalSieben(pPl^) else dec(pPl^.pm30Offset,cSiebGroesse DIV 30); inc(pPL); inc(i); until i > Prim30idx; inc(cnt,Auszaehlen); inc(k); until k >=8; result := cnt; end;
var cnt:NativeUInt; BEGIN PrimListeErstellen; cnt := AlleSieben;
writeln('Limit : ',cSiebGroesse); writeln('Suchfeldlimit : ',(High(SiebMod30)+1)*30*cBitSize); writeln('Anzahl : ',cnt); end. |
Jetzt muss nur noch eine clevere Lösung gefunden werden, dass für die Quadrupelsuche zu nutzen.
Die Platzersparnis für das Sieb ist ja enorm.Wenn man auf Bytebene bleibt ist es ein Faktor 30 auf Bitebene nochmals 8.
Bei 1MB levelII-Cache bei mir könnten da jetzt um 240 Mio Zahlen passen.
AMD-Rechner vor Ryzen sind leider im Speicherzugriff relativ lahm.Je "näher" an der CPU, desto besser.
Für jede Streichprimzahl mir auch noch OffsetMod30[0..3] und den momentanen Index [0..3] zu merken bläht die Primzahlliste auf, aber die wird ja sequentiell und damit schnell gelesen.
Vielleicht wieder holt sich das Muster wieder alle 30 0R1(=1),1R1(=31),2R1(=61)..30R1(=901)
Mal schauen
Gruß Horst
Edit:
Das Bitfeld ist jetzt implementiert.Jetzt fehlt die Berechnung der Überträge.
Edit2->3:
Die Uebtraege habe ich immer noch nicht.Aber ich kämpfe mit einer kompakten Speicherung der Siebprimzahlen.
Jetzt 8 Byte,weil ich die Differenz zum Vorgänger speichere.pm30Delta als Byte entspricht maximal 959 als Abstand.
Das ist irgendwo bei 1e15 der Fall...
Delphi-Quelltext 1: 2: 3: 4: 5: 6: 7: 8:
| type tPrimMod30 =packed record pm30Offset : LongWord; pmNextSivbNr : word; pm30AktIdx, pm30Delta : byte; end; |
EDIT3: Es geht einfacher, siehe Programm.
Zuletzt bearbeitet von Horst_H am Di 05.09.17 18:28, insgesamt 2-mal bearbeitet
Für diesen Beitrag haben gedankt: Mathematiker
|
|
Mathematiker
Beiträge: 2622
Erhaltene Danke: 1447
Win 7, 8.1, 10
Delphi 5, 7, 10.1
|
Verfasst: So 03.09.17 22:03
Hallo,
ich habe die Primzahlvierlingssuche jetzt auf 2 Threads aufgeteilt.
Der erste siebt jeweils einen Bereich von 500000 mit 10 Millionen Primzahlen vor und ermittelt die "Kandidaten" für einen Vierling, d.h. in dem Zehner dürfen die auf 1,3,7,9 endenden Zahlen nicht ausgesiebt sein.
Die Kandidaten werden über eine threadsichere Stringliste an den 2.Thread übergeben.
Dieser testet nun die 4 Zahlen eines Kandidaten auf Primzahleigenschaft und entfernt den getesteten Kandidaten aus der Stringliste.
Es funktioniert eigentlich ganz gut.
Leider ist der 2.Thread deutlich langsamer als der erste. Damit es nicht zu einem Stau in der Stringliste kommt, prüfe ich im ersten Thread mit
Delphi-Quelltext 1: 2:
| while liste.count>50 do begin end; | die Listeneinträge und bremse den 1.Thread bei mehr als 50 Kandidaten.
Wichtig ist das vor allem, wenn die Berechnung abgebrochen wird, da dann erst alle Kandidaten abgearbeitet werden müssen, d.h. der Hauptthread wartet bis der 2. fertig ist.
Schön ist das nicht.
Meine Frage nun: Seht ihr eine andere Möglichkeit, beide Threads auf etwa gleiche Geschwindigkeit zu bringen? Gibt es eine andere Möglichkeit, den ersten Thread mit dem 2. zu synchonisieren?
Der angehängte Quelltext ist Lazarus-64-Bit. Diese Exen sind deutlich schneller als bei meinem Delphi 5.
Danke
Steffen
Nachtrag: Nach 19 min Laufzeit stürzt es jetzt auch noch ab. Warum? Ich weiß es nicht.
Nebenbei: Unsere kleine Gruppe, die erfolgreich nach Primzahlvierlingen sucht, benötigt noch Rechenleistung. Mitstreiter sind herzlich willkommen.
siehe: matheplanet.com/math...lps=1680408#v1680408
bzw. mathematikalpha.de/primzahlvierlinge
Einloggen, um Attachments anzusehen!
_________________ Töten im Krieg ist nach meiner Auffassung um nichts besser als gewöhnlicher Mord. Albert Einstein
Zuletzt bearbeitet von Mathematiker am Mo 04.09.17 19:56, insgesamt 1-mal bearbeitet
|
|
Delphi-Laie
Beiträge: 1600
Erhaltene Danke: 232
Delphi 2 - RAD-Studio 10.1 Berlin
|
Verfasst: Mo 04.09.17 01:48
Statt Polling mit der while-Schleife sind Botschaften auf der Sendeseite und deren ressourcenschonender Empfang mit Getmessage ratsam. So programmierte ich einige parallele Sortieralgorithmen.
Ansonsten kannst Du es mit vertauschter Reihenfolge der Threaderstellung und / oder mit Änderung ihrer versuchen.
Generell läßt sich Windows aber in bezug auf die Abarbeitung der Threads und Tasks kaum bis gar nicht hineinreden, und dementsprechend sind die Möglichkeit, da einzugreifen, äußerst begrenzt und nur "pseudo".
Zuletzt bearbeitet von Delphi-Laie am Mo 04.09.17 15:54, insgesamt 1-mal bearbeitet
|
|
pzktupel
Hält's aus hier
Beiträge: 129
Erhaltene Danke: 30
|
Verfasst: Mo 04.09.17 06:17
Hallo,
ich bin neu hier und habe mich mal hier angemeldet.
Durch eine Anfrage mit meinem kleinsten 1000-stelligen Primzahl-4-Tupel , hatte ich mir
2 Tage mühe gemacht um was zu programmieren, um bis n=999 alle zu finden.
Im Moment schafft mein Programm für jeden Exponenten in derselben Zeit ca. k= 300 000 000 000 / Stunde für 10^n + k abzusieben.... und das bei nur 1 Kern einer modernen CPU. Bei RYZEN wahrscheinlich bis 1 Billion pro Stunde.
Wollte ich nur mal erwähnen...
Gruß
|
|
cryptnex
Beiträge: 23
Erhaltene Danke: 5
|
Verfasst: Mo 04.09.17 11:55
pzktupel hat folgendes geschrieben : | Durch eine Anfrage mit meinem kleinsten 1000-stelligen Primzahl-4-Tupel , hatte ich mir
2 Tage mühe gemacht um was zu programmieren, um bis n=999 alle zu finden.
Wollte ich nur mal erwähnen... |
Würdest Du den Quellcode zum Studieren hier veröffentlichen? Ich führe nur ungern Programme aus unbekannten Quellen auf meinem Rechner aus.
Viele Grüße
|
|
pzktupel
Hält's aus hier
Beiträge: 129
Erhaltene Danke: 30
|
Verfasst: Mo 04.09.17 12:56
Ok, es ist Freebasiccode, mit C++ habe ich es nicht so.
Anmerkung aber: Das ist ein Ergebnis nach Jahren der Tüftellei. Auch ist dieser unübersichtlich.
Letzte Fassung:
1: 2: 3: 4: 5: 6: 7: 8: 9: 10: 11: 12: 13: 14: 15: 16: 17: 18: 19: 20: 21: 22: 23: 24: 25: 26: 27: 28: 29: 30: 31: 32: 33: 34: 35: 36: 37: 38: 39: 40: 41: 42: 43: 44: 45: 46: 47: 48: 49: 50: 51: 52: 53: 54: 55: 56: 57: 58: 59: 60: 61: 62: 63: 64: 65: 66: 67: 68: 69: 70: 71: 72: 73: 74: 75: 76: 77: 78: 79: 80: 81: 82: 83: 84: 85: 86: 87: 88: 89: 90: 91: 92: 93: 94: 95: 96: 97: 98: 99: 100: 101: 102: 103: 104: 105: 106: 107: 108: 109: 110: 111: 112: 113: 114: 115: 116: 117: 118: 119: 120: 121: 122: 123: 124: 125: 126: 127: 128: 129: 130: 131: 132: 133: 134: 135: 136: 137: 138: 139: 140: 141: 142: 143: 144: 145: 146: 147: 148: 149: 150: 151: 152: 153: 154: 155: 156: 157: 158: 159: 160: 161: 162: 163: 164: 165:
| #INCLUDE "windows.bi" #INCLUDE "vbcompat.bi"
DIM AS INTEGER i,j,C,COUNT,XX,TT,A,PP,QZ,U,RR ,CON , n ,maxP , STC DIM AS LONGINT Start,ST,STT DIM AS STRING FF,FF1,FF2,FF3 DIM AS BYTE zaehler FF=CURDIR FF2=FF+"\primes.txt" FF3=FF+"\QuadR.txt" //Quelldatei der 4Tupel 101 bis 9699690 , welche ab 23 erst Teiler haben, sind 36855 Stück
REDIM Q(36855) AS INTEGER REDIM QC(36855) AS INTEGER INPUT " - 10^n+.... - Exponent n";n INPUT "Start Milliarden";Start:Start=Start*1000000000 INPUT "Stop Milliarde ";ST:St=St*1000000000 INPUT "Primzahlen ins Sieb ab 100000 bis 80 M";maxP maxP=INT (maxP/log(maxP)) // Abschätzung der Anzahl für maxP RR=1 FOR J=1 TO n RR=(RR*10) MOD 96996900 // 10 ^ n MOD 96996900 = 10 Block a 9699690 wegen Rasterversatz NEXT J Start= (Start \ 96996900) * 96996900 - RR <- siehe hier OPEN FF3 FOR INPUT AS #1 //einlesen der 4Tupel FOR C=1 TO 36855 INPUT #1,A:Q(C)=A+8 NEXT C CLOSE #1
REDIM P(maxP) AS INTEGER // einlesen der Primteiler bis maxP OPEN FF2 FOR INPUT AS #2 FOR C=1 TO maxP INPUT #2,A P(C)=A NEXT C CLOSE #2 REM Erstellung der RESTTABELLE REDIM R(maxP) AS INTEGER FOR I=1 TO maxP // ermitteln der Resttabelle aller Primteiler für 10^n+Startwert, für Verschiebungsstart des Blockes nötig RR=1 FOR J=1 TO n RR=(RR*10) MOD P(I) NEXT J R(i)=(RR+Start) MOD P(I) // Der Restwert für 10^n+Start MOD Primteiler NEXT I REDIM R2(maxP) AS INTEGER // ersten 3000 Teiler Restwert für 96996900 Primteiler für die Unterblockverschiebung 9699690 FOR I=1 TO 3000 R2(I)=9699690 MOD P(I) NEXT I FOR I=3001 TO maxP // ab 3001 Teiler Restwert für 96996900 Primteiler für die Unterblockverschiebung 96996900 R2(I)=96996900 MOD P(I) NEXT I REM ENDE REM Erstellung der RESTTABELLE 2 - 9699690 MOD P REM BEGIN SIEBEN CLS B0: STC=INT(Start/10^11) // Schnitt volle 100 Mrd Blöcke für paralelles PRPing mit PFGW.exe FF1=FF+"\"+STR(n)+"."+STR(STC) // Datei PFGW Style erstellung
IF FILEEXISTS(FF1) THEN OPEN FF1 FOR APPEND AS #9 ELSE OPEN FF1 FOR OUTPUT AS #9 PRINT #9,"ABC 10^";n;"+$a & 10^";n;"+$a+2 & 10^";n;"+$a+6 & 10^";n;"+$a+8" BEGIN: LOCATE 1,1:PRINT "10^";n;"+";Start REDIM F(96996900) AS BYTE
FOR zaehler=0 TO 9 // 10 mal 9699690 Block FOR I=1 TO 36855 //Konstantes Feld Q - Feld für jeden Blockdurchlauf in QC kopiert, weil im QC Feld hier reduziert wird QC(I)=Q(I) NEXT I
count=36855 FOR i=1 TO 3000 U=1 XX=R(I) PP=P(i) WHILE u<count // und das hier ist das Meisterstück ! REST (Primzahl bzgl Block ) MOD Primzahl soll hier 0,2,6, oder 8 sein, damit der 4er ungültig wird bei einer Bedingung QZ=(QC(u) + XX ) MOD PP IF QZ>8 THEN GOTO WEITER IF QZ MOD 2 = 1 OR QZ=4 THEN GOTO WEITER // <10 und ungerade fallen auch raus QC(u)=QC(COUNT):COUNT-=1:u-=1 // Primärsieb ! Hole das letzte 4Tupel ins erste, wenn das 1. rausfällt, reduziere Anzahl der 4er im Block um 1 WEITER: u+=1 WEND R(I)=(R(I)+R2(I)) MOD PP // Errechne neuen Rest bei Blockverschiebung 9699690 NEXT i
// kopiere die übrigen in ein BYTE Feld der Größe 96996900 für weiteres Sieben, weil hier erstmals die 96996900 / Primteiler die Anzahl der Restkandidaten unterschreitet F() ist das neue Siebfeld, der Endblock als Raster
FOR j=1 TO count F(QC(j)+zaehler*9699690)=1 F(QC(j)-2+zaehler*9699690)=1 F(QC(j)-6+zaehler*9699690)=1 F(QC(j)-8+zaehler*9699690)=1 NEXT j
NEXT zaehler
// Siebe alle Primteiler ab der 3001. bis maxP wie gewohnt durch FOR i=3001 TO maxP U=1 RR=P(i)-R(i) IF RR MOD 2= 0 THEN RR+=P(i) FOR j=RR TO 96996900 STEP P(i)*2 F(j)=0 NEXT j R(I)=(R(I)+R2(I)) MOD P(i) // Auch hier die Restblock bzgl Primteiler bei Verschiebung Block 96996900 NEXT i // Ausgabe aller übrigen 4er ab 41 in 30er Block aussuchen und in Datei fürs PRPing schreiben FOR i=41 TO 96996900 STEP 30 IF F(i)=1 AND F(i+2)=1 AND F(i+6)=1 AND F(i+8 )=1 THEN PRINT #9,I+Start NEXT i STC=INT(Start/10^11) Start+=96996900 IF START>ST THEN CLOSE #9:STOP // Prüft ob Datei abgeschlossen wird bei volle 100 Mrd IF INT(Start/10^11)>STC THEN CLOSE #9:GOTO B0 GOTO BEGIN |
Zusatz: Die Zahl 36855 ist die Anzahl aller 4 Tupel von 41+30k bis 96996900 , die Teiler ab 23 erst haben.
Diese werden aus einer fertigen Datei eingelesen.
Primzahlen bis 80 Mio werden auch eingelesen.
--------------
Baller gerade an 665 herum. Nach 8 Stunden 10^665+ ( 2 000 000 000 000 ) fast abgeschlossen, ohne Fund. Denke heute Abend wirds bei 2-5 Billionen was sein.
Moderiert von Th69: Code-Tags hinzugefügt
Zuletzt bearbeitet von pzktupel am Mo 04.09.17 16:00, insgesamt 1-mal bearbeitet
Für diesen Beitrag haben gedankt: cryptnex, Gammatester
|
|
|