Actual source code: pack.c
2: #include <../src/dm/impls/composite/packimpl.h>
3: #include <petsc/private/isimpl.h>
4: #include <petsc/private/glvisviewerimpl.h>
5: #include <petscds.h>
7: /*@C
8: DMCompositeSetCoupling - Sets user provided routines that compute the coupling between the
9: separate components `DM` in a `DMCOMPOSITE` to build the correct matrix nonzero structure.
11: Logically Collective; No Fortran Support
13: Input Parameters:
14: + dm - the composite object
15: - formcouplelocations - routine to set the nonzero locations in the matrix
17: Level: advanced
19: Note:
20: See `DMSetApplicationContext()` and `DMGetApplicationContext()` for how to get user information into
21: this routine
23: .seealso: `DMCOMPOSITE`, `DM`
24: @*/
25: PetscErrorCode DMCompositeSetCoupling(DM dm, PetscErrorCode (*FormCoupleLocations)(DM, Mat, PetscInt *, PetscInt *, PetscInt, PetscInt, PetscInt, PetscInt))
26: {
27: DM_Composite *com = (DM_Composite *)dm->data;
28: PetscBool flg;
30: PetscFunctionBegin;
31: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
32: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
33: com->FormCoupleLocations = FormCoupleLocations;
34: PetscFunctionReturn(PETSC_SUCCESS);
35: }
37: PetscErrorCode DMDestroy_Composite(DM dm)
38: {
39: struct DMCompositeLink *next, *prev;
40: DM_Composite *com = (DM_Composite *)dm->data;
42: PetscFunctionBegin;
43: next = com->next;
44: while (next) {
45: prev = next;
46: next = next->next;
47: PetscCall(DMDestroy(&prev->dm));
48: PetscCall(PetscFree(prev->grstarts));
49: PetscCall(PetscFree(prev));
50: }
51: PetscCall(PetscObjectComposeFunction((PetscObject)dm, "DMSetUpGLVisViewer_C", NULL));
52: /* This was originally freed in DMDestroy(), but that prevents reference counting of backend objects */
53: PetscCall(PetscFree(com));
54: PetscFunctionReturn(PETSC_SUCCESS);
55: }
57: PetscErrorCode DMView_Composite(DM dm, PetscViewer v)
58: {
59: PetscBool iascii;
60: DM_Composite *com = (DM_Composite *)dm->data;
62: PetscFunctionBegin;
63: PetscCall(PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &iascii));
64: if (iascii) {
65: struct DMCompositeLink *lnk = com->next;
66: PetscInt i;
68: PetscCall(PetscViewerASCIIPrintf(v, "DM (%s)\n", ((PetscObject)dm)->prefix ? ((PetscObject)dm)->prefix : "no prefix"));
69: PetscCall(PetscViewerASCIIPrintf(v, " contains %" PetscInt_FMT " DMs\n", com->nDM));
70: PetscCall(PetscViewerASCIIPushTab(v));
71: for (i = 0; lnk; lnk = lnk->next, i++) {
72: PetscCall(PetscViewerASCIIPrintf(v, "Link %" PetscInt_FMT ": DM of type %s\n", i, ((PetscObject)lnk->dm)->type_name));
73: PetscCall(PetscViewerASCIIPushTab(v));
74: PetscCall(DMView(lnk->dm, v));
75: PetscCall(PetscViewerASCIIPopTab(v));
76: }
77: PetscCall(PetscViewerASCIIPopTab(v));
78: }
79: PetscFunctionReturn(PETSC_SUCCESS);
80: }
82: /* --------------------------------------------------------------------------------------*/
83: PetscErrorCode DMSetUp_Composite(DM dm)
84: {
85: PetscInt nprev = 0;
86: PetscMPIInt rank, size;
87: DM_Composite *com = (DM_Composite *)dm->data;
88: struct DMCompositeLink *next = com->next;
89: PetscLayout map;
91: PetscFunctionBegin;
92: PetscCheck(!com->setup, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Packer has already been setup");
93: PetscCall(PetscLayoutCreate(PetscObjectComm((PetscObject)dm), &map));
94: PetscCall(PetscLayoutSetLocalSize(map, com->n));
95: PetscCall(PetscLayoutSetSize(map, PETSC_DETERMINE));
96: PetscCall(PetscLayoutSetBlockSize(map, 1));
97: PetscCall(PetscLayoutSetUp(map));
98: PetscCall(PetscLayoutGetSize(map, &com->N));
99: PetscCall(PetscLayoutGetRange(map, &com->rstart, NULL));
100: PetscCall(PetscLayoutDestroy(&map));
102: /* now set the rstart for each linked vector */
103: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
104: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)dm), &size));
105: while (next) {
106: next->rstart = nprev;
107: nprev += next->n;
108: next->grstart = com->rstart + next->rstart;
109: PetscCall(PetscMalloc1(size, &next->grstarts));
110: PetscCallMPI(MPI_Allgather(&next->grstart, 1, MPIU_INT, next->grstarts, 1, MPIU_INT, PetscObjectComm((PetscObject)dm)));
111: next = next->next;
112: }
113: com->setup = PETSC_TRUE;
114: PetscFunctionReturn(PETSC_SUCCESS);
115: }
117: /* ----------------------------------------------------------------------------------*/
119: /*@
120: DMCompositeGetNumberDM - Gets the number of `DM` objects in the `DMCOMPOSITE`
121: representation.
123: Not Collective
125: Input Parameter:
126: . dm - the `DMCOMPOSITE` object
128: Output Parameter:
129: . nDM - the number of `DM`
131: Level: beginner
133: .seealso: `DMCOMPOSITE`, `DM`
134: @*/
135: PetscErrorCode DMCompositeGetNumberDM(DM dm, PetscInt *nDM)
136: {
137: DM_Composite *com = (DM_Composite *)dm->data;
138: PetscBool flg;
140: PetscFunctionBegin;
142: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
143: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
144: *nDM = com->nDM;
145: PetscFunctionReturn(PETSC_SUCCESS);
146: }
148: /*@C
149: DMCompositeGetAccess - Allows one to access the individual packed vectors in their global
150: representation.
152: Collective
154: Input Parameters:
155: + dm - the `DMCOMPOSITE` object
156: - gvec - the global vector
158: Output Parameter:
159: . Vec* ... - the packed parallel vectors, NULL for those that are not needed
161: Level: advanced
163: Note:
164: Use `DMCompositeRestoreAccess()` to return the vectors when you no longer need them
166: Fortran Note:
167: Fortran callers must use numbered versions of this routine, e.g., DMCompositeGetAccess4(dm,gvec,vec1,vec2,vec3,vec4)
168: or use the alternative interface `DMCompositeGetAccessArray()`.
170: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeGetEntries()`, `DMCompositeScatter()`
171: @*/
172: PetscErrorCode DMCompositeGetAccess(DM dm, Vec gvec, ...)
173: {
174: va_list Argp;
175: struct DMCompositeLink *next;
176: DM_Composite *com = (DM_Composite *)dm->data;
177: PetscInt readonly;
178: PetscBool flg;
180: PetscFunctionBegin;
183: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
184: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
185: next = com->next;
186: if (!com->setup) PetscCall(DMSetUp(dm));
188: PetscCall(VecLockGet(gvec, &readonly));
189: /* loop over packed objects, handling one at at time */
190: va_start(Argp, gvec);
191: while (next) {
192: Vec *vec;
193: vec = va_arg(Argp, Vec *);
194: if (vec) {
195: PetscCall(DMGetGlobalVector(next->dm, vec));
196: if (readonly) {
197: const PetscScalar *array;
198: PetscCall(VecGetArrayRead(gvec, &array));
199: PetscCall(VecPlaceArray(*vec, array + next->rstart));
200: PetscCall(VecLockReadPush(*vec));
201: PetscCall(VecRestoreArrayRead(gvec, &array));
202: } else {
203: PetscScalar *array;
204: PetscCall(VecGetArray(gvec, &array));
205: PetscCall(VecPlaceArray(*vec, array + next->rstart));
206: PetscCall(VecRestoreArray(gvec, &array));
207: }
208: }
209: next = next->next;
210: }
211: va_end(Argp);
212: PetscFunctionReturn(PETSC_SUCCESS);
213: }
215: /*@C
216: DMCompositeGetAccessArray - Allows one to access the individual packed vectors in their global
217: representation.
219: Collective
221: Input Parameters:
222: + dm - the `DMCOMPOSITE`
223: . pvec - packed vector
224: . nwanted - number of vectors wanted
225: - wanted - sorted array of vectors wanted, or NULL to get all vectors
227: Output Parameter:
228: . vecs - array of requested global vectors (must be allocated)
230: Level: advanced
232: Note:
233: Use `DMCompositeRestoreAccessArray()` to return the vectors when you no longer need them
235: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeGetAccess()`, `DMCompositeGetEntries()`, `DMCompositeScatter()`, `DMCompositeGather()`
236: @*/
237: PetscErrorCode DMCompositeGetAccessArray(DM dm, Vec pvec, PetscInt nwanted, const PetscInt *wanted, Vec *vecs)
238: {
239: struct DMCompositeLink *link;
240: PetscInt i, wnum;
241: DM_Composite *com = (DM_Composite *)dm->data;
242: PetscInt readonly;
243: PetscBool flg;
245: PetscFunctionBegin;
248: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
249: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
250: if (!com->setup) PetscCall(DMSetUp(dm));
252: PetscCall(VecLockGet(pvec, &readonly));
253: for (i = 0, wnum = 0, link = com->next; link && wnum < nwanted; i++, link = link->next) {
254: if (!wanted || i == wanted[wnum]) {
255: Vec v;
256: PetscCall(DMGetGlobalVector(link->dm, &v));
257: if (readonly) {
258: const PetscScalar *array;
259: PetscCall(VecGetArrayRead(pvec, &array));
260: PetscCall(VecPlaceArray(v, array + link->rstart));
261: PetscCall(VecLockReadPush(v));
262: PetscCall(VecRestoreArrayRead(pvec, &array));
263: } else {
264: PetscScalar *array;
265: PetscCall(VecGetArray(pvec, &array));
266: PetscCall(VecPlaceArray(v, array + link->rstart));
267: PetscCall(VecRestoreArray(pvec, &array));
268: }
269: vecs[wnum++] = v;
270: }
271: }
272: PetscFunctionReturn(PETSC_SUCCESS);
273: }
275: /*@C
276: DMCompositeGetLocalAccessArray - Allows one to access the individual
277: packed vectors in their local representation.
279: Collective
281: Input Parameters:
282: + dm - the `DMCOMPOSITE`
283: . pvec - packed vector
284: . nwanted - number of vectors wanted
285: - wanted - sorted array of vectors wanted, or NULL to get all vectors
287: Output Parameter:
288: . vecs - array of requested local vectors (must be allocated)
290: Level: advanced
292: Note:
293: Use `DMCompositeRestoreLocalAccessArray()` to return the vectors
294: when you no longer need them.
296: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeRestoreLocalAccessArray()`, `DMCompositeGetAccess()`,
297: `DMCompositeGetEntries()`, `DMCompositeScatter()`, `DMCompositeGather()`
298: @*/
299: PetscErrorCode DMCompositeGetLocalAccessArray(DM dm, Vec pvec, PetscInt nwanted, const PetscInt *wanted, Vec *vecs)
300: {
301: struct DMCompositeLink *link;
302: PetscInt i, wnum;
303: DM_Composite *com = (DM_Composite *)dm->data;
304: PetscInt readonly;
305: PetscInt nlocal = 0;
306: PetscBool flg;
308: PetscFunctionBegin;
311: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
312: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
313: if (!com->setup) PetscCall(DMSetUp(dm));
315: PetscCall(VecLockGet(pvec, &readonly));
316: for (i = 0, wnum = 0, link = com->next; link && wnum < nwanted; i++, link = link->next) {
317: if (!wanted || i == wanted[wnum]) {
318: Vec v;
319: PetscCall(DMGetLocalVector(link->dm, &v));
320: if (readonly) {
321: const PetscScalar *array;
322: PetscCall(VecGetArrayRead(pvec, &array));
323: PetscCall(VecPlaceArray(v, array + nlocal));
324: // this method does not make sense. The local vectors are not updated with a global-to-local and the user can not do it because it is locked
325: PetscCall(VecLockReadPush(v));
326: PetscCall(VecRestoreArrayRead(pvec, &array));
327: } else {
328: PetscScalar *array;
329: PetscCall(VecGetArray(pvec, &array));
330: PetscCall(VecPlaceArray(v, array + nlocal));
331: PetscCall(VecRestoreArray(pvec, &array));
332: }
333: vecs[wnum++] = v;
334: }
336: nlocal += link->nlocal;
337: }
339: PetscFunctionReturn(PETSC_SUCCESS);
340: }
342: /*@C
343: DMCompositeRestoreAccess - Returns the vectors obtained with `DMCompositeGetAccess()`
344: representation.
346: Collective
348: Input Parameters:
349: + dm - the `DMCOMPOSITE` object
350: . gvec - the global vector
351: - Vec* ... - the individual parallel vectors, NULL for those that are not needed
353: Level: advanced
355: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
356: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeScatter()`,
357: `DMCompositeRestoreAccess()`, `DMCompositeGetAccess()`
358: @*/
359: PetscErrorCode DMCompositeRestoreAccess(DM dm, Vec gvec, ...)
360: {
361: va_list Argp;
362: struct DMCompositeLink *next;
363: DM_Composite *com = (DM_Composite *)dm->data;
364: PetscInt readonly;
365: PetscBool flg;
367: PetscFunctionBegin;
370: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
371: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
372: next = com->next;
373: if (!com->setup) PetscCall(DMSetUp(dm));
375: PetscCall(VecLockGet(gvec, &readonly));
376: /* loop over packed objects, handling one at at time */
377: va_start(Argp, gvec);
378: while (next) {
379: Vec *vec;
380: vec = va_arg(Argp, Vec *);
381: if (vec) {
382: PetscCall(VecResetArray(*vec));
383: if (readonly) PetscCall(VecLockReadPop(*vec));
384: PetscCall(DMRestoreGlobalVector(next->dm, vec));
385: }
386: next = next->next;
387: }
388: va_end(Argp);
389: PetscFunctionReturn(PETSC_SUCCESS);
390: }
392: /*@C
393: DMCompositeRestoreAccessArray - Returns the vectors obtained with `DMCompositeGetAccessArray()`
395: Collective
397: Input Parameters:
398: + dm - the `DMCOMPOSITE` object
399: . pvec - packed vector
400: . nwanted - number of vectors wanted
401: . wanted - sorted array of vectors wanted, or NULL to get all vectors
402: - vecs - array of global vectors to return
404: Level: advanced
406: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeRestoreAccess()`, `DMCompositeRestoreEntries()`, `DMCompositeScatter()`, `DMCompositeGather()`
407: @*/
408: PetscErrorCode DMCompositeRestoreAccessArray(DM dm, Vec pvec, PetscInt nwanted, const PetscInt *wanted, Vec *vecs)
409: {
410: struct DMCompositeLink *link;
411: PetscInt i, wnum;
412: DM_Composite *com = (DM_Composite *)dm->data;
413: PetscInt readonly;
414: PetscBool flg;
416: PetscFunctionBegin;
419: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
420: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
421: if (!com->setup) PetscCall(DMSetUp(dm));
423: PetscCall(VecLockGet(pvec, &readonly));
424: for (i = 0, wnum = 0, link = com->next; link && wnum < nwanted; i++, link = link->next) {
425: if (!wanted || i == wanted[wnum]) {
426: PetscCall(VecResetArray(vecs[wnum]));
427: if (readonly) PetscCall(VecLockReadPop(vecs[wnum]));
428: PetscCall(DMRestoreGlobalVector(link->dm, &vecs[wnum]));
429: wnum++;
430: }
431: }
432: PetscFunctionReturn(PETSC_SUCCESS);
433: }
435: /*@C
436: DMCompositeRestoreLocalAccessArray - Returns the vectors obtained with `DMCompositeGetLocalAccessArray()`.
438: Collective
440: Input Parameters:
441: + dm - the `DMCOMPOSITE` object
442: . pvec - packed vector
443: . nwanted - number of vectors wanted
444: . wanted - sorted array of vectors wanted, or NULL to restore all vectors
445: - vecs - array of local vectors to return
447: Level: advanced
449: Note:
450: nwanted and wanted must match the values given to `DMCompositeGetLocalAccessArray()`
451: otherwise the call will fail.
453: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeGetLocalAccessArray()`, `DMCompositeRestoreAccessArray()`,
454: `DMCompositeRestoreAccess()`, `DMCompositeRestoreEntries()`,
455: `DMCompositeScatter()`, `DMCompositeGather()`
456: @*/
457: PetscErrorCode DMCompositeRestoreLocalAccessArray(DM dm, Vec pvec, PetscInt nwanted, const PetscInt *wanted, Vec *vecs)
458: {
459: struct DMCompositeLink *link;
460: PetscInt i, wnum;
461: DM_Composite *com = (DM_Composite *)dm->data;
462: PetscInt readonly;
463: PetscBool flg;
465: PetscFunctionBegin;
468: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
469: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
470: if (!com->setup) PetscCall(DMSetUp(dm));
472: PetscCall(VecLockGet(pvec, &readonly));
473: for (i = 0, wnum = 0, link = com->next; link && wnum < nwanted; i++, link = link->next) {
474: if (!wanted || i == wanted[wnum]) {
475: PetscCall(VecResetArray(vecs[wnum]));
476: if (readonly) PetscCall(VecLockReadPop(vecs[wnum]));
477: PetscCall(DMRestoreLocalVector(link->dm, &vecs[wnum]));
478: wnum++;
479: }
480: }
481: PetscFunctionReturn(PETSC_SUCCESS);
482: }
484: /*@C
485: DMCompositeScatter - Scatters from a global packed vector into its individual local vectors
487: Collective
489: Input Parameters:
490: + dm - the `DMCOMPOSITE` object
491: . gvec - the global vector
492: - Vec ... - the individual sequential vectors, NULL for those that are not needed
494: Level: advanced
496: Note:
497: `DMCompositeScatterArray()` is a non-variadic alternative that is often more convenient for library callers and is
498: accessible from Fortran.
500: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
501: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
502: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
503: `DMCompositeScatterArray()`
504: @*/
505: PetscErrorCode DMCompositeScatter(DM dm, Vec gvec, ...)
506: {
507: va_list Argp;
508: struct DMCompositeLink *next;
509: PETSC_UNUSED PetscInt cnt;
510: DM_Composite *com = (DM_Composite *)dm->data;
511: PetscBool flg;
513: PetscFunctionBegin;
516: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
517: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
518: if (!com->setup) PetscCall(DMSetUp(dm));
520: /* loop over packed objects, handling one at at time */
521: va_start(Argp, gvec);
522: for (cnt = 3, next = com->next; next; cnt++, next = next->next) {
523: Vec local;
524: local = va_arg(Argp, Vec);
525: if (local) {
526: Vec global;
527: const PetscScalar *array;
529: PetscCall(DMGetGlobalVector(next->dm, &global));
530: PetscCall(VecGetArrayRead(gvec, &array));
531: PetscCall(VecPlaceArray(global, array + next->rstart));
532: PetscCall(DMGlobalToLocalBegin(next->dm, global, INSERT_VALUES, local));
533: PetscCall(DMGlobalToLocalEnd(next->dm, global, INSERT_VALUES, local));
534: PetscCall(VecRestoreArrayRead(gvec, &array));
535: PetscCall(VecResetArray(global));
536: PetscCall(DMRestoreGlobalVector(next->dm, &global));
537: }
538: }
539: va_end(Argp);
540: PetscFunctionReturn(PETSC_SUCCESS);
541: }
543: /*@
544: DMCompositeScatterArray - Scatters from a global packed vector into its individual local vectors
546: Collective
548: Input Parameters:
549: + dm - the `DMCOMPOSITE` object
550: . gvec - the global vector
551: - lvecs - array of local vectors, NULL for any that are not needed
553: Level: advanced
555: Note:
556: This is a non-variadic alternative to `DMCompositeScatter()`
558: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`
559: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
560: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
561: @*/
562: PetscErrorCode DMCompositeScatterArray(DM dm, Vec gvec, Vec *lvecs)
563: {
564: struct DMCompositeLink *next;
565: PetscInt i;
566: DM_Composite *com = (DM_Composite *)dm->data;
567: PetscBool flg;
569: PetscFunctionBegin;
572: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
573: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
574: if (!com->setup) PetscCall(DMSetUp(dm));
576: /* loop over packed objects, handling one at at time */
577: for (i = 0, next = com->next; next; next = next->next, i++) {
578: if (lvecs[i]) {
579: Vec global;
580: const PetscScalar *array;
582: PetscCall(DMGetGlobalVector(next->dm, &global));
583: PetscCall(VecGetArrayRead(gvec, &array));
584: PetscCall(VecPlaceArray(global, (PetscScalar *)array + next->rstart));
585: PetscCall(DMGlobalToLocalBegin(next->dm, global, INSERT_VALUES, lvecs[i]));
586: PetscCall(DMGlobalToLocalEnd(next->dm, global, INSERT_VALUES, lvecs[i]));
587: PetscCall(VecRestoreArrayRead(gvec, &array));
588: PetscCall(VecResetArray(global));
589: PetscCall(DMRestoreGlobalVector(next->dm, &global));
590: }
591: }
592: PetscFunctionReturn(PETSC_SUCCESS);
593: }
595: /*@C
596: DMCompositeGather - Gathers into a global packed vector from its individual local vectors
598: Collective
600: Input Parameters:
601: + dm - the `DMCOMPOSITE` object
602: . gvec - the global vector
603: . imode - `INSERT_VALUES` or `ADD_VALUES`
604: - Vec ... - the individual sequential vectors, NULL for any that are not needed
606: Level: advanced
608: Fortran Note:
609: Fortran users should use `DMCompositeGatherArray()`
611: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
612: `DMCompositeScatter()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
613: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
614: @*/
615: PetscErrorCode DMCompositeGather(DM dm, InsertMode imode, Vec gvec, ...)
616: {
617: va_list Argp;
618: struct DMCompositeLink *next;
619: DM_Composite *com = (DM_Composite *)dm->data;
620: PETSC_UNUSED PetscInt cnt;
621: PetscBool flg;
623: PetscFunctionBegin;
626: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
627: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
628: if (!com->setup) PetscCall(DMSetUp(dm));
630: /* loop over packed objects, handling one at at time */
631: va_start(Argp, gvec);
632: for (cnt = 3, next = com->next; next; cnt++, next = next->next) {
633: Vec local;
634: local = va_arg(Argp, Vec);
635: if (local) {
636: PetscScalar *array;
637: Vec global;
639: PetscCall(DMGetGlobalVector(next->dm, &global));
640: PetscCall(VecGetArray(gvec, &array));
641: PetscCall(VecPlaceArray(global, array + next->rstart));
642: PetscCall(DMLocalToGlobalBegin(next->dm, local, imode, global));
643: PetscCall(DMLocalToGlobalEnd(next->dm, local, imode, global));
644: PetscCall(VecRestoreArray(gvec, &array));
645: PetscCall(VecResetArray(global));
646: PetscCall(DMRestoreGlobalVector(next->dm, &global));
647: }
648: }
649: va_end(Argp);
650: PetscFunctionReturn(PETSC_SUCCESS);
651: }
653: /*@
654: DMCompositeGatherArray - Gathers into a global packed vector from its individual local vectors
656: Collective
658: Input Parameters:
659: + dm - the `DMCOMPOSITE` object
660: . gvec - the global vector
661: . imode - `INSERT_VALUES` or `ADD_VALUES`
662: - lvecs - the individual sequential vectors, NULL for any that are not needed
664: Level: advanced
666: Note:
667: This is a non-variadic alternative to `DMCompositeGather()`.
669: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
670: `DMCompositeScatter()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
671: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`,
672: @*/
673: PetscErrorCode DMCompositeGatherArray(DM dm, InsertMode imode, Vec gvec, Vec *lvecs)
674: {
675: struct DMCompositeLink *next;
676: DM_Composite *com = (DM_Composite *)dm->data;
677: PetscInt i;
678: PetscBool flg;
680: PetscFunctionBegin;
683: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
684: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
685: if (!com->setup) PetscCall(DMSetUp(dm));
687: /* loop over packed objects, handling one at at time */
688: for (next = com->next, i = 0; next; next = next->next, i++) {
689: if (lvecs[i]) {
690: PetscScalar *array;
691: Vec global;
693: PetscCall(DMGetGlobalVector(next->dm, &global));
694: PetscCall(VecGetArray(gvec, &array));
695: PetscCall(VecPlaceArray(global, array + next->rstart));
696: PetscCall(DMLocalToGlobalBegin(next->dm, lvecs[i], imode, global));
697: PetscCall(DMLocalToGlobalEnd(next->dm, lvecs[i], imode, global));
698: PetscCall(VecRestoreArray(gvec, &array));
699: PetscCall(VecResetArray(global));
700: PetscCall(DMRestoreGlobalVector(next->dm, &global));
701: }
702: }
703: PetscFunctionReturn(PETSC_SUCCESS);
704: }
706: /*@
707: DMCompositeAddDM - adds a `DM` vector to a `DMCOMPOSITE`
709: Collective
711: Input Parameters:
712: + dmc - the `DMCOMPOSITE` object
713: - dm - the `DM` object
715: Level: advanced
717: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeGather()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
718: `DMCompositeScatter()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
719: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
720: @*/
721: PetscErrorCode DMCompositeAddDM(DM dmc, DM dm)
722: {
723: PetscInt n, nlocal;
724: struct DMCompositeLink *mine, *next;
725: Vec global, local;
726: DM_Composite *com = (DM_Composite *)dmc->data;
727: PetscBool flg;
729: PetscFunctionBegin;
732: PetscCall(PetscObjectTypeCompare((PetscObject)dmc, DMCOMPOSITE, &flg));
733: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
734: next = com->next;
735: PetscCheck(!com->setup, PetscObjectComm((PetscObject)dmc), PETSC_ERR_ARG_WRONGSTATE, "Cannot add a DM once you have used the DMComposite");
737: /* create new link */
738: PetscCall(PetscNew(&mine));
739: PetscCall(PetscObjectReference((PetscObject)dm));
740: PetscCall(DMGetGlobalVector(dm, &global));
741: PetscCall(VecGetLocalSize(global, &n));
742: PetscCall(DMRestoreGlobalVector(dm, &global));
743: PetscCall(DMGetLocalVector(dm, &local));
744: PetscCall(VecGetSize(local, &nlocal));
745: PetscCall(DMRestoreLocalVector(dm, &local));
747: mine->n = n;
748: mine->nlocal = nlocal;
749: mine->dm = dm;
750: mine->next = NULL;
751: com->n += n;
752: com->nghost += nlocal;
754: /* add to end of list */
755: if (!next) com->next = mine;
756: else {
757: while (next->next) next = next->next;
758: next->next = mine;
759: }
760: com->nDM++;
761: com->nmine++;
762: PetscFunctionReturn(PETSC_SUCCESS);
763: }
765: #include <petscdraw.h>
766: PETSC_EXTERN PetscErrorCode VecView_MPI(Vec, PetscViewer);
767: PetscErrorCode VecView_DMComposite(Vec gvec, PetscViewer viewer)
768: {
769: DM dm;
770: struct DMCompositeLink *next;
771: PetscBool isdraw;
772: DM_Composite *com;
774: PetscFunctionBegin;
775: PetscCall(VecGetDM(gvec, &dm));
776: PetscCheck(dm, PetscObjectComm((PetscObject)gvec), PETSC_ERR_ARG_WRONG, "Vector not generated from a DMComposite");
777: com = (DM_Composite *)dm->data;
778: next = com->next;
780: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERDRAW, &isdraw));
781: if (!isdraw) {
782: /* do I really want to call this? */
783: PetscCall(VecView_MPI(gvec, viewer));
784: } else {
785: PetscInt cnt = 0;
787: /* loop over packed objects, handling one at at time */
788: while (next) {
789: Vec vec;
790: const PetscScalar *array;
791: PetscInt bs;
793: /* Should use VecGetSubVector() eventually, but would need to forward the DM for that to work */
794: PetscCall(DMGetGlobalVector(next->dm, &vec));
795: PetscCall(VecGetArrayRead(gvec, &array));
796: PetscCall(VecPlaceArray(vec, (PetscScalar *)array + next->rstart));
797: PetscCall(VecRestoreArrayRead(gvec, &array));
798: PetscCall(VecView(vec, viewer));
799: PetscCall(VecResetArray(vec));
800: PetscCall(VecGetBlockSize(vec, &bs));
801: PetscCall(DMRestoreGlobalVector(next->dm, &vec));
802: PetscCall(PetscViewerDrawBaseAdd(viewer, bs));
803: cnt += bs;
804: next = next->next;
805: }
806: PetscCall(PetscViewerDrawBaseAdd(viewer, -cnt));
807: }
808: PetscFunctionReturn(PETSC_SUCCESS);
809: }
811: PetscErrorCode DMCreateGlobalVector_Composite(DM dm, Vec *gvec)
812: {
813: DM_Composite *com = (DM_Composite *)dm->data;
815: PetscFunctionBegin;
817: PetscCall(DMSetFromOptions(dm));
818: PetscCall(DMSetUp(dm));
819: PetscCall(VecCreate(PetscObjectComm((PetscObject)dm), gvec));
820: PetscCall(VecSetType(*gvec, dm->vectype));
821: PetscCall(VecSetSizes(*gvec, com->n, com->N));
822: PetscCall(VecSetDM(*gvec, dm));
823: PetscCall(VecSetOperation(*gvec, VECOP_VIEW, (void (*)(void))VecView_DMComposite));
824: PetscFunctionReturn(PETSC_SUCCESS);
825: }
827: PetscErrorCode DMCreateLocalVector_Composite(DM dm, Vec *lvec)
828: {
829: DM_Composite *com = (DM_Composite *)dm->data;
831: PetscFunctionBegin;
833: if (!com->setup) {
834: PetscCall(DMSetFromOptions(dm));
835: PetscCall(DMSetUp(dm));
836: }
837: PetscCall(VecCreate(PETSC_COMM_SELF, lvec));
838: PetscCall(VecSetType(*lvec, dm->vectype));
839: PetscCall(VecSetSizes(*lvec, com->nghost, PETSC_DECIDE));
840: PetscCall(VecSetDM(*lvec, dm));
841: PetscFunctionReturn(PETSC_SUCCESS);
842: }
844: /*@C
845: DMCompositeGetISLocalToGlobalMappings - gets an `ISLocalToGlobalMapping` for each `DM` in the `DMCOMPOSITE`, maps to the composite global space
847: Collective; No Fortran Support
849: Input Parameter:
850: . dm - the `DMCOMPOSITE` object
852: Output Parameter:
853: . ltogs - the individual mappings for each packed vector. Note that this includes
854: all the ghost points that individual ghosted `DMDA` may have.
856: Level: advanced
858: Note:
859: Each entry of ltogs should be destroyed with `ISLocalToGlobalMappingDestroy()`, the ltogs array should be freed with `PetscFree()`.
861: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
862: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetAccess()`, `DMCompositeScatter()`,
863: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
864: @*/
865: PetscErrorCode DMCompositeGetISLocalToGlobalMappings(DM dm, ISLocalToGlobalMapping **ltogs)
866: {
867: PetscInt i, *idx, n, cnt;
868: struct DMCompositeLink *next;
869: PetscMPIInt rank;
870: DM_Composite *com = (DM_Composite *)dm->data;
871: PetscBool flg;
873: PetscFunctionBegin;
875: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
876: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
877: PetscCall(DMSetUp(dm));
878: PetscCall(PetscMalloc1(com->nDM, ltogs));
879: next = com->next;
880: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
882: /* loop over packed objects, handling one at at time */
883: cnt = 0;
884: while (next) {
885: ISLocalToGlobalMapping ltog;
886: PetscMPIInt size;
887: const PetscInt *suboff, *indices;
888: Vec global;
890: /* Get sub-DM global indices for each local dof */
891: PetscCall(DMGetLocalToGlobalMapping(next->dm, <og));
892: PetscCall(ISLocalToGlobalMappingGetSize(ltog, &n));
893: PetscCall(ISLocalToGlobalMappingGetIndices(ltog, &indices));
894: PetscCall(PetscMalloc1(n, &idx));
896: /* Get the offsets for the sub-DM global vector */
897: PetscCall(DMGetGlobalVector(next->dm, &global));
898: PetscCall(VecGetOwnershipRanges(global, &suboff));
899: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)global), &size));
901: /* Shift the sub-DM definition of the global space to the composite global space */
902: for (i = 0; i < n; i++) {
903: PetscInt subi = indices[i], lo = 0, hi = size, t;
904: /* There's no consensus on what a negative index means,
905: except for skipping when setting the values in vectors and matrices */
906: if (subi < 0) {
907: idx[i] = subi - next->grstarts[rank];
908: continue;
909: }
910: /* Binary search to find which rank owns subi */
911: while (hi - lo > 1) {
912: t = lo + (hi - lo) / 2;
913: if (suboff[t] > subi) hi = t;
914: else lo = t;
915: }
916: idx[i] = subi - suboff[lo] + next->grstarts[lo];
917: }
918: PetscCall(ISLocalToGlobalMappingRestoreIndices(ltog, &indices));
919: PetscCall(ISLocalToGlobalMappingCreate(PetscObjectComm((PetscObject)dm), 1, n, idx, PETSC_OWN_POINTER, &(*ltogs)[cnt]));
920: PetscCall(DMRestoreGlobalVector(next->dm, &global));
921: next = next->next;
922: cnt++;
923: }
924: PetscFunctionReturn(PETSC_SUCCESS);
925: }
927: /*@C
928: DMCompositeGetLocalISs - Gets index sets for each component of a composite local vector
930: Not Collective; No Fortran Support
932: Input Parameter:
933: . dm - the `DMCOMPOSITE`
935: Output Parameter:
936: . is - array of serial index sets for each each component of the `DMCOMPOSITE`
938: Level: intermediate
940: Notes:
941: At present, a composite local vector does not normally exist. This function is used to provide index sets for
942: `MatGetLocalSubMatrix()`. In the future, the scatters for each entry in the `DMCOMPOSITE` may be be merged into a single
943: scatter to a composite local vector. The user should not typically need to know which is being done.
945: To get the composite global indices at all local points (including ghosts), use `DMCompositeGetISLocalToGlobalMappings()`.
947: To get index sets for pieces of the composite global vector, use `DMCompositeGetGlobalISs()`.
949: Each returned `IS` should be destroyed with `ISDestroy()`, the array should be freed with `PetscFree()`.
951: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeGetGlobalISs()`, `DMCompositeGetISLocalToGlobalMappings()`, `MatGetLocalSubMatrix()`, `MatCreateLocalRef()`
952: @*/
953: PetscErrorCode DMCompositeGetLocalISs(DM dm, IS **is)
954: {
955: DM_Composite *com = (DM_Composite *)dm->data;
956: struct DMCompositeLink *link;
957: PetscInt cnt, start;
958: PetscBool flg;
960: PetscFunctionBegin;
963: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
964: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
965: PetscCall(PetscMalloc1(com->nmine, is));
966: for (cnt = 0, start = 0, link = com->next; link; start += link->nlocal, cnt++, link = link->next) {
967: PetscInt bs;
968: PetscCall(ISCreateStride(PETSC_COMM_SELF, link->nlocal, start, 1, &(*is)[cnt]));
969: PetscCall(DMGetBlockSize(link->dm, &bs));
970: PetscCall(ISSetBlockSize((*is)[cnt], bs));
971: }
972: PetscFunctionReturn(PETSC_SUCCESS);
973: }
975: /*@C
976: DMCompositeGetGlobalISs - Gets the index sets for each composed object in a `DMCOMPOSITE`
978: Collective
980: Input Parameter:
981: . dm - the `DMCOMPOSITE` object
983: Output Parameter:
984: . is - the array of index sets
986: Level: advanced
988: Notes:
989: The is entries should be destroyed with `ISDestroy()`, the is array should be freed with `PetscFree()`
991: These could be used to extract a subset of vector entries for a "multi-physics" preconditioner
993: Use `DMCompositeGetLocalISs()` for index sets in the packed local numbering, and
994: `DMCompositeGetISLocalToGlobalMappings()` for to map local sub-`DM` (including ghost) indices to packed global
995: indices.
997: Fortran Note:
998: The output argument 'is' must be an allocated array of sufficient length, which can be learned using `DMCompositeGetNumberDM()`.
1000: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
1001: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetAccess()`, `DMCompositeScatter()`,
1002: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
1003: @*/
1004: PetscErrorCode DMCompositeGetGlobalISs(DM dm, IS *is[])
1005: {
1006: PetscInt cnt = 0;
1007: struct DMCompositeLink *next;
1008: PetscMPIInt rank;
1009: DM_Composite *com = (DM_Composite *)dm->data;
1010: PetscBool flg;
1012: PetscFunctionBegin;
1014: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1015: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1016: PetscCheck(dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Must call DMSetUp() before");
1017: PetscCall(PetscMalloc1(com->nDM, is));
1018: next = com->next;
1019: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
1021: /* loop over packed objects, handling one at at time */
1022: while (next) {
1023: PetscDS prob;
1025: PetscCall(ISCreateStride(PetscObjectComm((PetscObject)dm), next->n, next->grstart, 1, &(*is)[cnt]));
1026: PetscCall(DMGetDS(dm, &prob));
1027: if (prob) {
1028: MatNullSpace space;
1029: Mat pmat;
1030: PetscObject disc;
1031: PetscInt Nf;
1033: PetscCall(PetscDSGetNumFields(prob, &Nf));
1034: if (cnt < Nf) {
1035: PetscCall(PetscDSGetDiscretization(prob, cnt, &disc));
1036: PetscCall(PetscObjectQuery(disc, "nullspace", (PetscObject *)&space));
1037: if (space) PetscCall(PetscObjectCompose((PetscObject)(*is)[cnt], "nullspace", (PetscObject)space));
1038: PetscCall(PetscObjectQuery(disc, "nearnullspace", (PetscObject *)&space));
1039: if (space) PetscCall(PetscObjectCompose((PetscObject)(*is)[cnt], "nearnullspace", (PetscObject)space));
1040: PetscCall(PetscObjectQuery(disc, "pmat", (PetscObject *)&pmat));
1041: if (pmat) PetscCall(PetscObjectCompose((PetscObject)(*is)[cnt], "pmat", (PetscObject)pmat));
1042: }
1043: }
1044: cnt++;
1045: next = next->next;
1046: }
1047: PetscFunctionReturn(PETSC_SUCCESS);
1048: }
1050: PetscErrorCode DMCreateFieldIS_Composite(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
1051: {
1052: PetscInt nDM;
1053: DM *dms;
1054: PetscInt i;
1056: PetscFunctionBegin;
1057: PetscCall(DMCompositeGetNumberDM(dm, &nDM));
1058: if (numFields) *numFields = nDM;
1059: PetscCall(DMCompositeGetGlobalISs(dm, fields));
1060: if (fieldNames) {
1061: PetscCall(PetscMalloc1(nDM, &dms));
1062: PetscCall(PetscMalloc1(nDM, fieldNames));
1063: PetscCall(DMCompositeGetEntriesArray(dm, dms));
1064: for (i = 0; i < nDM; i++) {
1065: char buf[256];
1066: const char *splitname;
1068: /* Split naming precedence: object name, prefix, number */
1069: splitname = ((PetscObject)dm)->name;
1070: if (!splitname) {
1071: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)dms[i], &splitname));
1072: if (splitname) {
1073: size_t len;
1074: PetscCall(PetscStrncpy(buf, splitname, sizeof(buf)));
1075: buf[sizeof(buf) - 1] = 0;
1076: PetscCall(PetscStrlen(buf, &len));
1077: if (buf[len - 1] == '_') buf[len - 1] = 0; /* Remove trailing underscore if it was used */
1078: splitname = buf;
1079: }
1080: }
1081: if (!splitname) {
1082: PetscCall(PetscSNPrintf(buf, sizeof(buf), "%" PetscInt_FMT, i));
1083: splitname = buf;
1084: }
1085: PetscCall(PetscStrallocpy(splitname, &(*fieldNames)[i]));
1086: }
1087: PetscCall(PetscFree(dms));
1088: }
1089: PetscFunctionReturn(PETSC_SUCCESS);
1090: }
1092: /*
1093: This could take over from DMCreateFieldIS(), as it is more general,
1094: making DMCreateFieldIS() a special case -- calling with dmlist == NULL;
1095: At this point it's probably best to be less intrusive, however.
1096: */
1097: PetscErrorCode DMCreateFieldDecomposition_Composite(DM dm, PetscInt *len, char ***namelist, IS **islist, DM **dmlist)
1098: {
1099: PetscInt nDM;
1100: PetscInt i;
1102: PetscFunctionBegin;
1103: PetscCall(DMCreateFieldIS_Composite(dm, len, namelist, islist));
1104: if (dmlist) {
1105: PetscCall(DMCompositeGetNumberDM(dm, &nDM));
1106: PetscCall(PetscMalloc1(nDM, dmlist));
1107: PetscCall(DMCompositeGetEntriesArray(dm, *dmlist));
1108: for (i = 0; i < nDM; i++) PetscCall(PetscObjectReference((PetscObject)((*dmlist)[i])));
1109: }
1110: PetscFunctionReturn(PETSC_SUCCESS);
1111: }
1113: /* -------------------------------------------------------------------------------------*/
1114: /*@C
1115: DMCompositeGetLocalVectors - Gets local vectors for each part of a `DMCOMPOSITE`
1116: Use `DMCompositeRestoreLocalVectors()` to return them.
1118: Not Collective; No Fortran Support
1120: Input Parameter:
1121: . dm - the `DMCOMPOSITE` object
1123: Output Parameter:
1124: . Vec ... - the individual sequential `Vec`s
1126: Level: advanced
1128: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
1129: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
1130: `DMCompositeRestoreLocalVectors()`, `DMCompositeScatter()`, `DMCompositeGetEntries()`
1131: @*/
1132: PetscErrorCode DMCompositeGetLocalVectors(DM dm, ...)
1133: {
1134: va_list Argp;
1135: struct DMCompositeLink *next;
1136: DM_Composite *com = (DM_Composite *)dm->data;
1137: PetscBool flg;
1139: PetscFunctionBegin;
1141: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1142: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1143: next = com->next;
1144: /* loop over packed objects, handling one at at time */
1145: va_start(Argp, dm);
1146: while (next) {
1147: Vec *vec;
1148: vec = va_arg(Argp, Vec *);
1149: if (vec) PetscCall(DMGetLocalVector(next->dm, vec));
1150: next = next->next;
1151: }
1152: va_end(Argp);
1153: PetscFunctionReturn(PETSC_SUCCESS);
1154: }
1156: /*@C
1157: DMCompositeRestoreLocalVectors - Restores local vectors for each part of a `DMCOMPOSITE`
1159: Not Collective; No Fortran Support
1161: Input Parameter:
1162: . dm - the `DMCOMPOSITE` object
1164: Output Parameter:
1165: . Vec ... - the individual sequential `Vec`
1167: Level: advanced
1169: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
1170: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
1171: `DMCompositeGetLocalVectors()`, `DMCompositeScatter()`, `DMCompositeGetEntries()`
1172: @*/
1173: PetscErrorCode DMCompositeRestoreLocalVectors(DM dm, ...)
1174: {
1175: va_list Argp;
1176: struct DMCompositeLink *next;
1177: DM_Composite *com = (DM_Composite *)dm->data;
1178: PetscBool flg;
1180: PetscFunctionBegin;
1182: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1183: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1184: next = com->next;
1185: /* loop over packed objects, handling one at at time */
1186: va_start(Argp, dm);
1187: while (next) {
1188: Vec *vec;
1189: vec = va_arg(Argp, Vec *);
1190: if (vec) PetscCall(DMRestoreLocalVector(next->dm, vec));
1191: next = next->next;
1192: }
1193: va_end(Argp);
1194: PetscFunctionReturn(PETSC_SUCCESS);
1195: }
1197: /* -------------------------------------------------------------------------------------*/
1198: /*@C
1199: DMCompositeGetEntries - Gets the `DM` for each entry in a `DMCOMPOSITE`.
1201: Not Collective
1203: Input Parameter:
1204: . dm - the `DMCOMPOSITE` object
1206: Output Parameter:
1207: . DM ... - the individual entries `DM`
1209: Level: advanced
1211: Fortran Note:
1212: Available as `DMCompositeGetEntries()` for one output `DM`, DMCompositeGetEntries2() for 2, etc
1214: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`, `DMCompositeGetEntriesArray()`
1215: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
1216: `DMCompositeRestoreLocalVectors()`, `DMCompositeGetLocalVectors()`, `DMCompositeScatter()`,
1217: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`
1218: @*/
1219: PetscErrorCode DMCompositeGetEntries(DM dm, ...)
1220: {
1221: va_list Argp;
1222: struct DMCompositeLink *next;
1223: DM_Composite *com = (DM_Composite *)dm->data;
1224: PetscBool flg;
1226: PetscFunctionBegin;
1228: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1229: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1230: next = com->next;
1231: /* loop over packed objects, handling one at at time */
1232: va_start(Argp, dm);
1233: while (next) {
1234: DM *dmn;
1235: dmn = va_arg(Argp, DM *);
1236: if (dmn) *dmn = next->dm;
1237: next = next->next;
1238: }
1239: va_end(Argp);
1240: PetscFunctionReturn(PETSC_SUCCESS);
1241: }
1243: /*@C
1244: DMCompositeGetEntriesArray - Gets the DM for each entry in a `DMCOMPOSITE`
1246: Not Collective
1248: Input Parameter:
1249: . dm - the `DMCOMPOSITE` object
1251: Output Parameter:
1252: . dms - array of sufficient length (see `DMCompositeGetNumberDM()`) to hold the individual `DM`
1254: Level: advanced
1256: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`, `DMCompositeGetEntries()`
1257: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
1258: `DMCompositeRestoreLocalVectors()`, `DMCompositeGetLocalVectors()`, `DMCompositeScatter()`,
1259: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`
1260: @*/
1261: PetscErrorCode DMCompositeGetEntriesArray(DM dm, DM dms[])
1262: {
1263: struct DMCompositeLink *next;
1264: DM_Composite *com = (DM_Composite *)dm->data;
1265: PetscInt i;
1266: PetscBool flg;
1268: PetscFunctionBegin;
1270: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1271: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1272: /* loop over packed objects, handling one at at time */
1273: for (next = com->next, i = 0; next; next = next->next, i++) dms[i] = next->dm;
1274: PetscFunctionReturn(PETSC_SUCCESS);
1275: }
1277: typedef struct {
1278: DM dm;
1279: PetscViewer *subv;
1280: Vec *vecs;
1281: } GLVisViewerCtx;
1283: static PetscErrorCode DestroyGLVisViewerCtx_Private(void *vctx)
1284: {
1285: GLVisViewerCtx *ctx = (GLVisViewerCtx *)vctx;
1286: PetscInt i, n;
1288: PetscFunctionBegin;
1289: PetscCall(DMCompositeGetNumberDM(ctx->dm, &n));
1290: for (i = 0; i < n; i++) PetscCall(PetscViewerDestroy(&ctx->subv[i]));
1291: PetscCall(PetscFree2(ctx->subv, ctx->vecs));
1292: PetscCall(DMDestroy(&ctx->dm));
1293: PetscCall(PetscFree(ctx));
1294: PetscFunctionReturn(PETSC_SUCCESS);
1295: }
1297: static PetscErrorCode DMCompositeSampleGLVisFields_Private(PetscObject oX, PetscInt nf, PetscObject oXfield[], void *vctx)
1298: {
1299: Vec X = (Vec)oX;
1300: GLVisViewerCtx *ctx = (GLVisViewerCtx *)vctx;
1301: PetscInt i, n, cumf;
1303: PetscFunctionBegin;
1304: PetscCall(DMCompositeGetNumberDM(ctx->dm, &n));
1305: PetscCall(DMCompositeGetAccessArray(ctx->dm, X, n, NULL, ctx->vecs));
1306: for (i = 0, cumf = 0; i < n; i++) {
1307: PetscErrorCode (*g2l)(PetscObject, PetscInt, PetscObject[], void *);
1308: void *fctx;
1309: PetscInt nfi;
1311: PetscCall(PetscViewerGLVisGetFields_Private(ctx->subv[i], &nfi, NULL, NULL, &g2l, NULL, &fctx));
1312: if (!nfi) continue;
1313: if (g2l) PetscCall((*g2l)((PetscObject)ctx->vecs[i], nfi, oXfield + cumf, fctx));
1314: else PetscCall(VecCopy(ctx->vecs[i], (Vec)(oXfield[cumf])));
1315: cumf += nfi;
1316: }
1317: PetscCall(DMCompositeRestoreAccessArray(ctx->dm, X, n, NULL, ctx->vecs));
1318: PetscFunctionReturn(PETSC_SUCCESS);
1319: }
1321: static PetscErrorCode DMSetUpGLVisViewer_Composite(PetscObject odm, PetscViewer viewer)
1322: {
1323: DM dm = (DM)odm, *dms;
1324: Vec *Ufds;
1325: GLVisViewerCtx *ctx;
1326: PetscInt i, n, tnf, *sdim;
1327: char **fecs;
1329: PetscFunctionBegin;
1330: PetscCall(PetscNew(&ctx));
1331: PetscCall(PetscObjectReference((PetscObject)dm));
1332: ctx->dm = dm;
1333: PetscCall(DMCompositeGetNumberDM(dm, &n));
1334: PetscCall(PetscMalloc1(n, &dms));
1335: PetscCall(DMCompositeGetEntriesArray(dm, dms));
1336: PetscCall(PetscMalloc2(n, &ctx->subv, n, &ctx->vecs));
1337: for (i = 0, tnf = 0; i < n; i++) {
1338: PetscInt nf;
1340: PetscCall(PetscViewerCreate(PetscObjectComm(odm), &ctx->subv[i]));
1341: PetscCall(PetscViewerSetType(ctx->subv[i], PETSCVIEWERGLVIS));
1342: PetscCall(PetscViewerGLVisSetDM_Private(ctx->subv[i], (PetscObject)dms[i]));
1343: PetscCall(PetscViewerGLVisGetFields_Private(ctx->subv[i], &nf, NULL, NULL, NULL, NULL, NULL));
1344: tnf += nf;
1345: }
1346: PetscCall(PetscFree(dms));
1347: PetscCall(PetscMalloc3(tnf, &fecs, tnf, &sdim, tnf, &Ufds));
1348: for (i = 0, tnf = 0; i < n; i++) {
1349: PetscInt *sd, nf, f;
1350: const char **fec;
1351: Vec *Uf;
1353: PetscCall(PetscViewerGLVisGetFields_Private(ctx->subv[i], &nf, &fec, &sd, NULL, (PetscObject **)&Uf, NULL));
1354: for (f = 0; f < nf; f++) {
1355: PetscCall(PetscStrallocpy(fec[f], &fecs[tnf + f]));
1356: Ufds[tnf + f] = Uf[f];
1357: sdim[tnf + f] = sd[f];
1358: }
1359: tnf += nf;
1360: }
1361: PetscCall(PetscViewerGLVisSetFields(viewer, tnf, (const char **)fecs, sdim, DMCompositeSampleGLVisFields_Private, (PetscObject *)Ufds, ctx, DestroyGLVisViewerCtx_Private));
1362: for (i = 0; i < tnf; i++) PetscCall(PetscFree(fecs[i]));
1363: PetscCall(PetscFree3(fecs, sdim, Ufds));
1364: PetscFunctionReturn(PETSC_SUCCESS);
1365: }
1367: PetscErrorCode DMRefine_Composite(DM dmi, MPI_Comm comm, DM *fine)
1368: {
1369: struct DMCompositeLink *next;
1370: DM_Composite *com = (DM_Composite *)dmi->data;
1371: DM dm;
1373: PetscFunctionBegin;
1375: if (comm == MPI_COMM_NULL) PetscCall(PetscObjectGetComm((PetscObject)dmi, &comm));
1376: PetscCall(DMSetUp(dmi));
1377: next = com->next;
1378: PetscCall(DMCompositeCreate(comm, fine));
1380: /* loop over packed objects, handling one at at time */
1381: while (next) {
1382: PetscCall(DMRefine(next->dm, comm, &dm));
1383: PetscCall(DMCompositeAddDM(*fine, dm));
1384: PetscCall(PetscObjectDereference((PetscObject)dm));
1385: next = next->next;
1386: }
1387: PetscFunctionReturn(PETSC_SUCCESS);
1388: }
1390: PetscErrorCode DMCoarsen_Composite(DM dmi, MPI_Comm comm, DM *fine)
1391: {
1392: struct DMCompositeLink *next;
1393: DM_Composite *com = (DM_Composite *)dmi->data;
1394: DM dm;
1396: PetscFunctionBegin;
1398: PetscCall(DMSetUp(dmi));
1399: if (comm == MPI_COMM_NULL) PetscCall(PetscObjectGetComm((PetscObject)dmi, &comm));
1400: next = com->next;
1401: PetscCall(DMCompositeCreate(comm, fine));
1403: /* loop over packed objects, handling one at at time */
1404: while (next) {
1405: PetscCall(DMCoarsen(next->dm, comm, &dm));
1406: PetscCall(DMCompositeAddDM(*fine, dm));
1407: PetscCall(PetscObjectDereference((PetscObject)dm));
1408: next = next->next;
1409: }
1410: PetscFunctionReturn(PETSC_SUCCESS);
1411: }
1413: PetscErrorCode DMCreateInterpolation_Composite(DM coarse, DM fine, Mat *A, Vec *v)
1414: {
1415: PetscInt m, n, M, N, nDM, i;
1416: struct DMCompositeLink *nextc;
1417: struct DMCompositeLink *nextf;
1418: Vec gcoarse, gfine, *vecs;
1419: DM_Composite *comcoarse = (DM_Composite *)coarse->data;
1420: DM_Composite *comfine = (DM_Composite *)fine->data;
1421: Mat *mats;
1423: PetscFunctionBegin;
1426: PetscCall(DMSetUp(coarse));
1427: PetscCall(DMSetUp(fine));
1428: /* use global vectors only for determining matrix layout */
1429: PetscCall(DMGetGlobalVector(coarse, &gcoarse));
1430: PetscCall(DMGetGlobalVector(fine, &gfine));
1431: PetscCall(VecGetLocalSize(gcoarse, &n));
1432: PetscCall(VecGetLocalSize(gfine, &m));
1433: PetscCall(VecGetSize(gcoarse, &N));
1434: PetscCall(VecGetSize(gfine, &M));
1435: PetscCall(DMRestoreGlobalVector(coarse, &gcoarse));
1436: PetscCall(DMRestoreGlobalVector(fine, &gfine));
1438: nDM = comfine->nDM;
1439: PetscCheck(nDM == comcoarse->nDM, PetscObjectComm((PetscObject)fine), PETSC_ERR_ARG_INCOMP, "Fine DMComposite has %" PetscInt_FMT " entries, but coarse has %" PetscInt_FMT, nDM, comcoarse->nDM);
1440: PetscCall(PetscCalloc1(nDM * nDM, &mats));
1441: if (v) PetscCall(PetscCalloc1(nDM, &vecs));
1443: /* loop over packed objects, handling one at at time */
1444: for (nextc = comcoarse->next, nextf = comfine->next, i = 0; nextc; nextc = nextc->next, nextf = nextf->next, i++) {
1445: if (!v) PetscCall(DMCreateInterpolation(nextc->dm, nextf->dm, &mats[i * nDM + i], NULL));
1446: else PetscCall(DMCreateInterpolation(nextc->dm, nextf->dm, &mats[i * nDM + i], &vecs[i]));
1447: }
1448: PetscCall(MatCreateNest(PetscObjectComm((PetscObject)fine), nDM, NULL, nDM, NULL, mats, A));
1449: if (v) PetscCall(VecCreateNest(PetscObjectComm((PetscObject)fine), nDM, NULL, vecs, v));
1450: for (i = 0; i < nDM * nDM; i++) PetscCall(MatDestroy(&mats[i]));
1451: PetscCall(PetscFree(mats));
1452: if (v) {
1453: for (i = 0; i < nDM; i++) PetscCall(VecDestroy(&vecs[i]));
1454: PetscCall(PetscFree(vecs));
1455: }
1456: PetscFunctionReturn(PETSC_SUCCESS);
1457: }
1459: static PetscErrorCode DMGetLocalToGlobalMapping_Composite(DM dm)
1460: {
1461: DM_Composite *com = (DM_Composite *)dm->data;
1462: ISLocalToGlobalMapping *ltogs;
1463: PetscInt i;
1465: PetscFunctionBegin;
1466: /* Set the ISLocalToGlobalMapping on the new matrix */
1467: PetscCall(DMCompositeGetISLocalToGlobalMappings(dm, <ogs));
1468: PetscCall(ISLocalToGlobalMappingConcatenate(PetscObjectComm((PetscObject)dm), com->nDM, ltogs, &dm->ltogmap));
1469: for (i = 0; i < com->nDM; i++) PetscCall(ISLocalToGlobalMappingDestroy(<ogs[i]));
1470: PetscCall(PetscFree(ltogs));
1471: PetscFunctionReturn(PETSC_SUCCESS);
1472: }
1474: PetscErrorCode DMCreateColoring_Composite(DM dm, ISColoringType ctype, ISColoring *coloring)
1475: {
1476: PetscInt n, i, cnt;
1477: ISColoringValue *colors;
1478: PetscBool dense = PETSC_FALSE;
1479: ISColoringValue maxcol = 0;
1480: DM_Composite *com = (DM_Composite *)dm->data;
1482: PetscFunctionBegin;
1484: PetscCheck(ctype != IS_COLORING_LOCAL, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Only global coloring supported");
1485: if (ctype == IS_COLORING_GLOBAL) {
1486: n = com->n;
1487: } else SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "Unknown ISColoringType");
1488: PetscCall(PetscMalloc1(n, &colors)); /* freed in ISColoringDestroy() */
1490: PetscCall(PetscOptionsGetBool(((PetscObject)dm)->options, ((PetscObject)dm)->prefix, "-dmcomposite_dense_jacobian", &dense, NULL));
1491: if (dense) {
1492: for (i = 0; i < n; i++) colors[i] = (ISColoringValue)(com->rstart + i);
1493: maxcol = com->N;
1494: } else {
1495: struct DMCompositeLink *next = com->next;
1496: PetscMPIInt rank;
1498: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
1499: cnt = 0;
1500: while (next) {
1501: ISColoring lcoloring;
1503: PetscCall(DMCreateColoring(next->dm, IS_COLORING_GLOBAL, &lcoloring));
1504: for (i = 0; i < lcoloring->N; i++) colors[cnt++] = maxcol + lcoloring->colors[i];
1505: maxcol += lcoloring->n;
1506: PetscCall(ISColoringDestroy(&lcoloring));
1507: next = next->next;
1508: }
1509: }
1510: PetscCall(ISColoringCreate(PetscObjectComm((PetscObject)dm), maxcol, n, colors, PETSC_OWN_POINTER, coloring));
1511: PetscFunctionReturn(PETSC_SUCCESS);
1512: }
1514: PetscErrorCode DMGlobalToLocalBegin_Composite(DM dm, Vec gvec, InsertMode mode, Vec lvec)
1515: {
1516: struct DMCompositeLink *next;
1517: PetscScalar *garray, *larray;
1518: DM_Composite *com = (DM_Composite *)dm->data;
1520: PetscFunctionBegin;
1524: if (!com->setup) PetscCall(DMSetUp(dm));
1526: PetscCall(VecGetArray(gvec, &garray));
1527: PetscCall(VecGetArray(lvec, &larray));
1529: /* loop over packed objects, handling one at at time */
1530: next = com->next;
1531: while (next) {
1532: Vec local, global;
1533: PetscInt N;
1535: PetscCall(DMGetGlobalVector(next->dm, &global));
1536: PetscCall(VecGetLocalSize(global, &N));
1537: PetscCall(VecPlaceArray(global, garray));
1538: PetscCall(DMGetLocalVector(next->dm, &local));
1539: PetscCall(VecPlaceArray(local, larray));
1540: PetscCall(DMGlobalToLocalBegin(next->dm, global, mode, local));
1541: PetscCall(DMGlobalToLocalEnd(next->dm, global, mode, local));
1542: PetscCall(VecResetArray(global));
1543: PetscCall(VecResetArray(local));
1544: PetscCall(DMRestoreGlobalVector(next->dm, &global));
1545: PetscCall(DMRestoreLocalVector(next->dm, &local));
1547: larray += next->nlocal;
1548: garray += next->n;
1549: next = next->next;
1550: }
1552: PetscCall(VecRestoreArray(gvec, NULL));
1553: PetscCall(VecRestoreArray(lvec, NULL));
1554: PetscFunctionReturn(PETSC_SUCCESS);
1555: }
1557: PetscErrorCode DMGlobalToLocalEnd_Composite(DM dm, Vec gvec, InsertMode mode, Vec lvec)
1558: {
1559: PetscFunctionBegin;
1563: PetscFunctionReturn(PETSC_SUCCESS);
1564: }
1566: PetscErrorCode DMLocalToGlobalBegin_Composite(DM dm, Vec lvec, InsertMode mode, Vec gvec)
1567: {
1568: struct DMCompositeLink *next;
1569: PetscScalar *larray, *garray;
1570: DM_Composite *com = (DM_Composite *)dm->data;
1572: PetscFunctionBegin;
1577: if (!com->setup) PetscCall(DMSetUp(dm));
1579: PetscCall(VecGetArray(lvec, &larray));
1580: PetscCall(VecGetArray(gvec, &garray));
1582: /* loop over packed objects, handling one at at time */
1583: next = com->next;
1584: while (next) {
1585: Vec global, local;
1587: PetscCall(DMGetLocalVector(next->dm, &local));
1588: PetscCall(VecPlaceArray(local, larray));
1589: PetscCall(DMGetGlobalVector(next->dm, &global));
1590: PetscCall(VecPlaceArray(global, garray));
1591: PetscCall(DMLocalToGlobalBegin(next->dm, local, mode, global));
1592: PetscCall(DMLocalToGlobalEnd(next->dm, local, mode, global));
1593: PetscCall(VecResetArray(local));
1594: PetscCall(VecResetArray(global));
1595: PetscCall(DMRestoreGlobalVector(next->dm, &global));
1596: PetscCall(DMRestoreLocalVector(next->dm, &local));
1598: garray += next->n;
1599: larray += next->nlocal;
1600: next = next->next;
1601: }
1603: PetscCall(VecRestoreArray(gvec, NULL));
1604: PetscCall(VecRestoreArray(lvec, NULL));
1606: PetscFunctionReturn(PETSC_SUCCESS);
1607: }
1609: PetscErrorCode DMLocalToGlobalEnd_Composite(DM dm, Vec lvec, InsertMode mode, Vec gvec)
1610: {
1611: PetscFunctionBegin;
1615: PetscFunctionReturn(PETSC_SUCCESS);
1616: }
1618: PetscErrorCode DMLocalToLocalBegin_Composite(DM dm, Vec vec1, InsertMode mode, Vec vec2)
1619: {
1620: struct DMCompositeLink *next;
1621: PetscScalar *array1, *array2;
1622: DM_Composite *com = (DM_Composite *)dm->data;
1624: PetscFunctionBegin;
1629: if (!com->setup) PetscCall(DMSetUp(dm));
1631: PetscCall(VecGetArray(vec1, &array1));
1632: PetscCall(VecGetArray(vec2, &array2));
1634: /* loop over packed objects, handling one at at time */
1635: next = com->next;
1636: while (next) {
1637: Vec local1, local2;
1639: PetscCall(DMGetLocalVector(next->dm, &local1));
1640: PetscCall(VecPlaceArray(local1, array1));
1641: PetscCall(DMGetLocalVector(next->dm, &local2));
1642: PetscCall(VecPlaceArray(local2, array2));
1643: PetscCall(DMLocalToLocalBegin(next->dm, local1, mode, local2));
1644: PetscCall(DMLocalToLocalEnd(next->dm, local1, mode, local2));
1645: PetscCall(VecResetArray(local2));
1646: PetscCall(DMRestoreLocalVector(next->dm, &local2));
1647: PetscCall(VecResetArray(local1));
1648: PetscCall(DMRestoreLocalVector(next->dm, &local1));
1650: array1 += next->nlocal;
1651: array2 += next->nlocal;
1652: next = next->next;
1653: }
1655: PetscCall(VecRestoreArray(vec1, NULL));
1656: PetscCall(VecRestoreArray(vec2, NULL));
1658: PetscFunctionReturn(PETSC_SUCCESS);
1659: }
1661: PetscErrorCode DMLocalToLocalEnd_Composite(DM dm, Vec lvec, InsertMode mode, Vec gvec)
1662: {
1663: PetscFunctionBegin;
1667: PetscFunctionReturn(PETSC_SUCCESS);
1668: }
1670: /*MC
1671: DMCOMPOSITE = "composite" - A `DM` object that is used to manage data for a collection of `DM`
1673: Level: intermediate
1675: .seealso: `DMType`, `DM`, `DMDACreate()`, `DMCreate()`, `DMSetType()`, `DMCompositeCreate()`
1676: M*/
1678: PETSC_EXTERN PetscErrorCode DMCreate_Composite(DM p)
1679: {
1680: DM_Composite *com;
1682: PetscFunctionBegin;
1683: PetscCall(PetscNew(&com));
1684: p->data = com;
1685: com->n = 0;
1686: com->nghost = 0;
1687: com->next = NULL;
1688: com->nDM = 0;
1690: p->ops->createglobalvector = DMCreateGlobalVector_Composite;
1691: p->ops->createlocalvector = DMCreateLocalVector_Composite;
1692: p->ops->getlocaltoglobalmapping = DMGetLocalToGlobalMapping_Composite;
1693: p->ops->createfieldis = DMCreateFieldIS_Composite;
1694: p->ops->createfielddecomposition = DMCreateFieldDecomposition_Composite;
1695: p->ops->refine = DMRefine_Composite;
1696: p->ops->coarsen = DMCoarsen_Composite;
1697: p->ops->createinterpolation = DMCreateInterpolation_Composite;
1698: p->ops->creatematrix = DMCreateMatrix_Composite;
1699: p->ops->getcoloring = DMCreateColoring_Composite;
1700: p->ops->globaltolocalbegin = DMGlobalToLocalBegin_Composite;
1701: p->ops->globaltolocalend = DMGlobalToLocalEnd_Composite;
1702: p->ops->localtoglobalbegin = DMLocalToGlobalBegin_Composite;
1703: p->ops->localtoglobalend = DMLocalToGlobalEnd_Composite;
1704: p->ops->localtolocalbegin = DMLocalToLocalBegin_Composite;
1705: p->ops->localtolocalend = DMLocalToLocalEnd_Composite;
1706: p->ops->destroy = DMDestroy_Composite;
1707: p->ops->view = DMView_Composite;
1708: p->ops->setup = DMSetUp_Composite;
1710: PetscCall(PetscObjectComposeFunction((PetscObject)p, "DMSetUpGLVisViewer_C", DMSetUpGLVisViewer_Composite));
1711: PetscFunctionReturn(PETSC_SUCCESS);
1712: }
1714: /*@
1715: DMCompositeCreate - Creates a `DMCOMPOSITE`, used to generate "composite"
1716: vectors made up of several subvectors.
1718: Collective
1720: Input Parameter:
1721: . comm - the processors that will share the global vector
1723: Output Parameter:
1724: . packer - the `DMCOMPOSITE` object
1726: Level: advanced
1728: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCompositeScatter()`, `DMCOMPOSITE`, `DMCreate()`
1729: `DMCompositeGather()`, `DMCreateGlobalVector()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`
1730: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
1731: @*/
1732: PetscErrorCode DMCompositeCreate(MPI_Comm comm, DM *packer)
1733: {
1734: PetscFunctionBegin;
1736: PetscCall(DMCreate(comm, packer));
1737: PetscCall(DMSetType(*packer, DMCOMPOSITE));
1738: PetscFunctionReturn(PETSC_SUCCESS);
1739: }