#include<stdlib.h>
#include <stdio.h>
#include<time.h>
#include <omp.h>
#define NUM_THREADS 2
#define Nu 1000  // Lato della matrice quadrata
int main() {
	// clock_t ticks1, ticks2;
	omp_set_num_threads(NUM_THREADS);
	int i,j,k;
	double *A, *B, *C, *D, tmp, start_time, end_time, ok_time, bad_time;
	A = (double *)malloc(Nu*Nu*sizeof(double)); // Riserva la memoria per
	B = (double *)malloc(Nu*Nu*sizeof(double)); 
	C = (double *)malloc(Nu*Nu*sizeof(double));
	D = (double *)malloc(Nu*Nu*sizeof(double));	
	printf(" Matrix line lenght (bytes): %d \n",Nu*sizeof(double));
	
	
// Inizializza le matrici A e B

start_time=omp_get_wtime(); // fa partire il timer
#pragma omp parallel for 	
	for(j=0;j<Nu;j++)
	{
		for(i=0;i<Nu;i++)
		{			
			*(A+(j*Nu+i))=i+j;
			*(B+(j*Nu+i))=i-j;
		}
	}
	
	
	
// Prepara la matrice trasposta  D ( Trasposta di B )
#pragma omp  parallel for 
	for (k=0; k<Nu; k++)
	{
		for (j=0; j<Nu; j++)
		{
			*(D+(k*Nu+j))= *(B+(j*Nu+k));
		}
	}

//#pragma omp parallel for
	for (i=0; i<Nu; i++)
	{
		for (j=0; j<Nu; j++)
		{
			*(C+(i*Nu+j))=0.0;
#pragma omp parallel for
			for(k=0;k<Nu ;k++)
			{
				*(C+(i*Nu+j))  += *(A+(i*Nu+k)) *  *(D+(j*Nu+k)); /* C(i,j) = sum(over k) A(i,k) * B(k,j) */
			}
		}
	}

	end_time=omp_get_wtime();
	ok_time=end_time-start_time;
	printf("Tempo ok : %f  \n",ok_time);
	start_time=omp_get_wtime();
	
#pragma omp parallel for 	
	for(j=0;j<Nu;j++){
		for(i=0;i<Nu;i++){			
			*(A+(i*Nu+j))=i+j;
			*(B+(i*Nu+j))=i-j;
		}
	}
	
#pragma omp parallel for	
	for (i=0; i<Nu; i++)
	{
		for (j=0; j<Nu; j++)
		{
			*(C+(i*Nu+j))=0.0;
			for(k=0;k<Nu ;k++)
			{
				*(C+(i*Nu+j))  += *(A+(i*Nu+k)) *  *(B+(k*Nu+j)); /* C(i,j) = sum(over k) A(i,k) * B(k,j) */
			}
			
		}	
	}
		end_time=omp_get_wtime();
		bad_time=end_time-start_time;
		printf("Tempo bad : %f  \n",bad_time);
		printf("Rapporto : %f , Per una matrice quadrata di lato %d\n",bad_time/ok_time, Nu );
	
	
}