1 const char help[] = "Test clearing stale AMR data (example contributed by Berend van Wachem)";
2
3 #include <petscdmplex.h>
4 #include <petscdmforest.h>
5
CloneDMWithNewSection(DM OriginalDM,DM * NewDM,PetscInt NFields)6 PetscErrorCode CloneDMWithNewSection(DM OriginalDM, DM *NewDM, PetscInt NFields)
7 {
8 PetscSection section;
9 PetscInt *NumComp, *NumDof;
10
11 PetscFunctionBegin;
12 PetscCall(DMClone(OriginalDM, NewDM));
13 PetscCall(DMPlexDistributeSetDefault(*NewDM, PETSC_FALSE));
14 PetscCall(DMClearDS(*NewDM));
15 PetscCall(PetscCalloc2(1, &NumComp, 3, &NumDof));
16 NumComp[0] = 1;
17 NumDof[2] = NFields;
18 PetscCall(DMSetNumFields(*NewDM, 1));
19 PetscCall(DMSetFromOptions(*NewDM));
20 PetscCall(DMPlexCreateSection(*NewDM, NULL, NumComp, NumDof, 0, NULL, NULL, NULL, NULL, §ion));
21 PetscCall(DMSetLocalSection(*NewDM, section));
22 PetscCall(PetscFree2(NumComp, NumDof));
23 PetscCall(PetscSectionDestroy(§ion));
24 PetscFE fe;
25 PetscCall(PetscFECreateDefault(PETSC_COMM_WORLD, 2, 1, PETSC_FALSE, NULL, PETSC_DEFAULT, &fe));
26 PetscCall(DMSetField(*NewDM, 0, NULL, (PetscObject)fe));
27 PetscCall(PetscFEDestroy(&fe));
28 PetscCall(DMCreateDS(*NewDM));
29 PetscFunctionReturn(PETSC_SUCCESS);
30 }
31
main(int argc,char ** argv)32 int main(int argc, char **argv)
33 {
34 PetscInt dim = 2;
35 PetscInt cells_per_dir[] = {1, 1};
36 PetscReal dir_min[] = {0.0, 0.0};
37 PetscReal dir_max[] = {1.0, 1.0};
38 DMBoundaryType bcs[] = {DM_BOUNDARY_NONE, DM_BOUNDARY_NONE};
39 DM forest, plex, NewDM;
40 Vec NewDMVecGlobal, NewDMVecLocal;
41
42 PetscFunctionBeginUser;
43 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
44 PetscCall(DMCreate(PETSC_COMM_WORLD, &forest));
45 PetscCall(DMSetType(forest, DMP4EST));
46 {
47 DM dm_base;
48 PetscCall(DMPlexCreateBoxMesh(PETSC_COMM_WORLD, dim, /* simplex */ PETSC_FALSE, cells_per_dir, dir_min, dir_max, bcs, /* interpolate */ PETSC_TRUE, 0, PETSC_TRUE, &dm_base));
49 PetscCall(DMSetFromOptions(dm_base));
50 PetscCall(DMViewFromOptions(dm_base, NULL, "-dm_base_view"));
51 PetscCall(DMCopyFields(dm_base, PETSC_DETERMINE, PETSC_DETERMINE, forest));
52 PetscCall(DMForestSetBaseDM(forest, dm_base));
53 PetscCall(DMDestroy(&dm_base));
54 }
55 PetscCall(DMSetFromOptions(forest));
56 PetscCall(DMSetUp(forest));
57
58 PetscCall(DMViewFromOptions(forest, NULL, "-dm_forest_view"));
59
60 PetscCall(DMConvert(forest, DMPLEX, &plex));
61
62 PetscInt numFields = 2;
63 PetscInt numComp[2] = {1, 1};
64 PetscInt numDof[6] = {0};
65 for (PetscInt i = 0; i < numFields; i++) numDof[i * (dim + 1) + dim] = 1;
66
67 PetscCall(DMSetNumFields(plex, numFields));
68 PetscCall(DMSetNumFields(forest, numFields));
69
70 PetscSection section;
71 PetscCall(DMPlexCreateSection(plex, NULL, numComp, numDof, 0, NULL, NULL, NULL, NULL, §ion));
72
73 const char *names[] = {"field 0", "field 1"};
74 for (PetscInt i = 0; i < numFields; i++) PetscCall(PetscSectionSetFieldName(section, i, names[i]));
75 PetscCall(DMSetLocalSection(forest, section));
76 PetscCall(DMSetLocalSection(plex, section));
77 PetscCall(PetscSectionDestroy(§ion));
78
79 PetscFE fe;
80 PetscCall(PetscFECreateDefault(PETSC_COMM_WORLD, dim, 1, PETSC_FALSE, NULL, PETSC_DEFAULT, &fe));
81 for (PetscInt i = 0; i < numFields; i++) {
82 PetscCall(DMSetField(plex, i, NULL, (PetscObject)fe));
83 PetscCall(DMSetField(forest, i, NULL, (PetscObject)fe));
84 }
85 PetscCall(PetscFEDestroy(&fe));
86
87 PetscCall(DMCreateDS(plex));
88 PetscCall(DMCreateDS(forest));
89
90 /* Make another DM, based on the layout of the previous DM, but with a different number of fields */
91 PetscCall(CloneDMWithNewSection(plex, &NewDM, 1));
92 PetscCall(DMCreateDS(NewDM));
93 PetscCall(DMCreateGlobalVector(NewDM, &NewDMVecGlobal));
94 PetscCall(DMCreateLocalVector(NewDM, &NewDMVecLocal));
95 PetscCall(VecSet(NewDMVecGlobal, 3.141592));
96 PetscCall(DMGlobalToLocal(NewDM, NewDMVecGlobal, INSERT_VALUES, NewDMVecLocal));
97 PetscCall(VecDestroy(&NewDMVecLocal));
98 PetscCall(VecDestroy(&NewDMVecGlobal));
99 PetscCall(DMClearDS(NewDM));
100 PetscCall(DMDestroy(&NewDM));
101
102 Vec g_vec, l_vec;
103 PetscCall(DMCreateGlobalVector(plex, &g_vec));
104 PetscCall(VecSet(g_vec, 1.0));
105 PetscCall(DMCreateLocalVector(plex, &l_vec));
106 PetscCall(DMGlobalToLocal(plex, g_vec, INSERT_VALUES, l_vec));
107 PetscCall(VecViewFromOptions(l_vec, NULL, "-local_vec_view"));
108 PetscCall(VecDestroy(&l_vec));
109 PetscCall(VecDestroy(&g_vec));
110
111 PetscCall(DMViewFromOptions(forest, NULL, "-dm_plex_view"));
112 PetscCall(DMDestroy(&plex));
113 PetscCall(DMDestroy(&forest));
114 PetscCall(PetscFinalize());
115 return 0;
116 }
117
118 /*TEST
119
120 test:
121 suffix: 0
122 requires: p4est
123 args: -dm_forest_initial_refinement 1 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash
124 output_file: output/empty.out
125
126 TEST*/
127