[petsc-users] DMPlex from periodic Gmsh and ghost cells

[Thibault Bridel-Bertomeu]

Dear all,

I am playing around with creating DMPlex from a periodic Gmsh mesh (still
in the same finite volume code that solves the Euler equations) because I
want to run some classical periodic test cases from the literature (like
convecting a vortex and so on).
I got from the mailing list and the examples that I gotta use
-dm_plex_gmsh_periodic for the reader and the DM creator to do their job
properly and read the $Periodic section from the Gmsh mesh. That works
fine, I can write the DM in HDF5 format afterwards and the periodicity is
clearly there (i can see in Paraview that opposite sides are "attached" to
each other).
Now, in the context of the finite volume solver, I wrote my one right-hand
side term computing routine, and I rely heavily on routines such as
DMPlexGetFaceFields. Those routines, in turn rely on the ghost cells and
... I am having an issue with DMPlexConstructGhostCells. When I use
-dm_plex_gmsh_periodic, DMPlexCreateFromFile works fine,
DMSetBasicAdjacency(dm, true, false) works fine, DMPlexDistribute works
fine, but when I get to DMPlexConstructGhostCells, I get

 "DM has boundary face 2048 with 2 support cells"

I understand the DM begin periodic, even the "boundary" faces have two
neighbors because of the topological folding, so ... yeah, the ghost cells
cannot really exist. But then how can I proceed ? Because I need that
"ghost" label for all the routines I call from the RHS computing routine ...

I attach to this email the periodic Gmsh mesh I am playing with, and the
routine from my code where I create the DM and so on:

subroutine initmesh

PetscErrorCode :: ierr
DM :: dmDist, dmGhost
integer :: overlap
PetscViewer :: vtkViewer

call PetscLogEventBegin(timer_initmesh, ierr); CHKERRA(ierr)

! Number of neighbours taken into account in MP communications(1 - Order 1;
2 - Order 2)
overlap = 1

call PetscPrintf(PETSC_COMM_WORLD, "Initializing mesh...\n", ierr) ;
CHKERRA(ierr)

! Force DMPlex to use gmsh marker
! call PetscOptionsSetValue(PETSC_NULL_OPTIONS, "-dm_plex_gmsh_use_marker",
"true", ierr); CHKERRA(ierr)

! Read mesh from file name 'meshname'
call DMPlexCreateFromFile(PETSC_COMM_WORLD, meshname, PETSC_TRUE, dm,
ierr); CHKERRA(ierr)

! Distribute on processors
! Start with connectivity
call DMSetBasicAdjacency(dm, PETSC_TRUE, PETSC_FALSE, ierr) ; CHKERRA(ierr)

! Distribute on processors
call DMPlexDistribute(dm, overlap, PETSC_NULL_SF, dmDist, ierr) ;
CHKERRA(ierr)

! Security check
if (dmDist /= PETSC_NULL_DM) then
! Destroy previous dm
call DMDestroy(dm, ierr) ; CHKERRA(ierr)
! Replace with dmDist
dm = dmDist
end if

! Finalize setup of the object
call DMSetFromOptions(dm, ierr) ; CHKERRA(ierr)

! Boundary condition with ghost cells
call DMPlexConstructGhostCells(dm, PETSC_NULL_CHARACTER,
PETSC_NULL_INTEGER, dmGhost, ierr); CHKERRA(ierr)

! Security check
if (dmGhost /= PETSC_NULL_DM) then
! Destroy previous dm
call DMDestroy(dm, ierr) ; CHKERRA(ierr)
! Replace with dmGhost
dm = dmGhost
end if

if (debug) then
! Show in terminal
call PetscPrintf(PETSC_COMM_WORLD, ":: [DEBUG] Visualizing DM in console
::\n", ierr); CHKERRA(ierr)
call DMView(dm, PETSC_VIEWER_STDOUT_WORLD, ierr) ; CHKERRA(ierr)
! VTK viewer
call PetscViewerCreate(PETSC_COMM_WORLD, vtkViewer, ierr) ; CHKERRA(ierr)
call PetscViewerSetType(vtkViewer, PETSCVIEWERHDF5, ierr) ; CHKERRA(ierr)
call PetscViewerFileSetMode(vtkViewer, FILE_MODE_WRITE, ierr) ;
CHKERRA(ierr)
call PetscViewerFileSetName(vtkViewer, "debug_initmesh.h5", ierr) ;
CHKERRA(ierr)
call DMView(dm, vtkViewer, ierr) ; CHKERRA(ierr)
call PetscViewerDestroy(vtkViewer, ierr) ; CHKERRA(ierr)
end if

call PetscPrintf(PETSC_COMM_WORLD, "Done !\n", ierr) ; CHKERRA(ierr)

call PetscLogEventEnd(timer_initmesh, ierr); CHKERRA(ierr)

end subroutine initmesh

Thank you very much for your help !

Best regards,

Thibault Bridel-Bertomeu
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20210518/30438f89/attachment-0001.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: vortex.msh
Type: application/octet-stream
Size: 69414 bytes
Desc: not available
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20210518/30438f89/attachment-0001.obj>