[petsc-users] Periodic BC in petsc DMPlex / firedrake
Stefano Zampini
stefano.zampini at gmail.com
Thu Oct 25 08:36:15 CDT 2018
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
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20181025/e997dff5/attachment-0001.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: periodic.png
Type: image/png
Size: 255676 bytes
Desc: not available
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20181025/e997dff5/attachment-0001.png>
More information about the petsc-users
mailing list