// Suma_Matrices.cu : Defines the entry point for the console application.
//

#include <stdio.h>

/********************************************/
/********************************************/
// Función Kernel que se ejecuta en el Device.
__global__ void Suma_Matrices(float *c,float *a,float *b, 
							  int fil,int col)
{
  int idx = blockIdx.x * blockDim.x + threadIdx.x;  
  // Escribe el identificador de la posicion en "y"
  int idy = blockIdx.y* blockDim.y + threadIdx.y; 

  // Posicion del indice en el vector
  int index=idx*col+idy;
  
  // Complete the condition on idx
  if (idx<fil && idy<col){
	  c[index] = a[index] + b[index];
  }	
}

int div_up(int a,int b){
	return(a/b + (a%b == 0 ? 0:1));
}

/********************************************/
/********************************************/
// Código principal que se ejecuta en el Host
int main(void){
	float *a_h,*b_h,*c_h; //Punteros a matrices en el Host 
	float *a_d,*b_d,*c_d;  //Punteros a matrices en el Device
	int fil = 16;  //Número de filas
	int col = 16;  //Número de columnas
	int N=fil*col;  //Número de elementos de la matriz

	
	size_t size=N * sizeof(float);

	a_h = (float *)malloc(size); // Pedimos memoria en el Host
	b_h = (float *)malloc(size);
	c_h = (float *)malloc(size);//También se puede con cudaMallocHost
	
	//Inicializamos las matrices a,b en el Host
	for (int i=0; i<fil; i++){
		for(int j=0;j<col;j++){
			a_h[i*col+j] = (float)i;
			b_h[i*col+j] = (float)(i+1);	
		}
	}

	printf("\nMatriz a:\n");
	for (int i=0; i<fil; i++){
		for(int j=0;j<col;j++){
			printf("%.2f ", a_h[i*col+j]);
		}
		printf("\n");
	}
	printf("\n\nMatriz b:\n");
	for (int i=0; i<fil; i++){
		for(int j=0;j<col;j++){
			printf("%.2f ", b_h[i*col+j]);
		}
		printf("\n");
	}
	
	cudaMalloc((void **) &a_d,size);   // Pedimos memoria en el Device
	/*HERE*/((void **) &b_d,size);
	cudaMalloc((void **) &c_d,size);
	
	//Pasamos las matrices a y b del Host al Device
	/*HERE*/(a_d, a_h, size, cudaMemcpyHostToDevice);
	cudaMemcpy(b_d, b_h, size, cudaMemcpyHostToDevice);
	
	//Realizamos el culculo en el Device
	dim3 block_size(4,4);
	dim3 n_blocks(div_up(fil,block_size.x),div_up(col,block_size.y)) ;
	Suma_Matrices <<< /*HERE*/, block_size >>> (c_d,a_d,b_d,fil,col);
		
	//Pasamos el resultado del Device al Host
	/*HERE*/(c_h, c_d, size, cudaMemcpyDeviceToHost);
    
    
	//Resultado
	printf("\n\nMatriz c:\n");
	for (int i=0; i<fil; i++){
		for(int j=0;j<col;j++){
			printf("%.2f ", c_h[i*col+j]);
		}
		printf("\n");
	}
	
		
	// Liberamos la memoria del Host
	free(a_h); 
	free(b_h); 
	free(c_h); 

	// Liberamos la memoria del Device
	/*HERE*/(a_d); 
	cudaFree(b_d); 
	cudaFree(c_d); 
	return(0);
}


