1 #include <petsc/private/dmpleximpl.h> /*I "petscdmplex.h" I*/ 2 3 /*@ 4 DMPlexGetPointLocal - get location of point data in local `Vec` 5 6 Not Collective 7 8 Input Parameters: 9 + dm - `DM` defining the topological space 10 - point - topological point 11 12 Output Parameters: 13 + start - start of point data 14 - end - end of point data 15 16 Level: intermediate 17 18 Note: 19 This is a half open interval [start, end) 20 21 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMPlexGetPointLocalField()`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexPointLocalRead()`, `DMPlexPointLocalRead()`, `DMPlexPointLocalRef()` 22 @*/ 23 PetscErrorCode DMPlexGetPointLocal(DM dm, PetscInt point, PetscInt *start, PetscInt *end) 24 { 25 PetscInt s, e; 26 27 PetscFunctionBegin; 28 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 29 if (start) PetscValidIntPointer(start, 3); 30 if (end) PetscValidIntPointer(end, 4); 31 PetscCall(DMGetLocalOffset_Private(dm, point, &s, &e)); 32 if (start) *start = s; 33 if (end) *end = e; 34 PetscFunctionReturn(PETSC_SUCCESS); 35 } 36 37 /*@ 38 DMPlexPointLocalRead - return read access to a point in local array 39 40 Not Collective 41 42 Input Parameters: 43 + dm - `DM` defining topological space 44 . point - topological point 45 - array - array to index into 46 47 Output Parameter: 48 . ptr - address of read reference to point data, type generic so user can place in structure 49 50 Level: intermediate 51 52 Note: 53 A common usage when data sizes are known statically: 54 .vb 55 const struct { PetscScalar foo,bar,baz; } *ptr; 56 DMPlexPointLocalRead(dm,point,array,&ptr); 57 x = 2*ptr->foo + 3*ptr->bar + 5*ptr->baz; 58 .ve 59 60 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRead()` 61 @*/ 62 PetscErrorCode DMPlexPointLocalRead(DM dm, PetscInt point, const PetscScalar *array, void *ptr) 63 { 64 PetscInt start, end; 65 66 PetscFunctionBegin; 67 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 68 PetscValidScalarPointer(array, 3); 69 PetscValidPointer(ptr, 4); 70 PetscCall(DMGetLocalOffset_Private(dm, point, &start, &end)); 71 *(const PetscScalar **)ptr = (start < end) ? array + start : NULL; 72 PetscFunctionReturn(PETSC_SUCCESS); 73 } 74 75 /*@ 76 DMPlexPointLocalRef - return read/write access to a point in local array 77 78 Not Collective 79 80 Input Parameters: 81 + dm - `DM` defining topological space 82 . point - topological point 83 - array - array to index into 84 85 Output Parameter: 86 . ptr - address of reference to point data, type generic so user can place in structure 87 88 Level: intermediate 89 90 Note: 91 A common usage when data sizes are known statically: 92 .vb 93 struct { PetscScalar foo,bar,baz; } *ptr; 94 DMPlexPointLocalRef(dm,point,array,&ptr); 95 ptr->foo = 2; ptr->bar = 3; ptr->baz = 5; 96 .ve 97 98 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRef()` 99 @*/ 100 PetscErrorCode DMPlexPointLocalRef(DM dm, PetscInt point, PetscScalar *array, void *ptr) 101 { 102 PetscInt start, end; 103 104 PetscFunctionBegin; 105 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 106 PetscValidScalarPointer(array, 3); 107 PetscValidPointer(ptr, 4); 108 PetscCall(DMGetLocalOffset_Private(dm, point, &start, &end)); 109 *(PetscScalar **)ptr = (start < end) ? array + start : NULL; 110 PetscFunctionReturn(PETSC_SUCCESS); 111 } 112 113 /*@ 114 DMPlexGetPointLocalField - get location of point field data in local Vec 115 116 Not Collective 117 118 Input Parameters: 119 + dm - `DM` defining the topological space 120 . point - topological point 121 - field - the field number 122 123 Output Parameters: 124 + start - start of point data 125 - end - end of point data 126 127 Level: intermediate 128 129 Note: 130 This is a half open interval [start, end) 131 132 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMPlexGetPointLocal()`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexPointLocalRead()`, `DMPlexPointLocalRead()`, `DMPlexPointLocalRef()` 133 @*/ 134 PetscErrorCode DMPlexGetPointLocalField(DM dm, PetscInt point, PetscInt field, PetscInt *start, PetscInt *end) 135 { 136 PetscInt s, e; 137 138 PetscFunctionBegin; 139 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 140 if (start) PetscValidIntPointer(start, 4); 141 if (end) PetscValidIntPointer(end, 5); 142 PetscCall(DMGetLocalFieldOffset_Private(dm, point, field, &s, &e)); 143 if (start) *start = s; 144 if (end) *end = e; 145 PetscFunctionReturn(PETSC_SUCCESS); 146 } 147 148 /*@ 149 DMPlexPointLocalFieldRead - return read access to a field on a point in local array 150 151 Not Collective 152 153 Input Parameters: 154 + dm - `DM` defining topological space 155 . point - topological point 156 . field - field number 157 - array - array to index into 158 159 Output Parameter: 160 . ptr - address of read reference to point data, type generic so user can place in structure 161 162 Level: intermediate 163 164 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRef()` 165 @*/ 166 PetscErrorCode DMPlexPointLocalFieldRead(DM dm, PetscInt point, PetscInt field, const PetscScalar *array, void *ptr) 167 { 168 PetscInt start, end; 169 170 PetscFunctionBegin; 171 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 172 PetscValidScalarPointer(array, 4); 173 PetscValidPointer(ptr, 5); 174 PetscCall(DMGetLocalFieldOffset_Private(dm, point, field, &start, &end)); 175 *(const PetscScalar **)ptr = array + start; 176 PetscFunctionReturn(PETSC_SUCCESS); 177 } 178 179 /*@ 180 DMPlexPointLocalFieldRef - return read/write access to a field on a point in local array 181 182 Not Collective 183 184 Input Parameters: 185 + dm - `DM` defining topological space 186 . point - topological point 187 . field - field number 188 - array - array to index into 189 190 Output Parameter: 191 . ptr - address of reference to point data, type generic so user can place in structure 192 193 Level: intermediate 194 195 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRef()` 196 @*/ 197 PetscErrorCode DMPlexPointLocalFieldRef(DM dm, PetscInt point, PetscInt field, PetscScalar *array, void *ptr) 198 { 199 PetscInt start, end; 200 201 PetscFunctionBegin; 202 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 203 PetscValidScalarPointer(array, 4); 204 PetscValidPointer(ptr, 5); 205 PetscCall(DMGetLocalFieldOffset_Private(dm, point, field, &start, &end)); 206 *(PetscScalar **)ptr = array + start; 207 PetscFunctionReturn(PETSC_SUCCESS); 208 } 209 210 /*@ 211 DMPlexGetPointGlobal - get location of point data in global Vec 212 213 Not Collective 214 215 Input Parameters: 216 + dm - `DM` defining the topological space 217 - point - topological point 218 219 Output Parameters: 220 + start - start of point data; returns -(globalStart+1) if point is not owned 221 - end - end of point data; returns -(globalEnd+1) if point is not owned 222 223 Level: intermediate 224 225 Note: 226 This is a half open interval [start, end) 227 228 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMPlexGetPointGlobalField()`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexPointGlobalRead()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRead()`, `DMPlexPointGlobalRef()` 229 @*/ 230 PetscErrorCode DMPlexGetPointGlobal(DM dm, PetscInt point, PetscInt *start, PetscInt *end) 231 { 232 PetscInt s, e; 233 234 PetscFunctionBegin; 235 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 236 if (start) PetscValidIntPointer(start, 3); 237 if (end) PetscValidIntPointer(end, 4); 238 PetscCall(DMGetGlobalOffset_Private(dm, point, &s, &e)); 239 if (start) *start = s; 240 if (end) *end = e; 241 PetscFunctionReturn(PETSC_SUCCESS); 242 } 243 244 /*@ 245 DMPlexPointGlobalRead - return read access to a point in global array 246 247 Not Collective 248 249 Input Parameters: 250 + dm - `DM` defining topological space 251 . point - topological point 252 - array - array to index into 253 254 Output Parameter: 255 . ptr - address of read reference to point data, type generic so user can place in structure; returns NULL if global point is not owned 256 257 Level: intermediate 258 259 Note: 260 A common usage when data sizes are known statically: 261 .vb 262 const struct { PetscScalar foo,bar,baz; } *ptr; 263 DMPlexPointGlobalRead(dm,point,array,&ptr); 264 x = 2*ptr->foo + 3*ptr->bar + 5*ptr->baz; 265 .ve 266 267 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointGlobal()`, `DMPlexPointLocalRead()`, `DMPlexPointGlobalRef()` 268 @*/ 269 PetscErrorCode DMPlexPointGlobalRead(DM dm, PetscInt point, const PetscScalar *array, const void *ptr) 270 { 271 PetscInt start, end; 272 273 PetscFunctionBegin; 274 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 275 PetscValidScalarPointer(array, 3); 276 PetscValidPointer(ptr, 4); 277 PetscCall(DMGetGlobalOffset_Private(dm, point, &start, &end)); 278 *(const PetscScalar **)ptr = (start < end) ? array + start - dm->map->rstart : NULL; 279 PetscFunctionReturn(PETSC_SUCCESS); 280 } 281 282 /*@ 283 DMPlexPointGlobalRef - return read/write access to a point in global array 284 285 Not Collective 286 287 Input Parameters: 288 + dm - `DM` defining topological space 289 . point - topological point 290 - array - array to index into 291 292 Output Parameter: 293 . ptr - address of reference to point data, type generic so user can place in structure; returns NULL if global point is not owned 294 295 Level: intermediate 296 297 Note: 298 A common usage when data sizes are known statically: 299 .vb 300 struct { PetscScalar foo,bar,baz; } *ptr; 301 DMPlexPointGlobalRef(dm,point,array,&ptr); 302 ptr->foo = 2; ptr->bar = 3; ptr->baz = 5; 303 .ve 304 305 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointGlobal()`, `DMPlexPointLocalRef()`, `DMPlexPointGlobalRead()` 306 @*/ 307 PetscErrorCode DMPlexPointGlobalRef(DM dm, PetscInt point, PetscScalar *array, void *ptr) 308 { 309 PetscInt start, end; 310 311 PetscFunctionBegin; 312 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 313 PetscValidScalarPointer(array, 3); 314 PetscValidPointer(ptr, 4); 315 PetscCall(DMGetGlobalOffset_Private(dm, point, &start, &end)); 316 *(PetscScalar **)ptr = (start < end) ? array + start - dm->map->rstart : NULL; 317 PetscFunctionReturn(PETSC_SUCCESS); 318 } 319 320 /*@ 321 DMPlexGetPointGlobalField - get location of point field data in global `Vec` 322 323 Not Collective 324 325 Input Parameters: 326 + dm - `DM` defining the topological space 327 . point - topological point 328 - field - the field number 329 330 Output Parameters: 331 + start - start of point data; returns -(globalStart+1) if point is not owned 332 - end - end of point data; returns -(globalEnd+1) if point is not owned 333 334 Level: intermediate 335 336 Note: 337 This is a half open interval [start, end) 338 339 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMPlexGetPointGlobal()`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexPointGlobalRead()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRead()`, `DMPlexPointGlobalRef()` 340 @*/ 341 PetscErrorCode DMPlexGetPointGlobalField(DM dm, PetscInt point, PetscInt field, PetscInt *start, PetscInt *end) 342 { 343 PetscInt s, e; 344 345 PetscFunctionBegin; 346 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 347 if (start) PetscValidIntPointer(start, 4); 348 if (end) PetscValidIntPointer(end, 5); 349 PetscCall(DMGetGlobalFieldOffset_Private(dm, point, field, &s, &e)); 350 if (start) *start = s; 351 if (end) *end = e; 352 PetscFunctionReturn(PETSC_SUCCESS); 353 } 354 355 /*@ 356 DMPlexPointGlobalFieldRead - return read access to a field on a point in global array 357 358 Not Collective 359 360 Input Parameters: 361 + dm - `DM` defining topological space 362 . point - topological point 363 . field - field number 364 - array - array to index into 365 366 Output Parameter: 367 . ptr - address of read reference to point data, type generic so user can place in structure; returns NULL if global point is not owned 368 369 Level: intermediate 370 371 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointGlobal()`, `DMPlexPointLocalRead()`, `DMPlexPointGlobalRef()` 372 @*/ 373 PetscErrorCode DMPlexPointGlobalFieldRead(DM dm, PetscInt point, PetscInt field, const PetscScalar *array, void *ptr) 374 { 375 PetscInt start, end; 376 377 PetscFunctionBegin; 378 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 379 PetscValidScalarPointer(array, 4); 380 PetscValidPointer(ptr, 5); 381 PetscCall(DMGetGlobalFieldOffset_Private(dm, point, field, &start, &end)); 382 *(const PetscScalar **)ptr = (start < end) ? array + start - dm->map->rstart : NULL; 383 PetscFunctionReturn(PETSC_SUCCESS); 384 } 385 386 /*@ 387 DMPlexPointGlobalFieldRef - return read/write access to a field on a point in global array 388 389 Not Collective 390 391 Input Parameters: 392 + dm - `DM` defining topological space 393 . point - topological point 394 . field - field number 395 - array - array to index into 396 397 Output Parameter: 398 . ptr - address of reference to point data, type generic so user can place in structure; returns NULL if global point is not owned 399 400 Level: intermediate 401 402 .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointGlobal()`, `DMPlexPointLocalRef()`, `DMPlexPointGlobalRead()` 403 @*/ 404 PetscErrorCode DMPlexPointGlobalFieldRef(DM dm, PetscInt point, PetscInt field, PetscScalar *array, void *ptr) 405 { 406 PetscInt start, end; 407 408 PetscFunctionBegin; 409 PetscValidHeaderSpecific(dm, DM_CLASSID, 1); 410 PetscValidScalarPointer(array, 4); 411 PetscValidPointer(ptr, 5); 412 PetscCall(DMGetGlobalFieldOffset_Private(dm, point, field, &start, &end)); 413 *(PetscScalar **)ptr = (start < end) ? array + start - dm->map->rstart : NULL; 414 PetscFunctionReturn(PETSC_SUCCESS); 415 } 416