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: Note: This is a half open interval [start, end)
18: Level: intermediate
20: .seealso: `DMPlexGetPointLocalField()`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexPointLocalRead()`, `DMPlexPointLocalRead()`, `DMPlexPointLocalRef()`
21: @*/
22: PetscErrorCode DMPlexGetPointLocal(DM dm, PetscInt point, PetscInt *start, PetscInt *end)
23: {
24: PetscInt s, e;
29: DMGetLocalOffset_Private(dm, point, &s, &e);
30: if (start) *start = s;
31: if (end) *end = e;
32: return 0;
33: }
35: /*@
36: DMPlexPointLocalRead - return read access to a point in local array
38: Not Collective
40: Input Parameters:
41: + dm - DM defining topological space
42: . point - topological point
43: - array - array to index into
45: Output Parameter:
46: . ptr - address of read reference to point data, type generic so user can place in structure
48: Level: intermediate
50: Note:
51: A common usage when data sizes are known statically:
53: $ const struct { PetscScalar foo,bar,baz; } *ptr;
54: $ DMPlexPointLocalRead(dm,point,array,&ptr);
55: $ x = 2*ptr->foo + 3*ptr->bar + 5*ptr->baz;
57: .seealso: `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRead()`
58: @*/
59: PetscErrorCode DMPlexPointLocalRead(DM dm, PetscInt point, const PetscScalar *array, void *ptr)
60: {
61: PetscInt start, end;
66: DMGetLocalOffset_Private(dm, point, &start, &end);
67: *(const PetscScalar **)ptr = (start < end) ? array + start : NULL;
68: return 0;
69: }
71: /*@
72: DMPlexPointLocalRef - return read/write access to a point in local array
74: Not Collective
76: Input Parameters:
77: + dm - DM defining topological space
78: . point - topological point
79: - array - array to index into
81: Output Parameter:
82: . ptr - address of reference to point data, type generic so user can place in structure
84: Level: intermediate
86: Note:
87: A common usage when data sizes are known statically:
89: $ struct { PetscScalar foo,bar,baz; } *ptr;
90: $ DMPlexPointLocalRef(dm,point,array,&ptr);
91: $ ptr->foo = 2; ptr->bar = 3; ptr->baz = 5;
93: .seealso: `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRef()`
94: @*/
95: PetscErrorCode DMPlexPointLocalRef(DM dm, PetscInt point, PetscScalar *array, void *ptr)
96: {
97: PetscInt start, end;
102: DMGetLocalOffset_Private(dm, point, &start, &end);
103: *(PetscScalar **)ptr = (start < end) ? array + start : NULL;
104: return 0;
105: }
107: /*@
108: DMPlexGetPointLocalField - get location of point field data in local Vec
110: Not Collective
112: Input Parameters:
113: + dm - DM defining the topological space
114: . point - topological point
115: - field - the field number
117: Output Parameters:
118: + start - start of point data
119: - end - end of point data
121: Note: This is a half open interval [start, end)
123: Level: intermediate
125: .seealso: `DMPlexGetPointLocal()`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexPointLocalRead()`, `DMPlexPointLocalRead()`, `DMPlexPointLocalRef()`
126: @*/
127: PetscErrorCode DMPlexGetPointLocalField(DM dm, PetscInt point, PetscInt field, PetscInt *start, PetscInt *end)
128: {
129: PetscInt s, e;
134: DMGetLocalFieldOffset_Private(dm, point, field, &s, &e);
135: if (start) *start = s;
136: if (end) *end = e;
137: return 0;
138: }
140: /*@
141: DMPlexPointLocalFieldRead - return read access to a field on a point in local array
143: Not Collective
145: Input Parameters:
146: + dm - DM defining topological space
147: . point - topological point
148: . field - field number
149: - array - array to index into
151: Output Parameter:
152: . ptr - address of read reference to point data, type generic so user can place in structure
154: Level: intermediate
156: .seealso: `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRef()`
157: @*/
158: PetscErrorCode DMPlexPointLocalFieldRead(DM dm, PetscInt point, PetscInt field, const PetscScalar *array, void *ptr)
159: {
160: PetscInt start, end;
165: DMGetLocalFieldOffset_Private(dm, point, field, &start, &end);
166: *(const PetscScalar **)ptr = array + start;
167: return 0;
168: }
170: /*@
171: DMPlexPointLocalFieldRef - return read/write access to a field on a point in local array
173: Not Collective
175: Input Parameters:
176: + dm - DM defining topological space
177: . point - topological point
178: . field - field number
179: - array - array to index into
181: Output Parameter:
182: . ptr - address of reference to point data, type generic so user can place in structure
184: Level: intermediate
186: .seealso: `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRef()`
187: @*/
188: PetscErrorCode DMPlexPointLocalFieldRef(DM dm, PetscInt point, PetscInt field, PetscScalar *array, void *ptr)
189: {
190: PetscInt start, end;
195: DMGetLocalFieldOffset_Private(dm, point, field, &start, &end);
196: *(PetscScalar **)ptr = array + start;
197: return 0;
198: }
200: /*@
201: DMPlexGetPointGlobal - get location of point data in global Vec
203: Not Collective
205: Input Parameters:
206: + dm - DM defining the topological space
207: - point - topological point
209: Output Parameters:
210: + start - start of point data; returns -(globalStart+1) if point is not owned
211: - end - end of point data; returns -(globalEnd+1) if point is not owned
213: Note: This is a half open interval [start, end)
215: Level: intermediate
217: .seealso: `DMPlexGetPointGlobalField()`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexPointGlobalRead()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRead()`, `DMPlexPointGlobalRef()`
218: @*/
219: PetscErrorCode DMPlexGetPointGlobal(DM dm, PetscInt point, PetscInt *start, PetscInt *end)
220: {
221: PetscInt s, e;
226: DMGetGlobalOffset_Private(dm, point, &s, &e);
227: if (start) *start = s;
228: if (end) *end = e;
229: return 0;
230: }
232: /*@
233: DMPlexPointGlobalRead - return read access to a point in global array
235: Not Collective
237: Input Parameters:
238: + dm - DM defining topological space
239: . point - topological point
240: - array - array to index into
242: Output Parameter:
243: . ptr - address of read reference to point data, type generic so user can place in structure; returns NULL if global point is not owned
245: Level: intermediate
247: Note:
248: A common usage when data sizes are known statically:
250: $ const struct { PetscScalar foo,bar,baz; } *ptr;
251: $ DMPlexPointGlobalRead(dm,point,array,&ptr);
252: $ x = 2*ptr->foo + 3*ptr->bar + 5*ptr->baz;
254: .seealso: `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointGlobal()`, `DMPlexPointLocalRead()`, `DMPlexPointGlobalRef()`
255: @*/
256: PetscErrorCode DMPlexPointGlobalRead(DM dm, PetscInt point, const PetscScalar *array, const void *ptr)
257: {
258: PetscInt start, end;
263: DMGetGlobalOffset_Private(dm, point, &start, &end);
264: *(const PetscScalar **)ptr = (start < end) ? array + start - dm->map->rstart : NULL;
265: return 0;
266: }
268: /*@
269: DMPlexPointGlobalRef - return read/write access to a point in global array
271: Not Collective
273: Input Parameters:
274: + dm - DM defining topological space
275: . point - topological point
276: - array - array to index into
278: Output Parameter:
279: . ptr - address of reference to point data, type generic so user can place in structure; returns NULL if global point is not owned
281: Level: intermediate
283: Note:
284: A common usage when data sizes are known statically:
286: $ struct { PetscScalar foo,bar,baz; } *ptr;
287: $ DMPlexPointGlobalRef(dm,point,array,&ptr);
288: $ ptr->foo = 2; ptr->bar = 3; ptr->baz = 5;
290: .seealso: `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointGlobal()`, `DMPlexPointLocalRef()`, `DMPlexPointGlobalRead()`
291: @*/
292: PetscErrorCode DMPlexPointGlobalRef(DM dm, PetscInt point, PetscScalar *array, void *ptr)
293: {
294: PetscInt start, end;
299: DMGetGlobalOffset_Private(dm, point, &start, &end);
300: *(PetscScalar **)ptr = (start < end) ? array + start - dm->map->rstart : NULL;
301: return 0;
302: }
304: /*@
305: DMPlexGetPointGlobalField - get location of point field data in global Vec
307: Not Collective
309: Input Parameters:
310: + dm - DM defining the topological space
311: . point - topological point
312: - field - the field number
314: Output Parameters:
315: + start - start of point data; returns -(globalStart+1) if point is not owned
316: - end - end of point data; returns -(globalEnd+1) if point is not owned
318: Note: This is a half open interval [start, end)
320: Level: intermediate
322: .seealso: `DMPlexGetPointGlobal()`, `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexPointGlobalRead()`, `DMPlexGetPointLocal()`, `DMPlexPointGlobalRead()`, `DMPlexPointGlobalRef()`
323: @*/
324: PetscErrorCode DMPlexGetPointGlobalField(DM dm, PetscInt point, PetscInt field, PetscInt *start, PetscInt *end)
325: {
326: PetscInt s, e;
331: DMGetGlobalFieldOffset_Private(dm, point, field, &s, &e);
332: if (start) *start = s;
333: if (end) *end = e;
334: return 0;
335: }
337: /*@
338: DMPlexPointGlobalFieldRead - return read access to a field on a point in global array
340: Not Collective
342: Input Parameters:
343: + dm - DM defining topological space
344: . point - topological point
345: . field - field number
346: - array - array to index into
348: Output Parameter:
349: . ptr - address of read reference to point data, type generic so user can place in structure; returns NULL if global point is not owned
351: Level: intermediate
353: .seealso: `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointGlobal()`, `DMPlexPointLocalRead()`, `DMPlexPointGlobalRef()`
354: @*/
355: PetscErrorCode DMPlexPointGlobalFieldRead(DM dm, PetscInt point, PetscInt field, const PetscScalar *array, void *ptr)
356: {
357: PetscInt start, end;
362: DMGetGlobalFieldOffset_Private(dm, point, field, &start, &end);
363: *(const PetscScalar **)ptr = (start < end) ? array + start - dm->map->rstart : NULL;
364: return 0;
365: }
367: /*@
368: DMPlexPointGlobalFieldRef - return read/write access to a field on a point in global array
370: Not Collective
372: Input Parameters:
373: + dm - DM defining topological space
374: . point - topological point
375: . field - field number
376: - array - array to index into
378: Output Parameter:
379: . ptr - address of reference to point data, type generic so user can place in structure; returns NULL if global point is not owned
381: Level: intermediate
383: .seealso: `DMGetLocalSection()`, `PetscSectionGetOffset()`, `PetscSectionGetDof()`, `DMPlexGetPointGlobal()`, `DMPlexPointLocalRef()`, `DMPlexPointGlobalRead()`
384: @*/
385: PetscErrorCode DMPlexPointGlobalFieldRef(DM dm, PetscInt point, PetscInt field, PetscScalar *array, void *ptr)
386: {
387: PetscInt start, end;
392: DMGetGlobalFieldOffset_Private(dm, point, field, &start, &end);
393: *(PetscScalar **)ptr = (start < end) ? array + start - dm->map->rstart : NULL;
394: return 0;
395: }