[petsc-users] Periodic BC in petsc DMPlex / firedrake
Matthew Knepley
knepley at gmail.com
Thu Oct 25 09:43:57 CDT 2018
Good catch Stefano.
Matt
On Thu, Oct 25, 2018 at 9:36 AM Stefano Zampini <stefano.zampini at gmail.com>
wrote:
> Maybe this is a fix
>
> diff --git a/src/dm/impls/plex/plexvtu.c b/src/dm/impls/plex/plexvtu.c
> index acdea12c2f..1a8bbada6a 100644
> --- a/src/dm/impls/plex/plexvtu.c
> +++ b/src/dm/impls/plex/plexvtu.c
> @@ -465,10 +465,11 @@ PetscErrorCode DMPlexVTKWriteAll_VTU(DM
> dm,PetscViewer viewer)
> if ((closure[v] >= vStart) && (closure[v] < vEnd)) {
> PetscScalar *xpoint;
>
> - ierr =
> DMPlexPointLocalRead(dm,v,x,&xpoint);CHKERRQ(ierr);
> + ierr =
> DMPlexPointLocalRead(dm,closure[v],x,&xpoint);CHKERRQ(ierr);
> y[cnt + off++] = xpoint[i];
> }
> }
> + cnt += off;
> ierr = DMPlexRestoreTransitiveClosure(dm, c, PETSC_TRUE,
> &closureSize, &closure);CHKERRQ(ierr);
> }
> }
>
> Max, does this fix your problem? If you confirm, I'll fix this in the
> maint branch
>
> If I run the below command line with the patch and with snes tutorials
> ex12 I get the nice picture attached
> $ ./ex12 -quiet -run_type test -interpolate 1 -bc_type dirichlet
> -petscspace_degree 1 -vec_view vtk:test.vtu:vtk_vtu -f
> ${PETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh
> -dm_plex_gmsh_periodic -x_periodicity periodic -y_periodicity periodic
> -dm_refine 4
>
> Il giorno gio 25 ott 2018 alle ore 15:11 Stefano Zampini <
> stefano.zampini at gmail.com> ha scritto:
>
>> Matt,
>>
>> you can reproduce it via
>>
>> $ valgrind ./ex12 -quiet -run_type test -interpolate 1 -bc_type dirichlet
>> -petscspace_degree 1 -vec_view vtk:test.vtu:vtk_vtu -f
>> ${PETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh
>> -dm_plex_gmsh_periodic
>>
>> Long time ago I added support for viewing meshes with periodic vertices
>> in the VTK_VTU viewer, but I did not fix the part that writes fields
>>
>>
>> Il giorno mer 24 ott 2018 alle ore 21:04 Matthew Knepley <
>> knepley at gmail.com> ha scritto:
>>
>>> On Wed, Oct 24, 2018 at 11:36 AM Maximilian Hartig <
>>> imilian.hartig at gmail.com> wrote:
>>>
>>>>
>>>>
>>>> On 24. Oct 2018, at 12:49, Matthew Knepley <knepley at gmail.com> wrote:
>>>>
>>>> On Wed, Oct 24, 2018 at 6:29 AM Lawrence Mitchell <wence at gmx.li> wrote:
>>>>
>>>>> Hi Max,
>>>>>
>>>>> (I'm cc'ing in the petsc-users mailing list which may have more
>>>>> advice, if you are using PETSc you should definitely subscribe!
>>>>>
>>>>> > On 24 Oct 2018, at 09:27, Maximilian Hartig <
>>>>> imilian.hartig at gmail.com> wrote:
>>>>> >
>>>>> > Hello Lawrence,
>>>>> >
>>>>> > sorry to message you out of the blue. My name is Max and I found
>>>>> your post on GitHub (
>>>>> https://github.com/firedrakeproject/firedrake/issues/1246 ) on DMPlex
>>>>> being able to read periodic gmsh files. I am currently trying to do just
>>>>> that (creating a periodic DMPlex mesh with gmsh) in the context of my PhD
>>>>> work. So far I haven’ t found any documentation on the periodic BC’s with
>>>>> DMPlex and gmsh in the official petsc documentation.
>>>>> > I was wondering whether you’d be so kind as to point me in a general
>>>>> direction concerning how to achieve this. You seem experienced in using
>>>>> petsc and I would greatly appreciate your help.
>>>>>
>>>>>
>>>>> I think the answer is "it depends". If you're just using DMPlex
>>>>> directly and all the of the functionality with PetscDS, then I /think/ that
>>>>> reading periodic meshes via gmsh (assuming you're using the appropriate
>>>>> gmsh mesh format [v2]) "just works".
>>>>>
>>>>
>>>> There are two phases here: topological and geometric. DMPlex represents
>>>> the periodic topological entity directly. For example, a circle is just a
>>>> segment with one end hooked to the other. Vertices are not duplicated, or
>>>> mapped to each other. This makes topology simple and easy to implement.
>>>> However, then geometry is more complicated. What Plex does is allow
>>>> coordinates to be represented by a discontinuous field taking values on
>>>> cells, in addition to vertices. In our circle example, each cells near the
>>>> cut will have 2 coordinates, one for each vertex, but they will not agree
>>>> across the cut. If you define a periodic domain, then Plex can construct
>>>> this coordinate field automatically using DMPlexLocalize(). These DG
>>>> coordinates are then used by the integration routines.
>>>>
>>>>
>>>> Ok, I think I understand the concept. DMPlex reads information about
>>>> both topology and coordinates from the .msh file. Creating a periodic mesh
>>>> in gmsh then should allow DMPlex to identify the periodic boundaries as the
>>>> “cut” and build the mesh topology accordingly. Coordinate information is
>>>> handled separately.
>>>> That means, as Lawrence suggested, building periodic meshes in gmsh and
>>>> reading them in to petsc’s DMPlex should indeed “just work”. (From the
>>>> user perspective). The only extra step is to call DMLocalizeCoordinates()
>>>> after DMPlexReadFromFile(). Sorry to reiterate, I am just trying to make
>>>> sense of this.
>>>>
>>>>
>>>>
>>>>> From my side, the issue is to do with mapping that coordinate field
>>>>> into one that I understand (within Firedrake). You may not have this
>>>>> problem.
>>>>>
>>>>
>>>> Firedrake uses its own coordinate mapping and integration routines, so
>>>> they must manage the second part independently. I hope to get change this
>>>> slightly soon by making the Firedrake representation a DMField, so that it
>>>> looks the same to Plex.
>>>>
>>>> Thanks,
>>>>
>>>> Matt
>>>>
>>>>
>>>>> Thanks,
>>>>>
>>>>> Lawrence
>>>>
>>>>
>>>>
>>>> --
>>>> What most experimenters take for granted before they begin their
>>>> experiments is infinitely more interesting than any results to which their
>>>> experiments lead.
>>>> -- Norbert Wiener
>>>>
>>>> https://www.cse.buffalo.edu/~knepley/
>>>>
>>>>
>>>> To read periodic meshes from GMSH, you need to use the option
>>>> -dm_plex_gmsh_periodic and DMPlexCreateFromFile
>>>>
>>>>
>>>> Ahh, thanks. I was missing the option " -dm_plex_gmsh_periodic “. But
>>>> using this option I now generate a segmentation fault error when calling
>>>> VecView() on the solution vector with vtk and hdf5 viewers. Any suggestions?
>>>>
>>>
>>> Small example? VTK is deprecated. HDF5 should work, although it will
>>> require you to have proper coordinates I think. We have to
>>> think about what you mean. If its for a checkpoint, there is no problem,
>>> but for viz, those programs do not understand periodicity. Thus I embed it
>>> in a higher dimensional space.
>>>
>>> Matt
>>>
>>>> See src/dm/impls/plex/examples/tests/ex1.c. An example runs
>>>>
>>>> $ ./ex1 -filename
>>>> ${PETSC_DIR}/share/petsc/datafiles/meshes/cube_periodic_bin.msh
>>>> -dm_plex_gmsh_periodic -dm_view ::ascii_info_detail -interpolate -test_shape
>>>>
>>>> generating periodic meshes in gmsh may be tricky, Lisandro for sure may
>>>> advice.
>>>>
>>>>
>>>> Thanks,
>>>> Max
>>>>
>>>>
>>>
>>> --
>>> What most experimenters take for granted before they begin their
>>> experiments is infinitely more interesting than any results to which their
>>> experiments lead.
>>> -- Norbert Wiener
>>>
>>> https://www.cse.buffalo.edu/~knepley/
>>> <http://www.cse.buffalo.edu/~knepley/>
>>>
>>
>>
>> --
>> Stefano
>>
>
>
> --
> Stefano
>
--
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which their
experiments lead.
-- Norbert Wiener
https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20181025/8fb15752/attachment.html>
More information about the petsc-users
mailing list