Programmierwettbewerb: Primzahlzerlegung

[i_hasser]

Hi

Ziel ist es für eine sehr große Zahl möglichst schnell alle Darstellungsmöglichkeiten durch a*b zu finden, wo bei a und b eine Primzahl sein muss.

Also zb. 21=3*7

Solche Zerlegungen spielen zb. beim Knacken von Verschlüsselungen eine wichtige Rolle.
Sprachen usw. sind alle erlaubt, sagen wir außer reinem Assembler (inline Assembler ist natürlich auch erlaubt). Betriebssysteme geht auch alles, aber ich denk es wäre erstmal besser wir bleiben bei x86 (-16, -32 und -64 ;)). Lösungen auf anderen Architekturen sind natürlich auch willkommen.

Und damit wir nicht in das Dilemma mit Zahlenlibs kommen denk ich beschränken wir das erstmal auf ausgewählte 64bit Integer und sehen wie die Zeiten dafür aussehen. Zum Testen gibts hier eine Liste mit gültigen Zahlen (die sich auch tatsächlich in a*b für Primzahlen zerlegen lassen):

21
(kommt noch mehr ;))


Viel Spaß!

[mj]

So, ich liefere gleich mal den ersten Beitrag, erstmal unkommentiert. Die Kommentare mach ich grade, wer es also noch nicht versteht, etwas Geduld bitte ;)

Prim.java (http://home.in.tum.de/~jungowsk/sw/java/prime/Prim.java)
Start.java (http://home.in.tum.de/~jungowsk/sw/java/prime/Start.java)

Das Prinzip dahinter ist ganz simpel: Es wird zunächst ein Feld erzeugt mit allen Primzahlen zw. 0 und 400. Dies ist variabel, so dass diese Lösung grundsätzlich auch für größere Zahlen geeignet ist.
Anschließend wird die gegebene Zahl (in unserem Fall 21) der Reihe nach, von der kleinsten Zahl beginnend durch jede Zahl die in dem Feld gespeichert ist geteilt und festgestellt, ob das Ergebnis auch eine Primzahl ist. Ist dem der Fall, dann haben wir das Ergebnis - fertig.

Edit: So, mittlerweile hab ich auch die Kommentare hinzugefügt, mehr gibt's da eigentlich nicht zu sagen. Für kleine Zahlen ist diese Lösung sehr schnell, mal sehen was ich ändern muss wenn größere Zahlen ins Spiel kommen. Natürlich ist diese Lösung as-is höchstens für Zahlen des Datentyps long zulässig.

[[ab]noname]

Mein unoptimierter C++ Anfang:


#include <iostream>

using namespace std;
typedef __int64 i64;

long IsPrime(i64 nr) {
for(i64 d = 2; d*d <= nr+1; d++) {
if (!(nr % d)) {
return false;
}
}
return true;
}

void Search2PrimFacts(i64 Number) {
i64 a;

for(i64 i = 3; i <= Number; i += 2) {
if (IsPrime(i)) {
if (Number % i == 0) {
a = Number / i;
if (IsPrime(a)) {
printf("%I64d\n", a);
printf("%I64d\n", i);
break;
}
}
}
}
}

int main() {
Search2PrimFacts(21);
}

[JKuehl]

ich werde die tage mal probieren unseren lenstra mit elliptischen kurven algo von maple auf c zu porten.

evtl vorher auch die pollard p-1 methode, weil die wesentliche einfacher zu implementieren ist.

[mj]

Ich denk, ich werd mich heut Abend mal dran machen, das Sieb des Eratosthenes in Java implementieren. Das sollte die Erstellung des Arrays mit Primzahlen stark beschleunigen und auch für größere Zahlen fitmachen ;)

@[ab]noname:
Dein Code compiliert nicht.

[[ab]noname]

Original geschrieben von D'Espice
...
@[ab]noname:
Dein Code compiliert nicht.
Compiler?
Ist mit VC++ 6.0 geschrieben. Ich weiss jetzt nicht aus dem Stehgreif, ob der Datentyp __int64 auch bei anderen Compilern verfügbar ist. Sonst setz bei dem Typedef einfach mal long ein - dann ist natürlich bei 32bit Schluss.

[mj]

Compiler gcc 2.95.4 unter Debian Linux. Funktioniert im übrigen auch nicht, wenn ich __int64 durch long ersetzt im typedef. Die genaue Fehlermeldung:

martin@debian:~/sw$ g++ prime.cpp -o prime
prime.cpp: In function `void Search2PrimFacts(long int)':
prime.cpp:23: implicit declaration of function `int printf(...)'

Ich installiere grade gcc-3.3, mal sehen ob's damit geht.

[[ab]noname]

aso, da wird wohl im #include die cstdlib gebraucht...

[mj]

Nope, auch ein

#include <cstdlib>

hilft nicht, gleicher Fehler.

[[ab]noname]

omg ;D

dann nimm mal das (halt nur 32bit):


#include <iostream>

using namespace std;
typedef long i64;

long IsPrime(i64 nr) {
for(i64 d = 2; d*d <= nr+1; d++) {
if (!(nr % d)) {
return false;
}
}
return true;
}

void Search2PrimFacts(i64 Number) {
i64 a;

for(i64 i = 3; i <= Number; i += 2) {
if (IsPrime(i)) {
if (Number % i == 0) {
a = Number / i;
if (IsPrime(a)) {
//printf("%I64d\n", a);
//printf("%I64d\n", i);
cout << a << " " << i << endl;
break;
}
}
}
}
}

int main() {
Search2PrimFacts(21);
}


Das nächste mal code ich besser gleich mit linux, dann geht es wenigstens überall...

EDIT:
Hab inzwischen mein gutes Debian laufen. Also habs mir schon gedacht, der GCC kommt nicht mit den Formatangaben für int64 zurecht. War also wieder mal ein propr. MS Feature. Der 2. Code compiliert mit dem GCC 3.x bei mir.

[mj]

So, jetzt hab ich das ganze mittels BigInteger implementiert. Sieht jetzt zwar deutlich leichter und strukturierter aus, ist aber langsamer. Das liegt nicht nur daran, dass
Java's BigInteger relativ lahm ist (Strukturbedingt) sondern auch daran, dass eine Implementierung eines Primzahlsiebs nicht nur unnötig kompliziert, sondern auch unnötig geworden ist.

Der Quellcode kann hier (http://home.in.tum.de/~jungowsk/sw/java/bigprime/) gefunden werden, derzeit läuft grade ein Testlauf auf einem Pentium 4 / 1.8 GHz und Debian Linux mit der kleinsten der von JKuehl geposteten Zahlen.

Edit: HAHA! Im Moment in dem ich es schreibe, bringt er auch schon ein Ergebnis, ich werd jetzt auch mal die anderen Zahlen durchlaufen lassen. An der nächstgrößeren Zahl rechnet er jetzt schon eine Ewigkeit rum...



Ergebnisse Intel Pentium4 (1.8 GHz):

outdated, siehe unten


Ergebnisse AMD Opteron 144 (1.8 GHz):

outdated, siehe unten

[i_hasser]

Mal sehen ob ich morgen Zeit hab (im Moment hab ich leider garkeine :() - dann werd ich euch mit meinem C-Hickhack beglücken *chatt*

[mj]

Hmm... der zweite Durchlauf läuft immer noch, mittlerweile schon gut über 40 Minuten :-/
Ich glaub, da ist noch einiges an Optimierungsarbeit zu leisten. Werd's trotzdem mal über Nacht auf'm Opteron laufen lassen, allerdings nicht die zweitkleinste sondern gleich die größte der von JKuehl genannten Zahlen *buck*

Edit:
HOLY F***
Vergesst die Ergebnisse oben, durch minimales Optimieren bin ich jetzt auf folgendes Ergebnis gekommen (gemessen auf einem AMD Opteron 144 (1.8 GHz) mit 1024MB DDR333 Speicher unter Debian Linux, Java 1.4.2):

1949565409282837: 0m6.785s
312133288776418481: 2m11.723s
11081452861688528621: 20m48.564s

Der optimierte Quellcode findet sich hier (http://home.in.tum.de/~jungowsk/sw/java/bigprimeopt). Und für diejenigen, die Probleme habe das ganze mit javac zu compilieren, hab ich auch gleich noch die .class Dateien beigefügt.

[BoMbY]

Tagchen,

wo ich gerade das hier sehe:

for(i64 i = 3; i <= Number; i += 2)

Ihr braucht übrigens nur soviel hier zu testen:

for(i64 i = 3; i <= sqrt(Number); i += 2)

weil mathematisch beweisbar:

sqrt(pq) > p oder sqrt(pq) > q

heißt im Prinzip, einer der beiden Primfaktoren ist auf jedenfall kleiner als die Wurzel des Produkts.

Ich würde also sowas machen (für max. zwei Primfaktoren):

void Search2PrimFacts(i64 Number) {
i64 i = 3; // naja, die zwei ist ja fast wirklich egal
bool found = false;

while ( i <= round(sqrt(Number)) && !found)
{
if (isPrime(i))
{
if (Number % i == 0)
{
printf("%I64d*%I64d=%i64d\n", i, Number / i, Number);
found = true;
}
}
i += 2;
}
}

m.f.g.
BoMbY

[BoMbY]

Wobei ich gerade Überlege: Ist die statistische Wahrscheinlichkeit nicht höher das die Gesuchte Primzahl p irgendwo zwischen sqrt(pq)/2 und sqrt(pq) liegt, als zwischen 3 und sqrt(pq)/2? Wenn dem so ist würde ich nicht von 3 hochzählen sondern von sqrt(pq) runter. Wobei, je höher die getestete Zahl, desto länger dauert es ja auch das Modulo zu berechnen...

Das müsste ich mal genau ausrechnen/testen.

m.f.g.
BoMbY

[[ab]noname]

Original geschrieben von BoMbY
...
void Search2PrimFacts(i64 Number) {
i64 i = 3; // naja, die zwei ist ja fast wirklich egal
bool found = false;

while ( i <= round(sqrt(Number)) && !found)
{
if (isPrime(i))
{
if (Number % i == 0)
{
printf("%I64d*%I64d=%i64d\n", i, Number / i, Number);
found = true;
}
}
i += 2;
}
}

Das zählen bis zur Quadratwurzel ist eigentlich schon in meiner IsPrime Funktion integriert und vergleich auch nicht das Quadrat, sondern das Produkt, da Wurzelziehen langsamer ist, als ein Produkt zu berechnen:
- i*i <= Number ist äquivalent zu i <= sqrt(Number), jedoch schneller!

Dennoch ist die Frage nach dem Sinn, da mein ursprünglicher Algorithmus bei dem ersten gefundenen Ergebnis abgebrochen hat. Theoretisch bleibt dann der Aufwand sowieso unterhalb der SQRT Grenze, wenn es denn stimmt, dass bis zu sqrt(Zahl) je 2 Primzahlen als Faktoren gefunden wurden.
Hätteste du nähere Infos zu dem Beweis? In der inneren Schleife reicht ein Test bis zur Quadratwurzel, jedoch bin ich mir in der äußeren dessen nicht sicher.

Original geschrieben von BoMbY
Wobei ich gerade Überlege: Ist die statistische Wahrscheinlichkeit nicht höher das die Gesuchte Primzahl p irgendwo zwischen sqrt(pq)/2 und sqrt(pq) liegt, als zwischen 3 und sqrt(pq)/2? Wenn dem so ist würde ich nicht von 3 hochzählen sondern von sqrt(pq) runter. Wobei, je höher die getestete Zahl, desto länger dauert es ja auch das Modulo zu berechnen...

Das müsste ich mal genau ausrechnen/testen.

m.f.g.
BoMbY
Ich komme eben aus der Schule und hab in Programmieren etwas weitergebastelt und bin dabei auch auf die Idee gekommen, dass man von einer "Mitte" ausgehend die beiden Faktoren suchen kann. Einen Geschwindigkeitsvorteil hab ich allerdings noch nicht erreicht.


Die 1. o.g. Zahl wird in 288s zerlegt, das ist weitaus langsamer als D'Espice Verfahren, also heissts weiter optimieren.

[BoMbY]

Also wenn ich nicht ganz dumm bin, ist das hier der Beweis (V=oder):

Annahme: sqrt(pq) > p V sqrt(pq) > q

Zwischenschritt: pq > p² V pq > q²

Schluss: q > p V p > q

Sprich also q ist größer p oder p ist größer als q <- und davon können wir ja wohl ausgehen, oder? Okay, ich war jetzt zu faul den Fall (gleich =) zu berücksichtigen, aber der wird auch abgedeckt.

m.f.g.
BoMbY

Edit: PS: Klar, das ist nicht schneller wenn man einen Primfaktor findet, nur wenn man keinen findet!

[Sargnagel]

Original geschrieben von [ab]noname
Das zählen bis zur Quadratwurzel ist eigentlich schon in meiner IsPrime Funktion integriert und vergleich auch nicht das Quadrat, sondern das Produkt, da Wurzelziehen langsamer ist, als ein Produkt zu berechnen:
- i*i <= Number ist äquivalent zu i <= sqrt(Number), jedoch schneller!

Sicher ist es schneller zu quadrieren, als die Wurzel zu ziehen. Aber da Number in der while-Schleife nicht verändert wird, braucht die Berechnung nur einmal pro Aufruf der Funktion Search2PrimFacts() durchgeführt zu werden. Dann sollte sich das wieder relativieren, es sei denn Du rufst besagte Funktion sehr oft (!) auf.
Theoretisch sollte der Compiler den Quellcode dahingehend optimieren, daß die angesprochene Berechnung auch nur einmal durchgeführt wird, aber zur Sicherheit solltest Du das mal austesten, indem Du round(sqrt(Number)) vor der while-Schleife berechnest.


Und die Boolsche Variable brauchst Du auch nicht. Einfach mit break die Schleife beenden. Dadurch sparst Du pro Durchlauf einen Test (!found).

void Search2PrimFacts(i64 Number) {
i64 i = 3, // naja, die zwei ist ja fast wirklich egal
tmp;

tmp = round(sqrt(Number));

while ( i <= tmp)
{
if (isPrime(i))
{
if (Number % i == 0)
{
printf("%I64d*%I64d=%i64d\n", i, Number / i, Number);
break;
}
}
i += 2;
}
}

[Shaft99]

Nur mal so ne generelle Frage,
brecht ihr nach der ersten Möglichkeit ab, oder lasst
ihr alle Möglichkeiten ausrechnen ???
Wenn ich aber Zeiten von 6s sehe, denke ich, dass nach der
1. gefundenen Möglichkeit abgebrochen wird.
Interessant wird IMHO aber erst wenn man alle sucht,
das dauert dann wirklich mal nen paar Stunden für so Zahlen wie oben genannt.

[[ab]noname]

Das ist korrekt...

Der Algorithmus ist dennoch erheblich zu langsam. Interessant ist die Frage woran es hängt.

D'Espice liest erstmal die Primzahlen in ein Array so wie ich das verstanden habe - berechnen muss man sie aber dennoch. Bei mir wird halt pro Primzahl getestet, wie es sich mit dem 2. Faktor der bei einer Division entsteht verhält. - kein Plan, warum das so hinkt :]

[[ab]noname]

Original geschrieben von Sargnagel
Sicher ist es schneller zu quadrieren, als die Wurzel zu ziehen. Aber da Number in der while-Schleife nicht verändert wird, braucht die Berechnung nur einmal pro Aufruf der Funktion Search2PrimFacts() durchgeführt zu werden. Dann sollte sich das wieder relativieren, es sei denn Du rufst besagte Funktion sehr oft (!) auf.
Theoretisch sollte der Compiler den Quellcode dahingehend optimieren, daß die angesprochene Berechnung auch nur einmal durchgeführt wird, aber zur Sicherheit solltest Du das mal austesten, indem Du round(sqrt(Number)) vor der while-Schleife berechnest.
[/EDIT]
In der Funktion IsPrime sollte das schon einen Geschwindigkeitsvorteil bringen. In der äußeren Schleife hast du dagegen Recht.

[Sargnagel]

Original geschrieben von [ab]noname
In der Funktion IsPrime sollte das schon einen Geschwindigkeitsvorteil bringen. In der äußeren Schleife hast du dagegen Recht.
Hmm .. ich habe gerade Deinen Sourcecode ein paar Postings weiter oben entdeckt. Da hast Du das ja schon ganz anders geregelt. Jetzt sehe ich Deinen vorherigen Kommentar in einem ganz anderen Lichte.
Ich lasse mal meine grauen Zellen ackern, was da eventuell noch zu optimieren ist ... *kopfkratz



for(i64 d = 2; d*d <= nr+1; d++)

Das schreit geradezu danach, optimiert zu werden! Vergleich doch nicht das Quadrat von d mit nr+1, sondern d mit der Wurzel aus nr+1 oder was auch immer in Frage käme ... so genau habe ich nun nicht darüber nachgedacht. Die Wurzel natürlich nur 1x vor der Schleife berechnen! So sparst Du irre viele Quadrierungen!

[[ab]noname]

Selbst n-faches Quadrieren ist bei Zahlen dieser Größenordnung schneller als die Wurzel ziehen.

long IsPrime(i64 nr) {
for(i64 d = 2; d <= sqrt(nr)+1; d++) {
if (!(nr % d)) {
return false;
}
}
return true;
}
erreicht:
98326469 * 19827473 = 1949565409282837
t[s] = 330


long IsPrime(i64 nr) {
for(i64 d = 2; d*d <= nr+1; d++) {
if (!(nr % d)) {
return false;
}
}
return true;
}
erreicht:
98326469 * 19827473 = 1949565409282837
t[s] = 288

[Sargnagel]

Original geschrieben von [ab]noname
Selbst n-faches Quadrieren ist bei Zahlen dieser Größenordnung schneller als die Wurzel ziehen.
Da ist doch was faul ... nämlich der Compiler. Wie ich oben schon schrieb, führe die Berechnung der Wurzel doch mal vor der Schleife aus. Ich habe das ungute Gefühl, das für jeden Durchlauf der while-Schleife die Wurzelberechnung beim Test der Schleifenbedingungen erneut ausgeführt wird, was natürlich ein Performancefresser ist!


bool IsPrime(i64 nr) {
i64 tmp = ((i64)sqrt(nr))+1;

for(i64 d = 2; d <= tmp; d++) {
if (!(nr % d)) {
return false;
}
}
return true;
}

Benchmarke doch bitte mal mit obiger Funktion.

Und die Funktion sollte doch sicherlich bool und nicht long zurückgeben?

[[ab]noname]

Hab die Wurzel ausgelagert, 330s... Das optimiert auch jeder anständige Compiler.

Ich bin fest der Überzeugung, dass ein FP Sqrt viel zu viel Zeit braucht, als ne (vergleichsweise) simple Multiplikation, gleich, dass sie vielfach vorkommt.

Habs auch erst nicht geglaubt... :]

Code:

#include <iostream>
#include <time.h>
#include <math.h>

using namespace std;
typedef __int64 i64;

bool IsPrime(i64 nr) {
for(i64 d = 2; d*d <= nr; d++) {
if (!(nr % d)) {
return false;
}
}
return true;
}
//Unoptimiert:
// 288sec:
// 98326469 * 19827473 = 1949565409282837
void Search2PrimFacts(i64 Number) {
i64 a;

for(i64 i = 3; (i*i) <= Number; i += 2) {
if (IsPrime(i)) {
if (Number % i == 0) {
a = Number / i;
if (IsPrime(a)) {
printf("%I64d * %I64d = %I64d\n", a, i, Number);
break;
}
}
}
}
}

int main() {
time_t t1 = time(NULL);
Search2PrimFacts(1949565409282837);
cout << "t[s] = " << time(NULL) - t1;
}

[mj]

Original geschrieben von [ab]noname
D'Espice liest erstmal die Primzahlen in ein Array so wie ich das verstanden habe - berechnen muss man sie aber dennoch. Bei mir wird halt pro Primzahl getestet, wie es sich mit dem 2. Faktor der bei einer Division entsteht verhält. - kein Plan, warum das so hinkt :]
Das war die erste Fassung, ja. Hat Vorteile, wenn man mehrere Zahlen hintereinander testen will weil dann nur einmal das Array durchlaufen wird was mit einer Komplexität von O(n) geschieht.
Die zweite Version, auf BigInteger basierend ist bei einzelnen Zahlen viel schneller da erst getestet wird, ob die Zahl modulo Schleifenzähler gleich Null ist. Ist dem nicht so, wird der Schleifenzähler erhöht und weitergemacht, ist dem so wird getestet, ob der Schleifenzähler eine Primzahl ist (BigInteger.isProbablePrime()) und ist auch das erfüllt, dann wird geschaut ob die Zahl geteilt durch Schleifenzähler auch eine Primzahl ist - wenn das alles zutrifft (schön der Reihe nach), dann haben wir ein Ergebnis.

[Sargnagel]

Original geschrieben von [ab]noname
Hab die Wurzel ausgelagert, 330s... Das optimiert auch jeder anständige Compiler.

Ich bin fest der Überzeugung, dass ein FP Sqrt viel zu viel Zeit braucht, als ne (vergleichsweise) simple Multiplikation, gleich, dass sie vielfach vorkommt.

Habs auch erst nicht geglaubt, aber die Code zeigts doch? :]
Jaja, okay. Ich seh's ja ein. Hatte aus den Augen verloren, daß isPrime() andauernd aufgerufen wird. *kippt Asche über's eigene Haupt* Da es sich um eine Integer-Wurzel handelt, muß es doch haufenweise Optimierungen für eine solche Berechnung geben! Vielleicht sogar mit Hilfe von MMX. Was sagen denn die AMD bzw. Intel TechDocs? Steht da zufällig schon eine Musterlösung drin?

[mj]

Naja nein, nicht wirklich. Das ziehen einer Wurzel ist einfach extrem aufwändig.

@[ab]noname:
Für welche Zahl braucht dein Algorithmus 330s? Die erste der von JKuehl genannten? Wenn ja, dann musst du noch sehr viel optimieren, ich bin da mit Java bereits bei unter 7s angekommen ;)
Du kannst deinen Algorithmus extrem beschleunigen, indem du die Reihenfolge bei der Abfrage umdrehst. Du prüfst derzeit zunächst für jede beliebige Zahl ob sie Prim ist und erst dann, ob sie Modulo Null ergibt - dreh das um! Es reicht wenn du nur die Zahlen auf prim prüfst, die wirklich Modulo gerechnet auch Null ergeben. Das spart dir unglaublich viel Zeit ein, du wirst erstaunt sein.

[[ab]noname]

Hrhr der Modulo Trick ;D

Jetzt hab ichs auf t = 0 bis 1sec für die 1. Zahl - Yippieeh! thx D'Espice


#include <iostream>
#include <time.h>
#include <math.h>

using namespace std;
typedef __int64 i64;

bool IsPrime(i64 nr) {
for(i64 d = 2; d*d <= nr; d++) {
if (!(nr % d)) {
return false;
}
}
return true;
}

void Search2PrimFacts(i64 Number) {
i64 a;

for(i64 i = 3; (i*i) <= Number; i += 2) {
if (Number % i == 0) {
if (IsPrime(i)) {
a = Number / i;
if (IsPrime(a)) {
printf("%I64d * %I64d = %I64d\n", a, i, Number);
break;
}
}
}
}
}

int main() {
time_t t1 = time(NULL);
Search2PrimFacts(1949565409282837);
cout << "t[s] = " << time(NULL) - t1;
}

[Sargnagel]

Original geschrieben von D'Espice
Naja nein, nicht wirklich. Das ziehen einer Wurzel ist einfach extrem aufwändig.

Jo, hab's gerade gefunden:
http://www.azillionmonkeys.com/qed/sqroot.html

Eine einfache Quadrierung sollte in der Tat schneller sein, als die optimierteste Assembler-Version der Quadratwurzelberechnung von obiger Seite. :(

[mj]

Original geschrieben von [ab]noname
Hrhr der Modulo Trick ;D

Jetzt hab ichs auf t = 0 bis 1sec für die 1. Zahl - Yippieeh! thx D'Espice
Na siehste, geht doch ;D
Jetzt muss ich wieder kontern, meine 7s wurden unterboten. Leider wird das wohl kaum möglich sein, da BigInteger einfach extrem lahm ist. Aber ich werd heut Nacht überlegen, ob sich da noch was drehen und optimieren lässt.
Was ich mir überlegt habe ist ein Feld anzulegen, welches alle Zahlen enthält die bereits als negativ geprüft wurden. Für jede neue Zahl kann dann zunächst geprüft werden, ob sie ein Vielfaches einer der bereits gespeicherten Zahlen ist. So kann man sich nochmals einige hundertausend Rechnungen sparen.

[[ab]noname]

Wie kann man jetzt genauer messen, weil ich bin unter der Messgenauigkeit mit meinem aktuellen C++ Algorithmus?

[mj]

Unter Linux kann ich das genau messen, ich werd gleich versuchen ob ich deinen Quellcode mit gcc compilieren kann und wenn ja, sag ich dir die genaue Zeit auf dem Opteron für den Algorithmus ;)

EDIT:
Mir fällt grad ein, wird eh nix bringen. Denn gcc kann mit __int64 nix anfangen und long ist zu klein für diese Zahl. Vergessen wir die Linux-Idee also gleich wieder ;)

