From c90c4291b87087c3a8de72df68e7190deedbb6f2 Mon Sep 17 00:00:00 2001
From: Christoph Lehmann <christoph.lehmann@ufz.de>
Date: Mon, 25 Apr 2022 08:31:20 +0200
Subject: [PATCH] [MeL] Surface normal computation also for lines

---
 MeshLib/Elements/Utils.h | 49 ++++++++++++++++++++++++++++++----------
 1 file changed, 37 insertions(+), 12 deletions(-)

diff --git a/MeshLib/Elements/Utils.h b/MeshLib/Elements/Utils.h
index 6fa00f580cb..f4eb2d25c4c 100644
--- a/MeshLib/Elements/Utils.h
+++ b/MeshLib/Elements/Utils.h
@@ -46,19 +46,44 @@ inline Eigen::Vector3d calculateNormalizedSurfaceNormal(
     MeshLib::Element const& bulk_element)
 {
     Eigen::Vector3d surface_element_normal;
-    if (surface_element.getDimension() < 2)
-    {
-        auto const bulk_element_normal =
-            MeshLib::FaceRule::getSurfaceNormal(bulk_element);
-        auto const& v0 = surface_element.getNode(0)->asEigenVector3d();
-        auto const& v1 = surface_element.getNode(1)->asEigenVector3d();
-        Eigen::Vector3d const edge_vector = v1 - v0;
-        surface_element_normal = bulk_element_normal.cross(edge_vector);
-    }
-    else
+
+    switch (surface_element.getDimension())
     {
-        surface_element_normal =
-            MeshLib::FaceRule::getSurfaceNormal(surface_element);
+        case 2:
+            surface_element_normal =
+                MeshLib::FaceRule::getSurfaceNormal(surface_element);
+            break;
+        case 1:
+        {
+            auto const bulk_element_normal =
+                MeshLib::FaceRule::getSurfaceNormal(bulk_element);
+            auto const& v0 = surface_element.getNode(0)->asEigenVector3d();
+            auto const& v1 = surface_element.getNode(1)->asEigenVector3d();
+            Eigen::Vector3d const edge_vector = v1 - v0;
+            surface_element_normal = -bulk_element_normal.cross(edge_vector);
+            break;
+        }
+        case 0:
+        {
+            assert(surface_element.getCellType() == CellType::POINT1);
+            assert(bulk_element.getCellType() == CellType::LINE2 ||
+                   bulk_element.getCellType() == CellType::LINE3);
+
+            auto const& x = surface_element.getNode(0)->asEigenVector3d();
+
+            // The start of the line element.
+            auto const& a = bulk_element.getNode(0)->asEigenVector3d();
+
+            // The end of the line element is the 2nd base node of the line,
+            // which is the 2nd node.
+            auto const& b = bulk_element.getNode(1)->asEigenVector3d();
+
+            // x coincides either with the start or with the end of the line.
+            // a + b - 2 * x evaluates to a - b or to b - a, respectively.
+            // The formula assumes that the line is perfectly straight, even in
+            // the LINE3 case.
+            surface_element_normal = a + b - 2 * x;
+        }
     }
 
     surface_element_normal.normalize();
-- 
GitLab