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 a 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: }
337: PetscFunctionReturn(PETSC_SUCCESS);
338: }
340: /*@C
341: DMCompositeRestoreAccess - Returns the vectors obtained with `DMCompositeGetAccess()`
342: representation.
344: Collective
346: Input Parameters:
347: + dm - the `DMCOMPOSITE` object
348: . gvec - the global vector
349: - ... - the individual parallel vectors, `NULL` for those that are not needed
351: Level: advanced
353: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
354: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeScatter()`,
355: `DMCompositeGetAccess()`
356: @*/
357: PetscErrorCode DMCompositeRestoreAccess(DM dm, Vec gvec, ...)
358: {
359: va_list Argp;
360: struct DMCompositeLink *next;
361: DM_Composite *com = (DM_Composite *)dm->data;
362: PetscInt readonly;
363: PetscBool flg;
365: PetscFunctionBegin;
368: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
369: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
370: next = com->next;
371: if (!com->setup) PetscCall(DMSetUp(dm));
373: PetscCall(VecLockGet(gvec, &readonly));
374: /* loop over packed objects, handling one at a time */
375: va_start(Argp, gvec);
376: while (next) {
377: Vec *vec;
378: vec = va_arg(Argp, Vec *);
379: if (vec) {
380: PetscCall(VecResetArray(*vec));
381: if (readonly) PetscCall(VecLockReadPop(*vec));
382: PetscCall(DMRestoreGlobalVector(next->dm, vec));
383: }
384: next = next->next;
385: }
386: va_end(Argp);
387: PetscFunctionReturn(PETSC_SUCCESS);
388: }
390: /*@C
391: DMCompositeRestoreAccessArray - Returns the vectors obtained with `DMCompositeGetAccessArray()`
393: Collective
395: Input Parameters:
396: + dm - the `DMCOMPOSITE` object
397: . pvec - packed vector
398: . nwanted - number of vectors wanted
399: . wanted - sorted array of vectors wanted, or NULL to get all vectors
400: - vecs - array of global vectors to return
402: Level: advanced
404: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeRestoreAccess()`, `DMCompositeRestoreEntries()`, `DMCompositeScatter()`, `DMCompositeGather()`
405: @*/
406: PetscErrorCode DMCompositeRestoreAccessArray(DM dm, Vec pvec, PetscInt nwanted, const PetscInt *wanted, Vec *vecs)
407: {
408: struct DMCompositeLink *link;
409: PetscInt i, wnum;
410: DM_Composite *com = (DM_Composite *)dm->data;
411: PetscInt readonly;
412: PetscBool flg;
414: PetscFunctionBegin;
417: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
418: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
419: if (!com->setup) PetscCall(DMSetUp(dm));
421: PetscCall(VecLockGet(pvec, &readonly));
422: for (i = 0, wnum = 0, link = com->next; link && wnum < nwanted; i++, link = link->next) {
423: if (!wanted || i == wanted[wnum]) {
424: PetscCall(VecResetArray(vecs[wnum]));
425: if (readonly) PetscCall(VecLockReadPop(vecs[wnum]));
426: PetscCall(DMRestoreGlobalVector(link->dm, &vecs[wnum]));
427: wnum++;
428: }
429: }
430: PetscFunctionReturn(PETSC_SUCCESS);
431: }
433: /*@C
434: DMCompositeRestoreLocalAccessArray - Returns the vectors obtained with `DMCompositeGetLocalAccessArray()`.
436: Collective
438: Input Parameters:
439: + dm - the `DMCOMPOSITE` object
440: . pvec - packed vector
441: . nwanted - number of vectors wanted
442: . wanted - sorted array of vectors wanted, or NULL to restore all vectors
443: - vecs - array of local vectors to return
445: Level: advanced
447: Note:
448: nwanted and wanted must match the values given to `DMCompositeGetLocalAccessArray()`
449: otherwise the call will fail.
451: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeGetLocalAccessArray()`, `DMCompositeRestoreAccessArray()`,
452: `DMCompositeRestoreAccess()`, `DMCompositeRestoreEntries()`,
453: `DMCompositeScatter()`, `DMCompositeGather()`
454: @*/
455: PetscErrorCode DMCompositeRestoreLocalAccessArray(DM dm, Vec pvec, PetscInt nwanted, const PetscInt *wanted, Vec *vecs)
456: {
457: struct DMCompositeLink *link;
458: PetscInt i, wnum;
459: DM_Composite *com = (DM_Composite *)dm->data;
460: PetscInt readonly;
461: PetscBool flg;
463: PetscFunctionBegin;
466: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
467: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
468: if (!com->setup) PetscCall(DMSetUp(dm));
470: PetscCall(VecLockGet(pvec, &readonly));
471: for (i = 0, wnum = 0, link = com->next; link && wnum < nwanted; i++, link = link->next) {
472: if (!wanted || i == wanted[wnum]) {
473: PetscCall(VecResetArray(vecs[wnum]));
474: if (readonly) PetscCall(VecLockReadPop(vecs[wnum]));
475: PetscCall(DMRestoreLocalVector(link->dm, &vecs[wnum]));
476: wnum++;
477: }
478: }
479: PetscFunctionReturn(PETSC_SUCCESS);
480: }
482: /*@C
483: DMCompositeScatter - Scatters from a global packed vector into its individual local vectors
485: Collective
487: Input Parameters:
488: + dm - the `DMCOMPOSITE` object
489: . gvec - the global vector
490: - ... - the individual sequential vectors, `NULL` for those that are not needed
492: Level: advanced
494: Note:
495: `DMCompositeScatterArray()` is a non-variadic alternative that is often more convenient for library callers and is
496: accessible from Fortran.
498: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
499: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
500: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
501: `DMCompositeScatterArray()`
502: @*/
503: PetscErrorCode DMCompositeScatter(DM dm, Vec gvec, ...)
504: {
505: va_list Argp;
506: struct DMCompositeLink *next;
507: PETSC_UNUSED PetscInt cnt;
508: DM_Composite *com = (DM_Composite *)dm->data;
509: PetscBool flg;
511: PetscFunctionBegin;
514: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
515: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
516: if (!com->setup) PetscCall(DMSetUp(dm));
518: /* loop over packed objects, handling one at a time */
519: va_start(Argp, gvec);
520: for (cnt = 3, next = com->next; next; cnt++, next = next->next) {
521: Vec local;
522: local = va_arg(Argp, Vec);
523: if (local) {
524: Vec global;
525: const PetscScalar *array;
527: PetscCall(DMGetGlobalVector(next->dm, &global));
528: PetscCall(VecGetArrayRead(gvec, &array));
529: PetscCall(VecPlaceArray(global, array + next->rstart));
530: PetscCall(DMGlobalToLocalBegin(next->dm, global, INSERT_VALUES, local));
531: PetscCall(DMGlobalToLocalEnd(next->dm, global, INSERT_VALUES, local));
532: PetscCall(VecRestoreArrayRead(gvec, &array));
533: PetscCall(VecResetArray(global));
534: PetscCall(DMRestoreGlobalVector(next->dm, &global));
535: }
536: }
537: va_end(Argp);
538: PetscFunctionReturn(PETSC_SUCCESS);
539: }
541: /*@
542: DMCompositeScatterArray - Scatters from a global packed vector into its individual local vectors
544: Collective
546: Input Parameters:
547: + dm - the `DMCOMPOSITE` object
548: . gvec - the global vector
549: - lvecs - array of local vectors, NULL for any that are not needed
551: Level: advanced
553: Note:
554: This is a non-variadic alternative to `DMCompositeScatter()`
556: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`
557: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
558: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
559: @*/
560: PetscErrorCode DMCompositeScatterArray(DM dm, Vec gvec, Vec *lvecs)
561: {
562: struct DMCompositeLink *next;
563: PetscInt i;
564: DM_Composite *com = (DM_Composite *)dm->data;
565: PetscBool flg;
567: PetscFunctionBegin;
570: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
571: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
572: if (!com->setup) PetscCall(DMSetUp(dm));
574: /* loop over packed objects, handling one at a time */
575: for (i = 0, next = com->next; next; next = next->next, i++) {
576: if (lvecs[i]) {
577: Vec global;
578: const PetscScalar *array;
580: PetscCall(DMGetGlobalVector(next->dm, &global));
581: PetscCall(VecGetArrayRead(gvec, &array));
582: PetscCall(VecPlaceArray(global, (PetscScalar *)array + next->rstart));
583: PetscCall(DMGlobalToLocalBegin(next->dm, global, INSERT_VALUES, lvecs[i]));
584: PetscCall(DMGlobalToLocalEnd(next->dm, global, INSERT_VALUES, lvecs[i]));
585: PetscCall(VecRestoreArrayRead(gvec, &array));
586: PetscCall(VecResetArray(global));
587: PetscCall(DMRestoreGlobalVector(next->dm, &global));
588: }
589: }
590: PetscFunctionReturn(PETSC_SUCCESS);
591: }
593: /*@C
594: DMCompositeGather - Gathers into a global packed vector from its individual local vectors
596: Collective
598: Input Parameters:
599: + dm - the `DMCOMPOSITE` object
600: . imode - `INSERT_VALUES` or `ADD_VALUES`
601: . gvec - the global vector
602: - ... - the individual sequential vectors, `NULL` for any that are not needed
604: Level: advanced
606: Fortran Notes:
607: Fortran users should use `DMCompositeGatherArray()`
609: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
610: `DMCompositeScatter()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
611: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
612: @*/
613: PetscErrorCode DMCompositeGather(DM dm, InsertMode imode, Vec gvec, ...)
614: {
615: va_list Argp;
616: struct DMCompositeLink *next;
617: DM_Composite *com = (DM_Composite *)dm->data;
618: PETSC_UNUSED PetscInt cnt;
619: PetscBool flg;
621: PetscFunctionBegin;
624: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
625: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
626: if (!com->setup) PetscCall(DMSetUp(dm));
628: /* loop over packed objects, handling one at a time */
629: va_start(Argp, gvec);
630: for (cnt = 3, next = com->next; next; cnt++, next = next->next) {
631: Vec local;
632: local = va_arg(Argp, Vec);
633: if (local) {
634: PetscScalar *array;
635: Vec global;
637: PetscCall(DMGetGlobalVector(next->dm, &global));
638: PetscCall(VecGetArray(gvec, &array));
639: PetscCall(VecPlaceArray(global, array + next->rstart));
640: PetscCall(DMLocalToGlobalBegin(next->dm, local, imode, global));
641: PetscCall(DMLocalToGlobalEnd(next->dm, local, imode, global));
642: PetscCall(VecRestoreArray(gvec, &array));
643: PetscCall(VecResetArray(global));
644: PetscCall(DMRestoreGlobalVector(next->dm, &global));
645: }
646: }
647: va_end(Argp);
648: PetscFunctionReturn(PETSC_SUCCESS);
649: }
651: /*@
652: DMCompositeGatherArray - Gathers into a global packed vector from its individual local vectors
654: Collective
656: Input Parameters:
657: + dm - the `DMCOMPOSITE` object
658: . gvec - the global vector
659: . imode - `INSERT_VALUES` or `ADD_VALUES`
660: - lvecs - the individual sequential vectors, NULL for any that are not needed
662: Level: advanced
664: Note:
665: This is a non-variadic alternative to `DMCompositeGather()`.
667: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
668: `DMCompositeScatter()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
669: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`,
670: @*/
671: PetscErrorCode DMCompositeGatherArray(DM dm, InsertMode imode, Vec gvec, Vec *lvecs)
672: {
673: struct DMCompositeLink *next;
674: DM_Composite *com = (DM_Composite *)dm->data;
675: PetscInt i;
676: PetscBool flg;
678: PetscFunctionBegin;
681: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
682: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
683: if (!com->setup) PetscCall(DMSetUp(dm));
685: /* loop over packed objects, handling one at a time */
686: for (next = com->next, i = 0; next; next = next->next, i++) {
687: if (lvecs[i]) {
688: PetscScalar *array;
689: Vec global;
691: PetscCall(DMGetGlobalVector(next->dm, &global));
692: PetscCall(VecGetArray(gvec, &array));
693: PetscCall(VecPlaceArray(global, array + next->rstart));
694: PetscCall(DMLocalToGlobalBegin(next->dm, lvecs[i], imode, global));
695: PetscCall(DMLocalToGlobalEnd(next->dm, lvecs[i], imode, global));
696: PetscCall(VecRestoreArray(gvec, &array));
697: PetscCall(VecResetArray(global));
698: PetscCall(DMRestoreGlobalVector(next->dm, &global));
699: }
700: }
701: PetscFunctionReturn(PETSC_SUCCESS);
702: }
704: /*@
705: DMCompositeAddDM - adds a `DM` vector to a `DMCOMPOSITE`
707: Collective
709: Input Parameters:
710: + dmc - the `DMCOMPOSITE` object
711: - dm - the `DM` object
713: Level: advanced
715: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeGather()`, `DMCreateGlobalVector()`,
716: `DMCompositeScatter()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
717: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
718: @*/
719: PetscErrorCode DMCompositeAddDM(DM dmc, DM dm)
720: {
721: PetscInt n, nlocal;
722: struct DMCompositeLink *mine, *next;
723: Vec global, local;
724: DM_Composite *com = (DM_Composite *)dmc->data;
725: PetscBool flg;
727: PetscFunctionBegin;
730: PetscCall(PetscObjectTypeCompare((PetscObject)dmc, DMCOMPOSITE, &flg));
731: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
732: next = com->next;
733: PetscCheck(!com->setup, PetscObjectComm((PetscObject)dmc), PETSC_ERR_ARG_WRONGSTATE, "Cannot add a DM once you have used the DMComposite");
735: /* create new link */
736: PetscCall(PetscNew(&mine));
737: PetscCall(PetscObjectReference((PetscObject)dm));
738: PetscCall(DMGetGlobalVector(dm, &global));
739: PetscCall(VecGetLocalSize(global, &n));
740: PetscCall(DMRestoreGlobalVector(dm, &global));
741: PetscCall(DMGetLocalVector(dm, &local));
742: PetscCall(VecGetSize(local, &nlocal));
743: PetscCall(DMRestoreLocalVector(dm, &local));
745: mine->n = n;
746: mine->nlocal = nlocal;
747: mine->dm = dm;
748: mine->next = NULL;
749: com->n += n;
750: com->nghost += nlocal;
752: /* add to end of list */
753: if (!next) com->next = mine;
754: else {
755: while (next->next) next = next->next;
756: next->next = mine;
757: }
758: com->nDM++;
759: com->nmine++;
760: PetscFunctionReturn(PETSC_SUCCESS);
761: }
763: #include <petscdraw.h>
764: PETSC_EXTERN PetscErrorCode VecView_MPI(Vec, PetscViewer);
765: static PetscErrorCode VecView_DMComposite(Vec gvec, PetscViewer viewer)
766: {
767: DM dm;
768: struct DMCompositeLink *next;
769: PetscBool isdraw;
770: DM_Composite *com;
772: PetscFunctionBegin;
773: PetscCall(VecGetDM(gvec, &dm));
774: PetscCheck(dm, PetscObjectComm((PetscObject)gvec), PETSC_ERR_ARG_WRONG, "Vector not generated from a DMComposite");
775: com = (DM_Composite *)dm->data;
776: next = com->next;
778: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERDRAW, &isdraw));
779: if (!isdraw) {
780: /* do I really want to call this? */
781: PetscCall(VecView_MPI(gvec, viewer));
782: } else {
783: PetscInt cnt = 0;
785: /* loop over packed objects, handling one at a time */
786: while (next) {
787: Vec vec;
788: const PetscScalar *array;
789: PetscInt bs;
791: /* Should use VecGetSubVector() eventually, but would need to forward the DM for that to work */
792: PetscCall(DMGetGlobalVector(next->dm, &vec));
793: PetscCall(VecGetArrayRead(gvec, &array));
794: PetscCall(VecPlaceArray(vec, (PetscScalar *)array + next->rstart));
795: PetscCall(VecRestoreArrayRead(gvec, &array));
796: PetscCall(VecView(vec, viewer));
797: PetscCall(VecResetArray(vec));
798: PetscCall(VecGetBlockSize(vec, &bs));
799: PetscCall(DMRestoreGlobalVector(next->dm, &vec));
800: PetscCall(PetscViewerDrawBaseAdd(viewer, bs));
801: cnt += bs;
802: next = next->next;
803: }
804: PetscCall(PetscViewerDrawBaseAdd(viewer, -cnt));
805: }
806: PetscFunctionReturn(PETSC_SUCCESS);
807: }
809: static PetscErrorCode DMCreateGlobalVector_Composite(DM dm, Vec *gvec)
810: {
811: DM_Composite *com = (DM_Composite *)dm->data;
813: PetscFunctionBegin;
815: PetscCall(DMSetFromOptions(dm));
816: PetscCall(DMSetUp(dm));
817: PetscCall(VecCreate(PetscObjectComm((PetscObject)dm), gvec));
818: PetscCall(VecSetType(*gvec, dm->vectype));
819: PetscCall(VecSetSizes(*gvec, com->n, com->N));
820: PetscCall(VecSetDM(*gvec, dm));
821: PetscCall(VecSetOperation(*gvec, VECOP_VIEW, (void (*)(void))VecView_DMComposite));
822: PetscFunctionReturn(PETSC_SUCCESS);
823: }
825: static PetscErrorCode DMCreateLocalVector_Composite(DM dm, Vec *lvec)
826: {
827: DM_Composite *com = (DM_Composite *)dm->data;
829: PetscFunctionBegin;
831: if (!com->setup) {
832: PetscCall(DMSetFromOptions(dm));
833: PetscCall(DMSetUp(dm));
834: }
835: PetscCall(VecCreate(PETSC_COMM_SELF, lvec));
836: PetscCall(VecSetType(*lvec, dm->vectype));
837: PetscCall(VecSetSizes(*lvec, com->nghost, PETSC_DECIDE));
838: PetscCall(VecSetDM(*lvec, dm));
839: PetscFunctionReturn(PETSC_SUCCESS);
840: }
842: /*@C
843: DMCompositeGetISLocalToGlobalMappings - gets an `ISLocalToGlobalMapping` for each `DM` in the `DMCOMPOSITE`, maps to the composite global space
845: Collective; No Fortran Support
847: Input Parameter:
848: . dm - the `DMCOMPOSITE` object
850: Output Parameter:
851: . ltogs - the individual mappings for each packed vector. Note that this includes
852: all the ghost points that individual ghosted `DMDA` may have.
854: Level: advanced
856: Note:
857: Each entry of ltogs should be destroyed with `ISLocalToGlobalMappingDestroy()`, the ltogs array should be freed with `PetscFree()`.
859: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
860: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetAccess()`, `DMCompositeScatter()`,
861: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
862: @*/
863: PetscErrorCode DMCompositeGetISLocalToGlobalMappings(DM dm, ISLocalToGlobalMapping **ltogs)
864: {
865: PetscInt i, *idx, n, cnt;
866: struct DMCompositeLink *next;
867: PetscMPIInt rank;
868: DM_Composite *com = (DM_Composite *)dm->data;
869: PetscBool flg;
871: PetscFunctionBegin;
873: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
874: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
875: PetscCall(DMSetUp(dm));
876: PetscCall(PetscMalloc1(com->nDM, ltogs));
877: next = com->next;
878: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
880: /* loop over packed objects, handling one at a time */
881: cnt = 0;
882: while (next) {
883: ISLocalToGlobalMapping ltog;
884: PetscMPIInt size;
885: const PetscInt *suboff, *indices;
886: Vec global;
888: /* Get sub-DM global indices for each local dof */
889: PetscCall(DMGetLocalToGlobalMapping(next->dm, <og));
890: PetscCall(ISLocalToGlobalMappingGetSize(ltog, &n));
891: PetscCall(ISLocalToGlobalMappingGetIndices(ltog, &indices));
892: PetscCall(PetscMalloc1(n, &idx));
894: /* Get the offsets for the sub-DM global vector */
895: PetscCall(DMGetGlobalVector(next->dm, &global));
896: PetscCall(VecGetOwnershipRanges(global, &suboff));
897: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)global), &size));
899: /* Shift the sub-DM definition of the global space to the composite global space */
900: for (i = 0; i < n; i++) {
901: PetscInt subi = indices[i], lo = 0, hi = size, t;
902: /* There's no consensus on what a negative index means,
903: except for skipping when setting the values in vectors and matrices */
904: if (subi < 0) {
905: idx[i] = subi - next->grstarts[rank];
906: continue;
907: }
908: /* Binary search to find which rank owns subi */
909: while (hi - lo > 1) {
910: t = lo + (hi - lo) / 2;
911: if (suboff[t] > subi) hi = t;
912: else lo = t;
913: }
914: idx[i] = subi - suboff[lo] + next->grstarts[lo];
915: }
916: PetscCall(ISLocalToGlobalMappingRestoreIndices(ltog, &indices));
917: PetscCall(ISLocalToGlobalMappingCreate(PetscObjectComm((PetscObject)dm), 1, n, idx, PETSC_OWN_POINTER, &(*ltogs)[cnt]));
918: PetscCall(DMRestoreGlobalVector(next->dm, &global));
919: next = next->next;
920: cnt++;
921: }
922: PetscFunctionReturn(PETSC_SUCCESS);
923: }
925: /*@C
926: DMCompositeGetLocalISs - Gets index sets for each component of a composite local vector
928: Not Collective; No Fortran Support
930: Input Parameter:
931: . dm - the `DMCOMPOSITE`
933: Output Parameter:
934: . is - array of serial index sets for each component of the `DMCOMPOSITE`
936: Level: intermediate
938: Notes:
939: At present, a composite local vector does not normally exist. This function is used to provide index sets for
940: `MatGetLocalSubMatrix()`. In the future, the scatters for each entry in the `DMCOMPOSITE` may be merged into a single
941: scatter to a composite local vector. The user should not typically need to know which is being done.
943: To get the composite global indices at all local points (including ghosts), use `DMCompositeGetISLocalToGlobalMappings()`.
945: To get index sets for pieces of the composite global vector, use `DMCompositeGetGlobalISs()`.
947: Each returned `IS` should be destroyed with `ISDestroy()`, the array should be freed with `PetscFree()`.
949: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeGetGlobalISs()`, `DMCompositeGetISLocalToGlobalMappings()`, `MatGetLocalSubMatrix()`, `MatCreateLocalRef()`
950: @*/
951: PetscErrorCode DMCompositeGetLocalISs(DM dm, IS **is)
952: {
953: DM_Composite *com = (DM_Composite *)dm->data;
954: struct DMCompositeLink *link;
955: PetscInt cnt, start;
956: PetscBool flg;
958: PetscFunctionBegin;
960: PetscAssertPointer(is, 2);
961: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
962: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
963: PetscCall(PetscMalloc1(com->nmine, is));
964: for (cnt = 0, start = 0, link = com->next; link; start += link->nlocal, cnt++, link = link->next) {
965: PetscInt bs;
966: PetscCall(ISCreateStride(PETSC_COMM_SELF, link->nlocal, start, 1, &(*is)[cnt]));
967: PetscCall(DMGetBlockSize(link->dm, &bs));
968: PetscCall(ISSetBlockSize((*is)[cnt], bs));
969: }
970: PetscFunctionReturn(PETSC_SUCCESS);
971: }
973: /*@C
974: DMCompositeGetGlobalISs - Gets the index sets for each composed object in a `DMCOMPOSITE`
976: Collective
978: Input Parameter:
979: . dm - the `DMCOMPOSITE` object
981: Output Parameter:
982: . is - the array of index sets
984: Level: advanced
986: Notes:
987: The is entries should be destroyed with `ISDestroy()`, the is array should be freed with `PetscFree()`
989: These could be used to extract a subset of vector entries for a "multi-physics" preconditioner
991: Use `DMCompositeGetLocalISs()` for index sets in the packed local numbering, and
992: `DMCompositeGetISLocalToGlobalMappings()` for to map local sub-`DM` (including ghost) indices to packed global
993: indices.
995: Fortran Notes:
996: The output argument 'is' must be an allocated array of sufficient length, which can be learned using `DMCompositeGetNumberDM()`.
998: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
999: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetAccess()`, `DMCompositeScatter()`,
1000: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
1001: @*/
1002: PetscErrorCode DMCompositeGetGlobalISs(DM dm, IS *is[])
1003: {
1004: PetscInt cnt = 0;
1005: struct DMCompositeLink *next;
1006: PetscMPIInt rank;
1007: DM_Composite *com = (DM_Composite *)dm->data;
1008: PetscBool flg;
1010: PetscFunctionBegin;
1012: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1013: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1014: PetscCheck(dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Must call DMSetUp() before");
1015: PetscCall(PetscMalloc1(com->nDM, is));
1016: next = com->next;
1017: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
1019: /* loop over packed objects, handling one at a time */
1020: while (next) {
1021: PetscDS prob;
1023: PetscCall(ISCreateStride(PetscObjectComm((PetscObject)dm), next->n, next->grstart, 1, &(*is)[cnt]));
1024: PetscCall(DMGetDS(dm, &prob));
1025: if (prob) {
1026: MatNullSpace space;
1027: Mat pmat;
1028: PetscObject disc;
1029: PetscInt Nf;
1031: PetscCall(PetscDSGetNumFields(prob, &Nf));
1032: if (cnt < Nf) {
1033: PetscCall(PetscDSGetDiscretization(prob, cnt, &disc));
1034: PetscCall(PetscObjectQuery(disc, "nullspace", (PetscObject *)&space));
1035: if (space) PetscCall(PetscObjectCompose((PetscObject)(*is)[cnt], "nullspace", (PetscObject)space));
1036: PetscCall(PetscObjectQuery(disc, "nearnullspace", (PetscObject *)&space));
1037: if (space) PetscCall(PetscObjectCompose((PetscObject)(*is)[cnt], "nearnullspace", (PetscObject)space));
1038: PetscCall(PetscObjectQuery(disc, "pmat", (PetscObject *)&pmat));
1039: if (pmat) PetscCall(PetscObjectCompose((PetscObject)(*is)[cnt], "pmat", (PetscObject)pmat));
1040: }
1041: }
1042: cnt++;
1043: next = next->next;
1044: }
1045: PetscFunctionReturn(PETSC_SUCCESS);
1046: }
1048: static PetscErrorCode DMCreateFieldIS_Composite(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
1049: {
1050: PetscInt nDM;
1051: DM *dms;
1052: PetscInt i;
1054: PetscFunctionBegin;
1055: PetscCall(DMCompositeGetNumberDM(dm, &nDM));
1056: if (numFields) *numFields = nDM;
1057: PetscCall(DMCompositeGetGlobalISs(dm, fields));
1058: if (fieldNames) {
1059: PetscCall(PetscMalloc1(nDM, &dms));
1060: PetscCall(PetscMalloc1(nDM, fieldNames));
1061: PetscCall(DMCompositeGetEntriesArray(dm, dms));
1062: for (i = 0; i < nDM; i++) {
1063: char buf[256];
1064: const char *splitname;
1066: /* Split naming precedence: object name, prefix, number */
1067: splitname = ((PetscObject)dm)->name;
1068: if (!splitname) {
1069: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)dms[i], &splitname));
1070: if (splitname) {
1071: size_t len;
1072: PetscCall(PetscStrncpy(buf, splitname, sizeof(buf)));
1073: buf[sizeof(buf) - 1] = 0;
1074: PetscCall(PetscStrlen(buf, &len));
1075: if (buf[len - 1] == '_') buf[len - 1] = 0; /* Remove trailing underscore if it was used */
1076: splitname = buf;
1077: }
1078: }
1079: if (!splitname) {
1080: PetscCall(PetscSNPrintf(buf, sizeof(buf), "%" PetscInt_FMT, i));
1081: splitname = buf;
1082: }
1083: PetscCall(PetscStrallocpy(splitname, &(*fieldNames)[i]));
1084: }
1085: PetscCall(PetscFree(dms));
1086: }
1087: PetscFunctionReturn(PETSC_SUCCESS);
1088: }
1090: /*
1091: This could take over from DMCreateFieldIS(), as it is more general,
1092: making DMCreateFieldIS() a special case -- calling with dmlist == NULL;
1093: At this point it's probably best to be less intrusive, however.
1094: */
1095: static PetscErrorCode DMCreateFieldDecomposition_Composite(DM dm, PetscInt *len, char ***namelist, IS **islist, DM **dmlist)
1096: {
1097: PetscInt nDM;
1098: PetscInt i;
1100: PetscFunctionBegin;
1101: PetscCall(DMCreateFieldIS_Composite(dm, len, namelist, islist));
1102: if (dmlist) {
1103: PetscCall(DMCompositeGetNumberDM(dm, &nDM));
1104: PetscCall(PetscMalloc1(nDM, dmlist));
1105: PetscCall(DMCompositeGetEntriesArray(dm, *dmlist));
1106: for (i = 0; i < nDM; i++) PetscCall(PetscObjectReference((PetscObject)((*dmlist)[i])));
1107: }
1108: PetscFunctionReturn(PETSC_SUCCESS);
1109: }
1111: /* -------------------------------------------------------------------------------------*/
1112: /*@C
1113: DMCompositeGetLocalVectors - Gets local vectors for each part of a `DMCOMPOSITE`
1114: Use `DMCompositeRestoreLocalVectors()` to return them.
1116: Not Collective; No Fortran Support
1118: Input Parameter:
1119: . dm - the `DMCOMPOSITE` object
1121: Output Parameter:
1122: . ... - the individual sequential `Vec`s
1124: Level: advanced
1126: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
1127: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
1128: `DMCompositeRestoreLocalVectors()`, `DMCompositeScatter()`, `DMCompositeGetEntries()`
1129: @*/
1130: PetscErrorCode DMCompositeGetLocalVectors(DM dm, ...)
1131: {
1132: va_list Argp;
1133: struct DMCompositeLink *next;
1134: DM_Composite *com = (DM_Composite *)dm->data;
1135: PetscBool flg;
1137: PetscFunctionBegin;
1139: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1140: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1141: next = com->next;
1142: /* loop over packed objects, handling one at a time */
1143: va_start(Argp, dm);
1144: while (next) {
1145: Vec *vec;
1146: vec = va_arg(Argp, Vec *);
1147: if (vec) PetscCall(DMGetLocalVector(next->dm, vec));
1148: next = next->next;
1149: }
1150: va_end(Argp);
1151: PetscFunctionReturn(PETSC_SUCCESS);
1152: }
1154: /*@C
1155: DMCompositeRestoreLocalVectors - Restores local vectors for each part of a `DMCOMPOSITE`
1157: Not Collective; No Fortran Support
1159: Input Parameter:
1160: . dm - the `DMCOMPOSITE` object
1162: Output Parameter:
1163: . ... - the individual sequential `Vec`
1165: Level: advanced
1167: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
1168: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
1169: `DMCompositeGetLocalVectors()`, `DMCompositeScatter()`, `DMCompositeGetEntries()`
1170: @*/
1171: PetscErrorCode DMCompositeRestoreLocalVectors(DM dm, ...)
1172: {
1173: va_list Argp;
1174: struct DMCompositeLink *next;
1175: DM_Composite *com = (DM_Composite *)dm->data;
1176: PetscBool flg;
1178: PetscFunctionBegin;
1180: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1181: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1182: next = com->next;
1183: /* loop over packed objects, handling one at a time */
1184: va_start(Argp, dm);
1185: while (next) {
1186: Vec *vec;
1187: vec = va_arg(Argp, Vec *);
1188: if (vec) PetscCall(DMRestoreLocalVector(next->dm, vec));
1189: next = next->next;
1190: }
1191: va_end(Argp);
1192: PetscFunctionReturn(PETSC_SUCCESS);
1193: }
1195: /* -------------------------------------------------------------------------------------*/
1196: /*@C
1197: DMCompositeGetEntries - Gets the `DM` for each entry in a `DMCOMPOSITE`.
1199: Not Collective
1201: Input Parameter:
1202: . dm - the `DMCOMPOSITE` object
1204: Output Parameter:
1205: . ... - the individual entries `DM`
1207: Level: advanced
1209: Fortran Notes:
1210: Available as `DMCompositeGetEntries()` for one output `DM`, DMCompositeGetEntries2() for 2, etc
1212: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`, `DMCompositeGetEntriesArray()`
1213: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
1214: `DMCompositeRestoreLocalVectors()`, `DMCompositeGetLocalVectors()`, `DMCompositeScatter()`
1215: @*/
1216: PetscErrorCode DMCompositeGetEntries(DM dm, ...)
1217: {
1218: va_list Argp;
1219: struct DMCompositeLink *next;
1220: DM_Composite *com = (DM_Composite *)dm->data;
1221: PetscBool flg;
1223: PetscFunctionBegin;
1225: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1226: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1227: next = com->next;
1228: /* loop over packed objects, handling one at a time */
1229: va_start(Argp, dm);
1230: while (next) {
1231: DM *dmn;
1232: dmn = va_arg(Argp, DM *);
1233: if (dmn) *dmn = next->dm;
1234: next = next->next;
1235: }
1236: va_end(Argp);
1237: PetscFunctionReturn(PETSC_SUCCESS);
1238: }
1240: /*@C
1241: DMCompositeGetEntriesArray - Gets the DM for each entry in a `DMCOMPOSITE`
1243: Not Collective
1245: Input Parameter:
1246: . dm - the `DMCOMPOSITE` object
1248: Output Parameter:
1249: . dms - array of sufficient length (see `DMCompositeGetNumberDM()`) to hold the individual `DM`
1251: Level: advanced
1253: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`, `DMCompositeGetEntries()`
1254: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
1255: `DMCompositeRestoreLocalVectors()`, `DMCompositeGetLocalVectors()`, `DMCompositeScatter()`
1256: @*/
1257: PetscErrorCode DMCompositeGetEntriesArray(DM dm, DM dms[])
1258: {
1259: struct DMCompositeLink *next;
1260: DM_Composite *com = (DM_Composite *)dm->data;
1261: PetscInt i;
1262: PetscBool flg;
1264: PetscFunctionBegin;
1266: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1267: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1268: /* loop over packed objects, handling one at a time */
1269: for (next = com->next, i = 0; next; next = next->next, i++) dms[i] = next->dm;
1270: PetscFunctionReturn(PETSC_SUCCESS);
1271: }
1273: typedef struct {
1274: DM dm;
1275: PetscViewer *subv;
1276: Vec *vecs;
1277: } GLVisViewerCtx;
1279: static PetscErrorCode DestroyGLVisViewerCtx_Private(void *vctx)
1280: {
1281: GLVisViewerCtx *ctx = (GLVisViewerCtx *)vctx;
1282: PetscInt i, n;
1284: PetscFunctionBegin;
1285: PetscCall(DMCompositeGetNumberDM(ctx->dm, &n));
1286: for (i = 0; i < n; i++) PetscCall(PetscViewerDestroy(&ctx->subv[i]));
1287: PetscCall(PetscFree2(ctx->subv, ctx->vecs));
1288: PetscCall(DMDestroy(&ctx->dm));
1289: PetscCall(PetscFree(ctx));
1290: PetscFunctionReturn(PETSC_SUCCESS);
1291: }
1293: static PetscErrorCode DMCompositeSampleGLVisFields_Private(PetscObject oX, PetscInt nf, PetscObject oXfield[], void *vctx)
1294: {
1295: Vec X = (Vec)oX;
1296: GLVisViewerCtx *ctx = (GLVisViewerCtx *)vctx;
1297: PetscInt i, n, cumf;
1299: PetscFunctionBegin;
1300: PetscCall(DMCompositeGetNumberDM(ctx->dm, &n));
1301: PetscCall(DMCompositeGetAccessArray(ctx->dm, X, n, NULL, ctx->vecs));
1302: for (i = 0, cumf = 0; i < n; i++) {
1303: PetscErrorCode (*g2l)(PetscObject, PetscInt, PetscObject[], void *);
1304: void *fctx;
1305: PetscInt nfi;
1307: PetscCall(PetscViewerGLVisGetFields_Internal(ctx->subv[i], &nfi, NULL, NULL, &g2l, NULL, &fctx));
1308: if (!nfi) continue;
1309: if (g2l) PetscCall((*g2l)((PetscObject)ctx->vecs[i], nfi, oXfield + cumf, fctx));
1310: else PetscCall(VecCopy(ctx->vecs[i], (Vec)oXfield[cumf]));
1311: cumf += nfi;
1312: }
1313: PetscCall(DMCompositeRestoreAccessArray(ctx->dm, X, n, NULL, ctx->vecs));
1314: PetscFunctionReturn(PETSC_SUCCESS);
1315: }
1317: static PetscErrorCode DMSetUpGLVisViewer_Composite(PetscObject odm, PetscViewer viewer)
1318: {
1319: DM dm = (DM)odm, *dms;
1320: Vec *Ufds;
1321: GLVisViewerCtx *ctx;
1322: PetscInt i, n, tnf, *sdim;
1323: char **fecs;
1325: PetscFunctionBegin;
1326: PetscCall(PetscNew(&ctx));
1327: PetscCall(PetscObjectReference((PetscObject)dm));
1328: ctx->dm = dm;
1329: PetscCall(DMCompositeGetNumberDM(dm, &n));
1330: PetscCall(PetscMalloc1(n, &dms));
1331: PetscCall(DMCompositeGetEntriesArray(dm, dms));
1332: PetscCall(PetscMalloc2(n, &ctx->subv, n, &ctx->vecs));
1333: for (i = 0, tnf = 0; i < n; i++) {
1334: PetscInt nf;
1336: PetscCall(PetscViewerCreate(PetscObjectComm(odm), &ctx->subv[i]));
1337: PetscCall(PetscViewerSetType(ctx->subv[i], PETSCVIEWERGLVIS));
1338: PetscCall(PetscViewerGLVisSetDM_Internal(ctx->subv[i], (PetscObject)dms[i]));
1339: PetscCall(PetscViewerGLVisGetFields_Internal(ctx->subv[i], &nf, NULL, NULL, NULL, NULL, NULL));
1340: tnf += nf;
1341: }
1342: PetscCall(PetscFree(dms));
1343: PetscCall(PetscMalloc3(tnf, &fecs, tnf, &sdim, tnf, &Ufds));
1344: for (i = 0, tnf = 0; i < n; i++) {
1345: PetscInt *sd, nf, f;
1346: const char **fec;
1347: Vec *Uf;
1349: PetscCall(PetscViewerGLVisGetFields_Internal(ctx->subv[i], &nf, &fec, &sd, NULL, (PetscObject **)&Uf, NULL));
1350: for (f = 0; f < nf; f++) {
1351: PetscCall(PetscStrallocpy(fec[f], &fecs[tnf + f]));
1352: Ufds[tnf + f] = Uf[f];
1353: sdim[tnf + f] = sd[f];
1354: }
1355: tnf += nf;
1356: }
1357: PetscCall(PetscViewerGLVisSetFields(viewer, tnf, (const char **)fecs, sdim, DMCompositeSampleGLVisFields_Private, (PetscObject *)Ufds, ctx, DestroyGLVisViewerCtx_Private));
1358: for (i = 0; i < tnf; i++) PetscCall(PetscFree(fecs[i]));
1359: PetscCall(PetscFree3(fecs, sdim, Ufds));
1360: PetscFunctionReturn(PETSC_SUCCESS);
1361: }
1363: static PetscErrorCode DMRefine_Composite(DM dmi, MPI_Comm comm, DM *fine)
1364: {
1365: struct DMCompositeLink *next;
1366: DM_Composite *com = (DM_Composite *)dmi->data;
1367: DM dm;
1369: PetscFunctionBegin;
1371: if (comm == MPI_COMM_NULL) PetscCall(PetscObjectGetComm((PetscObject)dmi, &comm));
1372: PetscCall(DMSetUp(dmi));
1373: next = com->next;
1374: PetscCall(DMCompositeCreate(comm, fine));
1376: /* loop over packed objects, handling one at a time */
1377: while (next) {
1378: PetscCall(DMRefine(next->dm, comm, &dm));
1379: PetscCall(DMCompositeAddDM(*fine, dm));
1380: PetscCall(PetscObjectDereference((PetscObject)dm));
1381: next = next->next;
1382: }
1383: PetscFunctionReturn(PETSC_SUCCESS);
1384: }
1386: static PetscErrorCode DMCoarsen_Composite(DM dmi, MPI_Comm comm, DM *fine)
1387: {
1388: struct DMCompositeLink *next;
1389: DM_Composite *com = (DM_Composite *)dmi->data;
1390: DM dm;
1392: PetscFunctionBegin;
1394: PetscCall(DMSetUp(dmi));
1395: if (comm == MPI_COMM_NULL) PetscCall(PetscObjectGetComm((PetscObject)dmi, &comm));
1396: next = com->next;
1397: PetscCall(DMCompositeCreate(comm, fine));
1399: /* loop over packed objects, handling one at a time */
1400: while (next) {
1401: PetscCall(DMCoarsen(next->dm, comm, &dm));
1402: PetscCall(DMCompositeAddDM(*fine, dm));
1403: PetscCall(PetscObjectDereference((PetscObject)dm));
1404: next = next->next;
1405: }
1406: PetscFunctionReturn(PETSC_SUCCESS);
1407: }
1409: static PetscErrorCode DMCreateInterpolation_Composite(DM coarse, DM fine, Mat *A, Vec *v)
1410: {
1411: PetscInt m, n, M, N, nDM, i;
1412: struct DMCompositeLink *nextc;
1413: struct DMCompositeLink *nextf;
1414: Vec gcoarse, gfine, *vecs;
1415: DM_Composite *comcoarse = (DM_Composite *)coarse->data;
1416: DM_Composite *comfine = (DM_Composite *)fine->data;
1417: Mat *mats;
1419: PetscFunctionBegin;
1422: PetscCall(DMSetUp(coarse));
1423: PetscCall(DMSetUp(fine));
1424: /* use global vectors only for determining matrix layout */
1425: PetscCall(DMGetGlobalVector(coarse, &gcoarse));
1426: PetscCall(DMGetGlobalVector(fine, &gfine));
1427: PetscCall(VecGetLocalSize(gcoarse, &n));
1428: PetscCall(VecGetLocalSize(gfine, &m));
1429: PetscCall(VecGetSize(gcoarse, &N));
1430: PetscCall(VecGetSize(gfine, &M));
1431: PetscCall(DMRestoreGlobalVector(coarse, &gcoarse));
1432: PetscCall(DMRestoreGlobalVector(fine, &gfine));
1434: nDM = comfine->nDM;
1435: 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);
1436: PetscCall(PetscCalloc1(nDM * nDM, &mats));
1437: if (v) PetscCall(PetscCalloc1(nDM, &vecs));
1439: /* loop over packed objects, handling one at a time */
1440: for (nextc = comcoarse->next, nextf = comfine->next, i = 0; nextc; nextc = nextc->next, nextf = nextf->next, i++) {
1441: if (!v) PetscCall(DMCreateInterpolation(nextc->dm, nextf->dm, &mats[i * nDM + i], NULL));
1442: else PetscCall(DMCreateInterpolation(nextc->dm, nextf->dm, &mats[i * nDM + i], &vecs[i]));
1443: }
1444: PetscCall(MatCreateNest(PetscObjectComm((PetscObject)fine), nDM, NULL, nDM, NULL, mats, A));
1445: if (v) PetscCall(VecCreateNest(PetscObjectComm((PetscObject)fine), nDM, NULL, vecs, v));
1446: for (i = 0; i < nDM * nDM; i++) PetscCall(MatDestroy(&mats[i]));
1447: PetscCall(PetscFree(mats));
1448: if (v) {
1449: for (i = 0; i < nDM; i++) PetscCall(VecDestroy(&vecs[i]));
1450: PetscCall(PetscFree(vecs));
1451: }
1452: PetscFunctionReturn(PETSC_SUCCESS);
1453: }
1455: static PetscErrorCode DMGetLocalToGlobalMapping_Composite(DM dm)
1456: {
1457: DM_Composite *com = (DM_Composite *)dm->data;
1458: ISLocalToGlobalMapping *ltogs;
1459: PetscInt i;
1461: PetscFunctionBegin;
1462: /* Set the ISLocalToGlobalMapping on the new matrix */
1463: PetscCall(DMCompositeGetISLocalToGlobalMappings(dm, <ogs));
1464: PetscCall(ISLocalToGlobalMappingConcatenate(PetscObjectComm((PetscObject)dm), com->nDM, ltogs, &dm->ltogmap));
1465: for (i = 0; i < com->nDM; i++) PetscCall(ISLocalToGlobalMappingDestroy(<ogs[i]));
1466: PetscCall(PetscFree(ltogs));
1467: PetscFunctionReturn(PETSC_SUCCESS);
1468: }
1470: static PetscErrorCode DMCreateColoring_Composite(DM dm, ISColoringType ctype, ISColoring *coloring)
1471: {
1472: PetscInt n, i, cnt;
1473: ISColoringValue *colors;
1474: PetscBool dense = PETSC_FALSE;
1475: ISColoringValue maxcol = 0;
1476: DM_Composite *com = (DM_Composite *)dm->data;
1478: PetscFunctionBegin;
1480: PetscCheck(ctype != IS_COLORING_LOCAL, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Only global coloring supported");
1481: if (ctype == IS_COLORING_GLOBAL) {
1482: n = com->n;
1483: } else SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "Unknown ISColoringType");
1484: PetscCall(PetscMalloc1(n, &colors)); /* freed in ISColoringDestroy() */
1486: PetscCall(PetscOptionsGetBool(((PetscObject)dm)->options, ((PetscObject)dm)->prefix, "-dmcomposite_dense_jacobian", &dense, NULL));
1487: if (dense) {
1488: for (i = 0; i < n; i++) colors[i] = (ISColoringValue)(com->rstart + i);
1489: maxcol = com->N;
1490: } else {
1491: struct DMCompositeLink *next = com->next;
1492: PetscMPIInt rank;
1494: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
1495: cnt = 0;
1496: while (next) {
1497: ISColoring lcoloring;
1499: PetscCall(DMCreateColoring(next->dm, IS_COLORING_GLOBAL, &lcoloring));
1500: for (i = 0; i < lcoloring->N; i++) colors[cnt++] = maxcol + lcoloring->colors[i];
1501: maxcol += lcoloring->n;
1502: PetscCall(ISColoringDestroy(&lcoloring));
1503: next = next->next;
1504: }
1505: }
1506: PetscCall(ISColoringCreate(PetscObjectComm((PetscObject)dm), maxcol, n, colors, PETSC_OWN_POINTER, coloring));
1507: PetscFunctionReturn(PETSC_SUCCESS);
1508: }
1510: static PetscErrorCode DMGlobalToLocalBegin_Composite(DM dm, Vec gvec, InsertMode mode, Vec lvec)
1511: {
1512: struct DMCompositeLink *next;
1513: PetscScalar *garray, *larray;
1514: DM_Composite *com = (DM_Composite *)dm->data;
1516: PetscFunctionBegin;
1520: if (!com->setup) PetscCall(DMSetUp(dm));
1522: PetscCall(VecGetArray(gvec, &garray));
1523: PetscCall(VecGetArray(lvec, &larray));
1525: /* loop over packed objects, handling one at a time */
1526: next = com->next;
1527: while (next) {
1528: Vec local, global;
1529: PetscInt N;
1531: PetscCall(DMGetGlobalVector(next->dm, &global));
1532: PetscCall(VecGetLocalSize(global, &N));
1533: PetscCall(VecPlaceArray(global, garray));
1534: PetscCall(DMGetLocalVector(next->dm, &local));
1535: PetscCall(VecPlaceArray(local, larray));
1536: PetscCall(DMGlobalToLocalBegin(next->dm, global, mode, local));
1537: PetscCall(DMGlobalToLocalEnd(next->dm, global, mode, local));
1538: PetscCall(VecResetArray(global));
1539: PetscCall(VecResetArray(local));
1540: PetscCall(DMRestoreGlobalVector(next->dm, &global));
1541: PetscCall(DMRestoreLocalVector(next->dm, &local));
1543: larray += next->nlocal;
1544: garray += next->n;
1545: next = next->next;
1546: }
1548: PetscCall(VecRestoreArray(gvec, NULL));
1549: PetscCall(VecRestoreArray(lvec, NULL));
1550: PetscFunctionReturn(PETSC_SUCCESS);
1551: }
1553: static PetscErrorCode DMGlobalToLocalEnd_Composite(DM dm, Vec gvec, InsertMode mode, Vec lvec)
1554: {
1555: PetscFunctionBegin;
1559: PetscFunctionReturn(PETSC_SUCCESS);
1560: }
1562: static PetscErrorCode DMLocalToGlobalBegin_Composite(DM dm, Vec lvec, InsertMode mode, Vec gvec)
1563: {
1564: struct DMCompositeLink *next;
1565: PetscScalar *larray, *garray;
1566: DM_Composite *com = (DM_Composite *)dm->data;
1568: PetscFunctionBegin;
1573: if (!com->setup) PetscCall(DMSetUp(dm));
1575: PetscCall(VecGetArray(lvec, &larray));
1576: PetscCall(VecGetArray(gvec, &garray));
1578: /* loop over packed objects, handling one at a time */
1579: next = com->next;
1580: while (next) {
1581: Vec global, local;
1583: PetscCall(DMGetLocalVector(next->dm, &local));
1584: PetscCall(VecPlaceArray(local, larray));
1585: PetscCall(DMGetGlobalVector(next->dm, &global));
1586: PetscCall(VecPlaceArray(global, garray));
1587: PetscCall(DMLocalToGlobalBegin(next->dm, local, mode, global));
1588: PetscCall(DMLocalToGlobalEnd(next->dm, local, mode, global));
1589: PetscCall(VecResetArray(local));
1590: PetscCall(VecResetArray(global));
1591: PetscCall(DMRestoreGlobalVector(next->dm, &global));
1592: PetscCall(DMRestoreLocalVector(next->dm, &local));
1594: garray += next->n;
1595: larray += next->nlocal;
1596: next = next->next;
1597: }
1599: PetscCall(VecRestoreArray(gvec, NULL));
1600: PetscCall(VecRestoreArray(lvec, NULL));
1601: PetscFunctionReturn(PETSC_SUCCESS);
1602: }
1604: static PetscErrorCode DMLocalToGlobalEnd_Composite(DM dm, Vec lvec, InsertMode mode, Vec gvec)
1605: {
1606: PetscFunctionBegin;
1610: PetscFunctionReturn(PETSC_SUCCESS);
1611: }
1613: static PetscErrorCode DMLocalToLocalBegin_Composite(DM dm, Vec vec1, InsertMode mode, Vec vec2)
1614: {
1615: struct DMCompositeLink *next;
1616: PetscScalar *array1, *array2;
1617: DM_Composite *com = (DM_Composite *)dm->data;
1619: PetscFunctionBegin;
1624: if (!com->setup) PetscCall(DMSetUp(dm));
1626: PetscCall(VecGetArray(vec1, &array1));
1627: PetscCall(VecGetArray(vec2, &array2));
1629: /* loop over packed objects, handling one at a time */
1630: next = com->next;
1631: while (next) {
1632: Vec local1, local2;
1634: PetscCall(DMGetLocalVector(next->dm, &local1));
1635: PetscCall(VecPlaceArray(local1, array1));
1636: PetscCall(DMGetLocalVector(next->dm, &local2));
1637: PetscCall(VecPlaceArray(local2, array2));
1638: PetscCall(DMLocalToLocalBegin(next->dm, local1, mode, local2));
1639: PetscCall(DMLocalToLocalEnd(next->dm, local1, mode, local2));
1640: PetscCall(VecResetArray(local2));
1641: PetscCall(DMRestoreLocalVector(next->dm, &local2));
1642: PetscCall(VecResetArray(local1));
1643: PetscCall(DMRestoreLocalVector(next->dm, &local1));
1645: array1 += next->nlocal;
1646: array2 += next->nlocal;
1647: next = next->next;
1648: }
1650: PetscCall(VecRestoreArray(vec1, NULL));
1651: PetscCall(VecRestoreArray(vec2, NULL));
1652: PetscFunctionReturn(PETSC_SUCCESS);
1653: }
1655: static PetscErrorCode DMLocalToLocalEnd_Composite(DM dm, Vec lvec, InsertMode mode, Vec gvec)
1656: {
1657: PetscFunctionBegin;
1661: PetscFunctionReturn(PETSC_SUCCESS);
1662: }
1664: /*MC
1665: DMCOMPOSITE = "composite" - A `DM` object that is used to manage data for a collection of `DM`
1667: Level: intermediate
1669: .seealso: `DMType`, `DM`, `DMDACreate()`, `DMCreate()`, `DMSetType()`, `DMCompositeCreate()`
1670: M*/
1672: PETSC_EXTERN PetscErrorCode DMCreate_Composite(DM p)
1673: {
1674: DM_Composite *com;
1676: PetscFunctionBegin;
1677: PetscCall(PetscNew(&com));
1678: p->data = com;
1679: com->n = 0;
1680: com->nghost = 0;
1681: com->next = NULL;
1682: com->nDM = 0;
1684: p->ops->createglobalvector = DMCreateGlobalVector_Composite;
1685: p->ops->createlocalvector = DMCreateLocalVector_Composite;
1686: p->ops->getlocaltoglobalmapping = DMGetLocalToGlobalMapping_Composite;
1687: p->ops->createfieldis = DMCreateFieldIS_Composite;
1688: p->ops->createfielddecomposition = DMCreateFieldDecomposition_Composite;
1689: p->ops->refine = DMRefine_Composite;
1690: p->ops->coarsen = DMCoarsen_Composite;
1691: p->ops->createinterpolation = DMCreateInterpolation_Composite;
1692: p->ops->creatematrix = DMCreateMatrix_Composite;
1693: p->ops->getcoloring = DMCreateColoring_Composite;
1694: p->ops->globaltolocalbegin = DMGlobalToLocalBegin_Composite;
1695: p->ops->globaltolocalend = DMGlobalToLocalEnd_Composite;
1696: p->ops->localtoglobalbegin = DMLocalToGlobalBegin_Composite;
1697: p->ops->localtoglobalend = DMLocalToGlobalEnd_Composite;
1698: p->ops->localtolocalbegin = DMLocalToLocalBegin_Composite;
1699: p->ops->localtolocalend = DMLocalToLocalEnd_Composite;
1700: p->ops->destroy = DMDestroy_Composite;
1701: p->ops->view = DMView_Composite;
1702: p->ops->setup = DMSetUp_Composite;
1704: PetscCall(PetscObjectComposeFunction((PetscObject)p, "DMSetUpGLVisViewer_C", DMSetUpGLVisViewer_Composite));
1705: PetscFunctionReturn(PETSC_SUCCESS);
1706: }
1708: /*@
1709: DMCompositeCreate - Creates a `DMCOMPOSITE`, used to generate "composite"
1710: vectors made up of several subvectors.
1712: Collective
1714: Input Parameter:
1715: . comm - the processors that will share the global vector
1717: Output Parameter:
1718: . packer - the `DMCOMPOSITE` object
1720: Level: advanced
1722: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCompositeScatter()`, `DMCreate()`
1723: `DMCompositeGather()`, `DMCreateGlobalVector()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`
1724: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
1725: @*/
1726: PetscErrorCode DMCompositeCreate(MPI_Comm comm, DM *packer)
1727: {
1728: PetscFunctionBegin;
1729: PetscAssertPointer(packer, 2);
1730: PetscCall(DMCreate(comm, packer));
1731: PetscCall(DMSetType(*packer, DMCOMPOSITE));
1732: PetscFunctionReturn(PETSC_SUCCESS);
1733: }