Actual source code: ex52.c

  1: const char help[] = "Test clearing stale AMR data (example contributed by Berend van Wachem)";

  3: #include <petscdmplex.h>
  4: #include <petscdmforest.h>

  6: PetscErrorCode CloneDMWithNewSection(DM OriginalDM, DM *NewDM, PetscInt NFields)
  7: {
  8:   PetscSection section;
  9:   PetscInt    *NumComp, *NumDof;

 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, &section));
 21:   PetscCall(DMSetLocalSection(*NewDM, section));
 22:   PetscCall(PetscFree2(NumComp, NumDof));
 23:   PetscCall(PetscSectionDestroy(&section));
 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: }

 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;

 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));

 58:   PetscCall(DMViewFromOptions(forest, NULL, "-dm_forest_view"));

 60:   PetscCall(DMConvert(forest, DMPLEX, &plex));

 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;

 67:   PetscCall(DMSetNumFields(plex, numFields));
 68:   PetscCall(DMSetNumFields(forest, numFields));

 70:   PetscSection section;
 71:   PetscCall(DMPlexCreateSection(plex, NULL, numComp, numDof, 0, NULL, NULL, NULL, NULL, &section));

 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(&section));

 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));

 87:   PetscCall(DMCreateDS(plex));
 88:   PetscCall(DMCreateDS(forest));

 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));

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));

111:   PetscCall(DMViewFromOptions(forest, NULL, "-dm_plex_view"));
112:   PetscCall(DMDestroy(&plex));
113:   PetscCall(DMDestroy(&forest));
114:   PetscCall(PetscFinalize());
115:   return 0;
116: }

118: /*TEST

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

125: TEST*/