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:
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~~~