From aeb1546f5d801e89a7b1143f1a504c3d6804abd2 Mon Sep 17 00:00:00 2001 From: "Samuel F. Potter" Date: Mon, 19 Jul 2021 10:08:07 -0400 Subject: [PATCH 1/7] enable use of insert_points in meshpy.tet.tetrahedralize (wip) --- meshpy/tet.py | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/meshpy/tet.py b/meshpy/tet.py index 2e86641..808ebce 100644 --- a/meshpy/tet.py +++ b/meshpy/tet.py @@ -141,7 +141,7 @@ def __init__(self, switches, **kwargs): setattr(self, k, v) -def tetrahedralize(mesh_info, options): +def tetrahedralize(mesh_info, options, insert_points=None): mesh = MeshInfo() # restore "C" locale--otherwise tetgen might mis-parse stuff like "a0.01" @@ -155,7 +155,7 @@ def tetrahedralize(mesh_info, options): locale.setlocale(locale.LC_NUMERIC, "C") try: - internals.tetrahedralize(options, mesh_info, mesh) + internals.tetrahedralize(options, mesh_info, mesh, insert_points) finally: # restore previous locale if we've changed it if have_locale: @@ -171,6 +171,8 @@ def build(mesh_info, options=Options("pq"), verbose=False, options.quiet = 1 if insert_points is not None: + if not isinstance(insert_points, MeshInfo): + raise ValueError('insert_points should be an instance of MeshInfo') options.insertaddpoints = 1 if attributes: @@ -183,4 +185,4 @@ def build(mesh_info, options=Options("pq"), verbose=False, if diagnose: options.diagnose = 1 - return tetrahedralize(mesh_info, options) + return tetrahedralize(mesh_info, options, insert_points) From e1652ea1e800d256d9c974196a733a3c71377bf0 Mon Sep 17 00:00:00 2001 From: "Samuel F. Potter" Date: Mon, 19 Jul 2021 10:37:52 -0400 Subject: [PATCH 2/7] add test_insert_points.py example --- examples/test_insert_points.py | 44 ++++++++++++++++++++++++++++++++++ 1 file changed, 44 insertions(+) create mode 100644 examples/test_insert_points.py diff --git a/examples/test_insert_points.py b/examples/test_insert_points.py new file mode 100644 index 0000000..9edeaed --- /dev/null +++ b/examples/test_insert_points.py @@ -0,0 +1,44 @@ +import logging +import numpy as np + +from meshpy.tet import MeshInfo, Options, build + +if __name__ == '__main__': + logger = logging.getLogger('test_insert_points.py') + + points = [(0, 0, 0), + (0, 0, 1), + (1, 0, 0), + (1, 0, 1), + (1, 1, 0), + (1, 1, 1), + (0, 1, 0), + (0, 1, 1)] + + facets = [(0, 1, 3, 2), + (2, 3, 5, 4), + (4, 5, 7, 6), + (6, 7, 1, 0), + (0, 2, 4, 6), + (1, 3, 5, 7)] + + mesh_info = MeshInfo() + mesh_info.set_points(points) + mesh_info.set_facets(facets) + + options = Options(switches='', plc=True, verbose=True, quiet=False) + + # Insert an interior point of the cube as a constrained point + interior_point = (0.33, 0.7, 0.91) + + insert_points_mesh_info = MeshInfo() + insert_points_mesh_info.set_points([interior_point]) + + mesh = build(mesh_info, options, max_volume=0.1, + insert_points=insert_points_mesh_info) + + mesh_points = np.array(mesh.points) + min_dist = np.sqrt(np.sum((interior_point - mesh_points))**2).min() + + if min_dist > 0: + logger.error('tetrahedron mesh does not contain contrained point') From 023ce61beea0d3e9cab1f86c5beaa2ef156e0aa0 Mon Sep 17 00:00:00 2001 From: "Samuel F. Potter" Date: Mon, 19 Jul 2021 10:37:52 -0400 Subject: [PATCH 3/7] add test_insert_points.py example --- examples/test_insert_points.py | 44 ++++++++++++++++++++++++++++++++++ 1 file changed, 44 insertions(+) create mode 100644 examples/test_insert_points.py diff --git a/examples/test_insert_points.py b/examples/test_insert_points.py new file mode 100644 index 0000000..e6abbe3 --- /dev/null +++ b/examples/test_insert_points.py @@ -0,0 +1,44 @@ +import logging +import numpy as np + +from meshpy.tet import MeshInfo, Options, build + +if __name__ == '__main__': + logger = logging.getLogger('test_insert_points.py') + + points = [(0, 0, 0), + (0, 0, 1), + (1, 0, 0), + (1, 0, 1), + (1, 1, 0), + (1, 1, 1), + (0, 1, 0), + (0, 1, 1)] + + facets = [(0, 1, 3, 2), + (2, 3, 5, 4), + (4, 5, 7, 6), + (6, 7, 1, 0), + (0, 2, 4, 6), + (1, 3, 5, 7)] + + mesh_info = MeshInfo() + mesh_info.set_points(points) + mesh_info.set_facets(facets) + + options = Options(switches='', plc=True, verbose=True, quiet=False) + + # Insert an interior point of the cube as a constrained point + interior_point = (0.33, 0.7, 0.91) + + insert_points_mesh_info = MeshInfo() + insert_points_mesh_info.set_points([interior_point]) + + mesh = build(mesh_info, options, max_volume=0.1, + insert_points=insert_points_mesh_info) + + mesh_points = np.array(mesh.points) + min_dist = np.sqrt(np.sum((interior_point - mesh_points), axis=1)**2).min() + + if min_dist > 0: + logger.error('tetrahedron mesh does not contain contrained point') From d637c6a70e18d5ea062d9e906ba7cb8b91563988 Mon Sep 17 00:00:00 2001 From: "Samuel F. Potter" Date: Mon, 19 Jul 2021 10:55:17 -0400 Subject: [PATCH 4/7] don't use Options's keywords --- examples/test_insert_points.py | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/examples/test_insert_points.py b/examples/test_insert_points.py index e6abbe3..b574539 100644 --- a/examples/test_insert_points.py +++ b/examples/test_insert_points.py @@ -26,15 +26,13 @@ mesh_info.set_points(points) mesh_info.set_facets(facets) - options = Options(switches='', plc=True, verbose=True, quiet=False) - # Insert an interior point of the cube as a constrained point interior_point = (0.33, 0.7, 0.91) insert_points_mesh_info = MeshInfo() insert_points_mesh_info.set_points([interior_point]) - mesh = build(mesh_info, options, max_volume=0.1, + mesh = build(mesh_info, max_volume=0.25, insert_points=insert_points_mesh_info) mesh_points = np.array(mesh.points) From 4ce94047e2f674b13ddd9cfbbd84c45b57051ec4 Mon Sep 17 00:00:00 2001 From: "Samuel F. Potter" Date: Mon, 19 Jul 2021 12:22:47 -0400 Subject: [PATCH 5/7] add test for insert_points --- test/test_meshpy.py | 50 +++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 50 insertions(+) diff --git a/test/test_meshpy.py b/test/test_meshpy.py index 52c59d7..187a4c7 100644 --- a/test/test_meshpy.py +++ b/test/test_meshpy.py @@ -162,6 +162,56 @@ def test_tetgen_points(): #mesh.write_vtk("test.vtk") + +def test_tetgen_insert_points(): + import numpy as np + + from meshpy.tet import MeshInfo, Options, build + + points = [(0, 0, 0), + (0, 0, 1), + (1, 0, 0), + (1, 0, 1), + (1, 1, 0), + (1, 1, 1), + (0, 1, 0), + (0, 1, 1)] + + facets = [(0, 1, 3, 2), + (2, 3, 5, 4), + (4, 5, 7, 6), + (6, 7, 1, 0), + (0, 2, 4, 6), + (1, 3, 5, 7)] + + mesh_info = MeshInfo() + mesh_info.set_points(points) + mesh_info.set_facets(facets) + + interior_point = (0.33, 0.7, 0.91) + + insert_points_mesh_info = MeshInfo() + insert_points_mesh_info.set_points([interior_point]) + + mesh = build(mesh_info, max_volume=0.25, + insert_points=insert_points_mesh_info) + + mesh_points = np.array(mesh.points) + + # because of the max_volume==0.25 constraint, TetGen should have + # inserted some points in addition to the 8 cube vertices and the + # constranied interior point + assert mesh_points.shape[0] > 9 + + interior_point_index = \ + np.where((mesh_points == np.array(interior_point)).all(1))[0][0] + + mesh_elements = np.array(mesh.elements, dtype=int) + + # make sure that there cells including the constrained point have + # been added to the mesh + assert (mesh_elements == interior_point_index).any() + # }}} From 04e71d53b6daad155d9a6d0f6e9b24db39ef71de Mon Sep 17 00:00:00 2001 From: "Samuel F. Potter" Date: Mon, 19 Jul 2021 12:47:28 -0400 Subject: [PATCH 6/7] fix flake8 failures --- examples/test_insert_points.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/examples/test_insert_points.py b/examples/test_insert_points.py index b574539..cd35ea1 100644 --- a/examples/test_insert_points.py +++ b/examples/test_insert_points.py @@ -1,10 +1,10 @@ import logging import numpy as np -from meshpy.tet import MeshInfo, Options, build +from meshpy.tet import MeshInfo, build -if __name__ == '__main__': - logger = logging.getLogger('test_insert_points.py') +if __name__ == "__main__": + logger = logging.getLogger("test_insert_points.py") points = [(0, 0, 0), (0, 0, 1), @@ -39,4 +39,4 @@ min_dist = np.sqrt(np.sum((interior_point - mesh_points), axis=1)**2).min() if min_dist > 0: - logger.error('tetrahedron mesh does not contain contrained point') + logger.error("tetrahedron mesh does not contain contrained point") From 2ce10c0bd1d9c1d5dd1cdf06f11257b929978f81 Mon Sep 17 00:00:00 2001 From: "Samuel F. Potter" Date: Mon, 19 Jul 2021 16:09:28 -0400 Subject: [PATCH 7/7] fix flake8 failure in test/test_meshpy.py --- test/test_meshpy.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/test_meshpy.py b/test/test_meshpy.py index 187a4c7..e1c1651 100644 --- a/test/test_meshpy.py +++ b/test/test_meshpy.py @@ -166,7 +166,7 @@ def test_tetgen_points(): def test_tetgen_insert_points(): import numpy as np - from meshpy.tet import MeshInfo, Options, build + from meshpy.tet import MeshInfo, build points = [(0, 0, 0), (0, 0, 1),