[[ab]noname]

Wird Zeit, sich mal damit zu befassen: http://www.swox.com/gmp/

Würde das ganze lieber unter dem GCC laufen lassen, da man so mehr Systeme fürs Benchmarking nutzen kann.

[Sargnagel]

Original geschrieben von D'Espice
EDIT:
Mir fällt grad ein, wird eh nix bringen. Denn gcc kann mit __int64 nix anfangen und long ist zu klein für diese Zahl. Vergessen wir die Linux-Idee also gleich wieder ;)
Kein Problem:


long long number;
[...]
Search2PrimFacts(1949565409282837LL);


@GMP-Library:
Das ist ganz easy und die Library ist auch schön schnell. :)

[Sargnagel]

So, hier ist meine Version mit der GMP-Library. Da ist noch nix entsprechend optimiert. Es müßte noch die Anzahl der Aufrufe der mpz_init()-Funktion in der IsPrime()-Funktion drastisch verringert werden und was weiß ich noch alles, um mehr Performance zu erreichen, aber dafür habe ich leider keine Zeit mehr.

#include <stdio.h>
#include <stdbool.h>
#include <time.h>
#include <math.h>
#include <gmp.h>

bool IsPrime(mpz_t nr) {
mpz_t d, tmp;

mpz_init(d);
mpz_init(tmp);
mpz_set_ui(d, 2);
mpz_pow_ui(tmp, d, 2);
//printf("IsPrime\n");
while(mpz_cmp(tmp, nr) <= 0) {
//gmp_printf("d=%Zd, tmp=%Zd\n", d, tmp);
if(mpz_divisible_p(nr, d))
return false;
mpz_add_ui(d, d, 1);
mpz_pow_ui(tmp, d, 2);
}
return true;
}

void Search2PrimeFacts(mpz_t number) {
mpz_t a, i, tmp;

mpz_init(a);
mpz_init(i);
mpz_init(tmp);
mpz_set_ui(i, 3);
mpz_sqrt(tmp, number);

while(mpz_cmp(i, tmp) <= 0) {
if(mpz_divisible_p(number, i)) {
if(IsPrime(i)) {
mpz_cdiv_q(a, number, i);
if(IsPrime(a)) {
gmp_printf("%Zd * %Zd = %Zd\n", a, i, number);
break;
}
}
}
mpz_add_ui(i, i, 2);
}
}

int main(void) {
float millisecs;
clock_t start, end;

mpz_t number;

mpz_init(number);

mpz_set_str(number, "312133288776418481", 10);

start = clock();
Search2PrimeFacts(number);
end = clock();
millisecs = (float)(end - start)/((float)CLOCKS_PER_SEC/1000);
printf("Elapsed time in milliseconds: %f\n", millisecs);

return 0;
}

Ein paar print-Befehle zum Debuggen habe ich mal nicht gelöscht.

Viel Spaß damit. ;)


Ich habe mal die Search2PrimeFacts() Funktion leicht abgewandelt in Anlehnung an BoMbY's Beweis. Die Berechnungsdauer für "312133288776418481" wurde um etwas mehr als die Hälfte gesenkt.

[DarkAvenger]

Ich hoffe, es ist bekannt, daß Primzahltest in NP (und sogar coNP) liegt und die gängigen Alg auch nicht darüberhinaus kommen? Ein Bekannter hat mir zwar gesagt, daß es mittlerweile einen O(N) Algorithmus gibt, aber der soll wohl erst bei Zahlen mit einigen Dutzend (oder waren es hundert?) Nullen, schneller sein, als die bisher bekannten...aufgrund einer groooooßen Konstate vor dem "n".

Ist es nicht stupide, irgendwelche Algorithmen zu implementieren, ohne sich schlau zu machen, was es für Alg gibt bzw. ohne großartig Theorie zu verstehen, dumm draufloszucoden? Etwa das Mersenne Prime Projekt liefert interessante und vor allem schnelle Algorithmen....da zu optimieren wäre interessanter und auch von bleibenden Wert.

Naja, mehr als Beschäftigungstherapie wird das so nicht werden... :]

[[ab]noname]

Na dann auf und lass uns was sehen... ;)

[BoMbY]

Original geschrieben von DarkAvenger
[...]Naja, mehr als Beschäftigungstherapie wird das so nicht werden... :]

Morgen,

wenn man bei den Basics anfängt und sich dann langsam zu den guten bekannten Algorithmen vorarbeitet, kann man vieleicht sogar verstehen was die machen.

Also ich finde den Ansatz hier nicht schlecht, klar gibt's gerade auf diesem Gebiet tausende die schon 100x weiter sind - das Thema ist ja nicht gerade neu ...

m.f.g.
BoMbY

[DarkAvenger]

Ich fände es eher interessaner, wenn man Gebiete betrachten würde, wo noch wenig Gescheites bei rum gekommen ist. Da könnte evtl dem ein oder anderen doch etwas interessanteres gelingen. Etwa beim SAT Problem. Auch wenn das ein oft untersuchtes Problem ist, so gibt es dazu keinen effizienten Algorithmus (zumal es NP vollständig ist ;)).

[i_hasser]

Original geschrieben von DarkAvenger
Ist es nicht stupide, irgendwelche Algorithmen zu implementieren, ohne sich schlau zu machen, was es für Alg gibt bzw. ohne großartig Theorie zu verstehen, dumm draufloszucoden? Etwa das Mersenne Prime Projekt liefert interessante und vor allem schnelle Algorithmen....da zu optimieren wäre interessanter und auch von bleibenden Wert.

Naja, mehr als Beschäftigungstherapie wird das so nicht werden... :]

Tut mir ja leid, aber nicht jeder hier ist Prof. für Mathe und zufälligerweise auch noch für die richtigen Spezialbereiche. Selbst denken macht schlau, klar kann man sich die wirklich guten Algs angucken und verstehen - aber dahinter steckt so extrem viel geistige Leistung von so vielen Leuten die auf dem Gebiet wirklich Ahnung haben... da finden wir mal eben beim P3D Prog-Wettbewerb nicht einen schnelleren Algorithmus, und wenn doch wäre das mit einem Lotto-Gewinn gleichzusetzen.

... jeden Tag zerbrechen sich über genau diese Algorithmen mehr Leute den Kopf als es aktive Member im ganzen Progg-Forum gibt.

[DarkAvenger]

Richtig, genau deshlab sollte man sich ein Gebiet aussuchen, wo noch nicht so viel gearbeitet wurde, oder was von der Aufgabenstellung komplexer ist. Etwa im Bereich KI oder Bilderkennung oder so. Das würde mehr Sinn machen, da man eher so brauchbaren Ergebnissen käme. (Natürlich wird da auch viel geforscht, aber letztendendlich ist man da nicht sehr weit...) Alles andere ist -wie gesagt- nur Beschäftigungstherapie, bzw man lernt evt seinen compiler besser kennen...

Wie schon gesagt, falls Interess besteht, kann ich das SAT Problem mal formulieren und auch eine Beispiel Implementation mitliefern (so kann jeder seinen Alg gegenchecken) und auch einen Problemgenerator dafür. Das Problem ist recht einfach formuliert, aber sehr schwer gelöst (die Schritte sind einfach, aber man muß sehr sehr viele machen)...

Wer Primzahlen lösen ohne FFT anfängt, kann es auch gleich sein lassen...

[Sargnagel]

Original geschrieben von DarkAvenger
Richtig, genau deshlab sollte man sich ein Gebiet aussuchen, wo noch nicht so viel gearbeitet wurde, oder was von der Aufgabenstellung komplexer ist.
Ähem ... wir machen hier keine Grundlagenforschung. ;) Und wir wollen die Ergebnisse dieser Wettbewerbe auch sicherlich nicht an Nature oder Science zur Veröffentlichung schicken ... zumal sie schon öffentlich sind ... *chatt*
Wir wollen lediglich unsere kleinen Grauen Zellen nicht einrosten lassen. Und, wie intel_hasser schon schrieb, sind wir keine Mathe-Profs.
Original geschrieben von DarkAvenger
Wie schon gesagt, falls Interess besteht, kann ich das SAT Problem mal formulieren und auch eine Beispiel Implementation mitliefern (so kann jeder seinen Alg gegenchecken) und auch einen Problemgenerator dafür. Das Problem ist recht einfach formuliert, aber sehr schwer gelöst (die Schritte sind einfach, aber man muß sehr sehr viele machen)...

Mach einfach einen neuen Programmierwettbewerbthread auf. Da werden sich sicherlich Leute finden, die das interessiert.
Original geschrieben von DarkAvenger

Wer Primzahlen lösen ohne FFT anfängt, kann es auch gleich sein lassen...
Na, dann rate mal, womit die GMP-Library intern rechnet ... :]

[DarkAvenger]

Ach, was solls, ich werde heir einfach mal da SAT Problem formulieren, evtl wird das jemanden motivieren etwas in der Richtung zu probieren:


SAT=satifyability, Es geht um das Problem, ob eine boolsche FOrmel erfüllbar ist. Man betrachtet die sog. KNF (konjunktive Normalform; ist nicht weiter wichtig... ;))

Eine SAT Formel ist folgendermaßen aufgebaut:

Man hat "Literale". Das sind im Prinzip boolsche Variablen, die entweder 1 oder 0 sein können. (hier xnum, wobei num eine natürlich Zahl).

EIn negiertes Literal wird durch ein "-" angedeutet.

Die Literale werden per "OR" zu "Klauseln" verbunden, hier wird das "OR" mit "+" abgeküzt, Kaluseln werden durch Klammern abgeschlossen.

Die Klauseln wiederum werden mit "AND" verknüpft (hier "*"). Fertig ist die SAT Formel.

Beispiele:

(x1) offenabr erfüllbar, wenn x1=1

(x1)*(x2) dito, wenn x1=1 und x2=1

(-x1)*(x2) dito, wenn x1=0 und x2=1


Heir reicht erst einmal die Betrachtung von 3-SAT, dh. die Klauseln bestehen aus 3 Litaralen.

Beispiele:

(x1+x2+ -x2), IMMER erfüllbar (warum? x2 kann nur 0 oder 1 sein, aber x2 und -x2 kommen in der Klausel vor, also ist diese immer erfüllt wegen der OR Verknüpfung zwischen den Literalen)

(x1+x2+x2)*(-x1+x2+x2)*(x1+ -x2+ -x2)*(-x1+ -x2 -x2) NICHT erfüllbar! (Man probiere alle Kombinationen durch: x1=0, x2=0, x1=1 x2=0, x1=0 x2=1, x1=1 x2=1)


Das wäre auch schon das Problem! Man muß also eine erfüllbar "Belegung" der Literale finden, um zu zeigen, daß die (3-)SAT Formel erfüllbar ist.

In der theor. Informatik wurde gezeigt, daß 3-SAT schon NP-vollständig ist, dh. es wird mit großer Wahrscheinlichkeit keinen effizienten Alg geben.


Um also zu zeigen, daß eine SAT Formel erfüllbar ist, muß man eine Belegung finden, dh. mit geschickter Wahl des Alg kann man versuchen die (eine) richtige Belegung möglichst schnell zu raten. Nur i.a. wird man da schon sehr sehr viele Konbinationen prüfen müssen, da man (soweit bekannt) es an der Formel nicht direkt ansehen kann, wie denn die erfüllende Belegung auszusehen hat.

Um dagegen zu zeigen, daß eine SAT Formel unerfüllbar ist (streng genommen das coSAT Problem), muß man zeigen, daß ALLE Belegungen die Formel nicht erfüllen. Wenn wir n verschiedene Literale haben, gibt es 2^n Kombinationen, wodurch schon klar sein werden sollte, warum i.a. die Laufzeit so mies ist.


Wer sich also hieran die Zähne ausbeißen möchte, viel Vergnügen! Wer also die entsprechenden Formelgenratoren und einen Beispielalg haben möchte, so möge der hier posten oder mir eine email schreiben.


noch eine Anmerkung: Die Laufzeiten bei 3-SAT sind nach experimentellen Resultaten am längesten, wenn man etwa #Klauseln=4,25* #Literale wählt.


Warum SAT so wichtig ist? Ihr erinnert auch bestimmt noch an den Pentium FPU bug? Eine CPU/FPU ist ja letztendlich auch nur einen boolsche Formel. Nur hatte man damals nicht die Rechenpower zu überprüfen, ob jenes Formel auch wirklich unerfüllbar ist (das wäre zu zeigen, wenn ich es richtig verstanden habe). Das wurde nicht gemacht und deshalb gab es so lustige Probleme mit der FPU....

