Actual source code: handle.hip.cpp

  1: /*
  2:  Management of HIPBLAS and HIPSOLVER handles

  4:  Unlike CUDA, hipSOLVER is just for dense matrices so there is
  5:  no distinguishing being dense and sparse.  Also, hipSOLVER is
  6:  very immature so we often have to do the mapping between roc and
  7:  cuda manually.
  8:  */

 10: #include <petscsys.h>
 11: #include <petsc/private/petscimpl.h>
 12: #include <petscdevice.h>

 14: static hipblasHandle_t     hipblasv2handle   = NULL;
 15: static hipsolverHandle_t   hipsolverhandle = NULL;

 17: /*
 18:    Destroys the HIPBLAS handle.
 19:    This function is intended and registered for PetscFinalize - do not call manually!
 20:  */
 21: static PetscErrorCode PetscHIPBLASDestroyHandle()
 22: {
 23:   hipblasStatus_t cberr;

 26:   if (hipblasv2handle) {
 27:     cberr          = hipblasDestroy(hipblasv2handle);CHKERRHIPBLAS(cberr);
 28:     hipblasv2handle = NULL;  /* Ensures proper reinitialization */
 29:   }
 30:   return(0);
 31: }

 33: /*
 34:     Initializing the hipBLAS handle can take 1/2 a second therefore
 35:     initialize in PetscInitialize() before being timing so it does
 36:     not distort the -log_view information
 37: */
 38: PetscErrorCode PetscHIPBLASInitializeHandle(void)
 39: {
 41:   hipblasStatus_t cberr;

 44:   if (!hipblasv2handle) {
 45:     cberr = hipblasCreate(&hipblasv2handle);CHKERRHIPBLAS(cberr);
 46:     /* Make sure that the handle will be destroyed properly */
 47:     PetscRegisterFinalize(PetscHIPBLASDestroyHandle);
 48:   }
 49:   return(0);
 50: }

 52: PetscErrorCode PetscHIPBLASGetHandle(hipblasHandle_t *handle)
 53: {

 58:   if (!hipblasv2handle) {PetscHIPBLASInitializeHandle();}
 59:   *handle = hipblasv2handle;
 60:   return(0);
 61: }

 63: /* hipsolver */
 64: static PetscErrorCode PetscHIPSOLVERDestroyHandle()
 65: {
 66:   hipsolverStatus_t  cerr;

 69:   if (hipsolverhandle) {
 70:     cerr             = hipsolverDestroy(hipsolverhandle);CHKERRHIPSOLVER(cerr);
 71:     hipsolverhandle = NULL;  /* Ensures proper reinitialization */
 72:   }
 73:   return(0);
 74: }

 76: PetscErrorCode PetscHIPSOLVERInitializeHandle(void)
 77: {
 78:   PetscErrorCode    ierr;
 79:   hipsolverStatus_t  cerr;

 82:   if (!hipsolverhandle) {
 83:     cerr = hipsolverCreate(&hipsolverhandle);CHKERRHIPSOLVER(cerr);
 84:     PetscRegisterFinalize(PetscHIPSOLVERDestroyHandle);
 85:   }
 86:   return(0);
 87: }

 89: PetscErrorCode PetscHIPSOLVERGetHandle(hipsolverHandle_t *handle)
 90: {
 91:   PetscErrorCode     ierr;

 95:   if (!hipsolverhandle) {PetscHIPSOLVERInitializeHandle();}
 96:   *handle = hipsolverhandle;
 97:   return(0);
 98: }