Actual source code: checkptr.c
1: #include <petsc/private/petscimpl.h>
3: #if defined(PETSC_HAVE_CUDA)
4: #include <cuda_runtime.h>
5: #endif
7: #if defined(PETSC_HAVE_HIP)
8: #include <hip/hip_runtime.h>
9: #endif
11: static PetscInt petsc_checkpointer_intensity = 1;
13: /*@
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.
18: Not Collective
20: Input Parameter:
21: . intensity - how much to check pointers for validity
23: Options Database:
24: . -check_pointer_intensity - intensity (0, 1, or 2)
26: Level: advanced
29: @*/
31: {
34: switch (intensity) {
35: case 0:
36: case 1:
37: case 2:
38: petsc_checkpointer_intensity = intensity;
39: break;
40: default: SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Intensity %D not in 0,1,2",intensity);
41: }
42: return(0);
43: }
45: /* ---------------------------------------------------------------------------------------*/
47: #if defined(PETSC_HAVE_SETJMP_H)
48: #include <setjmp.h>
49: static jmp_buf PetscSegvJumpBuf;
50: static PetscBool PetscSegvJumpBuf_set;
52: /*@C
53: PetscSignalSegvCheckPointerOrMpi - To be called from a signal handler for SIGSEGV. If the signal was received
55: with no effect. This function is called automatically by PetscSignalHandlerDefault().
57: Not Collective
59: Level: developer
61: .seealso: PetscPushSignalHandler()
62: @*/
63: void PetscSignalSegvCheckPointerOrMpi(void)
64: {
65: if (PetscSegvJumpBuf_set) longjmp(PetscSegvJumpBuf,1);
66: }
68: /*@C
71: Not Collective
73: Input Parameters:
74: + ptr - the pointer
75: - dtype - the type of data the pointer is suppose to point to
77: Level: developer
80: @*/
82: {
84: if (PETSC_RUNNING_ON_VALGRIND) return PETSC_TRUE;
85: if (!ptr) return PETSC_FALSE;
86: if (petsc_checkpointer_intensity < 1) return PETSC_TRUE;
88: #if PetscDefined(USE_DEBUG)
89: /* Skip the verbose check if we are inside a hot function. */
90: if (petscstack.hotdepth > 0 && petsc_checkpointer_intensity < 2) return PETSC_TRUE;
91: #endif
93: PetscSegvJumpBuf_set = PETSC_TRUE;
95: if (setjmp(PetscSegvJumpBuf)) {
96: /* A segv was triggered in the code below hence we return with an error code */
97: PetscSegvJumpBuf_set = PETSC_FALSE;
98: return PETSC_FALSE;
99: } else {
100: switch (dtype) {
101: case PETSC_INT:{
102: PETSC_UNUSED PetscInt x = (PetscInt)*(volatile PetscInt*)ptr;
103: break;
104: }
105: #if defined(PETSC_USE_COMPLEX)
106: case PETSC_SCALAR:{ /* C++ is seriously dysfunctional with volatile std::complex. */
107: #if defined(PETSC_USE_CXXCOMPLEX)
108: PetscReal xreal = ((volatile PetscReal*)ptr)[0],ximag = ((volatile PetscReal*)ptr)[1];
109: PETSC_UNUSED volatile PetscScalar x = xreal + PETSC_i*ximag;
110: #else
111: PETSC_UNUSED PetscScalar x = *(volatile PetscScalar*)ptr;
112: #endif
113: break;
114: }
115: #endif
116: case PETSC_REAL:{
117: PETSC_UNUSED PetscReal x = *(volatile PetscReal*)ptr;
118: break;
119: }
120: case PETSC_BOOL:{
121: PETSC_UNUSED PetscBool x = *(volatile PetscBool*)ptr;
122: break;
123: }
124: case PETSC_ENUM:{
125: PETSC_UNUSED PetscEnum x = *(volatile PetscEnum*)ptr;
126: break;
127: }
128: case PETSC_CHAR:{
129: PETSC_UNUSED char x = *(volatile char*)ptr;
130: break;
131: }
132: case PETSC_OBJECT:{
133: PETSC_UNUSED volatile PetscClassId classid = ((PetscObject)ptr)->classid;
134: break;
135: }
136: default:;
137: }
138: }
139: PetscSegvJumpBuf_set = PETSC_FALSE;
140: return PETSC_TRUE;
141: }
143: #define PetscMPICUPMAwarnessCheckFunction \
144: PetscBool PetscMPICUPMAwarenessCheck(void) \
145: { \
146: cupmError_t cerr=cupmSuccess; \
147: int ierr,hbuf[2]={1,0},*dbuf=NULL; \
148: PetscBool awareness=PETSC_FALSE; \
149: cerr = cupmMalloc((void**)&dbuf,sizeof(int)*2);if (cerr != cupmSuccess) return PETSC_FALSE; \
150: cerr = cupmMemcpy(dbuf,hbuf,sizeof(int)*2,cupmMemcpyHostToDevice);if (cerr != cupmSuccess) return PETSC_FALSE; \
151: PetscSegvJumpBuf_set = PETSC_TRUE; \
152: if (setjmp(PetscSegvJumpBuf)) { \
153: /* If a segv was triggered in the MPI_Allreduce below, it is very likely due to the MPI is not GPU-aware */ \
154: awareness = PETSC_FALSE; \
155: } else { \
156: MPI_Allreduce(dbuf,dbuf+1,1,MPI_INT,MPI_SUM,PETSC_COMM_SELF); \
157: if (!ierr) awareness = PETSC_TRUE; \
158: } \
159: PetscSegvJumpBuf_set = PETSC_FALSE; \
160: cerr = cupmFree(dbuf);if (cerr != cupmSuccess) return PETSC_FALSE; \
161: return awareness; \
162: }
164: #if defined(PETSC_HAVE_CUDA)
165: #define cupmError_t cudaError_t
166: #define cupmMalloc cudaMalloc
167: #define cupmMemcpy cudaMemcpy
168: #define cupmFree cudaFree
169: #define cupmSuccess cudaSuccess
170: #define cupmMemcpyHostToDevice cudaMemcpyHostToDevice
171: #define PetscMPICUPMAwarenessCheck PetscMPICUDAAwarenessCheck
172: PetscMPICUPMAwarnessCheckFunction
173: #endif
175: #if defined(PETSC_HAVE_HIP)
176: #define cupmError_t hipError_t
177: #define cupmMalloc hipMalloc
178: #define cupmMemcpy hipMemcpy
179: #define cupmFree hipFree
180: #define cupmSuccess hipSuccess
181: #define cupmMemcpyHostToDevice hipMemcpyHostToDevice
182: #define PetscMPICUPMAwarenessCheck PetscMPIHIPAwarenessCheck
183: PetscMPICUPMAwarnessCheckFunction
184: #endif
186: #else
187: void PetscSignalSegvCheckPointerOrMpi(void)
188: {
189: return;
190: }
193: {
194: if (!ptr) return PETSC_FALSE;
195: return PETSC_TRUE;
196: }
198: #if defined (PETSC_HAVE_CUDA)
199: PetscBool PetscMPICUDAAwarenessCheck(void)
200: {
201: /* If no setjmp (rare), return true and let users code run (and segfault if they should) */
202: return PETSC_TRUE;
203: }
204: #endif
206: #if defined (PETSC_HAVE_HIP)
207: PetscBool PetscMPIHIPAwarenessCheck(void)
208: {
209: /* If no setjmp (rare), return true and let users code run (and segfault if they should) */
210: return PETSC_TRUE;
211: }
212: #endif
214: #endif