Goldbach Vermutung



  • Das hier läuft schon deutlich fixer:

    BigUnsigned mulmod(BigUnsigned a, BigUnsigned b, BigUnsigned mod)
    {
    	/////////////
    	QueryPerformanceCounter((LARGE_INTEGER*)&startCount);
    	/////////////
    	BigUnsigned x = 0;
    	BigUnsigned y = a % mod;
    
    	while (b > 0)
    	{
    		if (((uint32_t)b.getBlock(0)) & 1)
    		{
    			x += y;
    			x %= mod;			
    		}
    
    		y <<= 1;
    		y %= mod;
    
    		b >>= 1;		
    	}
    	/////////////
    	QueryPerformanceCounter((LARGE_INTEGER*)&lastCount);
    	sumTime_Function += (lastCount - startCount);
    	/////////////
    	return x;
    }
    


  • Erhard Henkes schrieb:

    Das hier läuft schon deutlich fixer:

    Cool-O-Matik! 😃

    Gibt wahrscheinlich noch ne Reihe andere Stellen im Code, die davon profitieren würden. Ist nur ein bissl blöd dass das der Code-Lesbarkeit schadet. Das ist für mich eher so C-Hacker-Stil von jemandem der grad von Assembler kommt 🙄



  • Deine Variante ist 8% schneller, daher baue ich sie gerne ein. Insgesamt ein großer Schritt. 👍

    Also neuer Sourcecode: http://codepad.org/ReYuzGjk (nur mulmod)



  • ... auch an anderen Stellen verändert: http://codepad.org/f5zbHebm

    @Finnegan: Danke! 👍

    mulmod vorher: 70% jetzt: 45 % (bei 10^54)

    Es lohnt sich also, weiter an dem Interface BigInteger und mulmod zu tunen.



  • Erhard Henkes schrieb:

    mulmod vorher: 70% jetzt: 45 % (bei 10^54)

    Gern. Schön wenn man so ne Funktion hat, die im Profiler leuchtet wie ein Weihnachtsbaum. Da ist das Optimieren leicht.
    Leider hat man auch oft mit Programmen zu tun wo alle Funktionen irgendwie gleich "langsam" sind, das ist immer ein Krampf.
    ...aber bei so ner schönen fetten spanischen Galleone schraub ich mir glatt mein Holzbein an und leg' ne Überstunde ein 😃



  • Vielleicht könnte man auch eine eigene MyBigInteger-Klasse auf Basis uint64_t anstelle uint32_t basteln, nur mit den Funktionen, die man hier braucht. Dann könnte man das auf Geschwindigkeit genau dafür trimmen.



  • Was ist das Ziel des ganzen hier? Willst du so weit kommen, wie du nur kannst? Dann würde ich jetzt eher aufhören, zu mikrooptimieren, sondern mir die BigInt-Implementierung ansehen. BigInts gibt's viele, aber richtig gute nicht. Wenn du's selber machen willst, dann würde ich sie auf doppelter Wortbreite basieren lassen und ebenfalls im Zweierkomplement darstellen.

    Die Addition muss mit ziemlicher Sicherheit in Assembler gemacht werden, da es dort Mnemonics wie ADC gibt, die einem das explizite Carryen sparen (sofern der Compiler es nicht schon so generiert).

    Danach Multiplikation: da kannst du mit Karazuba oder gleich mit Toom-Cook anfangen. Oder du machst vorher etwas Mathe, lädst dir die Fast Fourier Transform runter und implementierst Schönhage-Strassen. Die Algorithmen lassen sich auch allesamt parallelisieren + cache-optimieren.



  • Hallo

    Erhard Henkes schrieb:

    Ich sehe bezüglich weiterer Optimierungen folgende Felder:
    ...

    Wenn ich den Code richtig verstehe, beginnt Ihr bei jeder Zahl wieder vorn mit den Test potentieller Primzahlpärchen. Bei den großen Zahlen wie 10^20 oder sogar 10^240 sind die gefundenen großen Primzahlen doch viel zu wertvoll, um sie gleich wieder wegzuwerfen.

    1.) Wenn ich z.B. für x1 = 13 + (x1-13) gefunden habe, weil x1-13 eine Primzahl ist, dann kann ich für alle Primzahlen größer 13 aus meiner Primzahltabelle weiter PrimePairs bilden und die aus der Liste der zu prüfenden Zahlen streichen.
    Also:
    x3 = (x1+4) = 17 + (X1-13)
    x4 = (x1+6) = 19 + (X1-13)
    u.s.w.
    2.) für x2 = (x1+2) passt das noch nicht, aber auch hier gibt es eine Optimierung:
    beim Prüfen von x1 habe ich für x1-3, x1-5, x1-7, x1-11 schon festgestellt, dass das keine Primzahlen sind. Für x2 verschiebt sich das um 2, somit kann ich mir ein paar Primzahltests gleich sparen: von x2-5, x2-7, x2-13 weiß ich schon, dass das keine Primzahlen sind. Hierfür kann man sich intern eine Tabelle halten für alle ungeraden Zahlen von endNumber-3 (Index 0) bis startNumber-n mit jeweils eine Flag für "Prime", "not Prime", oder "Unbestimmt" (evtl. noch ein "wird gerade berechnet" fürs parallelisieren). Natürlich muss die Tabelle irgendwie begrenzt werden, aber ich schätze mal, dass man die meisten PrimPairs findet, bevor die Tabelle unhandlich groß wird.



  • Erhard Henkes schrieb:

    Vielleicht könnte man auch eine eigene MyBigInteger-Klasse auf Basis uint64_t anstelle uint32_t basteln, nur mit den Funktionen, die man hier braucht. Dann könnte man das auf Geschwindigkeit genau dafür trimmen.

    Ich würde so frech sein und behaupten, dass der Compiler das schon macht soweit er kann. Zumindest bei dem was ist schon in von Visual C++ generiertem Maschinencode gesehen habe, würde ich sogar davon ausgehen, dass er, wenn z.B. die entsprechende Architektur aktiviert ist (bspw. -arch:avx oder -arch:avx2 - SSE2 ist per default an), er sogar schon entsprechende Vektor-Instruktionen erzeugt... aber probieren geht über studieren. Hab zwar jetzt keine Zeit mehr dafür, aber ein kurzes Mini-Testprogramm und dann im Debug-Modus mit Code-Optimierungen und mit Breakpoint auf die Operation auf "Show Disassembly" gehen, sollte Aufschluss geben (xmm/ymm-Register und vpadd*-Instruktionen wären z.B. bei Addition ein Hinweis... wichtig bei so einem Test: Variablen mit denen man sowas testet sollten nicht durch irgendwelche Konstanten/Literale gefüllt werden, sonst rechnet der das Ergebnis einfach beim Kompilieren aus und baut es direkt in die EXE ein [Stichw. Konstantenfaltung] 😃 ... ich hab mir da immer mit dem volatile -Keyword ausgeholfen wenn ich wirkliche Instruktionen sehen wollte).

    Falls man es wirklich selbst implementieren will, könnte es allerdings mehr bringen, dem CPU-Cache besser entgegenzukommen. Da muss man allerdings einiges an Hirnschmalz reinstecken, auch wenn die Grundregel simpel ist: Möglichst linear auf den Speicher zugreifen und keine wilden Pointer-Sprünge quer durch den RAM. Daher habe ich auch vorhin bemerkt dass die BigInts alle auf dem Heap liegen (theoretisch kreuz- und quer verteilt). Je nach deren grösse macht es Sinn die in einem Array auf dem Stack zu haben, dann sind z.B. die lokalen BigInts auf denen man gerade arbeitet im (Stack-)RAM "nah" beieinander und man kann besser von den CPU Cache-Leveln profitieren. Die std::strings haben wegen sowas z.B. die SSO (Small String Optimization). Kurze Strings liegen auf dem Stack, erst wenn sie länger werden wird der Heap genommen. Eine richtig gute BigInt-Library würde das wahrscheinlich ähnlich machen.

    Finnegan

    P.S.: Noch eine einfachere Möglichkeit, mit der man noch ein bisschen rausquetschen könnte: Whole Program Optimization und Link Time Code Generation. Wenn die cpp-Dateien der BigInt-Lib Teil des Projekts sind kann er dann auch aus diesen Funktionen inlinen oder die Nutzung der CPU-Register für das ganze Programm optimieren. Es gibt Fälle wo sowas durchaus nochmal um die 10% bringen kann, wenn mn Glück hat



  • Bei den großen Zahlen wie 10^20 oder sogar 10^240 sind die gefundenen großen Primzahlen doch viel zu wertvoll, um sie gleich wieder wegzuwerfen.

    Ja, das ist bei großen Zahlen völlig richtig.

    Was ist das Ziel des ganzen hier? Willst du so weit kommen, wie du nur kannst? Dann würde ich jetzt eher aufhören, zu mikrooptimieren, sondern mir die BigInt-Implementierung ansehen.

    Verbesserte Abläufe, mulmod und BigUnsigned sind die Brennpunkte. So etwas in der Art habe ich bisher noch nicht unternommen. Ich finde es interessant zu sehen, was man heute auf dem eigenen PC untersuchen kann und welche Optimierungsstrategien dabei wichtig sind.

    Eine Widerlegung durch Finden einer Anti-Goldbach-Zahl erwarte ich nicht mehr, aber genau das versucht das Programm. 😃 😉



  • Ich habe den Vorschlag von DJohn, hohe Primzahlen zu speichern und vor dem Miller-Rabin per Lookup zu checken, mal auf die Schnelle eingefügt:
    Den Erfolg kann man durch Aktivierung von // cout << " LH "; in isPrime(...) überprüfen.

    #define NOMINMAX     // due to conflict with max()
    #include <windows.h> // warning: no C++ standard!
    
    #include <iostream>
    #include <iomanip>
    #include <cstdint>
    #include <cstdlib>
    #include <cassert>
    #include <algorithm>
    #include <cmath>
    #include <vector>
    #include <unordered_set>
    #include <fstream>
    #include <thread>
    #include <mutex>
    #include "BigIntegerLibrary.hh" // https://mattmccutchen.net/bigint/
    
    using namespace std;
    
    ///#define MILLER_RABIN_LIAR_CHECK
    #define FUNCTION_TIME_ANALYSIS
    
    const BigUnsigned startNumber = stringToBigUnsigned("1000000000000000000000000000000000000000000000000000000000000");
    const BigUnsigned endNumber = stringToBigUnsigned("1000000000000000000000000000000000000000000000000000000000100");
    const int numberOfThreads = 11; // count of even numbers from start to end = (end - start)/2 + 1
    
    BigUnsigned hiPrimes[100000];   // calculated primes (high values)
    uint64_t countHiNumbersToBeTested = 0;
    uint64_t countSuccessfulDivisionPrecheck = 0;
    uint64_t countSavedMR = 0;
    uint64_t countMR = 0;
    
    uint64_t frequency;
    
    #ifdef FUNCTION_TIME_ANALYSIS
    	uint64_t startCount=0, lastCount=0;
    	uint64_t sumTime_Function = 0;
    #endif
    
    void setConsole()
    {
    	_CONSOLE_SCREEN_BUFFER_INFO info;
    	HANDLE handle = GetStdHandle(STD_OUTPUT_HANDLE);
    	GetConsoleScreenBufferInfo(handle, &info);
    	COORD c;
    	c.X = std::max<SHORT>(info.dwSize.X, 150);
    	c.Y = std::max<SHORT>(info.dwSize.Y, 1000);
    	SetConsoleScreenBufferSize(handle, c);
    
    	SMALL_RECT RECT;
    	RECT.Left = 0;
    	RECT.Top = 0;
    	RECT.Bottom = std::max(info.srWindow.Bottom - info.srWindow.Top, 40 - 1);
    	RECT.Right = std::max(c.X - 1, 100 - 1);
    	SetConsoleWindowInfo(handle, TRUE, &RECT);
    }
    
    void textcolor(unsigned short color = 15)
    {
    	SetConsoleTextAttribute(::GetStdHandle(STD_OUTPUT_HANDLE), color);
    }
    
    void wait()
    {
    	cout << "Press any key to continue." << endl;
    	cin.clear();
    	cin.ignore(numeric_limits<streamsize>::max(), '\n');
    	cin.get();
    }
    
    mutex mutex0;
    mutex mutex1;
    mutex mutex2;
    mutex mutex3;
    mutex mutex4;
    
    unordered_set<uint64_t> setLiar;
    uint64_t lastLiar = 0;
    
    uint64_t convertBigToU64(BigUnsigned num)
    {
    	BigUnsigned ZweiHochVierUndSechzig = stringToBigUnsigned("18446744073709551616");
    	assert(num < ZweiHochVierUndSechzig);
    	uint64_t jHi = num.getBlock(1);
    	uint64_t jLo = num.getBlock(0);
    	return (jHi << 32 | jLo);
    }
    
    void convertU64ToBig(uint64_t num, BigUnsigned& b)
    {
    	uint32_t numLo = (uint32_t)(num & 0xFFFFFFFF);
    	uint32_t numHi = (uint32_t)((num >> 32) & 0xFFFFFFFF);
    	b.setBlock(0, numLo);
    	b.setBlock(1, numHi);
    }
    
    BigUnsigned mulmod(BigUnsigned a, BigUnsigned b, const BigUnsigned& mod)
    {
    #ifdef FUNCTION_TIME_ANALYSIS
    	QueryPerformanceCounter((LARGE_INTEGER*)&startCount);
    #endif
    
    	BigUnsigned x = 0;
    	BigUnsigned& y = a;
    
    	while (b > 0)
    	{
    		if (b.getBit(0))
    		{
    			x += y;
    			x %= mod;
    		}
    		y <<= 1;
    		y %= mod;
    		b >>= 1;
    	}
    
    #ifdef FUNCTION_TIME_ANALYSIS
    	QueryPerformanceCounter((LARGE_INTEGER*)&lastCount);
    	sumTime_Function += (lastCount - startCount);
    #endif
    
    	return x;
    }
    
    BigUnsigned powmod_Volkard(BigUnsigned base, BigUnsigned exp, BigUnsigned modul)
    {
    	//rechnet 'base hoch exp mod modul'
    	BigUnsigned a1 = base, z1 = exp, x = 1, m = modul;
    	while (z1 != 0)
    	{
    		while (!z1.getBit(0))
    		{
    			z1 >>= 1;
    			a1 = mulmod(a1, a1, m);
    		}
    		--z1;
    		x = mulmod(x, a1, m);
    	}
    	return x;
    }
    
    bool IsProbablyPrime_MillerRabinOptimized(BigUnsigned n)
    {
        mutex4.lock();
        countMR++;
        mutex4.unlock();
    
    	/*
    	example:
    	step: "calculate s and d"
    	n = 1001
    	n-1 = 1000
    	s = 3
    	d = 125
    	1000 / 2^3 = 125
    
    	step: "test with base a"
    	ad = a^d mod n  (cf. function powmod)
    	...
    	*/
    
    	if (n == 2)
    		return true;
    	// "s" is the logarithm to base base 2 of the biggest number 2^s dividing n-1
    	// d = (n-1)/(2^s)
    	// calculate d and s
    	BigUnsigned d = n - 1;
    	BigUnsigned s = 0;
    
    	while ((d & 1) == 0)
    	{
    		++s; d >>= 1;
    	}
    
    	// We do not use random base a = {2 ... n-1}, but fix it to a = 2
    	BigUnsigned a = 2; // base ("witness")
    	BigUnsigned ad = (powmod_Volkard(a%n, d, n)); // (a^d) mod n
    
    	if (ad == 1)
    		return true; // a^d mod n == 1
    	if ((s > 0) && (ad == (n - 1)))
    		return true; // a^d mod n == -1 (tests (a^d)^(2^0))
    
    					 // test for r = {0 ... s-1}
    	for (BigUnsigned r = 1; r < s; ++r)
    	{
    		// ad is currently ad^(2^(r-1)), thus we square it to get ad^(2^r))
    		ad = mulmod(ad, ad, n);
    		if (ad == (n - 1))
    			return true; // tests (a^d)^(2^(r+1)))
    	}
    	return 0; // false, composite
    }
    
    bool IsPrimeDivisionTest(BigUnsigned number)
    {
    	if (number<2)    return false;
    	if (number == 2)   return true;
    	if (number % 2 == 0) return false;
    	for (BigUnsigned i = 3; i*i <= number; i += 2)
    	{
    		if (number%i == 0) return false;
    	}
    	return true;
    }
    
    bool IsPrime(vector<bool>& primes, BigUnsigned number, BigUnsigned primeLimit, unordered_set<uint64_t>& setLiar)
    {
    	// lookup from bitset (in the RAM) - low primes
    	if (number <= primeLimit)
    		return primes[number.toUnsignedLong()];
    
    	// lookup from bitset (in the RAM) - high primes
    	if ((endNumber - number) <= 100000)
    	{
    		uint32_t element = (endNumber - number).toUnsignedLong();
    		if (hiPrimes[element] == 1)
    		{
    			// cout << " LH Prime: endN - " << element << " ";
    			mutex3.lock();
    			countSavedMR++;
    			mutex3.unlock();
    			return true;
    		}
    		if (hiPrimes[element] == 0)
    		{
    			// cout << " LH No-Prime: endN - " << element << " ";
    			mutex3.lock();
    			countSavedMR++;
    			mutex3.unlock();
    			return false;
    		}
    	}
    
    	// division test with first prime numbers
    	const uint32_t DIVISOR_PRECHECK_MAX = 1999993;
    	for (uint32_t divisor = 3; divisor <= DIVISOR_PRECHECK_MAX; ++divisor)
    	{
    		if (primes[divisor])
    		{
    			if (number % divisor == 0)
    			{
    				if ((endNumber - number) <= 100000)
    				{
    					uint32_t element = (endNumber - number).toUnsignedLong();
    					hiPrimes[element] = 0;
    				}
    				mutex2.lock();
    				countSuccessfulDivisionPrecheck++;
    				mutex2.unlock();
    				return false;
    			}
    		}
    	}
    
    	// Miller-Rabin, not perfect, but fast
    	if (!IsProbablyPrime_MillerRabinOptimized(number))
    	{
    		if ((endNumber - number) <= 100000)
    		{
    			uint32_t element = (endNumber - number).toUnsignedLong();
    			hiPrimes[element] = 0;
    		}
    		return false;
    	}
    	else
    	{
    		if (lastLiar == 0)
    		{
    			return true;
    		}
    
    		BigUnsigned bigLastLiar;
    		BigUnsigned& refBigLastLiar = bigLastLiar;
    		convertU64ToBig(lastLiar, refBigLastLiar);
    
    		if (number > bigLastLiar)
    		{
    			if ((endNumber - number) <= 100000)
    			{
    				uint32_t element = (endNumber - number).toUnsignedLong();
    				hiPrimes[element] = 1;
    			}
    			return true;
    		}
    		else
    		{
    			if (setLiar.find(convertBigToU64(number)) != setLiar.cend())
    			{
    				// cout << "ALARM: Miller-Rabin-Liar found: " << number << "  ";
    				// wait();
    				return false; // function tells the truth by failure list (0 ... 2^64), above lastLiar there will be "Prime Hell" again.
    			}
    			return true;
    		}
    	}
    }
    
    void PrimePairFound(vector<bool>& primes, BigUnsigned i, const BigUnsigned primeLimit, unordered_set<uint64_t>& setLiar, bool* retVal)
    {
    	*retVal = false;
    
    	BigUnsigned grenze = (i >> 1);
    	BigUnsigned a = 0;
    	for (a = 3; a <= grenze; a += 2)
    	{
    	    bool retValLo = IsPrime(primes, a, primeLimit, setLiar);
    	    if (retValLo == false)
    	    {
    	        continue;
    	    }
    
            mutex1.lock();
    	    countHiNumbersToBeTested++;
    	    mutex1.unlock();
    
    	    bool retValHi = IsPrime(primes, i - a, primeLimit, setLiar);
    
    		if (retValLo && retValHi)
    		{
    			*retVal = true;
    			break;
    		}
    	}
    
    	mutex0.lock();
    	cout << i << "  " << setw(5) << a << endl;
    	mutex0.unlock();
    }
    
    int main()
    {
    	setConsole(); // windows.h // big console for large numbers! ^^
    	textcolor(14); // windows.h // 10 = light green, 14 = light yellow
    
    	QueryPerformanceFrequency((LARGE_INTEGER*)&frequency);
    	cout << "cpu frequency: " << frequency << " Hz" << endl << endl;
    
    	uint64_t const primeLimit = 2000000; //200000000000;
    	cout << "We generate vector<bool>(primes) up to: " << primeLimit << endl; // (this amount depends on the capability of your RAM)
    
    	vector<bool> primes(primeLimit); // calculated primes (low values)
    
        // calculated primes (high values)
    	for (int element = 0; element < 100000; ++element)
    	{
    		if ((element & 1))
    			hiPrimes[element] = 2; // not determined
    		else
    			hiPrimes[element] = 0; // no prime
    	}
    	// endNumber (index 0) to endnumber-100000 (index 100000)
    	// 0: no prime
    	// 1: (probably) prime
    	// 2: not determined
    	// 3: just being calculated
    
    	time_t startTime, endTime;       // measuring execution time
    	uint64_t global_count = 0;       // total number of primes
    	time(&startTime);
    
    	// initialize the number array
    	primes[0] = false;
    	primes[1] = false;
    	for (uint64_t i = 2; i<primeLimit + 1; ++i)
    	{
    		primes[i] = true;
    	}
    
    	// sieving loop
    	uint64_t iMax = (uint64_t)sqrt((double)primeLimit) + 1;
    	for (uint64_t i = 2; i<iMax; ++i)
    	{
    		if (primes[i])
    		{
    			for (uint64_t j = 2 * i; j<primeLimit + 1; j += i)
    			{
    				primes[j] = false;
    			}
    		}
    	}
    
    	time(&endTime);
    
    	// count the number of primes
    	for (uint64_t i = 0; i<primeLimit + 1; ++i)
    	{
    		if (primes[i])
    			global_count++;
    	}
    
    	cout << "In " << difftime(endTime, startTime) << " s we found " << global_count
    		<< " prime numbers (2 to " << primeLimit << " )." << endl << endl;
    
    #ifdef MILLER_RABIN_LIAR_CHECK
    	// Prepare Miller-Rabin "Liar-Liste" and detect Goldbach prime pairs with Miller-Rabin test
    	ifstream fileLiar("F:\\Daten\\Miller_Rabin_Liar\\2-SPRP-2-to-64.txt", fstream::in); // https://miller-rabin.appspot.com/
    	cout << "Miller-Rabin Liar file starts to be read from file." << endl;
    	uint64_t count = 0;
    	uint64_t u64LiarNumber = 0;
    	while (!fileLiar.eof())
    	{
    		fileLiar >> u64LiarNumber;
    		// cout << u64LiarNumber << "\t";
    		setLiar.insert(u64LiarNumber);
    		count++;
    	}
    	fileLiar.close();
    	lastLiar = u64LiarNumber;
    	cout << "Miller-Rabin Liar file with " << count << " liars is read." << endl << "Last liar: " << u64LiarNumber << endl << endl;
    #endif // MILLER_RABIN_LIAR_CHECK
    
    	unordered_set<uint64_t>& refSetLiar = setLiar;
    
    #ifdef MILLER_RABIN_LIAR_CHECK
    	cout << "Now Goldbach first prime number pair will be detected:" << endl;
    #else
    	cout << "Now Goldbach probable (no liar list available) first prime number pair will be detected:" << endl;
    #endif
    
    	uint64_t zeit1, zeit2, zeit2_old, delta;
    	QueryPerformanceCounter((LARGE_INTEGER*)&zeit1);
    	zeit2 = zeit1;
    	uint64_t sumTime = 0;
    
    	for (BigUnsigned i = startNumber; i <= endNumber; i += (2 * numberOfThreads))
    	{
    		bool retVal[numberOfThreads];
    		vector<thread> t;
    
    		for (int x = 0; x<numberOfThreads; ++x)
    		{
    			retVal[x] = false;
    			if ((i+2*x) <= endNumber)
    				t.emplace_back(PrimePairFound, std::ref(primes), i+2*x, BigUnsigned((unsigned long)primeLimit), std::ref(refSetLiar), retVal + x);
    		}
    
    		for (int x = 0; x < numberOfThreads; ++x)
    		{
    			if ((i + 2 * x) <= endNumber)
    			{
    				t[x].join(); // we have to wait for the result of the thread!
    				if (retVal[x] == false)
    					cout << "Counterevidence found for this number: " << i + 2 * x << endl;
    			}
    		}
    
    		zeit2_old = zeit2;
    		QueryPerformanceCounter((LARGE_INTEGER*)&zeit2);
    		delta = (zeit2 - zeit2_old);
    		sumTime += delta;
    		cout << "time: " << 1000 * delta / frequency << " ms" << endl << endl;
    	}
    
    	uint64_t totalTimeMS = 1000 * sumTime / frequency;
    	cout << "total time: " << totalTimeMS << " ms" << endl << endl;
    
    #ifdef FUNCTION_TIME_ANALYSIS
    	uint64_t totalTimeAtFunctionMS = 1000 * sumTime_Function / frequency;
    	cout << "total time at mulmod: " << totalTimeAtFunctionMS << " ms (" << 100 * totalTimeAtFunctionMS / totalTimeMS << " %)" << endl << endl;
    #endif
    
        cout << "High Numbers to be tested: " << countHiNumbersToBeTested << endl;
        cout << setprecision(3);
        cout << "Number of successful Division Prechecks: " << countSuccessfulDivisionPrecheck << " (" << 100 * (double)countSuccessfulDivisionPrecheck/(double)countHiNumbersToBeTested << " %)" << endl;
        cout << "Number of MR tests: " << countMR << " (" << 100 * (double)countMR/(double)countHiNumbersToBeTested << " %)" << endl;
        cout << "Number of saved MR tests by high prime lookup: " << countSavedMR << " (" << 100 * (double)countSavedMR/(double)countHiNumbersToBeTested << " %)" << endl;
    
    	wait();
    }
    

    `from 10^60 to 10^60 + 100:

    High Numbers to be tested: 1954

    Number of successful Division Prechecks: 547 (28 😵

    Number of MR tests: 104 (5.32 😵

    Number of saved MR tests by high prime lookup: 1303 (66.7 😵

    `

    `from 2 * 10^60 to 2 * 10^60 + 1000:

    High Numbers to be tested: 46574

    Number of successful Division Prechecks: 1951 (4.19 😵

    Number of MR tests: 752 (1.61 😵

    Number of saved MR tests by high prime lookup: 43871 (94.2 😵

    `

    Die Idee des "high prime lookup" spart gerade beim Testen vieler Zahlen in einem hohen Zahlenbereich eine Vielzahl an Miller-Rabin-Tests.

    Thanks to DJohn. 👍



  • Für weitere Optimierungen nehme ich den Goldbach-Test (Geradzahl, einfaches Primzahlenpaar) bei 2*10^120 bis 2*10^120 + 1000 als Basis.
    `

    total time: 426900 ms

    total time at mulmod: 248637 ms (58 😵

    High Numbers to be tested: 98055

    Number of successful Division Prechecks: 2867 (2.92 😵

    Number of MR tests: 788 (0.804 😵

    Number of saved MR tests by high prime lookup: 94400 (96.3 %)`
    komplett: http://pastebin.com/ziNH4sna

    Zum eventuellen Herumschrauben an der freien BigInteger Library habe ich diese auf eine einfache h/cpp-Kombination konzentriert und auf das Wesentliche reduziert: http://henkessoft.de/Sonstiges/MyBigUnsigned.rar
    Ich bin allerdings nicht sicher, ob das Sinn macht.



  • Erhard Henkes schrieb:

    Die Idee des "high prime lookup" spart gerade beim Testen vieler Zahlen in einem hohen Zahlenbereich eine Vielzahl an Miller-Rabin-Tests.

    Thanks to DJohn. 👍

    Auch wenn ich mich wieder unbeliebt mache, schau mal nach "fensterweise" auf der ersten Seite dieses Threads. Ich fürchte, dahin geht es zunächst. Wir werden festestellen, daß bis zu grotesken Größen die beiden äußeren Fester reichen; evenuelle Zweifler kann man dann ja noch extra untersuchen.

    Ich denke gerade darüber nach, statt der äußeren Fenster nur ein einziges inneres Fenster um sqrt(n) zu führen. Dann könnte man bis 2^128 untersuchen obwohl man nur isPrime(UInt64) hat.

    An einem schnellen isPrime(UInt64) arbeite ich seit einer Woche, alse genau die Gegenrichtung von Dir.



  • Mich überrascht bei diesem Thema (noch immer unbewiesene Vermutung des Herren Goldbach) die Vielfalt der Möglichkeiten bei der Herangehensweise und dass man mit einem leistungsstarken PC in wirklich beeindruckende Zahlenbereiche vordringen kann. Ich bin gespannt auf Volkards Neuerungen.



  • Der Satz ist zwar unbewiesen, aber man kann sehr einfach beweisen, dass es verdammt unwahrscheinlich wäre, wenn er nicht wahr wäre 😃



  • Erhard Henkes schrieb:

    [...] und dass man mit einem leistungsstarken PC in wirklich beeindruckende Zahlenbereiche vordringen kann. Ich bin gespannt auf Volkards Neuerungen.

    Wobei man bei der ganzen Sache nicht vergessen sollte, das man zwar im Bereich 10^120 einzelne Zahlen in ein paar Millisekunden testen kann aber nicht alle Zahlen bis dahin. Selbst wenn du die gesamte Rechenkapazität der Erde nutzen würdest, kriegst du bis zum Weltuntergang nicht alle Zahlen getestet. Nur so als Vergleich: Das Universum ist ca. 4*10^17 Sekunden alt.



  • sebi707: Ja, das ist richtig. Egal, wie toll ein Algorithmus oder Rechner ist, einfach einige hundert Zehnerpotenzen weiter, und schon ist Schluss. Alle Zahlen bis 2^64 zu checken erfordert umfangreiches verteiltes Rechnen. Bei Collatz habe ich da mal vor ein paar Jahren Teil genommen. Da gabs 20000er Pakete. http://www.ericr.nl/wondrous/

    Hier ist ein interessanter Artikel zum Berechnen der Goldbach Vermutung (schon 15 Jahre alt): http://www.ams.org/journals/mcom/2001-70-236/S0025-5718-00-01290-4/S0025-5718-00-01290-4.pdf



  • Ich bin mit dem aktuellen Programm wieder etwas weiter ins "Zahlen-Universum" vorgestoßen und habe den Bereich von 10^360 bis 10^360 + 100 durchgerechnet:

    `total time: 9302703 ms

    total time at mulmod: 8612644 ms (92 😵

    High Numbers to be tested: 28593

    Number of successful Division Prechecks: 8801 (30.8 😵

    Number of MR tests: 758 (2.65 😵

    Number of saved MR tests by high prime lookup: 19034 (66.6 %)`

    Alle Ergebnisse: http://pastebin.com/Ftgq6eMt

    Die höchste Spanne nach unten zur nächsten Primzahl ist dort 34259.



  • Bengo schrieb:

    Der Satz ist zwar unbewiesen, aber man kann sehr einfach beweisen, dass es verdammt unwahrscheinlich wäre, wenn er nicht wahr wäre 😃

    Kann man einem, der in Mathe nicht sonderlich bewandert ist, halbwegs einfach erklären wie dieser "Beweis der Unwahrscheinlichkeit" aussieht? Bzw. was da jetzt überhaupt mit unwahrscheinlich gemeint ist?
    Würde mich nämlich interessieren.


  • Mod

    hustbaer schrieb:

    Bengo schrieb:

    Der Satz ist zwar unbewiesen, aber man kann sehr einfach beweisen, dass es verdammt unwahrscheinlich wäre, wenn er nicht wahr wäre 😃

    Kann man einem, der in Mathe nicht sonderlich bewandert ist, halbwegs einfach erklären wie dieser "Beweis der Unwahrscheinlichkeit" aussieht? Bzw. was da jetzt überhaupt mit unwahrscheinlich gemeint ist?
    Würde mich nämlich interessieren.

    Das Argument beruht auf der Primzahlhäufigkeit. Wenn man zufällig eine Zahl X wählt, ist diese im Durchschnitt zu 1/log(X) eine Primzahl.

    Wenn wir nun eine große Zahl N nehmen, die sich als Summe aus M und N-M darstellen lässt, dann ist die Wahrscheinlichkeit, dass M und N-M beides Primzahlen sind 1/(log(M)*log(N-M)). Wenn wir nun alle Möglichkeiten für M betrachten (also M zwischen 3 und N/2), dann ist die durchschnittliche Anzahl von Primzahlpaaren M und N-M die Summe über alle M von obiger Wahrscheinlichkeit. Das ist für sehr große N ungefähr N/(2*log(N)^2). Die Anzahl von im Schnitt zu erwartenden Primzahlpaaren, um N darzustellen, wird also immer größer, wenn N größer wird. Da wir wissen, dass für kleine N stets solche Paare existieren, ist nicht zu erwarten, dass wir bei größeren N irgendwann kein solches Paar finden werden. Im Gegenteil sollte die Anzahl der Darstellungen als Summe eines Primzahlpaares sogar immer größer werden.

    Das wird natürlich einem strengen mathematischen Beweis nicht annähernd gerecht, aber es vermittelt eine starke Vorstellung davon, wieso die Vermutung wahrscheinlich richtig ist. Die Anzahl möglicher Zahlenpaare, die in der Summe eine bestimmte Zahl ergeben, nimmt schneller zu als die Wahrscheinlichkeit, dass beide diese Zahlen Primzahlen sind, abnimmt.


Anmelden zum Antworten