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