Autor |
Beitrag |
Th69
Beiträge: 4784
Erhaltene Danke: 1055
Win10
C#, C++ (VS 2017/19/22)
|
Verfasst: Di 05.09.17 16:50
Zu deinem Edit: man sollte diese Befehle allerdings nicht benutzen: Suspending Thread Execution
MSDN hat folgendes geschrieben: | The SuspendThread function is not intended to be used for thread synchronization because it does not control the point in the code at which the thread's execution is suspended. This function is primarily designed for use by debuggers. |
Für diesen Beitrag haben gedankt: Delphi-Laie
|
|
Gammatester
Beiträge: 328
Erhaltene Danke: 101
|
Verfasst: Di 05.09.17 20:55
Zur Abwechselung mal wieder etwas zu Quadrupeln:
Mit der Pinch-Liste der Carmichael-Zahlen bis 10^16 (Link: ftp.gwdg.de/pub/math...ael/carmichael-16.gz , 4066449 Bytes) habe ich folgende Tabelle berechnet. Sie enthält Fermat-Pseudoprimzahl-Quadrupel mit einer Carmichael-Zahl. (Hinweis: Zitat von de.wikipedia.org/wiki/Carmichael-Zahl Eine natürliche Zahl heißt Carmichael-Zahl, benannt nach dem Mathematiker Robert Daniel Carmichael, wenn sie eine fermatsche Pseudoprimzahl bezüglich aller zu ihr teilerfremden Basen ist.)
Interessanterweise tritt nur die Konstellation x-8, x-6, x-2, x auf, wobei x die Carmichael-Zahl ist und x-8, x-6, x-2 Fermat-Pseudoprimzahlen zur Basis 2.
Quelltext 1: 2: 3: 4: 5: 6: 7: 8: 9: 10: 11: 12: 13: 14: 15: 16: 17: 18: 19: 20:
| x Faktorisierung 1033669 7 13 37 307 8719309 19 37 79 157 9558334369 67 109 199 6577 107023329865 5 19 23 103 307 1549 317868427009 7 19 109 769 28513 1826363517265 5 7 53 79 157 163 487 2687463022465 5 13 17 19 67 1153 1657 3718714500865 5 7 13 19 397 769 1409 144692588218465 5 7 23 89 439 1237 3719 448933306317769 7 13 607 1213 6700249 637070091543865 5 7 107 397 8963 47807 842033901223945 5 53 109 313 937 99397 859152343812049 7 283 2179 7129 27919 1114610512354465 5 7 17 823 35117 64817 1206125685237889 13 59 73 929 1033 22447 2190392507242369 19 37 109 157 199 433 2113 4921398310999105 5 13 17 19 67 89 193 353 577 7248575534632705 5 7 13 17 37 97 433 523 1153 7459104727653085 5 13 37 67 199 277 829 1013 |
Ich habe schon auf math.stackexchange.com/questions/2417565/ eine entprechende Frage zu dieser Beobachtung gestellt, aber noch keine Antwort erhalten.
Für das Prim-Vierlings-Projekt haben diese Zahlen keine direkte Relevanz, die obigen Carmichael-Zahlen würden durch das Sieben ausgeschlossen werden. Es zeigt aber, daß ein Basis-2-Fermattest nicht ausreichend ist (wahrscheinlich gilt ähnliches für Basis-3-Fermattest, daß habe ich aber noch nicht durchgerechnet).
Für diesen Beitrag haben gedankt: cryptnex, Horst_H, Mathematiker
|
|
Delphi-Laie
Beiträge: 1600
Erhaltene Danke: 232
Delphi 2 - RAD-Studio 10.1 Berlin
|
Verfasst: Di 05.09.17 21:01
Th69 hat folgendes geschrieben : | Zu deinem Edit: man sollte diese Befehle allerdings nicht benutzen: Suspending Thread Execution
MSDN hat folgendes geschrieben: | The SuspendThread function is not intended to be used for thread synchronization because it does not control the point in the code at which the thread's execution is suspended. This function is primarily designed for use by debuggers. |
|
Danke! Tue ich ohnehin nicht, weil ich die beiden Befehle bis heute nicht verstehe. Ich wollte sie nur der Vollständigkeit halber als Edit erwähnen, weil ich vorher noch falsch behauptete, daß man in die Threadablaufsteuerung kaum eingreifen kann.
|
|
Horst_H
Beiträge: 1653
Erhaltene Danke: 243
WIN10,PuppyLinux
FreePascal,Lazarus
|
Verfasst: Mi 06.09.17 08:00
Hallo,
@ Gammatester
Ich hätte nicht ein so große Häufigkeit erwartet.( 6 von 1e15..1e16 == 6 aus 9e15 Zahlen )
Letztendlich ist der Fermattest, wie dass sieben, nur eine schnelle Vorauswahl, um die richtig aufwendigen Tests darauf los zu lassen.
Jetzt wäre doch die Frage, ob ein Kandidat, der sieben und Fermat überstanden hat, bei einem starkes Primtest durchgefallen ist.Das weiß sicher Mathematiker, dass wäre ja gefühlt eine Mini-Nadel im Heuhaufen.
Gruß Horst
|
|
OlafSt
Beiträge: 486
Erhaltene Danke: 99
Win7, Win81, Win10
Tokyo, VS2017
|
Verfasst: Mi 06.09.17 11:08
Delphi-Laie hat folgendes geschrieben : | Danke! Tue ich ohnehin nicht, weil ich die beiden Befehle bis heute nicht verstehe. Ich wollte sie nur der Vollständigkeit halber als Edit erwähnen, weil ich vorher noch falsch behauptete, daß man in die Threadablaufsteuerung kaum eingreifen kann. |
Nun, Suspend hält den Thread an. Der Thread bleibt augenblicklich stehen, egal was er auch immer gerade tun mag. Resume läßt einen mit Suspend angehaltenen Thread weiterlaufen.
Das Problem an Suspend ist, das der Thread irgendwo anhält. Das kann auch genauso gut mitten in einem inc(MyLongInt); sein, wo das Lo-Word bereits um eins hochgezählt ist, das Hi-Word aber noch nicht bearbeitet wurde. Man hat also überhaupt keine Kontrolle darüber, wie und wo der Thread zum Stillstand kommt.
Darum wurden Suspend und Resume auch mit deprecated markiert, so das man seine eigenen Anhalte- und Weiterlauf-Funktionen implementiert.
Dies nur als Hintergrund-Infos. Multithreading ist gar nicht so kompliziert, wie immer gesagt wird Man muß nur ein paar extra-Regeln beachten. Oh, und nachdem ich gestern Delphi Tokyo installieren und benutzen durfte, ist die Multithread-Unterstützung dort inzwischen ähnlich gut wie in C#. Wenn ich auch nur im Ansatz verstehen würde, was ihr da macht, würde ich da einsteigen und mal was mehrkerniges basteln Einfach nur um zu sehen, ob es was bringt.
_________________ Lies, was da steht. Denk dann drüber nach. Dann erst fragen.
Für diesen Beitrag haben gedankt: Delphi-Laie
|
|
Delphi-Laie
Beiträge: 1600
Erhaltene Danke: 232
Delphi 2 - RAD-Studio 10.1 Berlin
|
Verfasst: Mi 06.09.17 11:25
OlafSt hat folgendes geschrieben : | Delphi-Laie hat folgendes geschrieben : | Danke! Tue ich ohnehin nicht, weil ich die beiden Befehle bis heute nicht verstehe. Ich wollte sie nur der Vollständigkeit halber als Edit erwähnen, weil ich vorher noch falsch behauptete, daß man in die Threadablaufsteuerung kaum eingreifen kann. |
Nun, Suspend hält den Thread an. Der Thread bleibt augenblicklich stehen, egal was er auch immer gerade tun mag. Resume läßt einen mit Suspend angehaltenen Thread weiterlaufen. |
Soweit war ich auch schon, konnte das aber nie bis zu funktionierendem, gar produktivem Code weiterentwickeln. Was ich am wenigsten verstand: Wie soll sich ein angehaltener Thread "selbst am Schopfe aus dem Sumpfe ziehend" mit Suspend ein Weiterlaufen verpassen? Muß wohl von außerhalb passieren.
OlafSt hat folgendes geschrieben : | Multithreading ist gar nicht so kompliziert, wie immer gesagt wird Man muß nur ein paar extra-Regeln beachten. |
Naja, man muß sich schon damit beschäftigen, hinterher ist man immer schlauer, und was man weiß, ist einfach (alles Binsensweisheiten).
OlafSt hat folgendes geschrieben : | Wenn ich auch nur im Ansatz verstehen würde, was ihr da macht, würde ich da einsteigen und mal was mehrkerniges basteln Einfach nur um zu sehen, ob es was bringt. |
Nun, da gibt es eine einfache Faustregel, auf die ich selbst stieß, die aber eigentlich auch offensichtlich ist: Dinge, deren Abarbeitungsreihenfolge vertauscht werden kann oder gar egal ist (im Notfalle sogar ausprobierbar), können prinzipiell auch gleichzeitig bzw. simultan geschehen.
Einfaches Beispiel aus der Mathematik, eher noch dem Rechnen der unteren Schulklassen; erinnern wir uns daran, wie das abläuft:
- Bei der schriftlichen Addition, Subtraktion und Division (die eigentlich auch nur eine "gestaffelte" Subtraktion ist), ist die Abarbeitungsreihenfolge festgelegt (bei den ersteren von hinten / klein nach vorn / groß, bei letzterem umgekehrt), weil die Zwischenergebnisse benötigt und benutzt, also weiterverarbeitet werden. Sieht eher nicht nach Parallelisierbarkeit aus.
- Bei der schriftlichen Multiplikation hingegen sind die Teilmultiplikationen völlig unabhängig voneinander, keine Zwischenergebnisse werden für eine andere Teilmultiplikation benötigt. Also können diese Schritte auch parallisiert werden. Auch die Addition der erhaltenen Summanden kann "balanciert", also parallel erfolgen (wobei aber jede diese Addition intern dann wieder nichtparallel erfolgt).
Zuletzt bearbeitet von Delphi-Laie am Mi 06.09.17 21:08, insgesamt 2-mal bearbeitet
|
|
OlafSt
Beiträge: 486
Erhaltene Danke: 99
Win7, Win81, Win10
Tokyo, VS2017
|
Verfasst: Mi 06.09.17 14:42
Delphi-Laie hat folgendes geschrieben : |
Soweit war ich auch schon, konnte das aber nie bis zu funktionierendem, gar produktiven Code weiterentwickeln. Was ich am wenigsten verstand: Wie soll sich ein angehaltener Thread "selbst am Schopfe aus dem Sumpfe ziehend" mit Suspend ein Weiterlaufen verpassen? Muß wohl von außerhalb passieren.
|
Absolut richtig erkannt. Ein mit Suspend "eingeschläferter" Thread kann nur per Resume wieder zum Leben erweckt werden - oder eben gekillt werden. Durch ein Suspend bekommt der Thread keinerlei Rechenzeit mehr seitens des Betriebssystems zugewiesen, er kann sich also auch nicht "selbst an den Haaren aus dem Sumpf ziehen". Somit kann ein Suspended Thread auch nicht korrekt "aufgeräumt" werden, wenn das erstellende Programm terminiert.
_________________ Lies, was da steht. Denk dann drüber nach. Dann erst fragen.
Für diesen Beitrag haben gedankt: Delphi-Laie
|
|
Gammatester
Beiträge: 328
Erhaltene Danke: 101
|
Verfasst: Do 07.09.17 17:36
Hier neue Ergebnisse zum diesem Beitrag:
Nachdem ich für mehrere Basen gerechnet und keine neuen Quadrupel mit Carmichael-Zahlen befunden habe, stellte sich heraus, daß die x-8, x-6, x-2 alles Primzahlen sind (im Sinne von mp_is_pprime). Die Tabelle für die Basis b erhält man, wenn man aus Liste alle Zahlen streicht, die durch b teilbar sind.
Es bleibt weiterhin die offene Frage, warum die Carmichael-Zahl immer am Ende eines Quadrupels steht.
Gruß Gammatester
|
|
pzktupel
Hält's aus hier
Beiträge: 129
Erhaltene Danke: 30
|
Verfasst: Di 12.09.17 23:02
Gammatester hat folgendes geschrieben : | Zur Abwechselung mal wieder etwas zu Quadrupeln:
Mit der Pinch-Liste der Carmichael-Zahlen bis 10^16 (Link: ftp.gwdg.de/pub/math...ael/carmichael-16.gz , 4066449 Bytes) habe ich folgende Tabelle berechnet. Sie enthält Fermat-Pseudoprimzahl-Quadrupel mit einer Carmichael-Zahl. (Hinweis: Zitat von de.wikipedia.org/wiki/Carmichael-Zahl Eine natürliche Zahl heißt Carmichael-Zahl, benannt nach dem Mathematiker Robert Daniel Carmichael, wenn sie eine fermatsche Pseudoprimzahl bezüglich aller zu ihr teilerfremden Basen ist.)
Interessanterweise tritt nur die Konstellation x-8, x-6, x-2, x auf, wobei x die Carmichael-Zahl ist und x-8, x-6, x-2 Fermat-Pseudoprimzahlen zur Basis 2.
Quelltext 1: 2: 3: 4: 5: 6: 7: 8: 9: 10: 11: 12: 13: 14: 15: 16: 17: 18: 19: 20:
| x Faktorisierung 1033669 7 13 37 307 8719309 19 37 79 157 9558334369 67 109 199 6577 107023329865 5 19 23 103 307 1549 317868427009 7 19 109 769 28513 1826363517265 5 7 53 79 157 163 487 2687463022465 5 13 17 19 67 1153 1657 3718714500865 5 7 13 19 397 769 1409 144692588218465 5 7 23 89 439 1237 3719 448933306317769 7 13 607 1213 6700249 637070091543865 5 7 107 397 8963 47807 842033901223945 5 53 109 313 937 99397 859152343812049 7 283 2179 7129 27919 1114610512354465 5 7 17 823 35117 64817 1206125685237889 13 59 73 929 1033 22447 2190392507242369 19 37 109 157 199 433 2113 4921398310999105 5 13 17 19 67 89 193 353 577 7248575534632705 5 7 13 17 37 97 433 523 1153 7459104727653085 5 13 37 67 199 277 829 1013 |
Ich habe schon auf math.stackexchange.com/questions/2417565/ eine entprechende Frage zu dieser Beobachtung gestellt, aber noch keine Antwort erhalten.
Für das Prim-Vierlings-Projekt haben diese Zahlen keine direkte Relevanz, die obigen Carmichael-Zahlen würden durch das Sieben ausgeschlossen werden. Es zeigt aber, daß ein Basis-2-Fermattest nicht ausreichend ist (wahrscheinlich gilt ähnliches für Basis-3-Fermattest, daß habe ich aber noch nicht durchgerechnet). |
Interessant. Die obige Zahl erfüllt für alle Basen den Fermattest, jedoch fallen die beim verschärften Test durch und die Teiler liegen sofort auf der Hand.
(1033669-1)/4=258417
2^258147 MOD 1033669 = 276606. 276607=17x53x307 , 307 ist Teiler.
3^258147 MOD 1033669 = 238538. 238539=3x7x37x307 ggT(238539,1033669)=7x37x307
|
|
Gammatester
Beiträge: 328
Erhaltene Danke: 101
|
Verfasst: Di 12.09.17 23:51
|
|
pzktupel
Hält's aus hier
Beiträge: 129
Erhaltene Danke: 30
|
Verfasst: Mi 13.09.17 05:59
Moment, ich meine, das für alle Basen a von 2 bis 1000 erstmal erfüllt ist:
a^1033669=a MOD 1033669.
Quelle Derive: VECTOR(MOD(a^1033669,1033669),a,2,1000)
Ausgabe: [2,3,4,5,6,7,8.....]
|
|
Gammatester
Beiträge: 328
Erhaltene Danke: 101
|
Verfasst: Mi 13.09.17 10:12
pzktupel hat folgendes geschrieben : | Moment, ich meine, das für alle Basen a von 2 bis 1000 erstmal erfüllt ist:
a^1033669=a MOD 1033669.
Quelle Derive: VECTOR(MOD(a^1033669,1033669),a,2,1000)
Ausgabe: [2,3,4,5,6,7,8.....] |
Das ist trivial für Carmichael-Zahlen, hier gilt ja a^(n-1) = 1 mod n ( de.wikipedia.org/wik...hael-Zahl#Definition ) und damit a^n = a mod n.
Normalerweise ist n eine starke Pseudo-Primzahl zur Basis a, wenn mit n=d*2^s +1 gilt (vgl. de.wikipedia.org/wik...tarke_Pseudoprimzahl )
1. a^d = 1 mod, oder
2. a^(d*2^r) = -1 mod n für ein r mit 0 <= r < s
Ich weiß jetzt nicht, was Du unter 'verschärften Test' verstehst. Ich hatte angenommen, Du meinst einen starken Pseudo-Primzahltest.
Wie Du schon geschrieben hast, hat man zB 2^258417 mod 1033669 = 276606 <> 1, und da 2^(2^1*258417) mod 1033669 = 986049 <> -1 gilt, ist 1033669 keine starke Pseudoprimzahl zur Basis 2. Aber wegen 9^258417 mod 1033669 = 1 ist es eine starke Pseudoprimzahl zur Basis 9.
|
|
pzktupel
Hält's aus hier
Beiträge: 129
Erhaltene Danke: 30
|
Verfasst: Mi 13.09.17 11:24
Für alle Interessenten.
Habe eben das kleinste 100-stellige Primzahl - Quadrupel - Pärchen mit Abstand 30 gefunden.
Es lautet:
10^99+5'294'137'569'927'811 +d,d=0,2,6,8,30,32,36,38
LG
Norman
|
|
Horst_H
Beiträge: 1653
Erhaltene Danke: 243
WIN10,PuppyLinux
FreePascal,Lazarus
|
Verfasst: Do 14.09.17 21:10
Hallo,
anbei meine neueste Version, die gmp nutzt.Schafft einen Bereich von etwa 310e6 Zahlen in einer Stunde.
954 lag bei 7,5e12.Das dauert etwas ( 26 Stunden mit einer älteren Version 288e6/h )
Da ich mit Lazarus unter Linux 64 Bit arbeite, brauche ich keine gmp.dll für Windows.
Dafür habe ich nur eine alte 32-Bit Version von 2014, was sehr langsam war.
Wenn jemand eine sehr neue, am liebsten Win64 Version zur Verfügung stellen könnte, wäre das sicher hilfreich,
Gruß Horst
Einloggen, um Attachments anzusehen!
Für diesen Beitrag haben gedankt: Mathematiker
|
|
pzktupel
Hält's aus hier
Beiträge: 129
Erhaltene Danke: 30
|
Verfasst: Do 14.09.17 23:21
Horst_H hat folgendes geschrieben : | Hallo,
anbei meine neueste Version, die gmp nutzt.Schafft einen Bereich von etwa 310e6 Zahlen in einer Stunde.
954 lag bei 7,5e12.Das dauert etwas ( 26 Stunden mit einer älteren Version 288e6/h )
Da ich mit Lazarus unter Linux 64 Bit arbeite, brauche ich keine gmp.dll für Windows.
Dafür habe ich nur eine alte 32-Bit Version von 2014, was sehr langsam war.
Wenn jemand eine sehr neue, am liebsten Win64 Version zur Verfügung stellen könnte, wäre das sicher hilfreich,
Gruß Horst |
Horst , du meinst doch immer e9,oder ?
|
|
Horst_H
Beiträge: 1653
Erhaltene Danke: 243
WIN10,PuppyLinux
FreePascal,Lazarus
|
Verfasst: Fr 15.09.17 06:58
Hallo,
ja natürlich e9 bei Stunden, das wäre ja sonst zu traurig.
7.5e12 in 26 Stunden sind ja 288e9
Gruß Horst
|
|
Horst_H
Beiträge: 1653
Erhaltene Danke: 243
WIN10,PuppyLinux
FreePascal,Lazarus
|
Verfasst: Di 19.09.17 16:47
Hallo,
ich habe mal eine Version ohne Threads und mit etwas verbessertem Sieb.
Ich passe jetzt die Anzahl der verwendeten Siebprimzahlen an.
Zudem nutze ich jetzt Siebnr für große Siebprimzahlen.2e9 streicht im besten Fall nach 70 Sieben und dann nach 6.
Das bedeutet, eine kleine Ersparnis.
620Mio werden bei 998 in 773 ms gesiebt, also ein Bereich 0,8e9/s
bei 430 Stellen sind es schon 620e6/0,5 = 1,24e9/s.
Dort werden die Kandidaten auch 10-mal schneller getestet.
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:
| Siebgroesse 2586584 Byte // jedes Bit entspricht einem k aus 30*k+11-> 620.780.160 Start-Potenz 998 Start-Offset 0 //Zu Beginne werden alle Siebprimzahlen "ersiebt" Anzahl Siebprimzahlen 98222284 bis 1999999973 //Die Anzahl der verwendeten Siebprimzahlen Maximaler SiebprimIndex: 97930762 //inclusive Primzahlen ersieben und Übertrag berechnen für SiebprimIndex Siebzeit 00:00:21.674 Testzeit 00:00:18.886 Testzeit pro Kandidat 15.218 ms Siebzeit 00:00:00.773 Testzeit 00:00:17.744 Testzeit pro Kandidat 14.961 ms Siebzeit 00:00:00.771 Testzeit 00:00:16.967 Testzeit pro Kandidat 14.589 ms Siebzeit 00:00:00.795 Testzeit 00:00:17.689 Testzeit pro Kandidat 14.231 ms //Loesung gefunden, deshalb etwas weniger Zeit 998 1970797081 00:01:35.309 Gesamtzeit: 00:01:35.310 ---------- Start-Potenz 598 Start-Offset 0 Maximaler SiebprimIndex: 45959098 Siebzeit 00:00:04.617 //hier nur inclusive Übertrag berechnen für SiebprimIndex Testzeit 00:00:04.902 Testzeit pro Kandidat 3.494 ms Siebzeit 00:00:00.588 Testzeit 00:00:04.888 Testzeit pro Kandidat 3.409 ms Siebzeit 00:00:00.596 Testzeit 00:00:04.954 Testzeit pro Kandidat 3.452 ms Siebzeit 00:00:00.589 Testzeit 00:00:05.056 Testzeit pro Kandidat 3.451 ms Siebzeit 00:00:00.587 Testzeit 00:00:04.622 Testzeit pro Kandidat 3.470 ms Siebzeit 00:00:00.585 Testzeit 00:00:04.791 Testzeit pro Kandidat 3.417 ms Siebzeit 00:00:00.583 Testzeit 00:00:04.730 Testzeit pro Kandidat 3.433 ms Siebzeit 00:00:00.583 Testzeit 00:00:04.923 Testzeit pro Kandidat 3.435 ms Siebzeit 00:00:00.585 Testzeit 00:00:04.949 Testzeit pro Kandidat 3.420 ms Siebzeit 00:00:00.584 Testzeit 00:00:04.952 Testzeit pro Kandidat 3.441 ms Siebzeit 00:00:00.584 Testzeit 00:00:00.050 Testzeit pro Kandidat 0.035 ms //Loesung gefunden, deshalb weniger Zeit 598 6011086471 00:00:59.367
Start-Potenz 430 Start-Offset 0 Maximaler SiebprimIndex: 28413753 Siebzeit 00:00:02.713 Testzeit 00:00:02.320 Testzeit pro Kandidat 1.456 ms Siebzeit 00:00:00.504 Testzeit 00:00:02.337 Testzeit pro Kandidat 1.457 ms Siebzeit 00:00:00.505 Testzeit 00:00:02.189 Testzeit pro Kandidat 1.416 ms Siebzeit 00:00:00.504 Testzeit 00:00:02.340 Testzeit pro Kandidat 1.485 ms Siebzeit 00:00:00.503 Testzeit 00:00:02.184 Testzeit pro Kandidat 1.418 ms Siebzeit 00:00:00.502 Testzeit 00:00:02.262 Testzeit pro Kandidat 1.468 ms Siebzeit 00:00:00.503 Testzeit 00:00:02.274 Testzeit pro Kandidat 1.421 ms Siebzeit 00:00:00.513 Testzeit 00:00:02.339 Testzeit pro Kandidat 1.480 ms Siebzeit 00:00:00.503 Testzeit 00:00:02.206 Testzeit pro Kandidat 1.424 ms Siebzeit 00:00:00.502 Testzeit 00:00:02.370 Testzeit pro Kandidat 1.470 ms Siebzeit 00:00:00.504 Testzeit 00:00:02.140 Testzeit pro Kandidat 1.427 ms Siebzeit 00:00:00.503 Testzeit 00:00:01.572 Testzeit pro Kandidat 0.994 ms //Loesung gefunden, deshalb weniger Zeit 430 6969624301 00:00:34.836 |
Um alle Siebprimzahlen einmal nur abzuklappern, ohne was anderes zu tun ausser p-> p +deltaP ,dauert schon 0,1 Sekunden für die 98.222.284 Primzahlen.Irgendwie lahm.
Gruß Horst
Edit:
Unter wine als Win32 Programm braucht das Programm bei n= 998 für die ersten Kandidaten für den Test eines Kandidaten sage und schreibe 263 ms statt 15,2 ms unter Linux64 mit neuester gmp.so.
32-Bit siebt ein wenig schneller.Das hilft nur nichts.
Einloggen, um Attachments anzusehen!
Für diesen Beitrag haben gedankt: Mathematiker
|
|
Horst_H
Beiträge: 1653
Erhaltene Danke: 243
WIN10,PuppyLinux
FreePascal,Lazarus
|
Verfasst: Fr 03.11.17 10:46
Hallo,
das neueste Kommandozeilenmonster, das nur sie möglichen Quadrupel ersiebt und dann mittels pfgw64 aus sourceforge.net/projects/openpfgw/files testet.
Es wird kein gmp genutzt, um die Startposition im Quadsieb zu berechnen, was unter Linux doppelt so schnell wäre, aber die alte Windows gmp.dll ist viel langsamer, als reines Pascal mit klitzekleinem ASM 64 Bit..
Es werden die Siebprimzahlen immer wieder neu erzeugt und die Startposition im QuadrupelSieb bestimmt.Bei einem Limit von 50e9 müsste man mehr als 2e9 speichern.Geht nicht.
Funktioniert nur noch für 64-Bit.
//Man muss entsprechend den zu testenden Exponenten cZahlBereich und cMaxQuadSiebPrim anpassen.
Wird jetzt automatisch angepasst, wie man an den Laufzeiten in Funde.txt erkennt.
Unter Wine64 läuft das Programm auch, man muss nur pfgw64 im selben Verzeichnis haben ( reicht ja der Link/Verknüpfung )
Aufruf mit Parametern Stellenzahl/Exponent und gegebenfalls Offset.
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: 313: 314: 315: 316: 317: 318: 319: 320: 321: 322: 323: 324: 325: 326: 327: 328: 329: 330: 331: 332: 333: 334: 335: 336: 337: 338: 339: 340: 341: 342: 343: 344: 345: 346: 347: 348: 349: 350: 351: 352: 353: 354: 355: 356: 357: 358: 359: 360: 361: 362: 363: 364: 365: 366: 367: 368: 369: 370: 371: 372: 373: 374: 375: 376: 377: 378: 379: 380: 381: 382: 383: 384: 385: 386: 387: 388: 389: 390: 391: 392: 393: 394: 395: 396: 397: 398: 399: 400: 401: 402: 403: 404: 405: 406: 407: 408: 409: 410: 411: 412: 413: 414: 415: 416: 417: 418: 419: 420: 421: 422: 423: 424: 425: 426: 427: 428: 429: 430: 431: 432: 433: 434: 435: 436: 437: 438: 439: 440: 441: 442: 443: 444: 445: 446: 447: 448: 449: 450: 451: 452: 453: 454: 455: 456: 457: 458: 459: 460: 461: 462: 463: 464: 465: 466: 467: 468: 469: 470: 471: 472: 473: 474: 475: 476: 477: 478: 479: 480: 481: 482: 483: 484: 485: 486: 487: 488: 489: 490: 491: 492: 493: 494: 495: 496: 497: 498: 499: 500: 501: 502: 503: 504: 505: 506: 507: 508: 509: 510: 511: 512: 513: 514: 515: 516: 517: 518: 519: 520: 521: 522: 523: 524: 525: 526: 527: 528: 529: 530: 531: 532: 533: 534: 535: 536: 537: 538: 539: 540: 541: 542: 543: 544: 545: 546: 547: 548: 549: 550: 551: 552: 553: 554: 555: 556: 557: 558: 559: 560: 561: 562: 563: 564: 565: 566: 567: 568: 569: 570: 571: 572: 573: 574: 575: 576: 577: 578: 579: 580: 581: 582: 583: 584: 585: 586: 587: 588: 589: 590: 591: 592: 593: 594: 595: 596: 597: 598: 599: 600: 601: 602: 603: 604: 605: 606: 607: 608: 609: 610: 611: 612: 613: 614: 615: 616: 617: 618: 619: 620: 621: 622: 623: 624: 625: 626: 627: 628: 629: 630: 631: 632: 633: 634: 635: 636: 637: 638: 639: 640: 641: 642: 643: 644: 645: 646: 647: 648: 649: 650: 651: 652: 653: 654: 655: 656: 657: 658: 659: 660: 661: 662: 663: 664: 665: 666: 667: 668: 669: 670: 671: 672: 673: 674: 675: 676: 677: 678: 679: 680: 681: 682: 683: 684: 685: 686: 687: 688: 689: 690: 691: 692: 693: 694: 695: 696: 697: 698: 699: 700: 701: 702: 703: 704: 705: 706: 707: 708: 709: 710: 711: 712: 713: 714: 715: 716: 717: 718: 719: 720: 721: 722: 723: 724: 725: 726: 727: 728: 729: 730: 731: 732: 733: 734: 735: 736: 737: 738: 739: 740: 741: 742: 743: 744: 745: 746: 747: 748: 749: 750: 751: 752: 753: 754: 755: 756: 757: 758: 759: 760: 761: 762: 763: 764: 765: 766: 767: 768: 769: 770: 771: 772: 773: 774: 775: 776: 777: 778: 779: 780: 781: 782: 783: 784: 785: 786: 787: 788: 789: 790: 791: 792: 793: 794: 795: 796: 797: 798: 799: 800: 801: 802: 803: 804: 805: 806: 807: 808: 809: 810: 811: 812: 813: 814: 815: 816: 817: 818: 819: 820: 821: 822: 823: 824: 825: 826: 827: 828: 829: 830: 831: 832: 833: 834: 835: 836: 837: 838: 839: 840: 841: 842: 843: 844: 845: 846: 847: 848: 849: 850: 851: 852: 853: 854: 855: 856: 857: 858: 859: 860: 861: 862: 863: 864: 865: 866: 867: 868: 869: 870: 871: 872: 873: 874: 875: 876: 877: 878: 879: 880: 881: 882: 883: 884: 885: 886: 887: 888: 889: 890: 891: 892: 893: 894: 895: 896: 897: 898: 899: 900: 901: 902: 903: 904: 905: 906: 907: 908: 909: 910: 911: 912: 913: 914: 915: 916: 917: 918: 919: 920: 921: 922: 923: 924: 925: 926: 927: 928: 929: 930: 931: 932: 933: 934: 935: 936: 937: 938: 939: 940: 941: 942: 943: 944: 945: 946: 947: 948: 949: 950: 951: 952: 953: 954: 955: 956: 957: 958: 959: 960: 961: 962: 963: 964: 965: 966: 967: 968: 969: 970: 971: 972: 973: 974: 975: 976: 977: 978: 979: 980: 981: 982: 983: 984: 985: 986: 987: 988: 989: 990: 991: 992: 993: 994: 995: 996: 997: 998: 999: 1000: 1001: 1002: 1003: 1004: 1005: 1006: 1007: 1008: 1009: 1010: 1011: 1012: 1013: 1014: 1015: 1016: 1017: 1018: 1019: 1020: 1021: 1022: 1023: 1024: 1025: 1026: 1027: 1028: 1029: 1030: 1031: 1032: 1033: 1034: 1035: 1036: 1037: 1038: 1039: 1040: 1041: 1042:
| program QuadSieb; {$IFDEF DELPHI} {$APPTYPE CONSOLE} {$ENDIF}
{$IFDEF FPC} {$MODE Delphi} {$OPTIMIZATION ON,PEEPHOLE,CSE} {$CODEALIGN proc=32} {$ALIGN 8} {$ENDIF}
{$DEFINE SUCHE} uses SysUtils,process;
const PWGW_Log_Name = 'pfgw.log';
cE6 = 1000*1000; cE9 = 1000*cE6;
cMAXQuadSiebBereich = 500*cE9; cMaxQuadSiebPrim = 100*cE9; cDeltaBitLimit_II = 8 * 2*1024*1024; cElemZahlBereich = 64*30; type tSieveElem = NativeUint; tpSieveElem = pNativeUint; tSiebMod30 = array of tSieveElem;
tDelta30_64 = array[0..3] of Uint64; tpDelta30_64 =^tDelta30_64;
tMulOfs = record moMul,moOfs: byte; end; tMulPosArr = array[0..3] of tMulOfs;
tSegmentedQuadSievePrime = record sqDelta :tDelta30_64; sqPrim :Uint64; sqOfsIdx:Uint64; end; tsqFeld = array[0..40030] of tSegmentedQuadSievePrime; tpSegQuadSvPrim = ^tSegmentedQuadSievePrime; const cElemBitCnt = SizeOf(tSieveElem)*8; cAndMask = cElemBitCnt-1;
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);
cPrimFeldLimit = 78498; Max = 2*3*5*7*11*13; StartPrimNr = 5; Max_2 = (Max+1) shr 1;
AnzahlWifeldinSiebFeld = 2; AnzahlProSieb = AnzahlWifeldinSiebFeld*Max; Max_2_2 = AnzahlProSieb shr 1;
type tPrimrec = Record Uebertrag :Uint64; Primzahl :Uint64; end; tWiFeld = array [0..Max_2-1] of boolean; tSieb = array [0..Max_2_2-1] of boolean;
tprimFeld = array [0..cPrimFeldLimit] of tPrimrec;
var Sieb : ^tSieb; WiFeld : ^tWiFeld; PrimFeld : ^tprimFeld; StartZeit : TDAteTime; WiFeLaenge : NativeUint; gblAktPrimPos, gblMaxPrimPos : NativeUInt;
OberesLimit, AnzahlSiebPrim, OberGrenzeAlt, OberGrenzeNeu, SiebStartZahl, SiebStartPrim : Uint64;
QuadSieb : tSiebMod30;
myTextBuf: array[0..4096*15] of byte;
PFGW_Surfilename, PFGW_Tstfilename :AnsiString; PFGW_Dat, PFGW_Log, Funde_Dat : textFile;
sqFeld : tsqFeld; RepLaenge, gblExpo, gblBase, gblExpo2, gblFact2 :Uint64;
gblSiebMod, gblSiebOfs, gblSiebOfsDiv10, gblStartOffset, gblCompPot : Int64;
gblPrimZahlLimit, gblSiebElemCnt , gblSiebBitCnt, gblSiebNumbCnt, gblQuadSiebBereich, gblMaxQuadSiebPrim : Uint64;
Time1:TDateTime; gblSiebNr, gblsqFeldIdx : longWord;
startCnt, endCnt, gblDatCnt, ChkLidx :LongWord;
gefunden,RepFertig : Boolean;
abbruch : boolean;
procedure PowModAsm64_Init(Expo:NativeInt); var i,DeziDgtCnt : NativeUint; begin DeziDgtCnt := trunc(ln(High(Uint64) DIV cMaxQuadSiebPrim)/ln(10)); i := 19; repeat gblExpo2 := Expo DIV i; gblFact2 := Expo MOD i; If gblFact2 <=DeziDgtCnt then break; dec(i); until i <= 9; gblBase := 1; For i := i downto 1 do gblBase := gblBase*10; i := gblFact2; gblFact2 := 1; For i := i downto 1 do gblFact2 := gblFact2*10; end;
function PowMod(base,exponent,divisor:Uint64):LongWord; begin result := 1; IF (divisor = 0) OR (exponent = 0) then EXIT; base := base MOD divisor; IF base = 0 then Begin result := 0; EXIT; end; repeat IF ODD(Exponent) then result := (result*base) MOD divisor; base := sqr(base) MOD divisor; Exponent:= Exponent shr 1; until Exponent = 0; end;
procedure DateiTesten; var PFGW_Process : TProcess; begin PFGW_Process := Tprocess.create(NIl); PFGW_Process.Executable := 'pfgw64'; PFGW_Process.options:=PFGW_Process.options+[poWaitOnExit]; PFGW_Process.Parameters.clear; PFGW_Process.Parameters.Add('-f0'); PFGW_Process.Parameters.Add('-k'); PFGW_Process.Parameters.Add('-Cquiet'); PFGW_Process.Parameters.Add(PFGW_Tstfilename); PFGW_Process.execute; PFGW_Process.free; end;
procedure NextQuadSieb; begin RepFertig := false; inc(gblDatCnt); inc(gblSiebOfs,gblSiebNumbCnt); inc(gblSiebOfsDiv10,gblSiebNumbCnt DIV 10); end;
procedure ChangeStartCnt(Cnt:NativeInt); var Factor : double; begin IF cnt < 2 then cnt := 2; gblExpo := Cnt-1; gblSiebMod := PowMod(10,gblExpo,cElemZahlBereich); gblSiebOfsDiv10 := gblSiebMod DIV 10; gblStartOffset := gblStartOffset-gblStartOffset MOD cElemZahlBereich; gblSiebOfsDiv10 := (gblSiebMod+gblStartOffset) DIV 10; gblSiebOfs := gblStartOffset+gblSiebMod+11-20;
Factor := Cnt/1050; IF Factor> 1 then Factor := 1; IF Factor < 0.7 then Begin gblQuadSiebBereich := trunc(exp((3.1-2*Factor)*ln(Factor))*cMAXQuadSiebBereich)+60060; gblMaxQuadSiebPrim := trunc(exp((3.7-Factor)*ln(Factor))*cMaxQuadSiebPrim); end else Begin gblQuadSiebBereich := cMAXQuadSiebBereich; gblMaxQuadSiebPrim := trunc(Factor*cMaxQuadSiebPrim); end;
if gblMaxQuadSiebPrim <60060 then gblMaxQuadSiebPrim := 60060; gblPrimZahlLimit := trunc(sqrt(gblMaxQuadSiebPrim))+1;
gblSiebNumbCnt:= (gblQuadSiebBereich DIV cElemZahlBereich+1)*cElemZahlBereich; gblSiebBitCnt := gblSiebNumbCnt DIV 30; gblSiebElemCnt := gblSiebBitCnt DIV cElemBitCnt; RepFertig := false; writeln(Format('Digits %d QuadSiebGroesze %d MaxSiebPrim %d', [cnt,gblQuadSiebBereich,gblMaxQuadSiebPrim])); end;
procedure KandidatenEinsammeln; var kandPos, MyOffset : Int64; Elem : tSieveElem; i : LongWord; s: String; Begin AssignFile(PFGW_Dat,PFGW_Tstfilename); ReWrite(PFGW_Dat); s := Format('10^%d+$a',[gblExpo]); s := Format('ABC %s&%s+2&%s+6&%s+8',[s,s,s,s]); writeln(PFGW_Dat,s); s:= ''; ChkLIdx := 0; MyOffset := gblSiebOfs; For i := 0 to gblSiebElemCnt-1 do Begin Elem := NOT QuadSieb[i]; if Elem <> 0 then Begin kandPos := MyOffset; repeat if (Elem AND 255)= 0 then Begin Elem := Elem shr 8; inc(kandPos,30*8); continue; end; IF Boolean(Elem AND 1) then Begin writeln(PFGW_DAT,kandpos); inc(ChkLidx); if (ChkLidx AND 4095)=0 then flush(PFGW_Dat); end; inc(kandPos,30); Elem:= Elem shr 1; until Elem = 0; end; inc(MyOffset,(30*cElemBitCnt)); end; CloseFile(PFGW_Dat); end;
function P30toP(p30:NativeUint):NativeUint; Begin result := (p30 shr 3)*30+cIdxToMod30[p30 AND 7]; end;
function KleinePrimListeErsieben:NativeInt; var EraSieb : array of byte; i,p: NativeInt; Begin setlength(EraSieb,gblPrimZahlLimit+210); p := 2; i := p*p; repeat while i <= High(EraSieb) do Begin EraSieb[i] := 1; inc(i,p); end; repeat inc(p); until EraSieb[p]=0; i:= p*p; until i > gblPrimZahlLimit; i := 0; For p := 2 to High(EraSieb) do IF EraSieb[p] = 0 then with PrimFeld^[i] do Begin Primzahl := p; Uebertrag:= (p*p)SHR 1; inc(i); end; dec(i); setlength(EraSieb,0); result:= i; end;
procedure Wifeld_Aufbauen; var j,J_2,k, AnzahlSiebPrim, PosInSieb : integer; begin WiFeld^[0] := true; WiFeld^[1] := true; WiFeLaenge := 1; PosInSieb := 1; AnzahlSiebPrim := 0; while AnzahlSiebPrim < StartPrimNr do Begin if WiFeld^[PosInSieb] then begin J_2 := succ(PosInSieb shl 1); for j := 1 to J_2-1 do move(WiFeld^[0],WiFeld^[WiFeLaenge*j],WiFeLaenge); K := PosInSieb; while k < Max_2 do begin WiFeld^[k]:= false; inc(k,j_2); end; WiFeLaenge := WiFeLaenge*J_2; inc(AnzahlSiebPrim); end; inc(PosInSieb); end; end;
procedure SiebFeldaufbauen; var i : integer; begin for i := 0 to AnzahlWifeldinSiebFeld-1 do move(Wifeld^[0],Sieb^[i*WiFeLaenge],WiFeLaenge); end;
procedure NaechstesKleinPrimSieb; var MomTestPrim:^tPrimrec; pS : ^tSieb; u,p,i : Int64;
begin SiebFeldaufbauen; pS := @Sieb^[0]; MomTestprim := @Primfeld^[StartPrimNr+1]; i := gblAktPrimPos-(StartPrimNr+1); for i := i downto 0 do begin with MomTestPrim^ do begin p := Primzahl; u := Uebertrag; end; while u < MAX_2_2 do begin pS^[u]:= false; inc(u,p); end; MomTestPrim^.Uebertrag := u-Max_2_2; inc(MomTestPrim); end; i := Max_2_2*gblSiebNr; IF gblAktPrimPos < gblMaxPrimPos then Begin while MomTestPrim^.Uebertrag < i do Begin with MomTestPrim^ do begin p := Primzahl; u := Uebertrag; u := (u+Max_2_2)-i; end; repeat pS^[u]:= false; inc(u,p); until u >= MAX_2_2; MomTestPrim^.Uebertrag := u-Max_2_2; inc(MomTestPrim); inc(gblAktPrimPos); end; end; inc(gblSiebNr); end;
procedure KleineQuadPrimSieben(p,sqIdx:LongInt); var pSieb :pByte; pDelta30 : tpDelta30_64; LastRep,tmp,Offset30:Int64; i,idx : NativeInt; Begin pSieb :=@QuadSieb[0]; LastRep := RepLaenge; RepLaenge := LastRep*p; For i := 1 to p do Begin IF (i+1)*LastRep <= gblSiebBitCnt then move(pSieb[0],pSieb[i*LastRep DIV 8],LastRep DIV 8) else Begin move(pSieb[0],pSieb[i*LastRep DIV 8],(gblSiebBitCnt-i*LastRep) DIV 8); RepLaenge := gblSiebBitCnt; BREAK; end; end;
with sqFeld[sqIdx] do Begin pDelta30 := @sqDelta[0]; idx := sqOfsIdx AND 3; Offset30 := sqOfsIdx shr 2; end;
repeat tmp := Offset30 DIV 8; pSieb[tmp] := pSieb[tmp] OR (1 shl (Offset30 AND 7)); inc(Offset30,pDelta30[idx]); idx:= (idx+1) AND 3; until (Offset30 >= RepLaenge); end;
procedure EinmalQuadSiebLmt(sqF:tpSegQuadSvPrim;Limit:NativeUInt); label UndTschuesz; var pSieb : tpSieveElem; pDelta30 : tpDelta30_64; Offset30 : tSieveElem; mytmp,Idx: NativeInt; k, Eins : tSieveElem; begin with sqF^ do Begin idx := sqOfsIdx; Offset30 := idx shr 2; idx := idx AND 3; if Offset30>= Limit then EXIT; pDelta30 := @sqDelta[0]; end; pSieb := @QuadSieb[0]; Eins := 1;
while Idx <> 0 do Begin mytmp := Offset30 DIV cElemBitCnt; pSieb[mytmp]:= (Eins shl (Offset30 AND cAndMask)) or pSieb[mytmp]; inc(Offset30,pDelta30^[Idx]); Idx := (Idx + 1) AND 3; if Offset30> Limit then GOTO UndTschuesz; end;
k := (Limit-Offset30) DIV sqF^.sqPrim; while k> 0 do Begin mytmp := Offset30 DIV cElemBitCnt; pSieb[mytmp] := pSieb[mytmp] OR (Eins shl (Offset30 AND cAndMask)); inc(Offset30,pDelta30^[0]);
mytmp := Offset30 DIV cElemBitCnt; pSieb[mytmp]:= pSieb[mytmp] OR (Eins shl (Offset30 AND cAndMask)); inc(Offset30,pDelta30^[1]);
mytmp := Offset30 DIV cElemBitCnt; pSieb[mytmp]:= pSieb[mytmp] OR (Eins shl (Offset30 AND cAndMask)); inc(Offset30,pDelta30^[2]);
mytmp := Offset30 DIV cElemBitCnt; pSieb[mytmp]:= pSieb[mytmp] OR (Eins shl (Offset30 AND cAndMask)); inc(Offset30,pDelta30^[3]);
dec(k); end; while Offset30 < Limit do Begin mytmp := Offset30 DIV cElemBitCnt; pSieb[mytmp]:= pSieb[mytmp] OR (Eins shl (Offset30 AND cAndMask)); inc(Offset30,pDelta30^[Idx]); Idx := (Idx + 1) AND 3; end; UndTschuesz: sqF^.sqOfsIdx := Offset30 shl 2 +idx; end;
procedure Check; Begin write(#13,AnzahlSiebPrim:10,gblCompPot:13,FormatDateTime(' HH:NN:SS.ZZZ', StartZeit - now)); IF gblCompPot>= cE9 then inc(gblCompPot,cE9) else Begin gblCompPot := 4*gblCompPot; if gblCompPot > cE9 then gblCompPot := cE9; end; end;
procedure QuadrupelSiebSieben(idx:LongInt); const KleinePrim : array[0..21] of byte = (7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,97); var pSqFeld,pSqFeld0 : tpSegQuadSvPrim; Limit,deltaLmt: NativeInt; i :longInt; begin pSqFeld := @sqFeld[idx]; If RepFertig then Begin idx := gblsqFeldIdx-idx; IF pSqFeld^.sqPrim > cDeltaBitLimit_II then Begin Limit := gblSiebBitCnt; For i := idx downto 0 do Begin EinmalQuadSiebLmt(pSqFeld,Limit); inc(pSqFeld); end; end else Begin deltaLmt := cDeltaBitLimit_II; pSqFeld0 := pSqFeld; repeat IF Limit>= gblSiebBitCnt then Break; For i := idx downto 0 do Begin EinmalQuadSiebLmt(pSqFeld,Limit); inc(pSqFeld); end; pSqFeld := pSqFeld0; inc(Limit,deltaLmt); until false; Limit := gblSiebBitCnt; For i := idx downto 0 do Begin EinmalQuadSiebLmt(pSqFeld,Limit); inc(pSqFeld); end; dec(pSqFeld); end; end else Begin RepLaenge:=8; QuadSieb[0] := 0; repeat KleineQuadPrimSieben(KleinePrim[idx],Idx); inc(Idx); until RepLaenge = gblSiebBitCnt; RepFertig := true; QuadrupelSiebSieben(Idx); end; end;
procedure MulListErstellen(p30:Uint64;pDelta30:tpDelta30_64); var pMulPos : ^tMulOfs; begin pMulPos := @cMulOfs[p30 AND 7][0]; p30 := p30 shr 3; with pMulPos^ do pDelta30^[0] := p30*moMul+moOfs; pDelta30^[2] := pDelta30^[0]; inc(pMulPos); with pMulPos^ do pDelta30^[1] := p30*moMul+moOfs; inc(pMulPos,2); with pMulPos^ do pDelta30^[3] := p30*moMul+moOfs; end; {$IFDEF CPUX64} {$IFDEF WIN64} {$ASMMODE INTEL} function PowModAsm64_(base,exponent,divisor:NativeUint):NativeUint;assembler;
Label SchonSchlusz,Gerade,NochMalVonVorn; asm push r9 push r10 push r11 xor r11,r11 inc r11 mov r9,r11 test rdx,rdx jz SchonSchlusz test r8,r8 jz SchonSchlusz
mov r10,rdx xor rdx,rdx mov rax,rcx div r8 mov rcx,rdx
NochMalVonVorn: test r9,r10 jz Gerade xor rdx,rdx mov rax,r11 mul rcx div r8 mov r11,rdx Gerade: xor rdx,rdx mov rax,rcx mul rcx div r8 mov rcx,rdx shr r10,1 jnz NochMalVonVorn
SchonSchlusz: mov rax,r11
pop r11 pop r10 pop r9 end; {$ELSE} function PowModAsm64_(base,exponent,divisor:NativeUint):NativeUint;assembler;
asm pushq %r9 xorq %r8,%r8 inc %r8 movq %r8,%r9 testq %rsi,%rsi jz .Lj242 testq %rdx,%rdx jz .Lj242 pushq %rcx movq %rdx,%rcx cmpq %rdi,%rdx ja .Lj241 movq %rdi,%rax xorq %rdx,%rdx divq %rcx movq %rdx,%rdi .Lj241: xorq %rdx,%rdx testq %r9,%rsi je .Lj245 movq %rdi,%rax mulq %r8 divq %rcx movq %rdx,%r8 .Lj245: xorq %rdx,%rdx movq %rdi,%rax mulq %rdi divq %rcx shrq $1,%rsi movq %rdx,%rdi
.Lj243: jne .Lj241 popq %rcx .Lj242: movq %r8,%rax popq %r9 end; {$ENDIF} {$ENDIF}
function PowModAsm64_30(i:Uint64):Uint64; var p3: NativeUint; begin p3 := 3*i; if p3 >= gblFact2 then result := PowModAsm64_(gblBase,gblExpo2,p3)* gblFact2 else result := PowModAsm64_(gblBase,gblExpo2,p3)*(gblFact2 MOD p3); result:= ((result+gblSiebOfsDiv10) DIV 3) MOD i; end;
procedure QuadSiebCalc(sqIdx:NativeInt;p:Uint64); var pSqFeld : tpSegQuadSvPrim; pDelta30 : tpDelta30_64; p30, rest : Int64; idx : NativeInt; Begin rest:= PowModAsm64_30(p); pSqFeld := @sqFeld[sqIdx]; pSqFeld^.sqPrim := p; p30 := p DIV 30; p30 := cMod30ToIDx[p-30*p30]+p30 shl 3; pDelta30 := @pSqFeld^.sqDelta; MulListErstellen(p30,pDelta30); idx := 0; p30 := pDelta30^[3] shr 1; while p30 < rest do Begin inc(p30,pDelta30^[idx]); idx:= (idx+1) AND 3; end; pSqFeld^.sqOfsIdx := (p30-rest) shl 2 +idx; end;
function KleinesPrimSiebAbarbeiten(OG : INT64):NativeInt; var Siebzeiger : pBoolean; p : Uint64; i : NativeInt; begin result := gblsqFeldIdx; p := SiebStartPrim+1; Siebzeiger := @Sieb^[0]; for i := 1 to OG do begin inc(p,2); if Siebzeiger[i] then Begin QuadSiebCalc(result,p); inc(result); end; end; p := result-gblsqFeldIdx; gblsqFeldIdx := result-1; result := p; end;
procedure PrimSiebInit; begin new(WiFeld); new(sieb); new(primfeld);
Wifeld_Aufbauen; gblMaxPrimPos := KleinePrimListeErsieben;
OberesLimit := (gblMaxQuadSiebPrim-1)shr 1; gblSiebNr := 1; SiebStartZahl := 1; ObergrenzeNeu := (SiebStartZahl-1) shr 1; AnzahlSiebPrim := 0; RepLaenge := 8;
gblsqFeldIdx :=Low(tsqFeld); gblAktPrimPos :=StartPrimNr+1; NaechstesKleinPrimSieb; Sieb[( 1-1)shr 1] := false; Sieb[( 7-1)shr 1] := true; Sieb[(11-1)shr 1] := true; Sieb[(13-1)shr 1] := true; end;
procedure PrimSiebEnd; Begin dispose(Primfeld); dispose(Sieb); dispose(WiFeld); end;
procedure KomplettSieben; type tmyChecks = (gesiebt,fertig, nochAusgeben); var cnt: NativeInt; myChecks : set of tmyChecks; Begin myChecks :=[nochAusgeben]; StartZeit := time; PrimSiebInit; gblCompPot := 10; cnt := 0; PowModAsm64_Init(gblExpo-1); repeat OberGrenzeAlt := ObergrenzeNeu; SiebStartPrim := OberGrenzeAlt shl 1; inc(OberGrenzeNeu,Max_2_2); if OberesLimit < OberGrenzeNeu then Begin OberGrenzeNeu := OberesLimit; include(myChecks,fertig); end; inc(AnzahlSiebPrim,KleinesPrimSiebAbarbeiten(OberGrenzeNeu-ObergrenzeAlt)); if gblsqFeldIdx > 35410 then begin QuadrupelSiebSieben(0); include(myChecks,gesiebt); gblsqFeldIdx:= Low(tsqFeld); end; inc(cnt); NaechstesKleinPrimSieb; until fertig in myChecks; If Not(Gesiebt in myChecks) then QuadrupelSiebSieben(0); PrimSiebEnd; end; procedure CheckLog; var l : NativeInt; s,s1: AnsiString; Begin gefunden := false; Reset(PFGW_Log); s1:= ''; while NOT(EOF(PFGW_Log)) do Begin readln(PFGW_Log,s); l := length(s1); IF l > 0 then IF (s[1] <> '1') AND (s1[l] = '8') then Begin setlength(s1,l-2); writeln(' Fund ',s1); s := s1+FormatDatetime(' HH:NN:SS.ZZZ',Now-Time1); writeln(s); {$I-} Append(Funde_Dat); IF IoResult<> 0 then Rewrite(Funde_Dat); {$I+} writeln(Funde_Dat,s); CloseFile(Funde_Dat); gefunden := true; BREAK; end; s1:= s; end; CloseFile(PFGW_Log); Erase(PFGW_Log); If Not Gefunden then Erase(PFGW_Dat);
end; procedure Abarbeiten; var Time0 : TDateTime; begin AssignFile(PFGW_Log,PWGW_Log_Name); Rewrite(PFGW_Log); CloseFile(PFGW_Log); Time0 := now; system.SetTextBuf(PFGW_Dat,myTextBuf); system.SetTextBuf(PFGW_Log,myTextBuf); gblDatCnt := gblSiebOfs DIV gblSiebNumbCnt; {$IFDEF SUCHE} repeat {$ENDIF} ChangeStartCnt(startcnt); gblDatCnt := gblSiebOfs DIV gblSiebNumbCnt+1; time1 := now; repeat writeln('Suche bis ',gblSiebOfs); PFGW_Surfilename:= GetTempDir+Format('Quad_%.4d_%.3d',[StartCnt,gblDatCnt]); PFGW_Tstfilename:= PFGW_Surfilename+'.tst'; KomplettSieben; KandidatenEinsammeln; DateiTesten; CheckLog; If not gefunden then NextQuadSieb; until gefunden or abbruch; if gefunden then {$IFDEF SUCHE} inc(startCnt); IF startCnt > EndCnt then Begin dec(startCnt); abbruch := true; end else Begin gblStartOffset:=0; end; until abbruch; {$ENDIF} writeln('Gesamtzeit: ' + FormatDateTime('HH:NN:SS.ZZZ', time0 - now)); writeln('----------'); end;
procedure Button1Click; begin abbruch := NOT abbruch; IF abbruch then else Abarbeiten; end; function MainInit:boolean; Begin If ParamCount <1 then Begin Writeln('Zu wenig Parmameter'); Writeln('1. Startpotenz 10^ Startpotenz-1 '); Writeln('2. Offset 10^Startpotenz-1 +Offset default= 0)'); result := false; EXIT; end; IF NOT (tryStrtoInt(paramStr(1),longInt(StartCnt))) then StartCnt :=10; If paramCount = 2 then Begin writeln(paramStr(2)); IF NOT (tryStrtoInt64(paramStr(2),gblStartOffset)) then gblStartOffset := 0; end; RepFertig := false; Abbruch:= true; result := true; try endcnt:= startcnt+50; ChangeStartCnt(EndCnt); setlength(QuadSieb,gblSiebElemCnt+32); RepFertig := false; Abbruch:= true; result := true; AssignFile(Funde_Dat,'Funde.txt'); finally end; end;
procedure MainExit; Begin setlength(QuadSieb,0); end;
Begin If MainInit then Button1Click; MainExit; end. |
Ausgabe nun in eine Datei Funde.txt
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:
| 10^99+349781731 00:00:00.542 10^100+1053594241 00:00:01.114 10^101+37648981 00:00:00.573 10^102+514950241 00:00:00.604 10^103+29093071 00:00:00.734 10^104+951859441 00:00:01.300 10^105+659805571 00:00:01.317 10^106+1408385221 00:00:02.059 10^107+453185221 00:00:00.729 10^108+307174981 00:00:00.735 10^109+214038511 00:00:00.754 10^110+2323457011 00:00:03.121 10^111+281187091 00:00:00.833 10^112+584332741 00:00:00.821 10^113+487243951 00:00:00.848 10^114+1137967651 00:00:01.776 10^115+690074251 00:00:00.925 10^116+743896741 00:00:00.967 10^117+210840661 00:00:00.990 10^118+1125848521 00:00:02.055 10^119+2048632891 00:00:03.215 10^120+1718199841 00:00:02.211 10^121+424562761 00:00:01.142 10^122+109586071 00:00:01.172 10^123+3177022021 00:00:03.675 10^124+1871400811 00:00:02.616 10^125+2491588621 00:00:03.953 10^126+140357941 00:00:01.318 10^127+3246722431 00:00:04.077 10^128+18602761 00:00:01.379 10^129+2748589231 00:00:04.395 10^130+19879291 00:00:01.470 10^131+1053097021 00:00:01.524 10^132+2691761761 00:00:03.154 10^133+2633126251 00:00:03.281 10^134+493019851 00:00:01.653 10^135+4269273571 00:00:05.422 10^136+318447541 00:00:01.833 10^137+4015289641 00:00:05.835 10^138+1959210661 00:00:03.985 10^139+469899331 00:00:02.254 10^140+1501977541 00:00:02.445 10^141+3663136111 00:00:06.741 10^142+717523051 00:00:02.220 10^143+6693235021 00:00:09.056 10^144+800150731 00:00:02.357 10^145+131678371 00:00:02.334 10^146+563382661 00:00:02.426 10^147+7503579361 00:00:09.925 10^148+303068311 00:00:02.519 10^149+105012451 00:00:02.591 . 10^304+18543699601 00:01:55.553 |
Gruß Horst
EDIT:
Automatische Anpassung der Quadrupel Siebgröße und maximaler Suchprimzahl.
EDIT2: Anpassung angepasst, pfgw64 für Windows langsamer als für Linux, vielleicht wegen wine
Für diesen Beitrag haben gedankt: Mathematiker
|
|
|