Um die Rechenleistung zu verdeutlichen: (mit worst-case #Klauseln): Auf meinen Athlon-XP@2,2 GHZ, braucht man für eine 100 Literale Formel ein paar Sekunden um die Unerfüllbarkeit zu zeigen. Für eine 500 Literale Formel schon fast eine Stunde...

Dann kann man sich ausmalen, wie übel es sein wird die boolsche Formel einer CPU durchzutesten....


Wrum ich denke, daß sich bei diesem Problem eher was sinnvolles ergibt? s.o.: Man kann durch entsprechende "Heuristik" versuchen den ALg intelligent zu steuern um die erfüllende Belegung zu finden. Hierbei könnte auch einem Nicht-Wissenschaftler etwas gutes gelingen. Um Unerfüllbarkeit schnell nachzuweisen, müßte man schon andere Betrachtungen anstellen...

[mj]

@DarkAvenger:
Bool und Logik erklärt in einem Beitrag - reife Leistung ;)
Prinzipiell hast du ja recht, du kannst das 3-SAT Problem aber auf ein 2-SAT reduzieren indem du Abfragen vorschaltest. So geschehen hab ich mein 3-SAT Problem auf drei 1-SAT Probleme zurückgeführt und statt 2min nur noch 6s Laufzeit.

Solcher Tricks bedarf es um schnell zu sein, leider bin ich jetzt mit JAVA bereits so ziemlich am Ende angelangt. Das Problem ist, dass die BigInteger Arithmetik extrem lahm ist, Optimierungsansätze kann ich jetzt nur noch bei Object.isProbablePrime() und der eigentlichen Arithmetik (speziell Modulo) finden.

[DarkAvenger]

@Sargnagel

Naja, GMP habe ich mir nicht direkt angesehen, aber macht Sinn es zu benutzen. ;)

Hmm, jetzt habe ich das hier gepostet. Egal, kann man ja immer noch splitten, zumal ich meine Bedenken habe, daß jemand sich dem Problem annehmen wird, zumal mich wohl einige für nen Klugscheißer oder ähnl halten werden.

Ich finde es nur schade, wenn fähige Programmierer ohne große theor. Gedanken zu machen in die Tasten hauen. zumal ich denke, daß die Anwesenden nicht blöd sind, könnten Sie halt sinnvoller ihre gruen Zellen einsetzen. ;)

[DarkAvenger]

@D.Espice

Leider kannst du 3-Sat nciht so ohne weiteres auf 2 SAT zurückführen. Denn 2-sat ist aus P (es gibt eff. Alg). Du wirst im Endeffekt exponentiel viele "Zurückführungen" haben, darum ist ein effiz. Alg weit weit weg. ;)

Ich habe hier ganz einfache Formeln angegeben. bei komplexeren dürfte es klar werden, was ich meine.

[Sargnagel]

Original geschrieben von DarkAvenger
@Sargnagel

Naja, GMP habe ich mir nicht direkt angesehen, aber macht Sinn es zu benutzen. ;)

Yep! Nur leider ist der Performanceeinbruch gewaltig. Da muß noch einiges optimiert werden.
Original geschrieben von DarkAvenger

Hmm, jetzt habe ich das hier gepostet. Egal, kann man ja immer noch splitten, zumal ich meine Bedenken habe, daß jemand sich dem Problem annehmen wird, zumal mich wohl einige für nen Klugscheißer oder ähnl halten werden.
Nimm's als Kompliment. ;D
Also, ich finde das SAT-Problem schon interessant, aber ich habe dafür überhaupt keine Zeit. Obigen Algorithmus habe ich auch nur kurz für die GMP-Lib umgeschrieben ... :(
Original geschrieben von DarkAvenger

Ich finde es nur schade, wenn fähige Programmierer ohne große theor. Gedanken zu machen in die Tasten hauen. zumal ich denke, daß die Anwesenden nicht blöd sind, könnten Sie halt sinnvoller ihre gruen Zellen einsetzen. ;)
Sicherlich könnte man "wichtigere" Probleme versuchen zu lösen, aber dafür fehlen mir und den meisten anderen wohl entscheidende Kenntnisse der höheren Mathematik und vor allem die nötige Zeit! Ich bin "nur" Zellbiologe und Zeit wächst nicht auf Bäumen ... aber wer weiß, was eine Synthese aus Quantenphysik und Gentechnik zu Leisten im Stande ist. *chatt*

Und nun zurück zu den Primzahlen.

[DarkAvenger]

Bzw. um mal konkret zu verdeutlichen (im Prinzip arbeiten fast alle Sat solver so, auch der "fast parallel sat solver" von Max Böhm, den ich als Bspalg hätte, allerding als seqentielle Variante):

Du hast deine SAT Formel und generierst 2 neue SAT Formeln (i.a. auch 3-SAT, auch wenn einige 2-SAT Klausln vorkommen!!!) indem du dir ein Literal aussuchst (Heuristik!) und es auf 0 oder 1 setzt. Mit den beiden neuen Formeln verfährt man rekursiv fort (-> exponentiell viele Formeln). Es reicht ja, wenn eine der beiden Formeln erfüllbar (wenn man weiß, daß die Formel erfüllbar ist...ansonsten muß man immer beide testen)ist, dh. (wieder Heuristik) im seq. Fall sollte man hier aussuchen, welche man zuerst testet. So könnte man versuchen nahe polynomieller Laufzeit die richtige Belegung zu finden.

Der fast parallel sat solver wählt nun so, daß die Teilprobleme möglichst gleich groß sind und jedes Teilproblem klein wie möglcih (-> Literal auswählen, was am häufigsten vorkommt). Diese Heuristik ist für den parallelen Einsatz optimiert, dh. eines der Teilprobleme wird an einen anderen Rechner geschickt und dort weitergetestet. Dh. im seq. Fall ließe sich der fast parallel sat solver durch geschicktere Wahl unterbieten.

Dieser Alg ist nun schon fast 10 Jahre alt, aber war Jahre lang ungeschlagen. Und man sieht, daß man trotz einfacher Heurtistik (mehr ist es eigentlich wirklich nicht! Der Kniff sind noch die Datenstrukturen im Programm selbst) diverse Wettbewerbe gewinnen kann.

Bsp für das Aufteilen:

Formel: (x1+x2+x3)*(-x1+x2-x3)

Wähle x1:

x1=1: (1+x2+x3)*(0+x2-x3)=(x2-x3)
Wahl x2:
x2=1: (1-x3)=1 -> Belegung gefunden!
x2=0: (-x3) => x3=0 -> Belegung gefunden!

(andere Formel: )
x1=0: (0+x2+x3)*(1+x2-x3)=(x2+x3)
Wahl x2:
x2=1: (1+x3)=1 -> Belegung gefunden!
x2=0: (x3) => x3=1 -> Belegung gefunden!


Wenn man das sich als Baum aufmalt, sieht man zum einen, daß dieser exponentiell wächst (->Binärbaum) um man findet die gültigen Belegungen indem man den Baum "travesiert" von der Wurzel (Ursprungsformel) zu einem Blatt, wo zum Schluß eine 1 steht. Im obigen Fall gibt es nur 2 Fälle, welche die Formel nicht erfüllen, nämlich x1=1,x2=0,x3=1 und x1=0,x2=0,x3=0.

So ein Haufen Theorie für heute. ;)

[DarkAvenger]

@Sargnagel

Das mit GMP hat mit der Größe des Problems zu tun. z.b. O(N) ist nicht immer schneller als O(N^2), es hängt vom Vorfaktor ab. Um den "Killeralg" zu schreiben, müßte man verzweigen:

etwa Alg1 hat 2*N^2, Alg2 hat 10000*N an Laufzeit.

Man benchmark seine Kiste und berechnet daraus, bis zu welchem N Alg1 verwendet wird und danach ALg2. Das ist der Unterschied zwischen theor. Betrachtungen und Implementationen auf Rechnern. ;) Euer einfacher Alg ist schnell für "kleine" Zahlen, nur bei Zahlen mit einigen Dutzend Dezimalstellen, wird FFT nötig. Mathematich ist der FFT Alg "fast immer" (wegen Limesbilding) schneller, nur auf dem Computer interessiert die Unendlichkeit nicht so sehr, darum muß man etwas Engineering betreiben. ;)

[Shaft99]

So hab jetzt auch mal was fertig:
C++ und mit GMP

#include <iostream>
#include <gmpxx.h>

int main() {
mpz_class Zahl;
Zahl = "1949565409282837";
mpz_t Wurzel;
mpz_init(Wurzel);
mpz_sqrt(Wurzel,Zahl.get_mpz_t());
for (mpz_class div(Wurzel);div>1;--div)
{
if (Zahl % div == 0)
{
mpz_class division(Zahl/div);
if ( ((mpz_probab_prime_p(div.get_mpz_t(),10)==1) || (mpz_probab_prime_p(div.get_mpz_t(),10)==2)) && ((mpz_probab_prime_p(division.get_mpz_t(),10)==1) || (mpz_probab_prime_p(division.get_mpz_t(),10)==2)))
{
std::cout << div << "*" << Zahl/div << std::endl;
};
};
};
return 0;
}

Jetzt noch ne Frage:
Wieviel Sekunden dauert das jetzt?

$ time ./faktor
19827473*98326469

real 0m9.777s
user 0m9.304s
sys 0m0.026s


Und warum dauert das in real 9.7s und sys nur 0.026s ??
Wartet der auf den Speicher oder an was liegts?

[BoMbY]

Also ich glaube, die einzige Möglichkeit wirklich einen schnellen Algorithmus für die Primfaktorzeregung zu bekommen, wäre es schaffen einen der Primfaktoren näher einzugrenzen. Also mein sqrt(pq) gibt ja z.B. eine Obergrenze an, wenn man jetzt z.B. noch eine Berechnung finden würde, mit der man eine recht hohe Untergrenze angeben kann (für den gleichen Faktor), dann könnte man da echt was reißen ...

Aber leider ist mir bis jetzt nichts eingefallen um einen der beiden Faktoren besser zu umzingeln.

m.f.g.
BoMbY

[Sargnagel]

Ich habe mal meine Version mit GMP upgedatet (http://www.planet3dnow.de/vbulletin/showthread.php3?s=&postid=1620781#post1620781). Jetzt läufts viel schneller. :)

Hier mal die Zeiten auf meinem Athlon 1.33 GHz, 266MHz DDR-SDRAM:
98326469 * 19827473 = 1949565409282837
Elapsed time in milliseconds: 1151

783216463 * 398527487 = 312133288776418481
Elapsed time in milliseconds: 22822

3981527477 * 2783216473 = 11081452861688528621
Elapsed time in milliseconds: 158948

Original geschrieben von DarkAvenger
Das mit GMP hat mit der Größe des Problems zu tun. z.b. O(N) ist nicht immer schneller als O(N^2), es hängt vom Vorfaktor ab.
Jo, da hast Du recht. Ich nehme mal an, das die GMP schwer auf wirklich große Zahlen optimiert ist. Was wir hier bisher betrachten sind ja nur recht kleine Zahlen.

Original geschrieben von DarkAvenger
Mathematich ist der FFT Alg "fast immer" (wegen Limesbilding) schneller, nur auf dem Computer interessiert die Unendlichkeit nicht so sehr
*chatt* Wie war das noch gleich mit der Berechenbarkeit? *buck* ;D

Original geschrieben von BoMbY
Aber leider ist mir bis jetzt nichts eingefallen um einen der beiden Faktoren besser zu umzingeln.
Da gibt es doch sicherlich schon Literatur zu?!
Übrigends, Deinen Beweis habe ich mal ins Programm einfließen lassen. Bei großen Zahlen ist die einmalige Berechnung der Wurzel schneller, als immer wieder zu quadrieren.

[Shaft99]

Hab mir mal die Zeitmessung von Sargnagel "geklaut".

Mein neuer Code (C++ und GMP):


#include <iostream>
#include <gmpxx.h>
#include <time.h>

int main() {
float millisecs;
clock_t start, end;
mpz_class Zahl;
Zahl = "312133288776418481";
mpz_t Wurzel;
mpz_init(Wurzel);
start = clock();
mpz_sqrt(Wurzel,Zahl.get_mpz_t());
for (mpz_class div(Wurzel);div>1;--div)
{
if (mpz_divisible_p(Zahl.get_mpz_t(),div.get_mpz_t()))
{
mpz_class division(Zahl/div);
if ( (mpz_probab_prime_p(div.get_mpz_t(),10)!=0) && (mpz_probab_prime_p(division.get_mpz_t(),10)!=0))
{
std::cout << div << "*" << Zahl/div << " = " << Zahl << std::endl;
end = clock();
break;
};
};
};
millisecs = (float)(end - start)/((float)CLOCKS_PER_SEC/1000);
std::cout << "Elapsed time in milliseconds: " << millisecs << std::endl;
return 0;
}


Meine Ergebnisse:
398527487*783216463=312133288776418481
Elapsed time in milliseconds: 11720
Mit Sargnagel's Code:
783216463 * 398527487 = 312133288776418481
Elapsed time in milliseconds: 18590.000000

Auf nem 2600+ mit 2087 MHz und 512 MB
mit GCC 3.3.2 und Flags: -O3
Vielleicht holt der neue 3.4.0 noch nen paar ms ;D
gleich mal probieren

[mj]

:-/ Ohne vernünftige Mathe-Bibliothek kann ich in Java von sowas nur träumen :(

[Sargnagel]

Original geschrieben von Shaft99
Meine Ergebnisse:
398527487*783216463=312133288776418481
Elapsed time in milliseconds: 11720
Mit Sargnagel's Code:
783216463 * 398527487 = 312133288776418481
Elapsed time in milliseconds: 18590.000000

Sehr schön. Aber jetzt besteht das Programm nur noch aus vorgefertigten Funktionen. Sinn und Zweck dieser Übung ist es eigentlich, selber die nötigen Funktionen zu implementieren. Natürlich soll man das Rad nicht neu erfinden, aber bei dieser Aufgabe sollte zumindest die Funktion, die testet, ob eine Zahl eine Primzahl ist, selber implementiert worden sein, denn gerade das ist doch der Reiz an der Sache. *buck*

Und wie schnell ist das Programm bei 11081452861688528621 und 6864244079350478607287? Ich hatte meines mal mit der letzten Zahl gefüttert, aber nach 35 min Laufzeit hatte ich abgebrochen. :-X

@D'Espice:
Gibt es auch nix ordentliches bei sourceforge.net? Ich habe mit Google eine Lib von IBM gefunden:
http://www-124.ibm.com/developerworks/projects/mathlib


Hier ist noch eine große Übersicht zu Java-Math-Libs:
http://math.nist.gov/javanumerics/

[mj]

Ich hab mir die Bibliotheken mal durchgeschaut, aber leider sind alle, die halbwegs brauchbar wären (wie z.B. die Intel Math Kernel Library) nicht für Java. Naja, was solls. Langsam verliere ich eh Interesse, irgendwie war die Aufgabe zu leicht zu lösen, werd mich jetzt anderen Projekten widmen.

[Sonic]

is ja lustig. so ein programm muss ich für mein info praktikum auch coden, aber da isses egal wie langsam rechnet. hauptsache es tut...

also mein prof hat sich folgendes ausgedacht:


#include <iostream>

int main()
{
int n;
do
{
cout<<"Dieses Programm zerlegt eine ganze Zahl n in Primfaktoren.\n"
"Es muss n>1 gelten!!\n"
"Bitte geben Sie n ein: ";
cin>>n;
}
while(n<=1);
cout<<n<<" =";
bool multiply_symbol=false;
for(int faktor=2;n>1;faktor=faktor+(faktor%2+1))
{
int anz=0;
while(n%faktor==0)
{
n=n/faktor;
anz=anz+1;
}
if(anz>0)
{
cout<<" ";
if(multiply_symbol) cout<<"* ";
cout<<faktor<<"^"<<anz;
multiply_symbol=true;
}
}
cout<<endl;
}



vielleicht könnt ihr damit ja was anfangen... ich muss erstma weiter lernen...

[mj]

Sonic:
Du sollst hier schon Eigenentwicklungen posten, nicht etwas, was sich dein Professor ausgedacht hat ;)

[Sonic]

ja ja, ruhig an... erstma in die problemstellung einlesen...

sollte ja nur n kleiner anstoß sein. ihr seid eh schon weiter, glaub ich... aber ich hab noch ein interessantes aufgabenblatt zum thema "Monte Carlo Algorithmus" gefunden. das werd ich erst einmal durcharbeiten...

[DarkAvenger]

Monte Carlo ist probalistisch, ich weiß nicht, ob das der Aufgabenstellung entspräche... Einen RP Algorithmus wüßte ich auch...

[Shaft99]

Original geschrieben von Sargnagel
Sehr schön. Aber jetzt besteht das Programm nur noch aus vorgefertigten Funktionen. Sinn und Zweck dieser Übung ist es eigentlich, selber die nötigen Funktionen zu implementieren. Natürlich soll man das Rad nicht neu erfinden, aber bei dieser Aufgabe sollte zumindest die Funktion, die testet, ob eine Zahl eine Primzahl ist, selber implementiert worden sein, denn gerade das ist doch der Reiz an der Sache. *buck*

Und wie schnell ist das Programm bei 11081452861688528621 und 6864244079350478607287? Ich hatte meines mal mit der letzten Zahl gefüttert, aber nach 35 min Laufzeit hatte ich abgebrochen. :-X


Das stimmt, aber meine selbst geschriebene isprime() war einfach zu lahm für (richtig)
große Zahlen. Und selbst wenn ich wöllte, könnt ich keinen Alg. schreiben, der
für solche Zahlen auch nur halbwegs akzeptable Zeiten produziert. Kein Mathe-Student usw ...
Aber eigentlich mach ich hier ja mit, um meine Programmierkenntnisse zu erweitern und
da ich bisher recht wenig mit libs gemacht hab, war das für mich ne gute Übung.
Wir könnten ja ne separate Wertung für Code mit Libs machen ;D


2783216473*3981527477 = 11081452861688528621
Elapsed time in milliseconds: 40590

und bei 6864... rechnet er noch

[Sargnagel]

Original geschrieben von Shaft99
Das stimmt, aber meine selbst geschriebene isprime() war einfach zu lahm für (richtig)
große Zahlen. Und selbst wenn ich wöllte, könnt ich keinen Alg. schreiben, der
für solche Zahlen auch nur halbwegs akzeptable Zeiten produziert. Kein Mathe-Student usw ...

Tja, mir fällt leider auch nichts ein. Aber das spukt wohl die nächsten Wochen noch in meinem Kopfe herum ... *kopfkratz
Mathe-Student bin ich ebenfalls nicht.
Original geschrieben von Shaft99

Aber eigentlich mach ich hier ja mit, um meine Programmierkenntnisse zu erweitern und
da ich bisher recht wenig mit libs gemacht hab, war das für mich ne gute Übung.
Wir könnten ja ne separate Wertung für Code mit Libs machen ;D

Zum Üben ist das hier optimal. Und auch Libraries wollen richtig verwendet werden. Aber man kann Libraries bestimmt als Referenz verwenden, um die Qualität der eigenen Algorithmen zu testen.
Original geschrieben von Shaft99

2783216473*3981527477 = 11081452861688528621
Elapsed time in milliseconds: 40590

und bei 6864... rechnet er noch
Vielen Dank! Das sieht doch gut aus. Ich bin mal gespannt, ob der mit der 6864... noch in vertretbarer Zeit fertig wird ...

[EDIT]
Original geschrieben von DarkAvenger

Monte Carlo ist probalistisch, ich weiß nicht, ob das der Aufgabenstellung entspräche... Einen RP Algorithmus wüßte ich auch...
Hmmm ... da frage ich mich, ob man mit einem genetischen/evolutionären Algorithmus ebenfalls gute Zeiten erreichen könnte?

[Shaft99]

Original geschrieben von Sargnagel
Vielen Dank! Das sieht doch gut aus. Ich bin mal gespannt, ob der mit der 6864... noch in vertretbarer Zeit fertig wird ...

Heut jedenfalls nicht mehr, muss morgen früh aufstehen.
Ergebnisse dann morgen Abend.

[Marnem]

Ich hab mal ne dumme Frage:
Wie löst ihr das Problem mit Zahlen, die zwar keine Primzahlen sind, sich aber auch nicht in 2 Primzahlen zerlegen lassen? 27 ist so ein Beispiel. Es lässt sich in 3*9 zerlegen, aber 9 ist ja auch keine Primzahl?

[mj]

Original geschrieben von Marnem
Ich hab mal ne dumme Frage:
Wie löst ihr das Problem mit Zahlen, die zwar keine Primzahlen sind, sich aber auch nicht in 2 Primzahlen zerlegen lassen? 27 ist so ein Beispiel. Es lässt sich in 3*9 zerlegen, aber 9 ist ja auch keine Primzahl?
Naja, dann wird eben ausgegeben, dass keine Lösung gefunden werden kann ;)

