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