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 return (size_t)index / PETSC_BITS_PER_BYTE; 39 } 40 41 static inline char PetscBTMask_Internal(PetscInt index) { 42 return 1 << index % PETSC_BITS_PER_BYTE; 43 } 44 45 static inline size_t PetscBTLength(PetscInt m) { 46 return (size_t)m / PETSC_BITS_PER_BYTE + 1; 47 } 48 49 static inline PetscErrorCode PetscBTMemzero(PetscInt m, PetscBT array) { 50 return PetscArrayzero(array, PetscBTLength(m)); 51 } 52 53 static inline PetscErrorCode PetscBTDestroy(PetscBT *array) { 54 return (*array) ? PetscFree(*array) : 0; 55 } 56 57 static inline PetscErrorCode PetscBTCreate(PetscInt m, PetscBT *array) { 58 return PetscCalloc1(PetscBTLength(m), array); 59 } 60 61 static inline char PetscBTLookup(PetscBT array, PetscInt index) { 62 return array[PetscBTIndex_Internal(index)] & PetscBTMask_Internal(index); 63 } 64 65 static inline PetscErrorCode PetscBTSet(PetscBT array, PetscInt index) { 66 PetscFunctionBegin; 67 array[PetscBTIndex_Internal(index)] |= PetscBTMask_Internal(index); 68 PetscFunctionReturn(0); 69 } 70 71 static inline PetscErrorCode PetscBTNegate(PetscBT array, PetscInt index) { 72 PetscFunctionBegin; 73 array[PetscBTIndex_Internal(index)] ^= PetscBTMask_Internal(index); 74 PetscFunctionReturn(0); 75 } 76 77 static inline PetscErrorCode PetscBTClear(PetscBT array, PetscInt index) { 78 PetscFunctionBegin; 79 array[PetscBTIndex_Internal(index)] &= (char)~PetscBTMask_Internal(index); 80 PetscFunctionReturn(0); 81 } 82 83 static inline char PetscBTLookupSet(PetscBT array, PetscInt index) { 84 const char ret = PetscBTLookup(array, index); 85 PetscCallContinue(PetscBTSet(array, index)); 86 return ret; 87 } 88 89 static inline char PetscBTLookupClear(PetscBT array, PetscInt index) { 90 const char ret = PetscBTLookup(array, index); 91 PetscCallContinue(PetscBTClear(array, index)); 92 return ret; 93 } 94 95 static inline PetscErrorCode PetscBTView(PetscInt m, const PetscBT bt, PetscViewer viewer) { 96 PetscFunctionBegin; 97 if (m < 1) PetscFunctionReturn(0); 98 if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PETSC_COMM_SELF, &viewer)); 99 PetscCall(PetscViewerASCIIPushSynchronized(viewer)); 100 for (PetscInt i = 0; i < m; ++i) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "%" PetscInt_FMT " %d\n", i, (int)PetscBTLookup(bt, i))); 101 PetscCall(PetscViewerFlush(viewer)); 102 PetscCall(PetscViewerASCIIPopSynchronized(viewer)); 103 PetscFunctionReturn(0); 104 } 105 #endif /* PETSCBT_H */ 106