[i_hasser]

Ha! - jetzt hab ich endlich auch mal Zeit ;D

Gute Ideen sind mir auch schon ein paar gekommen, also macht euch auf was gefasst 8)

PS Die 2. Kanne Tee ist auch schon unterwegs *chatt*

[Sargnagel]

Original geschrieben von intel_hasser
Ha! - jetzt hab ich endlich auch mal Zeit ;D

Gute Ideen sind mir auch schon ein paar gekommen, also macht euch auf was gefasst 8)

PS Die 2. Kanne Tee ist auch schon unterwegs *chatt*
Yes! Endlich wieder völlig glasklare bitweise Operationen! Her damit! :] *buck* *chatt*

[i_hasser]

So, der allererste alg steht.

4.1 Sekunden Initialisierung und 0.019 Sekunden Suchzeit ;D (hab erstmal mit der 1949565409282837 begonnen)

Der Alg sucht momentan nur einen kleineren Bereich von Primzahlen ab, findet also bestimmte Teiler nicht.
Momentan sucht er im Bereich 19.49e6 bis 100e6 - ideal wäre natürlich 2 bis sqrt(1949565409282837)~44.15e6.

Nichtsdestotrotz hat er im Bereich ein Primzahlenpaar gefunden, nämlich 19827473*98326469.

Damit der Alg den gesamten Bereich durchsucht muss ich mir noch was ausdenken, aber sooo extrem schwierig dürfte das net werden - nur ein bisschen langsamer.


@Sarnagel

Freu dich auf den Code, der ist wiedermal haarsträubend ohne Ende *chatt*

[i_hasser]

So, jetzt werden alle Bereiche durchsucht. 4.19 Sekunden Initialisierung nach wie vor und 1.3 Sekunden Suchzeit. Jetzt lass ich ihn mal nach allen möglichen Kombinationen suchen (bisher nur nach der 1.).


So, es gibt für diese Zahl nur die eine Zerlegung. Jetzt hat er aber immerhin schon 6 Sekunden zum Suchen gebraucht (kann noch brachial optimiert werden ;)).

[i_hasser]

*lol* der Code war mit -O0 compiliert :]

Momentan ist aber eine Ecke drinnen :(

[i_hasser]

ne, also manchmal zweifel ich echt an mir. Fehler gefunden - hab statt for(initialcode;bedingungzumweiterendurchlauf;zählercode;) ein for(initialcode;zählercode;bedingungzumweiterendurchlauf) geschrieben, und jetzt wo ich diesen Text schreibe hab ich es anfangs auch wieder vertauscht.

*chatt*

Das gibt natürlich schönste Bugs.


Also für die 1. Zahl dauerts jetzt 4.7 Sekunden. Ich lass ihn mal die 2. berechnen...


€: 45 Sekunden für die 2. Zahl, nach 20 Sekunden bin ich über die erste (und einzige) Möglichkeit gestolpert.


Update: So langsam bekommt der Alg seinen Schliff (muss ihn nur noch von C++ nach C portieren *chatt*).

Zahl 1: 3.7 Sekunden, Zahl 2 dauert komischerweise immernoch so lange (jetzt sogar 46 Sekunden) - da werd ich mal den Profiler drübergrasen lassen.


Update: Keine Ahnung wieso Zahl1 vorhin schneller war, die Optimierung hab ich nämlich mittendrinn versaut. Jetzt braucht er für Zahl2 noch 42 Sekunden... hätte mir eigentlich mehr davon versprochen.

Heute mach ich da nimmer weiter, gute nacht allerseits *chatt*

[Marnem]

So, ich versuch mich auch an der Aufgabe.
Die Versionen von [ab]noname kann ich nicht kompilieren (Linux, gcc) und daher nicht nachvollziehen, wie er auf diese Zeiten kommt.
Momentan arbeite ich mit der "primitivsten" Lösung, heißt ein Array mit den Primzahlen anzulegen und die zu testende Zahl an den Primzahlen auszutesten. Dabei nutze ich das Wissen, daß der größere der beiden Faktoren von 1949565409282837 98326469 ist, berechne also nur die Primzahlen bis 100000000.
Trotzdem bekomme ich Zeiten von ettlichen Minuten.
Warum ist [ab]nonames Version so viel schneller? Mit seiner IsPrime Methode ists sogar noch langsamer als mit meiner (meine nutzt die Tabelle mit den Primzahlen und muß daher deutlich weniger Zahlen testen).

Es gibt 5761454 Primzahlen zwischen 1 und 100000000
Die Zahl 1949565409282837 lässt sich in die folgenden Primfaktoren zerlegen: 19827473 * 98326469

real 9m39.457s
user 9m5.410s
sys 0m0.327s


#include <iostream>
#include <math.h>

using namespace std;

typedef long long i64;

int isPrime(i64 primes[],i64 pos,i64 x) {
for (i64 i = 0; i < pos; i++) { //für alle primzahlen in primes
i64 a = primes[i];
if ((a*a) > x) {break;}
else if (!(x % a)) {return 0;}
}
return 1;
}

i64 fillArrays(i64 primes[],i64 arraygroesse,i64 maxwert) {
i64 primepos = 0;
for (i64 i = 3; i <= maxwert;i=i+2) { //gerade Zahlen können nie Primzahlen, allerdings kann nicht sinnvoll mit geraden
//Zahlen getestet werden
if (isPrime(primes,primepos,i)) {
primes[primepos] = i;
primepos++;
if (primepos == arraygroesse) return -1;
}
if (!(primepos % 100000)) printf("Array gefuellt bis %llu mit %llu\n",primepos,primes[primepos -1]);
printf("Es gibt %llu Primzahlen zwischen 1 und %llu\n",primepos,maxwert);
return primepos;
}

i64 zerlegungsmoeglichkeiten(i64 zerlegungen[],int groesse,i64 primes[],i64 lastprimepos,i64 nonprime) {
i64 x = (nonprime +1) / 2;
i64 y = 0;
int pos = 0;
for (i64 i = 0; (i <= lastprimepos) && ((y=primes[i]) < x); i++) {
if (nonprime % y) {}
else {
i64 z = nonprime / y;
if ((y <= z) && isPrime(primes,lastprimepos,z)) {
zerlegungen[pos] = y;
zerlegungen[pos+1] = (nonprime / y);
pos=pos+2;
if (pos >= (groesse - 1)) return -1;
}
}
}
return pos;
}

void printzerlegungen(i64 zerlegungen[], int lastpos, i64 nonprime) {
printf("Die Zahl %llu lässt sich in die folgenden Primfaktoren zerlegen: ", nonprime);
for (int i=0;i < lastpos;i= i+2) {
printf("%llu * %llu ",zerlegungen[i],zerlegungen[i+1]);
}
printf("\n");
}

int main() {
i64 arraygroesseprimes = 6000000;
i64 testnr = 1949565409282837ll; //zu testende Nummer, das ll nicht vergessen!! nicht mit geraden Zahlen testen !!
i64 maxwert = (100000000); //bis zu diesem Wert werden alle Primzahlen berechnet
i64 primes[arraygroesseprimes];
i64 lastprimepos = 0;

lastprimepos = fillArrays(primes,arraygroesseprimes,maxwert); //ohne "nonprime"
if (lastprimepos == -1) {printf("Array zu klein\n"); return 0;}
lastprimepos--;
i64 zerlegungen[100];
int lastpos = zerlegungsmoeglichkeiten(zerlegungen,100,primes, lastprimepos,testnr);
if (lastpos == -1) {
printf("Array zerlegungen zu klein\n"); // wenn es mehr als 50 zerlegungen gibt
return 0;
}
if (lastpos != 0) {
lastpos--;
printzerlegungen(zerlegungen,lastpos,testnr);
}
else printf("Die Zahl %i lässt sich nicht in Primfaktoren zerlegen.\n",testnr);
return 0;
}

[i_hasser]

#include <iostream>
#include <time.h>
#include <math.h>

using namespace std;
typedef unsigned long long int i64;

bool IsPrime(i64 nr) {
for(i64 d = 2; d*d <= nr; d++) {
if (!(nr % d)) {
return false;
}
}
return true;
}

void Search2PrimFacts(i64 Number) {
i64 a;

for(i64 i = 3; (i*i) <= Number; i += 2) {
if (Number % i == 0) {
if (IsPrime(i)) {
a = Number / i;
if (IsPrime(a)) {
printf("%Ld * %Ld = %Ld\n", a, i, Number);
break;
}
}
}
}
}

int main() {
clock_t t1 = clock();
Search2PrimFacts(1949565409282837);
cout << "t[s] = " << ((float)(clock() - t1))/CLOCKS_PER_SEC;
}



So rennt das auch unter Linux. Allerdings zweifle ich irgendwie an dem Algorithmus, viel zu einfach ;)
Mal sehen wie der mit größeren Zahlen skaliert...



EDIT: Bin durchgestiegen, ist vom Ansatz her ähnlich wie meiner. Allerdings dürfte der bei großen Zahlen deutlich schlechter skalieren ;)

EDIT2: Bei der 2. Zahl braucht er ~11 Sekunden. Allerdings sind das schon alles besondere Zahlen, weil die genau eine Zerlegung in sehr große Primzahlen haben ;)
Ich werd mal ein paar andere Zahlen durch meinen Alg laufen lassen und die Ergebnisse und Teiler hier posten.

[Shaft99]

So neue Ergebnisse:
(Aber alles GMP-Funktionen, sonst wirds viel zu langsam)

73981527493*92783216459 = 6864244079350478607287
Elapsed time in milliseconds: 1113390

Eigentlich schneller als ich dachte, gestern hab ich die gleiche Zahl
nach 80 ! Minuten abgebrochen.

[Sargnagel]

Original geschrieben von Shaft99
Eigentlich schneller als ich dachte, gestern hab ich die gleiche Zahl
nach 80 ! Minuten abgebrochen.
Danke für die Messung! Da bin ich froh, daß ich die Berechnung bei mir rechtzeitig abgebrochen hatte. Wer weiß, wie lange das noch hätte dauern können ... :P :-X

@intel_hasser:
Zeig doch mal Deinen Quellcode oder ist er Dir noch nicht kryptisch genug? *buck*

[DarkAvenger]

#include <iostream>
#include <time.h>
#include <math.h>

using namespace std;



Igitt, sollte das nicht eher


#include <iostream>
#include <ctime>
#include <cmath>

using namespace std;


sein? :)

[Shaft99]

Ich frag mich auch immer bei:

#include <iostream>

printf("...");


Aber jedem das seine *zweifel*

[i_hasser]

Hab das Programm jetzt auf C portiert - jetzt is'ses kryptisch genug ;D

Mir ist inzwischen eine neue Idee gekommen, die deutlich schneller sein sollte - da mach ich jetzt weiter. Den Alg gibts erstmal weitestgehend unkommentiert:

main.c

#ifdef HAVE_CONFIG_H
#include <config.h>
#endif

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>

#include <largenum.h>
#include <primeapi.h>

int
main(int argc, char *argv[])
{
// ---- Predefines ----

// Timing Vars
clock_t time1, time2;
clock_t timei, timei1;
clock_t timee, timee1;

// - the value to test -
//int64u alg_testval=1949565409282837;
int64u alg_testval=1949565409282838;
//int64u alg_testval=312133288776418481;
//int64u alg_testval=11081452861688528621;

// ---- Start Initialisation: ----
printf("\ninitialization...");
fflush(stdout);

int64u alg_sqrttestval=(int64u) sqrt( (double) alg_testval) +1;

// search bounds - TODO: Chg this
int64u alg_searchframe_high=100000000;
int64u alg_searchframe_low=2;

// Init Prime API
time1=clock();
prime_init(alg_searchframe_high);
time2=clock();
timei1=time2-time1;

timei=timei1;

// ---- Check Number ----
printf("\nchecking number...");
time1=clock();

// print some statistics
printf("\n search area: %f to %f", (double) alg_searchframe_low, (double) alg_sqrttestval);
fflush(stdout);

// start search: for every prime in searchrange do...
int64u alg_count, alg_count_old;
int32u i,rep_count=0;

for(alg_count=alg_searchframe_low;alg_count=prime_nextprime(alg_count);alg_count!=0) {
if(alg_count>alg_sqrttestval) break;

rep_count++;

// mkmodulo to determine if devisor fits...
if( alg_testval%((int64u) alg_count)==0) {
// if it does, chk if div eval is prime
i=prime_isprime(alg_testval/((int64u) alg_count));

if(i != PRIMECHK_NO_PRIME)
{
if(i==PRIMECHK_UNKNOW) {
printf("\nunknown state of prime %Ld", alg_testval);
exit(1);
}

printf("\n devisor found: %Ld * %Ld", alg_count, alg_testval/((int64u) alg_count));
}
}

if( (rep_count&0x3FFFF)==0) {
printf("\n %Ld of %Ld", alg_count, alg_sqrttestval);
fflush(stdout);
}

alg_count_old=alg_count;
}
if(alg_count>alg_sqrttestval) printf("\nsearchend reached");
if(alg_count==0) printf("\nout of prime: %Ld", alg_count_old);

time2=clock();
timee1=time2-time1;

timee=timee1;

// ---- Print Closeup Statistics ----

printf("\n---- time statistics ----");
printf("\ninitialisation time: %.3f sec",(float) timei/CLOCKS_PER_SEC);
printf("\n prime api: %.3f sec",(float) timei1/CLOCKS_PER_SEC);
printf("\nexecution time: %.3f secs", (float) timee/CLOCKS_PER_SEC);
printf("\n check devisors: %.3f secs", (float) timee1/CLOCKS_PER_SEC);


printf("\n");
return 0;
}


largenum.h

#ifndef _LARGENUM_H_
#define _LARGENUM_H_

// Only some typedefs
typedef unsigned long long int int64u;
typedef unsigned int int32u;
typedef unsigned short int int16u;
typedef unsigned char int8u;


#endif


primeapi.h

#ifndef _PRIMEAPI_H_
#define _PRIMEAPI_H_

#include <math.h>
#include <largenum.h>

// Constants
static const int32u PRIMECHK_IS_PRIME=1;
static const int32u PRIMECHK_NO_PRIME=2;
static const int32u PRIMECHK_UNKNOW=3;


// ---- Global Vars ----

// Prime Matrix
int8u *prime_mat; // Prime Matrix Field
int64u prime_max; // Size of Prime Matrix

// Helpers for calcs
int32u prime_bit_sel[120]; // Bit Selector
int32u prime_mod_tab1[256]; // Modulo Table
int32u prime_mod_tab2[256]; // Modulo Table


// ---- Functions ----

// Init Functions
void prime_init(int64u upperbound);
void prime_init_mat(int64u upperbound);

// Data Access Functions

int32u prime_isprime(int64u number);
int32u prime_islongprime(int64u number);
int64u prime_nextprime(int64u startnumber);

#endif


primeapi.c

#include <primeapi.h>
#include <math.h>

// Init Function
// * initializes the whole prime api
void
prime_init(int64u upperbound)
{
prime_max=upperbound;
prime_init_mat(prime_max);


}

// Init Function for Prime Matrix
// * inits the prime_mat
void
prime_init_mat(int64u upperbound)
{
// Bench, Timing and Steering Values
unsigned int search_to=upperbound;

// Sieb-Matrix
prime_mat=(unsigned char *) malloc((search_to/30+1)*sizeof(unsigned char));

long int PtrCur;
long int i,j;

// other variables
long int PrimeCount=search_to*4/15;
long int SearchEnd=(unsigned long int) sqrt((double) search_to);
long int SearchVal=7;
long int NewSearchValFound=0;

long int ta,td;
long int tc,tb,te;

// PreCheck Primes
for(i=0;i<120;i++) prime_bit_sel[i]=10;

for(i=0;i<4;i++) {
prime_bit_sel[1+i*30]=0;
prime_bit_sel[7+i*30]=1;
prime_bit_sel[11+i*30]=2;
prime_bit_sel[13+i*30]=3;
prime_bit_sel[17+i*30]=4;
prime_bit_sel[19+i*30]=5;
prime_bit_sel[23+i*30]=6;
prime_bit_sel[29+i*30]=7;
}

for(i=0;i<256;i++) prime_mod_tab1[i]=i%30;
for(i=0;i<256;i++) prime_mod_tab2[i]=(i*256)%30;

// Run-initialization
for(i=0;i<=search_to/30;i++) prime_mat[i]=0xFF;
PrimeCount=search_to*4/15;
SearchEnd=(long int) sqrt((double) search_to);
SearchVal=7;
NewSearchValFound=0;

// Primzahlen filtern
while(SearchVal<=SearchEnd) {
// Test primes for SearchVal
for(PtrCur=SearchVal<<1;PtrCur<=search_to;PtrCur+=SearchVal) {
// this equals to modulo 30
te= prime_mod_tab1[PtrCur&0xFF]+prime_mod_tab2[(PtrCur&0xFF00)>>8]+
prime_mod_tab2[(PtrCur&0xFF0000)>>16]+prime_mod_tab2[PtrCur>>24];

if(prime_bit_sel[te]-10) {
// read out data for Value
tb=PtrCur/30;
td=1<<prime_bit_sel[te];
ta=prime_mat[tb];

if(ta&td) {
PrimeCount--;
prime_mat[tb]=ta^td;
}
}
}

// Get next SearchVal
NewSearchValFound=0;

for(i=SearchVal+1;i<=SearchEnd;i++) {
// this equals to mod 30
te= prime_mod_tab1[i&0xFF]+prime_mod_tab2[(i&0xFF00)>>8]+
prime_mod_tab2[(i&0xFF0000)>>16]+prime_mod_tab2[i>>24];

if(prime_bit_sel[te]-10) {
// read out data for Value
tb=i/30;
tc=prime_bit_sel[te];
ta=*(prime_mat+tb);

if((ta>>tc)&0x1) {
SearchVal=i;
NewSearchValFound=1;
break;
} // if
} // if
} // for

if(NewSearchValFound==0) break;

} // while

return;
}


