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