21
21
Mesh
22
22
"""
23
23
from __future__ import annotations
24
+ import meshio
24
25
import numpy
25
26
26
- import meshio
27
27
28
28
__all__ = ['Mesh' ]
29
29
@@ -68,29 +68,31 @@ class Mesh:
68
68
N_BC_triangles (dict[str, int]): the number of triangles on the boundary of the domain associated to each
69
69
boundary label in self.BC_labels. For example self.N_BC_triangles['my_label'] returns the number of boundary
70
70
triangles associated to lable 'my_label'. 'my_label' must be a key of self.BC_labels.
71
- vertices (numpy.ndarray): An (self.N_vertices x M) array containing the M coordinates of the self.N_vertices vertices that make up
72
- the mesh. M specifies the geometric dimension of the mesh, such that the mesh describes an M-dimensional
73
- domain.
71
+ vertices (numpy.ndarray): An (self.N_vertices x M) array containing the M coordinates of the self.N_vertices
72
+ vertices that make up the mesh. M specifies the geometric dimension of the mesh, such that the mesh
73
+ describes an M-dimensional domain.
74
74
tets (numpy.ndarray): An (self.N_tets x 4) array containing the 4 indices of the vertices of the self.N_tets
75
75
tetrahedra that make up the mesh.
76
76
BC_triangles (dict[str, numpy.ndarray]):
77
77
78
78
Example:
79
79
An element of this class can be initialized in the following way
80
80
81
- >>> import edg_acoustics
82
- >>> BC_labels = {'CNRBC': 12, 'slip': 11, 'impedance': 33}
83
- >>> filename = "path_to_my_mesh_file/mesh_filename.msh"
84
- >>> mesh = edg_acoustics.Mesh(filename, BC_labels)
81
+ >>> import edg_acoustics
82
+ >>> BC_labels = {'slip': 11, 'impedance1': 13, 'impedance2': 14, 'impedance3': 15}
83
+ >>> filename = "../data/tests/mesh/CoarseMesh.msh"
84
+ >>> mesh = edg_acoustics.Mesh(filename, BC_labels)
85
+ <BLANKLINE>
86
+ >>> mesh.N_BC_triangles
87
+ {'slip': 5347, 'impedance1': 400, 'impedance2': 3576, 'impedance3': 3294}
88
+
85
89
"""
86
90
87
91
def __init__ (self , filename : str , BC_labels : dict [str , int ]):
88
92
# Init from file
89
93
self .__init_from_file (filename , BC_labels )
90
94
91
95
def __init_from_file (self , filename : str , BC_labels : dict [str , int ]):
92
- print ("Init from filename" )
93
-
94
96
# Load mesh data from mesh file
95
97
mesh_data = meshio .read (filename )
96
98
@@ -101,16 +103,16 @@ def __init_from_file(self, filename: str, BC_labels: dict[str, int]):
101
103
self .vertices = mesh_data .points
102
104
103
105
# Check if all labels provided as input exist in the mesh data and vice versa, if not, raise error
104
- BC_labels_in_mesh = sorted (numpy .unique (mesh_data .cell_data_dict ['gmsh:physical' ]['triangle' ])) # get all labels in the mesh, sort for faster comparison below
106
+ BC_labels_in_mesh = sorted (numpy .unique (mesh_data .cell_data_dict ['gmsh:physical' ]['triangle' ])) # get labels in the mesh, sort for faster comparison below
105
107
BC_labels_in_input = sorted (BC_labels .values ()) # get all labels specified in input
106
108
if not BC_labels_in_input == BC_labels_in_mesh :
107
109
raise ValueError (
108
110
"[edg_acoustics.Mesh] All BC labels must be present in the mesh and all labels in the mesh must be "
109
111
"present in BC_labels." )
110
112
111
113
# Read the boundary triangles and their definitions separating them into the different boundary condition labels
112
- self .N_BC_triangles = dict ()
113
- self .BC_triangles = dict ()
114
+ self .N_BC_triangles = {}
115
+ self .BC_triangles = {}
114
116
for BC_label in BC_labels :
115
117
triangles_have_label = (mesh_data .cell_data_dict ['gmsh:physical' ]['triangle' ] == BC_labels [BC_label ]) # array with bools specifying if triangle has BC_label or not
116
118
self .N_BC_triangles [BC_label ] = triangles_have_label .sum () # number of triangles with label BC_label
@@ -123,6 +125,23 @@ def __init_from_file(self, filename: str, BC_labels: dict[str, int]):
123
125
# Compute the mesh connectivity
124
126
self .EToE , self .EToF = self .__compute_element_connectivity (self .tets )
125
127
128
+ # Operators --------------------------------------------------------------------------------------------------------
129
+ def __eq__ (self , other ):
130
+ if isinstance (other , type (self )):
131
+ # If self and other are mesh objects then check if all properties contain the same data
132
+ are_equal = (self .N_vertices == other .N_vertices
133
+ and self .N_tets == other .N_tets
134
+ and self .N_BC_triangles == other .N_BC_triangles
135
+ and numpy .array_equal (self .vertices , other .vertices )
136
+ and numpy .array_equal (self .tets , other .tets )
137
+ and all (numpy .array_equal (item , other .BC_triangles [key ]) for key , item in self .BC_triangles .items ()))
138
+ else :
139
+ # If they are of different types, then they are not the same
140
+ are_equal = False
141
+
142
+ return are_equal
143
+
144
+ # ------------------------------------------------------------------------------------------------------------------
126
145
127
146
# Static methods ---------------------------------------------------------------------------------------------------
128
147
@staticmethod
@@ -156,7 +175,7 @@ def __compute_element_connectivity(tets: numpy.ndarray):
156
175
# Get information on the number of faces, tets, and vertices
157
176
N_faces_per_tet = 4 # number of faces per element
158
177
N_tets = tets .shape [0 ] # number of elements in the mesh
159
- N_vertices = tets .max () + 1 # number of vertices in the mesh, +1 because indexing starts at 0
178
+ # N_vertices = tets.max() + 1 # number of vertices in the mesh, +1 because indexing starts at 0
160
179
161
180
# Create a unique identifier for each face based on the vertices that make up the face
162
181
# the order of the vertices does not matter.
@@ -201,10 +220,8 @@ def __compute_element_connectivity(tets: numpy.ndarray):
201
220
# We now sort the face_ids so that we have the identical faces next to each other
202
221
face_ids_sort_indices = numpy .argsort (face_ids ) # get the ordering that sorts the face_ids
203
222
face_ids = face_ids [face_ids_sort_indices ] # sort the face ids
204
- face_vertices = face_vertices [face_ids_sort_indices , :] # reorder the faces so that their corresponding
205
- # face_ids are sorted
206
- face_vertices_idx = face_vertices_idx [face_ids_sort_indices ] # reorder the face indices so that their
207
- # corresponding face_ids are sorted
223
+ face_vertices = face_vertices [face_ids_sort_indices , :] # reorder the faces so that their corresponding face_ids are sorted
224
+ face_vertices_idx = face_vertices_idx [face_ids_sort_indices ] # reorder the face indices so that their corresponding face_ids are sorted
208
225
209
226
# Find the indices of face_ids of all interior faces, i.e., that are shared by two elements
210
227
# i.e., faces that appear twice (one time for each element)
@@ -257,5 +274,3 @@ def EToV(self):
257
274
tetrahedra that make up the mesh. It returns the value in self.tets, since it is the same data."""
258
275
return self .tets
259
276
# ------------------------------------------------------------------------------------------------------------------
260
-
261
-
0 commit comments