// IsPrime
// * Checks if given Number is a prime and returns that
int32u
prime_isprime(int64u number)
{

// if out of range - perform a long range check -> SLOW
if(number > prime_max) return prime_islongprime(number);

// return for hardcoded primes
switch(number) {
case 2: return PRIMECHK_IS_PRIME;
case 3: return PRIMECHK_IS_PRIME;
case 5: return PRIMECHK_IS_PRIME;
}

// Check for hardcoded primes and 2n cases
if(number&1==0) return PRIMECHK_NO_PRIME;

// compute bitposition in prime_mat entry
int32u bitpos;
bitpos=prime_mod_tab1[number&0xFF]+prime_mod_tab2[(number&0xFF00)>>8]+prime_mod_tab2[(number&0xFF0000)>>16]+prime_mod_tab2[number>>24];

// bailout if excluded by hardwired primes (could be removed?)
if(prime_bit_sel[bitpos]==10) return PRIMECHK_NO_PRIME;

// compute bytepos
int32u bytepos; // the byte in prime_mat
bytepos=number/30;

// read data for that value out of prime_mat
int32u bitsel; // specs the absolute bit in given byte
int8u bitmask; // specs the bitmask to easy access prime_mat
int8u valdata; // dataset in prime_mat

bitmask=1<<prime_bit_sel[bitpos];
valdata=prime_mat[bytepos];
if(valdata&&bitmask) return PRIMECHK_IS_PRIME;
else return PRIMECHK_NO_PRIME;
}


int32u
prime_islongprime(int64u number)
{
int64u rest, count, sqrtnumber=sqrt( (double) number);
for(count=2;count<=sqrtnumber;count=prime_nextprime(count)) {
if(count==0) return PRIMECHK_UNKNOW;
if(number%count==0) return PRIMECHK_NO_PRIME;
}

return PRIMECHK_IS_PRIME;
}

int64u
prime_nextprime(int64u startnumber)
{
int64u i=startnumber;

// fillup this block
while(i<=startnumber+30)
{
i++;
if(prime_isprime(i)==PRIMECHK_IS_PRIME) return i;
}
startnumber=i;

// chk following blocks
if(startnumber<prime_max-30)
{
// jump over complete blocks
int32u bytepos;
for(bytepos=(int32u) startnumber/30; *(prime_mat+bytepos) == 0; bytepos++);
startnumber=30*bytepos;
}

// only do this if out of prime_mat
while(1)
{
i++;
if(prime_isprime(i)==PRIMECHK_IS_PRIME) return i;
}

return 0;
}



Theoretisch ist der Alg schneller als [ab]nonames. Das Problem ist nur, dass die Zahlen dafür deutlich zu klein sind... ok, ich werd den Alg mal auf die 3. Zahl erweitern - dann könnte er schon schneller sein. Bei der 2. Zahl braucht meiner ja 42 Sekunden, [ab]nonames braucht bei mir 11 Sekunden - das ist ca. 4:1, bei der 1. Zahl braucht meiner ca. 8mal so lange.
Allerdings muss ich für die 3. Zahl das alles noch etwas erweitern, und der Ram-Verbrauch dürfte auf 100mb ansteigen.

Wenn das mit meiner neuen Idee klappt gehts allerdings rund, das ist ein ähnlicher Ansatz wie ich ihn bei meinem Primzahlalgs verwirklichen wollte, was ich dann aber doch wieder verworfen hab (da zuviel Aufwand). Hier ließe sich das recht einfach implementieren.

[[ab]noname]

Original geschrieben von Shaft99
Ich frag mich auch immer bei:

#include <iostream>

printf("...");


Aber jedem das seine *zweifel*
Wenn der Meister mir ne Möglichkeit nennt über den Ausgabestream eine 64bit Integer Zahl auszugeben bin ich gerne bereit, auf printf zu verzichten. Bis dahin gilt Dieter Nuhr ;D

iostream ist bei mir aus gewohnheit in Consolen Programmen drin, da ich das meist für div. Zwischenausgaben brauche, die ich dann wieder rausnehme.

[i_hasser]

So, er initialisiert gerade für die 3. Zahl (dauert ein bisschen - wer durch den Alg schon durchgestiegen ist dürfte erkannt haben, dass alle Primzahlen bis ~sqrt(3. Zahl) berechnet werden müssen).

... so, 1/3 der Berechnung ist durch... in ca. 120 Sekunden ;D

... die Hälfte...

... 2/3 - interessant ist, dass der Algorithmus jetzt immer schneller wird (ist aber auch so gewollt)...

So - durch ist er ;D
Eine Darstellung hat er auch gefunden (habs mir aber net notiert, weil das so schnell vorbeigescrollt ist ;)).

73 Sekunden für die Initialisierung und 278 Sekunden für das Heraussuchen der Zahl.


Der Aufwand dürfte für kleinere Zahlen (so bei 2^64) bei O<sqrt(n) liegen. Je größer n wird, desto größer ist sqrt(n)-O (der Alg arbeitet also immer schneller). Dummerweise liegt der Speicherverbrauch bei sqrt(n)/30 Bytes. Aber bevor das zu viel wird ist der 64bit Integer zuende ;)

Mal eine kleine Gegenüberstellung:

n=1.949565409282837e15, t=4.7s, w1=t/n ~ 2.41e-15
n=3.12133288776418481e17, t=42s, w2=t/n ~ 1.35e-16
n=1.1081452861688528621e19, t=351s, w3=t/n ~ 3.168e-17

Also der Aufwand verläuft tatsache unter-linear. [ab]nonames alg hat nach 99.9 Sekunden das richtige Ergebniss für die 3. Zahl gefunden, da siehts so aus:

n=1.949565409282837e15, t=0.5s, w4=t/n ~ 2.56e-16
n=3.12133288776418481e17, t=11s, w5=t/n ~ 3.524e-17
n=1.1081452861688528621e19, t=99.9s, w6=t/n ~ 9.0e-18

Warum der auch nicht-linear verläuft weis ich jetzt aus dem Stehgreif zwar nicht, aber vor allem bei der letzten Zahl sieht man, dass der Aufwand nicht so stark abnimmt.

Noch ein paar Werte:

w1/w2 = 17.85
w2/w3 = 4.261

w4/w5 = 7.26
w5/w6 = 3.92

Wie man sieht nimmt bei meinem die Arbeit bei größeren Zahlen scheller ab (irgendwie muss ich meinen Klopper von Alg ja jetzt noch schönreden ;)).

[mj]

GAH! Wenn ich das sehe krieg ich auch wieder Lust, meinen Algorithmus zu optimieren :( Was ihr mich wieder an Zeit kostet, das ist abartig! *motz* *chatt*

[[ab]noname]

Hatte die letzten Tage eher wenig Zeit. Jetzt muss ich mich langsam mal wieder ransetzen.

Meine 1. Lösung war ja eher trivialer Natur und relativ intuitiv. Bei intel_hassers Lösung sieht man schon, dass es krasser geht ;D (Kleine Frage, bist du Student, in Ausbildung oder ist das rein Hobby?)

Mal schaun, was ich noch zu Stande bekomme...

[i_hasser]

das ist abartig! ;)
das ist abartig! *motz* *chatt*

Starke Stimmungsschwankungen, was? *chatt* ;)

[mj]

Original geschrieben von intel_hasser
Starke Stimmungsschwankungen, was? *chatt* ;)
Immer ;)

[Shaft99]

Original geschrieben von [ab]noname
Wenn der Meister mir ne Möglichkeit nennt über den Ausgabestream eine 64bit Integer Zahl auszugeben bin ich gerne bereit, auf printf zu verzichten. Bis dahin gilt Dieter Nuhr ;D

Das soll jetzt nicht beleigend oder herablassend klingen, aber wie
wärs mit überladen vom operator<< für 64bit Integer ?
Klar ist das mehr Arbeit.
Mit meinem Beitrag waren aber eigentlich die Leute gemeint,
die sich in C++ immer noch mit printf rumquälen. ;D

[i_hasser]

Original geschrieben von [ab]noname
Hatte die letzten Tage eher wenig Zeit. Jetzt muss ich mich langsam mal wieder ransetzen.

Meine 1. Lösung war ja eher trivialer Natur und relativ intuitiv. Bei intel_hassers Lösung sieht man schon, dass es krasser geht ;D (Kleine Frage, bist du Student, in Ausbildung oder ist das rein Hobby?)

Mal schaun, was ich noch zu Stande bekomme...

Im Moment ist es Hobby ;)

[[ab]noname]

Original geschrieben von Shaft99
Das soll jetzt nicht beleigend oder herablassend klingen, aber wie
wärs mit überladen vom operator<< für 64bit Integer ?
Keine Angst klingts nicht ;D
Das Ding ist doch aber, dass es hier eigentlich um den Algorithmus geht und nicht um C++ spez. Dinge.
Ich benutz printf eigentlich überhaupt nicht, es sei denn man kommt damit einfacher zum Ziel...

[i_hasser]

Original geschrieben von [ab]noname
Hatte die letzten Tage eher wenig Zeit. Jetzt muss ich mich langsam mal wieder ransetzen.

Meine 1. Lösung war ja eher trivialer Natur und relativ intuitiv. Bei intel_hassers Lösung sieht man schon, dass es krasser geht ;D (Kleine Frage, bist du Student, in Ausbildung oder ist das rein Hobby?)

Mal schaun, was ich noch zu Stande bekomme...

Im Moment ist es Hobby ;)

Mir raucht hier gerade der Kopf, aber es lohnt sich - der neue Alg dürfte um Längen schneller sein als die bisherigen ;D

[i_hasser]

Puh, endlich!

Bin der schnellste *baeh*

Mein Alg braucht für die 1. Zahl 0.413 Sekunden. Für die 2. Zahl braucht er 5.24 Sekunden und mir fällt gerade auf, dass das ohne Optimierungen war :]

0.38 Sekunden für Zahl1 und 4.86 Sekunden für Zahl2. Die Einsparungen sind zwar nicht so groß wie ich gehofft hatte, aber reichen tuts allemal ;)

[mj]

Naja schaun mer mal, mir ist grade noch eine geniale Idee gekommen wie ich meinen Algorithmus... sagen wir mal, "vectorisieren" kann ;)
Ich bin schon eifrig am coden, hab grade nochmal Nachschub an Zigaretten und Kaffee besorgt und hoffe, dass es bis zwei Uhr früh hinhaut und ich trotz Java und sacklahmer BigInteger-Arithmetik noch ein wenig Performance rausholen kann.
Vielleicht probier ich das auch mal im ML zu schreiben, das könnte dank wunderbarer SMP/SMT Unterstützung nochmal einiges bringen.




@intel_hasser:
Quellcode bitte und Testplattform.

[DarkAvenger]

Ich habe nicht mehr so ganz im kopf, ob ihr das noch braucht, aber evtl habe ich einen Tip:

Eine schnelle Annäherung der Quadratwurzel für natürliche Zahlen.


Gesucht ist x mit x^2=y mit y \in N. Bekanntermaßen ist 2*ln(x)=ln(x^2)=ln(y) => ln(x)=0.5*ln(y)

Nur ist log bzw ln Berechnung nicht unbedingt schneller als sqrt, also noch ein Trick: Es gilt log_b(x)=ln(x)/ln(b). Wähle hier b=2, also lg2(x)=ln(x)/ln(2)

Damit habe ich nun lg2(x)=0.5*lg2(y). Nun was bedeutet lg2 anschaulisch? Es entspricht etwa der # der Binärstellen einer nat. Zahl! Um also die gewünschnte Wurzel zu bekommen ist zu rechnen:

2^(0.5*lg2(y), dh. wir bestimmen die # der Binärstellen von y (per binärer Suche mittles entsprechende shifts und gucken ob die Zahl =0 geworden ist...) und halbieren diese (hier überlegen, ob man ceil oder floor wählen sollte) und shiftet y also um dei ermittelte Hälfte nach rechts. Fertig!

Mathe hilft. ;)

[Sargnagel]

Original geschrieben von DarkAvenger
Mathe hilft. ;)
Wenn man's beherrscht! ;)

Ich bin mal gespannt, wie lange es dauert, bis intel_hasser mit einer hochoptimierten, äußerst bitoperatorhaltigen Funktion für diese Berechnung daherkommt ... *buck*

[i_hasser]

Ich denke die Quadratwurzel ist hier nicht so wirklich der bremsende Faktor - immerhin hat die FPU eigentlich so gut wie nichts zu tun, also eine willkommene Abwechslung für die CPU ;)

Leider sitze ich gerade nicht an meinem Rechner (bäh - Windoofs2000 mit IEh6... immerhin ein Duron), sonst würde ich gleich mal den Source posten. Platform ist ein Athlon XP TbredB 2.2GHz 200MHz FSB, 200MHz DC CL2.5 Ram auf einem Nforce 2 Brettel, das alles unter Linux 2.6.4.

Ich bin sogar am Überlegen die hälfte der Berechnungen auf die FPU auszulagern - den Modulo könnte man ziemlich einfach per FPU implementieren, möglicherweise könnte man das sogar noch einfacher umsetzen. Würde das klappen könnte man es auch mit SSE versuchen... da fällt mir nochwas ein - meinen Alg könnte man über ziemlich große Matritzen realisieren, und sogar SSE tauglich (ohne die Werte der Elemente ständig zu ändern).


Na da muss ich mal sehen, vielleicht kann man da auch einen Mix aus MMX und 3DNow machen *chatt*

[Sargnagel]

Original geschrieben von intel_hasser
da fällt mir nochwas ein - meinen Alg könnte man über ziemlich große Matritzen realisieren, und sogar SSE tauglich (ohne die Werte der Elemente ständig zu ändern).

Hmm ... da fällt mir wiederum was ein. Wenn es sich um große Matrizen handelt und die Operationen auf den einzelnen Elementen bestimmte Bedingungen erfüllen, könnte man die Berechnungen auch auf die Graphikkarte auslagern. *chatt*
(Sourcecode) http://sourceforge.net/projects/brook
(Doku) http://graphics.stanford.edu/projects/brookgpu/
(unbedingt Doku lesen, da dort auf die derzeitigen Beschränkungen von Brook eingegangen wird)

Zugegeben, das wird erst mit PCI-Express-Graphikkarten interessant, da die im Full-Duplex-Modus mit hohen Bandbreiten und sehr geringen Latenzzeiten operieren. Bisher läuft das, glaube ich, teilweise nur unter Windows, da das Visual Studio benötigt wird. :[

Selbst das GROMACS-Team (s. Folding@Home) hat ausprobiert, ob GROMACS auch auf die GPU umgesetzt werden kann. Leider mußte jede Menge Shader-Code per Hand geschrieben werden und die Performance blieb 50% hinter der aktueller CPUs. Wie gesagt, PCI-Express und Geforce 6800 und Nachfolger werden Abhilfe schaffen!

[mj]

Also ich hab's mittlerweile geschafft, die Zeiten nochmals etwas mehr als zu halbieren... für die zweite Zahl ist er jetzt von 142s auf 59s runter.

[mj]

Original geschrieben von DarkAvenger
Mathe hilft. ;)
Klar tut's das... leider völlig irrelevant wenn die Arithmetik für große Zahlen quasi nix kann. Vom ln kann ich bei BigInteger nur träumen, dafür hab ich einen eigenen Weg gefunden die Wurzel anzunähern und bin damit auch erstaunlich nah dran ;)

[BoMbY]

Ehm, für die Wurzel hab ich grad mal ne Berechnung nach dem Newton-Verfahren in Java implementiert (oder besser aus ner alten Pascal Unit geklaut):

public static final BigInteger TWO = new BigInteger("2");

// trunc(sqrt(a)) nach Newton
public static BigInteger sqrt(BigInteger a)
{
// a < 0 return 0
if (a.compareTo(BigInteger.ZERO) == -1) {
return BigInteger.ZERO;
}

// a < 2 return 1
if (a.compareTo(TWO) == -1) {
return BigInteger.ONE;
}

BigInteger x = a, y = BigInteger.ONE;

while (x.compareTo(y) >= 0) {
y = y.add(y);
x = x.divide(TWO);
} do {
y = y.add(x);
x = y.divide(TWO);
y = a.divide(x);
} while (y.compareTo(x) < 0);

return x;
}


Also die ist so schnell, ich konnte auf meinem P4HT 2,4GHz hier auf der Arbeit nur 0ms für die Wurzelberechnung von 11081452861688528621 messen.

Vieleicht hilft's ja dem ein oder anderem... :)

m.f.g.
BoMbY

[mj]

Hmm... aus irgendeinem Grund ist mein neuer Algorithmus nicht linear... bei der ersten Zahl braucht er doppelt so lang wie zuvor, bei der zweiten Zahl wieder nur halb so lang, bei der dritten Zahl nur 1/3 so lang wie zuvor... irgendwie verwirrt mich das *kopfkratz*

[BoMbY]

uhm, okay, die Funktion ist wirklich gut... :D

sqrt(22739118617244071047525087654196635660439990845070408495437407426675870448826942 7103772722301836586752914565712472315683808112009277027836001343973520161) = 15079495554309524397740177588608012694510007187696452637587281796104079812881 <=> 16ms

[BoMbY]

Original geschrieben von D'Espice
Hmm... aus irgendeinem Grund ist mein neuer Algorithmus nicht linear... bei der ersten Zahl braucht er doppelt so lang wie zuvor, bei der zweiten Zahl wieder nur halb so lang, bei der dritten Zahl nur 1/3 so lang wie zuvor... irgendwie verwirrt mich das *kopfkratz*

Hat das vieleicht was mit der Optimierung durch die Runtime und den Prozessor (durch den Cache oder so) zu tun?

m.f.g.
BoMbY

[DarkAvenger]

So ich habe mal auch mal etwas überlegt. Da mir mein beschriebener Alg zu grob war für Abschätzung der Wurzel, habe ich mir mal etwas überlegt. Ich denke, daß es recht schnell sein wird. Es ist unoptimiert, ich denke sogar, daß man es komplett in integer rechnen kann. Habe nicht soviel Ahnung von 64bit oder großen Zahlen, darum hier nur der Beispielalg. der die Wurzel immer ziemlich gut trifft. Die Idee dahinter ist ganz einfach und ist eine Verallgemeinerung meines vorherigen beschrieben Vf. Wenn ich x um eine Position nach links shifte, wird das Quadrat um 2 Positionen nach links geshiftet. Umegekehrt wenn meine Zahl y um einen Bruchteil "geschifted" wird, muß ich entsprechendes mit der Wurzel machen. Dfür berechne ich den Korrekturfaktor. Dieser ließe sich unter Einbeziehung von psq (der unkorrigierten Wurzel, im Prinzip mein oben angegebenes VF) ganz in Integern mit ein paar shifts realisieren, aber hatte keine Lust, etzt daran rumzudoktorn, muß auch was sinnvolles tun... Ebenso ließe sich das Abbruchkriterium der Schleife vorziehen (bis die Bruchteile von psq eigentlich keinen ganzzahligen Anteil mehr geben.)

