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 PetscErrorCode ierr; 12c4762a1bSJed Brown 13c4762a1bSJed Brown PetscFunctionBegin; 14c4762a1bSJed Brown ierr = DMCreateLabel(dm, "identity");CHKERRQ(ierr); 15c4762a1bSJed Brown ierr = DMPlexGetChart(dm, &pStart, &pEnd);CHKERRQ(ierr); 16c4762a1bSJed Brown for (p = pStart; p < pEnd; p++) {ierr = DMSetLabelValue(dm, "identity", p, p);CHKERRQ(ierr);} 17c4762a1bSJed Brown PetscFunctionReturn(0); 18c4762a1bSJed Brown } 19c4762a1bSJed Brown 20c4762a1bSJed Brown static PetscErrorCode CreateAdaptivityLabel(DM forest,DMLabel *adaptLabel) 21c4762a1bSJed Brown { 22c4762a1bSJed Brown DMLabel identLabel; 23c4762a1bSJed Brown PetscInt cStart, cEnd, c; 24c4762a1bSJed Brown PetscErrorCode ierr; 25c4762a1bSJed Brown 26c4762a1bSJed Brown PetscFunctionBegin; 27c4762a1bSJed Brown ierr = DMLabelCreate(PETSC_COMM_SELF,"adapt",adaptLabel);CHKERRQ(ierr); 28c4762a1bSJed Brown ierr = DMLabelSetDefaultValue(*adaptLabel,DM_ADAPT_COARSEN);CHKERRQ(ierr); 29c4762a1bSJed Brown ierr = DMGetLabel(forest,"identity",&identLabel);CHKERRQ(ierr); 30c4762a1bSJed Brown ierr = DMForestGetCellChart(forest,&cStart,&cEnd);CHKERRQ(ierr); 31c4762a1bSJed Brown for (c = cStart; c < cEnd; c++) { 32c4762a1bSJed Brown PetscInt basePoint; 33c4762a1bSJed Brown 34c4762a1bSJed Brown ierr = DMLabelGetValue(identLabel,c,&basePoint);CHKERRQ(ierr); 35c4762a1bSJed Brown if (!basePoint) {ierr = DMLabelSetValue(*adaptLabel,c,DM_ADAPT_REFINE);CHKERRQ(ierr);} 36c4762a1bSJed Brown } 37c4762a1bSJed Brown PetscFunctionReturn(0); 38c4762a1bSJed Brown } 39c4762a1bSJed Brown 40c4762a1bSJed Brown static PetscErrorCode LinearFunction(PetscInt dim,PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar u[], void *ctx) 41c4762a1bSJed Brown { 42c4762a1bSJed Brown PetscFunctionBeginUser; 43c4762a1bSJed Brown u[0] = (x[0] * 2.0 + 1.) + (x[1] * 20.0 + 10.) + ((dim == 3) ? (x[2] * 200.0 + 100.) : 0.); 44c4762a1bSJed Brown PetscFunctionReturn(0); 45c4762a1bSJed Brown } 46c4762a1bSJed Brown 47c4762a1bSJed Brown static PetscErrorCode MultiaffineFunction(PetscInt dim,PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar u[], void *ctx) 48c4762a1bSJed Brown { 49c4762a1bSJed Brown PetscFunctionBeginUser; 50c4762a1bSJed Brown u[0] = (x[0] * 1.0 + 2.0) * (x[1] * 3.0 - 4.0) * ((dim == 3) ? (x[2] * 5.0 + 6.0) : 1.); 51c4762a1bSJed Brown PetscFunctionReturn(0); 52c4762a1bSJed Brown } 53c4762a1bSJed Brown 54c4762a1bSJed Brown static PetscErrorCode CoordsFunction(PetscInt dim,PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar u[], void *ctx) 55c4762a1bSJed Brown { 56c4762a1bSJed Brown PetscInt f; 57c4762a1bSJed Brown 58c4762a1bSJed Brown PetscFunctionBeginUser; 59c4762a1bSJed Brown for (f=0;f<Nf;f++) u[f] = x[f]; 60c4762a1bSJed Brown PetscFunctionReturn(0); 61c4762a1bSJed Brown } 62c4762a1bSJed Brown 63c4762a1bSJed Brown typedef struct _bc_func_ctx 64c4762a1bSJed Brown { 65c4762a1bSJed Brown PetscErrorCode (*func) (PetscInt,PetscReal,const PetscReal [], PetscInt, PetscScalar [], void *); 66c4762a1bSJed Brown PetscInt dim; 67c4762a1bSJed Brown PetscInt Nf; 68c4762a1bSJed Brown void *ctx; 69c4762a1bSJed Brown } 70c4762a1bSJed Brown bc_func_ctx; 71c4762a1bSJed Brown 72c4762a1bSJed Brown static PetscErrorCode bc_func_fv (PetscReal time, const PetscReal *c, const PetscReal *n, const PetscScalar *xI, PetscScalar *xG, void *ctx) 73c4762a1bSJed Brown { 74c4762a1bSJed Brown bc_func_ctx *bcCtx; 75c4762a1bSJed Brown PetscErrorCode ierr; 76c4762a1bSJed Brown 77c4762a1bSJed Brown PetscFunctionBegin; 78c4762a1bSJed Brown bcCtx = (bc_func_ctx *) ctx; 79c4762a1bSJed Brown ierr = (bcCtx->func)(bcCtx->dim,time,c,bcCtx->Nf,xG,bcCtx->ctx);CHKERRQ(ierr); 80c4762a1bSJed Brown PetscFunctionReturn(0); 81c4762a1bSJed Brown } 82c4762a1bSJed Brown 83c4762a1bSJed Brown static PetscErrorCode IdentifyBadPoints (DM dm, Vec vec, PetscReal tol) 84c4762a1bSJed Brown { 85c4762a1bSJed Brown DM dmplex; 86c4762a1bSJed Brown PetscInt p, pStart, pEnd, maxDof; 87c4762a1bSJed Brown Vec vecLocal; 88c4762a1bSJed Brown DMLabel depthLabel; 89c4762a1bSJed Brown PetscSection section; 90c4762a1bSJed Brown PetscErrorCode ierr; 91c4762a1bSJed Brown 92c4762a1bSJed Brown PetscFunctionBegin; 93c4762a1bSJed Brown ierr = DMCreateLocalVector(dm, &vecLocal);CHKERRQ(ierr); 94c4762a1bSJed Brown ierr = DMGlobalToLocalBegin(dm, vec, INSERT_VALUES, vecLocal);CHKERRQ(ierr); 95c4762a1bSJed Brown ierr = DMGlobalToLocalEnd(dm, vec, INSERT_VALUES, vecLocal);CHKERRQ(ierr); 96c4762a1bSJed Brown ierr = DMConvert(dm ,DMPLEX, &dmplex);CHKERRQ(ierr); 97c4762a1bSJed Brown ierr = DMPlexGetChart(dmplex, &pStart, &pEnd);CHKERRQ(ierr); 98c4762a1bSJed Brown ierr = DMPlexGetDepthLabel(dmplex, &depthLabel);CHKERRQ(ierr); 99c4762a1bSJed Brown ierr = DMGetLocalSection(dmplex, §ion);CHKERRQ(ierr); 100c4762a1bSJed Brown ierr = PetscSectionGetMaxDof(section, &maxDof);CHKERRQ(ierr); 101c4762a1bSJed Brown for (p = pStart; p < pEnd; p++) { 102c4762a1bSJed Brown PetscInt s, c, cSize, parent, childID, numChildren; 103c4762a1bSJed Brown PetscInt cl, closureSize, *closure = NULL; 104c4762a1bSJed Brown PetscScalar *values = NULL; 105c4762a1bSJed Brown PetscBool bad = PETSC_FALSE; 106c4762a1bSJed Brown 107c4762a1bSJed Brown ierr = VecGetValuesSection(vecLocal, section, p, &values);CHKERRQ(ierr); 108c4762a1bSJed Brown ierr = PetscSectionGetDof(section, p, &cSize);CHKERRQ(ierr); 109c4762a1bSJed Brown for (c = 0; c < cSize; c++) { 110*2f613bf5SBarry Smith PetscReal absDiff = PetscAbsScalar(values[c]); 111c4762a1bSJed Brown if (absDiff > tol) {bad = PETSC_TRUE; break;} 112c4762a1bSJed Brown } 113c4762a1bSJed Brown if (!bad) continue; 114c4762a1bSJed Brown ierr = PetscPrintf(PETSC_COMM_SELF, "Bad point %D\n", p);CHKERRQ(ierr); 115c4762a1bSJed Brown ierr = DMLabelGetValue(depthLabel, p, &s);CHKERRQ(ierr); 116c4762a1bSJed Brown ierr = PetscPrintf(PETSC_COMM_SELF, " Depth %D\n", s);CHKERRQ(ierr); 117c4762a1bSJed Brown ierr = DMPlexGetTransitiveClosure(dmplex, p, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr); 118c4762a1bSJed Brown for (cl = 0; cl < closureSize; cl++) { 119c4762a1bSJed Brown PetscInt cp = closure[2 * cl]; 120c4762a1bSJed Brown ierr = DMPlexGetTreeParent(dmplex, cp, &parent, &childID);CHKERRQ(ierr); 121c4762a1bSJed Brown if (parent != cp) { 122c4762a1bSJed Brown ierr = PetscPrintf(PETSC_COMM_SELF, " Closure point %D (%D) child of %D (ID %D)\n", cl, cp, parent, childID);CHKERRQ(ierr); 123c4762a1bSJed Brown } 124c4762a1bSJed Brown ierr = DMPlexGetTreeChildren(dmplex, cp, &numChildren, NULL);CHKERRQ(ierr); 125c4762a1bSJed Brown if (numChildren) { 126c4762a1bSJed Brown ierr = PetscPrintf(PETSC_COMM_SELF, " Closure point %D (%D) is parent\n", cl, cp);CHKERRQ(ierr); 127c4762a1bSJed Brown } 128c4762a1bSJed Brown } 129c4762a1bSJed Brown ierr = DMPlexRestoreTransitiveClosure(dmplex, p, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr); 130c4762a1bSJed Brown for (c = 0; c < cSize; c++) { 131*2f613bf5SBarry Smith PetscReal absDiff = PetscAbsScalar(values[c]); 132c4762a1bSJed Brown if (absDiff > tol) { 133c4762a1bSJed Brown ierr = PetscPrintf(PETSC_COMM_SELF, " Bad dof %D\n", c);CHKERRQ(ierr); 134c4762a1bSJed Brown } 135c4762a1bSJed Brown } 136c4762a1bSJed Brown } 137c4762a1bSJed Brown ierr = DMDestroy(&dmplex);CHKERRQ(ierr); 138c4762a1bSJed Brown ierr = VecDestroy(&vecLocal);CHKERRQ(ierr); 139c4762a1bSJed Brown PetscFunctionReturn(0); 140c4762a1bSJed Brown } 141c4762a1bSJed Brown 142c4762a1bSJed Brown int main(int argc, char **argv) 143c4762a1bSJed Brown { 144c4762a1bSJed Brown MPI_Comm comm; 145c4762a1bSJed Brown DM base, preForest, postForest; 14630602db0SMatthew G. Knepley PetscInt dim, Nf = 1; 147c4762a1bSJed Brown PetscInt step, adaptSteps = 1; 148c4762a1bSJed Brown PetscInt preCount, postCount; 149c4762a1bSJed Brown Vec preVec, postVecTransfer, postVecExact; 150c4762a1bSJed Brown PetscErrorCode (*funcs[1]) (PetscInt,PetscReal,const PetscReal [],PetscInt,PetscScalar [], void *) = {MultiaffineFunction}; 151c4762a1bSJed Brown void *ctxs[1] = {NULL}; 152c4762a1bSJed Brown PetscReal diff, tol = PETSC_SMALL; 153c4762a1bSJed Brown PetscBool linear = PETSC_FALSE; 154c4762a1bSJed Brown PetscBool coords = PETSC_FALSE; 155c4762a1bSJed Brown PetscBool useFV = PETSC_FALSE; 156c4762a1bSJed Brown PetscBool conv = PETSC_FALSE; 157c4762a1bSJed Brown PetscBool transfer_from_base[2] = {PETSC_TRUE,PETSC_FALSE}; 158c4762a1bSJed Brown PetscBool use_bcs = PETSC_TRUE; 159c4762a1bSJed Brown bc_func_ctx bcCtx; 160c4762a1bSJed Brown DMLabel adaptLabel; 161c4762a1bSJed Brown PetscErrorCode ierr; 162c4762a1bSJed Brown 163c4762a1bSJed Brown ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr; 164c4762a1bSJed Brown comm = PETSC_COMM_WORLD; 165c4762a1bSJed Brown ierr = PetscOptionsBegin(comm, "", "DMForestTransferVec() Test Options", "DMFOREST");CHKERRQ(ierr); 166c4762a1bSJed Brown ierr = PetscOptionsBool("-linear","Transfer a simple linear function", "ex2.c", linear, &linear, NULL);CHKERRQ(ierr); 167c4762a1bSJed Brown ierr = PetscOptionsBool("-coords","Transfer a simple coordinate function", "ex2.c", coords, &coords, NULL);CHKERRQ(ierr); 168c4762a1bSJed Brown ierr = PetscOptionsBool("-use_fv","Use a finite volume approximation", "ex2.c", useFV, &useFV, NULL);CHKERRQ(ierr); 169c4762a1bSJed Brown ierr = PetscOptionsBool("-test_convert","Test conversion to DMPLEX",NULL,conv,&conv,NULL);CHKERRQ(ierr); 170c4762a1bSJed Brown ierr = PetscOptionsBool("-transfer_from_base","Transfer a vector from base DM to DMForest", "ex2.c", transfer_from_base[0], &transfer_from_base[0], NULL);CHKERRQ(ierr); 171c4762a1bSJed Brown transfer_from_base[1] = transfer_from_base[0]; 172c4762a1bSJed Brown ierr = 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);CHKERRQ(ierr); 173c4762a1bSJed Brown ierr = PetscOptionsBool("-use_bcs","Use dirichlet boundary conditions", "ex2.c", use_bcs, &use_bcs, NULL);CHKERRQ(ierr); 174c4762a1bSJed Brown ierr = PetscOptionsBoundedInt("-adapt_steps","Number of adaptivity steps", "ex2.c", adaptSteps, &adaptSteps, NULL,0);CHKERRQ(ierr); 175c4762a1bSJed Brown ierr = PetscOptionsEnd();CHKERRQ(ierr); 176c4762a1bSJed Brown 177c4762a1bSJed Brown tol = PetscMax(1.e-10,tol); /* XXX fix for quadruple precision -> why do I need to do this? */ 178c4762a1bSJed Brown 17930602db0SMatthew G. Knepley /* the base mesh */ 18030602db0SMatthew G. Knepley ierr = DMCreate(comm, &base);CHKERRQ(ierr); 18130602db0SMatthew G. Knepley ierr = DMSetType(base, DMPLEX);CHKERRQ(ierr); 18230602db0SMatthew G. Knepley ierr = DMSetFromOptions(base);CHKERRQ(ierr); 18330602db0SMatthew G. Knepley 18430602db0SMatthew G. Knepley ierr = AddIdentityLabel(base);CHKERRQ(ierr); 18530602db0SMatthew G. Knepley ierr = DMGetDimension(base, &dim);CHKERRQ(ierr); 186c4762a1bSJed Brown 187c4762a1bSJed Brown if (linear) { 188c4762a1bSJed Brown funcs[0] = LinearFunction; 189c4762a1bSJed Brown } 190c4762a1bSJed Brown if (coords) { 191c4762a1bSJed Brown funcs[0] = CoordsFunction; 192c4762a1bSJed Brown Nf = dim; 193c4762a1bSJed Brown } 194c4762a1bSJed Brown 195c4762a1bSJed Brown bcCtx.func = funcs[0]; 196c4762a1bSJed Brown bcCtx.dim = dim; 197c4762a1bSJed Brown bcCtx.Nf = Nf; 198c4762a1bSJed Brown bcCtx.ctx = NULL; 199c4762a1bSJed Brown 200c4762a1bSJed Brown if (useFV) { 201c4762a1bSJed Brown PetscFV fv; 202c4762a1bSJed Brown PetscLimiter limiter; 203c4762a1bSJed Brown DM baseFV; 204c4762a1bSJed Brown 205c4762a1bSJed Brown ierr = DMPlexConstructGhostCells(base,NULL,NULL,&baseFV);CHKERRQ(ierr); 20654fcfd0cSMatthew G. Knepley ierr = DMViewFromOptions(baseFV, NULL, "-fv_dm_view");CHKERRQ(ierr); 207c4762a1bSJed Brown ierr = DMDestroy(&base);CHKERRQ(ierr); 208c4762a1bSJed Brown base = baseFV; 209c4762a1bSJed Brown ierr = PetscFVCreate(comm, &fv);CHKERRQ(ierr); 210c4762a1bSJed Brown ierr = PetscFVSetSpatialDimension(fv,dim);CHKERRQ(ierr); 211c4762a1bSJed Brown ierr = PetscFVSetType(fv,PETSCFVLEASTSQUARES);CHKERRQ(ierr); 212c4762a1bSJed Brown ierr = PetscFVSetNumComponents(fv,Nf);CHKERRQ(ierr); 213c4762a1bSJed Brown ierr = PetscLimiterCreate(comm,&limiter);CHKERRQ(ierr); 214c4762a1bSJed Brown ierr = PetscLimiterSetType(limiter,PETSCLIMITERNONE);CHKERRQ(ierr); 215c4762a1bSJed Brown ierr = PetscFVSetLimiter(fv,limiter);CHKERRQ(ierr); 216c4762a1bSJed Brown ierr = PetscLimiterDestroy(&limiter);CHKERRQ(ierr); 217c4762a1bSJed Brown ierr = PetscFVSetFromOptions(fv);CHKERRQ(ierr); 218c4762a1bSJed Brown ierr = DMSetField(base,0,NULL,(PetscObject)fv);CHKERRQ(ierr); 219c4762a1bSJed Brown ierr = PetscFVDestroy(&fv);CHKERRQ(ierr); 220c4762a1bSJed Brown } else { 221c4762a1bSJed Brown PetscFE fe; 222c4762a1bSJed Brown 223c4762a1bSJed Brown ierr = PetscFECreateDefault(comm,dim,Nf,PETSC_FALSE,NULL,PETSC_DEFAULT,&fe);CHKERRQ(ierr); 224c4762a1bSJed Brown ierr = DMSetField(base,0,NULL,(PetscObject)fe);CHKERRQ(ierr); 225c4762a1bSJed Brown ierr = PetscFEDestroy(&fe);CHKERRQ(ierr); 226c4762a1bSJed Brown } 227c4762a1bSJed Brown ierr = DMCreateDS(base);CHKERRQ(ierr); 228c4762a1bSJed Brown 229c4762a1bSJed Brown if (use_bcs) { 230c4762a1bSJed Brown PetscInt ids[] = {1, 2, 3, 4, 5, 6}; 23145480ffeSMatthew G. Knepley DMLabel label; 232c4762a1bSJed Brown 23345480ffeSMatthew G. Knepley ierr = DMGetLabel(base, "marker", &label);CHKERRQ(ierr); 23445480ffeSMatthew G. Knepley ierr = 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);CHKERRQ(ierr); 235c4762a1bSJed Brown } 236c4762a1bSJed Brown ierr = DMViewFromOptions(base,NULL,"-dm_base_view");CHKERRQ(ierr); 237c4762a1bSJed Brown 238c4762a1bSJed Brown /* the pre adaptivity forest */ 239c4762a1bSJed Brown ierr = DMCreate(comm,&preForest);CHKERRQ(ierr); 240c4762a1bSJed Brown ierr = DMSetType(preForest,(dim == 2) ? DMP4EST : DMP8EST);CHKERRQ(ierr); 241c4762a1bSJed Brown ierr = DMCopyDisc(base,preForest);CHKERRQ(ierr); 242c4762a1bSJed Brown ierr = DMForestSetBaseDM(preForest,base);CHKERRQ(ierr); 243c4762a1bSJed Brown ierr = DMForestSetMinimumRefinement(preForest,0);CHKERRQ(ierr); 244c4762a1bSJed Brown ierr = DMForestSetInitialRefinement(preForest,1);CHKERRQ(ierr); 245c4762a1bSJed Brown ierr = DMSetFromOptions(preForest);CHKERRQ(ierr); 246c4762a1bSJed Brown ierr = DMSetUp(preForest);CHKERRQ(ierr); 247c4762a1bSJed Brown ierr = DMViewFromOptions(preForest,NULL,"-dm_pre_view");CHKERRQ(ierr); 248c4762a1bSJed Brown 249c4762a1bSJed Brown /* the pre adaptivity field */ 250c4762a1bSJed Brown ierr = DMCreateGlobalVector(preForest,&preVec);CHKERRQ(ierr); 251c4762a1bSJed Brown ierr = DMProjectFunction(preForest,0.,funcs,ctxs,INSERT_VALUES,preVec);CHKERRQ(ierr); 252c4762a1bSJed Brown ierr = VecViewFromOptions(preVec,NULL,"-vec_pre_view");CHKERRQ(ierr); 253c4762a1bSJed Brown 254c4762a1bSJed Brown /* communicate between base and pre adaptivity forest */ 255c4762a1bSJed Brown if (transfer_from_base[0]) { 256c4762a1bSJed Brown Vec baseVec, baseVecMapped; 257c4762a1bSJed Brown 258c4762a1bSJed Brown ierr = DMGetGlobalVector(base,&baseVec);CHKERRQ(ierr); 259c4762a1bSJed Brown ierr = DMProjectFunction(base,0.,funcs,ctxs,INSERT_VALUES,baseVec);CHKERRQ(ierr); 260c4762a1bSJed Brown ierr = PetscObjectSetName((PetscObject)baseVec,"Function Base");CHKERRQ(ierr); 261c4762a1bSJed Brown ierr = VecViewFromOptions(baseVec,NULL,"-vec_base_view");CHKERRQ(ierr); 262c4762a1bSJed Brown 263c4762a1bSJed Brown ierr = DMGetGlobalVector(preForest,&baseVecMapped);CHKERRQ(ierr); 264c4762a1bSJed Brown ierr = DMForestTransferVecFromBase(preForest,baseVec,baseVecMapped);CHKERRQ(ierr); 265c4762a1bSJed Brown ierr = VecViewFromOptions(baseVecMapped,NULL,"-vec_map_base_view");CHKERRQ(ierr); 266c4762a1bSJed Brown 267c4762a1bSJed Brown /* compare */ 268c4762a1bSJed Brown ierr = VecAXPY(baseVecMapped,-1.,preVec);CHKERRQ(ierr); 269c4762a1bSJed Brown ierr = VecViewFromOptions(baseVecMapped,NULL,"-vec_map_diff_view");CHKERRQ(ierr); 270c4762a1bSJed Brown ierr = VecNorm(baseVecMapped,NORM_2,&diff);CHKERRQ(ierr); 271c4762a1bSJed Brown 272c4762a1bSJed Brown /* output */ 273c4762a1bSJed Brown if (diff < tol) { 274c4762a1bSJed Brown ierr = PetscPrintf(comm,"DMForestTransferVecFromBase() passes.\n");CHKERRQ(ierr); 275c4762a1bSJed Brown } else { 276c4762a1bSJed Brown ierr = PetscPrintf(comm,"DMForestTransferVecFromBase() fails with error %g and tolerance %g\n",(double)diff,(double)tol);CHKERRQ(ierr); 277c4762a1bSJed Brown } 278c4762a1bSJed Brown 279c4762a1bSJed Brown ierr = DMRestoreGlobalVector(base,&baseVec);CHKERRQ(ierr); 280c4762a1bSJed Brown ierr = DMRestoreGlobalVector(preForest,&baseVecMapped);CHKERRQ(ierr); 281c4762a1bSJed Brown } 282c4762a1bSJed Brown 283c4762a1bSJed Brown for (step = 0; step < adaptSteps; ++step) { 284c4762a1bSJed Brown 285c4762a1bSJed Brown if (!transfer_from_base[1]) { 286c4762a1bSJed Brown ierr = PetscObjectGetReference((PetscObject)preForest,&preCount);CHKERRQ(ierr); 287c4762a1bSJed Brown } 288c4762a1bSJed Brown 289c4762a1bSJed Brown /* adapt */ 290c4762a1bSJed Brown ierr = CreateAdaptivityLabel(preForest,&adaptLabel);CHKERRQ(ierr); 291c4762a1bSJed Brown ierr = DMForestTemplate(preForest,comm,&postForest);CHKERRQ(ierr); 292c4762a1bSJed Brown if (step) { ierr = DMForestSetAdaptivityLabel(postForest,adaptLabel);CHKERRQ(ierr); } 293c4762a1bSJed Brown ierr = DMLabelDestroy(&adaptLabel);CHKERRQ(ierr); 294c4762a1bSJed Brown ierr = DMSetUp(postForest);CHKERRQ(ierr); 295c4762a1bSJed Brown ierr = DMViewFromOptions(postForest,NULL,"-dm_post_view");CHKERRQ(ierr); 296c4762a1bSJed Brown 297c4762a1bSJed Brown /* transfer */ 298c4762a1bSJed Brown ierr = DMCreateGlobalVector(postForest,&postVecTransfer);CHKERRQ(ierr); 299c4762a1bSJed Brown ierr = DMForestTransferVec(preForest,preVec,postForest,postVecTransfer,PETSC_TRUE,0.0);CHKERRQ(ierr); 300c4762a1bSJed Brown ierr = VecViewFromOptions(postVecTransfer,NULL,"-vec_post_transfer_view");CHKERRQ(ierr); 301c4762a1bSJed Brown 302c4762a1bSJed Brown /* the exact post adaptivity field */ 303c4762a1bSJed Brown ierr = DMCreateGlobalVector(postForest,&postVecExact);CHKERRQ(ierr); 304c4762a1bSJed Brown ierr = DMProjectFunction(postForest,0.,funcs,ctxs,INSERT_VALUES,postVecExact);CHKERRQ(ierr); 305c4762a1bSJed Brown ierr = VecViewFromOptions(postVecExact,NULL,"-vec_post_exact_view");CHKERRQ(ierr); 306c4762a1bSJed Brown 307c4762a1bSJed Brown /* compare */ 308c4762a1bSJed Brown ierr = VecAXPY(postVecExact,-1.,postVecTransfer);CHKERRQ(ierr); 309c4762a1bSJed Brown ierr = VecViewFromOptions(postVecExact,NULL,"-vec_diff_view");CHKERRQ(ierr); 310c4762a1bSJed Brown ierr = VecNorm(postVecExact,NORM_2,&diff);CHKERRQ(ierr); 311c4762a1bSJed Brown 312c4762a1bSJed Brown /* output */ 313c4762a1bSJed Brown if (diff < tol) { 314c4762a1bSJed Brown ierr = PetscPrintf(comm,"DMForestTransferVec() passes.\n");CHKERRQ(ierr); 315c4762a1bSJed Brown } else { 316c4762a1bSJed Brown ierr = PetscPrintf(comm,"DMForestTransferVec() fails with error %g and tolerance %g\n",(double)diff,(double)tol);CHKERRQ(ierr); 317c4762a1bSJed Brown ierr = IdentifyBadPoints(postForest, postVecExact, tol);CHKERRQ(ierr); 318c4762a1bSJed Brown } 319c4762a1bSJed Brown ierr = VecDestroy(&postVecExact);CHKERRQ(ierr); 320c4762a1bSJed Brown 321c4762a1bSJed Brown /* disconnect preForest from postForest if we don't test the transfer throughout the entire refinement process */ 322c4762a1bSJed Brown if (!transfer_from_base[1]) { 323c4762a1bSJed Brown ierr = DMForestSetAdaptivityForest(postForest,NULL);CHKERRQ(ierr); 324c4762a1bSJed Brown ierr = PetscObjectGetReference((PetscObject)preForest,&postCount);CHKERRQ(ierr); 325c4762a1bSJed Brown if (postCount != preCount) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Adaptation not memory neutral: reference count increase from %d to %d\n",preCount,postCount); 326c4762a1bSJed Brown } 327c4762a1bSJed Brown 328c4762a1bSJed Brown if (conv) { 329c4762a1bSJed Brown DM dmConv; 330c4762a1bSJed Brown 331c4762a1bSJed Brown ierr = DMConvert(postForest,DMPLEX,&dmConv);CHKERRQ(ierr); 332c4762a1bSJed Brown ierr = DMViewFromOptions(dmConv,NULL,"-dm_conv_view");CHKERRQ(ierr); 333c4762a1bSJed Brown ierr = DMPlexCheckCellShape(dmConv,PETSC_TRUE,PETSC_DETERMINE);CHKERRQ(ierr); 334c4762a1bSJed Brown ierr = DMDestroy(&dmConv);CHKERRQ(ierr); 335c4762a1bSJed Brown } 336c4762a1bSJed Brown 337c4762a1bSJed Brown ierr = VecDestroy(&preVec);CHKERRQ(ierr); 338c4762a1bSJed Brown ierr = DMDestroy(&preForest);CHKERRQ(ierr); 339c4762a1bSJed Brown 340c4762a1bSJed Brown preVec = postVecTransfer; 341c4762a1bSJed Brown preForest = postForest; 342c4762a1bSJed Brown } 343c4762a1bSJed Brown 344c4762a1bSJed Brown if (transfer_from_base[1]) { 345c4762a1bSJed Brown Vec baseVec, baseVecMapped; 346c4762a1bSJed Brown 347c4762a1bSJed Brown /* communicate between base and last adapted forest */ 348c4762a1bSJed Brown ierr = DMGetGlobalVector(base,&baseVec);CHKERRQ(ierr); 349c4762a1bSJed Brown ierr = DMProjectFunction(base,0.,funcs,ctxs,INSERT_VALUES,baseVec);CHKERRQ(ierr); 350c4762a1bSJed Brown ierr = PetscObjectSetName((PetscObject)baseVec,"Function Base");CHKERRQ(ierr); 351c4762a1bSJed Brown ierr = VecViewFromOptions(baseVec,NULL,"-vec_base_view");CHKERRQ(ierr); 352c4762a1bSJed Brown 353c4762a1bSJed Brown ierr = DMGetGlobalVector(preForest,&baseVecMapped);CHKERRQ(ierr); 354c4762a1bSJed Brown ierr = DMForestTransferVecFromBase(preForest,baseVec,baseVecMapped);CHKERRQ(ierr); 355c4762a1bSJed Brown ierr = VecViewFromOptions(baseVecMapped,NULL,"-vec_map_base_view");CHKERRQ(ierr); 356c4762a1bSJed Brown 357c4762a1bSJed Brown /* compare */ 358c4762a1bSJed Brown ierr = VecAXPY(baseVecMapped,-1.,preVec);CHKERRQ(ierr); 359c4762a1bSJed Brown ierr = VecViewFromOptions(baseVecMapped,NULL,"-vec_map_diff_view");CHKERRQ(ierr); 360c4762a1bSJed Brown ierr = VecNorm(baseVecMapped,NORM_2,&diff);CHKERRQ(ierr); 361c4762a1bSJed Brown 362c4762a1bSJed Brown /* output */ 363c4762a1bSJed Brown if (diff < tol) { 364c4762a1bSJed Brown ierr = PetscPrintf(comm,"DMForestTransferVecFromBase() passes.\n");CHKERRQ(ierr); 365c4762a1bSJed Brown } else { 366c4762a1bSJed Brown ierr = PetscPrintf(comm,"DMForestTransferVecFromBase() fails with error %g and tolerance %g\n",(double)diff,(double)tol);CHKERRQ(ierr); 367c4762a1bSJed Brown } 368c4762a1bSJed Brown 369c4762a1bSJed Brown ierr = DMRestoreGlobalVector(base,&baseVec);CHKERRQ(ierr); 370c4762a1bSJed Brown ierr = DMRestoreGlobalVector(preForest,&baseVecMapped);CHKERRQ(ierr); 371c4762a1bSJed Brown } 372c4762a1bSJed Brown 373c4762a1bSJed Brown /* cleanup */ 374c4762a1bSJed Brown ierr = VecDestroy(&preVec);CHKERRQ(ierr); 375c4762a1bSJed Brown ierr = DMDestroy(&preForest);CHKERRQ(ierr); 376c4762a1bSJed Brown ierr = DMDestroy(&base);CHKERRQ(ierr); 377c4762a1bSJed Brown ierr = PetscFinalize(); 378c4762a1bSJed Brown return ierr; 379c4762a1bSJed Brown } 380c4762a1bSJed Brown 381c4762a1bSJed Brown /*TEST 38230602db0SMatthew G. Knepley testset: 38330602db0SMatthew G. Knepley args: -dm_plex_simplex 0 -dm_plex_box_faces 3,3,3 -petscspace_type tensor 384c4762a1bSJed Brown 385c4762a1bSJed Brown test: 386c4762a1bSJed Brown output_file: output/ex2_2d.out 387c4762a1bSJed Brown suffix: p4est_2d 38830602db0SMatthew G. Knepley args: -petscspace_degree 2 389c4762a1bSJed Brown nsize: 3 390c4762a1bSJed Brown requires: p4est !single 391c4762a1bSJed Brown 392c4762a1bSJed Brown test: 393c4762a1bSJed Brown output_file: output/ex2_2d.out 394c4762a1bSJed Brown suffix: p4est_2d_deg4 39530602db0SMatthew G. Knepley args: -petscspace_degree 4 396c4762a1bSJed Brown requires: p4est !single 397c4762a1bSJed Brown 398c4762a1bSJed Brown test: 399c4762a1bSJed Brown output_file: output/ex2_2d.out 400c4762a1bSJed Brown suffix: p4est_2d_deg8 40130602db0SMatthew G. Knepley args: -petscspace_degree 8 402c4762a1bSJed Brown requires: p4est !single 403c4762a1bSJed Brown 404c4762a1bSJed Brown test: 405c4762a1bSJed Brown output_file: output/ex2_steps2.out 406c4762a1bSJed Brown suffix: p4est_2d_deg2_steps2 40730602db0SMatthew G. Knepley args: -petscspace_degree 2 -coords -adapt_steps 2 408c4762a1bSJed Brown nsize: 3 409c4762a1bSJed Brown requires: p4est !single 410c4762a1bSJed Brown 411c4762a1bSJed Brown test: 412c4762a1bSJed Brown output_file: output/ex2_steps3.out 413c4762a1bSJed Brown suffix: p4est_2d_deg3_steps3 41430602db0SMatthew G. Knepley args: -petscspace_degree 3 -coords -adapt_steps 3 -petscdualspace_lagrange_node_type equispaced -petscdualspace_lagrange_node_endpoints 1 415c4762a1bSJed Brown nsize: 3 416c4762a1bSJed Brown requires: p4est !single 417c4762a1bSJed Brown 418c4762a1bSJed Brown test: 419c4762a1bSJed Brown output_file: output/ex2_steps3.out 420c4762a1bSJed Brown suffix: p4est_2d_deg3_steps3_L2_periodic 42130602db0SMatthew 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 422c4762a1bSJed Brown nsize: 3 423c4762a1bSJed Brown requires: p4est !single 424c4762a1bSJed Brown 425c4762a1bSJed Brown test: 426c4762a1bSJed Brown output_file: output/ex2_steps3.out 427c4762a1bSJed Brown suffix: p4est_3d_deg2_steps3_L2_periodic 42830602db0SMatthew 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 429c4762a1bSJed Brown nsize: 3 430c4762a1bSJed Brown requires: p4est !single 431c4762a1bSJed Brown 432c4762a1bSJed Brown test: 433c4762a1bSJed Brown output_file: output/ex2_steps2.out 434c4762a1bSJed Brown suffix: p4est_3d_deg2_steps2 43530602db0SMatthew G. Knepley args: -dm_plex_dim 3 -petscspace_degree 2 -coords -adapt_steps 2 436c4762a1bSJed Brown nsize: 3 437c4762a1bSJed Brown requires: p4est !single 438c4762a1bSJed Brown 439c4762a1bSJed Brown test: 440c4762a1bSJed Brown output_file: output/ex2_steps3.out 441c4762a1bSJed Brown suffix: p4est_3d_deg3_steps3 44230602db0SMatthew G. Knepley args: -dm_plex_dim 3 -petscspace_degree 3 -coords -adapt_steps 3 -petscdualspace_lagrange_node_type equispaced -petscdualspace_lagrange_node_endpoints 1 443c4762a1bSJed Brown nsize: 3 444c4762a1bSJed Brown requires: p4est !single 445c4762a1bSJed Brown 446c4762a1bSJed Brown test: 447c4762a1bSJed Brown output_file: output/ex2_3d.out 448c4762a1bSJed Brown suffix: p4est_3d 44930602db0SMatthew G. Knepley args: -dm_plex_dim 3 -petscspace_degree 1 450c4762a1bSJed Brown nsize: 3 451c4762a1bSJed Brown requires: p4est !single 452c4762a1bSJed Brown 453c4762a1bSJed Brown test: 454c4762a1bSJed Brown output_file: output/ex2_3d.out 455c4762a1bSJed Brown suffix: p4est_3d_deg3 45630602db0SMatthew G. Knepley args: -dm_plex_dim 3 -petscspace_degree 3 457c4762a1bSJed Brown nsize: 3 458c4762a1bSJed Brown requires: p4est !single 459c4762a1bSJed Brown 460c4762a1bSJed Brown test: 461c4762a1bSJed Brown output_file: output/ex2_2d.out 462c4762a1bSJed Brown suffix: p4est_2d_deg2_coords 46330602db0SMatthew G. Knepley args: -petscspace_degree 2 -coords 464c4762a1bSJed Brown nsize: 3 465c4762a1bSJed Brown requires: p4est !single 466c4762a1bSJed Brown 467c4762a1bSJed Brown test: 468c4762a1bSJed Brown output_file: output/ex2_3d.out 469c4762a1bSJed Brown suffix: p4est_3d_deg2_coords 47030602db0SMatthew G. Knepley args: -dm_plex_dim 3 -petscspace_degree 2 -coords 471c4762a1bSJed Brown nsize: 3 472c4762a1bSJed Brown requires: p4est !single 473c4762a1bSJed Brown 474c4762a1bSJed Brown test: 475c4762a1bSJed Brown suffix: p4est_3d_nans 47630602db0SMatthew G. Knepley args: -dm_plex_dim 3 -dm_forest_partition_overlap 1 -test_convert -petscspace_degree 1 477c4762a1bSJed Brown nsize: 2 478c4762a1bSJed Brown requires: p4est !single 479c4762a1bSJed Brown 480c4762a1bSJed Brown test: 481c4762a1bSJed Brown TODO: not broken, but the 3D case below is broken, so I do not trust this one 482c4762a1bSJed Brown output_file: output/ex2_steps2.out 483c4762a1bSJed Brown suffix: p4est_2d_tfb_distributed_nc 48430602db0SMatthew G. Knepley args: -petscspace_degree 3 -dm_forest_maximum_refinement 2 -dm_p4est_refine_pattern hash -use_bcs 0 -coords -adapt_steps 2 -dm_distribute -petscpartitioner_type shell -petscpartitioner_shell_random 485c4762a1bSJed Brown nsize: 3 486c4762a1bSJed Brown requires: p4est !single 487c4762a1bSJed Brown 488c4762a1bSJed Brown test: 489c4762a1bSJed Brown TODO: broken 490c4762a1bSJed Brown output_file: output/ex2_steps2.out 491c4762a1bSJed Brown suffix: p4est_3d_tfb_distributed_nc 49230602db0SMatthew 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 -dm_distribute -petscpartitioner_type shell -petscpartitioner_shell_random 493c4762a1bSJed Brown nsize: 3 494c4762a1bSJed Brown requires: p4est !single 495c4762a1bSJed Brown 49630602db0SMatthew G. Knepley testset: 497012bc364SMatthew G. Knepley args: -petscspace_type tensor -dm_coord_space 0 -dm_plex_transform_type refine_tobox 49830602db0SMatthew G. Knepley 499c4762a1bSJed Brown test: 500c4762a1bSJed Brown TODO: broken 501c4762a1bSJed Brown output_file: output/ex2_3d.out 502c4762a1bSJed Brown suffix: p4est_3d_transfer_fails 503012bc364SMatthew 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 504c4762a1bSJed Brown requires: p4est !single 505c4762a1bSJed Brown 506c4762a1bSJed Brown test: 507c4762a1bSJed Brown TODO: broken 508c4762a1bSJed Brown output_file: output/ex2_steps2_notfb.out 509c4762a1bSJed Brown suffix: p4est_3d_transfer_fails_2 510012bc364SMatthew 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 511c4762a1bSJed Brown requires: p4est !single 512c4762a1bSJed Brown 513c4762a1bSJed Brown test: 514c4762a1bSJed Brown output_file: output/ex2_steps2.out 515c4762a1bSJed Brown suffix: p4est_3d_multi_transfer_s2t 516012bc364SMatthew 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 517012bc364SMatthew G. Knepley requires: p4est !single 518c4762a1bSJed Brown 519c4762a1bSJed Brown test: 520c4762a1bSJed Brown output_file: output/ex2_steps2.out 521c4762a1bSJed Brown suffix: p4est_3d_coords_transfer_s2t 522012bc364SMatthew 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 523012bc364SMatthew G. Knepley requires: p4est !single 52430602db0SMatthew G. Knepley 52530602db0SMatthew G. Knepley testset: 52630602db0SMatthew G. Knepley args: -dm_plex_simplex 0 -dm_plex_box_faces 3,3,3 52730602db0SMatthew G. Knepley 52830602db0SMatthew G. Knepley test: 52930602db0SMatthew G. Knepley output_file: output/ex2_2d_fv.out 53030602db0SMatthew G. Knepley suffix: p4est_2d_fv 53130602db0SMatthew G. Knepley args: -transfer_from_base 0 -use_fv -linear -dm_forest_partition_overlap 1 53230602db0SMatthew G. Knepley nsize: 3 53330602db0SMatthew G. Knepley requires: p4est !single 53430602db0SMatthew G. Knepley 53530602db0SMatthew G. Knepley test: 53630602db0SMatthew G. Knepley TODO: broken (codimension adjacency) 53730602db0SMatthew G. Knepley output_file: output/ex2_2d_fv.out 53830602db0SMatthew G. Knepley suffix: p4est_2d_fv_adjcodim 53930602db0SMatthew G. Knepley args: -transfer_from_base 0 -use_fv -linear -dm_forest_partition_overlap 1 -dm_forest_adjacency_codimension 1 54030602db0SMatthew G. Knepley nsize: 2 54130602db0SMatthew G. Knepley requires: p4est !single 54230602db0SMatthew G. Knepley 54330602db0SMatthew G. Knepley test: 54430602db0SMatthew G. Knepley TODO: broken (dimension adjacency) 54530602db0SMatthew G. Knepley output_file: output/ex2_2d_fv.out 54630602db0SMatthew G. Knepley suffix: p4est_2d_fv_adjdim 54730602db0SMatthew G. Knepley args: -transfer_from_base 0 -use_fv -linear -dm_forest_partition_overlap 1 -dm_forest_adjacency_dimension 1 54830602db0SMatthew G. Knepley nsize: 2 54930602db0SMatthew G. Knepley requires: p4est !single 55030602db0SMatthew G. Knepley 55130602db0SMatthew G. Knepley test: 55230602db0SMatthew G. Knepley output_file: output/ex2_2d_fv.out 55330602db0SMatthew G. Knepley suffix: p4est_2d_fv_zerocells 55430602db0SMatthew G. Knepley args: -transfer_from_base 0 -use_fv -linear -dm_forest_partition_overlap 1 55530602db0SMatthew G. Knepley nsize: 10 55630602db0SMatthew G. Knepley requires: p4est !single 55730602db0SMatthew G. Knepley 55830602db0SMatthew G. Knepley test: 55930602db0SMatthew G. Knepley output_file: output/ex2_3d_fv.out 56030602db0SMatthew G. Knepley suffix: p4est_3d_fv 56130602db0SMatthew G. Knepley args: -dm_plex_dim 3 -transfer_from_base 0 -use_fv -linear -dm_forest_partition_overlap 1 56230602db0SMatthew G. Knepley nsize: 3 563c4762a1bSJed Brown requires: p4est !single 564c4762a1bSJed Brown 565c4762a1bSJed Brown TEST*/ 566