Python rechnet uneffizient?

Wenn du dir nicht sicher bist, in welchem der anderen Foren du die Frage stellen sollst, dann bist du hier im Forum für allgemeine Fragen sicher richtig.
DerTürke

Hallo BlackJack,

ich sah mich nun gezwungen das Sieb des Eratosthenes in einem von mir veränderten Algorithmus in c++ zu erstellen, hierbei wollte ich auch kurz demonstrieren was ich mit dynamischer Liste ( in der STL ist dies ein TemplateClass welche vector<T> heisst )
ich habe Testhalber mal versucht Primzahlen bis 200.000 zu erstellen nach ein paar Minuten habe ich dann einfach abgebrochen, weil es einfach zu lange dauert. Ich habe wie du siehst von vornherein die geraden Zahlen gemieden, und dann habe ich alle 2
verschachtelete Schleifen welche die Vielfachen ermitteln und diese aus der List ( vecotr<int> ) entfernen somit habe ich mit jedem äußeren durchlauf immer weniger elemente in der List (vector<int>). So das zu letzt die Liste (vector<int>)
nur noch die Primzahlen enthält.
Und dennoch muss ich sagen auch wenn der algorithmus nun dem Sieb entspricht und auch einige Modifikationen daran gemacht wurden, muss ich sagen das mir die Assembler variante einfach besser gefällt da diese extrem schnell ist.

Code: Alles auswählen

#include<iostream>
#include<vector>

using namespace std;

int main()
{
	
	vector<int> m_Sieb = vector<int>(5000);

	// Explizites initialisieren des ersten Elements mit der ersten Primzahl
	m_Sieb[0] = 2;
	
	// Initialisiere vector(array) mit alle ungeraden Zahlen von 3 bis xxx
	for (int i = 1; i < 5000; m_Sieb[i] = i++ * 2 + 1){}

	vector<int>::iterator m_CurrentIndex = m_Sieb.begin();
	
	m_CurrentIndex++; ///erfordelich da die Vielfachen von 2 nicht enthalten sind 

	while (m_CurrentIndex < m_Sieb.end())
	{
		if (*m_CurrentIndex != 0)
		{
			/// Dann bitte finden und eliminieren
			vector<int>::iterator m_VielfacheIndex = m_CurrentIndex;
			while (++m_VielfacheIndex < m_Sieb.end())
			{
				
				if (*m_VielfacheIndex%*m_CurrentIndex == 0)
				{
					//*m_VielfacheIndex = 0;
					vector<int>::iterator tmpIterator = m_VielfacheIndex - 1;
					m_Sieb.erase(m_VielfacheIndex);
					m_VielfacheIndex = tmpIterator;
				}
			}
		}
		m_CurrentIndex++;
	}

	// *********************************** an diesem Punkt
	m_CurrentIndex = m_Sieb.begin();
	while (++m_CurrentIndex < m_Sieb.end())
	{
		cout << "Primzahl : " << *m_CurrentIndex << endl;
	}

	int m_Value;
	cin >> m_Value; /// Total sinnlos einfach nur um die Konsole nicht sofort zu schliessen;
	return 0;
}

ok die Zeile
if (*m_CurrentIndex != 0) ist unnötig da alle Felder bei dieser Zeile auf jeden Fall einen Wert haben, dies war ein Überbleibsel aus der Markierung mit 0 wo ich noch nicht die Liste (vector<int>) gelöscht hatte, aber diese unnötige Zeile sollte nicht allzuviel zusätzl. Taktzyklen erzeugen.


das sieb hat // *********************************** an diesem Punkt
nur noch eine Größe von 1229 elemente
dauer bis ca. 5 sekunden

Ich habe auch einfach mal etwas aus dem compilat (des c++ compilers ) rausgesucht.

Es geht schlichtweg nur um das initialisieren der Liste (vector<int>)

// Explizites initialisieren des ersten Elements mit der ersten Primzahl
m_Sieb[0] = 2;
00A16F0E push 0
00A16F10 lea ecx,[m_Sieb]
00A16F13 call std::vector<int,std::allocator<int> >::operator[] (0A115CDh)
00A16F18 mov dword ptr [eax],2

// Initialisiere vector(array) mit alle ungeraden Zahlen von 3 bis xxx
for (int i = 1; i < 5000; m_Sieb = i++ * 2 + 1){}
00A16F1E mov dword ptr [ebp-30h],1
00A16F25 jmp main+95h (0A16F45h)
00A16F27 mov eax,dword ptr [ebp-30h]
00A16F2A lea esi,[eax+eax+1]
00A16F2E mov ecx,dword ptr [ebp-30h]
00A16F31 push ecx
00A16F32 lea ecx,[m_Sieb]
00A16F35 call std::vector<int,std::allocator<int> >::operator[] (0A115CDh)
00A16F3A mov dword ptr [eax],esi
00A16F3C mov edx,dword ptr [ebp-30h]
00A16F3F add edx,1
00A16F42 mov dword ptr [ebp-30h],edx
00A16F45 cmp dword ptr [ebp-30h],1388h
00A16F4C jge main+0A0h (0A16F50h)

