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:   PetscFunctionBegin;

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

 30:   PetscFunctionReturn(PETSC_SUCCESS);
 31: }

 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;

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

 60:   PetscCall(DMViewFromOptions(forest, NULL, "-dm_forest_view"));
 61:   DM plex;

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

 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;

 70:   PetscCall(DMSetNumFields(plex, numFields));
 71:   PetscCall(DMSetNumFields(forest, numFields));

 73:   PetscSection section;
 74:   PetscCall(DMPlexCreateSection(plex, NULL, numComp, numDof, 0, NULL, NULL, NULL, NULL, &section));

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

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

 90:   PetscCall(DMCreateDS(plex));
 91:   PetscCall(DMCreateDS(forest));

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

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

114:   PetscCall(DMViewFromOptions(forest, NULL, "-dm_plex_view"));
115:   PetscCall(DMDestroy(&plex));
116:   PetscCall(DMDestroy(&forest));
117:   PetscCall(PetscFinalize());
118:   return 0;
119: }

121: /*TEST

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

128: TEST*/