Ich brauche einen Primzahltest, der möglichst schnell und effizient (und sicher) testet, ob eine (seeehr große) Zahl eine Primzahl ist.
Zur Vorgeschichte:
Eigentlich wollte ich das Ganze über das Sieb des Eratosthenes machen.. sprich alle Primzahlen bis zur Wurzel der Zahl berechnen und dann auf Teilbarkeit prüfen.
Problem ist folgendes:
Bis zur Größenordnung 10^300 gibt es etwa 10^294 Primzahlen. Angenommen die sind etwa gleich verteilt muss ich also bis zur Wurzel einer Zahl der Größe 10^300 eine Anzahl an 10^147 Primzahlen berechnen.
Wenn ich das Ganze auf einem 32bit Computer mache mit 4 Byte Zeiger (auf Datenstrukturen, die dann die jeweiligen Primzahlen enthalten), dann verbrauche ich allein für die Zeiger 4 * 10^147 Byte Speicher, was echt nicht geht^^
Deswegen versuche ich es nun mit dem sog. Miller-Rabin-Test.
Der liefert zwar nicht mit 100% Sicherheit Primzahlen, allerdings beschränkt sich der Fehler auf (1/4)^k bei k Durchläufen, was akzeptabel ist.
Nun hab ich das in Code umgesetzt und es scheint auch zu funktionieren, jedoch bin ich mir sicher, dass man das noch optimieren kann (siehe die Geschwindigkeit von DIESEM Test bei Zahlen mit mehr als 500 Stellen)
Nur weiß ich nicht WIE.
Ich glaube es hängt damit zusammen, wie ich meine zufällig gewählten Zahlen wähle
(Achja: xint sind vorzeichenbehaftete ints beliebiger Größe)
Code: Alles auswählen
// Modulare Exponentation
// Berechnet also (num^exp)%n
// Ich denke auch, dass es funktioniert ;)
xint mod_exp(xint num, xint exp, xint n)
{
xint res;
res = 1;
xint aktpot;
aktpot = num;
while ( exp > 0) {
if (exp % 2 == 1) {
res = res * aktpot;
}
aktpot = aktpot * aktpot;
aktpot = aktpot % n;
exp = exp / 2;
}
res = res % n;
return res;
}
bool miller_rabin(xint n)
{
xint a;
xint t;
xint alpha;
xint p;
bool pseudoprim;
// Zufällige Basis a wählen
srand(time(NULL)+4);
a = (rand()%123456)+2;
while (a < n) {
a *= (rand()%123457)+1;
}
while (a >= n) {
a -= n-1;
}
// Hier wird der Exponent berechnet
t = n - 1;
alpha = 0;
while ( t%2 == 0) {
t = t/2;
alpha = alpha + 1;
}
// Und nun das Ergebnis von (a^t)%n
p = 1;
p = mod_exp(a,t,n);
pseudoprim = false;
alpha = alpha - 1;
if ((p == n-1) || (p==1)) {
pseudoprim = true;
}
while ((!pseudoprim) && (alpha>0) && (p>1)) {
p = p * p;
p = p % n;
alpha = alpha - 1;
if (p == n-1) {
pseudoprim = true;
}
}
return pseudoprim;
}