From 1e92bee7fc3c0aecbdea719d772901e31cfa6973 Mon Sep 17 00:00:00 2001
From: Jonathan Youett <youett@math.fu-berlin.de>
Date: Wed, 15 Feb 2017 14:10:01 +0100
Subject: [PATCH] Fix case when the grids overlap a little

---
 dune/contact/projections/closestpointprojection.hh | 13 +++++++------
 1 file changed, 7 insertions(+), 6 deletions(-)

diff --git a/dune/contact/projections/closestpointprojection.hh b/dune/contact/projections/closestpointprojection.hh
index 08c6d2ff..f42dea44 100644
--- a/dune/contact/projections/closestpointprojection.hh
+++ b/dune/contact/projections/closestpointprojection.hh
@@ -65,7 +65,7 @@ void ClosestPointProjection<BoundaryPatchType>::project(const BoundaryPatchType&
                                                         const BoundaryPatchType& mortar,
                                                         const ctype couplingDist)
 {
-    ctype eps = 1e-5;
+    ctype eps = 1e-11;
 
     // initialize
     this->obsDirections_.resize(nonmortar.numVertices());
@@ -131,19 +131,20 @@ void ClosestPointProjection<BoundaryPatchType>::project(const BoundaryPatchType&
         // if we found a feasible closest point save the distance and
         if (bestTri != -1) {
             this->obstacles_[localIdx] = dist;
-
-            // if obstacle distance is very small, use minus the outer mortar normal as direction
-            if (dist<1e-5) {
+            // if there is overlap, fix orientation of the direction
+            if (dist<this->overlap_ and (projDirection*mortarSegments[bestTri].unitOuterNormal)>eps) {
+                this->obsDirections_[localIdx] = projDirection;
+                this->obsDirections_[localIdx] *= -1;
+            // if obstacle distance is very small, use the inner mortar normal as direction
+            } else if (dist<eps) {
                 this->obsDirections_[localIdx] = mortarSegments[bestTri].unitOuterNormal;
                 this->obsDirections_[localIdx] *= -1;
             } else
                 this->obsDirections_[localIdx] = projDirection;
 
         } else {
-
             this->obstacles_[localIdx] = std::numeric_limits<ctype>::max();
             this->obsDirections_[localIdx] = averagedNonmortarNormals[globalIdx];
-
         }
     }
 }
-- 
GitLab