0
$\begingroup$

I am using python API of Gmsh to generate a mesh for a rectangular domain. I am really new at this. My code looks like this,

import sys
import gmsh

gmsh.initialize()

gmsh.model.add("t1")

lc = 1e-1
gmsh.model.geo.addPoint(-1, 0, 0, lc, 1)

gmsh.model.geo.addPoint(1, 0, 0, lc, 2)
gmsh.model.geo.addPoint(1, 1, 0, lc, 3)

gmsh.model.geo.addPoint(-1, 1, 0, lc, 4)

gmsh.model.geo.addLine(1, 2, 1)
gmsh.model.geo.addLine(3, 2, 2)
gmsh.model.geo.addLine(3, 4, 3)
gmsh.model.geo.addLine(4, 1, 4)

gmsh.model.geo.addCurveLoop([4, 1, -2, 3], 1)

gmsh.model.geo.addPlaneSurface([1], 1)
gmsh.model.geo.synchronize()

gmsh.model.addPhysicalGroup(1, [1, 2, 4], 5)
gmsh.model.addPhysicalGroup(2, [1], name = "My surface")

gmsh.model.mesh.generate(2)

gmsh.write("t1.msh")


if '-nopopup' not in sys.argv:
    gmsh.fltk.run()

gmsh.finalize()

I want to impose a mesh regularity condition, that the lines joining the cell centres must be orthogonal to the edges. I cannot find how to do that, any reference or solution would be of help.

Thank you in advance.

Edit: I am primarily interested in triangular meshes and I use the centroids of the mesh elements as cell centres.

$\endgroup$
6
  • 1
    $\begingroup$ How do you define "cell center"? $\endgroup$ Commented Jan 16, 2023 at 16:44
  • $\begingroup$ I use the centroid of the triangular mesh elements as the cell centres. $\endgroup$
    – Mainak
    Commented Jan 16, 2023 at 16:54
  • 2
    $\begingroup$ I think it's unlikely the gmsh or any other unstructured mesh generator will produce a mesh with the regularity you're describing. Generally these are engineered to produce meshes with certain good quality measures, like avoidance of very large or small angles or bounded triangle aspect ratio. There are mesh improvement strategies that can take you closer to the condition you describe but again I don't think you can expect this out of the box. $\endgroup$ Commented Jan 16, 2023 at 19:12
  • 1
    $\begingroup$ I wonder if it is even possible to create such meshes. (But don't know the answer to this.) What happens if you take an arbitrary triangle and you subdivided it into four triangles by using the edge midpoints. Does the result satisfy your constraints? $\endgroup$ Commented Jan 16, 2023 at 20:28
  • 1
    $\begingroup$ Besides regular triangles and rectangles are other tesselations that admit the condition that you want to impose? $\endgroup$
    – nicoguaro
    Commented Jan 16, 2023 at 20:56

1 Answer 1

3
$\begingroup$

If you're prepared to be a little more flexible with your definition of cell "centres", it's possible to achieve the regularity constraints you seek via Delaunay Triangulations and Voronoi Tessellations, specifically Centroidal Voronoi Tessellations, or, even more generally, Centroidal Laguerre-Power Tessellations.

In such tessellations, edges in the (primal) triangulation and (dual) polygonal tessellation are orthogonal, with dual edges spanning between "circumcentres" in adjacent triangles. To form tessellations that are "centroidal" as well as orthogonal, an additional optimisation step is used to adjust the position and topology of the mesh and drive each circumcentre toward its associated centroid:

primal-dual mesh

(In practice, it's hard to force the mesh to be identically centroidal, but CVT optimisers generally do a pretty good job, and any remaining offsets typically manifest as the difference between 2nd- and 1st-order errors).

I'm not sure whether gmsh supports the generation of CVTs, but there are a couple of mesh generators that do. One is my library jigsaw, another option is Bruno Levy's geogram.

$\endgroup$
1
  • $\begingroup$ For my case, the orthogonality of the line joining the cell centres to the edges is paramount. The choice of the cell centre being the centroid is for convenience and can be modified with little effort. Thanks, for pointing me in the right direction! $\endgroup$
    – Mainak
    Commented Jan 18, 2023 at 9:08

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.