1 #include <petsc-private/pleximpl.h> /*I "petscdmplex.h" I*/ 2 3 #undef __FUNCT__ 4 #define __FUNCT__ "DMPlexGetLocalOffset_Private" 5 PETSC_STATIC_INLINE PetscErrorCode DMPlexGetLocalOffset_Private(DM dm,PetscInt point,PetscInt *offset) 6 { 7 PetscFunctionBegin; 8 #if defined(PETSC_USE_DEBUG) 9 { 10 PetscErrorCode ierr; 11 ierr = PetscSectionGetOffset(dm->defaultSection,point,offset);CHKERRQ(ierr); 12 } 13 #else 14 { 15 PetscSection s = dm->defaultSection; 16 *offset = s->atlasOff[point - s->atlasLayout.pStart]; 17 } 18 #endif 19 PetscFunctionReturn(0); 20 } 21 22 #undef __FUNCT__ 23 #define __FUNCT__ "DMPlexGetGlobalOffset_Private" 24 PETSC_STATIC_INLINE PetscErrorCode DMPlexGetGlobalOffset_Private(DM dm,PetscInt point,PetscInt *offset) 25 { 26 PetscFunctionBegin; 27 #if defined(PETSC_USE_DEBUG) 28 { 29 PetscErrorCode ierr; 30 PetscInt dof,cdof; 31 ierr = PetscSectionGetOffset(dm->defaultGlobalSection,point,offset);CHKERRQ(ierr); 32 ierr = PetscSectionGetDof(dm->defaultGlobalSection,point,&dof);CHKERRQ(ierr); 33 ierr = PetscSectionGetConstraintDof(dm->defaultGlobalSection,point,&cdof);CHKERRQ(ierr); 34 if (dof-cdof <= 0) *offset = -1; /* Indicates no data */ 35 } 36 #else 37 { 38 PetscSection s = dm->defaultGlobalSection; 39 PetscInt dof,cdof; 40 *offset = s->atlasOff[point - s->atlasLayout.pStart]; 41 dof = s->atlasDof[point - s->atlasLayout.pStart]; 42 cdof = s->bc ? s->bc->atlasDof[point - s->bc->atlasLayout.pStart] : 0; 43 if (dof-cdof <= 0) *offset = -1; 44 } 45 #endif 46 PetscFunctionReturn(0); 47 } 48 49 #undef __FUNCT__ 50 #define __FUNCT__ "DMPlexGetPointLocal" 51 /*@ 52 DMPlexGetPointLocal - get location of point data in local Vec 53 54 Not Collective 55 56 Input Arguments: 57 + dm - DM defining the topological space 58 - point - topological point 59 60 Output Arguments: 61 + start - start of point data 62 - end - end of point data 63 64 Level: intermediate 65 66 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexPointLocalRead(), DMPlexPointLocalRead(), DMPlexPointLocalRef() 67 @*/ 68 PetscErrorCode DMPlexGetPointLocal(DM dm,PetscInt point,PetscInt *start,PetscInt *end) 69 { 70 PetscErrorCode ierr; 71 PetscInt offset,dof; 72 73 PetscFunctionBegin; 74 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 75 ierr = PetscSectionGetOffset(dm->defaultSection,point,&offset);CHKERRQ(ierr); 76 ierr = PetscSectionGetDof(dm->defaultSection,point,&dof);CHKERRQ(ierr); 77 if (start) *start = offset; 78 if (end) *end = offset + dof; 79 PetscFunctionReturn(0); 80 } 81 82 #undef __FUNCT__ 83 #define __FUNCT__ "DMPlexPointLocalRead" 84 /*@ 85 DMPlexPointLocalRead - return read access to a point in local array 86 87 Not Collective 88 89 Input Arguments: 90 + dm - DM defining topological space 91 . point - topological point 92 - array - array to index into 93 94 Output Arguments: 95 . ptr - address of read reference to point data, type generic so user can place in structure 96 97 Level: intermediate 98 99 Note: 100 A common usage when data sizes are known statically: 101 102 $ const struct { PetscScalar foo,bar,baz; } *ptr; 103 $ DMPlexPointLocalRead(dm,point,array,&ptr); 104 $ x = 2*ptr->foo + 3*ptr->bar + 5*ptr->baz; 105 106 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointLocal(), DMPlexPointGlobalRead() 107 @*/ 108 PetscErrorCode DMPlexPointLocalRead(DM dm,PetscInt point,const PetscScalar *array,const void *ptr) 109 { 110 PetscErrorCode ierr; 111 PetscInt start; 112 113 PetscFunctionBegin; 114 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 115 PetscValidScalarPointer(array,3); 116 PetscValidPointer(ptr,4); 117 ierr = DMPlexGetLocalOffset_Private(dm,point,&start);CHKERRQ(ierr); 118 *(const PetscScalar **)ptr = array + start; 119 PetscFunctionReturn(0); 120 } 121 122 #undef __FUNCT__ 123 #define __FUNCT__ "DMPlexPointLocalRef" 124 /*@ 125 DMPlexPointLocalRef - return read/write access to a point in local array 126 127 Not Collective 128 129 Input Arguments: 130 + dm - DM defining topological space 131 . point - topological point 132 - array - array to index into 133 134 Output Arguments: 135 . ptr - address of reference to point data, type generic so user can place in structure 136 137 Level: intermediate 138 139 Note: 140 A common usage when data sizes are known statically: 141 142 $ struct { PetscScalar foo,bar,baz; } *ptr; 143 $ DMPlexPointLocalRef(dm,point,array,&ptr); 144 $ ptr->foo = 2; ptr->bar = 3; ptr->baz = 5; 145 146 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointLocal(), DMPlexPointGlobalRef() 147 @*/ 148 PetscErrorCode DMPlexPointLocalRef(DM dm,PetscInt point,PetscScalar *array,void *ptr) 149 { 150 PetscErrorCode ierr; 151 PetscInt start; 152 153 PetscFunctionBegin; 154 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 155 PetscValidScalarPointer(array,3); 156 PetscValidPointer(ptr,4); 157 ierr = DMPlexGetLocalOffset_Private(dm,point,&start);CHKERRQ(ierr); 158 *(PetscScalar **)ptr = array + start; 159 PetscFunctionReturn(0); 160 } 161 162 #undef __FUNCT__ 163 #define __FUNCT__ "DMPlexGetPointGlobal" 164 /*@ 165 DMPlexGetPointGlobal - get location of point data in global Vec 166 167 Not Collective 168 169 Input Arguments: 170 + dm - DM defining the topological space 171 - point - topological point 172 173 Output Arguments: 174 + start - start of point data; returns -(global_start+1) if point is not owned 175 - end - end of point data; returns -(global_end+1) if point is not owned 176 177 Level: intermediate 178 179 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexPointGlobalRead(), DMPlexGetPointLocal(), DMPlexPointGlobalRead(), DMPlexPointGlobalRef() 180 @*/ 181 PetscErrorCode DMPlexGetPointGlobal(DM dm,PetscInt point,PetscInt *start,PetscInt *end) 182 { 183 PetscErrorCode ierr; 184 PetscInt offset,dof; 185 186 PetscFunctionBegin; 187 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 188 ierr = PetscSectionGetOffset(dm->defaultGlobalSection,point,&offset);CHKERRQ(ierr); 189 ierr = PetscSectionGetDof(dm->defaultGlobalSection,point,&dof);CHKERRQ(ierr); 190 if (start) *start = offset; 191 if (end) *end = offset + dof; 192 PetscFunctionReturn(0); 193 } 194 195 #undef __FUNCT__ 196 #define __FUNCT__ "DMPlexPointGlobalRead" 197 /*@ 198 DMPlexPointGlobalRead - return read access to a point in global array 199 200 Not Collective 201 202 Input Arguments: 203 + dm - DM defining topological space 204 . point - topological point 205 - array - array to index into 206 207 Output Arguments: 208 . ptr - address of read reference to point data, type generic so user can place in structure; returns PETSC_NULL if global point is not owned 209 210 Level: intermediate 211 212 Note: 213 A common usage when data sizes are known statically: 214 215 $ const struct { PetscScalar foo,bar,baz; } *ptr; 216 $ DMPlexPointGlobalRead(dm,point,array,&ptr); 217 $ x = 2*ptr->foo + 3*ptr->bar + 5*ptr->baz; 218 219 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointGlobal(), DMPlexPointLocalRead(), DMPlexPointGlobalRef() 220 @*/ 221 PetscErrorCode DMPlexPointGlobalRead(DM dm,PetscInt point,const PetscScalar *array,const void *ptr) 222 { 223 PetscErrorCode ierr; 224 PetscInt start; 225 226 PetscFunctionBegin; 227 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 228 PetscValidScalarPointer(array,3); 229 PetscValidPointer(ptr,4); 230 ierr = DMPlexGetGlobalOffset_Private(dm,point,&start);CHKERRQ(ierr); 231 *(const PetscScalar **)ptr = (start >= 0) ? array + start - dm->map->rstart : PETSC_NULL; 232 PetscFunctionReturn(0); 233 } 234 235 #undef __FUNCT__ 236 #define __FUNCT__ "DMPlexPointGlobalRef" 237 /*@ 238 DMPlexPointGlobalRef - return read/write access to a point in global array 239 240 Not Collective 241 242 Input Arguments: 243 + dm - DM defining topological space 244 . point - topological point 245 - array - array to index into 246 247 Output Arguments: 248 . ptr - address of reference to point data, type generic so user can place in structure; returns PETSC_NULL if global point is not owned 249 250 Level: intermediate 251 252 Note: 253 A common usage when data sizes are known statically: 254 255 $ struct { PetscScalar foo,bar,baz; } *ptr; 256 $ DMPlexPointGlobalRef(dm,point,array,&ptr); 257 $ ptr->foo = 2; ptr->bar = 3; ptr->baz = 5; 258 259 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointGlobal(), DMPlexPointLocalRef(), DMPlexPointGlobalRead() 260 @*/ 261 PetscErrorCode DMPlexPointGlobalRef(DM dm,PetscInt point,PetscScalar *array,void *ptr) 262 { 263 PetscErrorCode ierr; 264 PetscInt start; 265 266 PetscFunctionBegin; 267 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 268 PetscValidScalarPointer(array,3); 269 PetscValidPointer(ptr,4); 270 ierr = DMPlexGetGlobalOffset_Private(dm,point,&start);CHKERRQ(ierr); 271 *(PetscScalar **)ptr = (start >= 0) ? array + start - dm->map->rstart : PETSC_NULL; 272 PetscFunctionReturn(0); 273 } 274