Discussion :

[yan]

Bonjour.
Pour ecrire un bout de code fractal, j'ai voulue utiliser std::complex. Seulement c'est beaucoup plus lent que d'ecrire le calcul à la main

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
 
double zx = xf + xRatio*j;
double zy = yf + yRatio*i;
double cx = zx;
double cy = zy;
 
while(buf[j]<itermax && zx*zx+zy*zy<4)
            {
		double xnew = zx * zx - zy * zy + cx;
		zy = (2 * zx * zy )+ cy;
                zx = xnew;
                ++(buf[j]);
 
            }

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
 
z = c = std::complex<double>( xf + xRatio*j, yf + yRatio*i) ;	
while(buf[j]<itermax && std::abs(z)<4)
            {
                z = z*z+c;
                ++(buf[j]);
 
            }

es ce normale? je développe sous visual c++ 2008 express

[yan]

J'ai trouvé cela :
http://www.ddj.com/linux-open-source/199702312?pgno=1

Application of std::complex<T> could result in some shorter code, but the computational performance, which I try to maximize, would be lost.

y as rien a faire??

[droggo]

Goe,

y as rien a faire??

Si, coder toi-même ce qu'il faut, c'est à dire pas grand chose au total, en fait.

En évitant de faire une classe si la vitesse est primordiale, ce sera toujours ça de gagné (structure + les fonctions ad-hoc est suffisant [Je vais me faire insulter par les puristes, mais il faut savoir ce qu'on veut]).

[yan]

En faite ce qui me choque c'est que c'est 5 à 10 fois plus lent avec les complex!!!!
je trouve cela enorme!!! non?
je m'attendais à quasiment le même temp

[yan]

y aurais peut être une début d'explication
http://readlist.com/lists/gcc.gnu.or...help/0/12.html

ce serait à cause du test

std::abs(z)<4

Faudra que je teste

[edit]
http://shootout.alioth.debian.org/gp...ng=gcc%20&id=3
C'est bien sensé donner de trés bon resultat
http://shootout.alioth.debian.org/gp...&lang=gpp&id=0

[Sylvain Togni]

Cela n'a rien à voir avec les complexes. La différence de temps viens du fait que les deux codes ne font pas la même chose : le premier compare la norme au carré, le second la valeur absolue, c-à-d fait appel à sqrt().

[yan]

Cela n'a rien à voir avec les complexes. La différence de temps viens du fait que les deux codes ne font pas la même chose : le premier compare la norme au carré, le second la valeur absolue, c-à-d fait appel à sqrt().

Ha oui merci, je me suis embrouillé entre abs et norm

Normalement j'avais aussi essayé avec std::norm et c'était pareil...
Je vais vérifier ce soir

[nicroman]

Pas tant que ça... c'est pas du tout le même code !

Regardons les codes:
Ton code:

Code C++ :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
 
double zx = xf + xRatio*j;
// 1 addition, 1 multiplication
double zy = yf + yRatio*i;
// 1 addition, 1 multiplication
double cx = zx;
// 0
double cy = zy;
// 0
 
while(buf[j]<itermax && zx*zx+zy*zy<4) 
{
       // 1 indirection (pas en release), 2 comparaisons, 2 multiplications, 3 additions
       double xnew = zx * zx - zy * zy + cx;
       // 2 multiplications, 1 addition  A noter qu'en mode release, le compilo ne fait pas les 2 multiplications (déjà effectuées pour le test de la boucle)
       zy = (2 * zx * zy )+ cy;
       // 2 multiplications, 1 addition
       zx = xnew;
       // 1 affectation
       ++(buf[j]);
       // 1 incrémentation (à noter qu'en mode release cette opération est fait sur un registre qui est affecté après la boucle)
}

Résultat, disons pour 10 itérations: 52 additions, 42 multiplications, 20 comparaisons, 1 affectation avec indirection

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
 
z = c = std::complex<double>( xf + xRatio*j, yf + yRatio*i) ;
// idem code précédant
while(buf[j]<itermax && std::abs(z)<4)
{
      // std::abs == sqrt(std::norm(z))
      // donc: 2 comparaisons, 2 multiplications, 3 additions, 1 square-root
      z = z*z+c;
      // le compilo n'a pas fait d'optimisation ici: 
      // donc: 4 multiplications (2 par opérande du z*z), 4 additions (2 par opérande du z*z et 2 pour le +c)
      // a rajouter 1 construction temporaire d'objet, et 2 affectations
      ++(buf[j]);
      // idem code précédant
}

D'une part la condition de sortie de la boucle n'est pas la même (equivalent à <16 dans le premier code), ce qui doit induire plus d'itérations.
Même à nombre d'itération identique (10) on aura:
72 additions, 62 multiplications, 20 comparaisons, 21 affectations, 10 créations d'objet temporaire, 10 square-root

[yan]

Pas tant que ça... c'est pas du tout le même code !

Merci, très intéressant.
Par contre

double xnew = zx * zx - zy * zy + cx;
// 2 multiplications, 1 addition A noter qu'en mode release, le compilo ne fait pas les 2 multiplications (déjà effectuées pour le test de la boucle)

moi je voie 2 multiplications et 2 addition (et 1 changement de signe )
ce qui ferait avec la seconde instruction :
4 multiplication et 3 addition

z = z*z+c;
// le compilo n'a pas fait d'optimisation ici:
// donc: 4 multiplications (2 par opérande du z*z), 4 additions (2 par opérande du z*z et 2 pour le +c)
// a rajouter 1 construction temporaire d'objet, et 2 affectations

pourquoi ne serait'il pas optimiser et comment tu trouve

4 multiplications (2 par opérande du z*z), 4 additions (2 par opérande du z*z et 2 pour le +c)

[yan]

voici un code de teste :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
#include <iostream>
#include <fstream>
#include <complex>
#include <ctime>
 
int main(int argc, char* argv[])
{
    const int w = 1000;
    const int h = 1000;
    const int iterMax = 1000;
    const double p[] ={-2.,-2.};
    const double size[] ={4.,4.};
 
    const double xRatio =size[0]/(w-1);
    const double yRatio =size[1]/(h-1);
 
    std::complex<double> z;
    std::complex<double> c;
 
 
    int * pBuffer = new int[w*h];
 
    int t1 = clock();
    for (int i=0;i<h;++i)
    {
        for (int j=0;j<w;++j)
        {
        int id = i*w+j;
        int & refRes = pBuffer[id];
        refRes =0;
 
        z = c = std::complex<double>( p[0] + xRatio*j, p[1] + yRatio*i) ;
        while(refRes<iterMax && std::norm(z)<4)
            {
            #if 1
                z = z*z+c;
            #else
                z*=z;
                z+=c;
            #endif
                ++refRes;
            }
        if (refRes==iterMax) refRes =0;
 
        }
    }
    std::cout<< (double)(clock() - t1)/CLOCKS_PER_SEC<<std::endl;
 
 
    std::ofstream fout("image.pgm");
    fout<<"P2"<<std::endl;
    fout<<w<<" "<<h<<std::endl;
    fout<<"256"<<std::endl;
    const double max =.1*iterMax;
    for (int i=0;i<w*h;++i)
    {
        fout<<(int)(255.*(pBuffer[i]>max ? 1. : (double)pBuffer[i]/max))<<" ";
    }
 
 
 
    system("pause");
 
	return 0;
}

methode 1 : z =z*z+c
methode 2 : z =*z; z+=c;

avec mingw :
methode 1 : a peu prés 7
methode 2 : a peu prés 6.5

avec visual 2005 :
methode 1 : a peu prés 12.5
methode 2 : a peu prés 3.5