Hier mein schnell impl Alg:


unsigned long long ha=287530;//input
unsigned char bit;
long psq;
unsigned long long mask;
double factor;
double curfrac;
unsigned char i;

for (bit=0;ha>>bit!=0;++bit);
--bit;


psq=1<<(bit/2);
mask=1<<(bit-1);


if (bit%2) factor=1.5;
else factor=1.0;

curfrac=0.25;
for (i=bit;i>=bit/2;--i)
{
if (ha&mask) factor+=curfrac;
curfrac/=2.0;
mask>>=1;
}

printf("Number: %Lu Bits: %u plain sq: %u factor: %f cor. sq: %f Sqrt: %f\n",ha,bit,psq,factor,factor*psq,sqrt((double)ha));



Jemand kann denn Allg ja mal so für "große" Zahlen umbauen und Zeit stoppen. Ich wette, ich bin sehr schnell, obwohl noch einiges an Optimierung möglich ist.

[i_hasser]

Kurze Zwischenfrage - gibts eine Funktion die den Nachkommateil einer Fließkommazahl zurückgibt?

[BoMbY]

Original geschrieben von intel_hasser
Kurze Zwischenfrage - gibts eine Funktion die den Nachkommateil einer Fließkommazahl zurückgibt?

Ich kenne keine, aber wie wär's mit:

zahl -= trunc(zahl)

?

m.f.g.
BoMbY

[i_hasser]

Mist!

Gibts noch größere Fließkommazahlen als long double? Ich denke mal nicht... der ist leider schon zu klein für die großen Zahlen.

[i_hasser]

Ist doch nicht zu klein, aber dummerweise ist der Athlon mit Integer schneller als mit long double, und gcc erzeugt den Code so, dass der Athlon den schlecht parallelisieren kann :(


Ich glaub der GCC will mich vereiern - bekomm ein segmentation fault, weil ich 128bit long doubles wollte :] - also beim gcc, mein prog erstellt er erst garnet.

[Sargnagel]

Original geschrieben von intel_hasser
Gibts noch größere Fließkommazahlen als long double?
Laut C99-Standard nicht. :(

[i_hasser]

Naja, es funktioniert jetzt, ist aber wie gesagt ein kleines bisschen langsamer als mit Integern :(

Aber ich arbeite noch dran - gerade schau ich mir die 3dnow Spezifikationen an, ob man da überhaupt mit doubles > 64bit hantieren kann (ist ja Vorraussetzung).

Und ich will mal das Verhältniss bei den Operationen ändern, von 1:1 für int/float zu 2:1 für int/float. Dann will ich noch ein paar Codeteile verlagern, so dass der Athlon besser parallelisieren kann.

[Sargnagel]

Soweit ich weiß, kann 3dnow nicht mit double precision floating point numbers umgehen. Aber Du kannst die 64bits da sicherlich hineinschreiben und dann mit wilden 3dnow/mmx Befehlen manipulieren. ;D
Oder hat sich mit der 3dnow extension etwas verändert?

Eine 224bit floating point library wurde hier (http://www.planet3dnow.de/vbulletin/showthread.php3?s=&postid=1061858&highlight=mandelbrot+floating+point+library#post1061858) schon einmal erwähnt.

[i_hasser]

Ich glaub ich muss mir einen Mac kaufen, AltiVec kennt 64bit Integer ;D

Alles was auf x86 läuft kennt die leider nicht, und die "normalen" doubles mit 48 oder 64bit reichen auch nicht aus, müssen schon welche mit mindestens 80bit sein.

Allerdings kennt MMX wenigstens einzelne 64bit integer, damit entfällt hoffentlich der Registerhickhack und die 64bit Operationen dürften damit nochmal schneller werden - ich implementiere es gleich mal ;D

[i_hasser]

Nein, MMX bringt auch nix. Habs erstmal in einem neuen Prog getestet (einfach 500 Millionen mal +=2) - wenn man MMX anstellt scheint er 64bit integer immer mit mmx zu berechnen.


Ich werd die parallelisierungen jetzt wieder rausnehmen und den Code posten.


EDIT: Nein, es geht doch was ;D - momentan ist der Fließkommavergleich arschlahm, aber die Berechnung an sich ist viel schneller als mit Integer...


EDIT2: So, also das Integer-Modulo braucht ca. 2.5mal so lange wie eine Float-Division.

[Sargnagel]

Original geschrieben von intel_hasser
momentan ist der Fließkommavergleich arschlahm
Poste doch bitte mal den entsprechenden Quellcodeabschnitt. Eventuell kann man daran noch was drehen. Ich suche derweil mal nach der WWW-Quelle, wo ich die Lösung vermute ...


Ach, mist. Das geht nur mit 32bit floats. :-X

float x[MAX_ENTRIES];
long temp;

...

temp = *((long *)&(x[i]));
if( temp + temp ) {

....

}

Was es bedeutet, kann hier (http://www.azillionmonkeys.com/qed/optimize.html#tricks) nachgelesen werden.

[i_hasser]

if(fdummy==trunc(fdummy))

So sah es bisher aus, fdummy ist zahl1/zahl2, das if ist dann für zahl1%zahl2=0 wahr. Ich hab auch schon modfl versucht, aber das ist noch langsamer.

[DarkAvenger]

@sargnagel

Wenn das dein code ist: gcc3.4 wird ihn dir um die Ohren hauen (zumindest meckern). Du solltest lieber unions benutzen.

[Sargnagel]

@DarkAvenger:
Kein Sorge, der Quellcode ist nicht auf meinem Mist gewachsen. ;D Der Link in meinem obigen Posting führt zur Quelle.

@intel_hasser:
Ist fdummy denn immer noch vom Typ double oder ist es int? Weil ein "==" ist bei Fließkommazahlen immer so eine Sache ... *kopfkratz

[i_hasser]

Sind alles long doubles - das mit dem == weis ich, hat so aber auf Anhieb funktioniert und ich habs mir gespart da noch was passenderes zu implementieren ;)

[i_hasser]

ich könnt mich erschießen, wozu der ganze Aufwand wenn fmodl ein Modulo für floats über die FPU implementiert? *buck*

EDIT: 4.6 Sekunden, und jetzt fängt das Optimieren erst an ;)

EDIT2: Manchmal frag ich mich echt wieso ich überhaupt mit Programmieren angefangen hab... der Athlon schaffts wieso auch immer nicht die Integer und Float Operationen parallel auszuführen, also ist es am schnellsten wenn ich alles per Float berechnen lasse :P

[PuckPoltergeist]

Original geschrieben von intel_hasser
EDIT2: Manchmal frag ich mich echt wieso ich überhaupt mit Programmieren angefangen hab... der Athlon schaffts wieso auch immer nicht die Integer und Float Operationen parallel auszuführen, also ist es am schnellsten wenn ich alles per Float berechnen lasse :P

Das liegt aber mit ziehmlicher Sicherheit am Compiler, nicht am Athlon.

[i_hasser]

So, nun gibts Code ;D

Inzwischen beträgt die Berechnunszeit weniger als 3 Sekunden für die 2. Zahl (2.98 *buck*). Allerdings steigt dann die Initialisierungszeit auf 0.7 Sekunden (die ist unabhängig von der Größe der Zahl) - wenn ich die Optimierung ein bisschen zurückfahre beträgt die Initialisierungszeit noch 0.1 Sekunden, dafür die Rechenzeit wieder 3.05 Sekunden :P

Der Code ist ziemlich versaut, weil ich da "mal eben so" alles mögliche umgeändert habe, und inzwischen sicherlich keine Zeile unangetastet geblieben ist. Den Primzahlkrempel lass ich mal raus, weil der wirklich nur bei der Initialisierung zum Einsatz kommt und die Primzahlen alle sehr klein sind (eigentlich könnte ich den ganzen Primzahlkrempel komplett rausnehmen und ein paar einfache Modulo Funktionen implementieren).

Die Erkennungsroutine ob es sich beim gefundenen Teiler um eine Primzahl handelt hab ich auch erstmal rausgenommen, die Erkennungszeit dafür ist aber sehr klein (müssen ja nur 2 Primzahlen getestet werden).

Erklärung zum Code kommt im nächsten Beitrag ;) - einige auskommentierte Stellen hab ich jetzt gelöscht, weil die nichts zur Sache tun. Ach ja - Compilerparameter: "-O4 -march=athlon-xp -mno-ieee-fp"

#ifdef HAVE_CONFIG_H
#include <config.h>
#endif

#include <stdio.h>
#include <stdlib.h>
#include <time.h>

#include <primeapi.h>
#include <largenum.h>

// ---- Vars ----
int64u *alg_diffsel; // Table for Jumps between numbers
int32u alg_diffsel_size; // Entries in Table
int32u alg_diffsel_area; // Area the Table spans


void
create_divsel_tab(int32u arg)
{
printf("\n creating diffsel table...");

int32u diffsel_primes=arg; // Exclude the first ... primes

// calc size
int32u i,j,k,l;
int32u area_size=1;


// get area_size
for(i=0;i<diffsel_primes;i++) area_size=area_size*prime_getprimeno(i+1);
printf(" area_size=%ld",area_size);
fflush(stdout);

// count the size of divsel tab
int32u diffsel_size=0;

for(i=0;i<area_size;i++)
{
for(j=0;j<diffsel_primes;j++) if(i%prime_getprimeno(j+1)==0) break;

// if its a prime
if(j==diffsel_primes) {
diffsel_size++;
//printf("\np: %ld",i);
}
}

printf(" diffsel_size=%ld",diffsel_size);

k=0; // last prime found...
l=0; // prime count
int64u *diffsel=malloc(sizeof(int64u)*diffsel_size);

for(i=0;i<area_size;i++)
{
for(j=0;j<diffsel_primes;j++) if(i%prime_getprimeno(j+1)==0) break;

// if its a prime
if(j==diffsel_primes) {
diffsel[l]=(i-k);
//printf("\ndiff (%ld->%ld): %ld",k,i,(i-k));
l++;
k=i;
}
}

printf(" r=%ld%%",(diffsel_size*100)/(area_size/2));

alg_diffsel=diffsel;
alg_diffsel_size=diffsel_size;
alg_diffsel_area=area_size;

}


int
main(int argc, char *argv[])
{
// ---- Vars ----
int32u i,j;

clock_t time1, time2;
clock_t timei, timei1, timei2;
clock_t timee, timee1;

// ---- Initialization ----

// Init Prime API - do we need that really?
time1=clock();
printf("\ninitialisation...");
printf("\n creating prime api...");
prime_init(1000000);
time2=clock();
timei1=time2-time1;

// -- create diffsel --
time1=clock();
create_divsel_tab(7);
time2=clock();
timei2=time2-time1;
timei=timei1+timei2;

// ---- Calc ----
time1=clock();
int64u alg_inumber=312133288776418481;

long double alg_fnumber=(long double) alg_inumber;
int64u alg_sqrtnumber=sqrtl ( (long double) alg_inumber );

printf("\n\ncalculation...");
int64u alg_fullblocks=alg_sqrtnumber / (int64u) alg_diffsel_area;
int64u alg_partblock=alg_sqrtnumber % (int64u) alg_diffsel_area;

printf("\n fullblocks=%ld partblock=%ld",alg_fullblocks,alg_partblock);
fflush(stdout);

// Helpers
int64u i64, j64;

// Number Pointer
int64u alg_divptr=0;

// Area Number Pointer
int32u areaptr_size=alg_diffsel_size;
int32u areaptr_area=alg_diffsel_area;

int64u areaptr_iarea=(int64u) areaptr_area;
long double areaptr_farea=(long double) areaptr_area;

int64u* areaptr_inum=malloc(sizeof(int64u)*areaptr_size);
long double fdummy,*areaptr_fnum=malloc(sizeof(long double)*areaptr_size);


long double areafptr_size=(long double) alg_diffsel_size;
long double areafptr_area=(long double) alg_diffsel_area;

// -- Int64
areaptr_inum[0]=(int64u) alg_diffsel[0];
for(i=1;i<areaptr_size;i++) areaptr_inum[i]=areaptr_inum[i-1]+(int32u)alg_diffsel[i];

// -- long double
areaptr_fnum[0]=(int64u) alg_diffsel[0];
for(i=1;i<areaptr_size;i++) areaptr_fnum[i]=areaptr_fnum[i-1]+(long double)alg_diffsel[i];

// -- Modulo eval vars
int64u modiev1;
long double modfev1;

// Calc complete Blocks
for(i=0;i<alg_fullblocks;i++)
{
// Compute - half parallelized
for(j=0;j<areaptr_size;j++)
{
// do float calculation
modfev1=fmodl(alg_fnumber,areaptr_fnum[j]);
if(modfev1==0.0) printf("\n fdiv found: %.0Lf * %.0Lf",areaptr_fnum[j+0],(alg_fnumber/areaptr_fnum[j]));
areaptr_fnum[j]+=areaptr_farea;

} // for
}

// Calc last uncomplete Block
j=0;
alg_divptr=alg_fullblocks*alg_diffsel_area;
for(j=0;;j++) {

alg_divptr+=(int64u) alg_diffsel[j];
if(alg_divptr>alg_partblock+alg_fullblocks*alg_diffsel_area) break;

// if divideable
if(alg_inumber%alg_divptr==0)
{
printf("\n idiv found: %Ld * %Ld",alg_divptr,alg_inumber/alg_divptr);

/*if(prime_isprime(alg_divptr)==PRIMECHK_IS_PRIME)
{
if(prime_isprime(alg_inumber/alg_divptr)==PRIMECHK_IS_PRIME)
{
}
}*/
}
}

time2=clock();
timee1=time2-time1;
timee=timee1;

// ---- print time statistics ----
printf("\n\n---- time statistics ----");
printf("\ninitialisation: %.3fs",(double)timei/CLOCKS_PER_SEC);
printf("\n prime init: %.3fs",(double)timei1/CLOCKS_PER_SEC);
printf("\n diffsel init: %.3fs",(double)timei2/CLOCKS_PER_SEC);
printf("\n\ncalculation: %.3fs",(double)timee/CLOCKS_PER_SEC);

printf("\n");
return 0;
}

[i_hasser]

Original geschrieben von PuckPoltergeist
Das liegt aber mit ziehmlicher Sicherheit am Compiler, nicht am Athlon.

Klar, hatte nur nicht so auf meine Formulierung geachtet ;)

Vielleicht bastle ich das noch in Assembler ;D


Also zum Alg. Im Grunde genommen ist das auch nur ein "einfacher" Modulo-Alg, der für jede wichtige (<- das ist wichtig ;)) Zahl testet, ob sich die gegebene Zahl dadurch teilen lässt. Das Modulo ist als Float implementiert, was auf dem Athlon deutlich mehr Speed bringt (mein K5 dagegen würde die Integer-Variante sicherlich bevorzugen - in Punkto Integer-IPC macht der meinen Athlon platt ;D).

Der entscheidende Punkt sind die wichtigen Zahlen - tja, welche Zahlen sind nun wichtig. Wer sich meinen Primzahlalgorithmus mal angesehen hat ist sicherlich das hier aufgefallen:


for(i=0;i<4;i++) {
prime_bit_sel[1+i*30]=0;
prime_bit_sel[7+i*30]=1;
prime_bit_sel[11+i*30]=2;
prime_bit_sel[13+i*30]=3;
prime_bit_sel[17+i*30]=4;
prime_bit_sel[19+i*30]=5;
prime_bit_sel[23+i*30]=6;
prime_bit_sel[29+i*30]=7;
}

(die Variablie hieß früher BitSelector)

Die Tabelle schließt alle Vielfachen von 2, 3 oder 5 aus. Das i ist dabei nicht wichtig (das liegt nun wieder am Primzahlalg), also dafür kann man 0 einsetzen.

Die Sache basiert darauf, dass sich die Muster für Vielfaches/kein Vielfaches bei einer bestimmten Menge an Ausgangszahlen (oben zb. 2, 3 und 5) stetig wiederholen - genau nach 2*3*5=30 Zahlen.

Also ich mach es mal am Beispiel von 2 und 3:

0: vielfaches
1: kein vielfaches
2: vielfaches (1*2)
3: vielfaches (1*3)
4: vielfaches (2*2)
5: kein vielfaches
-------------------
6: -> 0
7: -> 1
8: -> 2
...

So - mein erster Alg (hat sich da eigentlich schon jemand durchgefunden? Ich glaub das würd ich selbst nimmer schaffen :]) basierte darauf, nur für Primzahlen zu testen ob DiegroßegegebeneZahl%Primzahl==0. Dummerweise ist ein Modulo aller Zahlen schneller erledigt als die Bestimmung der nächsten Primzahl (desswegen ist der Alg langsamer).

Im Alg von [ab]noname waren Vielfache von 2 schon alle ausgeschlossen - bei dem Teiler hieß es sinngemäß for(i=1;i<...;i+=2).

Das erspart uns schonmal 50% aller Zahlen. Aber schön wäre es ja nun, wenn wir die Vielfachen von noch mehr Primzahlen ausschließen könnten. Tja, und genau das macht ja die oben gegebene Tabelle zb. schon für 2, 3 und 5 - das würde die zu testenden Zahlen auf 8/30~26.6% gegenüber allen Zahlen, oder eben auf 53.333% gegenüber jeder 2. Zahl reduzieren.

Im Alg ist das dynamisch realisiert, mit "create_divsel_tab(7);" bestimmt man zb. dass die 1. 7 Primzahlen ausgeschlossen werden (mir ist kein besserer Name eingefallen, und um später nicht mit divptr, der inzwischen auch wieder halb verbannt ist und weil mir der Name einfach net gefallen hat hab ich das später zu diff... geändert, nur der Funktionsname fängt noch mit div... an - ich sag doch, schweinischer Code *chatt*).

Das wären also 2, 3, 5, 7, 11, 13 und 17. Das reduziert die zu testenden Zahlen gegenüber jeder 2. Zahl auf 36% und die anzulegenden Tabellen werden noch relativ klein (bei create_divsel_tab(8); werden die Tabellen immerhin 19mal so groß und die zu testende Zahlenmenge wird nur auf >18/19tel reduziert).
Damit lassen sich also effektiv Modulo-Operationen sparen (das ist das Ziel des ganzen ;)).

Nun wieder zu der Tabelle oben -

0: vielfaches
1: kein vielfaches
2: vielfaches (1*2)
3: vielfaches (1*3)
4: vielfaches (2*2)
5: kein vielfaches

(btw: die 6 entspricht der 0 (und daher endet die Liste bei 5), weil die letzte Zahl in der Liste immer ein vielfaches (genauer aller ausgeschlossenen Primzahlen) ist)

Wichtig sind ja nur die Zahlen, die keine Vielfache sind (wenn wir die gegebene Zahl%2 testen erübrigt sich ein %4 ja). Also könnten wir die Tabelle auch so umschreiben:

