Actual source code: cudavecimpl.h

petsc-master 2020-12-02
Report Typos and Errors

  4: #include <petscvec.h>
  5: #include <petsccublas.h>
  6: #include <petsc/private/vecimpl.h>

  8: typedef struct {
  9:   PetscScalar  *GPUarray;           /* this always holds the GPU data */
 10:   PetscScalar  *GPUarray_allocated; /* if the array was allocated by PETSc this is its pointer */
 11:   cudaStream_t stream;              /* A stream for doing asynchronous data transfers */
 12: } Vec_CUDA;

 14: PETSC_INTERN PetscErrorCode VecCUDAGetArrays_Private(Vec,const PetscScalar**,const PetscScalar**,PetscOffloadMask*);
 15: PETSC_INTERN PetscErrorCode VecDotNorm2_SeqCUDA(Vec,Vec,PetscScalar*, PetscScalar*);
 16: PETSC_INTERN PetscErrorCode VecPointwiseDivide_SeqCUDA(Vec,Vec,Vec);
 17: PETSC_INTERN PetscErrorCode VecWAXPY_SeqCUDA(Vec,PetscScalar,Vec,Vec);
 18: PETSC_INTERN PetscErrorCode VecMDot_SeqCUDA(Vec,PetscInt,const Vec[],PetscScalar*);
 19: PETSC_EXTERN PetscErrorCode VecSet_SeqCUDA(Vec,PetscScalar);
 20: PETSC_INTERN PetscErrorCode VecMAXPY_SeqCUDA(Vec,PetscInt,const PetscScalar*,Vec*);
 21: PETSC_INTERN PetscErrorCode VecAXPBYPCZ_SeqCUDA(Vec,PetscScalar,PetscScalar,PetscScalar,Vec,Vec);
 22: PETSC_INTERN PetscErrorCode VecPointwiseMult_SeqCUDA(Vec,Vec,Vec);
 23: PETSC_INTERN PetscErrorCode VecPlaceArray_SeqCUDA(Vec,const PetscScalar*);
 24: PETSC_INTERN PetscErrorCode VecResetArray_SeqCUDA(Vec);
 25: PETSC_INTERN PetscErrorCode VecReplaceArray_SeqCUDA(Vec,const PetscScalar*);
 26: PETSC_INTERN PetscErrorCode VecDot_SeqCUDA(Vec,Vec,PetscScalar*);
 27: PETSC_INTERN PetscErrorCode VecTDot_SeqCUDA(Vec,Vec,PetscScalar*);
 28: PETSC_INTERN PetscErrorCode VecScale_SeqCUDA(Vec,PetscScalar);
 29: PETSC_EXTERN PetscErrorCode VecCopy_SeqCUDA(Vec,Vec);
 30: PETSC_INTERN PetscErrorCode VecSwap_SeqCUDA(Vec,Vec);
 31: PETSC_EXTERN PetscErrorCode VecAXPY_SeqCUDA(Vec,PetscScalar,Vec);
 32: PETSC_INTERN PetscErrorCode VecAXPBY_SeqCUDA(Vec,PetscScalar,PetscScalar,Vec);
 33: PETSC_INTERN PetscErrorCode VecDuplicate_SeqCUDA(Vec,Vec*);
 34: PETSC_INTERN PetscErrorCode VecConjugate_SeqCUDA(Vec xin);
 35: PETSC_INTERN PetscErrorCode VecNorm_SeqCUDA(Vec,NormType,PetscReal*);
 36: PETSC_INTERN PetscErrorCode VecCUDACopyToGPU(Vec);
 37: PETSC_INTERN PetscErrorCode VecCUDAAllocateCheck(Vec);
 38: PETSC_EXTERN PetscErrorCode VecCreate_SeqCUDA(Vec);
 39: PETSC_INTERN PetscErrorCode VecCreate_SeqCUDA_Private(Vec,const PetscScalar*);
 40: PETSC_INTERN PetscErrorCode VecCreate_MPICUDA(Vec);
 41: PETSC_INTERN PetscErrorCode VecCreate_MPICUDA_Private(Vec,PetscBool,PetscInt,const PetscScalar*);
 42: PETSC_INTERN PetscErrorCode VecCreate_CUDA(Vec);
 43: PETSC_INTERN PetscErrorCode VecDestroy_SeqCUDA(Vec);
 44: PETSC_INTERN PetscErrorCode VecDestroy_MPICUDA(Vec);
 45: PETSC_INTERN PetscErrorCode VecAYPX_SeqCUDA(Vec,PetscScalar,Vec);
 46: PETSC_INTERN PetscErrorCode VecSetRandom_SeqCUDA(Vec,PetscRandom);
 47: PETSC_INTERN PetscErrorCode VecGetLocalVector_SeqCUDA(Vec,Vec);
 48: PETSC_INTERN PetscErrorCode VecRestoreLocalVector_SeqCUDA(Vec,Vec);
 49: PETSC_INTERN PetscErrorCode VecGetArrayWrite_SeqCUDA(Vec,PetscScalar**);
 50: PETSC_INTERN PetscErrorCode VecGetArray_SeqCUDA(Vec,PetscScalar**);
 51: PETSC_INTERN PetscErrorCode VecRestoreArray_SeqCUDA(Vec,PetscScalar**);
 52: PETSC_INTERN PetscErrorCode VecGetArrayAndMemType_SeqCUDA(Vec,PetscScalar**,PetscMemType*);
 53: PETSC_INTERN PetscErrorCode VecRestoreArrayAndMemType_SeqCUDA(Vec,PetscScalar**);
 54: PETSC_INTERN PetscErrorCode VecCopy_SeqCUDA_Private(Vec xin,Vec yin);
 55: PETSC_INTERN PetscErrorCode VecSetRandom_SeqCUDA_Private(Vec xin,PetscRandom r);
 56: PETSC_INTERN PetscErrorCode VecDestroy_SeqCUDA_Private(Vec v);
 57: PETSC_INTERN PetscErrorCode VecResetArray_SeqCUDA_Private(Vec vin);

 59: /* complex single */
 60: #if defined(PETSC_USE_COMPLEX)
 61: #if defined(PETSC_USE_REAL_SINGLE)
 62: #define cublasXaxpy(a,b,c,d,e,f,g)               cublasCaxpy((a),(b),(cuComplex*)(c),(cuComplex*)(d),(e),(cuComplex*)(f),(g))
 63: #define cublasXscal(a,b,c,d,e)                   cublasCscal((a),(b),(cuComplex*)(c),(cuComplex*)(d),(e))
 64: #define cublasXdotu(a,b,c,d,e,f,g)               cublasCdotu((a),(b),(cuComplex*)(c),(d),(cuComplex*)(e),(f),(cuComplex*)(g))
 65: #define cublasXdot(a,b,c,d,e,f,g)                cublasCdotc((a),(b),(cuComplex*)(c),(d),(cuComplex*)(e),(f),(cuComplex*)(g))
 66: #define cublasXswap(a,b,c,d,e,f)                 cublasCswap((a),(b),(cuComplex*)(c),(d),(cuComplex*)(e),(f))
 67: #define cublasXnrm2(a,b,c,d,e)                   cublasScnrm2((a),(b),(cuComplex*)(c),(d),(e))
 68: #define cublasIXamax(a,b,c,d,e)                  cublasIcamax((a),(b),(cuComplex*)(c),(d),(e))
 69: #define cublasXasum(a,b,c,d,e)                   cublasScasum((a),(b),(cuComplex*)(c),(d),(e))
 70: #define cublasXgemv(a,b,c,d,e,f,g,h,i,j,k,l)     cublasCgemv((a),(b),(c),(d),(cuComplex*)(e),(cuComplex*)(f),(g),(cuComplex*)(h),(i),(cuComplex*)(j),(cuComplex*)(k),(l))
 71: #define cublasXgemm(a,b,c,d,e,f,g,h,i,j,k,l,m,n) cublasCgemm((a),(b),(c),(d),(e),(f),(cuComplex*)(g),(cuComplex*)(h),(i),(cuComplex*)(j),(k),(cuComplex*)(l),(cuComplex*)(m),(n))
 72: #define cublasXgeam(a,b,c,d,e,f,g,h,i,j,k,l,m)   cublasCgeam((a),(b),(c),(d),(e),(cuComplex*)(f),(cuComplex*)(g),(h),(cuComplex*)(i),(cuComplex*)(j),(k),(cuComplex*)(l),(m))
 73: #else /* complex double */
 74: #define cublasXaxpy(a,b,c,d,e,f,g)               cublasZaxpy((a),(b),(cuDoubleComplex*)(c),(cuDoubleComplex*)(d),(e),(cuDoubleComplex*)(f),(g))
 75: #define cublasXscal(a,b,c,d,e)                   cublasZscal((a),(b),(cuDoubleComplex*)(c),(cuDoubleComplex*)(d),(e))
 76: #define cublasXdotu(a,b,c,d,e,f,g)               cublasZdotu((a),(b),(cuDoubleComplex*)(c),(d),(cuDoubleComplex*)(e),(f),(cuDoubleComplex*)(g))
 77: #define cublasXdot(a,b,c,d,e,f,g)                cublasZdotc((a),(b),(cuDoubleComplex*)(c),(d),(cuDoubleComplex*)(e),(f),(cuDoubleComplex*)(g))
 78: #define cublasXswap(a,b,c,d,e,f)                 cublasZswap((a),(b),(cuDoubleComplex*)(c),(d),(cuDoubleComplex*)(e),(f))
 79: #define cublasXnrm2(a,b,c,d,e)                   cublasDznrm2((a),(b),(cuDoubleComplex*)(c),(d),(e))
 80: #define cublasIXamax(a,b,c,d,e)                  cublasIzamax((a),(b),(cuDoubleComplex*)(c),(d),(e))
 81: #define cublasXasum(a,b,c,d,e)                   cublasDzasum((a),(b),(cuDoubleComplex*)(c),(d),(e))
 82: #define cublasXgemv(a,b,c,d,e,f,g,h,i,j,k,l)     cublasZgemv((a),(b),(c),(d),(cuDoubleComplex*)(e),(cuDoubleComplex*)(f),(g),(cuDoubleComplex*)(h),(i),(cuDoubleComplex*)(j),(cuDoubleComplex*)(k),(l))
 83: #define cublasXgemm(a,b,c,d,e,f,g,h,i,j,k,l,m,n) cublasZgemm((a),(b),(c),(d),(e),(f),(cuDoubleComplex*)(g),(cuDoubleComplex*)(h),(i),(cuDoubleComplex*)(j),(k),(cuDoubleComplex*)(l),(cuDoubleComplex*)(m),(n))
 84: #define cublasXgeam(a,b,c,d,e,f,g,h,i,j,k,l,m)   cublasZgeam((a),(b),(c),(d),(e),(cuDoubleComplex*)(f),(cuDoubleComplex*)(g),(h),(cuDoubleComplex*)(i),(cuDoubleComplex*)(j),(k),(cuDoubleComplex*)(l),(m))
 85: #endif
 86: #else /* real single */
 87: #if defined(PETSC_USE_REAL_SINGLE)
 88: #define cublasXaxpy  cublasSaxpy
 89: #define cublasXscal  cublasSscal
 90: #define cublasXdotu  cublasSdot
 91: #define cublasXdot   cublasSdot
 92: #define cublasXswap  cublasSswap
 93: #define cublasXnrm2  cublasSnrm2
 94: #define cublasIXamax cublasIsamax
 95: #define cublasXasum  cublasSasum
 96: #define cublasXgemv  cublasSgemv
 97: #define cublasXgemm  cublasSgemm
 98: #define cublasXgeam  cublasSgeam
 99: #else /* real double */
100: #define cublasXaxpy  cublasDaxpy
101: #define cublasXscal  cublasDscal
102: #define cublasXdotu  cublasDdot
103: #define cublasXdot   cublasDdot
104: #define cublasXswap  cublasDswap
105: #define cublasXnrm2  cublasDnrm2
106: #define cublasIXamax cublasIdamax
107: #define cublasXasum  cublasDasum
108: #define cublasXgemv  cublasDgemv
109: #define cublasXgemm  cublasDgemm
110: #define cublasXgeam  cublasDgeam
111: #endif
112: #endif

114: #endif