From a1595453de919b73cb637d35f1dd8f59cc300b3d Mon Sep 17 00:00:00 2001 From: Scott Sibole Date: Fri, 22 Jul 2016 14:41:09 -0600 Subject: [PATCH] Changed tetrahedral mesh generation to make use of tetmesh flagging boundary nodes. --- src/pyCellAnalyst/CellMech.py | 15 ++++++--------- src/pyCellAnalyst/FEA_GUI.py | 1 - 2 files changed, 6 insertions(+), 10 deletions(-) diff --git a/src/pyCellAnalyst/CellMech.py b/src/pyCellAnalyst/CellMech.py index dfc1ecf..019ef93 100644 --- a/src/pyCellAnalyst/CellMech.py +++ b/src/pyCellAnalyst/CellMech.py @@ -172,6 +172,7 @@ def _readstls(self): massProps = vtk.vtkMassProperties() massProps.SetInputData(self.rsurfs[-1]) massProps.Update() + print("Generating tetrahedral mesh from {:s}".format(fname)) self._make3Dmesh( str(os.path.normpath(self._ref_dir + os.sep + fname)), 'MATERIAL', massProps.GetVolume()) @@ -189,6 +190,7 @@ def _readstls(self): massProps = vtk.vtkMassProperties() massProps.SetInputData(self.dsurfs[-1]) massProps.Update() + print("Generating tetrahedral mesh from {:s}".format(fname)) self._make3Dmesh( str(os.path.normpath(self._def_dir + os.sep + fname)), 'SPATIAL', massProps.GetVolume()) @@ -473,14 +475,14 @@ def _make3Dmesh(self, filename, frame, vConst): * dvols * daxes """ - vConst /= 30000.0 + vConst /= 50000.0 edgeSize = (vConst*12/np.sqrt(2)) ** (1./3.) m = mesh.Mesher(inputname=filename, outputname="tmp.vtu", facetAngle=30.0, facetDistance=0.1, edgeLength=edgeSize, - edgeRatio=1.5) + edgeRatio=1.3) m.makeMesh() gridReader = vtk.vtkXMLUnstructuredGridReader() @@ -490,13 +492,8 @@ def _make3Dmesh(self, filename, frame, vConst): nodes = vtk_to_numpy(vtkMesh.GetPoints().GetData()) elements = vtk_to_numpy(vtkMesh.GetCells().GetData()).reshape( (vtkMesh.GetNumberOfCells(), 5))[:,1:] - extractSurface = vtk.vtkDataSetSurfaceFilter() - extractSurface.SetInputData(vtkMesh) - extractSurface.Update() - faces = extractSurface.GetOutput() - N = faces.GetNumberOfPolys() - faces = vtk_to_numpy(faces.GetPolys().GetData()).reshape((N, 4))[:,1:] - s_nodes = np.unique(np.ravel(faces)) + vertex_type = vtk_to_numpy(vtkMesh.GetPointData().GetArray("Vertex Type")) + s_nodes = np.argwhere(vertex_type==1).ravel() n1 = nodes[elements[:, 0], :] n2 = nodes[elements[:, 1], :] diff --git a/src/pyCellAnalyst/FEA_GUI.py b/src/pyCellAnalyst/FEA_GUI.py index 58804e3..b7c0382 100644 --- a/src/pyCellAnalyst/FEA_GUI.py +++ b/src/pyCellAnalyst/FEA_GUI.py @@ -405,7 +405,6 @@ def performFEA(self): mesh = febio.MeshDef() # would be good to vectorize these - print self.data[filename]['elements'] for i, e in enumerate(self.data[filename]['elements']): mesh.elements.append(['tet4', i + 1] + e)