+1
+4

Sinn: Ausgehend vom 1. Block befindet sich die 1. Zahl die kein Vielfaches ist an der Position+1, davon ausgehend die nächste relevante Zahl an dieser Position+4.
Das reduziert die Tabelle schonmal auf 2 Einträge.

Als ich das parallelisiert hab, hab ich das wieder abgeändert. Da sieht die Tabelle nun so aus:

+1
+5

Das sind also nur die absoluten Positionen, ausgehend vom Anfang des Blockes. So sieht die areaptr_inum und areaptr_fnum aus (die beiden unterscheiden sich übrigens nur dadurch, dass inum int64 Werte enthält und fnum 96bit float Werte - inum ist inzwischen überflüssig, da ich die integer operationen ja komplett rausgenommen hab).

areaptr_size ist für den Fall hier 2 (2 Einräge in der letztendlichen Tabelle) und areaptr_area ist 6 (alle 6 Zahlen wiederholt sich das Muster).

Jetzt reicht es aus die Geschichte so hier zu testen:

offset=0*areaptr_area
if zutestendezahl % (offset+areaptr_fnum[0]) == 0 -> Teiler gefunden
if zutestendezahl % (offset+areaptr_fnum[1]) == 0 -> Teiler gefunden

offset=1*areaptr_area
if zutestendezahl % (offset+areaptr_fnum[0]) == 0 -> Teiler gefunden
if zutestendezahl % (offset+areaptr_fnum[1]) == 0 -> Teiler gefunden

offset=2*areaptr_area
if zutestendezahl % (offset+areaptr_fnum[0]) == 0 -> Teiler gefunden
if zutestendezahl % (offset+areaptr_fnum[1]) == 0 -> Teiler gefunden

...

Im Alg hab ich das noch weiter optimiert. Da gibts nicht 1 Offset sondern (für das Beispiel hier) 2 Offsets, auf jedes wird nach einem if-durchlauf areaptr_area draufaddiert.

Tja, und das wars eigentlich schon. Die prime_... Funktionen sind eigentlich geblieben und spielen eben nur beim Erstellen der Tabellen eine Rolle. An den Funktionen hat sich nicht viel geändert (nur unten drunter), neu ist nur prime_getprimeno(int) - gibt die angegebene Primzahl zurück, also für 1->2, 2->3, 3->5, 4->7 usw - hat man schnell nachgeschrieben.

[BoMbY]

So, meine noch recht simple Java-Lösung:

public class Factorize
{
/**
* The BigInteger constant two.
*/
public static final BigInteger TWO = new BigInteger("2");

/**
* Returns a BigInteger whose value is the truncated SuqareRoot of a.<br>
*
* @param a - value from which the SquareRoot should be computed.
* @return trunc(sqrt(a))</tt>
*/
public static BigInteger sqrt(BigInteger a)
{
// a < 0 return 0
if (a.compareTo(BigInteger.ZERO) == -1) {
return BigInteger.ZERO;
}

// a < 2 return 1
if (a.compareTo(TWO) == -1) {
return BigInteger.ONE;
}

BigInteger x = a, y = BigInteger.ONE;

while (x.compareTo(y) >= 0) {
y = y.add(y);
x = x.divide(TWO);
} do {
y = y.add(x);
x = y.divide(TWO);
y = a.divide(x);
} while (y.compareTo(x) < 0);

return x;
}

public static BigInteger[] twoFactors(BigInteger n)
{
BigInteger[] result = new BigInteger[2];
BigInteger i = sqrt(n);

if (i.mod(TWO).compareTo(BigInteger.ZERO) == 0)
i = i.subtract(BigInteger.ONE);

while (i.compareTo(TWO) == 1)
{
if (n.remainder(i).compareTo(BigInteger.ZERO) == 0)
{
result[0] = i;
result[1] = n.divide(i);
return result;
}

i = i.subtract(TWO);
}

result[0] = result[1] = BigInteger.ZERO;

return result;
}
}

Mit diesen Zeiten auf meinem AthlonXP 2600+:

twoFactors(1949565409282837) = 19827473 * 98326469 <=> 9578ms
twoFactors(312133288776418481) = 398527487 * 783216463 <=> 62438ms
twoFactors(11081452861688528621) = 2783216473 * 3981527477 <=> 214719ms

Naja, das reicht für heute, aber da werd ich wohl nochmals ans Reißbrett müssen... :)

m.f.g. und gute Nacht,
BoMbY

[BoMbY]

@intel_hasser: Warum nennst Du Dein Kind eigentlich nicht beim Namen? :) Ich überlege die ganze Zeit woher mir Deine Vorgehensweise bekannt vorkommt, und dann finde ich das hier:

http://home.t-online.de/home/arndt.bruenner/mathe/scripts/eratosthenes.htm

Also basiert Dein Alorithmus praktisch auf dem Sieb des Eratosthenes.

m.f.g.
BoMbY

[BoMbY]

Ich hab grad mal ne nette Seite mit ein paar bekannten Algorithmen zur Primfaktorzerlegung gefunden. Da gibt's echt ne Menge Theorie nachzuholen bei mir, stell ich gerade fest: http://mathworld.wolfram.com/PrimeFactorizationAlgorithms.html

m.f.g.
BoMbY

[i_hasser]

Original geschrieben von BoMbY
@intel_hasser: Warum nennst Du Dein Kind eigentlich nicht beim Namen? :) Ich überlege die ganze Zeit woher mir Deine Vorgehensweise bekannt vorkommt, und dann finde ich das hier:

http://home.t-online.de/home/arndt.bruenner/mathe/scripts/eratosthenes.htm

Also basiert Dein Alorithmus praktisch auf dem Sieb des Eratosthenes.

m.f.g.
BoMbY

Ja, das hab ich ja schon im Primzahlwettbewerb angewandt. Von der Sache her ist es das Sieb, allerdings immer nur ein Ausschnitt davon und desswegen doch wieder etwas anders ;).

[BoMbY]

Original geschrieben von intel_hasser
Ja, das hab ich ja schon im Primzahlwettbewerb angewandt. Von der Sache her ist es das Sieb, allerdings immer nur ein Ausschnitt davon und desswegen doch wieder etwas anders ;).

Wobei ich zugeben muss: Ich verstehe zwar das Prinzip des Siebs, aber ich hab keine Ahnung was Du mit Deinem C-Code da machst, auch nicht mit Deiner Beschreibung ein paar Posts darunter ... :)

m.f.g,
BoMbY

[i_hasser]

Nagut, dann gibts jetzt ein paar Kommentare ;)

Also erstmal zur create_divsel_tab Funktion. Die erstellt die Tabelle, die nach dem Beispiel so aussieht:

+1
+4

Die hab ich ziemlich schnell zusammengehackt und zwischendurch auch komplett umgekrämpelt, also nicht wundern wenn man da viel optimieren kann. Die Initialisierungszeit ist aber im Vergleich zur Rechenzeit so gering, dass es sich nicht lohnt :)


// ---- Vars ----
int64u *alg_diffsel; // Table for Jumps between numbers
int32u alg_diffsel_size; // Entries in Table
int32u alg_diffsel_area; // Area the Table spans

void
create_divsel_tab(int32u arg)
{
// Hier hab ich nur einen Namen geändert, desswegen die Zuweisung
// Die Variable sagt aus, dass bis zur arg.ten Primzahl "vorgefiltert" werden soll,
// Also die Tabelle erfasst wenn arg=4 die Primzahlen 2, 3, 5 und 7
int32u diffsel_primes=arg;

// Das sind ein paar Hilfsvariablen
int32u i,j,k,l;
int32u area_size=1;


// Hier wird die Größe des Blocks berechnet, der sich periodisch immer
// wiederholt. Für 2, 3, 5 und 7 wäre das 2*3*5*7=210
for(i=0;i<diffsel_primes;i++) area_size=area_size*prime_getprimeno(i+1);
printf(" area_size=%ld",area_size);
fflush(stdout);

// Tja, der Teil ist Opfer einer umstellung geworden und daher ziemlich
// undurchsichtig. Hier wird berechnet, wie viele Zahlen es in einem Block
// gibt, die nicht für das Beispiel durch 2, 3, 5 und 7 teilbar sind.
// diffsel_size enthält dann diese Zahl.
int32u diffsel_size=0;
for(i=0;i<area_size;i++)
{
for(j=0;j<diffsel_primes;j++) if(i%prime_getprimeno(j+1)==0) break;

// if its a prime
if(j==diffsel_primes) {
diffsel_size++;
//printf("\np: %ld",i);
}
}

printf(" diffsel_size=%ld",diffsel_size);

// Nun gehts ans erstellen der Tabelle
k=0;
l=0;
int64u *diffsel=malloc(sizeof(int64u)*diffsel_size);

// Und hier wirds wieder etwas komplifizierter. k weist in der Schleife immer auf die
// Position der vorheringen nicht durch im Bsp. 2, 3, 5 und 7 teilbaren Zahl. Angefangen bei der 0,
// daher da oben auch das k=0
for(i=0;i<area_size;i++)
{
// Hier testen wir die Teilbarkeit auf im Bsp. 2, 3, 5 und 7 - diffsel_primes sagt ja, dass
// wir alles bis zur diffsel_primes.ten Primzahl ausschließen sollen
// prime_getprimeno(int64u x) gibt die x.te Primzahl zurück.
for(j=0;j<diffsel_primes;j++) if(i%prime_getprimeno(j+1)==0) break;

// Hier ist ein bisschen Optimierung am Werk. Wurde die schleife mitten drinnen abgebrochen
// ist die Durchlaufbedingung j<diffsel_primes wahr. Wenn die schleife bis zum Ende durchgelaufen
// ist, ist diese Bedingung falsch (daher wurde dann die for-Schleife auch abgebrochen).
// Genau das testen wir hier - wurde die Schleife oben mittendrinnen abgebrochen wurde ein Teiler
// gefunden und die Zahl interessiert nicht weiter, wurde kein Teiler gefunden ist es eine nicht
// durch die gegebenen Teiler teilbare Zahl -> die wollen wir.
if(j==diffsel_primes) {
// Berechnung der Differenz zur Position der letzten nichtteilbaren Zahl
// Die Tragen wir dann in die Liste ein
diffsel[l]=(i-k);

// l sagt, wie viele Zahlen wir bisher gefunden haben
// und k zeigt nun auf die Position der gerade gefundenen Zahl, da wir ja für die nächste
// Zahl die wir finden wieder die Differenz zur Position dieser Zahl brauchen.
l++;
k=i;
}
}

printf(" r=%ld%%",(diffsel_size*100)/(area_size/2));

// Hier schreiben wir die gefundenen Werte in die globalen Variablen und dann gehts raus aus der Funktion :)
alg_diffsel=diffsel;
alg_diffsel_size=diffsel_size;
alg_diffsel_area=area_size;
}

[i_hasser]

Mir schwirrt da gerade wieder eine kleine Idee durch den Kopf...

Mal die Ansätze:

für a=x/y gilt a->0 wenn y->unendlich
und für b=x/y-x/(y+1) gilt b->0 wenn y/(y+1) -> 0, das wiederum für y->unendlich gilt.

Oder auf gut Deutsch: für sehr große y gilt x/y-x/(y+1) ~ x/(y+1)-x/(y+2)

Also wenn c=x/y-x/(y+1) und x1/y ~ trunc(x1/y), dann dürfte auch gelten (für sehr große Zahlen) x1/(y+1/c) ~ trunc( x1/(y+1/c) )


Wenn das hinhaut dürfte die Ersparnis vor allem bei großen Zahlen extrem sein (Faktor 10 bis 1000).

[BoMbY]

Danke, aber vieleicht sollte ich mich heute auf das Wesentliche konzentrieren:
http://www.smiley-channel.de/grafiken/smiley/rauchen/rauchen028.gif+http://www.smiley-channel.de/grafiken/smiley/trinken/trinken015.gif=http://www.smiley-channel.de/grafiken/smiley/zensiert/zensiert003.gif

;D

m.f.g.
BoMbY

[i_hasser]

Also je mehr ich über den Alg nachdenke, desto besser gefällt er mir. Das Potenzial hab ich mir inzwischen direkt an Zahlen angesehen - ich lag vorhin um Dimensionen daneben, damit könnte man den Aufwand auf 1/1e6 reduzieren, oder sogar noch besser.

Ich hab inzwischen mal mit einem Kumpel (der Physik studiert) drüber gesprochen, der meint auch da könnte sehr viel Potenzial drinnen stecken.

[Sargnagel]

Das ist ja alles schön und gut, aber kannst Du mir als Nicht-Mathematiker noch erklären, was Du mit ~ meinst? Bestimmt "ist ungefähr gleich" und nicht den Komplement-Operator, oder? ;) Und was meinst Du mit "x1"?

[i_hasser]

Ja klar, das ~ steht für ungefähr gleich ;)

x1 sei einfach als eine Zahl definiert für die die Bedingung zutrifft - mehr nicht. Man könnte natürlich auch einfach y entsprechend definieren.

[i_hasser]

Na endlich! Hab mich die letzten paar Stunden im Kreis gedreht, aber jetzt bin ich auf der Zielgeraden ;D

[i_hasser]

So... in bestimmten Intervallen funktioniert der Alg jetzt, aber das ist alles noch sehr wackelig. Aber der Aufwand wird immerhin schon auf 1/20 reduziert (im Vergleich dazu alle Zahlen durchzutesten) - und da ist noch wesentlich mehr drinnen.

PS Ich schlag mich hier gerade mit Abeitungen rum :]

[i_hasser]

Also ich hab jetzt keine Zeit (und Lust) mehr mich mit dem Alg zu beschäftigen, ich trag einfach mal zusammen was ich so gemacht hab - vielleicht hat ja von euch jemand Lust das weiterzumachen oder den kompletten Alg zu verbessern.

Also erstmal ein paar Namenssachen:
a ist die Zahl deren Primzerlegung wir suchen
i ist der Zähler, mit dem wir Modulo-Operationen durchtesten. Also der einfachste Alg wäre
for(i=2; (i<sqrt(a)) && (a&i!=0) ;i++);
in schönster C-Schreibweise.

So, das Modulo können wir ja durch eine Division ersetzen - gilt nämlich (a/i)=rtintl(a/i) (rintl ist ähnlich trunc, funzt aber besser) gäbe das Modulo 0.
Glücklicherweise unterstützt die FPU auch ein Modulo, für lange Doubles fmodl - also fmodl(a,i)=a%i und fmodl(a,1)=a/i-rintl(a/i) (ich hab im Alg weitestgehend fmodl genommen, weil es am saubersten und schnellsten ist - rintl hat Rundungsregeln, trunc auch).

Bei großen Zahlen (mit denen wir ja zur Genüge zu tun haben ;)) fällt eins ziemlich stark auf:
a/i=a/(i+1)+c, wobei c ganz stark gegen 0 geht (also man kanns fast vernachlässigen).

a/i-a/(i+1) entspricht dabei auch fast der Ableitung von x(i)=a/i an der Stelle i - x'(i)=-a/(b*b)

x'(i)-x'(i+1) ist bei größeren Zahlen nun fast 0, also ist dementsprechend x''(i) nahe 0 (ist ja nur ne andere Schreibweise für die Differenz und etwas genauer).

Damit man das nachvollziehen kann ein Beispiel:


i= 19826480 a/i= 98331393.6353219028314925
i= 19826481 a/i= 98331388.6757229888680740
i= 19826482 a/i= 98331383.7161245751995011
i= 19826483 a/i= 98331378.7565266618403257
i= 19826484 a/i= 98331373.7969292487759958
i= 19826485 a/i= 98331368.8373323360137874
i= 19826486 a/i= 98331363.8777359235537006
i= 19826487 a/i= 98331358.9181400113884592
i= 19826488 a/i= 98331353.9585445995253394
i= 19826489 a/i= 98331348.9989496879643411
i= 19826490 a/i= 98331344.0393552767054643
i= 19826491 a/i= 98331339.0797613657414331
i= 19826492 a/i= 98331334.1201679550795234
i= 19826493 a/i= 98331329.1605750447124592
i= 19826494 a/i= 98331324.2009826346475165


Wie man sieht ändert sich a/i immer um einen halbwegs konstanten Wert. Der Wert sieht so aus:


i= 19826480 fdelta= -4.9595991641139478
i= 19826481 fdelta= -4.9595986638134618
i= 19826482 fdelta= -4.9595981635130516
i= 19826483 fdelta= -4.9595976632127171
i= 19826484 fdelta= -4.9595971629124583
i= 19826485 fdelta= -4.9595966626122752
i= 19826486 fdelta= -4.9595961623121679
i= 19826487 fdelta= -4.9595956620121362
i= 19826488 fdelta= -4.9595951617121802
i= 19826489 fdelta= -4.9595946614122999
i= 19826490 fdelta= -4.9595941611124953
i= 19826491 fdelta= -4.9595936608127664
i= 19826492 fdelta= -4.9595931605131132
i= 19826493 fdelta= -4.9595926602135357
i= 19826494 fdelta= -4.9595921599140340
i= 19826495 fdelta= -4.9595916596146079
i= 19826496 fdelta= -4.9595911593152575
i= 19826497 fdelta= -4.9595906590159828
i= 19826498 fdelta= -4.9595901587167838
i= 19826499 fdelta= -4.9595896584176606


Und der Faktor ist wie man sieht so ziemlich konstant. Nicht wundern wenn das nicht die genauen Differenzen für die Werte oben sind, die hier ergeben sich aus der Ableitung.

