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