soviel zur optimalen umsetzung des compilers in Maschinensprache.

Ich glaube da musst du mir doch hoffentlich Recht geben, wenn ich sage das mit jeder Erkennung eines Vielfachen auch die Gesamtmenge der Liste(vector<int>) verringert wird, das das vom Sinn ( der heran gehensweise ) optimal sein müsste, darum ging es dir doch das
der Algorithmus im Vordergrund stand, mir ging es lediglich um die Geschwindigkeit.

Achso weil du ja nochmal nachgehakt hattest

Das Argument mit dem gegenläufigen Verhalten verstehe ich nicht? Was verhält sich da gegenläufig? `r10` wird in jedem Durchlauf erhöht und `rax` bleibt doch wohl hoffentlich immer gleich,
nämlich der Kandidat der auf prim getestet wird. Und für den Fall, dass er das *ist*, wird `r10` bis zur Primzahl und nur bis zur Wurzel hochgezählt, was unnötig ist.
Das Hochzählen um 1 ist auch nicht besonders glücklich, weil man die Anzahl der Divisionen für einen Kandidaten locker mal halbieren kann, wenn man erst auf Teilbarkeit durch 2 testet, und dann nur noch ungerade Testdivisionen rechnet.


hier nochmal die Sequenz

@Weiter:
xor rdx, rdx
MOVE_PARAM rax, 1
div r10
cmp rdx, 0
jz @NoPrime
inc r10
cmp r10, rax
jb @Weiter

kurze Erklärung:

divident ist im Registerpaar rdx:rax
divisor ist im Register r10
nach efolgter divison (div r10)
ist Quotient im Reigster RAX ( den ich dann als Prüfbedingung im cmp Befehl verwenden werde )
und Remainder (Rest) im Register RDX
jetzt prüfe ich ob es keinen Rest gegeben hat und falls die Bedingung zutrifft springe ich zu NoPrime
ansonsten erhöhe ich r10
und prüfe ob divisor < als Quotient ( denn die Vielfachen des Quotienten brauchen ja nicht mehr geprüft zu werden daher zähle ich nie bis n komplett hoch )
wenn das so ist dann springe auf Weiter, falls nein ist es eine Primzahl

Beispiel:
die Zahle 101 = Primzahl

zuvor
r10 = 2

1. Schleifendurchlauf
rdx = 0
rax = 101
-> div r10
rdx = 1
rax = 50
-> cmp rdx, 0
-> nein also nicht springen
-> inc r10
r10 = 3
-> cmp r10, rax (3 < 50)
-> jb @Weiter ( sprung wird ausgeführt )

2. Schleifen durchlauf
rdx = 0
rax = 101
-> div r10
rdx = 2
rax = 33
-> cmp rdx, 0
-> nein also nicht springen
-> inc r10
r10 = 4
-> cmp r10, rax (4 < 33)
-> jb @Weiter ( sprung wird ausgeführt )

3. Schleifen durchlauf
rdx = 0
rax = 101
-> div r10
rdx = 1
rax = 25
-> cmp rdx, 0
-> nein also nicht springen
-> inc r10
r10 = 5
-> cmp r10, rax (5 < 25)
-> jb @Weiter ( sprung wird ausgeführt )

4. Schleifen durchlauf
rdx = 0
rax = 101
-> div r10
rdx = 1
rax = 20
-> cmp rdx, 0
-> nein also nicht springen
-> inc r10
r10 = 6
-> cmp r10, rax (6 < 20)
-> jb @Weiter ( sprung wird ausgeführt )

5. Schleifen durchlauf
rdx = 0
rax = 101
-> div r10
rdx = 5
rax = 16
-> cmp rdx, 0
-> nein also nicht springen
-> inc r10
r10 = 7
-> cmp r10, rax (7 < 16)
-> jb @Weiter ( sprung wird ausgeführt )

6. Schleifen durchlauf
rdx = 0
rax = 101
-> div r10
rdx = 3
rax = 14
-> cmp rdx, 0
-> nein also nicht springen
-> inc r10
r10 = 8
-> cmp r10, rax (8 < 14)
-> jb @Weiter ( sprung wird ausgeführt )

7. Schleifen durchlauf
rdx = 0
rax = 101
-> div r10
rdx = 5
rax = 12
-> cmp rdx, 0
-> nein also nicht springen
-> inc r10
r10 = 9
-> cmp r10, rax (9 < 12)
-> jb @Weiter ( sprung wird ausgeführt )

