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


#define TOL  0.5E-9  /*tolerance*/ 
#define MAXIT 100000   /*nombre maximale d'iterations*/
//#define N 360       /*nombre de colonne de la matrice*/

double t1, texe;

#define CPUID __asm __emit 0fh __asm __emit 0a2h
#define RDTSC __asm __emit 0fh __asm __emit 031h

unsigned __int64 dtime(){
		unsigned cycles_low, cycles_high;

		__asm {
				pushad
				CPUID			
				RDTSC
				mov	cycles_high, edx			
				mov	cycles_low, eax	
				popad
				}
			return ((unsigned __int64)cycles_high << 32) | cycles_low;
	}


int main(int argc, char * argv[]){
	double error, resid;
	int i,j,nbiter;
	int iter = 0;
	int stop = 0;
	int N;
	float **u; 
	float **uold;

	
	// Nombre d'iterations.
	printf("Choisissez la taille de la matrice : ");
	scanf(" %d", &N); 
	printf("N : %d \n", N);
		
	// Nombre d'iterations.
	printf("Choisissez votre nombre d'iterations: ");
	scanf(" %d", &nbiter); 
	if( nbiter > MAXIT ) 
		nbiter = MAXIT;

	u = (float **)malloc(sizeof(float *)*N);
	u[0] = (float *)malloc(sizeof(float)*N*N);
	uold = (float **)malloc(sizeof(float *)*N);
	uold[0] = (float *)malloc(sizeof(float)*N*N);
	for(i=1;i<N;i++){
        u[i]=u[i-1]+N;
		uold[i]=uold[i-1]+N;}



	
	//initialisation de la matrice
	for(i = 0; i < N ; i++)
		for(j = 0; j < N ; j++)
			u[i][j] = 0.0;
	//initialisation de la ligne du haut	
	for(i=0 ; i<N ;i++)
		u[0][i] = 1.0;
	
	
	error = 10.0 * TOL;

	t1=dtime();

	while (iter < nbiter && stop!=1) 
	{
	
		error = 0.0;
		//sauvegarde de l'ancienne matrice
		for(i = 0 ; i < N; i++)
			for(j = 0 ; j < N ; j++)
				uold[i][j] = u[i][j]; 


		for(i = 1 ; i < N-1; i++)
			for(j = 1 ; j < N-1 ; j++){
				u[i][j] = 0.25 * (uold[i-1][j] + uold[i+1][j] + uold[i][j-1] + uold[i][j+1]);
				resid = u[i][j] - uold[i][j];
				error += resid * resid;
			}
		iter++;
		error = sqrt(error)/(N*N);
		//printf("iteration %d : erreur : %f \n",iter,error);
		
		if(error < TOL) {
			stop = 1;
			//printf("test erreur %f",testerreur);
		}		
		
	}
//	printf("matrice resultante\n");
		
		//affichage de la matrice
/*	for(i = 0;i<N;i++){
		for(j = 0;j<N;j++){
			printf(" %f ",u[i][j]);
		}
		printf("\n");
	}	
*/	
	texe=dtime()-t1;
	
	printf( "CPP =  %f \n ", (double)texe/(N*N));
	printf ("iterations = %d \n", iter-1);
}