So, also die Werte ändern sich fast garnicht. Die Änderung von fdelta (also fdelta'(i)) sieht so aus:


i= 19826480 fdeltadelta= 0.0000005003005238
i= 19826481 fdeltadelta= 0.0000005003004481
i= 19826482 fdeltadelta= 0.0000005003003724
i= 19826483 fdeltadelta= 0.0000005003002966
i= 19826484 fdeltadelta= 0.0000005003002209
i= 19826485 fdeltadelta= 0.0000005003001452
i= 19826486 fdeltadelta= 0.0000005003000695
i= 19826487 fdeltadelta= 0.0000005002999938
i= 19826488 fdeltadelta= 0.0000005002999181
i= 19826489 fdeltadelta= 0.0000005002998424
i= 19826490 fdeltadelta= 0.0000005002997667
i= 19826491 fdeltadelta= 0.0000005002996910
i= 19826492 fdeltadelta= 0.0000005002996153
i= 19826493 fdeltadelta= 0.0000005002995396
i= 19826494 fdeltadelta= 0.0000005002994639
i= 19826495 fdeltadelta= 0.0000005002993882
i= 19826496 fdeltadelta= 0.0000005002993125
i= 19826497 fdeltadelta= 0.0000005002992368
i= 19826498 fdeltadelta= 0.0000005002991611
i= 19826499 fdeltadelta= 0.0000005002990854


So - wie man sieht ist die Änderung wirklich sehr sehr klein. fdeltadelta ändert sich auch minnimalst, aber das können wir dann wirklich vernachlässigen, wir wollen ja nur eine Aussage in einem bestimmten kleinen Bereich treffen, wo die Änderungen keine Rolle spielen.

Tja - welche Aussage wollen wir eigentlich treffen? Interessant sind die Nachkommastellen von a/i -


i= 19826480 fract= 0.6353219028314925
i= 19826481 fract= 0.6757229888680740
i= 19826482 fract= 0.7161245751995011
i= 19826483 fract= 0.7565266618403257
i= 19826484 fract= 0.7969292487759958
i= 19826485 fract= 0.8373323360137874
i= 19826486 fract= 0.8777359235537006
i= 19826487 fract= 0.9181400113884592
i= 19826488 fract= 0.9585445995253394
i= 19826489 fract= 0.9989496879643411
i= 19826490 fract= 0.0393552767054643
i= 19826491 fract= 0.0797613657414331
i= 19826492 fract= 0.1201679550795234
i= 19826493 fract= 0.1605750447124592
i= 19826494 fract= 0.2009826346475165
i= 19826495 fract= 0.2413907248846954
i= 19826496 fract= 0.2817993154239957
i= 19826497 fract= 0.3222084062581416
i= 19826498 fract= 0.3626179973871331
i= 19826499 fract= 0.4030280888255220


Da sich ja a/i immer um einen nahezu konstanten Wert ändert, ändern sich die Nachkommastellen auch um einen Nahezu konstanten Wert.

Aus dieser Änderung können wir nun schonmal eine wichtige Vorhersage treffen:
Damit a%i==0 (also teilbar) muss ja gelten fmodl(a/i,1)=0 (also keine Nachkommastellen).

Da sich fract (also fmodl(a/i,1)) immer um einen (nahezu) konstanten Wert ändert können wir ja eigentlich darüber auch eine Vorhersage treffen, "in wievielen i es an der 0 vorbeigeht" - also das i bestimmen für das gilt fmodl(a/i,1)=0.9... und fmodl(a/(i+1),1)=0.0...

Alle i dazwischen brauchen wir garnicht erst zu betrachten, weil die Nachkommastelle bei der Division a/i ungleich 0 ist (wenn a/i=...,8 und wir wissen, dass sich die Nachkommastelle mit jedem i um -0.2 ändert brauchen wir i+1, i+2 und i+3 ja garnicht erst untersuchen).

Diese Vorhersage hab ich mittels

deltaprej=(1.0F-fract)/(1.0F+fdeltafract);
deltaprej=(1.0F-fract)/(1.0F+fdeltafract+fdeltadelta*deltaprej);


getroffen. Die 2. Zeile ist dabei erstmal uninteressant, die macht die Ergebnsse nur leicht genauer (nach der habichvergessenregel gilt, dass die Summe aller Ableitungen einer Funktion (also f'+f''+f'''...) gegen die Funktion+c selbst strebt).

Die 1. Zeile liefert uns schon eine gute Vorhersage (für bestimmte Zahlenbereiche), wo die Nachkommastelle wieder von .9 in .0 umschlägt.
Das sieht dann so aus:


i= 19826480 deltaprej= 9.0254899533891951
i= 19826481 deltaprej= 8.0255955394189947
i= 19826482 deltaprej= 7.0256887398161795
i= 19826483 deltaprej= 6.0257695546862856
i= 19826484 deltaprej= 5.0258379848550709
i= 19826485 deltaprej= 4.0258940306081040
i= 19826486 deltaprej= 3.0259376924110001
i= 19826487 deltaprej= 2.0259689709094323
i= 19826488 deltaprej= 1.0259878663889270
i= 19826489 deltaprej= 0.0259943793150640
i= 19826490 deltaprej= 23.7679018828794378
i= 19826491 deltaprej= 22.7681900670391922
i= 19826492 deltaprej= 21.7684658623179095
i= 19826493 deltaprej= 20.7687292693610986
i= 19826494 deltaprej= 19.7689802884543148
i= 19826495 deltaprej= 18.7692189200630710
i= 19826496 deltaprej= 17.7694451646528675
i= 19826497 deltaprej= 16.7696590228691747
i= 19826498 deltaprej= 15.7698604951774762
i= 19826499 deltaprej= 14.7700495816832717
i= 19826500 deltaprej= 13.7702262833919727
i= 19826501 deltaprej= 12.7703906002290627
i= 19826502 deltaprej= 11.7705425330199499
i= 19826503 deltaprej= 10.7706820822300594
i= 19826504 deltaprej= 9.7708092481448064
i= 19826505 deltaprej= 8.7709240312295856
i= 19826506 deltaprej= 7.7710264321297810
i= 19826507 deltaprej= 6.7711164513107700
i= 19826508 deltaprej= 5.7711940888779061
i= 19826509 deltaprej= 4.7712593458365547
i= 19826510 deltaprej= 3.7713122222920435
i= 19826511 deltaprej= 2.7713527187096959
i= 19826512 deltaprej= 1.7713808357348387
i= 19826513 deltaprej= 0.7713965738327786
i= 19826514 deltaprej= 24.5055110617278750
i= 19826515 deltaprej= 23.5058082948643504
i= 19826516 deltaprej= 22.5060931425722256
i= 19826517 deltaprej= 21.5063656054966727
i= 19826518 deltaprej= 20.5066256839230244
i= 19826519 deltaprej= 19.5068733784964337




Wie man sieht können wir schon bei 19826480 vorhersagen, dass der Umschlag der Nachkommastelle bei i+9 erfolgen wird - und das tut es auch. Bis hierhin hab ich den Alg schon implementiert, aber leider sind kleine Fehler drinnen.

Wir können aber noch eine (viel bessere) Vorhersage treffen: Der Umschlag erfolgte ja in 9.0254899533891951 - daraus können wir Ableiten, dass wir die .0 nicht genau treffen werden, weil i ja eine ganze Zahl sein soll - also ist i+9 schonmal kein Teiler. Und das wissen wir schon 9 Zahlen vorher, bzw. wie man am weiteren Verlauf von Deltaprej sieht in dem Bereich schon ~24 Berechnungen vorher.

So - nun ist nochetwas auffällig ;D

i= 19826489 deltaprej= 0.0259943793150640
i= 19826490 deltaprej= 23.7679018828794378

i= 19826513 deltaprej= 0.7713965738327786
i= 19826514 deltaprej= 24.5055110617278750

Die Differenz der Nachkommastellen sieht nämlich so aus:

0.23408... für i= 19826489/+1 und
0.26588... für i= 19826513/+1

Die geringe Änderung ist kein Zufall. Von jetzt an wird der Alg reine Theorie, hab also noch nix realisiert.
Über die Differenz können wir ja auch vorhersagen wo deltaprej .00 wird (wenigstens ungefähr, später sicher noch genauer).

Das bedeutet, wir können über mehrere solcher Blöcke vorhersagen, wann ungefähr der nächste, wichtige Block kommt. Der wichtige Block sieht so aus:

i= 19827449 deltaprej= 23.9964763724516645
i= 19827450 deltaprej= 22.9967639308655706
i= 19827451 deltaprej= 21.9970392477387904
i= 19827452 deltaprej= 20.9973023233480782
i= 19827453 deltaprej= 19.9975531581480316
i= 19827454 deltaprej= 18.9977917527711024
i= 19827455 deltaprej= 17.9980181073161366
i= 19827456 deltaprej= 16.9982322227713002
i= 19827457 deltaprej= 15.9984340990575428
i= 19827458 deltaprej= 14.9986237369851465
i= 19827459 deltaprej= 13.9988011370086471
i= 19827460 deltaprej= 12.9989662992268121
i= 19827461 deltaprej= 11.9991192246277832
i= 19827462 deltaprej= 10.9992599131324203
i= 19827463 deltaprej= 9.9993883655509703
i= 19827464 deltaprej= 8.9995045823379123
i= 19827465 deltaprej= 7.9996085637698240
i= 19827466 deltaprej= 6.9997003103011553
i= 19827467 deltaprej= 5.9997798225642378
i= 19827468 deltaprej= 4.9998471008356089
i= 19827469 deltaprej= 3.9999021457475821
i= 19827470 deltaprej= 2.9999449577545677
i= 19827471 deltaprej= 1.9999755371330626
i= 19827472 deltaprej= 0.9999938843374479


Auffällig ist, dass die Nachkommastellen gegen 0 bzw. 1 streben - also wird a/(i+c) hier vermutlich genau auf .0 fallen, und das tut es auch für i=19827473.
Das dies wahrscheinlich so sein wird sehen wir schon am Anfang des Blockes.


Soviel zu dem Alg. Ist mal was anderes, und vor allem nicht so leicht zu implementieren. Das "dumme" sind die vielen kleinen Abweichungen, die es zu umschiffen gibt. Problematisch wird es zb. wenn fract nahe 0 oder nahe 1 liegt, dann liefert deltaprej Sprungvorhersagen von über einer Millionen, die aber wegen den ganz geringen Abweichungen meilenweit daneben liegen.

Ich hab mir auch nicht umsonst gerade diesen Zahlenbereich für i herausgepickt ;)
Problematisch sind auch Bereiche wo fract irgendwo um die 0.5 liegt, dann liefert deltaprej Sprungvorhersagen von 1 oder 2, und dann sparen wir nix ein - da müsste man dann über die Differenz der Nachkommastellen von Deltaprej von verschiedenen Blöcken rangehen, aber das hab ich wie gesagt noch überhauptnicht implementiert, weil der 1. Teil nichtmal fehlerfrei läuft. Momentan spare ich ca. 25% aller Berechnungen ein, wie gesagt ohne die Vorhersage über fmodl(deltaprej,1).

Ein wichtiger Punkt ist auch noch, dass aufeinanderfolgende Blocke immer in etwa gleich groß sind. Also wenn deltaprej bei i gleich 28 ist, ist es für i+29 auch ~28 (ganz leicht kleiner -> Nachkommastellen). Damit kann man dann einen Sprung über zb. die nächsten 10 Blöcke realisieren.

[BoMbY]

Tagchen,

also ich hab grad mal nen Algorithmus nach Pollards Rho Methode in Java auf meinen AthlonXP 2600+ ausprobiert:

twoFactors(1949565409282837) = 19827473 * 98326469 <=> 140ms
twoFactors(312133288776418481) = 398527487 * 783216463 <=> 344ms
twoFactors(11081452861688528621) = 3981527477 * 2783216473 <=> 203ms
twoFactors(6864244079350478607287) = 73981527493 * 92783216459 <=> 4594ms
twoFactors(138773836209286437923293) = 292783216493 * 473981527601 <=> 9828ms
twoFactors(61798910261151363714617153) = 7292783216489 * 8473981527577 <=> 95906ms

Also, es sieht so aus als verstehe der Mann was von seinem Job... :) Aber vieleicht läßt sich daran ja noch was optimieren... ;)

m.f.g.
BoMbY

[i_hasser]

Da gibts einen Haken - der Alg testet nicht den gesamten Bereich, sondern hört auf wenn ein Teiler gefunden wurde. Damit hängen die Zeiten stark von den Teilern ab, und nicht von der größe der Zahl.

[BoMbY]

Ja, stimmt schon ...

Aber Du weißt nicht zufällig eine vernünftige Ableitung f'(x) für f(x)=trunc(x) ? :)

m.f.g.
BoMbY

[i_hasser]

Klar. Nämlich f'(x==trunc(x))=n.d. und f'(x!=trunc(x))=0 ;)

[BoMbY]

Weil, die Hochpunkte von f(x)=(trunc(n/x)*x)/n sollten die die Primzahlen p und q für n=p*q sein. Übrigens ist f(p)=1 und f(q)=1, und wie's ausschaut für alle x > 0, x <> p und x <> q ist f(x) < 1. Problem ist halt die Ableitung von f(x) um deren Nullstellen zu bestimmen... :)

[i_hasser]

Mist! Das läuft auf mein Verfahren oben hinaus :-/

Hab ich Sarnagel schon geschrieben, wenn man diese Vorhersage oft genug wiederholt bekommt man irgendwann nur noch 2 Positionen heraus, und das sind dann die beiden Teiler (bzw. wenn man es nur bis zur Wurzel n treibt eben nur einen einzigen).

Hats also doch schon jemand entdeckt.

Aber wenigstens bin ich selbst drauf gekommen *baeh*


Das werd ich mir auf jeden Fall mal anschauen, und in den nächsten Tagen einen Alg dazu schreiben.

[Sargnagel]

Original geschrieben von BoMbY
Weil, die Hochpunkte von f(x)=(trunc(n/x)*x)/n sollten die die Primzahlen p und q für n=p*q sein. Übrigens ist f(p)=1 und f(q)=1, und wie's ausschaut für alle x > 0, x <> p und x <> q ist f(x) < 1. Problem ist halt die Ableitung von f(x) um deren Nullstellen zu bestimmen... :)
Nur mal so gefragt, damit ich der Diskussion noch folgen kann: Ist das Problem die trunc()-Funktion? Kann man die überhaupt als simple Formel ausdrücken, d.h. nicht numerisch ???

[BoMbY]

Original geschrieben von intel_hasser
Hats also doch schon jemand entdeckt.

Also diese f(x)-Formel hab ich mir auf Basis Deines Beitrags zusammengebastelt. Ich hab den ganzen Nachmittag damit rumgerechnet, aber alles was ich sagen kann ist: "Hier steh ich nun ich armer Tor und bin so klug als wie zuvor". Also wenn man das ableiten, und dann die nullgesetzte Ableitung noch zu x auflösen könnte, dann hätten wir (oder besser Du) viele neue Freunde (und vermutlich genausoviele neue Feinde), weil die Public-Key-Kryptographie mit einem mal wertlos wäre... ;D

m.f.g.
BoMbY

[i_hasser]

Ach so, ich dachte du hättest das woanders so gefunden ;)

Also ich hab mir den Verlauf mal angesehen. So ad hoc kommt mir da wieder eine Methode über die Ableitung in den Sinn, die kann man ja nährungsweise berechnen mit f(x+1)-f(x). Und bei großen Zahlen ändert die sich so schnell auch nicht wesentlich.

Mir ist aber aufgefallen, dass nach deiner Gleichung der 80bit Float der x86 FPU schon arg klein wird.

Jetzt hab ich leider keine Zeit daran weiterzudenken, bis morgen ;)

[BoMbY]

Original geschrieben von Sargnagel
Nur mal so gefragt, damit ich der Diskussion noch folgen kann: Ist das Problem die trunc()-Funktion? Kann man die überhaupt als simple Formel ausdrücken, d.h. nicht numerisch ???

Das ist genau die Frage, welche ich auch nicht beantworten kann ... Genau bei dieser Funktion hört's bei mir im Moment auf ... Also jedenfalls was Ableitungen angeht, und damit rechnen ist auch irgendwie Mist. Wie kann man z.B. sowas wie "x*trunc(n/x) = n" zu x auflösen? Ich vermute mal gar nicht, oder? ...

[Sargnagel]

Original geschrieben von BoMbY
Wie kann man z.B. sowas wie "x*trunc(n/x) = n" zu x auflösen? Ich vermute mal gar nicht, oder? ...
Hmmm ...

x*trunc(n/x) = n
--> trunc(n/x) = n/x
--> tja, was ist das Pendant zu trunc()?
(also wie Quadrat und Quadratwurzel ;) )

Das löst das Problem aber nicht ... *chatt* ... da wären wir wieder beim Thor, der so klug war wie zuvor ... *kopfkratz

[i_hasser]

hmmm... also die Formel ist ja eigentlich schnell hergeleitet (ich kenn sowas, manchmal formt man um und formt um und landet doch wieder bei einer Gleichung, die man mit einem Schritt bekommen hätte).

Wenn n/p=q für p,q element N, dann muss ja n/p-trunc(n/p)=0 sein, weil n/p keine Nachkommastellen geben darf.
Und n/p-trunc(n/p)=0 ist ja gleich trunc(n/p)*p=n bzw. trunc(n/p)*p/n=1.

Also sind wir wieder da, wo wir angefangen haben :(

Trunc kann man meines Wissens nicht so einfach mathematisch ausdrücken, weil überhaupt "Stellen" bei einer rein mathematischen Betrachtung nicht existieren - ob man das nun mit Dezimalzahlen, Dualzahlen oder Oktalzahlen durchrechnet macht ja keinen Unterschied, das Ergebniss ist immer das selbe. Will man aber die 10. Stelle nach dem Komma auslöschen geht das eben nicht so einfach, weil die sich je nach Zahlensystem immer verschiebt.

Es müsste also schonmal die Basis 10 irgendwo auftauchen.

[BoMbY]

Tja, wenn wir sagen "trunc(n/x) = n/x" dann folgt daraus:

f(x) = (trunc(n/x)*x)/n

f(x) = ((n/x)*x)/n

f(x) = n/n

f(x) = 1

also Quatsch...

Aber nu bin ich erst mal pennen, bis morgen dann und gute Nacht ...

[i_hasser]

trunc(n/x)=n/x gilt ja nur für x=p||x=q

[i_hasser]

*lol*

Ich habs.

sin(pi*x)+sin(pi*n/x)==0 für n=zu zerlegende Zahl.

Eine Lösung dieser Gleichung ist das gesuchte Ergebnis... dummerweise ist die Gleichung symbolisch unlößbar *chatt*

Weis zufällig jemand wie man mathematisch x=0 und y=0 in einer Gleichung unterbringen kann?

[Marnem]

Original geschrieben von intel_hasser
*lol*

Ich habs.

sin(pi*x)+sin(pi*n/x)==0 für n=zu zerlegende Zahl.

Eine Lösung dieser Gleichung ist das gesuchte Ergebnis... dummerweise ist die Gleichung symbolisch unlößbar *chatt*

Weis zufällig jemand wie man mathematisch x=0 und y=0 in einer Gleichung unterbringen kann?

hmm, ich seh da kein y. Typo?

[i_hasser]

Wieso y? Das ist doch keine Funktion ;)

Du musst die Sache nur symbolisch nach x auflösen. Mehr ist es nicht.

[Blackspeed]

Was Primzahlen angeht, habt ihr keine Chance. Ein Japaner hat mit einem Programm 1,2 Billionen Nachkommastellen berechnen lassen.

[[ab]noname]

Bei Primzahlen gibts 2 Wege bei der Berechnung: entweder optimiert für riesen große Zahlen oder optimiert um ganz schnell kleinere zu berechnen.

Ob nun einer 1,2 bill. Stellen berechnet oder nicht ist irrelevant. Sein Algorithmus wird nicht erheblich schneller sein als standard Varianten. Er hatte einfach Glück, dass seine nächste Zahl eine war...

EDIT:
Ganz groß erreicht man z.B. schon über die Folge: p_(n+1) = 2 ^ p_(n) - 1 (p aus Primzahlen)
Nur sind die Sprünge so enorm, dass man nur einen Bruchteil der Primzahlen erfasst!

[i_hasser]

Manche Leute haben den Sinn von den Progg Wettbewerben anscheinend noch nicht begriffen :]

Abgesehen davon - was meinst du damit? Eine Primzahl hat keine Nachkommastellen, und 1.2 Billionen halte ich für unrealistisch - Milliarden heist im englischen "billion" ;)