8. Schleifen durchlauf
rdx = 0
rax = 101
-> div r10
rdx = 2
rax = 11
-> cmp rdx, 0
-> nein also nicht springen
-> inc r10
r10 = 10
-> cmp r10, rax (10 < 11)
-> jb @Weiter ( sprung wird ausgeführt )

9. Schleifen durchlauf
rdx = 0
rax = 101
-> div r10
rdx = 1
rax = 10
-> cmp rdx, 0
-> nein also nicht springen
-> inc r10
r10 = 11
-> cmp r10, rax (11 < 10
-> jb @Weiter ( sprung wird NICHT ausgeführt )

-> Ende der Prüfungen die Zahl 101 ist eine Primzahl ( dazu habe ich insgesamt 9 mal die Schleife durchlaufen und wenn man dann analysiert kommt folgende Formel heraus : (Wurzel aus N)-1 divisionen ( im max. Fall )

Ich hoffe jetzt ist es ein bisschen transparenter

Viele Grüße
DT
Sirius3
User
Beiträge: 17754
Registriert: Sonntag 21. Oktober 2012, 17:20

@DerTürke: Du hast es wirklich geschafft, Deinen langsamen Primzahlsucher noch um eine lineare Suche in einer Liste und n^2 Kopieroperationen zu kombinieren. Setze doch mal das Sieb in seiner ursprünglichen Form um und lass Deine »Optimierungen« weg.
Noch einen Nachtrag zu meinen Messungen:
Primzahlen zwischen 5000000000000 und 5000000000100:
- per is_prime: 20ms
- per Sieb: 14ms

Egal wie ich es versuch, das Sieb ist immer schneller.
Zuletzt geändert von Sirius3 am Sonntag 15. Mai 2016, 00:23, insgesamt 2-mal geändert.
DerTürke

@Sirius

Da kann ich jetzt wirklich nur noch passen
Schade hat viel Spass gemacht

Viele Grüße
DT
Sirius3
User
Beiträge: 17754
Registriert: Sonntag 21. Oktober 2012, 17:20

PS: man sollte dem Compiler sagen, dass er optimieren soll:
[codebox=asm file=Unbenannt.asm]
## m_Sieb[0] = 2;
movl $2, (%r14)
## for (int i = 1; i < 5000; m_Sieb = i++ * 2 + 1){}
leaq 20000(%r14), %r13
movl $4, %ecx
.align 4, 0x90
LBB0_2:
leal -1(%rcx), %eax
movl %eax, (%r14,%rcx,2)
addq $2, %rcx
cmpl $10002, %ecx
jne LBB0_2
[/code]
Da bleibt nicht viel Luft, zum händischen Rumschrauben.
DerTürke

@Sirius:

ok was ist dass jetzt ?
Sirius3
User
Beiträge: 17754
Registriert: Sonntag 21. Oktober 2012, 17:20

@DerTürke: die Assemblerausgabe eines handelsüblichen C++-Compilers.
DerTürke

@Sirius:
und was bedeuten diese movl, leaq, movl .align auf was für einem Prozessor werden diese Befehle ausgeführt ( ARM ? SPARC ? ) also es ist definitiv nicht Intel und AMD und auch der 6502/6510 kennt diese nicht es sei denn die Mnemonic's werden vom compile später in verständliche Opcode's umgesetzt. Wie kann ich den dem Compiler sagen dass er optimieren soll wo mach ich dass ?

Nenn mir doch mal bitte deinen C++ Compiler der dir diesen Code erstellt hat
Danke und Gruß
DerTürke

@Sirius:
@DerTürke: Du hast es wirklich geschafft, Deinen langsamen Primzahlsucher noch um eine lineare Suche in einer Liste und n^2 Kopieroperationen zu kombinieren. Setze doch mal das Sieb in seiner ursprünglichen Form um und lass Deine »Optimierungen« weg.
Noch einen Nachtrag zu meinen Messungen:
Primzahlen zwischen 5000000000000 und 5000000000100:
- per is_prime: 20ms
- per Sieb: 14ms

Egal wie ich es versuch, das Sieb ist immer schneller.
Ich glaube langsam du hasst mal absolut keine Ahnung Kollege und davon so glaube ich ziemlich reichlich und viel
DerTürke

@BlackJack
bestimmt:
Ich habe ja selber nicht dran geglaubt und mal meinen Algo. in C++ verfasst, ich muss sagen rasend schnell, vielleicht nur 0,5 sek langsamer als in Assembler:
Schaus dir an wenn du möchtest compiliere:

Code: Alles auswählen

	for (int i = 2; i < 2000000; i++)
	{
		bool bIsPrime = IsPrime(i);
		/*
		if (bIsPrime)
		{
			cout << "Primzahl : " << i << endl;
		}
		*/
	}
	cin >> m_Value; /// Total sinnlos einfach nur um die Konsole nicht sofort zu schliessen;
	
bool IsPrime(int pValue)
{
	bool retValue = true;
	int m_Divisor = 2;
	int m_Divident = pValue;
	div_t div_Result;
	div_Result.quot = m_Divident;
	while (m_Divisor < div_Result.quot)
	{
		div_Result = div(m_Divident, m_Divisor);
		if (div_Result.rem == 0)
		{
			retValue = false;
			break;
		}
		m_Divisor++;
	}
	return retValue;
}
Das Sieb hatte ich ja bereits vorher implementiert da habe ich die PZ von 2-100.000 errechnet in ca. 50 Sek.
Viele Grüße
DT
BlackJack

@DerTürke: Dir gefällt die Assemblervariante besser weil Du in C++ absichtlich einen so schlechten Algorithmus implementieren kannst, der *noch* ineffizienter ist? Und dann noch nicht mal dem Compiler sagen das er optimieren soll. :roll:

Natürlich ist das vom Sinn her nicht optimal. Klar wird die Gesamtmenge kleiner, aber zu einem viel zu hohen Preis von linearem suchen und unnötigem verschieben von Unmengen von Bytes.

Der Assemblerquelltext den Sirius3 für die Initialisierung gezeigt hat, ist für die x64-Architektur, also für die Prozessoren, für die Du auch geschrieben hast. Du scheinst die AT&T-Syntax anscheinend nur nicht zu kennen. In Intel-Syntax sieht das so aus, vielleicht kommt Dir das bekannter vor, es ist aber letztlich der selbe Code:

[codebox=asm file=Unbenannt.asm] mov DWORD PTR [R14], 2
lea R12, QWORD PTR [R14 + 20000]
mov EAX, 4
.align 16, 0x90
.LBB0_2:
lea ECX, DWORD PTR [RAX - 1]
mov DWORD PTR [R14 + 2*RAX], ECX
add RAX, 2
cmp EAX, 10002
jne .LBB0_2[/code]

Das produziert der Clang-Compiler.

Wenn Dein Algorithmus spürbar langsamer ist als Deine Assemblervariante, dann hast Du sehr wahrscheinlich auch dort vergessen dem Compiler zu sagen, dass Du optimierten Code möchtest.

Die Assemblervariante (ohne die Tests auf die ersten Primzahlen) könnte in C so aussehen:
[codebox=c file=Unbenannt.c]_Bool is_prime(uint_fast32_t candidate)
{
uint_fast32_t i;

for (i = 2; i < candidate; ++i) {
if (candidate % i == 0) return false;
if (candidate / i < i) break;
}
return candidate > 1;
}[/code]
Daraus macht der GCC das hier wenn man schnellen Code möchte:
[codebox=asm file=Unbenannt.asm]is_prime:
.LFB26:
.cfi_startproc
cmp rdi, 2
jbe .L2
xor eax, eax
test dil, 1
je .L3
cmp rdi, 3
mov ecx, 2
jne .L5
jmp .L2
.p2align 4,,10
.p2align 3
.L6:
xor edx, edx
mov rax, rdi
div rcx
test rdx, rdx
je .L8
cmp rax, rcx
jb .L2
.L5:
add rcx, 1
cmp rcx, rdi
jne .L6
.L2:
cmp rdi, 1
seta al
ret
.p2align 4,,10
.p2align 3
.L8:
xor eax, eax
.L3:
rep
ret[/code]

Und das hier wenn man möglichst kurzen Code möchte:
[codebox=asm file=Unbenannt.asm]is_prime:
.LFB12:
.cfi_startproc
mov ecx, 2
jmp .L2
.L5:
xor edx, edx
mov rax, rdi
div rcx
test rdx, rdx
je .L6
cmp rax, rcx
jb .L4
inc rcx
.L2:
cmp rcx, rdi
jb .L5
.L4:
cmp rdi, 1
seta al
ret
.L6:
xor eax, eax
ret[/code]

Die Schleifen sehen Deinem handgeschriebenen Assembler doch verdammt ähnlich. Ohne das man Assembler schreiben muss. Der C-Code ist schneller geschrieben, verständlicher, und wenn man mal einen anderen Prozessor vor sich hat, kann man ihn auch dafür übersetzen. Zum Beispiel für ARMv7:
[codebox=arm file=Unbenannt.txt]is_prime:
@ args = 0, pretend = 0, frame = 0
@ frame_needed = 0, uses_anonymous_args = 0
cmp r1, #0
it eq
cmpeq r0, #3
push {r3, r4, r5, lr}
bcc .L6
and r4, r0, #1
movs r5, #0
orrs r4, r4, r5
beq .L4
cmp r1, #0
it eq
cmpeq r0, #6
bcs .L6
cmp r1, #0
it eq
cmpeq r0, #3
beq .L6
movs r2, #3
movs r3, #0
bl __aeabi_uldivmod
orrs r3, r2, r3
ite eq
moveq r0, #1
movne r0, #0
pop {r3, r4, r5, pc}
.L6:
movs r0, #0
pop {r3, r4, r5, pc}
.L4:
movs r0, #1
pop {r3, r4, r5, pc}[/code]
(Hier war ich erst etwas stutzig, weil da kein Label für die Schleife ist, bis mir aufgefallen ist, dass das hier über das Linkregister und `push`/`pop` gemacht wird. :-))

So jetzt musst Du nur noch ein effizientes Sieb implementieren und nachmessen was Sirius3 schon gemessen hat…
Benutzeravatar
snafu
User
Beiträge: 6741
Registriert: Donnerstag 21. Februar 2008, 17:31
Wohnort: Gelsenkirchen

Stutzig macht mich ja, dass der GCC beim kompakten Code `INC` zum Erhöhen der Register-Wertes in der Schleife verwendet, dies aber in der performanten Variante über ein explizites `ADD` erledigt. Die GCC-Entwickler werden sicherlich ihre Gründe haben, warum sie das so implementiert haben, aber eigenartig finde ich es trotzdem.
BlackJack

@snafu: Der Grund ist ganz einfach: ``INC reg`` ist kompakter und ``ADD reg, 1`` ist auf den meisten x64-Prozessoren mindestens gleich schnell. Eine Ebene tiefer werden beide ziemlich wahrscheinlich zum gleichen, oder zumindest sehr ähnlichen Mikroinstruktionen. Ganz äquivalent sind die beiden von der Funktion her nicht, weil INC nicht alle Flags beeinflusst die von ADD beeinflusst werden.

Für Intel-Prozessoren generiert mir der GCC für den gezeigten Quelltext immer die ADD-Variante. Für AMD-Prozessoren dagegen INC. Das ist ja auch so ein Punkt den man bei handoptimiertem Assembler beachten muss: Verschiedene Prozessoren behandeln die gleichen Opcodes unterschiedlich. Und dann kommen noch Sachen wie Pipelining dazu, die dafür sorgen können, dass neben der Wahl des Opcodes auch die Reihenfolge der Abarbeitung eine Rolle spielen. Wenn man also Assembler von Hand schreibt, weil man glaubt man ist grundsätzlich besser als der Compiler, dann sollte man sich auch wirklich sicher sein, dass der handgeschriebene Code auch dann noch besser ist, wenn das C- oder C++-Gegenstück auf einem beliebigen Prozessor (der gleichen Architektur) speziell für *den* Prozessor optimiert wurde.
DerTürke

Hallo Jungs,

ich kanns einfach nicht sein lassen aber mein Algo ist nun mal schneller.
Diesmal habe ich wie von euch beanstandet das Sieb in seiner ursprünglichen Form in reinem C geschrieben und dennoch ist meine IsPrime Funktion schneller, tut mir echt leid aber das ist nun mal Tatsache ob ihr das hinnehmen wollt oder nicht

Code: Alles auswählen

#include<iostream>
#include<vector>
#include<stdlib.h>

using namespace std;

bool IsPrime(int pValue);

int main()
{	
	int n_Sieb[200000];
	int n_Primes[33860];
	int n_PrimesIndex = 1;	

	n_Sieb[0] = 2;
	n_Primes[0] = 2;
	for (int i = 1; i < 200000; i++)
	{
		n_Sieb[i] = i * 2 + 1;
	}

	for (int i = 1; i < 200000; i++)
	{
		if (n_Sieb[i] > 0)
		{
			n_Primes[n_PrimesIndex++] = n_Sieb[i];
			for (int j = i + 1; j < 200000; j++)
			{
				if (n_Sieb[j] % n_Sieb[i] == 0)
				{
					n_Sieb[j] = 0;
				}
			}
		}
	}
	
	/// Jetzt kommt mein Algorithmus
	n_PrimesIndex = 0;
	for (int i = 2; i < 400000; i++)
	{
		bool bIsPrime = IsPrime(i);				
		if (bIsPrime)
		{
			n_Primes[n_PrimesIndex++] = i;
		}		
	}
	return 0;
}


bool IsPrime(int pValue)
{
	bool retValue = true;
	int m_Divisor = 2;
	int m_Divident = pValue;
	div_t div_Result;
	div_Result.quot = m_Divident;
	while (m_Divisor < div_Result.quot)
	{
		div_Result = div(m_Divident, m_Divisor);
		if (div_Result.rem == 0)
		{
			retValue = false;
			break;
		}
		m_Divisor++;
	}
	return retValue;
}
und auch mit Code Optimierung sieht das compilat zwar nun viel besser als vorher aus aber dennoch ist meins kompakter:

Code: Alles auswählen

using namespace std;

bool IsPrime(int pValue);

int main()
{	
00007FF6D6001000  mov         qword ptr [rsp+8],rbx  
00007FF6D6001005  mov         qword ptr [rsp+10h],rsi  
00007FF6D600100A  push        rdi  
00007FF6D600100B  mov         eax,0E4630h  
00007FF6D6001010  call        __chkstk (07FF6D6001EA0h)  
00007FF6D6001015  sub         rsp,rax  
00007FF6D6001018  movdqa      xmm2,xmmword ptr [__xmm@00000003000000020000000100000000 (07FF6D6002220h)]  
	int n_Sieb[200000];
	int n_Primes[33860];
	int n_PrimesIndex = 1;	

	n_Sieb[0] = 2;
00007FF6D6001020  lea         rcx,[rsp+21134h]  
00007FF6D6001028  movdqa      xmm3,xmmword ptr [__xmm@00000001000000010000000100000001 (07FF6D6002210h)]  
	n_Primes[0] = 2;
	for (int i = 1; i < 200000; i++)
00007FF6D6001030  mov         edx,1  
00007FF6D6001035  mov         dword ptr [n_Sieb],2  
00007FF6D6001040  lea         eax,[rdx+4]  
00007FF6D6001043  movd        xmm0,edx  
00007FF6D6001047  pshufd      xmm0,xmm0,0  
00007FF6D600104C  lea         rcx,[rcx+20h]  
00007FF6D6001050  paddd       xmm0,xmm2  
00007FF6D6001054  movd        xmm1,eax  
00007FF6D6001058  pshufd      xmm1,xmm1,0  
00007FF6D600105D  paddd       xmm0,xmm0  
00007FF6D6001061  paddd       xmm1,xmm2  
	{
		n_Sieb[i] = i * 2 + 1;
00007FF6D6001065  paddd       xmm0,xmm3  
00007FF6D6001069  movdqu      xmmword ptr [rcx-20h],xmm0  
00007FF6D600106E  add         edx,8  
00007FF6D6001071  paddd       xmm1,xmm1  
00007FF6D6001075  paddd       xmm1,xmm3  
	{
		n_Sieb[i] = i * 2 + 1;
00007FF6D6001079  movdqu      xmmword ptr [rcx-10h],xmm1  
00007FF6D600107E  cmp         edx,30D39h  
00007FF6D6001084  jl          main+40h (07FF6D6001040h)  
	n_Primes[0] = 2;
	for (int i = 1; i < 200000; i++)
00007FF6D6001086  cmp         edx,30D40h  
00007FF6D600108C  jge         main+0C0h (07FF6D60010C0h)  
00007FF6D600108E  movsxd      rax,edx  
00007FF6D6001091  lea         rcx,[n_Sieb]  
00007FF6D6001099  lea         rcx,[rcx+rax*4]  
00007FF6D600109D  lea         eax,[rdx*2+1]  
00007FF6D60010A4  nop         dword ptr [rax]  
00007FF6D60010A8  nop         dword ptr [rax+rax]  
	{
		n_Sieb[i] = i * 2 + 1;
00007FF6D60010B0  mov         dword ptr [rcx],eax  
00007FF6D60010B2  lea         rcx,[rcx+4]  
00007FF6D60010B6  add         eax,2  
00007FF6D60010B9  cmp         eax,61A81h  
00007FF6D60010BE  jl          main+0B0h (07FF6D60010B0h)  
	}

	for (int i = 1; i < 200000; i++)
00007FF6D60010C0  mov         r9d,2  
00007FF6D60010C6  lea         r8,[rsp+21134h]  
00007FF6D60010CE  lea         r10,[rsp+24h]  
00007FF6D60010D3  xor         ebx,ebx  
00007FF6D60010D5  mov         r11d,30D3Fh  
00007FF6D60010DB  nop         dword ptr [rax+rax]  
	{
		if (n_Sieb[i] > 0)
00007FF6D60010E0  mov         eax,dword ptr [r8]  
	{
		if (n_Sieb[i] > 0)
00007FF6D60010E3  test        eax,eax  
00007FF6D60010E5  jle         main+122h (07FF6D6001122h)  
		{
			n_Primes[n_PrimesIndex++] = n_Sieb[i];
00007FF6D60010E7  add         r10,4  
			for (int j = i + 1; j < 200000; j++)
00007FF6D60010EB  mov         rcx,r9  
00007FF6D60010EE  cmp         r9,30D40h  
00007FF6D60010F5  jge         main+122h (07FF6D6001122h)  
00007FF6D60010F7  nop         word ptr [rax+rax]  
			{
				if (n_Sieb[j] % n_Sieb[i] == 0)
00007FF6D6001100  mov         eax,dword ptr n_Sieb[rcx*4]  
00007FF6D6001107  cdq  
00007FF6D6001108  idiv        eax,dword ptr [r8]  
00007FF6D600110B  test        edx,edx  
00007FF6D600110D  jne         main+116h (07FF6D6001116h)  
				{
					n_Sieb[j] = 0;
00007FF6D600110F  mov         dword ptr n_Sieb[rcx*4],ebx  
			for (int j = i + 1; j < 200000; j++)
00007FF6D6001116  inc         rcx  
00007FF6D6001119  cmp         rcx,30D40h  
00007FF6D6001120  jl          main+100h (07FF6D6001100h)  
	}

	for (int i = 1; i < 200000; i++)
00007FF6D6001122  inc         r9  
00007FF6D6001125  add         r8,4  
00007FF6D6001129  sub         r11,1  
00007FF6D600112D  jne         main+0E0h (07FF6D60010E0h)  
				}
			}
		}
	}
Sirius3
User
Beiträge: 17754
Registriert: Sonntag 21. Oktober 2012, 17:20

@DerTürke: Du hast immer noch kein Sieb implementiert. Der Trick ist ja, dass man gar keine Divisionen braucht. Die innere Schleife sollte ungefähr so aussehen:
[codebox=c file=Unbenannt.c] n_Primes[n_PrimesIndex++] = n_Sieb;
for (int j = i + i; j < 200000; j+=i)
{
n_Sieb[j] = 0;
}
[/code]

Damit wäre bewiesen, wer absolut keine Ahnung hat :evil:.
Sirius3
User
Beiträge: 17754
Registriert: Sonntag 21. Oktober 2012, 17:20

Eine weitere simple Optimierung betrifft das Auslassen von durch 3 teilbaren Zahlen. Dadurch erhöht sich die Geschwindigkeit folglich um 1/3 :wink: .
Bei den Primzahlen bis 5000 sieht das Rennen nun so aus:
- is_prime: 107.791 Takte
- Sieb: 12.869 Takte

@DerTürke: Absolut gesehen holt damit natürlich is_prime auf :twisted: .
BlackJack

Okay, dann nehmen wir doch mal ein ganz einfaches, wirklich „straight forward“ implementiertes Sieb in Python, ohne irgendwas speziell zu optimieren:

Code: Alles auswählen

#!/usr/bin/env python
# coding: utf8
from __future__ import absolute_import, division, print_function
from math import sqrt


def sieve(limit):
    is_prime = [True] * limit
    is_prime[0] = is_prime[1] = False
    for i in xrange(2, int(sqrt(limit)) + 1):
        if is_prime[i]:
            for j in xrange(i * i, limit, i):
                is_prime[j] = False
    return (i for i, t in enumerate(is_prime) if t)


def main():
    print(sum(sieve(2000000)))


if __name__ == '__main__':
    main()
Und lassen das gegen das folgende C-Programm antreten, welches den Algorithmus von DerTürke implementiert:
[codebox=c file=Unbenannt.c]#include <inttypes.h>
#include <stdbool.h>
#include <stdio.h>


static _Bool is_prime(uint_fast32_t candidate)
{
uint_fast32_t i;

for (i = 2; i < candidate; ++i) {
if (candidate % i == 0) return false;
if (candidate / i < i) break;
}
return candidate > 1;
}


int main(void)
{
uint_fast64_t sum = 0;
uint_fast32_t i;

for (i = 1; i < 2000000; ++i) {
if (is_prime(i)) {
sum += i;
}
}

printf("%"PRIuFAST64"\n", sum);
return 0;
}[/code]

Kompiliert mit GCC und Optimierungen wird die Testfunktion direkt in die Hauptfunktion eingefügt, so dass dort der Funktionsaufruf eingespart wird. Ich habe den generierten Quelltext mal mit sprechenderen Label-Namen versehen und ein wenig kommentiert:
[codebox=asm file=Unbenannt.asm] .intel_syntax noprefix
.section .rodata.str1.1,"aMS",@progbits,1
.format_string:
.string "%lu\n"

.section .text.startup,"ax",@progbits
.p2align 4
.globl main
main:
; candidate = 1
mov esi, 1
; sum = 0
xor edi, edi
; Align loop start to 16 byte boundary if max 10 bytes away
; or to 8 byte boundary
.p2align 4,,10
.p2align 3
.outer_loop:
; if candidate <= 2
cmp rsi, 2
jbe .L8

; if candidate is even but not <=2 it is not prime
test sil, 1
je .is_not_prime

; if candidate is three then it is prime,
; else i = 2 and run inner loop
cmp rsi, 3
mov ecx, 2
jne .inner_loop_inc
jmp .is_prime
.p2align 4,,10
.p2align 3
.inner_loop:
; div/mod candidate and i
xor edx, edx
mov rax, rsi
div rcx
; if remainder then not prime
test rdx, rdx
je .is_not_prime
; if result of division < i then candidate is prime
cmp rax, rcx
jb .is_prime
.inner_loop_inc:
; increase i and if not reached candidate then test i
add rcx, 1
cmp rcx, rsi
jne .inner_loop
.is_prime:
; sum += candidate
add rdi, rsi
.is_not_prime:
; increase candidate and test if upper limit is hit
add rsi, 1
cmp rsi, 2000000
jne .outer_loop

; print sum
mov rdx, rdi
mov esi, OFFSET FLAT:.format_string
mov edi, 1
xor eax, eax
call __printf_chk

; return 0
xor eax, eax
ret
.L8:
; if candidate is 1 then candidate = 2
; otherwise candidate stays 2
cmp rsi, 1
mov eax, 2
cmove rsi, rax
jmp .is_prime[/code]

Wie man sieht entspricht das ziemlich dem was DerTürke an Assembler gezeigt hat. Etwas schräg ist der erste Test und der Code bei `.L8`. Das hätte man als Mensch wahrscheinlich anders geschrieben. Letztlich fällt das bei der Laufzeit aber kaum ins Gewicht, denn es betrifft nur einen von zwei Millionen Schleifendurchläufen.

Interessant ist, dass der Compiler selbst gemerkt hat, das ungerade Kandidaten ausser der 2 nicht zu einer Addition auf die Summe führen und das deswegen mit einem billigen Bit-Test gesondert behandelt, um die teurere Division für jede zweite Zahl zu sparen.

Und das Ergebnis: Python ist fast doppelt so schnell.

Um mal einen Fusballtrainer zu zitieren: Schön spielen ist nichts, wenn das Resultat nicht stimmt. Man kann sich viel Mühe geben handgeklöppelten Assembler zu schreiben, und dafür sicher ein Fleissbienchen verdienen, aber wenn man dann *deutlich* von einer dynamisch interpetierten Programmiersprache überholt wird, wo sogar noch der Compilerlauf von Quelltext nach Bytecode für die VM mit in der Laufzeit enthalten ist, dann war das wohl nix. Aber dicke Töne spucken, andere hätten keine Ahnung. Süss. :twisted:
Benutzeravatar
kbr
User
Beiträge: 1487
Registriert: Mittwoch 15. Oktober 2008, 09:27

Für alle, die sich auf die schnelle weder C noch Assembler einlesen können, habe ich ein einfaches Sieb in Python umgesetzt und dem erwähnten 'algo' gegenübergestellt. Diesen habe ich zudem etwas modifiziert, damit die gleiche Datenstruktur zurück gegeben wird. Ansonsten habe ich den 'algo' unverändert gelassen. Es fällt sofort auf, dass dieser nicht schneller sein kann. Das Sieb berechnet die ersten 2000000 Primzahlen in unter 1 sec. Bei dem 'algo' wage ich so große Zahlen gar nicht erst.

Code: Alles auswählen

def sieve(n):
    primes = [True] * n
    for i in range(2, n):
        if primes[i]:
            for j in range(i+i, n, i):
                primes[j] = False
    return primes

def algo(n):
    primes = [True] * n
    for i in range(2, n):
        for j in range(2, i):
            if not i % j:
                primes[i] = False
                break
    return primes
BlackJack

@kbr: Du hast `algo` falsch umgesetzt. Du müsstest die innere Schleife abbrechen sobald ``i // j < j`` gilt. Das hat DerTürke gemacht um nur bis j nur bis sqrt(i) laufen zu lassen. Wobei sich das in Assembler anbietet weil man den Modulo-Wert und den Quotienten ”kostenlos” von der Division bekommt die man sowieso durchführt.
DerTürke

@Sirius:
Du hasst dir wohl meinen Code nicht angeschaut, ansonsten würdest du nicht behaupten ich hätte das Sieb ( was ja auch vom Schwierigkeitsgrad nicht sehr fordernd ist ) nicht implementiert, vielmehr habe ich auch ein 2 Array of int (also wirklich keine List kein vector oder der Gleichen) benutzt um auch die im Sieb ermittelten Primzahlen zu sammeln, das mit der Ahnung war von meiner Seite nicht so gemeint sorry falls ich dir damit zu Nahe getreten bin. Ich glaube ich bin hier sowieso fehl am Platz ich sollte mich am besten in Foren aufhalten wo man keine so große Antipathie gegen Assembler pflegt. Nochmal auch an BlackJack sorry wenn ich mit meinen Äußerungen euch zu Nahe getreten bin ich wünsche euch allen gutes gelingen und eine schöne Zeit.
Benutzeravatar
kbr
User
Beiträge: 1487
Registriert: Mittwoch 15. Oktober 2008, 09:27

@BlackJack: Stimmt, das hatte ich übersehen und mich schon gewundert. Diese Optimierung hatte ich ursprünglich auch angewandt, was am Ergebnis aber nichts ändert.
Antworten