Actual source code: pack.c
1: #include <../src/dm/impls/composite/packimpl.h>
2: #include <petsc/private/isimpl.h>
3: #include <petsc/private/glvisviewerimpl.h>
4: #include <petscds.h>
6: /*@C
7: DMCompositeSetCoupling - Sets user provided routines that compute the coupling between the
8: separate components `DM` in a `DMCOMPOSITE` to build the correct matrix nonzero structure.
10: Logically Collective; No Fortran Support
12: Input Parameters:
13: + dm - the composite object
14: - FormCoupleLocations - routine to set the nonzero locations in the matrix
16: Level: advanced
18: Note:
19: See `DMSetApplicationContext()` and `DMGetApplicationContext()` for how to get user information into
20: this routine
22: .seealso: `DMCOMPOSITE`, `DM`
23: @*/
24: PetscErrorCode DMCompositeSetCoupling(DM dm, PetscErrorCode (*FormCoupleLocations)(DM, Mat, PetscInt *, PetscInt *, PetscInt, PetscInt, PetscInt, PetscInt))
25: {
26: DM_Composite *com = (DM_Composite *)dm->data;
27: PetscBool flg;
29: PetscFunctionBegin;
30: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
31: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
32: com->FormCoupleLocations = FormCoupleLocations;
33: PetscFunctionReturn(PETSC_SUCCESS);
34: }
36: static PetscErrorCode DMDestroy_Composite(DM dm)
37: {
38: struct DMCompositeLink *next, *prev;
39: DM_Composite *com = (DM_Composite *)dm->data;
41: PetscFunctionBegin;
42: next = com->next;
43: while (next) {
44: prev = next;
45: next = next->next;
46: PetscCall(DMDestroy(&prev->dm));
47: PetscCall(PetscFree(prev->grstarts));
48: PetscCall(PetscFree(prev));
49: }
50: PetscCall(PetscObjectComposeFunction((PetscObject)dm, "DMSetUpGLVisViewer_C", NULL));
51: /* This was originally freed in DMDestroy(), but that prevents reference counting of backend objects */
52: PetscCall(PetscFree(com));
53: PetscFunctionReturn(PETSC_SUCCESS);
54: }
56: static PetscErrorCode DMView_Composite(DM dm, PetscViewer v)
57: {
58: PetscBool iascii;
59: DM_Composite *com = (DM_Composite *)dm->data;
61: PetscFunctionBegin;
62: PetscCall(PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &iascii));
63: if (iascii) {
64: struct DMCompositeLink *lnk = com->next;
65: PetscInt i;
67: PetscCall(PetscViewerASCIIPrintf(v, "DM (%s)\n", ((PetscObject)dm)->prefix ? ((PetscObject)dm)->prefix : "no prefix"));
68: PetscCall(PetscViewerASCIIPrintf(v, " contains %" PetscInt_FMT " DMs\n", com->nDM));
69: PetscCall(PetscViewerASCIIPushTab(v));
70: for (i = 0; lnk; lnk = lnk->next, i++) {
71: PetscCall(PetscViewerASCIIPrintf(v, "Link %" PetscInt_FMT ": DM of type %s\n", i, ((PetscObject)lnk->dm)->type_name));
72: PetscCall(PetscViewerASCIIPushTab(v));
73: PetscCall(DMView(lnk->dm, v));
74: PetscCall(PetscViewerASCIIPopTab(v));
75: }
76: PetscCall(PetscViewerASCIIPopTab(v));
77: }
78: PetscFunctionReturn(PETSC_SUCCESS);
79: }
81: /* --------------------------------------------------------------------------------------*/
82: static PetscErrorCode DMSetUp_Composite(DM dm)
83: {
84: PetscInt nprev = 0;
85: PetscMPIInt rank, size;
86: DM_Composite *com = (DM_Composite *)dm->data;
87: struct DMCompositeLink *next = com->next;
88: PetscLayout map;
90: PetscFunctionBegin;
91: PetscCheck(!com->setup, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Packer has already been setup");
92: PetscCall(PetscLayoutCreate(PetscObjectComm((PetscObject)dm), &map));
93: PetscCall(PetscLayoutSetLocalSize(map, com->n));
94: PetscCall(PetscLayoutSetSize(map, PETSC_DETERMINE));
95: PetscCall(PetscLayoutSetBlockSize(map, 1));
96: PetscCall(PetscLayoutSetUp(map));
97: PetscCall(PetscLayoutGetSize(map, &com->N));
98: PetscCall(PetscLayoutGetRange(map, &com->rstart, NULL));
99: PetscCall(PetscLayoutDestroy(&map));
101: /* now set the rstart for each linked vector */
102: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
103: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)dm), &size));
104: while (next) {
105: next->rstart = nprev;
106: nprev += next->n;
107: next->grstart = com->rstart + next->rstart;
108: PetscCall(PetscMalloc1(size, &next->grstarts));
109: PetscCallMPI(MPI_Allgather(&next->grstart, 1, MPIU_INT, next->grstarts, 1, MPIU_INT, PetscObjectComm((PetscObject)dm)));
110: next = next->next;
111: }
112: com->setup = PETSC_TRUE;
113: PetscFunctionReturn(PETSC_SUCCESS);
114: }
116: /* ----------------------------------------------------------------------------------*/
118: /*@
119: DMCompositeGetNumberDM - Gets the number of `DM` objects in the `DMCOMPOSITE`
120: representation.
122: Not Collective
124: Input Parameter:
125: . dm - the `DMCOMPOSITE` object
127: Output Parameter:
128: . nDM - the number of `DM`
130: Level: beginner
132: .seealso: `DMCOMPOSITE`, `DM`
133: @*/
134: PetscErrorCode DMCompositeGetNumberDM(DM dm, PetscInt *nDM)
135: {
136: DM_Composite *com = (DM_Composite *)dm->data;
137: PetscBool flg;
139: PetscFunctionBegin;
141: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
142: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
143: *nDM = com->nDM;
144: PetscFunctionReturn(PETSC_SUCCESS);
145: }
147: /*@C
148: DMCompositeGetAccess - Allows one to access the individual packed vectors in their global
149: representation.
151: Collective
153: Input Parameters:
154: + dm - the `DMCOMPOSITE` object
155: - gvec - the global vector
157: Output Parameter:
158: . ... - the packed parallel vectors, `NULL` for those that are not needed
160: Level: advanced
162: Note:
163: Use `DMCompositeRestoreAccess()` to return the vectors when you no longer need them
165: Fortran Notes:
166: Fortran callers must use numbered versions of this routine, e.g., DMCompositeGetAccess4(dm,gvec,vec1,vec2,vec3,vec4)
167: or use the alternative interface `DMCompositeGetAccessArray()`.
169: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeGetEntries()`, `DMCompositeScatter()`
170: @*/
171: PetscErrorCode DMCompositeGetAccess(DM dm, Vec gvec, ...)
172: {
173: va_list Argp;
174: struct DMCompositeLink *next;
175: DM_Composite *com = (DM_Composite *)dm->data;
176: PetscInt readonly;
177: PetscBool flg;
179: PetscFunctionBegin;
182: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
183: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
184: next = com->next;
185: if (!com->setup) PetscCall(DMSetUp(dm));
187: PetscCall(VecLockGet(gvec, &readonly));
188: /* loop over packed objects, handling one at at time */
189: va_start(Argp, gvec);
190: while (next) {
191: Vec *vec;
192: vec = va_arg(Argp, Vec *);
193: if (vec) {
194: PetscCall(DMGetGlobalVector(next->dm, vec));
195: if (readonly) {
196: const PetscScalar *array;
197: PetscCall(VecGetArrayRead(gvec, &array));
198: PetscCall(VecPlaceArray(*vec, array + next->rstart));
199: PetscCall(VecLockReadPush(*vec));
200: PetscCall(VecRestoreArrayRead(gvec, &array));
201: } else {
202: PetscScalar *array;
203: PetscCall(VecGetArray(gvec, &array));
204: PetscCall(VecPlaceArray(*vec, array + next->rstart));
205: PetscCall(VecRestoreArray(gvec, &array));
206: }
207: }
208: next = next->next;
209: }
210: va_end(Argp);
211: PetscFunctionReturn(PETSC_SUCCESS);
212: }
214: /*@C
215: DMCompositeGetAccessArray - Allows one to access the individual packed vectors in their global
216: representation.
218: Collective
220: Input Parameters:
221: + dm - the `DMCOMPOSITE`
222: . pvec - packed vector
223: . nwanted - number of vectors wanted
224: - wanted - sorted array of vectors wanted, or NULL to get all vectors
226: Output Parameter:
227: . vecs - array of requested global vectors (must be allocated)
229: Level: advanced
231: Note:
232: Use `DMCompositeRestoreAccessArray()` to return the vectors when you no longer need them
234: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeGetAccess()`, `DMCompositeGetEntries()`, `DMCompositeScatter()`, `DMCompositeGather()`
235: @*/
236: PetscErrorCode DMCompositeGetAccessArray(DM dm, Vec pvec, PetscInt nwanted, const PetscInt *wanted, Vec *vecs)
237: {
238: struct DMCompositeLink *link;
239: PetscInt i, wnum;
240: DM_Composite *com = (DM_Composite *)dm->data;
241: PetscInt readonly;
242: PetscBool flg;
244: PetscFunctionBegin;
247: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
248: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
249: if (!com->setup) PetscCall(DMSetUp(dm));
251: PetscCall(VecLockGet(pvec, &readonly));
252: for (i = 0, wnum = 0, link = com->next; link && wnum < nwanted; i++, link = link->next) {
253: if (!wanted || i == wanted[wnum]) {
254: Vec v;
255: PetscCall(DMGetGlobalVector(link->dm, &v));
256: if (readonly) {
257: const PetscScalar *array;
258: PetscCall(VecGetArrayRead(pvec, &array));
259: PetscCall(VecPlaceArray(v, array + link->rstart));
260: PetscCall(VecLockReadPush(v));
261: PetscCall(VecRestoreArrayRead(pvec, &array));
262: } else {
263: PetscScalar *array;
264: PetscCall(VecGetArray(pvec, &array));
265: PetscCall(VecPlaceArray(v, array + link->rstart));
266: PetscCall(VecRestoreArray(pvec, &array));
267: }
268: vecs[wnum++] = v;
269: }
270: }
271: PetscFunctionReturn(PETSC_SUCCESS);
272: }
274: /*@C
275: DMCompositeGetLocalAccessArray - Allows one to access the individual
276: packed vectors in their local representation.
278: Collective
280: Input Parameters:
281: + dm - the `DMCOMPOSITE`
282: . pvec - packed vector
283: . nwanted - number of vectors wanted
284: - wanted - sorted array of vectors wanted, or NULL to get all vectors
286: Output Parameter:
287: . vecs - array of requested local vectors (must be allocated)
289: Level: advanced
291: Note:
292: Use `DMCompositeRestoreLocalAccessArray()` to return the vectors
293: when you no longer need them.
295: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeRestoreLocalAccessArray()`, `DMCompositeGetAccess()`,
296: `DMCompositeGetEntries()`, `DMCompositeScatter()`, `DMCompositeGather()`
297: @*/
298: PetscErrorCode DMCompositeGetLocalAccessArray(DM dm, Vec pvec, PetscInt nwanted, const PetscInt *wanted, Vec *vecs)
299: {
300: struct DMCompositeLink *link;
301: PetscInt i, wnum;
302: DM_Composite *com = (DM_Composite *)dm->data;
303: PetscInt readonly;
304: PetscInt nlocal = 0;
305: PetscBool flg;
307: PetscFunctionBegin;
310: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
311: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
312: if (!com->setup) PetscCall(DMSetUp(dm));
314: PetscCall(VecLockGet(pvec, &readonly));
315: for (i = 0, wnum = 0, link = com->next; link && wnum < nwanted; i++, link = link->next) {
316: if (!wanted || i == wanted[wnum]) {
317: Vec v;
318: PetscCall(DMGetLocalVector(link->dm, &v));
319: if (readonly) {
320: const PetscScalar *array;
321: PetscCall(VecGetArrayRead(pvec, &array));
322: PetscCall(VecPlaceArray(v, array + nlocal));
323: // 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
324: PetscCall(VecLockReadPush(v));
325: PetscCall(VecRestoreArrayRead(pvec, &array));
326: } else {
327: PetscScalar *array;
328: PetscCall(VecGetArray(pvec, &array));
329: PetscCall(VecPlaceArray(v, array + nlocal));
330: PetscCall(VecRestoreArray(pvec, &array));
331: }
332: vecs[wnum++] = v;
333: }
335: nlocal += link->nlocal;
336: }
338: PetscFunctionReturn(PETSC_SUCCESS);
339: }
341: /*@C
342: DMCompositeRestoreAccess - Returns the vectors obtained with `DMCompositeGetAccess()`
343: representation.
345: Collective
347: Input Parameters:
348: + dm - the `DMCOMPOSITE` object
349: . gvec - the global vector
350: - ... - the individual parallel vectors, `NULL` for those that are not needed
352: Level: advanced
354: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
355: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeScatter()`,
356: `DMCompositeGetAccess()`
357: @*/
358: PetscErrorCode DMCompositeRestoreAccess(DM dm, Vec gvec, ...)
359: {
360: va_list Argp;
361: struct DMCompositeLink *next;
362: DM_Composite *com = (DM_Composite *)dm->data;
363: PetscInt readonly;
364: PetscBool flg;
366: PetscFunctionBegin;
369: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
370: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
371: next = com->next;
372: if (!com->setup) PetscCall(DMSetUp(dm));
374: PetscCall(VecLockGet(gvec, &readonly));
375: /* loop over packed objects, handling one at at time */
376: va_start(Argp, gvec);
377: while (next) {
378: Vec *vec;
379: vec = va_arg(Argp, Vec *);
380: if (vec) {
381: PetscCall(VecResetArray(*vec));
382: if (readonly) PetscCall(VecLockReadPop(*vec));
383: PetscCall(DMRestoreGlobalVector(next->dm, vec));
384: }
385: next = next->next;
386: }
387: va_end(Argp);
388: PetscFunctionReturn(PETSC_SUCCESS);
389: }
391: /*@C
392: DMCompositeRestoreAccessArray - Returns the vectors obtained with `DMCompositeGetAccessArray()`
394: Collective
396: Input Parameters:
397: + dm - the `DMCOMPOSITE` object
398: . pvec - packed vector
399: . nwanted - number of vectors wanted
400: . wanted - sorted array of vectors wanted, or NULL to get all vectors
401: - vecs - array of global vectors to return
403: Level: advanced
405: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeRestoreAccess()`, `DMCompositeRestoreEntries()`, `DMCompositeScatter()`, `DMCompositeGather()`
406: @*/
407: PetscErrorCode DMCompositeRestoreAccessArray(DM dm, Vec pvec, PetscInt nwanted, const PetscInt *wanted, Vec *vecs)
408: {
409: struct DMCompositeLink *link;
410: PetscInt i, wnum;
411: DM_Composite *com = (DM_Composite *)dm->data;
412: PetscInt readonly;
413: PetscBool flg;
415: PetscFunctionBegin;
418: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
419: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
420: if (!com->setup) PetscCall(DMSetUp(dm));
422: PetscCall(VecLockGet(pvec, &readonly));
423: for (i = 0, wnum = 0, link = com->next; link && wnum < nwanted; i++, link = link->next) {
424: if (!wanted || i == wanted[wnum]) {
425: PetscCall(VecResetArray(vecs[wnum]));
426: if (readonly) PetscCall(VecLockReadPop(vecs[wnum]));
427: PetscCall(DMRestoreGlobalVector(link->dm, &vecs[wnum]));
428: wnum++;
429: }
430: }
431: PetscFunctionReturn(PETSC_SUCCESS);
432: }
434: /*@C
435: DMCompositeRestoreLocalAccessArray - Returns the vectors obtained with `DMCompositeGetLocalAccessArray()`.
437: Collective
439: Input Parameters:
440: + dm - the `DMCOMPOSITE` object
441: . pvec - packed vector
442: . nwanted - number of vectors wanted
443: . wanted - sorted array of vectors wanted, or NULL to restore all vectors
444: - vecs - array of local vectors to return
446: Level: advanced
448: Note:
449: nwanted and wanted must match the values given to `DMCompositeGetLocalAccessArray()`
450: otherwise the call will fail.
452: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeGetLocalAccessArray()`, `DMCompositeRestoreAccessArray()`,
453: `DMCompositeRestoreAccess()`, `DMCompositeRestoreEntries()`,
454: `DMCompositeScatter()`, `DMCompositeGather()`
455: @*/
456: PetscErrorCode DMCompositeRestoreLocalAccessArray(DM dm, Vec pvec, PetscInt nwanted, const PetscInt *wanted, Vec *vecs)
457: {
458: struct DMCompositeLink *link;
459: PetscInt i, wnum;
460: DM_Composite *com = (DM_Composite *)dm->data;
461: PetscInt readonly;
462: PetscBool flg;
464: PetscFunctionBegin;
467: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
468: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
469: if (!com->setup) PetscCall(DMSetUp(dm));
471: PetscCall(VecLockGet(pvec, &readonly));
472: for (i = 0, wnum = 0, link = com->next; link && wnum < nwanted; i++, link = link->next) {
473: if (!wanted || i == wanted[wnum]) {
474: PetscCall(VecResetArray(vecs[wnum]));
475: if (readonly) PetscCall(VecLockReadPop(vecs[wnum]));
476: PetscCall(DMRestoreLocalVector(link->dm, &vecs[wnum]));
477: wnum++;
478: }
479: }
480: PetscFunctionReturn(PETSC_SUCCESS);
481: }
483: /*@C
484: DMCompositeScatter - Scatters from a global packed vector into its individual local vectors
486: Collective
488: Input Parameters:
489: + dm - the `DMCOMPOSITE` object
490: . gvec - the global vector
491: - ... - the individual sequential vectors, `NULL` for those that are not needed
493: Level: advanced
495: Note:
496: `DMCompositeScatterArray()` is a non-variadic alternative that is often more convenient for library callers and is
497: accessible from Fortran.
499: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
500: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
501: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
502: `DMCompositeScatterArray()`
503: @*/
504: PetscErrorCode DMCompositeScatter(DM dm, Vec gvec, ...)
505: {
506: va_list Argp;
507: struct DMCompositeLink *next;
508: PETSC_UNUSED PetscInt cnt;
509: DM_Composite *com = (DM_Composite *)dm->data;
510: PetscBool flg;
512: PetscFunctionBegin;
515: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
516: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
517: if (!com->setup) PetscCall(DMSetUp(dm));
519: /* loop over packed objects, handling one at at time */
520: va_start(Argp, gvec);
521: for (cnt = 3, next = com->next; next; cnt++, next = next->next) {
522: Vec local;
523: local = va_arg(Argp, Vec);
524: if (local) {
525: Vec global;
526: const PetscScalar *array;
528: PetscCall(DMGetGlobalVector(next->dm, &global));
529: PetscCall(VecGetArrayRead(gvec, &array));
530: PetscCall(VecPlaceArray(global, array + next->rstart));
531: PetscCall(DMGlobalToLocalBegin(next->dm, global, INSERT_VALUES, local));
532: PetscCall(DMGlobalToLocalEnd(next->dm, global, INSERT_VALUES, local));
533: PetscCall(VecRestoreArrayRead(gvec, &array));
534: PetscCall(VecResetArray(global));
535: PetscCall(DMRestoreGlobalVector(next->dm, &global));
536: }
537: }
538: va_end(Argp);
539: PetscFunctionReturn(PETSC_SUCCESS);
540: }
542: /*@
543: DMCompositeScatterArray - Scatters from a global packed vector into its individual local vectors
545: Collective
547: Input Parameters:
548: + dm - the `DMCOMPOSITE` object
549: . gvec - the global vector
550: - lvecs - array of local vectors, NULL for any that are not needed
552: Level: advanced
554: Note:
555: This is a non-variadic alternative to `DMCompositeScatter()`
557: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`
558: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
559: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
560: @*/
561: PetscErrorCode DMCompositeScatterArray(DM dm, Vec gvec, Vec *lvecs)
562: {
563: struct DMCompositeLink *next;
564: PetscInt i;
565: DM_Composite *com = (DM_Composite *)dm->data;
566: PetscBool flg;
568: PetscFunctionBegin;
571: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
572: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
573: if (!com->setup) PetscCall(DMSetUp(dm));
575: /* loop over packed objects, handling one at at time */
576: for (i = 0, next = com->next; next; next = next->next, i++) {
577: if (lvecs[i]) {
578: Vec global;
579: const PetscScalar *array;
581: PetscCall(DMGetGlobalVector(next->dm, &global));
582: PetscCall(VecGetArrayRead(gvec, &array));
583: PetscCall(VecPlaceArray(global, (PetscScalar *)array + next->rstart));
584: PetscCall(DMGlobalToLocalBegin(next->dm, global, INSERT_VALUES, lvecs[i]));
585: PetscCall(DMGlobalToLocalEnd(next->dm, global, INSERT_VALUES, lvecs[i]));
586: PetscCall(VecRestoreArrayRead(gvec, &array));
587: PetscCall(VecResetArray(global));
588: PetscCall(DMRestoreGlobalVector(next->dm, &global));
589: }
590: }
591: PetscFunctionReturn(PETSC_SUCCESS);
592: }
594: /*@C
595: DMCompositeGather - Gathers into a global packed vector from its individual local vectors
597: Collective
599: Input Parameters:
600: + dm - the `DMCOMPOSITE` object
601: . imode - `INSERT_VALUES` or `ADD_VALUES`
602: . gvec - the global vector
603: - ... - the individual sequential vectors, `NULL` for any that are not needed
605: Level: advanced
607: Fortran Notes:
608: Fortran users should use `DMCompositeGatherArray()`
610: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
611: `DMCompositeScatter()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
612: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
613: @*/
614: PetscErrorCode DMCompositeGather(DM dm, InsertMode imode, Vec gvec, ...)
615: {
616: va_list Argp;
617: struct DMCompositeLink *next;
618: DM_Composite *com = (DM_Composite *)dm->data;
619: PETSC_UNUSED PetscInt cnt;
620: PetscBool flg;
622: PetscFunctionBegin;
625: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
626: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
627: if (!com->setup) PetscCall(DMSetUp(dm));
629: /* loop over packed objects, handling one at at time */
630: va_start(Argp, gvec);
631: for (cnt = 3, next = com->next; next; cnt++, next = next->next) {
632: Vec local;
633: local = va_arg(Argp, Vec);
634: if (local) {
635: PetscScalar *array;
636: Vec global;
638: PetscCall(DMGetGlobalVector(next->dm, &global));
639: PetscCall(VecGetArray(gvec, &array));
640: PetscCall(VecPlaceArray(global, array + next->rstart));
641: PetscCall(DMLocalToGlobalBegin(next->dm, local, imode, global));
642: PetscCall(DMLocalToGlobalEnd(next->dm, local, imode, global));
643: PetscCall(VecRestoreArray(gvec, &array));
644: PetscCall(VecResetArray(global));
645: PetscCall(DMRestoreGlobalVector(next->dm, &global));
646: }
647: }
648: va_end(Argp);
649: PetscFunctionReturn(PETSC_SUCCESS);
650: }
652: /*@
653: DMCompositeGatherArray - Gathers into a global packed vector from its individual local vectors
655: Collective
657: Input Parameters:
658: + dm - the `DMCOMPOSITE` object
659: . gvec - the global vector
660: . imode - `INSERT_VALUES` or `ADD_VALUES`
661: - lvecs - the individual sequential vectors, NULL for any that are not needed
663: Level: advanced
665: Note:
666: This is a non-variadic alternative to `DMCompositeGather()`.
668: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
669: `DMCompositeScatter()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
670: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`,
671: @*/
672: PetscErrorCode DMCompositeGatherArray(DM dm, InsertMode imode, Vec gvec, Vec *lvecs)
673: {
674: struct DMCompositeLink *next;
675: DM_Composite *com = (DM_Composite *)dm->data;
676: PetscInt i;
677: PetscBool flg;
679: PetscFunctionBegin;
682: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
683: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
684: if (!com->setup) PetscCall(DMSetUp(dm));
686: /* loop over packed objects, handling one at at time */
687: for (next = com->next, i = 0; next; next = next->next, i++) {
688: if (lvecs[i]) {
689: PetscScalar *array;
690: Vec global;
692: PetscCall(DMGetGlobalVector(next->dm, &global));
693: PetscCall(VecGetArray(gvec, &array));
694: PetscCall(VecPlaceArray(global, array + next->rstart));
695: PetscCall(DMLocalToGlobalBegin(next->dm, lvecs[i], imode, global));
696: PetscCall(DMLocalToGlobalEnd(next->dm, lvecs[i], imode, global));
697: PetscCall(VecRestoreArray(gvec, &array));
698: PetscCall(VecResetArray(global));
699: PetscCall(DMRestoreGlobalVector(next->dm, &global));
700: }
701: }
702: PetscFunctionReturn(PETSC_SUCCESS);
703: }
705: /*@
706: DMCompositeAddDM - adds a `DM` vector to a `DMCOMPOSITE`
708: Collective
710: Input Parameters:
711: + dmc - the `DMCOMPOSITE` object
712: - dm - the `DM` object
714: Level: advanced
716: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeGather()`, `DMCreateGlobalVector()`,
717: `DMCompositeScatter()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
718: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
719: @*/
720: PetscErrorCode DMCompositeAddDM(DM dmc, DM dm)
721: {
722: PetscInt n, nlocal;
723: struct DMCompositeLink *mine, *next;
724: Vec global, local;
725: DM_Composite *com = (DM_Composite *)dmc->data;
726: PetscBool flg;
728: PetscFunctionBegin;
731: PetscCall(PetscObjectTypeCompare((PetscObject)dmc, DMCOMPOSITE, &flg));
732: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
733: next = com->next;
734: PetscCheck(!com->setup, PetscObjectComm((PetscObject)dmc), PETSC_ERR_ARG_WRONGSTATE, "Cannot add a DM once you have used the DMComposite");
736: /* create new link */
737: PetscCall(PetscNew(&mine));
738: PetscCall(PetscObjectReference((PetscObject)dm));
739: PetscCall(DMGetGlobalVector(dm, &global));
740: PetscCall(VecGetLocalSize(global, &n));
741: PetscCall(DMRestoreGlobalVector(dm, &global));
742: PetscCall(DMGetLocalVector(dm, &local));
743: PetscCall(VecGetSize(local, &nlocal));
744: PetscCall(DMRestoreLocalVector(dm, &local));
746: mine->n = n;
747: mine->nlocal = nlocal;
748: mine->dm = dm;
749: mine->next = NULL;
750: com->n += n;
751: com->nghost += nlocal;
753: /* add to end of list */
754: if (!next) com->next = mine;
755: else {
756: while (next->next) next = next->next;
757: next->next = mine;
758: }
759: com->nDM++;
760: com->nmine++;
761: PetscFunctionReturn(PETSC_SUCCESS);
762: }
764: #include <petscdraw.h>
765: PETSC_EXTERN PetscErrorCode VecView_MPI(Vec, PetscViewer);
766: static PetscErrorCode VecView_DMComposite(Vec gvec, PetscViewer viewer)
767: {
768: DM dm;
769: struct DMCompositeLink *next;
770: PetscBool isdraw;
771: DM_Composite *com;
773: PetscFunctionBegin;
774: PetscCall(VecGetDM(gvec, &dm));
775: PetscCheck(dm, PetscObjectComm((PetscObject)gvec), PETSC_ERR_ARG_WRONG, "Vector not generated from a DMComposite");
776: com = (DM_Composite *)dm->data;
777: next = com->next;
779: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERDRAW, &isdraw));
780: if (!isdraw) {
781: /* do I really want to call this? */
782: PetscCall(VecView_MPI(gvec, viewer));
783: } else {
784: PetscInt cnt = 0;
786: /* loop over packed objects, handling one at at time */
787: while (next) {
788: Vec vec;
789: const PetscScalar *array;
790: PetscInt bs;
792: /* Should use VecGetSubVector() eventually, but would need to forward the DM for that to work */
793: PetscCall(DMGetGlobalVector(next->dm, &vec));
794: PetscCall(VecGetArrayRead(gvec, &array));
795: PetscCall(VecPlaceArray(vec, (PetscScalar *)array + next->rstart));
796: PetscCall(VecRestoreArrayRead(gvec, &array));
797: PetscCall(VecView(vec, viewer));
798: PetscCall(VecResetArray(vec));
799: PetscCall(VecGetBlockSize(vec, &bs));
800: PetscCall(DMRestoreGlobalVector(next->dm, &vec));
801: PetscCall(PetscViewerDrawBaseAdd(viewer, bs));
802: cnt += bs;
803: next = next->next;
804: }
805: PetscCall(PetscViewerDrawBaseAdd(viewer, -cnt));
806: }
807: PetscFunctionReturn(PETSC_SUCCESS);
808: }
810: static PetscErrorCode DMCreateGlobalVector_Composite(DM dm, Vec *gvec)
811: {
812: DM_Composite *com = (DM_Composite *)dm->data;
814: PetscFunctionBegin;
816: PetscCall(DMSetFromOptions(dm));
817: PetscCall(DMSetUp(dm));
818: PetscCall(VecCreate(PetscObjectComm((PetscObject)dm), gvec));
819: PetscCall(VecSetType(*gvec, dm->vectype));
820: PetscCall(VecSetSizes(*gvec, com->n, com->N));
821: PetscCall(VecSetDM(*gvec, dm));
822: PetscCall(VecSetOperation(*gvec, VECOP_VIEW, (void (*)(void))VecView_DMComposite));
823: PetscFunctionReturn(PETSC_SUCCESS);
824: }
826: static PetscErrorCode DMCreateLocalVector_Composite(DM dm, Vec *lvec)
827: {
828: DM_Composite *com = (DM_Composite *)dm->data;
830: PetscFunctionBegin;
832: if (!com->setup) {
833: PetscCall(DMSetFromOptions(dm));
834: PetscCall(DMSetUp(dm));
835: }
836: PetscCall(VecCreate(PETSC_COMM_SELF, lvec));
837: PetscCall(VecSetType(*lvec, dm->vectype));
838: PetscCall(VecSetSizes(*lvec, com->nghost, PETSC_DECIDE));
839: PetscCall(VecSetDM(*lvec, dm));
840: PetscFunctionReturn(PETSC_SUCCESS);
841: }
843: /*@C
844: DMCompositeGetISLocalToGlobalMappings - gets an `ISLocalToGlobalMapping` for each `DM` in the `DMCOMPOSITE`, maps to the composite global space
846: Collective; No Fortran Support
848: Input Parameter:
849: . dm - the `DMCOMPOSITE` object
851: Output Parameter:
852: . ltogs - the individual mappings for each packed vector. Note that this includes
853: all the ghost points that individual ghosted `DMDA` may have.
855: Level: advanced
857: Note:
858: Each entry of ltogs should be destroyed with `ISLocalToGlobalMappingDestroy()`, the ltogs array should be freed with `PetscFree()`.
860: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
861: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetAccess()`, `DMCompositeScatter()`,
862: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
863: @*/
864: PetscErrorCode DMCompositeGetISLocalToGlobalMappings(DM dm, ISLocalToGlobalMapping **ltogs)
865: {
866: PetscInt i, *idx, n, cnt;
867: struct DMCompositeLink *next;
868: PetscMPIInt rank;
869: DM_Composite *com = (DM_Composite *)dm->data;
870: PetscBool flg;
872: PetscFunctionBegin;
874: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
875: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
876: PetscCall(DMSetUp(dm));
877: PetscCall(PetscMalloc1(com->nDM, ltogs));
878: next = com->next;
879: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
881: /* loop over packed objects, handling one at at time */
882: cnt = 0;
883: while (next) {
884: ISLocalToGlobalMapping ltog;
885: PetscMPIInt size;
886: const PetscInt *suboff, *indices;
887: Vec global;
889: /* Get sub-DM global indices for each local dof */
890: PetscCall(DMGetLocalToGlobalMapping(next->dm, <og));
891: PetscCall(ISLocalToGlobalMappingGetSize(ltog, &n));
892: PetscCall(ISLocalToGlobalMappingGetIndices(ltog, &indices));
893: PetscCall(PetscMalloc1(n, &idx));
895: /* Get the offsets for the sub-DM global vector */
896: PetscCall(DMGetGlobalVector(next->dm, &global));
897: PetscCall(VecGetOwnershipRanges(global, &suboff));
898: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)global), &size));
900: /* Shift the sub-DM definition of the global space to the composite global space */
901: for (i = 0; i < n; i++) {
902: PetscInt subi = indices[i], lo = 0, hi = size, t;
903: /* There's no consensus on what a negative index means,
904: except for skipping when setting the values in vectors and matrices */
905: if (subi < 0) {
906: idx[i] = subi - next->grstarts[rank];
907: continue;
908: }
909: /* Binary search to find which rank owns subi */
910: while (hi - lo > 1) {
911: t = lo + (hi - lo) / 2;
912: if (suboff[t] > subi) hi = t;
913: else lo = t;
914: }
915: idx[i] = subi - suboff[lo] + next->grstarts[lo];
916: }
917: PetscCall(ISLocalToGlobalMappingRestoreIndices(ltog, &indices));
918: PetscCall(ISLocalToGlobalMappingCreate(PetscObjectComm((PetscObject)dm), 1, n, idx, PETSC_OWN_POINTER, &(*ltogs)[cnt]));
919: PetscCall(DMRestoreGlobalVector(next->dm, &global));
920: next = next->next;
921: cnt++;
922: }
923: PetscFunctionReturn(PETSC_SUCCESS);
924: }
926: /*@C
927: DMCompositeGetLocalISs - Gets index sets for each component of a composite local vector
929: Not Collective; No Fortran Support
931: Input Parameter:
932: . dm - the `DMCOMPOSITE`
934: Output Parameter:
935: . is - array of serial index sets for each each component of the `DMCOMPOSITE`
937: Level: intermediate
939: Notes:
940: At present, a composite local vector does not normally exist. This function is used to provide index sets for
941: `MatGetLocalSubMatrix()`. In the future, the scatters for each entry in the `DMCOMPOSITE` may be be merged into a single
942: scatter to a composite local vector. The user should not typically need to know which is being done.
944: To get the composite global indices at all local points (including ghosts), use `DMCompositeGetISLocalToGlobalMappings()`.
946: To get index sets for pieces of the composite global vector, use `DMCompositeGetGlobalISs()`.
948: Each returned `IS` should be destroyed with `ISDestroy()`, the array should be freed with `PetscFree()`.
950: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeGetGlobalISs()`, `DMCompositeGetISLocalToGlobalMappings()`, `MatGetLocalSubMatrix()`, `MatCreateLocalRef()`
951: @*/
952: PetscErrorCode DMCompositeGetLocalISs(DM dm, IS **is)
953: {
954: DM_Composite *com = (DM_Composite *)dm->data;
955: struct DMCompositeLink *link;
956: PetscInt cnt, start;
957: PetscBool flg;
959: PetscFunctionBegin;
961: PetscAssertPointer(is, 2);
962: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
963: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
964: PetscCall(PetscMalloc1(com->nmine, is));
965: for (cnt = 0, start = 0, link = com->next; link; start += link->nlocal, cnt++, link = link->next) {
966: PetscInt bs;
967: PetscCall(ISCreateStride(PETSC_COMM_SELF, link->nlocal, start, 1, &(*is)[cnt]));
968: PetscCall(DMGetBlockSize(link->dm, &bs));
969: PetscCall(ISSetBlockSize((*is)[cnt], bs));
970: }
971: PetscFunctionReturn(PETSC_SUCCESS);
972: }
974: /*@C
975: DMCompositeGetGlobalISs - Gets the index sets for each composed object in a `DMCOMPOSITE`
977: Collective
979: Input Parameter:
980: . dm - the `DMCOMPOSITE` object
982: Output Parameter:
983: . is - the array of index sets
985: Level: advanced
987: Notes:
988: The is entries should be destroyed with `ISDestroy()`, the is array should be freed with `PetscFree()`
990: These could be used to extract a subset of vector entries for a "multi-physics" preconditioner
992: Use `DMCompositeGetLocalISs()` for index sets in the packed local numbering, and
993: `DMCompositeGetISLocalToGlobalMappings()` for to map local sub-`DM` (including ghost) indices to packed global
994: indices.
996: Fortran Notes:
997: The output argument 'is' must be an allocated array of sufficient length, which can be learned using `DMCompositeGetNumberDM()`.
999: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
1000: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetAccess()`, `DMCompositeScatter()`,
1001: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
1002: @*/
1003: PetscErrorCode DMCompositeGetGlobalISs(DM dm, IS *is[])
1004: {
1005: PetscInt cnt = 0;
1006: struct DMCompositeLink *next;
1007: PetscMPIInt rank;
1008: DM_Composite *com = (DM_Composite *)dm->data;
1009: PetscBool flg;
1011: PetscFunctionBegin;
1013: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1014: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1015: PetscCheck(dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Must call DMSetUp() before");
1016: PetscCall(PetscMalloc1(com->nDM, is));
1017: next = com->next;
1018: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
1020: /* loop over packed objects, handling one at at time */
1021: while (next) {
1022: PetscDS prob;
1024: PetscCall(ISCreateStride(PetscObjectComm((PetscObject)dm), next->n, next->grstart, 1, &(*is)[cnt]));
1025: PetscCall(DMGetDS(dm, &prob));
1026: if (prob) {
1027: MatNullSpace space;
1028: Mat pmat;
1029: PetscObject disc;
1030: PetscInt Nf;
1032: PetscCall(PetscDSGetNumFields(prob, &Nf));
1033: if (cnt < Nf) {
1034: PetscCall(PetscDSGetDiscretization(prob, cnt, &disc));
1035: PetscCall(PetscObjectQuery(disc, "nullspace", (PetscObject *)&space));
1036: if (space) PetscCall(PetscObjectCompose((PetscObject)(*is)[cnt], "nullspace", (PetscObject)space));
1037: PetscCall(PetscObjectQuery(disc, "nearnullspace", (PetscObject *)&space));
1038: if (space) PetscCall(PetscObjectCompose((PetscObject)(*is)[cnt], "nearnullspace", (PetscObject)space));
1039: PetscCall(PetscObjectQuery(disc, "pmat", (PetscObject *)&pmat));
1040: if (pmat) PetscCall(PetscObjectCompose((PetscObject)(*is)[cnt], "pmat", (PetscObject)pmat));
1041: }
1042: }
1043: cnt++;
1044: next = next->next;
1045: }
1046: PetscFunctionReturn(PETSC_SUCCESS);
1047: }
1049: static PetscErrorCode DMCreateFieldIS_Composite(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
1050: {
1051: PetscInt nDM;
1052: DM *dms;
1053: PetscInt i;
1055: PetscFunctionBegin;
1056: PetscCall(DMCompositeGetNumberDM(dm, &nDM));
1057: if (numFields) *numFields = nDM;
1058: PetscCall(DMCompositeGetGlobalISs(dm, fields));
1059: if (fieldNames) {
1060: PetscCall(PetscMalloc1(nDM, &dms));
1061: PetscCall(PetscMalloc1(nDM, fieldNames));
1062: PetscCall(DMCompositeGetEntriesArray(dm, dms));
1063: for (i = 0; i < nDM; i++) {
1064: char buf[256];
1065: const char *splitname;
1067: /* Split naming precedence: object name, prefix, number */
1068: splitname = ((PetscObject)dm)->name;
1069: if (!splitname) {
1070: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)dms[i], &splitname));
1071: if (splitname) {
1072: size_t len;
1073: PetscCall(PetscStrncpy(buf, splitname, sizeof(buf)));
1074: buf[sizeof(buf) - 1] = 0;
1075: PetscCall(PetscStrlen(buf, &len));
1076: if (buf[len - 1] == '_') buf[len - 1] = 0; /* Remove trailing underscore if it was used */
1077: splitname = buf;
1078: }
1079: }
1080: if (!splitname) {
1081: PetscCall(PetscSNPrintf(buf, sizeof(buf), "%" PetscInt_FMT, i));
1082: splitname = buf;
1083: }
1084: PetscCall(PetscStrallocpy(splitname, &(*fieldNames)[i]));
1085: }
1086: PetscCall(PetscFree(dms));
1087: }
1088: PetscFunctionReturn(PETSC_SUCCESS);
1089: }
1091: /*
1092: This could take over from DMCreateFieldIS(), as it is more general,
1093: making DMCreateFieldIS() a special case -- calling with dmlist == NULL;
1094: At this point it's probably best to be less intrusive, however.
1095: */
1096: static PetscErrorCode DMCreateFieldDecomposition_Composite(DM dm, PetscInt *len, char ***namelist, IS **islist, DM **dmlist)
1097: {
1098: PetscInt nDM;
1099: PetscInt i;
1101: PetscFunctionBegin;
1102: PetscCall(DMCreateFieldIS_Composite(dm, len, namelist, islist));
1103: if (dmlist) {
1104: PetscCall(DMCompositeGetNumberDM(dm, &nDM));
1105: PetscCall(PetscMalloc1(nDM, dmlist));
1106: PetscCall(DMCompositeGetEntriesArray(dm, *dmlist));
1107: for (i = 0; i < nDM; i++) PetscCall(PetscObjectReference((PetscObject)((*dmlist)[i])));
1108: }
1109: PetscFunctionReturn(PETSC_SUCCESS);
1110: }
1112: /* -------------------------------------------------------------------------------------*/
1113: /*@C
1114: DMCompositeGetLocalVectors - Gets local vectors for each part of a `DMCOMPOSITE`
1115: Use `DMCompositeRestoreLocalVectors()` to return them.
1117: Not Collective; No Fortran Support
1119: Input Parameter:
1120: . dm - the `DMCOMPOSITE` object
1122: Output Parameter:
1123: . ... - the individual sequential `Vec`s
1125: Level: advanced
1127: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
1128: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
1129: `DMCompositeRestoreLocalVectors()`, `DMCompositeScatter()`, `DMCompositeGetEntries()`
1130: @*/
1131: PetscErrorCode DMCompositeGetLocalVectors(DM dm, ...)
1132: {
1133: va_list Argp;
1134: struct DMCompositeLink *next;
1135: DM_Composite *com = (DM_Composite *)dm->data;
1136: PetscBool flg;
1138: PetscFunctionBegin;
1140: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1141: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1142: next = com->next;
1143: /* loop over packed objects, handling one at at time */
1144: va_start(Argp, dm);
1145: while (next) {
1146: Vec *vec;
1147: vec = va_arg(Argp, Vec *);
1148: if (vec) PetscCall(DMGetLocalVector(next->dm, vec));
1149: next = next->next;
1150: }
1151: va_end(Argp);
1152: PetscFunctionReturn(PETSC_SUCCESS);
1153: }
1155: /*@C
1156: DMCompositeRestoreLocalVectors - Restores local vectors for each part of a `DMCOMPOSITE`
1158: Not Collective; No Fortran Support
1160: Input Parameter:
1161: . dm - the `DMCOMPOSITE` object
1163: Output Parameter:
1164: . ... - the individual sequential `Vec`
1166: Level: advanced
1168: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
1169: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
1170: `DMCompositeGetLocalVectors()`, `DMCompositeScatter()`, `DMCompositeGetEntries()`
1171: @*/
1172: PetscErrorCode DMCompositeRestoreLocalVectors(DM dm, ...)
1173: {
1174: va_list Argp;
1175: struct DMCompositeLink *next;
1176: DM_Composite *com = (DM_Composite *)dm->data;
1177: PetscBool flg;
1179: PetscFunctionBegin;
1181: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1182: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1183: next = com->next;
1184: /* loop over packed objects, handling one at at time */
1185: va_start(Argp, dm);
1186: while (next) {
1187: Vec *vec;
1188: vec = va_arg(Argp, Vec *);
1189: if (vec) PetscCall(DMRestoreLocalVector(next->dm, vec));
1190: next = next->next;
1191: }
1192: va_end(Argp);
1193: PetscFunctionReturn(PETSC_SUCCESS);
1194: }
1196: /* -------------------------------------------------------------------------------------*/
1197: /*@C
1198: DMCompositeGetEntries - Gets the `DM` for each entry in a `DMCOMPOSITE`.
1200: Not Collective
1202: Input Parameter:
1203: . dm - the `DMCOMPOSITE` object
1205: Output Parameter:
1206: . ... - the individual entries `DM`
1208: Level: advanced
1210: Fortran Notes:
1211: Available as `DMCompositeGetEntries()` for one output `DM`, DMCompositeGetEntries2() for 2, etc
1213: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`, `DMCompositeGetEntriesArray()`
1214: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
1215: `DMCompositeRestoreLocalVectors()`, `DMCompositeGetLocalVectors()`, `DMCompositeScatter()`
1216: @*/
1217: PetscErrorCode DMCompositeGetEntries(DM dm, ...)
1218: {
1219: va_list Argp;
1220: struct DMCompositeLink *next;
1221: DM_Composite *com = (DM_Composite *)dm->data;
1222: PetscBool flg;
1224: PetscFunctionBegin;
1226: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1227: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1228: next = com->next;
1229: /* loop over packed objects, handling one at at time */
1230: va_start(Argp, dm);
1231: while (next) {
1232: DM *dmn;
1233: dmn = va_arg(Argp, DM *);
1234: if (dmn) *dmn = next->dm;
1235: next = next->next;
1236: }
1237: va_end(Argp);
1238: PetscFunctionReturn(PETSC_SUCCESS);
1239: }
1241: /*@C
1242: DMCompositeGetEntriesArray - Gets the DM for each entry in a `DMCOMPOSITE`
1244: Not Collective
1246: Input Parameter:
1247: . dm - the `DMCOMPOSITE` object
1249: Output Parameter:
1250: . dms - array of sufficient length (see `DMCompositeGetNumberDM()`) to hold the individual `DM`
1252: Level: advanced
1254: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`, `DMCompositeGetEntries()`
1255: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
1256: `DMCompositeRestoreLocalVectors()`, `DMCompositeGetLocalVectors()`, `DMCompositeScatter()`
1257: @*/
1258: PetscErrorCode DMCompositeGetEntriesArray(DM dm, DM dms[])
1259: {
1260: struct DMCompositeLink *next;
1261: DM_Composite *com = (DM_Composite *)dm->data;
1262: PetscInt i;
1263: PetscBool flg;
1265: PetscFunctionBegin;
1267: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1268: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1269: /* loop over packed objects, handling one at at time */
1270: for (next = com->next, i = 0; next; next = next->next, i++) dms[i] = next->dm;
1271: PetscFunctionReturn(PETSC_SUCCESS);
1272: }
1274: typedef struct {
1275: DM dm;
1276: PetscViewer *subv;
1277: Vec *vecs;
1278: } GLVisViewerCtx;
1280: static PetscErrorCode DestroyGLVisViewerCtx_Private(void *vctx)
1281: {
1282: GLVisViewerCtx *ctx = (GLVisViewerCtx *)vctx;
1283: PetscInt i, n;
1285: PetscFunctionBegin;
1286: PetscCall(DMCompositeGetNumberDM(ctx->dm, &n));
1287: for (i = 0; i < n; i++) PetscCall(PetscViewerDestroy(&ctx->subv[i]));
1288: PetscCall(PetscFree2(ctx->subv, ctx->vecs));
1289: PetscCall(DMDestroy(&ctx->dm));
1290: PetscCall(PetscFree(ctx));
1291: PetscFunctionReturn(PETSC_SUCCESS);
1292: }
1294: static PetscErrorCode DMCompositeSampleGLVisFields_Private(PetscObject oX, PetscInt nf, PetscObject oXfield[], void *vctx)
1295: {
1296: Vec X = (Vec)oX;
1297: GLVisViewerCtx *ctx = (GLVisViewerCtx *)vctx;
1298: PetscInt i, n, cumf;
1300: PetscFunctionBegin;
1301: PetscCall(DMCompositeGetNumberDM(ctx->dm, &n));
1302: PetscCall(DMCompositeGetAccessArray(ctx->dm, X, n, NULL, ctx->vecs));
1303: for (i = 0, cumf = 0; i < n; i++) {
1304: PetscErrorCode (*g2l)(PetscObject, PetscInt, PetscObject[], void *);
1305: void *fctx;
1306: PetscInt nfi;
1308: PetscCall(PetscViewerGLVisGetFields_Internal(ctx->subv[i], &nfi, NULL, NULL, &g2l, NULL, &fctx));
1309: if (!nfi) continue;
1310: if (g2l) PetscCall((*g2l)((PetscObject)ctx->vecs[i], nfi, oXfield + cumf, fctx));
1311: else PetscCall(VecCopy(ctx->vecs[i], (Vec)(oXfield[cumf])));
1312: cumf += nfi;
1313: }
1314: PetscCall(DMCompositeRestoreAccessArray(ctx->dm, X, n, NULL, ctx->vecs));
1315: PetscFunctionReturn(PETSC_SUCCESS);
1316: }
1318: static PetscErrorCode DMSetUpGLVisViewer_Composite(PetscObject odm, PetscViewer viewer)
1319: {
1320: DM dm = (DM)odm, *dms;
1321: Vec *Ufds;
1322: GLVisViewerCtx *ctx;
1323: PetscInt i, n, tnf, *sdim;
1324: char **fecs;
1326: PetscFunctionBegin;
1327: PetscCall(PetscNew(&ctx));
1328: PetscCall(PetscObjectReference((PetscObject)dm));
1329: ctx->dm = dm;
1330: PetscCall(DMCompositeGetNumberDM(dm, &n));
1331: PetscCall(PetscMalloc1(n, &dms));
1332: PetscCall(DMCompositeGetEntriesArray(dm, dms));
1333: PetscCall(PetscMalloc2(n, &ctx->subv, n, &ctx->vecs));
1334: for (i = 0, tnf = 0; i < n; i++) {
1335: PetscInt nf;
1337: PetscCall(PetscViewerCreate(PetscObjectComm(odm), &ctx->subv[i]));
1338: PetscCall(PetscViewerSetType(ctx->subv[i], PETSCVIEWERGLVIS));
1339: PetscCall(PetscViewerGLVisSetDM_Internal(ctx->subv[i], (PetscObject)dms[i]));
1340: PetscCall(PetscViewerGLVisGetFields_Internal(ctx->subv[i], &nf, NULL, NULL, NULL, NULL, NULL));
1341: tnf += nf;
1342: }
1343: PetscCall(PetscFree(dms));
1344: PetscCall(PetscMalloc3(tnf, &fecs, tnf, &sdim, tnf, &Ufds));
1345: for (i = 0, tnf = 0; i < n; i++) {
1346: PetscInt *sd, nf, f;
1347: const char **fec;
1348: Vec *Uf;
1350: PetscCall(PetscViewerGLVisGetFields_Internal(ctx->subv[i], &nf, &fec, &sd, NULL, (PetscObject **)&Uf, NULL));
1351: for (f = 0; f < nf; f++) {
1352: PetscCall(PetscStrallocpy(fec[f], &fecs[tnf + f]));
1353: Ufds[tnf + f] = Uf[f];
1354: sdim[tnf + f] = sd[f];
1355: }
1356: tnf += nf;
1357: }
1358: PetscCall(PetscViewerGLVisSetFields(viewer, tnf, (const char **)fecs, sdim, DMCompositeSampleGLVisFields_Private, (PetscObject *)Ufds, ctx, DestroyGLVisViewerCtx_Private));
1359: for (i = 0; i < tnf; i++) PetscCall(PetscFree(fecs[i]));
1360: PetscCall(PetscFree3(fecs, sdim, Ufds));
1361: PetscFunctionReturn(PETSC_SUCCESS);
1362: }
1364: static PetscErrorCode DMRefine_Composite(DM dmi, MPI_Comm comm, DM *fine)
1365: {
1366: struct DMCompositeLink *next;
1367: DM_Composite *com = (DM_Composite *)dmi->data;
1368: DM dm;
1370: PetscFunctionBegin;
1372: if (comm == MPI_COMM_NULL) PetscCall(PetscObjectGetComm((PetscObject)dmi, &comm));
1373: PetscCall(DMSetUp(dmi));
1374: next = com->next;
1375: PetscCall(DMCompositeCreate(comm, fine));
1377: /* loop over packed objects, handling one at at time */
1378: while (next) {
1379: PetscCall(DMRefine(next->dm, comm, &dm));
1380: PetscCall(DMCompositeAddDM(*fine, dm));
1381: PetscCall(PetscObjectDereference((PetscObject)dm));
1382: next = next->next;
1383: }
1384: PetscFunctionReturn(PETSC_SUCCESS);
1385: }
1387: static PetscErrorCode DMCoarsen_Composite(DM dmi, MPI_Comm comm, DM *fine)
1388: {
1389: struct DMCompositeLink *next;
1390: DM_Composite *com = (DM_Composite *)dmi->data;
1391: DM dm;
1393: PetscFunctionBegin;
1395: PetscCall(DMSetUp(dmi));
1396: if (comm == MPI_COMM_NULL) PetscCall(PetscObjectGetComm((PetscObject)dmi, &comm));
1397: next = com->next;
1398: PetscCall(DMCompositeCreate(comm, fine));
1400: /* loop over packed objects, handling one at at time */
1401: while (next) {
1402: PetscCall(DMCoarsen(next->dm, comm, &dm));
1403: PetscCall(DMCompositeAddDM(*fine, dm));
1404: PetscCall(PetscObjectDereference((PetscObject)dm));
1405: next = next->next;
1406: }
1407: PetscFunctionReturn(PETSC_SUCCESS);
1408: }
1410: static PetscErrorCode DMCreateInterpolation_Composite(DM coarse, DM fine, Mat *A, Vec *v)
1411: {
1412: PetscInt m, n, M, N, nDM, i;
1413: struct DMCompositeLink *nextc;
1414: struct DMCompositeLink *nextf;
1415: Vec gcoarse, gfine, *vecs;
1416: DM_Composite *comcoarse = (DM_Composite *)coarse->data;
1417: DM_Composite *comfine = (DM_Composite *)fine->data;
1418: Mat *mats;
1420: PetscFunctionBegin;
1423: PetscCall(DMSetUp(coarse));
1424: PetscCall(DMSetUp(fine));
1425: /* use global vectors only for determining matrix layout */
1426: PetscCall(DMGetGlobalVector(coarse, &gcoarse));
1427: PetscCall(DMGetGlobalVector(fine, &gfine));
1428: PetscCall(VecGetLocalSize(gcoarse, &n));
1429: PetscCall(VecGetLocalSize(gfine, &m));
1430: PetscCall(VecGetSize(gcoarse, &N));
1431: PetscCall(VecGetSize(gfine, &M));
1432: PetscCall(DMRestoreGlobalVector(coarse, &gcoarse));
1433: PetscCall(DMRestoreGlobalVector(fine, &gfine));
1435: nDM = comfine->nDM;
1436: 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);
1437: PetscCall(PetscCalloc1(nDM * nDM, &mats));
1438: if (v) PetscCall(PetscCalloc1(nDM, &vecs));
1440: /* loop over packed objects, handling one at at time */
1441: for (nextc = comcoarse->next, nextf = comfine->next, i = 0; nextc; nextc = nextc->next, nextf = nextf->next, i++) {
1442: if (!v) PetscCall(DMCreateInterpolation(nextc->dm, nextf->dm, &mats[i * nDM + i], NULL));
1443: else PetscCall(DMCreateInterpolation(nextc->dm, nextf->dm, &mats[i * nDM + i], &vecs[i]));
1444: }
1445: PetscCall(MatCreateNest(PetscObjectComm((PetscObject)fine), nDM, NULL, nDM, NULL, mats, A));
1446: if (v) PetscCall(VecCreateNest(PetscObjectComm((PetscObject)fine), nDM, NULL, vecs, v));
1447: for (i = 0; i < nDM * nDM; i++) PetscCall(MatDestroy(&mats[i]));
1448: PetscCall(PetscFree(mats));
1449: if (v) {
1450: for (i = 0; i < nDM; i++) PetscCall(VecDestroy(&vecs[i]));
1451: PetscCall(PetscFree(vecs));
1452: }
1453: PetscFunctionReturn(PETSC_SUCCESS);
1454: }
1456: static PetscErrorCode DMGetLocalToGlobalMapping_Composite(DM dm)
1457: {
1458: DM_Composite *com = (DM_Composite *)dm->data;
1459: ISLocalToGlobalMapping *ltogs;
1460: PetscInt i;
1462: PetscFunctionBegin;
1463: /* Set the ISLocalToGlobalMapping on the new matrix */
1464: PetscCall(DMCompositeGetISLocalToGlobalMappings(dm, <ogs));
1465: PetscCall(ISLocalToGlobalMappingConcatenate(PetscObjectComm((PetscObject)dm), com->nDM, ltogs, &dm->ltogmap));
1466: for (i = 0; i < com->nDM; i++) PetscCall(ISLocalToGlobalMappingDestroy(<ogs[i]));
1467: PetscCall(PetscFree(ltogs));
1468: PetscFunctionReturn(PETSC_SUCCESS);
1469: }
1471: static PetscErrorCode DMCreateColoring_Composite(DM dm, ISColoringType ctype, ISColoring *coloring)
1472: {
1473: PetscInt n, i, cnt;
1474: ISColoringValue *colors;
1475: PetscBool dense = PETSC_FALSE;
1476: ISColoringValue maxcol = 0;
1477: DM_Composite *com = (DM_Composite *)dm->data;
1479: PetscFunctionBegin;
1481: PetscCheck(ctype != IS_COLORING_LOCAL, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Only global coloring supported");
1482: if (ctype == IS_COLORING_GLOBAL) {
1483: n = com->n;
1484: } else SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "Unknown ISColoringType");
1485: PetscCall(PetscMalloc1(n, &colors)); /* freed in ISColoringDestroy() */
1487: PetscCall(PetscOptionsGetBool(((PetscObject)dm)->options, ((PetscObject)dm)->prefix, "-dmcomposite_dense_jacobian", &dense, NULL));
1488: if (dense) {
1489: for (i = 0; i < n; i++) colors[i] = (ISColoringValue)(com->rstart + i);
1490: maxcol = com->N;
1491: } else {
1492: struct DMCompositeLink *next = com->next;
1493: PetscMPIInt rank;
1495: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
1496: cnt = 0;
1497: while (next) {
1498: ISColoring lcoloring;
1500: PetscCall(DMCreateColoring(next->dm, IS_COLORING_GLOBAL, &lcoloring));
1501: for (i = 0; i < lcoloring->N; i++) colors[cnt++] = maxcol + lcoloring->colors[i];
1502: maxcol += lcoloring->n;
1503: PetscCall(ISColoringDestroy(&lcoloring));
1504: next = next->next;
1505: }
1506: }
1507: PetscCall(ISColoringCreate(PetscObjectComm((PetscObject)dm), maxcol, n, colors, PETSC_OWN_POINTER, coloring));
1508: PetscFunctionReturn(PETSC_SUCCESS);
1509: }
1511: static PetscErrorCode DMGlobalToLocalBegin_Composite(DM dm, Vec gvec, InsertMode mode, Vec lvec)
1512: {
1513: struct DMCompositeLink *next;
1514: PetscScalar *garray, *larray;
1515: DM_Composite *com = (DM_Composite *)dm->data;
1517: PetscFunctionBegin;
1521: if (!com->setup) PetscCall(DMSetUp(dm));
1523: PetscCall(VecGetArray(gvec, &garray));
1524: PetscCall(VecGetArray(lvec, &larray));
1526: /* loop over packed objects, handling one at at time */
1527: next = com->next;
1528: while (next) {
1529: Vec local, global;
1530: PetscInt N;
1532: PetscCall(DMGetGlobalVector(next->dm, &global));
1533: PetscCall(VecGetLocalSize(global, &N));
1534: PetscCall(VecPlaceArray(global, garray));
1535: PetscCall(DMGetLocalVector(next->dm, &local));
1536: PetscCall(VecPlaceArray(local, larray));
1537: PetscCall(DMGlobalToLocalBegin(next->dm, global, mode, local));
1538: PetscCall(DMGlobalToLocalEnd(next->dm, global, mode, local));
1539: PetscCall(VecResetArray(global));
1540: PetscCall(VecResetArray(local));
1541: PetscCall(DMRestoreGlobalVector(next->dm, &global));
1542: PetscCall(DMRestoreLocalVector(next->dm, &local));
1544: larray += next->nlocal;
1545: garray += next->n;
1546: next = next->next;
1547: }
1549: PetscCall(VecRestoreArray(gvec, NULL));
1550: PetscCall(VecRestoreArray(lvec, NULL));
1551: PetscFunctionReturn(PETSC_SUCCESS);
1552: }
1554: static PetscErrorCode DMGlobalToLocalEnd_Composite(DM dm, Vec gvec, InsertMode mode, Vec lvec)
1555: {
1556: PetscFunctionBegin;
1560: PetscFunctionReturn(PETSC_SUCCESS);
1561: }
1563: static PetscErrorCode DMLocalToGlobalBegin_Composite(DM dm, Vec lvec, InsertMode mode, Vec gvec)
1564: {
1565: struct DMCompositeLink *next;
1566: PetscScalar *larray, *garray;
1567: DM_Composite *com = (DM_Composite *)dm->data;
1569: PetscFunctionBegin;
1574: if (!com->setup) PetscCall(DMSetUp(dm));
1576: PetscCall(VecGetArray(lvec, &larray));
1577: PetscCall(VecGetArray(gvec, &garray));
1579: /* loop over packed objects, handling one at at time */
1580: next = com->next;
1581: while (next) {
1582: Vec global, local;
1584: PetscCall(DMGetLocalVector(next->dm, &local));
1585: PetscCall(VecPlaceArray(local, larray));
1586: PetscCall(DMGetGlobalVector(next->dm, &global));
1587: PetscCall(VecPlaceArray(global, garray));
1588: PetscCall(DMLocalToGlobalBegin(next->dm, local, mode, global));
1589: PetscCall(DMLocalToGlobalEnd(next->dm, local, mode, global));
1590: PetscCall(VecResetArray(local));
1591: PetscCall(VecResetArray(global));
1592: PetscCall(DMRestoreGlobalVector(next->dm, &global));
1593: PetscCall(DMRestoreLocalVector(next->dm, &local));
1595: garray += next->n;
1596: larray += next->nlocal;
1597: next = next->next;
1598: }
1600: PetscCall(VecRestoreArray(gvec, NULL));
1601: PetscCall(VecRestoreArray(lvec, NULL));
1603: PetscFunctionReturn(PETSC_SUCCESS);
1604: }
1606: static PetscErrorCode DMLocalToGlobalEnd_Composite(DM dm, Vec lvec, InsertMode mode, Vec gvec)
1607: {
1608: PetscFunctionBegin;
1612: PetscFunctionReturn(PETSC_SUCCESS);
1613: }
1615: static PetscErrorCode DMLocalToLocalBegin_Composite(DM dm, Vec vec1, InsertMode mode, Vec vec2)
1616: {
1617: struct DMCompositeLink *next;
1618: PetscScalar *array1, *array2;
1619: DM_Composite *com = (DM_Composite *)dm->data;
1621: PetscFunctionBegin;
1626: if (!com->setup) PetscCall(DMSetUp(dm));
1628: PetscCall(VecGetArray(vec1, &array1));
1629: PetscCall(VecGetArray(vec2, &array2));
1631: /* loop over packed objects, handling one at at time */
1632: next = com->next;
1633: while (next) {
1634: Vec local1, local2;
1636: PetscCall(DMGetLocalVector(next->dm, &local1));
1637: PetscCall(VecPlaceArray(local1, array1));
1638: PetscCall(DMGetLocalVector(next->dm, &local2));
1639: PetscCall(VecPlaceArray(local2, array2));
1640: PetscCall(DMLocalToLocalBegin(next->dm, local1, mode, local2));
1641: PetscCall(DMLocalToLocalEnd(next->dm, local1, mode, local2));
1642: PetscCall(VecResetArray(local2));
1643: PetscCall(DMRestoreLocalVector(next->dm, &local2));
1644: PetscCall(VecResetArray(local1));
1645: PetscCall(DMRestoreLocalVector(next->dm, &local1));
1647: array1 += next->nlocal;
1648: array2 += next->nlocal;
1649: next = next->next;
1650: }
1652: PetscCall(VecRestoreArray(vec1, NULL));
1653: PetscCall(VecRestoreArray(vec2, NULL));
1655: PetscFunctionReturn(PETSC_SUCCESS);
1656: }
1658: static PetscErrorCode DMLocalToLocalEnd_Composite(DM dm, Vec lvec, InsertMode mode, Vec gvec)
1659: {
1660: PetscFunctionBegin;
1664: PetscFunctionReturn(PETSC_SUCCESS);
1665: }
1667: /*MC
1668: DMCOMPOSITE = "composite" - A `DM` object that is used to manage data for a collection of `DM`
1670: Level: intermediate
1672: .seealso: `DMType`, `DM`, `DMDACreate()`, `DMCreate()`, `DMSetType()`, `DMCompositeCreate()`
1673: M*/
1675: PETSC_EXTERN PetscErrorCode DMCreate_Composite(DM p)
1676: {
1677: DM_Composite *com;
1679: PetscFunctionBegin;
1680: PetscCall(PetscNew(&com));
1681: p->data = com;
1682: com->n = 0;
1683: com->nghost = 0;
1684: com->next = NULL;
1685: com->nDM = 0;
1687: p->ops->createglobalvector = DMCreateGlobalVector_Composite;
1688: p->ops->createlocalvector = DMCreateLocalVector_Composite;
1689: p->ops->getlocaltoglobalmapping = DMGetLocalToGlobalMapping_Composite;
1690: p->ops->createfieldis = DMCreateFieldIS_Composite;
1691: p->ops->createfielddecomposition = DMCreateFieldDecomposition_Composite;
1692: p->ops->refine = DMRefine_Composite;
1693: p->ops->coarsen = DMCoarsen_Composite;
1694: p->ops->createinterpolation = DMCreateInterpolation_Composite;
1695: p->ops->creatematrix = DMCreateMatrix_Composite;
1696: p->ops->getcoloring = DMCreateColoring_Composite;
1697: p->ops->globaltolocalbegin = DMGlobalToLocalBegin_Composite;
1698: p->ops->globaltolocalend = DMGlobalToLocalEnd_Composite;
1699: p->ops->localtoglobalbegin = DMLocalToGlobalBegin_Composite;
1700: p->ops->localtoglobalend = DMLocalToGlobalEnd_Composite;
1701: p->ops->localtolocalbegin = DMLocalToLocalBegin_Composite;
1702: p->ops->localtolocalend = DMLocalToLocalEnd_Composite;
1703: p->ops->destroy = DMDestroy_Composite;
1704: p->ops->view = DMView_Composite;
1705: p->ops->setup = DMSetUp_Composite;
1707: PetscCall(PetscObjectComposeFunction((PetscObject)p, "DMSetUpGLVisViewer_C", DMSetUpGLVisViewer_Composite));
1708: PetscFunctionReturn(PETSC_SUCCESS);
1709: }
1711: /*@
1712: DMCompositeCreate - Creates a `DMCOMPOSITE`, used to generate "composite"
1713: vectors made up of several subvectors.
1715: Collective
1717: Input Parameter:
1718: . comm - the processors that will share the global vector
1720: Output Parameter:
1721: . packer - the `DMCOMPOSITE` object
1723: Level: advanced
1725: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCompositeScatter()`, `DMCreate()`
1726: `DMCompositeGather()`, `DMCreateGlobalVector()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`
1727: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
1728: @*/
1729: PetscErrorCode DMCompositeCreate(MPI_Comm comm, DM *packer)
1730: {
1731: PetscFunctionBegin;
1732: PetscAssertPointer(packer, 2);
1733: PetscCall(DMCreate(comm, packer));
1734: PetscCall(DMSetType(*packer, DMCOMPOSITE));
1735: PetscFunctionReturn(PETSC_SUCCESS);
1736: }