1 #ifndef PETSCBT_H 2 #define PETSCBT_H 3 4 #include <petscviewer.h> 5 6 /* SUBMANSEC = Sys */ 7 8 /*S 9 PetscBT - PETSc bitarrays, efficient storage of arrays of boolean values 10 11 Level: advanced 12 13 Notes: 14 The following routines do not have their own manual pages 15 16 .vb 17 PetscBTCreate(m,&bt) - creates a bit array with enough room to hold m values 18 PetscBTDestroy(&bt) - destroys the bit array 19 PetscBTMemzero(m,bt) - zeros the entire bit array (sets all values to false) 20 PetscBTSet(bt,index) - sets a particular entry as true 21 PetscBTClear(bt,index) - sets a particular entry as false 22 PetscBTLookup(bt,index) - returns the value 23 PetscBTLookupSet(bt,index) - returns the value and then sets it true 24 PetscBTLookupClear(bt,index) - returns the value and then sets it false 25 PetscBTLength(m) - returns number of bytes in array with m bits 26 PetscBTView(m,bt,viewer) - prints all the entries in a bit array 27 .ve 28 29 PETSc does not check error flags on `PetscBTLookup()`, `PetcBTLookupSet()`, `PetscBTLength()` because error checking 30 would cost hundreds more cycles then the operation. 31 32 S*/ 33 typedef char *PetscBT; 34 35 /* convert an index i to an index suitable for indexing a PetscBT, such that 36 * bt[PetscBTIndex(i)] returns the i'th value of the bt */ 37 static inline size_t PetscBTIndex_Internal(PetscInt index) 38 { 39 return (size_t)index / PETSC_BITS_PER_BYTE; 40 } 41 42 static inline char PetscBTMask_Internal(PetscInt index) 43 { 44 return (char)(1 << index % PETSC_BITS_PER_BYTE); 45 } 46 47 static inline size_t PetscBTLength(PetscInt m) 48 { 49 return (size_t)m / PETSC_BITS_PER_BYTE + 1; 50 } 51 52 static inline PetscErrorCode PetscBTMemzero(PetscInt m, PetscBT array) 53 { 54 return PetscArrayzero(array, PetscBTLength(m)); 55 } 56 57 static inline PetscErrorCode PetscBTDestroy(PetscBT *array) 58 { 59 return (*array) ? PetscFree(*array) : PETSC_SUCCESS; 60 } 61 62 static inline PetscErrorCode PetscBTCreate(PetscInt m, PetscBT *array) 63 { 64 return PetscCalloc1(PetscBTLength(m), array); 65 } 66 67 static inline char PetscBTLookup(PetscBT array, PetscInt index) 68 { 69 return array[PetscBTIndex_Internal(index)] & PetscBTMask_Internal(index); 70 } 71 72 static inline PetscErrorCode PetscBTSet(PetscBT array, PetscInt index) 73 { 74 PetscFunctionBegin; 75 array[PetscBTIndex_Internal(index)] |= PetscBTMask_Internal(index); 76 PetscFunctionReturn(PETSC_SUCCESS); 77 } 78 79 static inline PetscErrorCode PetscBTNegate(PetscBT array, PetscInt index) 80 { 81 PetscFunctionBegin; 82 array[PetscBTIndex_Internal(index)] ^= PetscBTMask_Internal(index); 83 PetscFunctionReturn(PETSC_SUCCESS); 84 } 85 86 static inline PetscErrorCode PetscBTClear(PetscBT array, PetscInt index) 87 { 88 PetscFunctionBegin; 89 array[PetscBTIndex_Internal(index)] &= (char)~PetscBTMask_Internal(index); 90 PetscFunctionReturn(PETSC_SUCCESS); 91 } 92 93 static inline char PetscBTLookupSet(PetscBT array, PetscInt index) 94 { 95 const char ret = PetscBTLookup(array, index); 96 PetscCallContinue(PetscBTSet(array, index)); 97 return ret; 98 } 99 100 static inline char PetscBTLookupClear(PetscBT array, PetscInt index) 101 { 102 const char ret = PetscBTLookup(array, index); 103 PetscCallContinue(PetscBTClear(array, index)); 104 return ret; 105 } 106 107 static inline PetscErrorCode PetscBTView(PetscInt m, const PetscBT bt, PetscViewer viewer) 108 { 109 PetscFunctionBegin; 110 if (m < 1) PetscFunctionReturn(PETSC_SUCCESS); 111 if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PETSC_COMM_SELF, &viewer)); 112 PetscCall(PetscViewerASCIIPushSynchronized(viewer)); 113 for (PetscInt i = 0; i < m; ++i) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "%" PetscInt_FMT " %d\n", i, (int)PetscBTLookup(bt, i))); 114 PetscCall(PetscViewerFlush(viewer)); 115 PetscCall(PetscViewerASCIIPopSynchronized(viewer)); 116 PetscFunctionReturn(PETSC_SUCCESS); 117 } 118 #endif /* PETSCBT_H */ 119