0

Currently, I hava a mesh: points and faces, and I want to optimized it with gmsh. So, the first problem is to input my mesh to gmsh.

The points is a nx3 data, which indicates the coordinate. The faces is also a mx3 data, which indicates the connection of points.

After learning the gmsh tutorial, my code to input points/faces to gmsh is:

nodes = list(range(1, len(points) + 1))
coords = []
for p in points:
    coords.extend(p)

gmsh.model.addDiscreteEntity(2, 0, [])
gmsh.model.mesh.addNodes(2, 0, nodes, coords)
elementType = 2
elementTags = list(range(len(faces)))
nodeTags = []
for t in faces:
    tt = [idx + 1 for idx in t]
    nodeTags.extend(tt)
gmsh.model.mesh.addElementsByType(0, elementType, elementTags, nodeTags)

However, the result is different between that loading from points/faces and that loading from stl file.

For loading mesh from points/faces, the generated mesh has an unexpected hole:

enter image description here

The code to reproduce the result is:

import gmsh
import sys
import os
import math
import numpy as np

# the function to obtain points/faces from gmsh
def GetPointFace():
    idx, points, _ = gmsh.model.mesh.getNodes()
    idx -= 1
    idx = idx.tolist()

    m = {}
    for i, v in enumerate(idx):
        m[v] = i

    points = np.asarray(points).reshape(-1, 3)

    elem_types, elem_tags, node_tags = gmsh.model.mesh.getElements()
    idx = elem_types.tolist().index(2)
    elem_type, elem_tag, node_tag = elem_types[idx], elem_tags[idx], node_tags[idx]

    num_nodes_per_cell = gmsh.model.mesh.getElementProperties(elem_type)[3]
    node_tags_reshaped = np.asarray(node_tag).reshape(-1, num_nodes_per_cell) - 1
    node_tags_sorted = node_tags_reshaped[np.argsort(elem_tag)]

    tris = []
    for i, tri in enumerate(node_tags_sorted):
        p0 = tri[0]
        p1 = tri[1]
        p2 = tri[2]
        tris.append([m[p0], m[p1], m[p2]])
    return points, tris

gmsh.initialize(sys.argv)

path = os.path.dirname(os.path.abspath(__file__))
gmsh.merge(os.path.join(path, 'aneurysm_data.stl'))

points, faces = GetPointFace() ## get the aneurysm data from stl, including points and the faces


###############################################################################################
## comment the following code would get the correct mesh
## input points/faces to gmsh
gmsh.clear()
nodes = list(range(1, len(points) + 1))
coords = []
for p in points:
    coords.extend(p)

gmsh.model.addDiscreteEntity(2, 0, [])
gmsh.model.mesh.addNodes(2, 0, nodes, coords)
elementType = 2
elementTags = list(range(len(faces)))
nodeTags = []
for t in faces:
    tt = [idx + 1 for idx in t]
    nodeTags.extend(tt)
gmsh.model.mesh.addElementsByType(0, elementType, elementTags, nodeTags)
##################################################################################################
gmsh.model.mesh.classifySurfaces(math.pi, True, True)
gmsh.model.mesh.createGeometry()

gmsh.model.geo.synchronize()

gmsh.option.setNumber('Mesh.Algorithm', 1)
gmsh.option.setNumber('Mesh.MeshSizeMin', 1)
gmsh.option.setNumber('Mesh.MeshSizeMax', 1)
gmsh.option.setNumber('Mesh.MeshSizeFactor', 1)

gmsh.model.mesh.generate(2)

## display the result
points, faces = GetPointFace()
from vedo import *
mesh = Mesh([points, faces]).c('y').bc('c')
show(mesh)

The test data is: https://gitlab.onelab.info/gmsh/gmsh/-/blob/master/examples/api/aneurysm_data.stl

Any suggestion is appreciated~~~

Qiang Zhang
  • 820
  • 8
  • 32

0 Answers0