cblas_wrappers/cuda_auxiliary.hpp

#ifndef CUDA_AUXILIARY_H_
#define CUDA_AUXILIARY_H_

#ifdef HAVE_LIBCUBLAS
#include "fupermod_cblas.h"
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <cublas_v2.h>
#include <cuda_runtime.h>

#define gama 0.6
#define PRINT 1

/* CUBLAS status type returns: 
typedef enum{
      CUBLAS_STATUS_SUCCESS          =0,
      CUBLAS_STATUS_NOT_INITIALIZED  =1,
      CUBLAS_STATUS_ALLOC_FAILED     =3,
      CUBLAS_STATUS_INVALID_VALUE    =7,
      CUBLAS_STATUS_ARCH_MISMATCH    =8,
      CUBLAS_STATUS_MAPPING_ERROR    =11,
      CUBLAS_STATUS_EXECUTION_FAILED =13,
      CUBLAS_STATUS_INTERNAL_ERROR   =14
} cublasStatus_t;
*/

#  define CUBLAS_SAFE_CALL( call ) {                                         \
    cublasStatus_t err = call;                                               \
    if( CUBLAS_STATUS_SUCCESS != err) {                                      \
        fprintf(stderr, "Cublas error %d in file '%s' in line %i.\n",        \
                err, __FILE__, __LINE__ );                                   \
        exit(EXIT_FAILURE);                                                  \
    } }

inline cublasOperation_t trans_cblas2cublas(enum CBLAS_TRANSPOSE trans) {
  return (cublasOperation_t)(trans - CblasNoTrans);
}

#ifdef ENABLE_BLAS_SP
#define CUBLAS_GEMM cublasSgemm
#else
#define CUBLAS_GEMM cublasDgemm
#endif

#define hA_ref(i)   (A+(i)*(offset_a))
#define hB_ref(i)   (B+(i)*(offset_b))
#define hC_ref(i,j) (C+(i)*(offset_c)+(j)*(mm)*(N))
#define dA_ref(i)   (dA+(i)*(offset_a))
#define dB_ref(i)   (dB+(i)*(offset_b))
#define dC_ref(i,j) (dC+(i)*(offset_c)+(j)*(mm)*(dN))

#define hA_rref(i)   ((i)*(offset_a))
#define hB_rref(i)   ((i)*(offset_b))
#define hC_rref(i,j) ((i)*(offset_c)+(j)*(mm)*(N))
#define dA_rref(i)   ((i)*(offset_a))
#define dB_rref(i)   ((i)*(offset_b))
#define dC_rref(i,j) ((i)*(offset_c)+(j)*(mm)*(dN))

typedef struct cublasxgemm_params {
  fupermod_float **dA, **dB, **dC;
  // try to avoid conflicts with IB when allocating huge pinned memory for C
  // it will only allocate pinned memory of the memory capacuties of devices
  // use stream callback function to schedule data transfers
  fupermod_float **pC;
  fupermod_float *workspace_WA, *workspace_WB;
  // somehow MPI_Bcast does not work on pinned memory, so workspace_WB is 
  fupermod_float *dma_WB;
  cublasHandle_t* handles;
  int idx_outer; 
  int num_iteration; // note: starting from 1, not 0
  int kernel_type;
  int num_device_use;
    
  // parameters from conf
  int device_id;
  int num_devices;
  int memory_mode;
  int algorithm;
  int is_overlap;
} cublasxgemm_params;

inline size_t getDeviceMemory(int device) {
  if (cudaGetDevice(&device) != cudaSuccess)
    return 0;
  struct cudaDeviceProp properties;
  if (cudaGetDeviceProperties(&properties, device) != cudaSuccess)
    return 0;
  return properties.totalGlobalMem; //bytes
}

// get the sum of an array
inline int summer(int* array, int size)
{
  if (size == 0) {
    return 0;
  } else {
    int i = 0;
    int sum = 0;
    for (i = 0; i < size; i++) {
      sum += *(array + i);
    }
    return sum;
  }
}

// get the address of a block from index. Plz refer to my matrix multiplicaiton kernel documentation
inline fupermod_float* addr(const fupermod_float* ptr, int i, int j, int* height, int* width, int N, int K, int flag)
{
  fupermod_float* address = (fupermod_float*)ptr;
  if (flag == 1) // A(i)
    address += summer(height, i) * K * K;
  else if (flag == 2) // B(i)
    address += summer(width, j) * K;
  else if (flag == 3) // C(i,j)
    address += summer(height, i) * K * N + summer(width, j) * K;
  else if (flag == 4) // A_(i)
    address += i * height[0] * K * K;
  else // C_(i)
    address += i * height[0] * width[0] * K * K;
    
  return address;
}

