From cb821ed5250e27e4c664f8d76fa03d55e2036408 Mon Sep 17 00:00:00 2001 From: David Seus <david.seus@ians.uni-stuttgart.de> Date: Mon, 20 May 2019 11:44:38 +0200 Subject: [PATCH] fix marking of vertical interfaces --- LDDsimulation/boundary_and_interface.py | 18 +++++++++++------- 1 file changed, 11 insertions(+), 7 deletions(-) diff --git a/LDDsimulation/boundary_and_interface.py b/LDDsimulation/boundary_and_interface.py index 12b3ec2..af44457 100644 --- a/LDDsimulation/boundary_and_interface.py +++ b/LDDsimulation/boundary_and_interface.py @@ -171,17 +171,21 @@ class BoundaryPart(df.SubDomain): ymin = min(p1[1], p2[1]) ymax = max(p1[1], p2[1]) # print(f"test if point {p} is on line segment between {p1} or {p2}") - # equality check - if np.absolute((p[0] - xmax)*(p[0] - xmin)) < tol and np.absolute((p[1] - ymax)*(p[1] - ymin)) < tol: + # check if p == p1 or p == p2 + if np.fabs((p[0] - xmax)*(p[0] - xmin)) < tol and np.fabs((p[1] - ymax)*(p[1] - ymin)) < tol: #print(f"point {p} is close to either {p1} or {p2}") return True # check there holds p1[0] < p[0] < p2[0]. If not, p cannot be on the line segment - - # the test ((p[0] - xmax)*(p[0] - xmin) < 0) might fail if p[0] and one of p1[0] - # or p2[0] are equal. In pp1 or pp2 is vertical. We need to still calculate - # the distance to the line segment in this case - if ((p[0] - xmax)*(p[0] - xmin) < tol) or np.absolute((p[0] - xmax)*(p[0] - xmin)) < tol: + # same needs to be done for p1[1] < p[1] < p2[1] + + # the test ((p[0] - xmax)*(p[0] - xmin) < 0) might fail if p[0] and one + # of p1[0] or p2[0] are equal. In this case the line pp1 or pp2 is + # vertical. We need to still calculate the distance to the line segment + # in this case. Same is true for ((p[1] - ymax)*(p[1] - ymin) < 0) + x_is_in_range = ((p[0] - xmax)*(p[0] - xmin) < tol) or np.fabs((p[0] - xmax)*(p[0] - xmin)) < tol + y_is_in_range = ((p[1] - ymax)*(p[1] - ymin) < tol) or np.fabs((p[1] - ymax)*(p[1] - ymin)) < tol + if x_is_in_range and y_is_in_range: # here we have a chance to actually lie on the line segment. segment_direction = (p2 - p1)/np.linalg.norm(p2 - p1) distance = (p - p1) - np.dot(segment_direction, p - p1)*segment_direction -- GitLab