Actual source code: ex16.c
2: static char help[] = "Tests DMComposite routines.\n\n";
4: #include <petscdmredundant.h>
5: #include <petscdm.h>
6: #include <petscdmda.h>
7: #include <petscdmcomposite.h>
8: #include <petscpf.h>
10: int main(int argc, char **argv)
11: {
12: PetscInt nredundant1 = 5, nredundant2 = 2, i;
13: ISLocalToGlobalMapping *ltog;
14: PetscMPIInt rank, size;
15: DM packer;
16: Vec global, local1, local2, redundant1, redundant2;
17: PF pf;
18: DM da1, da2, dmred1, dmred2;
19: PetscScalar *redundant1a, *redundant2a;
20: PetscViewer sviewer;
21: PetscBool gather_add = PETSC_FALSE;
23: PetscFunctionBeginUser;
24: PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
25: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
26: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
28: PetscCall(PetscOptionsGetBool(NULL, NULL, "-gather_add", &gather_add, NULL));
30: PetscCall(DMCompositeCreate(PETSC_COMM_WORLD, &packer));
32: PetscCall(DMRedundantCreate(PETSC_COMM_WORLD, 0, nredundant1, &dmred1));
33: PetscCall(DMCreateLocalVector(dmred1, &redundant1));
34: PetscCall(DMCompositeAddDM(packer, dmred1));
36: PetscCall(DMDACreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, 8, 1, 1, NULL, &da1));
37: PetscCall(DMSetFromOptions(da1));
38: PetscCall(DMSetUp(da1));
39: PetscCall(DMCreateLocalVector(da1, &local1));
40: PetscCall(DMCompositeAddDM(packer, da1));
42: PetscCall(DMRedundantCreate(PETSC_COMM_WORLD, 1 % size, nredundant2, &dmred2));
43: PetscCall(DMCreateLocalVector(dmred2, &redundant2));
44: PetscCall(DMCompositeAddDM(packer, dmred2));
46: PetscCall(DMDACreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, 6, 1, 1, NULL, &da2));
47: PetscCall(DMSetFromOptions(da2));
48: PetscCall(DMSetUp(da2));
49: PetscCall(DMCreateLocalVector(da2, &local2));
50: PetscCall(DMCompositeAddDM(packer, da2));
52: PetscCall(DMCreateGlobalVector(packer, &global));
53: PetscCall(PFCreate(PETSC_COMM_WORLD, 1, 1, &pf));
54: PetscCall(PFSetType(pf, PFIDENTITY, NULL));
55: PetscCall(PFApplyVec(pf, NULL, global));
56: PetscCall(PFDestroy(&pf));
57: PetscCall(VecView(global, PETSC_VIEWER_STDOUT_WORLD));
59: PetscCall(DMCompositeScatter(packer, global, redundant1, local1, redundant2, local2));
60: PetscCall(PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD));
61: PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD, "[%d] My part of redundant1 vector\n", rank));
62: PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer));
63: PetscCall(VecView(redundant1, sviewer));
64: PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer));
65: PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD));
66: PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD, "[%d] My part of da1 vector\n", rank));
67: PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer));
68: PetscCall(VecView(local1, sviewer));
69: PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer));
70: PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD));
71: PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD, "[%d] My part of redundant2 vector\n", rank));
72: PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer));
73: PetscCall(VecView(redundant2, sviewer));
74: PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer));
75: PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD));
76: PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD, "[%d] My part of da2 vector\n", rank));
77: PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer));
78: PetscCall(VecView(local2, sviewer));
79: PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer));
80: PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD));
81: PetscCall(PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD));
83: PetscCall(VecGetArray(redundant1, &redundant1a));
84: PetscCall(VecGetArray(redundant2, &redundant2a));
85: for (i = 0; i < nredundant1; i++) redundant1a[i] = (rank + 2) * i;
86: for (i = 0; i < nredundant2; i++) redundant2a[i] = (rank + 10) * i;
87: PetscCall(VecRestoreArray(redundant1, &redundant1a));
88: PetscCall(VecRestoreArray(redundant2, &redundant2a));
90: PetscCall(DMCompositeGather(packer, gather_add ? ADD_VALUES : INSERT_VALUES, global, redundant1, local1, redundant2, local2));
91: PetscCall(VecView(global, PETSC_VIEWER_STDOUT_WORLD));
93: /* get the global numbering for each subvector element */
94: PetscCall(DMCompositeGetISLocalToGlobalMappings(packer, <og));
96: PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "Local to global mapping of redundant1 vector\n"));
97: PetscCall(ISLocalToGlobalMappingView(ltog[0], PETSC_VIEWER_STDOUT_WORLD));
98: PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "Local to global mapping of local1 vector\n"));
99: PetscCall(ISLocalToGlobalMappingView(ltog[1], PETSC_VIEWER_STDOUT_WORLD));
100: PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "Local to global mapping of redundant2 vector\n"));
101: PetscCall(ISLocalToGlobalMappingView(ltog[2], PETSC_VIEWER_STDOUT_WORLD));
102: PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "Local to global mapping of local2 vector\n"));
103: PetscCall(ISLocalToGlobalMappingView(ltog[3], PETSC_VIEWER_STDOUT_WORLD));
105: for (i = 0; i < 4; i++) PetscCall(ISLocalToGlobalMappingDestroy(<og[i]));
106: PetscCall(PetscFree(ltog));
108: PetscCall(DMDestroy(&da1));
109: PetscCall(DMDestroy(&dmred1));
110: PetscCall(DMDestroy(&dmred2));
111: PetscCall(DMDestroy(&da2));
112: PetscCall(VecDestroy(&redundant1));
113: PetscCall(VecDestroy(&redundant2));
114: PetscCall(VecDestroy(&local1));
115: PetscCall(VecDestroy(&local2));
116: PetscCall(VecDestroy(&global));
117: PetscCall(DMDestroy(&packer));
118: PetscCall(PetscFinalize());
119: return 0;
120: }
122: /*TEST
124: build:
125: requires: !complex
127: test:
128: nsize: 3
130: test:
131: suffix: 2
132: nsize: 3
133: args: -gather_add
135: TEST*/