Generating periodic meshes by mirroring


As far as I am aware, we can only generate periodic meshes in Netgen for geometries of class CSGeometry. I did not find a function equivalent to the PeriodicSurfaces() function in the STLGeometry class. As a work around, I am now mirroring the mesh generated from the STL geometry. I have two concerns with this approach:
[li]Will Netgen/NGSolve know that a mesh generated in this way will be periodic? The documentation for the class Periodic mentions that the mesh needs to be periodic in order to map periodic degrees of freedom from minion boundaries to master boundaries. So I guess the function PeriodicSurfaces marks the master and minion surfaces, and this is not done when generating a periodic mesh by mirroring. [/li]
[li]On mirroring, all boundaries of the initial mesh seem to still exist as boundaries even though some of them are now inside the geometry and are no longer external surfaces. Is there a way to remove these boundaries by for instance merging duplicate nodes and surface elements in the mesh? If the boundaries are not somehow removed, I guess NGSolve will apply a zero gradient von Neumann boundary condition by default on these boundaries, which I do not want in most cases.[/li]

Thanks in advance!



Edit: I found a function IdentifyPeriodicBoundaries() in the netgen mesh class, but since mirrored surfaces will be listed in the same boundary, I am not sure if it is possible to use this function to specify the faces which are to be tested for periodicity of the mesh.

Hm yes, its quite hacky, but you can do this:

Mirror the mesh, as long as you define no integrators on the interior boundaries nothing bad happens. The neumann boundaries come from the boundary term when deriving the weak formulation which cancel out on interior boundaries.
Then you have to set the surface elements of one of the boundaries to a different index, something like this should work:

from netgen.csg import unit_cube
from netgen.meshing import *

mesh = unit_cube.GenerateMesh(maxh=0.2)
mesh = mesh.Mirror((0,0,0),(-1,0,0))

new_index = mesh.GetNFaceDescriptors() + 1
mesh.Add(FaceDescriptor(new_index, 1, 0, new_index))
mesh.SetBCName(new_index-1, "new_bc")

for el in mesh.Elements2D():
    if all([abs(mesh[p][0]+1) < 1e-10 for p in el.points]): # check somehow if on mirrored boundary
        # print("points = ", [mesh[p] for p in el.points])
        el.index = new_index

mapping = Trafo((-2, 0, 0))
mesh.IdentifyPeriodicBoundaries("front", "new_bc", mapping)

from ngsolve import *
from ngsolve.internal import *
viewoptions.drawidentified = 1
viewoptions.drawfilledtrigs = 0

mesh = Mesh(mesh)

Hello Christopher, thank you for your reply. I think I understand your code. I will try to implement it in my script and get back to you next week.

Is it possible to remove the boundaries that are inside the mesh by merging duplicated nodes/elements? I would also like to use the mesh in other programs like OpenFOAM, and in OpenFOAM I have to assign a boundary condition to every “external surface”. I guess one solution would be to copy all 3D elements and only the 2D elements I need to a new mesh object, but is there a simpler solution?

Hello Christopher,

I was able to implement everything from your reply in my script, and it seems to work fine. Thank you once again.

I tried to copy just the elements I want into a new mesh object and export the mesh to OpenFOAM, but I was not able to load the mesh in OpenFOAM. NGSolve had also given me a warning when I used the Draw(mesh) command that some “elements are not matching”. Do you have any suggestions?

I am attaching my script and a geometry file with this post


Attachment: test.stl