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: [](chapter_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;
30: DMGetLocalOffset_Private(dm, point, &s, &e);
31: if (start) *start = s;
32: if (end) *end = e;
33: return 0;
34: }
36: /*@
37: DMPlexPointLocalRead - return read access to a point in local array
39: Not Collective
41: Input Parameters:
42: + dm - `DM` defining topological space
43: . point - topological point
44: - array - array to index into
46: Output Parameter:
47: . ptr - address of read reference to point data, type generic so user can place in structure
49: Level: intermediate
51: Note:
52: A common usage when data sizes are known statically:
53: .vb
54: const struct { PetscScalar foo,bar,baz; } *ptr;
55: DMPlexPointLocalRead(dm,point,array,&ptr);
56: x = 2*ptr->foo + 3*ptr->bar + 5*ptr->baz;
57: .ve
59: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRead()`
60: @*/
61: PetscErrorCode DMPlexPointLocalRead(DM dm, PetscInt point, const PetscScalar *array, void *ptr)
62: {
63: PetscInt start, end;
68: DMGetLocalOffset_Private(dm, point, &start, &end);
69: *(const PetscScalar **)ptr = (start < end) ? array + start : NULL;
70: return 0;
71: }
73: /*@
74: DMPlexPointLocalRef - return read/write access to a point in local array
76: Not Collective
78: Input Parameters:
79: + dm - `DM` defining topological space
80: . point - topological point
81: - array - array to index into
83: Output Parameter:
84: . ptr - address of reference to point data, type generic so user can place in structure
86: Level: intermediate
88: Note:
89: A common usage when data sizes are known statically:
90: .vb
91: struct { PetscScalar foo,bar,baz; } *ptr;
92: DMPlexPointLocalRef(dm,point,array,&ptr);
93: ptr->foo = 2; ptr->bar = 3; ptr->baz = 5;
94: .ve
96: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRef()`
97: @*/
98: PetscErrorCode DMPlexPointLocalRef(DM dm, PetscInt point, PetscScalar *array, void *ptr)
99: {
100: PetscInt start, end;
105: DMGetLocalOffset_Private(dm, point, &start, &end);
106: *(PetscScalar **)ptr = (start < end) ? array + start : NULL;
107: return 0;
108: }
110: /*@
111: DMPlexGetPointLocalField - get location of point field data in local Vec
113: Not Collective
115: Input Parameters:
116: + dm - `DM` defining the topological space
117: . point - topological point
118: - field - the field number
120: Output Parameters:
121: + start - start of point data
122: - end - end of point data
124: Level: intermediate
126: Note:
127: This is a half open interval [start, end)
129: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMPlexGetPointLocal()`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexPointLocalRead()`, `DMPlexPointLocalRead()`, `DMPlexPointLocalRef()`
130: @*/
131: PetscErrorCode DMPlexGetPointLocalField(DM dm, PetscInt point, PetscInt field, PetscInt *start, PetscInt *end)
132: {
133: PetscInt s, e;
138: DMGetLocalFieldOffset_Private(dm, point, field, &s, &e);
139: if (start) *start = s;
140: if (end) *end = e;
141: return 0;
142: }
144: /*@
145: DMPlexPointLocalFieldRead - return read access to a field on a point in local array
147: Not Collective
149: Input Parameters:
150: + dm - `DM` defining topological space
151: . point - topological point
152: . field - field number
153: - array - array to index into
155: Output Parameter:
156: . ptr - address of read reference to point data, type generic so user can place in structure
158: Level: intermediate
160: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRef()`
161: @*/
162: PetscErrorCode DMPlexPointLocalFieldRead(DM dm, PetscInt point, PetscInt field, const PetscScalar *array, void *ptr)
163: {
164: PetscInt start, end;
169: DMGetLocalFieldOffset_Private(dm, point, field, &start, &end);
170: *(const PetscScalar **)ptr = array + start;
171: return 0;
172: }
174: /*@
175: DMPlexPointLocalFieldRef - return read/write access to a field on a point in local array
177: Not Collective
179: Input Parameters:
180: + dm - `DM` defining topological space
181: . point - topological point
182: . field - field number
183: - array - array to index into
185: Output Parameter:
186: . ptr - address of reference to point data, type generic so user can place in structure
188: Level: intermediate
190: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRef()`
191: @*/
192: PetscErrorCode DMPlexPointLocalFieldRef(DM dm, PetscInt point, PetscInt field, PetscScalar *array, void *ptr)
193: {
194: PetscInt start, end;
199: DMGetLocalFieldOffset_Private(dm, point, field, &start, &end);
200: *(PetscScalar **)ptr = array + start;
201: return 0;
202: }
204: /*@
205: DMPlexGetPointGlobal - get location of point data in global Vec
207: Not Collective
209: Input Parameters:
210: + dm - `DM` defining the topological space
211: - point - topological point
213: Output Parameters:
214: + start - start of point data; returns -(globalStart+1) if point is not owned
215: - end - end of point data; returns -(globalEnd+1) if point is not owned
217: Level: intermediate
219: Note:
220: This is a half open interval [start, end)
222: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMPlexGetPointGlobalField()`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexPointGlobalRead()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRead()`, `DMPlexPointGlobalRef()`
223: @*/
224: PetscErrorCode DMPlexGetPointGlobal(DM dm, PetscInt point, PetscInt *start, PetscInt *end)
225: {
226: PetscInt s, e;
231: DMGetGlobalOffset_Private(dm, point, &s, &e);
232: if (start) *start = s;
233: if (end) *end = e;
234: return 0;
235: }
237: /*@
238: DMPlexPointGlobalRead - return read access to a point in global array
240: Not Collective
242: Input Parameters:
243: + dm - `DM` defining topological space
244: . point - topological point
245: - array - array to index into
247: Output Parameter:
248: . ptr - address of read reference to point data, type generic so user can place in structure; returns NULL if global point is not owned
250: Level: intermediate
252: Note:
253: A common usage when data sizes are known statically:
254: .vb
255: const struct { PetscScalar foo,bar,baz; } *ptr;
256: DMPlexPointGlobalRead(dm,point,array,&ptr);
257: x = 2*ptr->foo + 3*ptr->bar + 5*ptr->baz;
258: .ve
260: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointGlobal()`, `DMPlexPointLocalRead()`, `DMPlexPointGlobalRef()`
261: @*/
262: PetscErrorCode DMPlexPointGlobalRead(DM dm, PetscInt point, const PetscScalar *array, const void *ptr)
263: {
264: PetscInt start, end;
269: DMGetGlobalOffset_Private(dm, point, &start, &end);
270: *(const PetscScalar **)ptr = (start < end) ? array + start - dm->map->rstart : NULL;
271: return 0;
272: }
274: /*@
275: DMPlexPointGlobalRef - return read/write access to a point in global array
277: Not Collective
279: Input Parameters:
280: + dm - `DM` defining topological space
281: . point - topological point
282: - array - array to index into
284: Output Parameter:
285: . ptr - address of reference to point data, type generic so user can place in structure; returns NULL if global point is not owned
287: Level: intermediate
289: Note:
290: A common usage when data sizes are known statically:
291: .vb
292: struct { PetscScalar foo,bar,baz; } *ptr;
293: DMPlexPointGlobalRef(dm,point,array,&ptr);
294: ptr->foo = 2; ptr->bar = 3; ptr->baz = 5;
295: .ve
297: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointGlobal()`, `DMPlexPointLocalRef()`, `DMPlexPointGlobalRead()`
298: @*/
299: PetscErrorCode DMPlexPointGlobalRef(DM dm, PetscInt point, PetscScalar *array, void *ptr)
300: {
301: PetscInt start, end;
306: DMGetGlobalOffset_Private(dm, point, &start, &end);
307: *(PetscScalar **)ptr = (start < end) ? array + start - dm->map->rstart : NULL;
308: return 0;
309: }
311: /*@
312: DMPlexGetPointGlobalField - get location of point field data in global `Vec`
314: Not Collective
316: Input Parameters:
317: + dm - `DM` defining the topological space
318: . point - topological point
319: - field - the field number
321: Output Parameters:
322: + start - start of point data; returns -(globalStart+1) if point is not owned
323: - end - end of point data; returns -(globalEnd+1) if point is not owned
325: Level: intermediate
327: Note:
328: This is a half open interval [start, end)
330: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMPlexGetPointGlobal()`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexPointGlobalRead()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRead()`, `DMPlexPointGlobalRef()`
331: @*/
332: PetscErrorCode DMPlexGetPointGlobalField(DM dm, PetscInt point, PetscInt field, PetscInt *start, PetscInt *end)
333: {
334: PetscInt s, e;
339: DMGetGlobalFieldOffset_Private(dm, point, field, &s, &e);
340: if (start) *start = s;
341: if (end) *end = e;
342: return 0;
343: }
345: /*@
346: DMPlexPointGlobalFieldRead - return read access to a field on a point in global array
348: Not Collective
350: Input Parameters:
351: + dm - `DM` defining topological space
352: . point - topological point
353: . field - field number
354: - array - array to index into
356: Output Parameter:
357: . ptr - address of read reference to point data, type generic so user can place in structure; returns NULL if global point is not owned
359: Level: intermediate
361: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointGlobal()`, `DMPlexPointLocalRead()`, `DMPlexPointGlobalRef()`
362: @*/
363: PetscErrorCode DMPlexPointGlobalFieldRead(DM dm, PetscInt point, PetscInt field, const PetscScalar *array, void *ptr)
364: {
365: PetscInt start, end;
370: DMGetGlobalFieldOffset_Private(dm, point, field, &start, &end);
371: *(const PetscScalar **)ptr = (start < end) ? array + start - dm->map->rstart : NULL;
372: return 0;
373: }
375: /*@
376: DMPlexPointGlobalFieldRef - return read/write access to a field on a point in global array
378: Not Collective
380: Input Parameters:
381: + dm - `DM` defining topological space
382: . point - topological point
383: . field - field number
384: - array - array to index into
386: Output Parameter:
387: . ptr - address of reference to point data, type generic so user can place in structure; returns NULL if global point is not owned
389: Level: intermediate
391: .seealso: [](chapter_unstructured), `DM`, `DMPLEX`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointGlobal()`, `DMPlexPointLocalRef()`, `DMPlexPointGlobalRead()`
392: @*/
393: PetscErrorCode DMPlexPointGlobalFieldRef(DM dm, PetscInt point, PetscInt field, PetscScalar *array, void *ptr)
394: {
395: PetscInt start, end;
400: DMGetGlobalFieldOffset_Private(dm, point, field, &start, &end);
401: *(PetscScalar **)ptr = (start < end) ? array + start - dm->map->rstart : NULL;
402: return 0;
403: }