CGAL users discussion list

[Mael Rouxel-Labbé <>]

Hello,

The following code will collect the faces of the arbitrarily-chosen "boundary" of one tetrahedral tile of the periodic mesh: https://gist.github.com/MaelRL/dd267ad58560b3ae97d1b69102651f49.

It is adapted from the example mesh_implicit_shape.cpp of the documentation of Periodic_3_mesh_3; the key point is:

      if(c3t3.triangulation().neighbor_offset(c, ni) == Tr::Offset())
        continue;

which uses the fact that we use a single representation for each unique element of the mesh and when we reach the "boundary" of the tile, we add an offset to the neighbor. See the documentation of Periodic_2_triangulation_2 (2D, for visual simplicity, it's the same mechanism): https://doc.cgal.org/latest/Periodic_2_triangulation_2/index.html#Chapter_2D_Periodic_Triangulations. There are some more comments in the code if needed.

Attached is a picture of the result. (Note that I didn't pay attention to the orientation in which the facets are output.)

Would this work for you?

Best,
Mael

On 17/03/2021 16:18, Matthew Russell ( via cgal-discuss Mailing List) wrote:

Hi Mael,

It was to be able to easily identify the boundaries in code, by
selecting boundary faces within an epsilon of x == 0 etc. If it's easy
to identify the faces on each boundary of a general periodic mesh
using CGAL and export that data along with a single instance of the
mesh, then that might also work.

I can imagine a naive algorithm that for each boundary face, it loops
over all boundary faces again looking for a counterpart separated by a
translation, labelling the faces accordingly. Is there a nicer way to
achieve this using CGAL functionality? I haven't come across any
examples so far.

Thanks,
Matt

On Wed, 17 Mar 2021 at 08:24, Mael Rouxel-Labbé
 wrote:
> Hello,
>
> Assuming you cannot modify your FEM formulation to simply use a single
> mesh as Iordan was describing, do you wish to have the facets to be
> aligned on the side of the cube to identify them as boundary, or because
> it is necessary in your FEM formulation that they geometrically lie on
> the cube's side ? If it's not a geometrical constraint, extracting a
> canonical block of tetrahedral elements with well-identified boundary
> faces that reproduces periodically in all three directions (basically,
> what you see in the drawing in the manual) is very easy to do.
>
> Best,
> Mael
>
> Does the FEM formulation require the boundary tets to have their
> boundary face fully aligned with a side of the cube? After all,
>
> On 17/03/2021 01:19, Matthew Russell ( via
> cgal-discuss Mailing List) wrote:
> > Dear Iordan,
> >
> > Thanks for this, and, as far as I'm concerned, no need to apologise.
> > It's an interesting way to look at the the problem and it makes a lot
> > of sense conceptually. However, in my experience, it is not how PBCs
> > are usually implemented in finite element solvers. Do you know of any
> > examples describing how to reformulate a periodic problem as a problem
> > on a torus or it's numerical implementation? The problem I'm solving
> > is Stokes flow in a periodic unit cell, with no-slip on the
> > non-periodic part of the boundary (the part that's kind-of
> > "inward-facing"; see `gmsh.jpg` in my first post).
> >
> > Thanks,
> > Matt
> >
> > On Tue, 16 Mar 2021 at 23:30, Iordan Iordanov  wrote:
> > > Dear Matthew, Monique, all,
> > >
> > > Sorry to appear out of nowhere like this, but I would like to offer another perspective, perhaps.
> > >
> > > Matthew, you mention that in your FEM problem you want to define periodic boundary conditions.
> > > The need for specifying periodic boundary conditions in FEM problems comes from the fact that
> > > traditional meshes do not have a way to represent the periodic medium over which the equation
> > > to be solved is defined. In other words, you need to discretize your domain so that the vertices and
> > > faces coincide neatly, so that you know how to compute interactions and value propagations
> > > periodically, from one side to the other.
> > >
> > > The thing is, with a periodic mesh, you don't need to specify any boundary conditions, because there
> > > is no boundary. The CGAL periodic 3D mesh package represents a domain that "wraps around" in
> > > 3 dimensions, as the name implies, so the only things you need to define are the equation itself and the
> > > initial conditions. In this setting, boundary conditions actually do not make sense -- if you try to apply
> > > boundary conditions, you would end up splitting a torus only to glue it back together.
> > >
> > > Why not just use the torus? :)
> > >
> > > I am guessing it will require some re-formulating of your problem, but you will get much better results
> > > in the end, and your computation will be more efficient. The reason is that, when you define a traditional
> > > mesh, you need to store the boundary elements (vertices/faces) twice. Moreover, you need to keep track
> > > of which vertices coincide and when not to count interactions twice, so it is a bit messy. With the periodic
> > > mesh, you really don't save anything extra, and all the elements are unique. Every element knows its
> > > neighbors so computations are much more straightforward conceptually.
> > >
> > > I hope this has some meaning for you, and, once more, sorry for butting in like this.
> > >
> > > Best of luck and best regards,
> > >    - Iordan
> > >
> > > On Wed, Mar 17, 2021 at 3:29 AM Matthew Russell  wrote:
> > > > Dear Monique,
> > > >
> > > > Thanks - I saw that example and have tried using a unit cube (defined
> > > > implicitly) as the second function in the vector but I don't think it
> > > > solves my problem. Is this what you had in mind? The resulting mesh is
> > > > similar to the example attached to my original email, except that it
> > > > can also include the region with the opposite sign (which I don't
> > > > require), depending on the chosen sign vectors.
> > > >
> > > > Apologies if I described the problem poorly! I would like to create a
> > > > collection of different meshes, one for each of a collection of
> > > > implicitly-defined regions, so I don't think it's a multiple domains
> > > > problem. On each mesh I will run a separate finite element simulation
> > > > with periodic boundary conditions. We can ignore the multiple meshes
> > > > here and consider just one. To make it easy to identify the periodic
> > > > boundaries in the FEM code, my first thought was to try to create
> > > > meshes where the cells at the boundaries lie on the facets of the
> > > > bounding cube as I mentioned earlier. It seems this isn't possible
> > > > with the periodic mesh package, but there may be a workaround using
> > > > the regular 3D mesh package as Mael suggested. I was wondering if
> > > > someone had any further hints on how to achieve this? It seems like a
> > > > difficult problem to manually place vertices on the intersection
> > > > between the bound cube's facets and the implicit region. I really
> > > > appreciate all your help with this!
> > > >
> > > > Thanks,
> > > > Matt
> > > >
> > > > On Tue, 16 Mar 2021 at 15:44, Monique Teillaud
> > > >
> > > >
> > > >
> > > > On Tue, 16 Mar 2021 at 15:44, Monique Teillaud
> > > >
> > > >  wrote:
> > > > > Dear Matthew,
> > > > >
> > > > > If I understand correctly (I am not sure) you would like to mesh multiple domains.
> > > > > Does this example help you?
> > > > > https://doc.cgal.org/latest/Periodic_3_mesh_3/index.html#title19
> > > > >
> > > > > Best,
> > > > > --
> > > > > Monique Teillaud
> > > > > https://members.loria.fr/Monique.Teillaud/
> > > > > INRIA Nancy - Grand Est, LORIA
> > > > > Institut National de Recherche en Informatique et Automatique
> > > > >
> > > > > ----- Le 16 Mar 21, à 15:48, cgal-discuss  a écrit :
> > > > >
> > > > > > Dear Monique,
> > > > > >
> > > > > > I agree that adding points on the facets of the bounding cube in the
> > > > > > right places would be a way to solve this. But this technique requires
> > > > > > a description of the regions on each facet that should be
> > > > > > triangulated. This is difficult because many of the surfaces I am
> > > > > > dealing with are defined implicitly and do not cover the entirety of
> > > > > > the cube's facets. Is there something in CGAL that would help do this?
> > > > > >
> > > > > > However, I suppose my actual issue is being able to accurately specify
> > > > > > boundary conditions. I was trying to make this easy by forcing the
> > > > > > boundaries of the mesh to align with bounding cube.
> > > > > >
> > > > > > With a triangulation generated by the periodic mesh package, is there
> > > > > > a way to colour or otherwise somehow mark the cell faces on each
> > > > > > boundary? I'm not sure this is a meaningful question given the
> > > > > > irregular nature of the boundaries!
> > > > > >
> > > > > > Sorry for the basic questions - I'm new to CGAL and still finding my way around!
> > > > > >
> > > > > > Many thanks,
> > > > > > Matt
> > > > > >
> > > > > >
> > > > > > On Tue, 16 Mar 2021 at 10:54, Monique Teillaud
> > > > > >  wrote:
> > > > > > > Dear Matthew,
> > > > > > >
> > > > > > > In fact, if you add points on the facets of the cube and you ensure that the
> > > > > > > facets are triangulated as mael suggests, then you are enforcing the
> > > > > > > periodicity through these points. You are not using the functionality that is
> > > > > > > offered by the periodic mesh package.
> > > > > > > In other words, this should rather be done with the non-periodic meshing
> > > > > > > package.
> > > > > > >
> > > > > > > (In fact I would be curious why you need the facets to be triangulated... Our
> > > > > > > algorithm is precisely designed to avoid adding such artificial points...)
> > > > > > >
> > > > > > > Best regards,
> > > > > > > --
> > > > > > > Monique Teillaud
> > > > > > > https://members.loria.fr/Monique.Teillaud/
> > > > > > > INRIA Nancy - Grand Est, LORIA
> > > > > > > Institut National de Recherche en Informatique et Automatique
> > > > > > >
> > > > > > > ----- Le 15 Mar 21, à 17:21, cgal-discuss  a écrit :
> > > > > > >
> > > > > > > Hi Mael,
> > > > > > >
> > > > > > > Thanks for your reply.
> > > > > > >
> > > > > > > It's a shame the periodic mesher can't do this currently. Would you be able to
> > > > > > > give some hints on how to implement your idea of using the standard 3D mesh
> > > > > > > package? I'm slightly wary about manually moving vertices around - won't this
> > > > > > > reduce the quality of the resulting mesh?
> > > > > > >
> > > > > > > Thanks,
> > > > > > > Matt
> > > > > > >
> > > > > > > On Mon, 15 Mar 2021, 08:04 Mael Rouxel-Labbé,
> > > > > > >  wrote:
> > > > > > > > Hello,
> > > > > > > >
> > > > > > > > This is not currently possible in a straightforward manner in the
> > > > > > > > periodic mesh generator, as there is no 3D constrained Delaunay
> > > > > > > > triangulations in CGAL yet.
> > > > > > > >
> > > > > > > > In a less straight forward manner, if you wish your elements to lie
> > > > > > > > exactly on the boundary of the periodic domain, you might be able to do
> > > > > > > > away with the periodic mesh generation algorithm: generate a 3D mesh of
> > > > > > > > the implicit function clipped by the periodic domain, and then ensure
> > > > > > > > that the vertices are identical on both sides.
> > > > > > > >
> > > > > > > > Best,
> > > > > > > > Mael
> > > > > > > >
> > > > > > > > On 15/03/2021 02:05, Matthew Russell ( via
> > > > > > > > cgal-discuss Mailing List) wrote:
> > > > > > > > > Dear all,
> > > > > > > > >
> > > > > > > > > I would like to generate periodic 3D meshes of regions bounded by
> > > > > > > > > implicit functions and so I was pleased to find examples such as [1]
> > > > > > > > > that do almost exactly what I want. The problem that I'm having is
> > > > > > > > > that cells in the meshes generated from this example are able to lie
> > > > > > > > > across the bounding unit cube, but to use the meshes in finite element
> > > > > > > > > simulations I need the faces of the cells to align with the boundaries
> > > > > > > > > of the cube.
> > > > > > > > >
> > > > > > > > > For example, the CGAL example [1] (with a different implciit function)
> > > > > > > > > generates meshes like [cgal.jpg], but I would like one like
> > > > > > > > > [gmsh.jpg].
> > > > > > > > >
> > > > > > > > > I haven't been able to find a way to do this in CGAL. Is it possible?
> > > > > > > > >
> > > > > > > > > Many thanks,
> > > > > > > > > Matt
> > > > > > > > >
> > > > > > > > > [1]
> > > > > > > > > https://doc.cgal.org/latest/Periodic_3_mesh_3/Periodic_3_mesh_3_2mesh_implicit_shape_8cpp-example.html
> > > > > > > > >
> > > > > > > > >                       
> > > > > > > > --
> > > > > > > > You are currently subscribed to cgal-discuss.
> > > > > > > > To unsubscribe or access the archives, go to
> > > > > > > > https://sympa.inria.fr/sympa/info/cgal-discuss
> > > > > > > >
> > > > > > > >
> > > > > > > >                     
> > > > > > > --
> > > > > > > You are currently subscribed to cgal-discuss.
> > > > > > > To unsubscribe or access the archives, go to
> > > > > > > https://sympa.inria.fr/sympa/info/cgal-discuss
> > > > > > >
> > > > > > >
> > > > > > > --
> > > > > > > You are currently subscribed to cgal-discuss.
> > > > > > > To unsubscribe or access the archives, go to
> > > > > > > https://sympa.inria.fr/sympa/info/cgal-discuss
> > > > > > >
> > > > > > >                   
> > > > > > --
> > > > > > You are currently subscribed to cgal-discuss.
> > > > > > To unsubscribe or access the archives, go to
> > > > > > https://sympa.inria.fr/sympa/info/cgal-discuss
> > > > > >                 
> > > > > --
> > > > > You are currently subscribed to cgal-discuss.
> > > > > To unsubscribe or access the archives, go to
> > > > > https://sympa.inria.fr/sympa/info/cgal-discuss
> > > > >
> > > > >
> > > > >               
> > > > --
> > > > You are currently subscribed to cgal-discuss.
> > > > To unsubscribe or access the archives, go to
> > > > https://sympa.inria.fr/sympa/info/cgal-discuss
> > > >
> > > >
> > > >             
> > > --
> > >     Iordan Iordanov, PhD
> > >     https://imiordanov.github.io/
> > >
> > >
> > >
> > > --
> > > You are currently subscribed to cgal-discuss.
> > > To unsubscribe or access the archives, go to
> > > https://sympa.inria.fr/sympa/info/cgal-discuss
> > >
> > >           
> >         
> --
> You are currently subscribed to cgal-discuss.
> To unsubscribe or access the archives, go to
> https://sympa.inria.fr/sympa/info/cgal-discuss
>
>
>

Attachment: P3M3_interface.jpg
Description: JPEG image