Actual source code: plexpoint.c
  1: #include <petsc/private/dmpleximpl.h>
  3: /*@
  4:   DMPlexGetPointLocal - get location of point data in local `Vec`
  6:   Not Collective
  8:   Input Parameters:
  9: + dm    - `DM` defining the topological space
 10: - point - topological point
 12:   Output Parameters:
 13: + start - start of point data
 14: - end   - end of point data
 16:   Level: intermediate
 18:   Note:
 19:   This is a half open interval [start, end)
 21: .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMPlexGetPointLocalField()`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexPointLocalRead()`, `DMPlexPointLocalRef()`
 22: @*/
 23: PetscErrorCode DMPlexGetPointLocal(DM dm, PetscInt point, PetscInt *start, PetscInt *end)
 24: {
 25:   PetscInt s, e;
 27:   PetscFunctionBegin;
 29:   if (start) PetscAssertPointer(start, 3);
 30:   if (end) PetscAssertPointer(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: }
 37: /*@
 38:   DMPlexPointLocalRead - return read access to a point in local array
 40:   Not Collective
 42:   Input Parameters:
 43: + dm    - `DM` defining topological space
 44: . point - topological point
 45: - array - array to index into
 47:   Output Parameter:
 48: . ptr - address of read reference to point data, type generic so user can place in structure
 50:   Level: intermediate
 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
 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;
 66:   PetscFunctionBegin;
 68:   PetscAssertPointer(array, 3);
 69:   PetscAssertPointer(ptr, 4);
 70:   PetscCall(DMGetLocalOffset_Private(dm, point, &start, &end));
 71:   *(const PetscScalar **)ptr = (start < end) ? array + start : NULL;
 72:   PetscFunctionReturn(PETSC_SUCCESS);
 73: }
 75: /*@
 76:   DMPlexPointLocalRef - return read/write access to a point in local array
 78:   Not Collective
 80:   Input Parameters:
 81: + dm    - `DM` defining topological space
 82: . point - topological point
 83: - array - array to index into
 85:   Output Parameter:
 86: . ptr - address of reference to point data, type generic so user can place in structure
 88:   Level: intermediate
 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
 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;
104:   PetscFunctionBegin;
106:   PetscAssertPointer(array, 3);
107:   PetscAssertPointer(ptr, 4);
108:   PetscCall(DMGetLocalOffset_Private(dm, point, &start, &end));
109:   *(PetscScalar **)ptr = (start < end) ? array + start : NULL;
110:   PetscFunctionReturn(PETSC_SUCCESS);
111: }
113: /*@
114:   DMPlexGetPointLocalField - get location of point field data in local Vec
116:   Not Collective
118:   Input Parameters:
119: + dm    - `DM` defining the topological space
120: . point - topological point
121: - field - the field number
123:   Output Parameters:
124: + start - start of point data
125: - end   - end of point data
127:   Level: intermediate
129:   Note:
130:   This is a half open interval [start, end)
132: .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMPlexGetPointLocal()`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexPointLocalRead()`, `DMPlexPointLocalRef()`
133: @*/
134: PetscErrorCode DMPlexGetPointLocalField(DM dm, PetscInt point, PetscInt field, PetscInt *start, PetscInt *end)
135: {
136:   PetscInt s, e;
138:   PetscFunctionBegin;
140:   if (start) PetscAssertPointer(start, 4);
141:   if (end) PetscAssertPointer(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: }
148: /*@
149:   DMPlexPointLocalFieldRead - return read access to a field on a point in local array
151:   Not Collective
153:   Input Parameters:
154: + dm    - `DM` defining topological space
155: . point - topological point
156: . field - field number
157: - array - array to index into
159:   Output Parameter:
160: . ptr - address of read reference to point data, type generic so user can place in structure
162:   Level: intermediate
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;
170:   PetscFunctionBegin;
172:   PetscAssertPointer(array, 4);
173:   PetscAssertPointer(ptr, 5);
174:   PetscCall(DMGetLocalFieldOffset_Private(dm, point, field, &start, &end));
175:   *(const PetscScalar **)ptr = array + start;
176:   PetscFunctionReturn(PETSC_SUCCESS);
177: }
179: /*@
180:   DMPlexPointLocalFieldRef - return read/write access to a field on a point in local array
182:   Not Collective
184:   Input Parameters:
185: + dm    - `DM` defining topological space
186: . point - topological point
187: . field - field number
188: - array - array to index into
190:   Output Parameter:
191: . ptr - address of reference to point data, type generic so user can place in structure
193:   Level: intermediate
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;
201:   PetscFunctionBegin;
203:   PetscAssertPointer(array, 4);
204:   PetscAssertPointer(ptr, 5);
205:   PetscCall(DMGetLocalFieldOffset_Private(dm, point, field, &start, &end));
206:   *(PetscScalar **)ptr = array + start;
207:   PetscFunctionReturn(PETSC_SUCCESS);
208: }
210: /*@
211:   DMPlexGetPointGlobal - get location of point data in global Vec
213:   Not Collective
215:   Input Parameters:
216: + dm    - `DM` defining the topological space
217: - point - topological point
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
223:   Level: intermediate
225:   Note:
226:   This is a half open interval [start, end)
228: .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMPlexGetPointGlobalField()`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexPointGlobalRead()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRef()`
229: @*/
230: PetscErrorCode DMPlexGetPointGlobal(DM dm, PetscInt point, PetscInt *start, PetscInt *end)
231: {
232:   PetscInt s, e;
234:   PetscFunctionBegin;
236:   if (start) PetscAssertPointer(start, 3);
237:   if (end) PetscAssertPointer(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: }
244: /*@
245:   DMPlexPointGlobalRead - return read access to a point in global array
247:   Not Collective
249:   Input Parameters:
250: + dm    - `DM` defining topological space
251: . point - topological point
252: - array - array to index into
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
257:   Level: intermediate
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
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;
273:   PetscFunctionBegin;
275:   PetscAssertPointer(array, 3);
276:   PetscAssertPointer(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: }
282: /*@
283:   DMPlexPointGlobalRef - return read/write access to a point in global array
285:   Not Collective
287:   Input Parameters:
288: + dm    - `DM` defining topological space
289: . point - topological point
290: - array - array to index into
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
295:   Level: intermediate
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
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;
311:   PetscFunctionBegin;
313:   PetscAssertPointer(array, 3);
314:   PetscAssertPointer(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: }
320: /*@
321:   DMPlexGetPointGlobalField - get location of point field data in global `Vec`
323:   Not Collective
325:   Input Parameters:
326: + dm    - `DM` defining the topological space
327: . point - topological point
328: - field - the field number
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
334:   Level: intermediate
336:   Note:
337:   This is a half open interval [start, end)
339: .seealso: [](ch_unstructured), `DM`, `DMPLEX`, `DMPlexGetPointGlobal()`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexPointGlobalRead()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRef()`
340: @*/
341: PetscErrorCode DMPlexGetPointGlobalField(DM dm, PetscInt point, PetscInt field, PetscInt *start, PetscInt *end)
342: {
343:   PetscInt s, e;
345:   PetscFunctionBegin;
347:   if (start) PetscAssertPointer(start, 4);
348:   if (end) PetscAssertPointer(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: }
355: /*@
356:   DMPlexPointGlobalFieldRead - return read access to a field on a point in global array
358:   Not Collective
360:   Input Parameters:
361: + dm    - `DM` defining topological space
362: . point - topological point
363: . field - field number
364: - array - array to index into
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
369:   Level: intermediate
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;
377:   PetscFunctionBegin;
379:   PetscAssertPointer(array, 4);
380:   PetscAssertPointer(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: }
386: /*@
387:   DMPlexPointGlobalFieldRef - return read/write access to a field on a point in global array
389:   Not Collective
391:   Input Parameters:
392: + dm    - `DM` defining topological space
393: . point - topological point
394: . field - field number
395: - array - array to index into
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
400:   Level: intermediate
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;
408:   PetscFunctionBegin;
410:   PetscAssertPointer(array, 4);
411:   PetscAssertPointer(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: }