1 #include <petsc/private/petscimpl.h> 2 3 #if defined(PETSC_HAVE_CUDA) 4 #include <cuda_runtime.h> 5 #endif 6 7 #if defined(PETSC_HAVE_HIP) 8 #include <hip/hip_runtime.h> 9 #endif 10 11 static PetscInt petsc_checkpointer_intensity = 1; 12 13 /*@ 14 PetscCheckPointerSetIntensity - An intense pointer check registers a signal handler and attempts to dereference to 15 confirm whether the address is valid. An intensity of 0 never uses signal handlers, 1 uses them when not in a "hot" 16 function, and intensity of 2 always uses a signal handler. 17 18 Not Collective 19 20 Input Parameter: 21 . intensity - how much to check pointers for validity 22 23 Options Database Key: 24 . -check_pointer_intensity - intensity (0, 1, or 2) 25 26 Level: advanced 27 28 .seealso: `PetscCheckPointer()`, `PetscFunctionBeginHot()` 29 @*/ 30 PetscErrorCode PetscCheckPointerSetIntensity(PetscInt intensity) { 31 PetscFunctionBegin; 32 switch (intensity) { 33 case 0: 34 case 1: 35 case 2: petsc_checkpointer_intensity = intensity; break; 36 default: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Intensity %" PetscInt_FMT " not in 0,1,2", intensity); 37 } 38 PetscFunctionReturn(0); 39 } 40 41 /* ---------------------------------------------------------------------------------------*/ 42 43 #if defined(PETSC_HAVE_SETJMP_H) 44 #include <setjmp.h> 45 static jmp_buf PetscSegvJumpBuf; 46 static PetscBool PetscSegvJumpBuf_set; 47 48 /*@C 49 PetscSignalSegvCheckPointerOrMpi - To be called from a signal handler for SIGSEGV. If the signal was received 50 while executing PetscCheckPointer()/PetscCheckMpiXxxAwareness(), this function longjmps back there, otherwise returns 51 with no effect. This function is called automatically by PetscSignalHandlerDefault(). 52 53 Not Collective 54 55 Level: developer 56 57 .seealso: `PetscPushSignalHandler()` 58 @*/ 59 void PetscSignalSegvCheckPointerOrMpi(void) { 60 if (PetscSegvJumpBuf_set) longjmp(PetscSegvJumpBuf, 1); 61 } 62 63 /*@C 64 PetscCheckPointer - Returns `PETSC_TRUE` if a pointer points to accessible data 65 66 Not Collective 67 68 Input Parameters: 69 + ptr - the pointer 70 - dtype - the type of data the pointer is suppose to point to 71 72 Level: developer 73 74 Note: 75 This is a non-standard PETSc function in that it returns the result as the return code and does not return an error code 76 77 .seealso: `PetscCheckPointerSetIntensity()` 78 @*/ 79 PetscBool PetscCheckPointer(const void *ptr, PetscDataType dtype) { 80 if (PETSC_RUNNING_ON_VALGRIND) return PETSC_TRUE; 81 if (!ptr) return PETSC_FALSE; 82 if (petsc_checkpointer_intensity < 1) return PETSC_TRUE; 83 84 #if PetscDefined(USE_DEBUG) 85 /* Skip the verbose check if we are inside a hot function. */ 86 if (petscstack.hotdepth > 0 && petsc_checkpointer_intensity < 2) return PETSC_TRUE; 87 #endif 88 89 PetscSegvJumpBuf_set = PETSC_TRUE; 90 91 if (setjmp(PetscSegvJumpBuf)) { 92 /* A segv was triggered in the code below hence we return with an error code */ 93 PetscSegvJumpBuf_set = PETSC_FALSE; 94 return PETSC_FALSE; 95 } else { 96 switch (dtype) { 97 case PETSC_INT: { 98 PETSC_UNUSED PetscInt x = (PetscInt) * (volatile PetscInt *)ptr; 99 break; 100 } 101 #if defined(PETSC_USE_COMPLEX) 102 case PETSC_SCALAR: { /* C++ is seriously dysfunctional with volatile std::complex. */ 103 #if defined(PETSC_USE_CXXCOMPLEX) 104 PetscReal xreal = ((volatile PetscReal *)ptr)[0], ximag = ((volatile PetscReal *)ptr)[1]; 105 PETSC_UNUSED volatile PetscScalar x = xreal + PETSC_i * ximag; 106 #else 107 PETSC_UNUSED PetscScalar x = *(volatile PetscScalar *)ptr; 108 #endif 109 break; 110 } 111 #endif 112 case PETSC_REAL: { 113 PETSC_UNUSED PetscReal x = *(volatile PetscReal *)ptr; 114 break; 115 } 116 case PETSC_BOOL: { 117 PETSC_UNUSED PetscBool x = *(volatile PetscBool *)ptr; 118 break; 119 } 120 case PETSC_ENUM: { 121 PETSC_UNUSED PetscEnum x = *(volatile PetscEnum *)ptr; 122 break; 123 } 124 case PETSC_CHAR: { 125 PETSC_UNUSED char x = *(volatile char *)ptr; 126 break; 127 } 128 case PETSC_OBJECT: { 129 PETSC_UNUSED volatile PetscClassId classid = ((PetscObject)ptr)->classid; 130 break; 131 } 132 default:; 133 } 134 } 135 PetscSegvJumpBuf_set = PETSC_FALSE; 136 return PETSC_TRUE; 137 } 138 139 #define PetscMPICUPMAwarnessCheckFunction \ 140 PetscBool PetscMPICUPMAwarenessCheck(void) { \ 141 cupmError_t cerr = cupmSuccess; \ 142 int ierr, hbuf[2] = {1, 0}, *dbuf = NULL; \ 143 PetscBool awareness = PETSC_FALSE; \ 144 cerr = cupmMalloc((void **)&dbuf, sizeof(int) * 2); \ 145 if (cerr != cupmSuccess) return PETSC_FALSE; \ 146 cerr = cupmMemcpy(dbuf, hbuf, sizeof(int) * 2, cupmMemcpyHostToDevice); \ 147 if (cerr != cupmSuccess) return PETSC_FALSE; \ 148 PetscSegvJumpBuf_set = PETSC_TRUE; \ 149 if (setjmp(PetscSegvJumpBuf)) { \ 150 /* If a segv was triggered in the MPI_Allreduce below, it is very likely due to the MPI is not GPU-aware */ \ 151 awareness = PETSC_FALSE; \ 152 } else { \ 153 ierr = MPI_Allreduce(dbuf, dbuf + 1, 1, MPI_INT, MPI_SUM, PETSC_COMM_SELF); \ 154 if (!ierr) awareness = PETSC_TRUE; \ 155 } \ 156 PetscSegvJumpBuf_set = PETSC_FALSE; \ 157 cerr = cupmFree(dbuf); \ 158 if (cerr != cupmSuccess) return PETSC_FALSE; \ 159 return awareness; \ 160 } 161 162 #if defined(PETSC_HAVE_CUDA) 163 #define cupmError_t cudaError_t 164 #define cupmMalloc cudaMalloc 165 #define cupmMemcpy cudaMemcpy 166 #define cupmFree cudaFree 167 #define cupmSuccess cudaSuccess 168 #define cupmMemcpyHostToDevice cudaMemcpyHostToDevice 169 #define PetscMPICUPMAwarenessCheck PetscMPICUDAAwarenessCheck 170 PetscMPICUPMAwarnessCheckFunction 171 #endif 172 173 #if defined(PETSC_HAVE_HIP) 174 #define cupmError_t hipError_t 175 #define cupmMalloc hipMalloc 176 #define cupmMemcpy hipMemcpy 177 #define cupmFree hipFree 178 #define cupmSuccess hipSuccess 179 #define cupmMemcpyHostToDevice hipMemcpyHostToDevice 180 #define PetscMPICUPMAwarenessCheck PetscMPIHIPAwarenessCheck 181 PetscMPICUPMAwarnessCheckFunction 182 #endif 183 184 #else 185 void PetscSignalSegvCheckPointerOrMpi(void) { 186 return; 187 } 188 189 PetscBool PetscCheckPointer(const void *ptr, PETSC_UNUSED PetscDataType dtype) { 190 if (!ptr) return PETSC_FALSE; 191 return PETSC_TRUE; 192 } 193 194 #if defined(PETSC_HAVE_CUDA) 195 PetscBool PetscMPICUDAAwarenessCheck(void) { 196 /* If no setjmp (rare), return true and let users code run (and segfault if they should) */ 197 return PETSC_TRUE; 198 } 199 #endif 200 201 #if defined(PETSC_HAVE_HIP) 202 PetscBool PetscMPIHIPAwarenessCheck(void) { 203 /* If no setjmp (rare), return true and let users code run (and segfault if they should) */ 204 return PETSC_TRUE; 205 } 206 #endif 207 208 #endif 209