xref: /petsc/src/dm/impls/forest/tests/ex2.c (revision b122ec5aa1bd4469eb4e0673542fb7de3f411254)
1c4762a1bSJed Brown static char help[] = "Create a mesh, refine and coarsen simultaneously, and transfer a field\n\n";
2c4762a1bSJed Brown 
3c4762a1bSJed Brown #include <petscds.h>
4c4762a1bSJed Brown #include <petscdmplex.h>
5c4762a1bSJed Brown #include <petscdmforest.h>
6c4762a1bSJed Brown #include <petscoptions.h>
7c4762a1bSJed Brown 
8c4762a1bSJed Brown static PetscErrorCode AddIdentityLabel(DM dm)
9c4762a1bSJed Brown {
10c4762a1bSJed Brown   PetscInt       pStart,pEnd,p;
11c4762a1bSJed Brown 
12c4762a1bSJed Brown   PetscFunctionBegin;
135f80ce2aSJacob Faibussowitsch   CHKERRQ(DMCreateLabel(dm, "identity"));
145f80ce2aSJacob Faibussowitsch   CHKERRQ(DMPlexGetChart(dm, &pStart, &pEnd));
155f80ce2aSJacob Faibussowitsch   for (p = pStart; p < pEnd; p++) CHKERRQ(DMSetLabelValue(dm, "identity", p, p));
16c4762a1bSJed Brown   PetscFunctionReturn(0);
17c4762a1bSJed Brown }
18c4762a1bSJed Brown 
19c4762a1bSJed Brown static PetscErrorCode CreateAdaptivityLabel(DM forest,DMLabel *adaptLabel)
20c4762a1bSJed Brown {
21c4762a1bSJed Brown   DMLabel        identLabel;
22c4762a1bSJed Brown   PetscInt       cStart, cEnd, c;
23c4762a1bSJed Brown 
24c4762a1bSJed Brown   PetscFunctionBegin;
255f80ce2aSJacob Faibussowitsch   CHKERRQ(DMLabelCreate(PETSC_COMM_SELF,"adapt",adaptLabel));
265f80ce2aSJacob Faibussowitsch   CHKERRQ(DMLabelSetDefaultValue(*adaptLabel,DM_ADAPT_COARSEN));
275f80ce2aSJacob Faibussowitsch   CHKERRQ(DMGetLabel(forest,"identity",&identLabel));
285f80ce2aSJacob Faibussowitsch   CHKERRQ(DMForestGetCellChart(forest,&cStart,&cEnd));
29c4762a1bSJed Brown   for (c = cStart; c < cEnd; c++) {
30c4762a1bSJed Brown     PetscInt basePoint;
31c4762a1bSJed Brown 
325f80ce2aSJacob Faibussowitsch     CHKERRQ(DMLabelGetValue(identLabel,c,&basePoint));
335f80ce2aSJacob Faibussowitsch     if (!basePoint) CHKERRQ(DMLabelSetValue(*adaptLabel,c,DM_ADAPT_REFINE));
34c4762a1bSJed Brown   }
35c4762a1bSJed Brown   PetscFunctionReturn(0);
36c4762a1bSJed Brown }
37c4762a1bSJed Brown 
38c4762a1bSJed Brown static PetscErrorCode LinearFunction(PetscInt dim,PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar u[], void *ctx)
39c4762a1bSJed Brown {
40c4762a1bSJed Brown   PetscFunctionBeginUser;
41c4762a1bSJed Brown   u[0] = (x[0] * 2.0 + 1.) + (x[1] * 20.0 + 10.) + ((dim == 3) ? (x[2] * 200.0 + 100.) : 0.);
42c4762a1bSJed Brown   PetscFunctionReturn(0);
43c4762a1bSJed Brown }
44c4762a1bSJed Brown 
45c4762a1bSJed Brown static PetscErrorCode MultiaffineFunction(PetscInt dim,PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar u[], void *ctx)
46c4762a1bSJed Brown {
47c4762a1bSJed Brown   PetscFunctionBeginUser;
48c4762a1bSJed Brown   u[0] = (x[0] * 1.0 + 2.0) * (x[1] * 3.0 - 4.0) * ((dim == 3) ? (x[2] * 5.0 + 6.0) : 1.);
49c4762a1bSJed Brown   PetscFunctionReturn(0);
50c4762a1bSJed Brown }
51c4762a1bSJed Brown 
52c4762a1bSJed Brown static PetscErrorCode CoordsFunction(PetscInt dim,PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar u[], void *ctx)
53c4762a1bSJed Brown {
54c4762a1bSJed Brown   PetscInt f;
55c4762a1bSJed Brown 
56c4762a1bSJed Brown   PetscFunctionBeginUser;
57c4762a1bSJed Brown   for (f=0;f<Nf;f++) u[f] = x[f];
58c4762a1bSJed Brown   PetscFunctionReturn(0);
59c4762a1bSJed Brown }
60c4762a1bSJed Brown 
61c4762a1bSJed Brown typedef struct _bc_func_ctx
62c4762a1bSJed Brown {
63c4762a1bSJed Brown   PetscErrorCode (*func) (PetscInt,PetscReal,const PetscReal [], PetscInt, PetscScalar [], void *);
64c4762a1bSJed Brown   PetscInt dim;
65c4762a1bSJed Brown   PetscInt Nf;
66c4762a1bSJed Brown   void *ctx;
67c4762a1bSJed Brown }
68c4762a1bSJed Brown bc_func_ctx;
69c4762a1bSJed Brown 
70c4762a1bSJed Brown static PetscErrorCode bc_func_fv (PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *xI, PetscScalar *xG, void *ctx)
71c4762a1bSJed Brown {
72c4762a1bSJed Brown   bc_func_ctx    *bcCtx;
73c4762a1bSJed Brown 
74c4762a1bSJed Brown   PetscFunctionBegin;
75c4762a1bSJed Brown   bcCtx = (bc_func_ctx *) ctx;
765f80ce2aSJacob Faibussowitsch   CHKERRQ((bcCtx->func)(bcCtx->dim,time,c,bcCtx->Nf,xG,bcCtx->ctx));
77c4762a1bSJed Brown   PetscFunctionReturn(0);
78c4762a1bSJed Brown }
79c4762a1bSJed Brown 
80c4762a1bSJed Brown static PetscErrorCode IdentifyBadPoints (DM dm, Vec vec, PetscReal tol)
81c4762a1bSJed Brown {
82c4762a1bSJed Brown   DM             dmplex;
83c4762a1bSJed Brown   PetscInt       p, pStart, pEnd, maxDof;
84c4762a1bSJed Brown   Vec            vecLocal;
85c4762a1bSJed Brown   DMLabel        depthLabel;
86c4762a1bSJed Brown   PetscSection   section;
87c4762a1bSJed Brown 
88c4762a1bSJed Brown   PetscFunctionBegin;
895f80ce2aSJacob Faibussowitsch   CHKERRQ(DMCreateLocalVector(dm, &vecLocal));
905f80ce2aSJacob Faibussowitsch   CHKERRQ(DMGlobalToLocalBegin(dm, vec, INSERT_VALUES, vecLocal));
915f80ce2aSJacob Faibussowitsch   CHKERRQ(DMGlobalToLocalEnd(dm, vec, INSERT_VALUES, vecLocal));
925f80ce2aSJacob Faibussowitsch   CHKERRQ(DMConvert(dm ,DMPLEX, &dmplex));
935f80ce2aSJacob Faibussowitsch   CHKERRQ(DMPlexGetChart(dmplex, &pStart, &pEnd));
945f80ce2aSJacob Faibussowitsch   CHKERRQ(DMPlexGetDepthLabel(dmplex, &depthLabel));
955f80ce2aSJacob Faibussowitsch   CHKERRQ(DMGetLocalSection(dmplex, &section));
965f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscSectionGetMaxDof(section, &maxDof));
97c4762a1bSJed Brown   for (p = pStart; p < pEnd; p++) {
98c4762a1bSJed Brown     PetscInt     s, c, cSize, parent, childID, numChildren;
99c4762a1bSJed Brown     PetscInt     cl, closureSize, *closure = NULL;
100c4762a1bSJed Brown     PetscScalar *values = NULL;
101c4762a1bSJed Brown     PetscBool    bad = PETSC_FALSE;
102c4762a1bSJed Brown 
1035f80ce2aSJacob Faibussowitsch     CHKERRQ(VecGetValuesSection(vecLocal, section, p, &values));
1045f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscSectionGetDof(section, p, &cSize));
105c4762a1bSJed Brown     for (c = 0; c < cSize; c++) {
1062f613bf5SBarry Smith       PetscReal absDiff = PetscAbsScalar(values[c]);
107c4762a1bSJed Brown       if (absDiff > tol) {bad = PETSC_TRUE; break;}
108c4762a1bSJed Brown     }
109c4762a1bSJed Brown     if (!bad) continue;
1105f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscPrintf(PETSC_COMM_SELF, "Bad point %D\n", p));
1115f80ce2aSJacob Faibussowitsch     CHKERRQ(DMLabelGetValue(depthLabel, p, &s));
1125f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscPrintf(PETSC_COMM_SELF, "  Depth %D\n", s));
1135f80ce2aSJacob Faibussowitsch     CHKERRQ(DMPlexGetTransitiveClosure(dmplex, p, PETSC_TRUE, &closureSize, &closure));
114c4762a1bSJed Brown     for (cl = 0; cl < closureSize; cl++) {
115c4762a1bSJed Brown       PetscInt cp = closure[2 * cl];
1165f80ce2aSJacob Faibussowitsch       CHKERRQ(DMPlexGetTreeParent(dmplex, cp, &parent, &childID));
117c4762a1bSJed Brown       if (parent != cp) {
1185f80ce2aSJacob Faibussowitsch         CHKERRQ(PetscPrintf(PETSC_COMM_SELF, "  Closure point %D (%D) child of %D (ID %D)\n", cl, cp, parent, childID));
119c4762a1bSJed Brown       }
1205f80ce2aSJacob Faibussowitsch       CHKERRQ(DMPlexGetTreeChildren(dmplex, cp, &numChildren, NULL));
121c4762a1bSJed Brown       if (numChildren) {
1225f80ce2aSJacob Faibussowitsch         CHKERRQ(PetscPrintf(PETSC_COMM_SELF, "  Closure point %D (%D) is parent\n", cl, cp));
123c4762a1bSJed Brown       }
124c4762a1bSJed Brown     }
1255f80ce2aSJacob Faibussowitsch     CHKERRQ(DMPlexRestoreTransitiveClosure(dmplex, p, PETSC_TRUE, &closureSize, &closure));
126c4762a1bSJed Brown     for (c = 0; c < cSize; c++) {
1272f613bf5SBarry Smith       PetscReal absDiff = PetscAbsScalar(values[c]);
128c4762a1bSJed Brown       if (absDiff > tol) {
1295f80ce2aSJacob Faibussowitsch         CHKERRQ(PetscPrintf(PETSC_COMM_SELF, "  Bad dof %D\n", c));
130c4762a1bSJed Brown       }
131c4762a1bSJed Brown     }
132c4762a1bSJed Brown   }
1335f80ce2aSJacob Faibussowitsch   CHKERRQ(DMDestroy(&dmplex));
1345f80ce2aSJacob Faibussowitsch   CHKERRQ(VecDestroy(&vecLocal));
135c4762a1bSJed Brown   PetscFunctionReturn(0);
136c4762a1bSJed Brown }
137c4762a1bSJed Brown 
138c4762a1bSJed Brown int main(int argc, char **argv)
139c4762a1bSJed Brown {
140c4762a1bSJed Brown   MPI_Comm       comm;
141c4762a1bSJed Brown   DM             base, preForest, postForest;
14230602db0SMatthew G. Knepley   PetscInt       dim, Nf = 1;
143c4762a1bSJed Brown   PetscInt       step, adaptSteps = 1;
144c4762a1bSJed Brown   PetscInt       preCount, postCount;
145c4762a1bSJed Brown   Vec            preVec, postVecTransfer, postVecExact;
146c4762a1bSJed Brown   PetscErrorCode (*funcs[1]) (PetscInt,PetscReal,const PetscReal [],PetscInt,PetscScalar [], void *) = {MultiaffineFunction};
147c4762a1bSJed Brown   void           *ctxs[1] = {NULL};
148c4762a1bSJed Brown   PetscReal      diff, tol = PETSC_SMALL;
149c4762a1bSJed Brown   PetscBool      linear = PETSC_FALSE;
150c4762a1bSJed Brown   PetscBool      coords = PETSC_FALSE;
151c4762a1bSJed Brown   PetscBool      useFV = PETSC_FALSE;
152c4762a1bSJed Brown   PetscBool      conv = PETSC_FALSE;
153c4762a1bSJed Brown   PetscBool      transfer_from_base[2] = {PETSC_TRUE,PETSC_FALSE};
154c4762a1bSJed Brown   PetscBool      use_bcs = PETSC_TRUE;
155c4762a1bSJed Brown   bc_func_ctx    bcCtx;
156c4762a1bSJed Brown   DMLabel        adaptLabel;
157c4762a1bSJed Brown   PetscErrorCode ierr;
158c4762a1bSJed Brown 
159*b122ec5aSJacob Faibussowitsch   CHKERRQ(PetscInitialize(&argc, &argv, NULL,help));
160c4762a1bSJed Brown   comm = PETSC_COMM_WORLD;
161c4762a1bSJed Brown   ierr = PetscOptionsBegin(comm, "", "DMForestTransferVec() Test Options", "DMFOREST");CHKERRQ(ierr);
1625f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsBool("-linear","Transfer a simple linear function", "ex2.c", linear, &linear, NULL));
1635f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsBool("-coords","Transfer a simple coordinate function", "ex2.c", coords, &coords, NULL));
1645f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsBool("-use_fv","Use a finite volume approximation", "ex2.c", useFV, &useFV, NULL));
1655f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsBool("-test_convert","Test conversion to DMPLEX",NULL,conv,&conv,NULL));
1665f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsBool("-transfer_from_base","Transfer a vector from base DM to DMForest", "ex2.c", transfer_from_base[0], &transfer_from_base[0], NULL));
167c4762a1bSJed Brown   transfer_from_base[1] = transfer_from_base[0];
1685f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsBool("-transfer_from_base_steps","Transfer a vector from base DM to the latest DMForest after the adaptivity steps", "ex2.c", transfer_from_base[1], &transfer_from_base[1], NULL));
1695f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsBool("-use_bcs","Use dirichlet boundary conditions", "ex2.c", use_bcs, &use_bcs, NULL));
1705f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsBoundedInt("-adapt_steps","Number of adaptivity steps", "ex2.c", adaptSteps, &adaptSteps, NULL,0));
171c4762a1bSJed Brown   ierr = PetscOptionsEnd();CHKERRQ(ierr);
172c4762a1bSJed Brown 
173c4762a1bSJed Brown   tol = PetscMax(1.e-10,tol); /* XXX fix for quadruple precision -> why do I need to do this? */
174c4762a1bSJed Brown 
17530602db0SMatthew G. Knepley   /* the base mesh */
1765f80ce2aSJacob Faibussowitsch   CHKERRQ(DMCreate(comm, &base));
1775f80ce2aSJacob Faibussowitsch   CHKERRQ(DMSetType(base, DMPLEX));
1785f80ce2aSJacob Faibussowitsch   CHKERRQ(DMSetFromOptions(base));
17930602db0SMatthew G. Knepley 
1805f80ce2aSJacob Faibussowitsch   CHKERRQ(AddIdentityLabel(base));
1815f80ce2aSJacob Faibussowitsch   CHKERRQ(DMGetDimension(base, &dim));
182c4762a1bSJed Brown 
183c4762a1bSJed Brown   if (linear) {
184c4762a1bSJed Brown     funcs[0] = LinearFunction;
185c4762a1bSJed Brown   }
186c4762a1bSJed Brown   if (coords) {
187c4762a1bSJed Brown     funcs[0] = CoordsFunction;
188c4762a1bSJed Brown     Nf = dim;
189c4762a1bSJed Brown   }
190c4762a1bSJed Brown 
191c4762a1bSJed Brown   bcCtx.func = funcs[0];
192c4762a1bSJed Brown   bcCtx.dim  = dim;
193c4762a1bSJed Brown   bcCtx.Nf   = Nf;
194c4762a1bSJed Brown   bcCtx.ctx  = NULL;
195c4762a1bSJed Brown 
196c4762a1bSJed Brown   if (useFV) {
197c4762a1bSJed Brown     PetscFV      fv;
198c4762a1bSJed Brown     PetscLimiter limiter;
199c4762a1bSJed Brown     DM           baseFV;
200c4762a1bSJed Brown 
2015f80ce2aSJacob Faibussowitsch     CHKERRQ(DMPlexConstructGhostCells(base,NULL,NULL,&baseFV));
2025f80ce2aSJacob Faibussowitsch     CHKERRQ(DMViewFromOptions(baseFV, NULL, "-fv_dm_view"));
2035f80ce2aSJacob Faibussowitsch     CHKERRQ(DMDestroy(&base));
204c4762a1bSJed Brown     base = baseFV;
2055f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscFVCreate(comm, &fv));
2065f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscFVSetSpatialDimension(fv,dim));
2075f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscFVSetType(fv,PETSCFVLEASTSQUARES));
2085f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscFVSetNumComponents(fv,Nf));
2095f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscLimiterCreate(comm,&limiter));
2105f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscLimiterSetType(limiter,PETSCLIMITERNONE));
2115f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscFVSetLimiter(fv,limiter));
2125f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscLimiterDestroy(&limiter));
2135f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscFVSetFromOptions(fv));
2145f80ce2aSJacob Faibussowitsch     CHKERRQ(DMSetField(base,0,NULL,(PetscObject)fv));
2155f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscFVDestroy(&fv));
216c4762a1bSJed Brown   } else {
217c4762a1bSJed Brown     PetscFE fe;
218c4762a1bSJed Brown 
2195f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscFECreateDefault(comm,dim,Nf,PETSC_FALSE,NULL,PETSC_DEFAULT,&fe));
2205f80ce2aSJacob Faibussowitsch     CHKERRQ(DMSetField(base,0,NULL,(PetscObject)fe));
2215f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscFEDestroy(&fe));
222c4762a1bSJed Brown   }
2235f80ce2aSJacob Faibussowitsch   CHKERRQ(DMCreateDS(base));
224c4762a1bSJed Brown 
225c4762a1bSJed Brown   if (use_bcs) {
226c4762a1bSJed Brown     PetscInt ids[] = {1, 2, 3, 4, 5, 6};
22745480ffeSMatthew G. Knepley     DMLabel  label;
228c4762a1bSJed Brown 
2295f80ce2aSJacob Faibussowitsch     CHKERRQ(DMGetLabel(base, "marker", &label));
2305f80ce2aSJacob Faibussowitsch     CHKERRQ(DMAddBoundary(base,DM_BC_ESSENTIAL, "bc", label, 2 * dim, ids, 0, 0, NULL, useFV ? (void(*)(void)) bc_func_fv : (void(*)(void)) funcs[0], NULL, useFV ? (void *) &bcCtx : NULL, NULL));
231c4762a1bSJed Brown   }
2325f80ce2aSJacob Faibussowitsch   CHKERRQ(DMViewFromOptions(base,NULL,"-dm_base_view"));
233c4762a1bSJed Brown 
234c4762a1bSJed Brown   /* the pre adaptivity forest */
2355f80ce2aSJacob Faibussowitsch   CHKERRQ(DMCreate(comm,&preForest));
2365f80ce2aSJacob Faibussowitsch   CHKERRQ(DMSetType(preForest,(dim == 2) ? DMP4EST : DMP8EST));
2375f80ce2aSJacob Faibussowitsch   CHKERRQ(DMCopyDisc(base,preForest));
2385f80ce2aSJacob Faibussowitsch   CHKERRQ(DMForestSetBaseDM(preForest,base));
2395f80ce2aSJacob Faibussowitsch   CHKERRQ(DMForestSetMinimumRefinement(preForest,0));
2405f80ce2aSJacob Faibussowitsch   CHKERRQ(DMForestSetInitialRefinement(preForest,1));
2415f80ce2aSJacob Faibussowitsch   CHKERRQ(DMSetFromOptions(preForest));
2425f80ce2aSJacob Faibussowitsch   CHKERRQ(DMSetUp(preForest));
2435f80ce2aSJacob Faibussowitsch   CHKERRQ(DMViewFromOptions(preForest,NULL,"-dm_pre_view"));
244c4762a1bSJed Brown 
245c4762a1bSJed Brown   /* the pre adaptivity field */
2465f80ce2aSJacob Faibussowitsch   CHKERRQ(DMCreateGlobalVector(preForest,&preVec));
2475f80ce2aSJacob Faibussowitsch   CHKERRQ(DMProjectFunction(preForest,0.,funcs,ctxs,INSERT_VALUES,preVec));
2485f80ce2aSJacob Faibussowitsch   CHKERRQ(VecViewFromOptions(preVec,NULL,"-vec_pre_view"));
249c4762a1bSJed Brown 
250c4762a1bSJed Brown   /* communicate between base and pre adaptivity forest */
251c4762a1bSJed Brown   if (transfer_from_base[0]) {
252c4762a1bSJed Brown     Vec baseVec, baseVecMapped;
253c4762a1bSJed Brown 
2545f80ce2aSJacob Faibussowitsch     CHKERRQ(DMGetGlobalVector(base,&baseVec));
2555f80ce2aSJacob Faibussowitsch     CHKERRQ(DMProjectFunction(base,0.,funcs,ctxs,INSERT_VALUES,baseVec));
2565f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscObjectSetName((PetscObject)baseVec,"Function Base"));
2575f80ce2aSJacob Faibussowitsch     CHKERRQ(VecViewFromOptions(baseVec,NULL,"-vec_base_view"));
258c4762a1bSJed Brown 
2595f80ce2aSJacob Faibussowitsch     CHKERRQ(DMGetGlobalVector(preForest,&baseVecMapped));
2605f80ce2aSJacob Faibussowitsch     CHKERRQ(DMForestTransferVecFromBase(preForest,baseVec,baseVecMapped));
2615f80ce2aSJacob Faibussowitsch     CHKERRQ(VecViewFromOptions(baseVecMapped,NULL,"-vec_map_base_view"));
262c4762a1bSJed Brown 
263c4762a1bSJed Brown     /* compare */
2645f80ce2aSJacob Faibussowitsch     CHKERRQ(VecAXPY(baseVecMapped,-1.,preVec));
2655f80ce2aSJacob Faibussowitsch     CHKERRQ(VecViewFromOptions(baseVecMapped,NULL,"-vec_map_diff_view"));
2665f80ce2aSJacob Faibussowitsch     CHKERRQ(VecNorm(baseVecMapped,NORM_2,&diff));
267c4762a1bSJed Brown 
268c4762a1bSJed Brown     /* output */
269c4762a1bSJed Brown     if (diff < tol) {
2705f80ce2aSJacob Faibussowitsch       CHKERRQ(PetscPrintf(comm,"DMForestTransferVecFromBase() passes.\n"));
271c4762a1bSJed Brown     } else {
2725f80ce2aSJacob Faibussowitsch       CHKERRQ(PetscPrintf(comm,"DMForestTransferVecFromBase() fails with error %g and tolerance %g\n",(double)diff,(double)tol));
273c4762a1bSJed Brown     }
274c4762a1bSJed Brown 
2755f80ce2aSJacob Faibussowitsch     CHKERRQ(DMRestoreGlobalVector(base,&baseVec));
2765f80ce2aSJacob Faibussowitsch     CHKERRQ(DMRestoreGlobalVector(preForest,&baseVecMapped));
277c4762a1bSJed Brown   }
278c4762a1bSJed Brown 
279c4762a1bSJed Brown   for (step = 0; step < adaptSteps; ++step) {
280c4762a1bSJed Brown 
281c4762a1bSJed Brown     if (!transfer_from_base[1]) {
2825f80ce2aSJacob Faibussowitsch       CHKERRQ(PetscObjectGetReference((PetscObject)preForest,&preCount));
283c4762a1bSJed Brown     }
284c4762a1bSJed Brown 
285c4762a1bSJed Brown     /* adapt */
2865f80ce2aSJacob Faibussowitsch     CHKERRQ(CreateAdaptivityLabel(preForest,&adaptLabel));
2875f80ce2aSJacob Faibussowitsch     CHKERRQ(DMForestTemplate(preForest,comm,&postForest));
2885f80ce2aSJacob Faibussowitsch     if (step) CHKERRQ(DMForestSetAdaptivityLabel(postForest,adaptLabel));
2895f80ce2aSJacob Faibussowitsch     CHKERRQ(DMLabelDestroy(&adaptLabel));
2905f80ce2aSJacob Faibussowitsch     CHKERRQ(DMSetUp(postForest));
2915f80ce2aSJacob Faibussowitsch     CHKERRQ(DMViewFromOptions(postForest,NULL,"-dm_post_view"));
292c4762a1bSJed Brown 
293c4762a1bSJed Brown     /* transfer */
2945f80ce2aSJacob Faibussowitsch     CHKERRQ(DMCreateGlobalVector(postForest,&postVecTransfer));
2955f80ce2aSJacob Faibussowitsch     CHKERRQ(DMForestTransferVec(preForest,preVec,postForest,postVecTransfer,PETSC_TRUE,0.0));
2965f80ce2aSJacob Faibussowitsch     CHKERRQ(VecViewFromOptions(postVecTransfer,NULL,"-vec_post_transfer_view"));
297c4762a1bSJed Brown 
298c4762a1bSJed Brown     /* the exact post adaptivity field */
2995f80ce2aSJacob Faibussowitsch     CHKERRQ(DMCreateGlobalVector(postForest,&postVecExact));
3005f80ce2aSJacob Faibussowitsch     CHKERRQ(DMProjectFunction(postForest,0.,funcs,ctxs,INSERT_VALUES,postVecExact));
3015f80ce2aSJacob Faibussowitsch     CHKERRQ(VecViewFromOptions(postVecExact,NULL,"-vec_post_exact_view"));
302c4762a1bSJed Brown 
303c4762a1bSJed Brown     /* compare */
3045f80ce2aSJacob Faibussowitsch     CHKERRQ(VecAXPY(postVecExact,-1.,postVecTransfer));
3055f80ce2aSJacob Faibussowitsch     CHKERRQ(VecViewFromOptions(postVecExact,NULL,"-vec_diff_view"));
3065f80ce2aSJacob Faibussowitsch     CHKERRQ(VecNorm(postVecExact,NORM_2,&diff));
307c4762a1bSJed Brown 
308c4762a1bSJed Brown     /* output */
309c4762a1bSJed Brown     if (diff < tol) {
3105f80ce2aSJacob Faibussowitsch       CHKERRQ(PetscPrintf(comm,"DMForestTransferVec() passes.\n"));
311c4762a1bSJed Brown     } else {
3125f80ce2aSJacob Faibussowitsch       CHKERRQ(PetscPrintf(comm,"DMForestTransferVec() fails with error %g and tolerance %g\n",(double)diff,(double)tol));
3135f80ce2aSJacob Faibussowitsch       CHKERRQ(IdentifyBadPoints(postForest, postVecExact, tol));
314c4762a1bSJed Brown     }
3155f80ce2aSJacob Faibussowitsch     CHKERRQ(VecDestroy(&postVecExact));
316c4762a1bSJed Brown 
317c4762a1bSJed Brown     /* disconnect preForest from postForest if we don't test the transfer throughout the entire refinement process */
318c4762a1bSJed Brown     if (!transfer_from_base[1]) {
3195f80ce2aSJacob Faibussowitsch       CHKERRQ(DMForestSetAdaptivityForest(postForest,NULL));
3205f80ce2aSJacob Faibussowitsch       CHKERRQ(PetscObjectGetReference((PetscObject)preForest,&postCount));
3212c71b3e2SJacob Faibussowitsch       PetscCheckFalse(postCount != preCount,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Adaptation not memory neutral: reference count increase from %d to %d",preCount,postCount);
322c4762a1bSJed Brown     }
323c4762a1bSJed Brown 
324c4762a1bSJed Brown     if (conv) {
325c4762a1bSJed Brown       DM dmConv;
326c4762a1bSJed Brown 
3275f80ce2aSJacob Faibussowitsch       CHKERRQ(DMConvert(postForest,DMPLEX,&dmConv));
3285f80ce2aSJacob Faibussowitsch       CHKERRQ(DMViewFromOptions(dmConv,NULL,"-dm_conv_view"));
3295f80ce2aSJacob Faibussowitsch       CHKERRQ(DMPlexCheckCellShape(dmConv,PETSC_TRUE,PETSC_DETERMINE));
3305f80ce2aSJacob Faibussowitsch       CHKERRQ(DMDestroy(&dmConv));
331c4762a1bSJed Brown     }
332c4762a1bSJed Brown 
3335f80ce2aSJacob Faibussowitsch     CHKERRQ(VecDestroy(&preVec));
3345f80ce2aSJacob Faibussowitsch     CHKERRQ(DMDestroy(&preForest));
335c4762a1bSJed Brown 
336c4762a1bSJed Brown     preVec    = postVecTransfer;
337c4762a1bSJed Brown     preForest = postForest;
338c4762a1bSJed Brown   }
339c4762a1bSJed Brown 
340c4762a1bSJed Brown   if (transfer_from_base[1]) {
341c4762a1bSJed Brown     Vec baseVec, baseVecMapped;
342c4762a1bSJed Brown 
343c4762a1bSJed Brown     /* communicate between base and last adapted forest */
3445f80ce2aSJacob Faibussowitsch     CHKERRQ(DMGetGlobalVector(base,&baseVec));
3455f80ce2aSJacob Faibussowitsch     CHKERRQ(DMProjectFunction(base,0.,funcs,ctxs,INSERT_VALUES,baseVec));
3465f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscObjectSetName((PetscObject)baseVec,"Function Base"));
3475f80ce2aSJacob Faibussowitsch     CHKERRQ(VecViewFromOptions(baseVec,NULL,"-vec_base_view"));
348c4762a1bSJed Brown 
3495f80ce2aSJacob Faibussowitsch     CHKERRQ(DMGetGlobalVector(preForest,&baseVecMapped));
3505f80ce2aSJacob Faibussowitsch     CHKERRQ(DMForestTransferVecFromBase(preForest,baseVec,baseVecMapped));
3515f80ce2aSJacob Faibussowitsch     CHKERRQ(VecViewFromOptions(baseVecMapped,NULL,"-vec_map_base_view"));
352c4762a1bSJed Brown 
353c4762a1bSJed Brown     /* compare */
3545f80ce2aSJacob Faibussowitsch     CHKERRQ(VecAXPY(baseVecMapped,-1.,preVec));
3555f80ce2aSJacob Faibussowitsch     CHKERRQ(VecViewFromOptions(baseVecMapped,NULL,"-vec_map_diff_view"));
3565f80ce2aSJacob Faibussowitsch     CHKERRQ(VecNorm(baseVecMapped,NORM_2,&diff));
357c4762a1bSJed Brown 
358c4762a1bSJed Brown     /* output */
359c4762a1bSJed Brown     if (diff < tol) {
3605f80ce2aSJacob Faibussowitsch       CHKERRQ(PetscPrintf(comm,"DMForestTransferVecFromBase() passes.\n"));
361c4762a1bSJed Brown     } else {
3625f80ce2aSJacob Faibussowitsch       CHKERRQ(PetscPrintf(comm,"DMForestTransferVecFromBase() fails with error %g and tolerance %g\n",(double)diff,(double)tol));
363c4762a1bSJed Brown     }
364c4762a1bSJed Brown 
3655f80ce2aSJacob Faibussowitsch     CHKERRQ(DMRestoreGlobalVector(base,&baseVec));
3665f80ce2aSJacob Faibussowitsch     CHKERRQ(DMRestoreGlobalVector(preForest,&baseVecMapped));
367c4762a1bSJed Brown   }
368c4762a1bSJed Brown 
369c4762a1bSJed Brown   /* cleanup */
3705f80ce2aSJacob Faibussowitsch   CHKERRQ(VecDestroy(&preVec));
3715f80ce2aSJacob Faibussowitsch   CHKERRQ(DMDestroy(&preForest));
3725f80ce2aSJacob Faibussowitsch   CHKERRQ(DMDestroy(&base));
373*b122ec5aSJacob Faibussowitsch   CHKERRQ(PetscFinalize());
374*b122ec5aSJacob Faibussowitsch   return 0;
375c4762a1bSJed Brown }
376c4762a1bSJed Brown 
377c4762a1bSJed Brown /*TEST
37830602db0SMatthew G. Knepley   testset:
37930602db0SMatthew G. Knepley     args: -dm_plex_simplex 0 -dm_plex_box_faces 3,3,3 -petscspace_type tensor
380c4762a1bSJed Brown 
381c4762a1bSJed Brown     test:
382c4762a1bSJed Brown       output_file: output/ex2_2d.out
383c4762a1bSJed Brown       suffix: p4est_2d
38430602db0SMatthew G. Knepley       args: -petscspace_degree 2
385c4762a1bSJed Brown       nsize: 3
386c4762a1bSJed Brown       requires: p4est !single
387c4762a1bSJed Brown 
388c4762a1bSJed Brown     test:
389c4762a1bSJed Brown       output_file: output/ex2_2d.out
390c4762a1bSJed Brown       suffix: p4est_2d_deg4
39130602db0SMatthew G. Knepley       args: -petscspace_degree 4
392c4762a1bSJed Brown       requires: p4est !single
393c4762a1bSJed Brown 
394c4762a1bSJed Brown     test:
395c4762a1bSJed Brown       output_file: output/ex2_2d.out
396c4762a1bSJed Brown       suffix: p4est_2d_deg8
39730602db0SMatthew G. Knepley       args: -petscspace_degree 8
398c4762a1bSJed Brown       requires: p4est !single
399c4762a1bSJed Brown 
400c4762a1bSJed Brown     test:
401c4762a1bSJed Brown       output_file: output/ex2_steps2.out
402c4762a1bSJed Brown       suffix: p4est_2d_deg2_steps2
40330602db0SMatthew G. Knepley       args: -petscspace_degree 2 -coords -adapt_steps 2
404c4762a1bSJed Brown       nsize: 3
405c4762a1bSJed Brown       requires: p4est !single
406c4762a1bSJed Brown 
407c4762a1bSJed Brown     test:
408c4762a1bSJed Brown       output_file: output/ex2_steps3.out
409c4762a1bSJed Brown       suffix: p4est_2d_deg3_steps3
41030602db0SMatthew G. Knepley       args: -petscspace_degree 3 -coords -adapt_steps 3 -petscdualspace_lagrange_node_type equispaced -petscdualspace_lagrange_node_endpoints 1
411c4762a1bSJed Brown       nsize: 3
412c4762a1bSJed Brown       requires: p4est !single
413c4762a1bSJed Brown 
414c4762a1bSJed Brown     test:
415c4762a1bSJed Brown       output_file: output/ex2_steps3.out
416c4762a1bSJed Brown       suffix: p4est_2d_deg3_steps3_L2_periodic
41730602db0SMatthew G. Knepley       args: -petscspace_degree 3 -petscdualspace_lagrange_continuity 0 -coords -adapt_steps 3 -dm_plex_box_bd periodic,periodic -use_bcs 0 -petscdualspace_lagrange_node_type equispaced
418c4762a1bSJed Brown       nsize: 3
419c4762a1bSJed Brown       requires: p4est !single
420c4762a1bSJed Brown 
421c4762a1bSJed Brown     test:
422c4762a1bSJed Brown       output_file: output/ex2_steps3.out
423c4762a1bSJed Brown       suffix: p4est_3d_deg2_steps3_L2_periodic
42430602db0SMatthew G. Knepley       args: -dm_plex_dim 3 -petscspace_degree 2 -petscdualspace_lagrange_continuity 0 -coords -adapt_steps 3 -dm_plex_box_bd periodic,periodic,periodic -use_bcs 0
425c4762a1bSJed Brown       nsize: 3
426c4762a1bSJed Brown       requires: p4est !single
427c4762a1bSJed Brown 
428c4762a1bSJed Brown     test:
429c4762a1bSJed Brown       output_file: output/ex2_steps2.out
430c4762a1bSJed Brown       suffix: p4est_3d_deg2_steps2
43130602db0SMatthew G. Knepley       args: -dm_plex_dim 3 -petscspace_degree 2 -coords -adapt_steps 2
432c4762a1bSJed Brown       nsize: 3
433c4762a1bSJed Brown       requires: p4est !single
434c4762a1bSJed Brown 
435c4762a1bSJed Brown     test:
436c4762a1bSJed Brown       output_file: output/ex2_steps3.out
437c4762a1bSJed Brown       suffix: p4est_3d_deg3_steps3
43830602db0SMatthew G. Knepley       args: -dm_plex_dim 3 -petscspace_degree 3 -coords -adapt_steps 3 -petscdualspace_lagrange_node_type equispaced -petscdualspace_lagrange_node_endpoints 1
439c4762a1bSJed Brown       nsize: 3
440c4762a1bSJed Brown       requires: p4est !single
441c4762a1bSJed Brown 
442c4762a1bSJed Brown     test:
443c4762a1bSJed Brown       output_file: output/ex2_3d.out
444c4762a1bSJed Brown       suffix: p4est_3d
44530602db0SMatthew G. Knepley       args: -dm_plex_dim 3 -petscspace_degree 1
446c4762a1bSJed Brown       nsize: 3
447c4762a1bSJed Brown       requires: p4est !single
448c4762a1bSJed Brown 
449c4762a1bSJed Brown     test:
450c4762a1bSJed Brown       output_file: output/ex2_3d.out
451c4762a1bSJed Brown       suffix: p4est_3d_deg3
45230602db0SMatthew G. Knepley       args: -dm_plex_dim 3 -petscspace_degree 3
453c4762a1bSJed Brown       nsize: 3
454c4762a1bSJed Brown       requires: p4est !single
455c4762a1bSJed Brown 
456c4762a1bSJed Brown     test:
457c4762a1bSJed Brown       output_file: output/ex2_2d.out
458c4762a1bSJed Brown       suffix: p4est_2d_deg2_coords
45930602db0SMatthew G. Knepley       args: -petscspace_degree 2 -coords
460c4762a1bSJed Brown       nsize: 3
461c4762a1bSJed Brown       requires: p4est !single
462c4762a1bSJed Brown 
463c4762a1bSJed Brown     test:
464c4762a1bSJed Brown       output_file: output/ex2_3d.out
465c4762a1bSJed Brown       suffix: p4est_3d_deg2_coords
46630602db0SMatthew G. Knepley       args: -dm_plex_dim 3 -petscspace_degree 2 -coords
467c4762a1bSJed Brown       nsize: 3
468c4762a1bSJed Brown       requires: p4est !single
469c4762a1bSJed Brown 
470c4762a1bSJed Brown     test:
471c4762a1bSJed Brown       suffix: p4est_3d_nans
47230602db0SMatthew G. Knepley       args: -dm_plex_dim 3 -dm_forest_partition_overlap 1 -test_convert -petscspace_degree 1
473c4762a1bSJed Brown       nsize: 2
474c4762a1bSJed Brown       requires: p4est !single
475c4762a1bSJed Brown 
476c4762a1bSJed Brown     test:
477c4762a1bSJed Brown       TODO: not broken, but the 3D case below is broken, so I do not trust this one
478c4762a1bSJed Brown       output_file: output/ex2_steps2.out
479c4762a1bSJed Brown       suffix: p4est_2d_tfb_distributed_nc
480e600fa54SMatthew G. Knepley       args: -petscspace_degree 3 -dm_forest_maximum_refinement 2 -dm_p4est_refine_pattern hash -use_bcs 0 -coords -adapt_steps 2 -petscpartitioner_type shell -petscpartitioner_shell_random
481c4762a1bSJed Brown       nsize: 3
482c4762a1bSJed Brown       requires: p4est !single
483c4762a1bSJed Brown 
484c4762a1bSJed Brown     test:
485c4762a1bSJed Brown       TODO: broken
486c4762a1bSJed Brown       output_file: output/ex2_steps2.out
487c4762a1bSJed Brown       suffix: p4est_3d_tfb_distributed_nc
488e600fa54SMatthew G. Knepley       args: -dm_plex_dim 3 -petscspace_degree 2 -dm_forest_maximum_refinement 2 -dm_p4est_refine_pattern hash -use_bcs 0 -coords -adapt_steps 2 -petscpartitioner_type shell -petscpartitioner_shell_random
489c4762a1bSJed Brown       nsize: 3
490c4762a1bSJed Brown       requires: p4est !single
491c4762a1bSJed Brown 
49230602db0SMatthew G. Knepley   testset:
493012bc364SMatthew G. Knepley     args: -petscspace_type tensor -dm_coord_space 0 -dm_plex_transform_type refine_tobox
49430602db0SMatthew G. Knepley 
495c4762a1bSJed Brown     test:
496c4762a1bSJed Brown       TODO: broken
497c4762a1bSJed Brown       output_file: output/ex2_3d.out
498c4762a1bSJed Brown       suffix: p4est_3d_transfer_fails
499012bc364SMatthew G. Knepley       args: -petscspace_degree 1 -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh -adapt_steps 1 -dm_forest_initial_refinement 1 -use_bcs 0 -dm_refine
500c4762a1bSJed Brown       requires: p4est !single
501c4762a1bSJed Brown 
502c4762a1bSJed Brown     test:
503c4762a1bSJed Brown       TODO: broken
504c4762a1bSJed Brown       output_file: output/ex2_steps2_notfb.out
505c4762a1bSJed Brown       suffix: p4est_3d_transfer_fails_2
506012bc364SMatthew G. Knepley       args: -petscspace_degree 1 -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh -adapt_steps 2 -dm_forest_initial_refinement 0 -transfer_from_base 0 -use_bcs 0 -dm_refine
507c4762a1bSJed Brown       requires: p4est !single
508c4762a1bSJed Brown 
509c4762a1bSJed Brown     test:
510c4762a1bSJed Brown       output_file: output/ex2_steps2.out
511c4762a1bSJed Brown       suffix: p4est_3d_multi_transfer_s2t
512012bc364SMatthew G. Knepley       args: -petscspace_degree 3 -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh -adapt_steps 2 -dm_forest_initial_refinement 1 -petscdualspace_lagrange_continuity 0 -use_bcs 0 -dm_refine 1
513012bc364SMatthew G. Knepley       requires: p4est !single
514c4762a1bSJed Brown 
515c4762a1bSJed Brown     test:
516c4762a1bSJed Brown       output_file: output/ex2_steps2.out
517c4762a1bSJed Brown       suffix: p4est_3d_coords_transfer_s2t
518012bc364SMatthew G. Knepley       args: -petscspace_degree 3 -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh -adapt_steps 2 -dm_forest_initial_refinement 1 -petscdualspace_lagrange_continuity 0 -coords -use_bcs 0 -dm_refine 1
519012bc364SMatthew G. Knepley       requires: p4est !single
52030602db0SMatthew G. Knepley 
52130602db0SMatthew G. Knepley   testset:
52230602db0SMatthew G. Knepley     args: -dm_plex_simplex 0 -dm_plex_box_faces 3,3,3
52330602db0SMatthew G. Knepley 
52430602db0SMatthew G. Knepley     test:
52530602db0SMatthew G. Knepley       output_file: output/ex2_2d_fv.out
52630602db0SMatthew G. Knepley       suffix: p4est_2d_fv
52730602db0SMatthew G. Knepley       args: -transfer_from_base 0 -use_fv -linear -dm_forest_partition_overlap 1
52830602db0SMatthew G. Knepley       nsize: 3
52930602db0SMatthew G. Knepley       requires: p4est !single
53030602db0SMatthew G. Knepley 
53130602db0SMatthew G. Knepley     test:
53230602db0SMatthew G. Knepley       TODO: broken (codimension adjacency)
53330602db0SMatthew G. Knepley       output_file: output/ex2_2d_fv.out
53430602db0SMatthew G. Knepley       suffix: p4est_2d_fv_adjcodim
53530602db0SMatthew G. Knepley       args: -transfer_from_base 0 -use_fv -linear -dm_forest_partition_overlap 1 -dm_forest_adjacency_codimension 1
53630602db0SMatthew G. Knepley       nsize: 2
53730602db0SMatthew G. Knepley       requires: p4est !single
53830602db0SMatthew G. Knepley 
53930602db0SMatthew G. Knepley     test:
54030602db0SMatthew G. Knepley       TODO: broken (dimension adjacency)
54130602db0SMatthew G. Knepley       output_file: output/ex2_2d_fv.out
54230602db0SMatthew G. Knepley       suffix: p4est_2d_fv_adjdim
54330602db0SMatthew G. Knepley       args: -transfer_from_base 0 -use_fv -linear -dm_forest_partition_overlap 1 -dm_forest_adjacency_dimension 1
54430602db0SMatthew G. Knepley       nsize: 2
54530602db0SMatthew G. Knepley       requires: p4est !single
54630602db0SMatthew G. Knepley 
54730602db0SMatthew G. Knepley     test:
54830602db0SMatthew G. Knepley       output_file: output/ex2_2d_fv.out
54930602db0SMatthew G. Knepley       suffix: p4est_2d_fv_zerocells
55030602db0SMatthew G. Knepley       args: -transfer_from_base 0 -use_fv -linear -dm_forest_partition_overlap 1
55130602db0SMatthew G. Knepley       nsize: 10
55230602db0SMatthew G. Knepley       requires: p4est !single
55330602db0SMatthew G. Knepley 
55430602db0SMatthew G. Knepley     test:
55530602db0SMatthew G. Knepley       output_file: output/ex2_3d_fv.out
55630602db0SMatthew G. Knepley       suffix: p4est_3d_fv
55730602db0SMatthew G. Knepley       args: -dm_plex_dim 3 -transfer_from_base 0 -use_fv -linear -dm_forest_partition_overlap 1
55830602db0SMatthew G. Knepley       nsize: 3
559c4762a1bSJed Brown       requires: p4est !single
560c4762a1bSJed Brown 
561c4762a1bSJed Brown TEST*/
562