1: /* 2: Management of CUBLAS and CUSOLVER handles 3: */ 5: #include <petscsys.h> 6: #include <petsc/private/petscimpl.h> 7: #include <petsccublas.h> 9: static cublasHandle_t cublasv2handle = NULL; 10: static cusolverDnHandle_t cusolverdnhandle = NULL; 12: /* 13: Destroys the CUBLAS handle. 14: This function is intended and registered for PetscFinalize - do not call manually! 15: */ 16: static PetscErrorCode PetscCUBLASDestroyHandle() 17: { 18: cublasStatus_t cberr; 21: if (cublasv2handle) { 22: cberr = cublasDestroy(cublasv2handle);CHKERRCUBLAS(cberr); 23: cublasv2handle = NULL; /* Ensures proper reinitialization */ 24: } 25: return(0); 26: } 28: /* 29: Initializing the cuBLAS handle can take 1/2 a second therefore 30: initialize in PetscInitialize() before being timing so it does 31: not distort the -log_view information 32: */ 33: PetscErrorCode PetscCUBLASInitializeHandle(void) 34: { 36: cublasStatus_t cberr; 39: if (!cublasv2handle) { 40: cberr = cublasCreate(&cublasv2handle);CHKERRCUBLAS(cberr); 41: /* Make sure that the handle will be destroyed properly */ 42: PetscRegisterFinalize(PetscCUBLASDestroyHandle); 43: } 44: return(0); 45: } 47: PetscErrorCode PetscCUBLASGetHandle(cublasHandle_t *handle) 48: { 53: if (!cublasv2handle) {PetscCUBLASInitializeHandle();} 54: *handle = cublasv2handle; 55: return(0); 56: } 58: /* cusolver */ 59: static PetscErrorCode PetscCUSOLVERDnDestroyHandle() 60: { 61: cusolverStatus_t cerr; 64: if (cusolverdnhandle) { 65: cerr = cusolverDnDestroy(cusolverdnhandle);CHKERRCUSOLVER(cerr); 66: cusolverdnhandle = NULL; /* Ensures proper reinitialization */ 67: } 68: return(0); 69: } 71: PetscErrorCode PetscCUSOLVERDnInitializeHandle(void) 72: { 73: PetscErrorCode ierr; 74: cusolverStatus_t cerr; 77: if (!cusolverdnhandle) { 78: cerr = cusolverDnCreate(&cusolverdnhandle);CHKERRCUSOLVER(cerr); 79: PetscRegisterFinalize(PetscCUSOLVERDnDestroyHandle); 80: } 81: return(0); 82: } 84: PetscErrorCode PetscCUSOLVERDnGetHandle(cusolverDnHandle_t *handle) 85: { 86: PetscErrorCode ierr; 90: if (!cusolverdnhandle) {PetscCUSOLVERDnInitializeHandle();} 91: *handle = cusolverdnhandle; 92: return(0); 93: }