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__ "DMPlexGetGlobalFieldOffset_Private" 71 PETSC_STATIC_INLINE PetscErrorCode DMPlexGetGlobalFieldOffset_Private(DM dm,PetscInt point,PetscInt field,PetscInt *offset) 72 { 73 PetscFunctionBegin; 74 #if defined(PETSC_USE_DEBUG) 75 { 76 PetscErrorCode ierr; 77 PetscInt loff,lfoff,dof,cdof; 78 ierr = PetscSectionGetOffset(dm->defaultGlobalSection,point,offset);CHKERRQ(ierr); 79 ierr = PetscSectionGetOffset(dm->defaultSection,point,&loff);CHKERRQ(ierr); 80 ierr = PetscSectionGetFieldOffset(dm->defaultSection,point,field,&lfoff);CHKERRQ(ierr); 81 ierr = PetscSectionGetFieldDof(dm->defaultSection,point,field,&dof);CHKERRQ(ierr); 82 ierr = PetscSectionGetFieldConstraintDof(dm->defaultSection,point,field,&cdof);CHKERRQ(ierr); 83 *offset += lfoff - loff; 84 if (dof-cdof <= 0) *offset = -1; /* Indicates no data */ 85 } 86 #else 87 { 88 PetscSection s = dm->defaultSection; 89 PetscSection fs = dm->defaultSection->field[field]; 90 PetscSection gs = dm->defaultGlobalSection; 91 PetscInt dof,cdof,loff,lfoff; 92 loff = s->atlasOff[point - s->pStart]; 93 lfoff = fs->atlasOff[point - s->pStart]; 94 dof = s->atlasDof[point - s->pStart]; 95 cdof = s->bc ? s->bc->atlasDof[point - s->bc->pStart] : 0; 96 *offset = gs->atlasOff[point - s->pStart] + lfoff - loff; 97 if (dof-cdof <= 0) *offset = -1; 98 } 99 #endif 100 PetscFunctionReturn(0); 101 } 102 103 #undef __FUNCT__ 104 #define __FUNCT__ "DMPlexGetPointLocal" 105 /*@ 106 DMPlexGetPointLocal - get location of point data in local Vec 107 108 Not Collective 109 110 Input Arguments: 111 + dm - DM defining the topological space 112 - point - topological point 113 114 Output Arguments: 115 + start - start of point data 116 - end - end of point data 117 118 Level: intermediate 119 120 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexPointLocalRead(), DMPlexPointLocalRead(), DMPlexPointLocalRef() 121 @*/ 122 PetscErrorCode DMPlexGetPointLocal(DM dm,PetscInt point,PetscInt *start,PetscInt *end) 123 { 124 PetscErrorCode ierr; 125 PetscInt offset,dof; 126 127 PetscFunctionBegin; 128 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 129 ierr = PetscSectionGetOffset(dm->defaultSection,point,&offset);CHKERRQ(ierr); 130 ierr = PetscSectionGetDof(dm->defaultSection,point,&dof);CHKERRQ(ierr); 131 if (start) *start = offset; 132 if (end) *end = offset + dof; 133 PetscFunctionReturn(0); 134 } 135 136 #undef __FUNCT__ 137 #define __FUNCT__ "DMPlexPointLocalRead" 138 /*@ 139 DMPlexPointLocalRead - return read access to a point in local array 140 141 Not Collective 142 143 Input Arguments: 144 + dm - DM defining topological space 145 . point - topological point 146 - array - array to index into 147 148 Output Arguments: 149 . ptr - address of read reference to point data, type generic so user can place in structure 150 151 Level: intermediate 152 153 Note: 154 A common usage when data sizes are known statically: 155 156 $ const struct { PetscScalar foo,bar,baz; } *ptr; 157 $ DMPlexPointLocalRead(dm,point,array,&ptr); 158 $ x = 2*ptr->foo + 3*ptr->bar + 5*ptr->baz; 159 160 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointLocal(), DMPlexPointGlobalRead() 161 @*/ 162 PetscErrorCode DMPlexPointLocalRead(DM dm,PetscInt point,const PetscScalar *array,const void *ptr) 163 { 164 PetscErrorCode ierr; 165 PetscInt start; 166 167 PetscFunctionBegin; 168 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 169 PetscValidScalarPointer(array,3); 170 PetscValidPointer(ptr,4); 171 ierr = DMPlexGetLocalOffset_Private(dm,point,&start);CHKERRQ(ierr); 172 *(const PetscScalar**)ptr = array + start; 173 PetscFunctionReturn(0); 174 } 175 176 #undef __FUNCT__ 177 #define __FUNCT__ "DMPlexPointLocalRef" 178 /*@ 179 DMPlexPointLocalRef - return read/write access to a point in local array 180 181 Not Collective 182 183 Input Arguments: 184 + dm - DM defining topological space 185 . point - topological point 186 - array - array to index into 187 188 Output Arguments: 189 . ptr - address of reference to point data, type generic so user can place in structure 190 191 Level: intermediate 192 193 Note: 194 A common usage when data sizes are known statically: 195 196 $ struct { PetscScalar foo,bar,baz; } *ptr; 197 $ DMPlexPointLocalRef(dm,point,array,&ptr); 198 $ ptr->foo = 2; ptr->bar = 3; ptr->baz = 5; 199 200 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointLocal(), DMPlexPointGlobalRef() 201 @*/ 202 PetscErrorCode DMPlexPointLocalRef(DM dm,PetscInt point,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 = DMPlexGetLocalOffset_Private(dm,point,&start);CHKERRQ(ierr); 212 *(PetscScalar**)ptr = array + start; 213 PetscFunctionReturn(0); 214 } 215 216 #undef __FUNCT__ 217 #define __FUNCT__ "DMPlexPointLocalFieldRead" 218 /*@ 219 DMPlexPointLocalFieldRead - return read access to a field on a point in local array 220 221 Not Collective 222 223 Input Arguments: 224 + dm - DM defining topological space 225 . point - topological point 226 . field - field number 227 - array - array to index into 228 229 Output Arguments: 230 . ptr - address of read reference to point data, type generic so user can place in structure 231 232 Level: intermediate 233 234 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointLocal(), DMPlexPointGlobalRef() 235 @*/ 236 PetscErrorCode DMPlexPointLocalFieldRead(DM dm,PetscInt point,PetscInt field,const PetscScalar *array,const void *ptr) 237 { 238 PetscErrorCode ierr; 239 PetscInt start; 240 241 PetscFunctionBegin; 242 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 243 PetscValidScalarPointer(array,3); 244 PetscValidPointer(ptr,4); 245 ierr = DMPlexGetLocalFieldOffset_Private(dm,point,field,&start);CHKERRQ(ierr); 246 *(const PetscScalar**)ptr = array + start; 247 PetscFunctionReturn(0); 248 } 249 250 #undef __FUNCT__ 251 #define __FUNCT__ "DMPlexPointLocalFieldRef" 252 /*@ 253 DMPlexPointLocalFieldRef - return read/write access to a field on a point in local array 254 255 Not Collective 256 257 Input Arguments: 258 + dm - DM defining topological space 259 . point - topological point 260 . field - field number 261 - array - array to index into 262 263 Output Arguments: 264 . ptr - address of reference to point data, type generic so user can place in structure 265 266 Level: intermediate 267 268 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointLocal(), DMPlexPointGlobalRef() 269 @*/ 270 PetscErrorCode DMPlexPointLocalFieldRef(DM dm,PetscInt point,PetscInt field,PetscScalar *array,void *ptr) 271 { 272 PetscErrorCode ierr; 273 PetscInt start; 274 275 PetscFunctionBegin; 276 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 277 PetscValidScalarPointer(array,3); 278 PetscValidPointer(ptr,4); 279 ierr = DMPlexGetLocalFieldOffset_Private(dm,point,field,&start);CHKERRQ(ierr); 280 *(PetscScalar**)ptr = array + start; 281 PetscFunctionReturn(0); 282 } 283 284 #undef __FUNCT__ 285 #define __FUNCT__ "DMPlexGetPointGlobal" 286 /*@ 287 DMPlexGetPointGlobal - get location of point data in global Vec 288 289 Not Collective 290 291 Input Arguments: 292 + dm - DM defining the topological space 293 - point - topological point 294 295 Output Arguments: 296 + start - start of point data; returns -(global_start+1) if point is not owned 297 - end - end of point data; returns -(global_end+1) if point is not owned 298 299 Level: intermediate 300 301 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexPointGlobalRead(), DMPlexGetPointLocal(), DMPlexPointGlobalRead(), DMPlexPointGlobalRef() 302 @*/ 303 PetscErrorCode DMPlexGetPointGlobal(DM dm,PetscInt point,PetscInt *start,PetscInt *end) 304 { 305 PetscErrorCode ierr; 306 PetscInt offset,dof; 307 308 PetscFunctionBegin; 309 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 310 ierr = PetscSectionGetOffset(dm->defaultGlobalSection,point,&offset);CHKERRQ(ierr); 311 ierr = PetscSectionGetDof(dm->defaultGlobalSection,point,&dof);CHKERRQ(ierr); 312 if (start) *start = offset; 313 if (end) *end = offset + dof; 314 PetscFunctionReturn(0); 315 } 316 317 #undef __FUNCT__ 318 #define __FUNCT__ "DMPlexPointGlobalRead" 319 /*@ 320 DMPlexPointGlobalRead - return read access to a point in global array 321 322 Not Collective 323 324 Input Arguments: 325 + dm - DM defining topological space 326 . point - topological point 327 - array - array to index into 328 329 Output Arguments: 330 . ptr - address of read reference to point data, type generic so user can place in structure; returns NULL if global point is not owned 331 332 Level: intermediate 333 334 Note: 335 A common usage when data sizes are known statically: 336 337 $ const struct { PetscScalar foo,bar,baz; } *ptr; 338 $ DMPlexPointGlobalRead(dm,point,array,&ptr); 339 $ x = 2*ptr->foo + 3*ptr->bar + 5*ptr->baz; 340 341 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointGlobal(), DMPlexPointLocalRead(), DMPlexPointGlobalRef() 342 @*/ 343 PetscErrorCode DMPlexPointGlobalRead(DM dm,PetscInt point,const PetscScalar *array,const void *ptr) 344 { 345 PetscErrorCode ierr; 346 PetscInt start; 347 348 PetscFunctionBegin; 349 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 350 PetscValidScalarPointer(array,3); 351 PetscValidPointer(ptr,4); 352 ierr = DMPlexGetGlobalOffset_Private(dm,point,&start);CHKERRQ(ierr); 353 *(const PetscScalar**)ptr = (start >= 0) ? array + start - dm->map->rstart : NULL; 354 PetscFunctionReturn(0); 355 } 356 357 #undef __FUNCT__ 358 #define __FUNCT__ "DMPlexPointGlobalRef" 359 /*@ 360 DMPlexPointGlobalRef - return read/write access to a point in global array 361 362 Not Collective 363 364 Input Arguments: 365 + dm - DM defining topological space 366 . point - topological point 367 - array - array to index into 368 369 Output Arguments: 370 . ptr - address of reference to point data, type generic so user can place in structure; returns NULL if global point is not owned 371 372 Level: intermediate 373 374 Note: 375 A common usage when data sizes are known statically: 376 377 $ struct { PetscScalar foo,bar,baz; } *ptr; 378 $ DMPlexPointGlobalRef(dm,point,array,&ptr); 379 $ ptr->foo = 2; ptr->bar = 3; ptr->baz = 5; 380 381 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointGlobal(), DMPlexPointLocalRef(), DMPlexPointGlobalRead() 382 @*/ 383 PetscErrorCode DMPlexPointGlobalRef(DM dm,PetscInt point,PetscScalar *array,void *ptr) 384 { 385 PetscErrorCode ierr; 386 PetscInt start; 387 388 PetscFunctionBegin; 389 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 390 PetscValidScalarPointer(array,3); 391 PetscValidPointer(ptr,4); 392 ierr = DMPlexGetGlobalOffset_Private(dm,point,&start);CHKERRQ(ierr); 393 *(PetscScalar**)ptr = (start >= 0) ? array + start - dm->map->rstart : NULL; 394 PetscFunctionReturn(0); 395 } 396 397 #undef __FUNCT__ 398 #define __FUNCT__ "DMPlexPointGlobalFieldRead" 399 /*@ 400 DMPlexPointGlobalFieldRead - return read access to a field on a point in global array 401 402 Not Collective 403 404 Input Arguments: 405 + dm - DM defining topological space 406 . point - topological point 407 . field - field number 408 - array - array to index into 409 410 Output Arguments: 411 . ptr - address of read reference to point data, type generic so user can place in structure; returns NULL if global point is not owned 412 413 Level: intermediate 414 415 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointGlobal(), DMPlexPointLocalRead(), DMPlexPointGlobalRef() 416 @*/ 417 PetscErrorCode DMPlexPointGlobalFieldRead(DM dm,PetscInt point,PetscInt field,const PetscScalar *array,const void *ptr) 418 { 419 PetscErrorCode ierr; 420 PetscInt start; 421 422 PetscFunctionBegin; 423 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 424 PetscValidScalarPointer(array,3); 425 PetscValidPointer(ptr,4); 426 ierr = DMPlexGetGlobalFieldOffset_Private(dm,point,field,&start);CHKERRQ(ierr); 427 *(const PetscScalar**)ptr = (start >= 0) ? array + start - dm->map->rstart : NULL; 428 PetscFunctionReturn(0); 429 } 430 431 #undef __FUNCT__ 432 #define __FUNCT__ "DMPlexPointGlobalFieldRef" 433 /*@ 434 DMPlexPointGlobalFieldRef - return read/write access to a field on a point in global array 435 436 Not Collective 437 438 Input Arguments: 439 + dm - DM defining topological space 440 . point - topological point 441 . field - field number 442 - array - array to index into 443 444 Output Arguments: 445 . ptr - address of reference to point data, type generic so user can place in structure; returns NULL if global point is not owned 446 447 Level: intermediate 448 449 .seealso: DMGetDefaultSection(), PetscSectionGetOffset(), PetscSectionGetDof(), DMPlexGetPointGlobal(), DMPlexPointLocalRef(), DMPlexPointGlobalRead() 450 @*/ 451 PetscErrorCode DMPlexPointGlobalFieldRef(DM dm,PetscInt point,PetscInt field,PetscScalar *array,void *ptr) 452 { 453 PetscErrorCode ierr; 454 PetscInt start; 455 456 PetscFunctionBegin; 457 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 458 PetscValidScalarPointer(array,3); 459 PetscValidPointer(ptr,4); 460 ierr = DMPlexGetGlobalFieldOffset_Private(dm,point,field,&start);CHKERRQ(ierr); 461 *(PetscScalar**)ptr = (start >= 0) ? array + start - dm->map->rstart : NULL; 462 PetscFunctionReturn(0); 463 } 464