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