// get the address of a block from index. Plz refer to my matrix multiplicaiton kernel documentation
inline fupermod_float* addr1D(const fupermod_float* ptr, int i, int j, int* height, int* width, int h_, int w_, int K, int flag)
{
  fupermod_float* address = (fupermod_float*)ptr;
  if (flag == 1) // A(i)
    address += summer(height, i) * K * K;
  else if (flag == 2) // B(i)
    address += summer(width, j) * K * K;
  else if (flag == 3) {// C(i,j) 
    int l;
    for (l = 0; l < j; l++) {
      address += summer(height, h_) * K * width[j] * K;
    }
    address += summer(height, i) * K * width[j] * K;
  }
  else if (flag == 4) // A_(i)
    address += i * height[0] * K * K;
  else // C_(i)
    address += i * height[0] * width[0] * K * K;
    
  return address;
}

// struct for transfer data between C and temporary pinned buffer for C
typedef struct transfer_t {
  cublasxgemm_params *params;
  fupermod_float *C;
  int i;
  int j;
  int buf;
  int *height;
  int *width;
  int N;
  int K;
  int id;
  int direction;
} transfer_t;

inline transfer_t *get_transfer(cublasxgemm_params *params, fupermod_float *C, int i, int j, int buf, int *h, int *w, int N, int K, int id, int direction) {
  transfer_t *t = (transfer_t *) malloc(sizeof(transfer_t));
  t->params = params;
  t->C = C;
  t->i = i;
  t->j = j;
  t->buf = buf;
  t->height = h;
  t->width = w;
  t->N = N;
  t->K = K;
  t->id = id;
  t->direction = direction;
  return t;
}
/*
void CUDART_CB transfer_cb(cudaStream_t stream, cudaError_t status, void *transfer) {
  if (status != cudaSuccess) {
    fprintf(stderr, "%s %s: call back function status error", __FILE__, __func__);
  }  
  transfer_t *t = (transfer_t *) transfer;
  cublasxgemm_params *params = t->params;
  fupermod_float *C = t->C;
  int i = t->i;
  int j = t->j;
  int buf = t->buf;
  int *height = t->height;
  int *width = t->width;
  int N = t->N;
  int K = t->K;
  int id = t->id;
  int direction = t->direction;
  
  size_t count = sizeof(fupermod_float) * width[j+id] * K;
  int l;
  if (direction == 0) {
    // C --> pC, followed by pC --> dC
    fupermod_float *src = addr(C,i,j+id,height,width,N,K,3);
    fupermod_float *des = addr(params->pC[id],buf,j,height,width,N,K,5);
    for (l = 0; l < height[i] * K; l++) {
      memcpy(des + l * width[0] * K, src + l * N, count);
    }
  } else {
    // pC --> C, after dC --> pC done
    fupermod_float *src = addr(params->pC[id],buf,j,height,width,N,K,5);
    fupermod_float *des = addr(C,i,j+id,height,width,N,K,3);
    for (l = 0; l < height[i] * K; l++) {
      memcpy(des + l * N, src + l * width[0] * K, count);
    }
  }
  free(t); t = NULL;
}*/

/*-------------------------------------------------------------------------------------------------------------------*/
void cuda_gemm_basic(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);

void cuda_gemm(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);

/*-------------------------------------------------------------------------------------------------------------------*/
void gemm_kernel_gpu_single_fit(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);

void gemm_kernel_gpu_single_exceed(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);

void gemm_kernel_gpu_single_exceed_0(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);

void gemm_kernel_gpu_single_exceed_1(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);

void gemm_kernel_gpu_single_exceed_2(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);

void gemm_kernel_gpu_single_exceed_3(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);
/*
void gemm_kernel_gpu_single_exceed_4(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);
*/          
/*-------------------------------------------------------------------------------------------------------------------*/
void gemm_kernel_gpu_server_fit(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);

void gemm_kernel_gpu_server_exceed(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);

void gemm_kernel_gpu_server_exceed_nd2_1(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);

void gemm_kernel_gpu_server_exceed_nd2_2(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);

void gemm_kernel_gpu_server_exceed_nd2_3(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);
          
void gemm_kernel_gpu_server_exceed_np2_1(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);

void gemm_kernel_gpu_server_exceed_np2_2(fupermod_gemm* gemm,
          const enum CBLAS_ORDER Order, const enum CBLAS_TRANSPOSE TransA, const enum CBLAS_TRANSPOSE TransB,
          const int M, const int N, const int K, const fupermod_float alpha, const fupermod_float *A, const int lda,
          const fupermod_float *B, const int ldb, const fupermod_float beta, fupermod_float *C, const int ldc);

/*-------------------------------------------------------------------------------------------------------------------*/
void matrix_split_single_fit(int M, int N, int K, int *h_, int *w_, int **height, int **width);
void matrix_split_single_exceed(int M, int N, int K, int mem_size, int *h_, int *w_, int **height, int **width);
void matrix_split_server_exceed(int M, int N, int K, int mem_size_device, int num_device_use, int *h_, int *w_, int **height, int **width);

/*-------------------------------------------------------------------------------------------------------------------*/
#endif          
#endif /* CUDA_AUXILIARY_H_ */