10bb0ac0dSMatthew G. Knepley static char help[] = "Test HDF5 input and output.\n\ 20bb0ac0dSMatthew G. Knepley This exposed a bug with sharing discretizations.\n\n\n"; 30bb0ac0dSMatthew G. Knepley 40bb0ac0dSMatthew G. Knepley #include <petscdmforest.h> 50bb0ac0dSMatthew G. Knepley #include <petscdmplex.h> 60bb0ac0dSMatthew G. Knepley #include <petscviewerhdf5.h> 70bb0ac0dSMatthew G. Knepley 80bb0ac0dSMatthew G. Knepley int main (int argc, char **argv) 90bb0ac0dSMatthew G. Knepley { 100bb0ac0dSMatthew G. Knepley 110bb0ac0dSMatthew G. Knepley DM base, forest, plex; 120bb0ac0dSMatthew G. Knepley Vec g, g2; 130bb0ac0dSMatthew G. Knepley PetscSection s; 140bb0ac0dSMatthew G. Knepley PetscViewer viewer; 150bb0ac0dSMatthew G. Knepley PetscReal diff; 160bb0ac0dSMatthew G. Knepley PetscInt min_refine = 2, overlap = 0; 170bb0ac0dSMatthew G. Knepley PetscInt vStart, vEnd, v; 180bb0ac0dSMatthew G. Knepley PetscErrorCode ierr; 190bb0ac0dSMatthew G. Knepley 200bb0ac0dSMatthew G. Knepley ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr; 210bb0ac0dSMatthew G. Knepley 22*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMCreate(PETSC_COMM_WORLD, &base)); 23*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMSetType(base, DMPLEX)); 24*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMSetFromOptions(base)); 250bb0ac0dSMatthew G. Knepley 26*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMCreate(PETSC_COMM_WORLD, &forest)); 27*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMSetType(forest, DMP4EST)); 28*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMForestSetBaseDM(forest, base)); 29*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMForestSetInitialRefinement(forest, min_refine)); 30*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMForestSetPartitionOverlap(forest, overlap)); 31*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMSetUp(forest)); 32*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMDestroy(&base)); 33*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMViewFromOptions(forest, NULL, "-dm_view")); 340bb0ac0dSMatthew G. Knepley 35*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMConvert(forest, DMPLEX, &plex)); 36*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMPlexGetDepthStratum(plex, 0, &vStart, &vEnd)); 37*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMDestroy(&plex)); 38*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscSectionCreate(PetscObjectComm((PetscObject) forest), &s)); 39*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscSectionSetChart(s, vStart, vEnd)); 40*5f80ce2aSJacob Faibussowitsch for (v = vStart; v < vEnd; ++v) CHKERRQ(PetscSectionSetDof(s, v, 1)); 41*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscSectionSetUp(s)); 42*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMSetLocalSection(forest, s)); 43*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscSectionDestroy(&s)); 440bb0ac0dSMatthew G. Knepley 45*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMCreateGlobalVector(forest, &g)); 46*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscObjectSetName((PetscObject) g, "g")); 47*5f80ce2aSJacob Faibussowitsch CHKERRQ(VecSet(g, 1.0)); 48*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerHDF5Open(PETSC_COMM_WORLD, "forest.h5", FILE_MODE_WRITE, &viewer)); 49*5f80ce2aSJacob Faibussowitsch CHKERRQ(VecView(g, viewer)); 50*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerDestroy(&viewer)); 510bb0ac0dSMatthew G. Knepley 52*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMCreateGlobalVector(forest, &g2)); 53*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscObjectSetName((PetscObject) g2, "g")); 54*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerHDF5Open(PETSC_COMM_WORLD, "forest.h5", FILE_MODE_READ, &viewer)); 55*5f80ce2aSJacob Faibussowitsch CHKERRQ(VecLoad(g2, viewer)); 56*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerDestroy(&viewer)); 570bb0ac0dSMatthew G. Knepley 580bb0ac0dSMatthew G. Knepley /* Check if the data is the same*/ 59*5f80ce2aSJacob Faibussowitsch CHKERRQ(VecAXPY(g2, -1.0, g)); 60*5f80ce2aSJacob Faibussowitsch CHKERRQ(VecNorm(g2, NORM_INFINITY, &diff)); 61*5f80ce2aSJacob Faibussowitsch if (diff > PETSC_MACHINE_EPSILON) CHKERRQ(PetscPrintf(PETSC_COMM_WORLD, "Check failed: %g\n", (double) diff)); 620bb0ac0dSMatthew G. Knepley 63*5f80ce2aSJacob Faibussowitsch CHKERRQ(VecDestroy(&g)); 64*5f80ce2aSJacob Faibussowitsch CHKERRQ(VecDestroy(&g2)); 65*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMDestroy(&forest)); 660bb0ac0dSMatthew G. Knepley ierr = PetscFinalize(); 670bb0ac0dSMatthew G. Knepley return ierr; 680bb0ac0dSMatthew G. Knepley } 690bb0ac0dSMatthew G. Knepley 700bb0ac0dSMatthew G. Knepley /*TEST 710bb0ac0dSMatthew G. Knepley 720bb0ac0dSMatthew G. Knepley build: 730bb0ac0dSMatthew G. Knepley requires: hdf5 p4est 740bb0ac0dSMatthew G. Knepley 750bb0ac0dSMatthew G. Knepley test: 7630602db0SMatthew G. Knepley args: -dm_plex_simplex 0 -dm_plex_box_faces 3,3 770bb0ac0dSMatthew G. Knepley 780bb0ac0dSMatthew G. Knepley TEST*/ 79