Diese Woche geht es gleich weiter mit der Lösung des siebten Projecteuler-Problems mit Python. Die guten Vorsätze von letzter Woche scheinen noch nicht abgelegt zu sein, deswegen beschäftigen wir uns in diesem Artikel wieder etwas mit Primzahlen und dem Sieb des Eratosthenes. Natürlich wird das Problem mit der Programmiersprache Python angegangen und gelöst!
Problemvorstellung
Das siebte Projecteuler-Problem ist ein sehr kurzes Problem! Zumindest was die textliche Beschreibung angeht. Laut dem Motto “In der Kürze liegt die Würze” lautet die Aufgabenstellung folgendermaßen:
Finde die 10001.Primzahl!
Also kann mit drei Worten, einer Zahl und zwei Satzzeichen etwas formuliert werden, was jedes Grüblerherz höher schlagen lässt.
Diejenigen Leser, die meine Serie “Projecteuler: Lösungen mit Python und C++” bisher mitverfolgt haben, auch ohne aktives Eingreifen in Form von Kommentaren im Blog oder auch ohne die Probleme selber gelöst zu haben, werden wissen, dass es bereits eine Aufgabe gab, die Primzahlen zum Thema hatte.
Im Rahmen der Lösung zum dritten Problem, wo es um eine Primfaktorenzerlegung ging, finden wir wichtige Werkzeuge auch zur Lösung des siebten Problems. Doch dazu später mehr, zuvor etwas Statistik:
Ganz aktuell gemessen, haben von insgesamt 112279 registrierten Projecteuler-Nutzern genau 57841 das Problem mit Bravour gelöst. Das entspricht ziemlich genau der Hälfte der Nutzer. Und wir sind dabei 
!
Problemeingrenzung
Die Schwierigkeiten bei dieser Fragestellung liegen auf der Hand:
- Zunächst ist es nett zu wissen, welche Zahlen denn nun überhaupt Primzahlen sind. Denn nur dann können sie auch gezählt werden. Also ist die erste Herausforderung das Filtern von Primzahlen.
- Da nicht alle Primzahlen bis zu einer bestimmten ganzen Zahl gesucht sind, sondern eine bestimmte Anzahl an Primzahlen, sind die gängigsten Werkzeuge erst mal nicht zu gebrauchen, oder besser: sind zu überdenken. Es müssen weiterhin alle Primzahlen gesucht und gezählt werden, ein Limit ist jedoch nicht gegeben.
Wie eingangs bereits angedeutet, können wir uns aus der Lösung zu Problem 3 eines ganz nützlichen Werkzeugs bedienen: dem Sieb des Eratosthenes.
Mit diesem Sieb lassen sich alle Primzahlen bis zu einer gegebenen ganzen Zahl N bestimmen. Wie im verlinkten Wikipedia-Artikel erläutert, müssen hierbei einfach alle Zahlen bis N “aufgeschrieben” werden. Dann werden mit der 2 beginnend alle Vielfachen wieder gestrichen. Hat man das für die 2 erledigt, so geht es in der Reihenfolge mit der nächsthöheren Zahl weiter, also der 3, dann 5, 7, 11, …
Das Sieb des Eratosthenes liefert uns also alle Primzahlen, die kleiner oder gleich N sind. Gesucht ist jedoch die 10001.Primzahl. Somit kommen wir mit dem Sieb des Eratosthenes nicht direkt zum Ziel.
Sollen nun im Zahlenbereich bis N noch keine 10001 Primzahlen identifiziert sein, so muss der Bereich erweitert und die Suche fortgesetzt werden. Wie das realisiert werden kann, steht unter der Überschrift:
Problemlösung
Der erste Teil der Lösung ist klar und bereits aus der Lösung zum Problem 3 bekannt:
def eratosthenes_vect(N):
primes = np.arange(3,N+1,2) # only uneven figures
primes = np.insert(primes,0,2,0) # ... and two
nUneven = np.size(primes) - 1.0
for i in np.arange(3,np.ceil(np.sqrt(N))+2,2):
primes[(i*i+1)/2-1:nUneven+1.0:i] = 0;
return np.trim_zeros(np.sort(primes))
Damit erhalten wir alle Primzahlen bis zur Zahl N.
Die einfachste Lösung wäre nun, einfach einen großes N vorzugeben und zu hoffen, dass es mindestens 10001 Primzahlen im Sieb hängen bleiben. Ist man in der Lage, ein gut passendes N abzuschätzen, so hat man mit dem Sieb des Eratosthenes eine sehr effiziente Lösungsmethode gefunden.
Diese Tatsache möchten wir im Folgenden ausnutzen und durch eine Erweiterung etwas absichern:
def eratosthenes_limits(StartN, EndN, Primes):
"""
Rechnet alle Primzahlen im Zahlenraum StartN, EndN.
Primes: alle Primzahlen bis StartN.
"""
if (StartN%2 == 0):
SearchField = np.arange(StartN+1,EndN+1,2)
limitleft = StartN+1
limitright = EndN
else:
SearchField = np.arange(StartN,EndN+1,2)
limitleft = StartN
limitright = EndN
for item in Primes:
itemQ = item*item
if itemQ <= limitright:
for i in range(0,SearchField.size,1):
if (SearchField[i]%item == 0):
SearchField[i] = 0
SearchField = np.trim_zeros(np.sort(SearchField))
limitleft = SearchField[0]
limitright = SearchField[-1]
return np.append(Primes,SearchField)
def solver(N,NEra):
Primes = eratosthenes_vect(NEra)
while (Primes.size < N):
Primes = eratosthenes_limits(NEra+1, NEra*2, Primes)
NEra = NEra*2
print("Suchbereich bis N=%i" %NEra)
return Primes[N-1]
Der gegebene Vorschlag ist wahrscheinlich nicht die eleganteste Methode, aber er funktioniert für das Problem recht ordentlich.
Die Funktion eratosthenes_limits kommt zum Einsatz, wenn der Zahlenbereich, den wir geschätzt haben, zu klein ist. Sie geht im Grunde genau nach dem Prinzip des Siebes von Eratosthenes vor. Es werden alle Zahlen in einem Bereich zwischen StartN und EndN auf Teilbarkeit durch die Primzahlen bis zur Zahl StartN untersucht. Da wir davon ausgehen, dass der zu Beginn geschätzte Zahlenbereich nicht sehr stark vergrößert werden muss, steigern wir in der solver-Funktion den Suchbereich um das doppelte, wenn nicht genug Primzahlen gefunden wurden.
Auf diese Weise kann es nicht vorkommen, dass die Funktion eratosthenes_limits ins Stocken gerät, weil Primzahlen aus dem aktuellen neuen Suchbereich benötigt werden, um alle Vielfachen zu streichen. Dafür müsste der Suchbereich schon quadratisch vergrößert werden: StartN = N, EndN=N*N.
Schätzen des Suchbereichs zu Beginn
Das Schätzen des richtigen Suchbereichs ist das A und O bei der vorgeschlagenen Lösung. Dazu lohnt sich ein Blick auf die Statistik. Betrachtet man bereits die Häufigkeit von Primzahlen im Bereich von 1-10 und von 20-30, so stellt man fest, dass die Anzahl abnimmt. Diese Tatsache lässt sich mathematisch beweisen, was wir hier aber nicht machen wollen.
In dem Zahlenbereich, der für das siebte Projecteuler-Problem interessant ist, treten laut der Analyse auf maerku69.ch Primzahlen mit einer Häufigkeit von etwas unter 10% auf.
Daraus schließen wir, dass eine Schätzung des Suchbereichs mit etwas mehr als dem zehnfachen von 10001 ganz ordentlich sein sollte.
Lösungsanalyse
Genug der trockenen Theorie, nun kommt die Analyse:
import time as t
import numpy as np
def eratosthenes_vect(N):
...
def eratosthenes_limits(StartN, EndN, Primes):
...
def solver(N,NEra):
...
N = 10001
N1 = 11*N
t0 = t.clock()
Prim = solver(N,N1)
tn = t.clock() - t0
print("Benoetigte Zeit: %15.6f" %tn)
print("Gefundene Primzahl: %i" %Prim)
Mit einem derart optimierten Suchbereich (11*10001) ist eine Erweiterung nicht nötig, die Lösung wird allein durch das Sieb des Eratosthenes gefunden (Ausgabe):
Suchbereich bis N=110011
Benoetigte Zeit: 0.320000
Gefundene Primzahl: 104743
Unterschätzt man den Suchbereich etwas (10*10001), dann dauert die Suche schon etwas länger:
Suchbereich bis N=200020
Benoetigte Zeit: 2.810000
Gefundene Primzahl: 104743
Hierbei wurde der Suchbereich ein mal vergrößert.
Wird dagegen der Suchbereich überschätzt (12*10001), dann dauert die Suche nur unwesentlich länger:
Suchbereich bis N=120012
Benoetigte Zeit: 0.350000
Gefundene Primzahl: 104743
Schlusswort
Primzahlen sind immer wieder gut für Knobel-Spielchen geeignet, so auch diesmal.
Das Lösen des siebte Problem hat mir persönlich viel Spass bereitet. Ich hoffe, dass es meinen Lesern nicht anders ging oder noch gehen wird und freue mich auf andere Lösungsvorschläge in den Kommentaren oder vielleicht als Reaktion im eigenen Blog.
Meinungen, Anregungen, Verbesserungen sind immer willkommen.
Schätzen des Suchbereichs zu Beginn
Das Schätzen des richtigen Suchbereichs ist das A und O bei der vorgeschlagenen Lösung. Dazu lohnt sich ein Blick auf die Statistik. Betrachtet man bereits die Häufigkeit von Primzahlen im Bereich von 1-10 und von 20-30, so stellt man fest, dass die Anzahl abnimmt. Diese Tatsache lässt sich mathematisch beweisen, was wir hier aber nicht machen wollen.
In dem Zahlenbereich, der für das siebte Projecteuler-Problem interessant ist, treten laut der Analyse auf maerku69.ch Primzahlen mit einer Häufigkeit von etwas unter 10% auf.
Daraus schließen wir, dass eine Schätzung des Suchbereichs mit etwas mehr als dem zehnfachen von 10001 ganz ordentlich sein sollte.
Lösungsanalyse
Genug der trockenen Theorie, nun kommt die Analyse:
import time as t
import numpy as np
def eratosthenes_vect(N):
...
def eratosthenes_limits(StartN, EndN, Primes):
...
def solver(N,NEra):
...
N = 10001
N1 = 11*N
t0 = t.clock()
Prim = solver(N,N1)
tn = t.clock() - t0
print("Benoetigte Zeit: %15.6f" %tn)
print("Gefundene Primzahl: %i" %Prim)
Mit einem derart optimierten Suchbereich (11*10001) ist eine Erweiterung nicht nötig, die Lösung wird allein durch das Sieb des Eratosthenes gefunden (Ausgabe):
Suchbereich bis N=110011
Benoetigte Zeit: 0.320000
Gefundene Primzahl: 104743
Unterschätzt man den Suchbereich etwas (10*10001), dann dauert die Suche schon etwas länger:
Suchbereich bis N=200020
Benoetigte Zeit: 2.810000
Gefundene Primzahl: 104743
Hierbei wurde der Suchbereich ein mal vergrößert.
Wird dagegen der Suchbereich überschätzt (12*10001), dann dauert die Suche nur unwesentlich länger:
Suchbereich bis N=120012
Benoetigte Zeit: 0.350000
Gefundene Primzahl: 104743
Schlusswort
Primzahlen sind immer wieder gut für Knobel-Spielchen geeignet, so auch diesmal.
Das Lösen des siebte Problem hat mir persönlich viel Spass bereitet. Ich hoffe, dass es meinen Lesern nicht anders ging oder noch gehen wird und freue mich auf andere Lösungsvorschläge in den Kommentaren oder vielleicht als Reaktion im eigenen Blog.
Meinungen, Anregungen, Verbesserungen sind immer willkommen.
-->Artikel aus der selben Kategorie:
Meinungsfreiheit für alle!
Es gibt noch keine Kommentare zu “Projecteuler: Lösung zu Problem 7”.
Jede Meinung ist willkommen!