static char help[] = "Second Order TVD Finite Volume Example.\n";


#include <petscdmplex.h>
#include <petscds.h>
#include <petscts.h>
#include <petscsf.h> /* For SplitFaces() */

#define DIM 2                   /* Geometric dimension */
#define ALEN(a) (sizeof(a)/sizeof((a)[0]))

#undef __FUNCT__
#define __FUNCT__ "main"
int main(int argc, char **argv)
{
    MPI_Comm          comm;
    DM             dm, dmDist = NULL;
    PetscErrorCode    ierr;
    char              filename[PETSC_MAX_PATH_LEN] = "cavity.exo";

    PetscSection   section;
    PetscViewer    viewer;

    PetscInt       i;
    PetscInt       dim = 2, numFields, numBC, overlap;
    PetscInt       numComp[3];
    PetscInt       numDof[12];
    PetscInt       *bcField;
    IS             *bcPointIS, regionIS, innerIS;
    const PetscInt *regions;

    Vec X, localX;
    PetscBool hasLabel;
    PetscInt numRegions;

    ierr = PetscInitialize(&argc, &argv, (char*) 0, help);CHKERRQ(ierr);
    comm = PETSC_COMM_WORLD;

    overlap = 2;

    /* 
     * Reading the grid. 
     */
    ierr = DMPlexCreateFromFile(comm, filename, PETSC_TRUE, &dm);CHKERRQ(ierr);

    ierr = DMCreateLabel(dm, "Face Sets");CHKERRQ(ierr);

    /* 
     * Distribute mesh over processes 
     */
    {
        ierr = DMPlexSetAdjacencyUseCone(dm, PETSC_TRUE);CHKERRQ(ierr);
        ierr = DMPlexSetAdjacencyUseClosure(dm, PETSC_FALSE);CHKERRQ(ierr);
        ierr = DMPlexDistribute(dm, overlap, NULL, &dmDist);CHKERRQ(ierr);
        if (dmDist) {
            ierr = DMDestroy(&dm);CHKERRQ(ierr);
            dm   = dmDist;
        }
    }

    {
        DM gdm;
        PetscInt numGhostCells;

        ierr = DMPlexConstructGhostCells(dm, "Face Sets", &numGhostCells, &gdm);CHKERRQ(ierr);
        ierr = DMDestroy(&dm);CHKERRQ(ierr);
        dm   = gdm;
        ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr);
    }

    ierr = PetscFinalize();
    return(0);
}