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: /*@
117: DMCompositeGetNumberDM - Gets the number of `DM` objects in the `DMCOMPOSITE`
118: representation.
120: Not Collective
122: Input Parameter:
123: . dm - the `DMCOMPOSITE` object
125: Output Parameter:
126: . nDM - the number of `DM`
128: Level: beginner
130: .seealso: `DMCOMPOSITE`, `DM`
131: @*/
132: PetscErrorCode DMCompositeGetNumberDM(DM dm, PetscInt *nDM)
133: {
134: DM_Composite *com = (DM_Composite *)dm->data;
135: PetscBool flg;
137: PetscFunctionBegin;
139: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
140: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
141: *nDM = com->nDM;
142: PetscFunctionReturn(PETSC_SUCCESS);
143: }
145: /*@C
146: DMCompositeGetAccess - Allows one to access the individual packed vectors in their global
147: representation.
149: Collective
151: Input Parameters:
152: + dm - the `DMCOMPOSITE` object
153: - gvec - the global vector
155: Output Parameter:
156: . ... - the packed parallel vectors, `NULL` for those that are not needed
158: Level: advanced
160: Note:
161: Use `DMCompositeRestoreAccess()` to return the vectors when you no longer need them
163: Fortran Notes:
164: Fortran callers must use numbered versions of this routine, e.g., DMCompositeGetAccess4(dm,gvec,vec1,vec2,vec3,vec4)
165: or use the alternative interface `DMCompositeGetAccessArray()`.
167: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeGetEntries()`, `DMCompositeScatter()`
168: @*/
169: PetscErrorCode DMCompositeGetAccess(DM dm, Vec gvec, ...)
170: {
171: va_list Argp;
172: struct DMCompositeLink *next;
173: DM_Composite *com = (DM_Composite *)dm->data;
174: PetscInt readonly;
175: PetscBool flg;
177: PetscFunctionBegin;
180: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
181: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
182: next = com->next;
183: if (!com->setup) PetscCall(DMSetUp(dm));
185: PetscCall(VecLockGet(gvec, &readonly));
186: /* loop over packed objects, handling one at a time */
187: va_start(Argp, gvec);
188: while (next) {
189: Vec *vec;
190: vec = va_arg(Argp, Vec *);
191: if (vec) {
192: PetscCall(DMGetGlobalVector(next->dm, vec));
193: if (readonly) {
194: const PetscScalar *array;
195: PetscCall(VecGetArrayRead(gvec, &array));
196: PetscCall(VecPlaceArray(*vec, array + next->rstart));
197: PetscCall(VecLockReadPush(*vec));
198: PetscCall(VecRestoreArrayRead(gvec, &array));
199: } else {
200: PetscScalar *array;
201: PetscCall(VecGetArray(gvec, &array));
202: PetscCall(VecPlaceArray(*vec, array + next->rstart));
203: PetscCall(VecRestoreArray(gvec, &array));
204: }
205: }
206: next = next->next;
207: }
208: va_end(Argp);
209: PetscFunctionReturn(PETSC_SUCCESS);
210: }
212: /*@
213: DMCompositeGetAccessArray - Allows one to access the individual packed vectors in their global
214: representation.
216: Collective
218: Input Parameters:
219: + dm - the `DMCOMPOSITE`
220: . pvec - packed vector
221: . nwanted - number of vectors wanted
222: - wanted - sorted array of vectors wanted, or `NULL` to get all vectors, length `nwanted`
224: Output Parameter:
225: . vecs - array of requested global vectors (must be previously allocated and of length `nwanted`)
227: Level: advanced
229: Note:
230: Use `DMCompositeRestoreAccessArray()` to return the vectors when you no longer need them
232: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeGetAccess()`, `DMCompositeGetEntries()`, `DMCompositeScatter()`, `DMCompositeGather()`
233: @*/
234: PetscErrorCode DMCompositeGetAccessArray(DM dm, Vec pvec, PetscInt nwanted, const PetscInt wanted[], Vec vecs[])
235: {
236: struct DMCompositeLink *link;
237: PetscInt i, wnum;
238: DM_Composite *com = (DM_Composite *)dm->data;
239: PetscInt readonly;
240: PetscBool flg;
242: PetscFunctionBegin;
245: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
246: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
247: if (!com->setup) PetscCall(DMSetUp(dm));
249: PetscCall(VecLockGet(pvec, &readonly));
250: for (i = 0, wnum = 0, link = com->next; link && wnum < nwanted; i++, link = link->next) {
251: if (!wanted || i == wanted[wnum]) {
252: Vec v;
253: PetscCall(DMGetGlobalVector(link->dm, &v));
254: if (readonly) {
255: const PetscScalar *array;
256: PetscCall(VecGetArrayRead(pvec, &array));
257: PetscCall(VecPlaceArray(v, array + link->rstart));
258: PetscCall(VecLockReadPush(v));
259: PetscCall(VecRestoreArrayRead(pvec, &array));
260: } else {
261: PetscScalar *array;
262: PetscCall(VecGetArray(pvec, &array));
263: PetscCall(VecPlaceArray(v, array + link->rstart));
264: PetscCall(VecRestoreArray(pvec, &array));
265: }
266: vecs[wnum++] = v;
267: }
268: }
269: PetscFunctionReturn(PETSC_SUCCESS);
270: }
272: /*@
273: DMCompositeGetLocalAccessArray - Allows one to access the individual
274: packed vectors in their local representation.
276: Collective
278: Input Parameters:
279: + dm - the `DMCOMPOSITE`
280: . pvec - packed vector
281: . nwanted - number of vectors wanted
282: - wanted - sorted array of vectors wanted, or `NULL` to get all vectors, length `nwanted`
284: Output Parameter:
285: . vecs - array of requested local vectors (must be allocated and of length `nwanted`)
287: Level: advanced
289: Note:
290: Use `DMCompositeRestoreLocalAccessArray()` to return the vectors
291: when you no longer need them.
293: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeRestoreLocalAccessArray()`, `DMCompositeGetAccess()`,
294: `DMCompositeGetEntries()`, `DMCompositeScatter()`, `DMCompositeGather()`
295: @*/
296: PetscErrorCode DMCompositeGetLocalAccessArray(DM dm, Vec pvec, PetscInt nwanted, const PetscInt wanted[], Vec vecs[])
297: {
298: struct DMCompositeLink *link;
299: PetscInt i, wnum;
300: DM_Composite *com = (DM_Composite *)dm->data;
301: PetscInt readonly;
302: PetscInt nlocal = 0;
303: PetscBool flg;
305: PetscFunctionBegin;
308: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
309: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
310: if (!com->setup) PetscCall(DMSetUp(dm));
312: PetscCall(VecLockGet(pvec, &readonly));
313: for (i = 0, wnum = 0, link = com->next; link && wnum < nwanted; i++, link = link->next) {
314: if (!wanted || i == wanted[wnum]) {
315: Vec v;
316: PetscCall(DMGetLocalVector(link->dm, &v));
317: if (readonly) {
318: const PetscScalar *array;
319: PetscCall(VecGetArrayRead(pvec, &array));
320: PetscCall(VecPlaceArray(v, array + nlocal));
321: // 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
322: PetscCall(VecLockReadPush(v));
323: PetscCall(VecRestoreArrayRead(pvec, &array));
324: } else {
325: PetscScalar *array;
326: PetscCall(VecGetArray(pvec, &array));
327: PetscCall(VecPlaceArray(v, array + nlocal));
328: PetscCall(VecRestoreArray(pvec, &array));
329: }
330: vecs[wnum++] = v;
331: }
333: nlocal += link->nlocal;
334: }
335: PetscFunctionReturn(PETSC_SUCCESS);
336: }
338: /*@C
339: DMCompositeRestoreAccess - Returns the vectors obtained with `DMCompositeGetAccess()`
340: representation.
342: Collective
344: Input Parameters:
345: + dm - the `DMCOMPOSITE` object
346: . gvec - the global vector
347: - ... - the individual parallel vectors, `NULL` for those that are not needed
349: Level: advanced
351: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
352: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeScatter()`,
353: `DMCompositeGetAccess()`
354: @*/
355: PetscErrorCode DMCompositeRestoreAccess(DM dm, Vec gvec, ...)
356: {
357: va_list Argp;
358: struct DMCompositeLink *next;
359: DM_Composite *com = (DM_Composite *)dm->data;
360: PetscInt readonly;
361: PetscBool flg;
363: PetscFunctionBegin;
366: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
367: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
368: next = com->next;
369: if (!com->setup) PetscCall(DMSetUp(dm));
371: PetscCall(VecLockGet(gvec, &readonly));
372: /* loop over packed objects, handling one at a time */
373: va_start(Argp, gvec);
374: while (next) {
375: Vec *vec;
376: vec = va_arg(Argp, Vec *);
377: if (vec) {
378: PetscCall(VecResetArray(*vec));
379: if (readonly) PetscCall(VecLockReadPop(*vec));
380: PetscCall(DMRestoreGlobalVector(next->dm, vec));
381: }
382: next = next->next;
383: }
384: va_end(Argp);
385: PetscFunctionReturn(PETSC_SUCCESS);
386: }
388: /*@
389: DMCompositeRestoreAccessArray - Returns the vectors obtained with `DMCompositeGetAccessArray()`
391: Collective
393: Input Parameters:
394: + dm - the `DMCOMPOSITE` object
395: . pvec - packed vector
396: . nwanted - number of vectors wanted
397: . wanted - sorted array of vectors wanted, or `NULL` to restore all vectors
398: - vecs - array of global vectors
400: Level: advanced
402: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeRestoreAccess()`, `DMCompositeRestoreEntries()`, `DMCompositeScatter()`, `DMCompositeGather()`
403: @*/
404: PetscErrorCode DMCompositeRestoreAccessArray(DM dm, Vec pvec, PetscInt nwanted, const PetscInt wanted[], Vec vecs[])
405: {
406: struct DMCompositeLink *link;
407: PetscInt i, wnum;
408: DM_Composite *com = (DM_Composite *)dm->data;
409: PetscInt readonly;
410: PetscBool flg;
412: PetscFunctionBegin;
415: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
416: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
417: if (!com->setup) PetscCall(DMSetUp(dm));
419: PetscCall(VecLockGet(pvec, &readonly));
420: for (i = 0, wnum = 0, link = com->next; link && wnum < nwanted; i++, link = link->next) {
421: if (!wanted || i == wanted[wnum]) {
422: PetscCall(VecResetArray(vecs[wnum]));
423: if (readonly) PetscCall(VecLockReadPop(vecs[wnum]));
424: PetscCall(DMRestoreGlobalVector(link->dm, &vecs[wnum]));
425: wnum++;
426: }
427: }
428: PetscFunctionReturn(PETSC_SUCCESS);
429: }
431: /*@
432: DMCompositeRestoreLocalAccessArray - Returns the vectors obtained with `DMCompositeGetLocalAccessArray()`.
434: Collective
436: Input Parameters:
437: + dm - the `DMCOMPOSITE` object
438: . pvec - packed vector
439: . nwanted - number of vectors wanted
440: . wanted - sorted array of vectors wanted, or `NULL` to restore all vectors
441: - vecs - array of local vectors
443: Level: advanced
445: Note:
446: `nwanted` and `wanted` must match the values given to `DMCompositeGetLocalAccessArray()`
447: otherwise the call will fail.
449: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeGetLocalAccessArray()`, `DMCompositeRestoreAccessArray()`,
450: `DMCompositeRestoreAccess()`, `DMCompositeRestoreEntries()`,
451: `DMCompositeScatter()`, `DMCompositeGather()`
452: @*/
453: PetscErrorCode DMCompositeRestoreLocalAccessArray(DM dm, Vec pvec, PetscInt nwanted, const PetscInt wanted[], Vec *vecs)
454: {
455: struct DMCompositeLink *link;
456: PetscInt i, wnum;
457: DM_Composite *com = (DM_Composite *)dm->data;
458: PetscInt readonly;
459: PetscBool flg;
461: PetscFunctionBegin;
464: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
465: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
466: if (!com->setup) PetscCall(DMSetUp(dm));
468: PetscCall(VecLockGet(pvec, &readonly));
469: for (i = 0, wnum = 0, link = com->next; link && wnum < nwanted; i++, link = link->next) {
470: if (!wanted || i == wanted[wnum]) {
471: PetscCall(VecResetArray(vecs[wnum]));
472: if (readonly) PetscCall(VecLockReadPop(vecs[wnum]));
473: PetscCall(DMRestoreLocalVector(link->dm, &vecs[wnum]));
474: wnum++;
475: }
476: }
477: PetscFunctionReturn(PETSC_SUCCESS);
478: }
480: /*@C
481: DMCompositeScatter - Scatters from a global packed vector into its individual local vectors
483: Collective
485: Input Parameters:
486: + dm - the `DMCOMPOSITE` object
487: . gvec - the global vector
488: - ... - the individual sequential vectors, `NULL` for those that are not needed
490: Level: advanced
492: Note:
493: `DMCompositeScatterArray()` is a non-variadic alternative that is often more convenient for library callers and is
494: accessible from Fortran.
496: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
497: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
498: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
499: `DMCompositeScatterArray()`
500: @*/
501: PetscErrorCode DMCompositeScatter(DM dm, Vec gvec, ...)
502: {
503: va_list Argp;
504: struct DMCompositeLink *next;
505: PETSC_UNUSED PetscInt cnt;
506: DM_Composite *com = (DM_Composite *)dm->data;
507: PetscBool flg;
509: PetscFunctionBegin;
512: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
513: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
514: if (!com->setup) PetscCall(DMSetUp(dm));
516: /* loop over packed objects, handling one at a time */
517: va_start(Argp, gvec);
518: for (cnt = 3, next = com->next; next; cnt++, next = next->next) {
519: Vec local;
520: local = va_arg(Argp, Vec);
521: if (local) {
522: Vec global;
523: const PetscScalar *array;
525: PetscCall(DMGetGlobalVector(next->dm, &global));
526: PetscCall(VecGetArrayRead(gvec, &array));
527: PetscCall(VecPlaceArray(global, array + next->rstart));
528: PetscCall(DMGlobalToLocalBegin(next->dm, global, INSERT_VALUES, local));
529: PetscCall(DMGlobalToLocalEnd(next->dm, global, INSERT_VALUES, local));
530: PetscCall(VecRestoreArrayRead(gvec, &array));
531: PetscCall(VecResetArray(global));
532: PetscCall(DMRestoreGlobalVector(next->dm, &global));
533: }
534: }
535: va_end(Argp);
536: PetscFunctionReturn(PETSC_SUCCESS);
537: }
539: /*@
540: DMCompositeScatterArray - Scatters from a global packed vector into its individual local vectors
542: Collective
544: Input Parameters:
545: + dm - the `DMCOMPOSITE` object
546: . gvec - the global vector
547: - lvecs - array of local vectors, NULL for any that are not needed
549: Level: advanced
551: Note:
552: This is a non-variadic alternative to `DMCompositeScatter()`
554: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`
555: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
556: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
557: @*/
558: PetscErrorCode DMCompositeScatterArray(DM dm, Vec gvec, Vec *lvecs)
559: {
560: struct DMCompositeLink *next;
561: PetscInt i;
562: DM_Composite *com = (DM_Composite *)dm->data;
563: PetscBool flg;
565: PetscFunctionBegin;
568: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
569: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
570: if (!com->setup) PetscCall(DMSetUp(dm));
572: /* loop over packed objects, handling one at a time */
573: for (i = 0, next = com->next; next; next = next->next, i++) {
574: if (lvecs[i]) {
575: Vec global;
576: const PetscScalar *array;
578: PetscCall(DMGetGlobalVector(next->dm, &global));
579: PetscCall(VecGetArrayRead(gvec, &array));
580: PetscCall(VecPlaceArray(global, (PetscScalar *)array + next->rstart));
581: PetscCall(DMGlobalToLocalBegin(next->dm, global, INSERT_VALUES, lvecs[i]));
582: PetscCall(DMGlobalToLocalEnd(next->dm, global, INSERT_VALUES, lvecs[i]));
583: PetscCall(VecRestoreArrayRead(gvec, &array));
584: PetscCall(VecResetArray(global));
585: PetscCall(DMRestoreGlobalVector(next->dm, &global));
586: }
587: }
588: PetscFunctionReturn(PETSC_SUCCESS);
589: }
591: /*@C
592: DMCompositeGather - Gathers into a global packed vector from its individual local vectors
594: Collective
596: Input Parameters:
597: + dm - the `DMCOMPOSITE` object
598: . imode - `INSERT_VALUES` or `ADD_VALUES`
599: . gvec - the global vector
600: - ... - the individual sequential vectors, `NULL` for any that are not needed
602: Level: advanced
604: Fortran Notes:
605: Fortran users should use `DMCompositeGatherArray()`
607: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
608: `DMCompositeScatter()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
609: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
610: @*/
611: PetscErrorCode DMCompositeGather(DM dm, InsertMode imode, Vec gvec, ...)
612: {
613: va_list Argp;
614: struct DMCompositeLink *next;
615: DM_Composite *com = (DM_Composite *)dm->data;
616: PETSC_UNUSED PetscInt cnt;
617: PetscBool flg;
619: PetscFunctionBegin;
622: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
623: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
624: if (!com->setup) PetscCall(DMSetUp(dm));
626: /* loop over packed objects, handling one at a time */
627: va_start(Argp, gvec);
628: for (cnt = 3, next = com->next; next; cnt++, next = next->next) {
629: Vec local;
630: local = va_arg(Argp, Vec);
631: if (local) {
632: PetscScalar *array;
633: Vec global;
635: PetscCall(DMGetGlobalVector(next->dm, &global));
636: PetscCall(VecGetArray(gvec, &array));
637: PetscCall(VecPlaceArray(global, array + next->rstart));
638: PetscCall(DMLocalToGlobalBegin(next->dm, local, imode, global));
639: PetscCall(DMLocalToGlobalEnd(next->dm, local, imode, global));
640: PetscCall(VecRestoreArray(gvec, &array));
641: PetscCall(VecResetArray(global));
642: PetscCall(DMRestoreGlobalVector(next->dm, &global));
643: }
644: }
645: va_end(Argp);
646: PetscFunctionReturn(PETSC_SUCCESS);
647: }
649: /*@
650: DMCompositeGatherArray - Gathers into a global packed vector from its individual local vectors
652: Collective
654: Input Parameters:
655: + dm - the `DMCOMPOSITE` object
656: . gvec - the global vector
657: . imode - `INSERT_VALUES` or `ADD_VALUES`
658: - lvecs - the individual sequential vectors, NULL for any that are not needed
660: Level: advanced
662: Note:
663: This is a non-variadic alternative to `DMCompositeGather()`.
665: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
666: `DMCompositeScatter()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
667: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`,
668: @*/
669: PetscErrorCode DMCompositeGatherArray(DM dm, InsertMode imode, Vec gvec, Vec *lvecs)
670: {
671: struct DMCompositeLink *next;
672: DM_Composite *com = (DM_Composite *)dm->data;
673: PetscInt i;
674: PetscBool flg;
676: PetscFunctionBegin;
679: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
680: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
681: if (!com->setup) PetscCall(DMSetUp(dm));
683: /* loop over packed objects, handling one at a time */
684: for (next = com->next, i = 0; next; next = next->next, i++) {
685: if (lvecs[i]) {
686: PetscScalar *array;
687: Vec global;
689: PetscCall(DMGetGlobalVector(next->dm, &global));
690: PetscCall(VecGetArray(gvec, &array));
691: PetscCall(VecPlaceArray(global, array + next->rstart));
692: PetscCall(DMLocalToGlobalBegin(next->dm, lvecs[i], imode, global));
693: PetscCall(DMLocalToGlobalEnd(next->dm, lvecs[i], imode, global));
694: PetscCall(VecRestoreArray(gvec, &array));
695: PetscCall(VecResetArray(global));
696: PetscCall(DMRestoreGlobalVector(next->dm, &global));
697: }
698: }
699: PetscFunctionReturn(PETSC_SUCCESS);
700: }
702: /*@
703: DMCompositeAddDM - adds a `DM` vector to a `DMCOMPOSITE`
705: Collective
707: Input Parameters:
708: + dmc - the `DMCOMPOSITE` object
709: - dm - the `DM` object
711: Level: advanced
713: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeGather()`, `DMCreateGlobalVector()`,
714: `DMCompositeScatter()`, `DMCompositeCreate()`, `DMCompositeGetISLocalToGlobalMappings()`, `DMCompositeGetAccess()`,
715: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
716: @*/
717: PetscErrorCode DMCompositeAddDM(DM dmc, DM dm)
718: {
719: PetscInt n, nlocal;
720: struct DMCompositeLink *mine, *next;
721: Vec global, local;
722: DM_Composite *com = (DM_Composite *)dmc->data;
723: PetscBool flg;
725: PetscFunctionBegin;
728: PetscCall(PetscObjectTypeCompare((PetscObject)dmc, DMCOMPOSITE, &flg));
729: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
730: next = com->next;
731: PetscCheck(!com->setup, PetscObjectComm((PetscObject)dmc), PETSC_ERR_ARG_WRONGSTATE, "Cannot add a DM once you have used the DMComposite");
733: /* create new link */
734: PetscCall(PetscNew(&mine));
735: PetscCall(PetscObjectReference((PetscObject)dm));
736: PetscCall(DMCreateGlobalVector(dm, &global)); // Not using DMGetGlobalVector() since it will stash the vector with a type decided by dm,
737: PetscCall(VecGetLocalSize(global, &n)); // while we may want to set dmc's vectype later, say via DMSetFromOptions(dmc), and we
738: PetscCall(VecDestroy(&global)); // want to propagate the type to dm.
739: PetscCall(DMCreateLocalVector(dm, &local)); // Not using DMGetLocalVector(), same reason as above.
740: PetscCall(VecGetSize(local, &nlocal));
741: PetscCall(VecDestroy(&local));
743: mine->n = n;
744: mine->nlocal = nlocal;
745: mine->dm = dm;
746: mine->next = NULL;
747: com->n += n;
748: com->nghost += nlocal;
750: /* add to end of list */
751: if (!next) com->next = mine;
752: else {
753: while (next->next) next = next->next;
754: next->next = mine;
755: }
756: com->nDM++;
757: com->nmine++;
758: PetscFunctionReturn(PETSC_SUCCESS);
759: }
761: #include <petscdraw.h>
762: PETSC_EXTERN PetscErrorCode VecView_MPI(Vec, PetscViewer);
763: static PetscErrorCode VecView_DMComposite(Vec gvec, PetscViewer viewer)
764: {
765: DM dm;
766: struct DMCompositeLink *next;
767: PetscBool isdraw;
768: DM_Composite *com;
770: PetscFunctionBegin;
771: PetscCall(VecGetDM(gvec, &dm));
772: PetscCheck(dm, PetscObjectComm((PetscObject)gvec), PETSC_ERR_ARG_WRONG, "Vector not generated from a DMComposite");
773: com = (DM_Composite *)dm->data;
774: next = com->next;
776: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERDRAW, &isdraw));
777: if (!isdraw) {
778: /* do I really want to call this? */
779: PetscCall(VecView_MPI(gvec, viewer));
780: } else {
781: PetscInt cnt = 0;
783: /* loop over packed objects, handling one at a time */
784: while (next) {
785: Vec vec;
786: const PetscScalar *array;
787: PetscInt bs;
789: /* Should use VecGetSubVector() eventually, but would need to forward the DM for that to work */
790: PetscCall(DMGetGlobalVector(next->dm, &vec));
791: PetscCall(VecGetArrayRead(gvec, &array));
792: PetscCall(VecPlaceArray(vec, (PetscScalar *)array + next->rstart));
793: PetscCall(VecRestoreArrayRead(gvec, &array));
794: PetscCall(VecView(vec, viewer));
795: PetscCall(VecResetArray(vec));
796: PetscCall(VecGetBlockSize(vec, &bs));
797: PetscCall(DMRestoreGlobalVector(next->dm, &vec));
798: PetscCall(PetscViewerDrawBaseAdd(viewer, bs));
799: cnt += bs;
800: next = next->next;
801: }
802: PetscCall(PetscViewerDrawBaseAdd(viewer, -cnt));
803: }
804: PetscFunctionReturn(PETSC_SUCCESS);
805: }
807: static PetscErrorCode DMCreateGlobalVector_Composite(DM dm, Vec *gvec)
808: {
809: DM_Composite *com = (DM_Composite *)dm->data;
811: PetscFunctionBegin;
813: PetscCall(DMSetFromOptions(dm));
814: PetscCall(DMSetUp(dm));
815: PetscCall(VecCreate(PetscObjectComm((PetscObject)dm), gvec));
816: PetscCall(VecSetType(*gvec, dm->vectype));
817: PetscCall(VecSetSizes(*gvec, com->n, com->N));
818: PetscCall(VecSetDM(*gvec, dm));
819: PetscCall(VecSetOperation(*gvec, VECOP_VIEW, (void (*)(void))VecView_DMComposite));
820: PetscFunctionReturn(PETSC_SUCCESS);
821: }
823: static PetscErrorCode DMCreateLocalVector_Composite(DM dm, Vec *lvec)
824: {
825: DM_Composite *com = (DM_Composite *)dm->data;
827: PetscFunctionBegin;
829: if (!com->setup) {
830: PetscCall(DMSetFromOptions(dm));
831: PetscCall(DMSetUp(dm));
832: }
833: PetscCall(VecCreate(PETSC_COMM_SELF, lvec));
834: PetscCall(VecSetType(*lvec, dm->vectype));
835: PetscCall(VecSetSizes(*lvec, com->nghost, PETSC_DECIDE));
836: PetscCall(VecSetDM(*lvec, dm));
837: PetscFunctionReturn(PETSC_SUCCESS);
838: }
840: /*@C
841: DMCompositeGetISLocalToGlobalMappings - gets an `ISLocalToGlobalMapping` for each `DM` in the `DMCOMPOSITE`, maps to the composite global space
843: Collective; No Fortran Support
845: Input Parameter:
846: . dm - the `DMCOMPOSITE` object
848: Output Parameter:
849: . ltogs - the individual mappings for each packed vector. Note that this includes
850: all the ghost points that individual ghosted `DMDA` may have.
852: Level: advanced
854: Note:
855: Each entry of `ltogs` should be destroyed with `ISLocalToGlobalMappingDestroy()`, `ltogs` should be freed with `PetscFree()`.
857: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
858: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetAccess()`, `DMCompositeScatter()`,
859: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
860: @*/
861: PetscErrorCode DMCompositeGetISLocalToGlobalMappings(DM dm, ISLocalToGlobalMapping *ltogs[])
862: {
863: PetscInt i, *idx, n, cnt;
864: struct DMCompositeLink *next;
865: PetscMPIInt rank;
866: DM_Composite *com = (DM_Composite *)dm->data;
867: PetscBool flg;
869: PetscFunctionBegin;
871: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
872: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
873: PetscCall(DMSetUp(dm));
874: PetscCall(PetscMalloc1(com->nDM, ltogs));
875: next = com->next;
876: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
878: /* loop over packed objects, handling one at a time */
879: cnt = 0;
880: while (next) {
881: ISLocalToGlobalMapping ltog;
882: PetscMPIInt size;
883: const PetscInt *suboff, *indices;
884: Vec global;
886: /* Get sub-DM global indices for each local dof */
887: PetscCall(DMGetLocalToGlobalMapping(next->dm, <og));
888: PetscCall(ISLocalToGlobalMappingGetSize(ltog, &n));
889: PetscCall(ISLocalToGlobalMappingGetIndices(ltog, &indices));
890: PetscCall(PetscMalloc1(n, &idx));
892: /* Get the offsets for the sub-DM global vector */
893: PetscCall(DMGetGlobalVector(next->dm, &global));
894: PetscCall(VecGetOwnershipRanges(global, &suboff));
895: PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)global), &size));
897: /* Shift the sub-DM definition of the global space to the composite global space */
898: for (i = 0; i < n; i++) {
899: PetscInt subi = indices[i], lo = 0, hi = size, t;
900: /* There's no consensus on what a negative index means,
901: except for skipping when setting the values in vectors and matrices */
902: if (subi < 0) {
903: idx[i] = subi - next->grstarts[rank];
904: continue;
905: }
906: /* Binary search to find which rank owns subi */
907: while (hi - lo > 1) {
908: t = lo + (hi - lo) / 2;
909: if (suboff[t] > subi) hi = t;
910: else lo = t;
911: }
912: idx[i] = subi - suboff[lo] + next->grstarts[lo];
913: }
914: PetscCall(ISLocalToGlobalMappingRestoreIndices(ltog, &indices));
915: PetscCall(ISLocalToGlobalMappingCreate(PetscObjectComm((PetscObject)dm), 1, n, idx, PETSC_OWN_POINTER, &(*ltogs)[cnt]));
916: PetscCall(DMRestoreGlobalVector(next->dm, &global));
917: next = next->next;
918: cnt++;
919: }
920: PetscFunctionReturn(PETSC_SUCCESS);
921: }
923: /*@C
924: DMCompositeGetLocalISs - Gets index sets for each component of a composite local vector
926: Not Collective; No Fortran Support
928: Input Parameter:
929: . dm - the `DMCOMPOSITE`
931: Output Parameter:
932: . is - array of serial index sets for each component of the `DMCOMPOSITE`
934: Level: intermediate
936: Notes:
937: At present, a composite local vector does not normally exist. This function is used to provide index sets for
938: `MatGetLocalSubMatrix()`. In the future, the scatters for each entry in the `DMCOMPOSITE` may be merged into a single
939: scatter to a composite local vector. The user should not typically need to know which is being done.
941: To get the composite global indices at all local points (including ghosts), use `DMCompositeGetISLocalToGlobalMappings()`.
943: To get index sets for pieces of the composite global vector, use `DMCompositeGetGlobalISs()`.
945: Each returned `IS` should be destroyed with `ISDestroy()`, the array should be freed with `PetscFree()`.
947: Fortran Note:
948: Pass in an array long enough to hold all the `IS`, see `DMCompositeGetNumberDM()`
950: .seealso: `DMCOMPOSITE`, `DM`, `DMCompositeGetGlobalISs()`, `DMCompositeGetISLocalToGlobalMappings()`, `MatGetLocalSubMatrix()`,
951: `MatCreateLocalRef()`, `DMCompositeGetNumberDM()`
952: @*/
953: PetscErrorCode DMCompositeGetLocalISs(DM dm, IS *is[])
954: {
955: DM_Composite *com = (DM_Composite *)dm->data;
956: struct DMCompositeLink *link;
957: PetscInt cnt, start;
958: PetscBool flg;
960: PetscFunctionBegin;
962: PetscAssertPointer(is, 2);
963: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
964: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
965: PetscCall(PetscMalloc1(com->nmine, is));
966: for (cnt = 0, start = 0, link = com->next; link; start += link->nlocal, cnt++, link = link->next) {
967: PetscInt bs;
968: PetscCall(ISCreateStride(PETSC_COMM_SELF, link->nlocal, start, 1, &(*is)[cnt]));
969: PetscCall(DMGetBlockSize(link->dm, &bs));
970: PetscCall(ISSetBlockSize((*is)[cnt], bs));
971: }
972: PetscFunctionReturn(PETSC_SUCCESS);
973: }
975: /*@C
976: DMCompositeGetGlobalISs - Gets the index sets for each composed object in a `DMCOMPOSITE`
978: Collective
980: Input Parameter:
981: . dm - the `DMCOMPOSITE` object
983: Output Parameter:
984: . is - the array of index sets
986: Level: advanced
988: Notes:
989: The `is` entries should be destroyed with `ISDestroy()`, `is` should be freed with `PetscFree()`
991: These could be used to extract a subset of vector entries for a "multi-physics" preconditioner
993: Use `DMCompositeGetLocalISs()` for index sets in the packed local numbering, and
994: `DMCompositeGetISLocalToGlobalMappings()` for to map local sub-`DM` (including ghost) indices to packed global
995: indices.
997: Fortran Notes:
998: The output argument 'is' must be an allocated array of sufficient length, which can be learned using `DMCompositeGetNumberDM()`.
1000: .seealso: `DMCOMPOSITE`, `DM`, `DMDestroy()`, `DMCompositeAddDM()`, `DMCreateGlobalVector()`,
1001: `DMCompositeGather()`, `DMCompositeCreate()`, `DMCompositeGetAccess()`, `DMCompositeScatter()`,
1002: `DMCompositeGetLocalVectors()`, `DMCompositeRestoreLocalVectors()`, `DMCompositeGetEntries()`
1003: @*/
1004: PetscErrorCode DMCompositeGetGlobalISs(DM dm, IS *is[])
1005: {
1006: PetscInt cnt = 0;
1007: struct DMCompositeLink *next;
1008: PetscMPIInt rank;
1009: DM_Composite *com = (DM_Composite *)dm->data;
1010: PetscBool flg;
1012: PetscFunctionBegin;
1014: PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMCOMPOSITE, &flg));
1015: PetscCheck(flg, PetscObjectComm((PetscObject)dm), PETSC_ERR_USER, "Not for type %s", ((PetscObject)dm)->type_name);
1016: PetscCheck(dm->setupcalled, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "Must call DMSetUp() before");
1017: PetscCall(PetscMalloc1(com->nDM, is));
1018: next = com->next;
1019: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
1021: /* loop over packed objects, handling one at a time */
1022: while (next) {
1023: PetscDS prob;
1025: PetscCall(ISCreateStride(PetscObjectComm((PetscObject)dm), next->n, next->grstart, 1, &(*is)[cnt]));
1026: PetscCall(DMGetDS(dm, &prob));
1027: if (prob) {
1028: MatNullSpace space;
1029: Mat pmat;
1030: PetscObject disc;
1031: PetscInt Nf;
1033: PetscCall(PetscDSGetNumFields(prob, &Nf));
1034: if (cnt < Nf) {
1035: PetscCall(PetscDSGetDiscretization(prob, cnt, &disc));
1036: PetscCall(PetscObjectQuery(disc, "nullspace", (PetscObject *)&space));
1037: if (space) PetscCall(PetscObjectCompose((PetscObject)(*is)[cnt], "nullspace", (PetscObject)space));
1038: PetscCall(PetscObjectQuery(disc, "nearnullspace", (PetscObject *)&space));
1039: if (space) PetscCall(PetscObjectCompose((PetscObject)(*is)[cnt], "nearnullspace", (PetscObject)space));
1040: PetscCall(PetscObjectQuery(disc, "pmat", (PetscObject *)&pmat));
1041: if (pmat) PetscCall(PetscObjectCompose((PetscObject)(*is)[cnt], "pmat", (PetscObject)pmat));
1042: }
1043: }
1044: cnt++;
1045: next = next->next;
1046: }
1047: PetscFunctionReturn(PETSC_SUCCESS);
1048: }
1050: static PetscErrorCode DMCreateFieldIS_Composite(DM dm, PetscInt *numFields, char ***fieldNames, IS **fields)
1051: {
1052: PetscInt nDM;
1053: DM *dms;
1054: PetscInt i;
1056: PetscFunctionBegin;
1057: PetscCall(DMCompositeGetNumberDM(dm, &nDM));
1058: if (numFields) *numFields = nDM;
1059: PetscCall(DMCompositeGetGlobalISs(dm, fields));
1060: if (fieldNames) {
1061: PetscCall(PetscMalloc1(nDM, &dms));
1062: PetscCall(PetscMalloc1(nDM, fieldNames));
1063: PetscCall(DMCompositeGetEntriesArray(dm, dms));
1064: for (i = 0; i < nDM; i++) {
1065: char buf[256];
1066: const char *splitname;
1068: /* Split naming precedence: object name, prefix, number */
1069: splitname = ((PetscObject)dm)->name;
1070: if (!splitname) {
1071: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)dms[i], &splitname));
1072: if (splitname) {
1073: size_t len;
1074: PetscCall(PetscStrncpy(buf, splitname, sizeof(buf)));
1075: buf[sizeof(buf) - 1] = 0;
1076: PetscCall(PetscStrlen(buf, &len));
1077: if (buf[len - 1] == '_') buf[len - 1] = 0; /* Remove trailing underscore if it was used */
1078: splitname = buf;
1079: }
1080: }
1081: if (!splitname) {
1082: PetscCall(PetscSNPrintf(buf, sizeof(buf), "%" PetscInt_FMT, i));
1083: splitname = buf;
1084: }
1085: PetscCall(PetscStrallocpy(splitname, &(*fieldNames)[i]));
1086: }
1087: PetscCall(PetscFree(dms));
1088: }
1089: PetscFunctionReturn(PETSC_SUCCESS);
1090: }
1092: /*
1093: This could take over from DMCreateFieldIS(), as it is more general,
1094: making DMCreateFieldIS() a special case -- calling with dmlist == NULL;
1095: At this point it's probably best to be less intrusive, however.
1096: */
1097: static PetscErrorCode DMCreateFieldDecomposition_Composite(DM dm, PetscInt *len, char ***namelist, IS **islist, DM **dmlist)
1098: {
1099: PetscInt nDM;
1100: PetscInt i;
1102: PetscFunctionBegin;
1103: PetscCall(DMCreateFieldIS_Composite(dm, len, namelist, islist));
1104: if (dmlist) {
1105: PetscCall(DMCompositeGetNumberDM(dm, &nDM));
1106: PetscCall(PetscMalloc1(nDM, dmlist));
1107: PetscCall(DMCompositeGetEntriesArray(dm, *dmlist));
1108: for (i = 0; i < nDM; i++) PetscCall(PetscObjectReference((PetscObject)((*dmlist)[i])));
1109: }
1110: PetscFunctionReturn(PETSC_SUCCESS);
1111: }
1113: /*@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 a 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 a 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: /*@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: Use `DMCompositeGetEntriesArray()`
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: /*@
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: PetscCall(ISColoringValueCast(com->N, &maxcol));
1489: for (i = 0; i < n; i++) PetscCall(ISColoringValueCast(com->rstart + i, colors + i));
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++) PetscCall(ISColoringValueCast(maxcol + lcoloring->colors[i], colors + cnt++));
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: }