From 7781656347df88ba4c732b2bd4eba37f50eb3cac Mon Sep 17 00:00:00 2001
From: Henry Weller <http://cfd.direct>
Date: Fri, 4 Nov 2016 17:13:04 +0000
Subject: [PATCH] pointConstraint: Added constrainDisplacement

which directly applies the constraint to the displacement without
external tensor ops.

Patch contributed by Mattijs Janssens
---
 .../pointConstraint/pointConstraint.H         |  7 ++-
 .../pointConstraint/pointConstraintI.H        | 46 ++++++++++++++++---
 2 files changed, 44 insertions(+), 9 deletions(-)

diff --git a/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointConstraint/pointConstraint.H b/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointConstraint/pointConstraint.H
index 9696efeed2..ae8235eec3 100644
--- a/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointConstraint/pointConstraint.H
+++ b/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointConstraint/pointConstraint.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -87,6 +87,9 @@ public:
         //- Return the accumulated unconstrained directions. Directions
         //  coded as first n rows of tensor.
         inline void unconstrainedDirections(label& n, tensor& vecs) const;
+
+        //- Constrain a displacement
+        inline vector constrainDisplacement(const vector& disp) const;
 };
 
 
@@ -103,7 +106,7 @@ inline pointConstraint transform(const tensor& tt, const pointConstraint& v);
 //- contiguous
 template<class T> bool contiguous();
 template<>
-inline bool contiguous<pointConstraint>()                {return true;}
+inline bool contiguous<pointConstraint>() {return true;}
 
 
 
diff --git a/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointConstraint/pointConstraintI.H b/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointConstraint/pointConstraintI.H
index 27aa1debf4..0a3d14d2dd 100644
--- a/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointConstraint/pointConstraintI.H
+++ b/src/OpenFOAM/fields/pointPatchFields/pointPatchField/pointConstraint/pointConstraintI.H
@@ -45,7 +45,7 @@ inline Foam::pointConstraint::pointConstraint(Istream& is)
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-void Foam::pointConstraint::applyConstraint(const vector& cd)
+inline void Foam::pointConstraint::applyConstraint(const vector& cd)
 {
     if (first() == 0)
     {
@@ -74,7 +74,7 @@ void Foam::pointConstraint::applyConstraint(const vector& cd)
 }
 
 
-void Foam::pointConstraint::combine(const pointConstraint& pc)
+inline void Foam::pointConstraint::combine(const pointConstraint& pc)
 {
     if (first() == 0)
     {
@@ -115,7 +115,7 @@ void Foam::pointConstraint::combine(const pointConstraint& pc)
 }
 
 
-Foam::tensor Foam::pointConstraint::constraintTransformation() const
+inline Foam::tensor Foam::pointConstraint::constraintTransformation() const
 {
     if (first() == 0)
     {
@@ -136,8 +136,11 @@ Foam::tensor Foam::pointConstraint::constraintTransformation() const
 }
 
 
-void Foam::pointConstraint::unconstrainedDirections(label& n, tensor& tt)
-const
+inline void Foam::pointConstraint::unconstrainedDirections
+(
+    label& n,
+    tensor& tt
+) const
 {
     n = 3-first();
 
@@ -179,7 +182,36 @@ const
 }
 
 
-void Foam::combineConstraintsEqOp::operator()
+inline Foam::vector Foam::pointConstraint::constrainDisplacement
+(
+    const vector& d
+) const
+{
+    vector cd;
+
+    if (first() == 0)
+    {
+        cd = d;
+    }
+    else if (first() == 1)
+    {
+        // Remove plane normal
+        cd = d-(d&second())*second();
+    }
+    else if (first() == 2)
+    {
+        // Keep line direction only
+        cd = (d&second())*second();
+    }
+    else
+    {
+        cd = Zero;
+    }
+    return cd;
+}
+
+
+inline void Foam::combineConstraintsEqOp::operator()
 (
     pointConstraint& x,
     const pointConstraint& y
@@ -189,7 +221,7 @@ void Foam::combineConstraintsEqOp::operator()
 }
 
 
-Foam::pointConstraint Foam::transform
+inline Foam::pointConstraint Foam::transform
 (
     const tensor& tt,
     const pointConstraint& v