Actual source code: checkptr.c
petsc-3.14.6 2021-03-30
1: #include <petsc/private/petscimpl.h>
2: #include <petscvalgrind.h>
4: #if defined(PETSC_HAVE_CUDA)
5: #include <cuda_runtime.h>
6: #endif
8: #if defined(PETSC_HAVE_HIP)
9: #include <hip/hip_runtime.h>
10: #endif
12: static PetscInt petsc_checkpointer_intensity = 1;
14: /*@
16: confirm whether the address is valid. An intensity of 0 never uses signal handlers, 1 uses them when not in a "hot"
17: function, and intensity of 2 always uses a signal handler.
19: Not Collective
21: Input Arguments:
22: . intensity - how much to check pointers for validity
24: Options Database:
25: . -check_pointer_intensity - intensity (0, 1, or 2)
27: Level: advanced
30: @*/
32: {
35: switch (intensity) {
36: case 0:
37: case 1:
38: case 2:
39: petsc_checkpointer_intensity = intensity;
40: break;
41: default: SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Intensity %D not in 0,1,2",intensity);
42: }
43: return(0);
44: }
46: /* ---------------------------------------------------------------------------------------*/
48: #if defined(PETSC_HAVE_SETJMP_H)
49: #include <setjmp.h>
50: static jmp_buf PetscSegvJumpBuf;
51: static PetscBool PetscSegvJumpBuf_set;
53: /*@C
54: PetscSignalSegvCheckPointerOrMpi - To be called from a signal handler for SIGSEGV. If the signal was received
56: with no effect. This function is called automatically by PetscSignalHandlerDefault().
58: Not Collective
60: Level: developer
62: .seealso: PetscPushSignalHandler()
63: @*/
64: void PetscSignalSegvCheckPointerOrMpi(void) {
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: /* Skip the verbose check if we are inside a hot function. */
89: if (petscstack && petscstack->hotdepth > 0 && petsc_checkpointer_intensity < 2) return PETSC_TRUE;
91: PetscSegvJumpBuf_set = PETSC_TRUE;
93: if (setjmp(PetscSegvJumpBuf)) {
94: /* A segv was triggered in the code below hence we return with an error code */
95: PetscSegvJumpBuf_set = PETSC_FALSE;
96: return PETSC_FALSE;
97: } else {
98: switch (dtype) {
99: case PETSC_INT:{
100: PETSC_UNUSED PetscInt x = (PetscInt)*(volatile PetscInt*)ptr;
101: break;
102: }
103: #if defined(PETSC_USE_COMPLEX)
104: case PETSC_SCALAR:{ /* C++ is seriously dysfunctional with volatile std::complex. */
105: #if defined(PETSC_USE_CXXCOMPLEX)
106: PetscReal xreal = ((volatile PetscReal*)ptr)[0],ximag = ((volatile PetscReal*)ptr)[1];
107: PETSC_UNUSED volatile PetscScalar x = xreal + PETSC_i*ximag;
108: #else
109: PETSC_UNUSED PetscScalar x = *(volatile PetscScalar*)ptr;
110: #endif
111: break;
112: }
113: #endif
114: case PETSC_REAL:{
115: PETSC_UNUSED PetscReal x = *(volatile PetscReal*)ptr;
116: break;
117: }
118: case PETSC_BOOL:{
119: PETSC_UNUSED PetscBool x = *(volatile PetscBool*)ptr;
120: break;
121: }
122: case PETSC_ENUM:{
123: PETSC_UNUSED PetscEnum x = *(volatile PetscEnum*)ptr;
124: break;
125: }
126: case PETSC_CHAR:{
127: PETSC_UNUSED char x = *(volatile char*)ptr;
128: break;
129: }
130: case PETSC_OBJECT:{
131: PETSC_UNUSED volatile PetscClassId classid = ((PetscObject)ptr)->classid;
132: break;
133: }
134: default:;
135: }
136: }
137: PetscSegvJumpBuf_set = PETSC_FALSE;
138: return PETSC_TRUE;
139: }
141: #define PetscMPICUPMAwarnessCheckFunction \
142: PetscBool PetscMPICUPMAwarenessCheck(void) \
143: { \
144: cupmError_t cerr=cupmSuccess; \
145: int ierr,hbuf[2]={1,0},*dbuf=NULL; \
146: PetscBool awareness=PETSC_FALSE; \
147: cerr = cupmMalloc((void**)&dbuf,sizeof(int)*2);if (cerr != cupmSuccess) return PETSC_FALSE; \
148: cerr = cupmMemcpy(dbuf,hbuf,sizeof(int)*2,cupmMemcpyHostToDevice);if (cerr != cupmSuccess) return PETSC_FALSE; \
149: PetscSegvJumpBuf_set = PETSC_TRUE; \
150: if (setjmp(PetscSegvJumpBuf)) { \
151: /* If a segv was triggered in the MPI_Allreduce below, it is very likely due to the MPI is not GPU-aware */ \
152: awareness = PETSC_FALSE; \
153: } else { \
154: MPI_Allreduce(dbuf,dbuf+1,1,MPI_INT,MPI_SUM,PETSC_COMM_SELF); \
155: if (!ierr) awareness = PETSC_TRUE; \
156: } \
157: PetscSegvJumpBuf_set = PETSC_FALSE; \
158: cerr = cupmFree(dbuf);if (cerr != cupmSuccess) return PETSC_FALSE; \
159: return awareness; \
160: }
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
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
184: #else
185: void PetscSignalSegvCheckPointerOrMpi(void) {
186: return;
187: }
190: {
191: if (!ptr) return PETSC_FALSE;
192: return PETSC_TRUE;
193: }
195: #if defined (PETSC_HAVE_CUDA)
196: PetscBool PetscMPICUDAAwarenessCheck(void)
197: {
198: /* If no setjmp (rare), return true and let users code run (and segfault if they should) */
199: return PETSC_TRUE;
200: }
201: #endif
203: #if defined (PETSC_HAVE_HIP)
204: PetscBool PetscMPIHIPAwarenessCheck(void)
205: {
206: /* If no setjmp (rare), return true and let users code run (and segfault if they should) */
207: return PETSC_TRUE;
208: }
209: #endif
211: #endif