     {
         printf("  v[%u] = %f\n", i, v[i]);
     }
-    
-    stat = cublasGetMatrix (dim, dim, sizeof(*block_pair), dBlockPair, dim, block_pair, dim);
-    if (stat != CUBLAS_STATUS_SUCCESS) {
-        printf ("data upload failed");
-        cudaFree (dBlockPair);
-        cublasDestroy(handle);
-        return EXIT_FAILURE;
-    }
-    
-    cudaFree (dBlockPair);
-    cublasDestroy(handle);
 
     free(v);
     

    
    if(x < dim && y < dim)
    {

    }
}

/* Compute the new block pair using computed Q matrix */
__global__ void UpdateBlockPair(float *block_pair, float *q, int dim)
{
    int x = threadIdx.x + blockIdx.x*blockDim.x;
    int y = threadIdx.y + blockIdx.y*blockDim.y;
    
    if(x < dim && y < dim)
    {

    }
}

/* Tridiagonalize first row & colum of block pair by computing Q using Householder vector */
/*void UpdateBlockPair(float *block_pair, float *v, float *q, int dim)
{
    //ComputeQ();
    //ComputeBlockPair();
}*/

// Reduce the matrix to a tridiagonal matrix via Householder transformations
int BlockPairReduction(float *q, float *column, float *block_pair, int dim)
{
    cublasStatus_t stat;
    cublasHandle_t handle;
    float *dBlockPair;
    
    cudaMalloc(&dBlockPair, (size_t)dim*dim*sizeof(*block_pair)); CUDA_CHECK;
    
    stat = cublasCreate(&handle);
    if (stat != CUBLAS_STATUS_SUCCESS) {
        printf ("CUBLAS initialization failed\n");
        return EXIT_FAILURE;
    }
    
    stat = cublasSetMatrix (dim, dim, sizeof(*block_pair), block_pair, dim, dBlockPair, dim);
    if (stat != CUBLAS_STATUS_SUCCESS) {
        printf ("data download failed");
        cudaFree (dBlockPair);
        cublasDestroy(handle);
        return EXIT_FAILURE;
    }

    /* HOUSEHOLDER TRANSFORMATION */

    float alpha = 0.0;

    // TODO: Put all GPU operations into one kernel

    ComputeHouseVec<<<grid, block>>>(gpuV, gpuColumn, alpha, r, dim);
    UpdateBlockPair(block_pair, gpuV, q, dim);

    // TODO: Retreive column!

    // Copy data back to CPU and free GPU memory
    cudaMemcpy(v, gpuV, (size_t)dim*sizeof(float), cudaMemcpyDeviceToHost); CUDA_CHECK;
