From e0dae182e2840148d3cef7aed49c162b941e98f3 Mon Sep 17 00:00:00 2001 From: Keith Roberts Date: Sun, 18 Apr 2021 16:17:56 -0300 Subject: [PATCH] fixing points in iterative laplacian2 (#210) * fixing points in iterative laplacian2 * update change log Co-authored-by: Keith Roberts --- README.md | 3 +++ SeismicMesh/geometry/utils.py | 19 ++++++++++++++++++- 2 files changed, 21 insertions(+), 1 deletion(-) diff --git a/README.md b/README.md index 2e59e93..912c327 100644 --- a/README.md +++ b/README.md @@ -823,6 +823,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ### Added - Smoothed sets (e.g., intersections, differences, and unions) - Conversion of velocity data from feet-second to meters-second +- Support for fixed points in iterative Laplacian mesh smoother. ### Improved - Simplified pybind11 build system. @@ -833,6 +834,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - Visuzlization of signed distance functions ### Fixed - Support for Python 3.9 +### Improved +- Fixed points in iterative Laplacian smooth ## [3.4.0]-2021-02-14 ### Added diff --git a/SeismicMesh/geometry/utils.py b/SeismicMesh/geometry/utils.py index b3615a7..600481b 100644 --- a/SeismicMesh/geometry/utils.py +++ b/SeismicMesh/geometry/utils.py @@ -547,7 +547,7 @@ def laplacian2_fixed_point(vertices, entities): return vertices_new, entities -def laplacian2(vertices, entities, max_iter=20, tol=0.01, verbose=1): +def laplacian2(vertices, entities, max_iter=20, tol=0.01, verbose=1, pfix=None): """Move vertices to the average position of their connected neighbors with the goal to hopefully improve geometric entity quality. @@ -559,6 +559,10 @@ def laplacian2(vertices, entities, max_iter=20, tol=0.01, verbose=1): :type max_iter: `int`, optional :param tol: iterations will cease when movement < tol :type tol: `float`, optional + :param verbose: display progress to the screen + :type verbose: `float`, optional + :param pfix: coordinates that you don't wish to move + :type pfix: array-like :return vertices: updated vertices of mesh :rtype: numpy.ndarray[`float` x dim] @@ -568,6 +572,12 @@ def laplacian2(vertices, entities, max_iter=20, tol=0.01, verbose=1): if vertices.ndim != 2: raise NotImplementedError("Laplacian smoothing only works in 2D for now") + def _closest_node(node, nodes): + nodes = np.asarray(nodes) + deltas = nodes - node + dist_2 = np.einsum("ij,ij->i", deltas, deltas) + return np.argmin(dist_2) + eps = np.finfo(float).eps n = len(vertices) @@ -580,6 +590,13 @@ def laplacian2(vertices, entities, max_iter=20, tol=0.01, verbose=1): ) bnd = get_boundary_vertices(entities) edge = get_edges(entities) + if pfix is not None: + ifix = [] + for fix in pfix: + ifix.append(_closest_node(fix, vertices)) + ifix = np.asarray(ifix) + bnd = np.concatenate((bnd, ifix)) + W = np.sum(S, 1) if np.any(W == 0): print("Invalid mesh. Disjoint vertices found. Returning", flush=True)