From f4db00045d7ea9412eabcd7e8cfd8ad0460061a8 Mon Sep 17 00:00:00 2001
From: Vuko Vukcevic <vuko.vukcevic@fsb.hr>
Date: Thu, 2 Mar 2017 09:58:12 +0100
Subject: [PATCH 1/2] Updates to dynamicRefineFvMesh

Ported protection for non-hex cells from Vanilla
---
 .../dynamicRefineFvMesh/dynamicRefineFvMesh.C | 146 +++++++++++++++++-
 .../dynamicRefineFvMesh/dynamicRefineFvMesh.H |  20 +++
 .../directTopoChange/directActions/hexRef8.C  |  10 ++
 3 files changed, 174 insertions(+), 2 deletions(-)

diff --git a/src/dynamicMesh/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.C b/src/dynamicMesh/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.C
index e7727a6f3..baac97b89 100644
--- a/src/dynamicMesh/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.C
+++ b/src/dynamicMesh/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.C
@@ -31,6 +31,7 @@ License
 #include "syncTools.H"
 #include "pointFields.H"
 #include "directTopoChange.H"
+#include "cellSet.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -59,7 +60,16 @@ label dynamicRefineFvMesh::count
         {
             n++;
         }
+
+        // debug also serves to get-around Clang compiler trying to optimise
+        // out this forAll loop under O3 optimisation
+
+        if (debug)
+        {
+            Info<< "n=" << n << endl;
+        }
     }
+
     return n;
 }
 
@@ -860,11 +870,78 @@ void dynamicRefineFvMesh::extendMarkedCells(PackedBoolList& markedCell) const
 }
 
 
+void Foam::dynamicRefineFvMesh::checkEightAnchorPoints
+(
+    PackedBoolList& protectedCell,
+    label& nProtected
+) const
+{
+    const labelList& cellLevel = meshCutter_.cellLevel();
+    const labelList& pointLevel = meshCutter_.pointLevel();
+
+    labelList nAnchorPoints(nCells(), 0);
+
+    forAll(pointLevel, pointI)
+    {
+        const labelList& pCells = pointCells(pointI);
+
+        forAll(pCells, pCellI)
+        {
+            label cellI = pCells[pCellI];
+
+            if (pointLevel[pointI] <= cellLevel[cellI])
+            {
+                // Check if cell has already 8 anchor points -> protect cell
+                if (nAnchorPoints[cellI] == 8)
+                {
+                    if (protectedCell.set(cellI, true))
+                    {
+                        nProtected++;
+                    }
+                }
+
+                if (!protectedCell[cellI])
+                {
+                    nAnchorPoints[cellI]++;
+                }
+            }
+        }
+    }
+
+
+    forAll(protectedCell, cellI)
+    {
+        if (!protectedCell[cellI] && nAnchorPoints[cellI] != 8)
+        {
+            protectedCell.set(cellI, true);
+            nProtected++;
+        }
+    }
+}
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 dynamicRefineFvMesh::dynamicRefineFvMesh(const IOobject& io)
 :
     dynamicFvMesh(io),
+    singleMotionUpdate_
+    (
+        IOdictionary
+        (
+            IOobject
+            (
+                "dynamicMeshDict",
+                time().constant(),
+                *this,
+                IOobject::MUST_READ,
+                IOobject::NO_WRITE,
+                false
+            )
+        ).subDict(typeName + "Coeffs")
+        .lookupOrDefault<Switch>("singleMotionUpdate", true)
+    ),
+    curTimeIndex_(-1),
     meshCutter_(*this),
     dumpLevel_(false),
     nRefinementIterations_(0),
@@ -984,12 +1061,63 @@ dynamicRefineFvMesh::dynamicRefineFvMesh(const IOobject& io)
                 nProtected++;
             }
         }
+
+        // Also protect any cells that are less than hex
+        forAll(cells(), cellI)
+        {
+            const cell& cFaces = cells()[cellI];
+
+            if (cFaces.size() < 6)
+            {
+                if (protectedCell_.set(cellI, 1))
+                {
+                    nProtected++;
+                }
+            }
+            else
+            {
+                forAll(cFaces, cFaceI)
+                {
+                    if (faces()[cFaces[cFaceI]].size() < 4)
+                    {
+                        if (protectedCell_.set(cellI, 1))
+                        {
+                            nProtected++;
+                        }
+                        break;
+                    }
+                }
+            }
+        }
+
+        // Check cells for 8 corner points
+        checkEightAnchorPoints(protectedCell_, nProtected);
     }
 
     if (returnReduce(nProtected, sumOp<label>()) == 0)
     {
         protectedCell_.clear();
     }
+    else
+    {
+
+        cellSet protectedCells(*this, "protectedCells", nProtected);
+        forAll(protectedCell_, cellI)
+        {
+            if (protectedCell_[cellI])
+            {
+                protectedCells.insert(cellI);
+            }
+        }
+
+        Info<< "Detected " << returnReduce(nProtected, sumOp<label>())
+            << " cells that are protected from refinement."
+            << " Writing these to cellSet "
+            << protectedCells.name()
+            << "." << endl;
+
+        protectedCells.write();
+    }
 }
 
 
@@ -1003,6 +1131,20 @@ dynamicRefineFvMesh::~dynamicRefineFvMesh()
 
 bool dynamicRefineFvMesh::update()
 {
+    // Handling multiple calls in a single time step
+    if
+    (
+        singleMotionUpdate_
+     && curTimeIndex_ == this->time().timeIndex()
+    )
+    {
+        // This is not the first call to update, simply return false
+        return false;
+    }
+
+    // Update local time index
+    curTimeIndex_ = this->time().timeIndex();
+
     // Re-read dictionary. Choosen since usually -small so trivial amount
     // of time compared to actual refinement. Also very useful to be able
     // to modify on-the-fly.
@@ -1071,9 +1213,9 @@ bool dynamicRefineFvMesh::update()
                 << exit(FatalError);
         }
 
-        word field(refineDict.lookup("field"));
+        const word fieldName(refineDict.lookup("field"));
 
-        const volScalarField& vFld = lookupObject<volScalarField>(field);
+        const volScalarField& vFld = lookupObject<volScalarField>(fieldName);
 
         const scalar lowerRefineLevel =
             readScalar(refineDict.lookup("lowerRefineLevel"));
diff --git a/src/dynamicMesh/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.H b/src/dynamicMesh/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.H
index 094e9f893..15c57ec5c 100644
--- a/src/dynamicMesh/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.H
+++ b/src/dynamicMesh/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.H
@@ -55,6 +55,20 @@ class dynamicRefineFvMesh
 :
     public dynamicFvMesh
 {
+
+    // Private data
+
+        // Helper variables to enable switching between a single and multiple
+        // mesh motion updates within a time step (if update() is called more
+        // than once in a single time step)
+
+            //- Switch for single motion update (true by default)
+            Switch singleMotionUpdate_;
+
+            //- Helper varaible: current time index
+            label curTimeIndex_;
+
+
 protected:
 
         //- Mesh cutting engine
@@ -148,6 +162,12 @@ protected:
             //- Extend markedCell with cell-face-cell.
             void extendMarkedCells(PackedBoolList& markedCell) const;
 
+            //- Check all cells have 8 anchor points
+            void checkEightAnchorPoints
+            (
+                PackedBoolList& protectedCell,
+                label& nProtected
+            ) const;
 
 private:
 
diff --git a/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/hexRef8.C b/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/hexRef8.C
index 697fa47e6..e72740a1c 100644
--- a/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/hexRef8.C
+++ b/src/dynamicMesh/dynamicMesh/directTopoChange/directTopoChange/directActions/hexRef8.C
@@ -353,6 +353,16 @@ void Foam::hexRef8::modFace
 // Bit complex way to determine the unrefined edge length.
 Foam::scalar Foam::hexRef8::getLevel0EdgeLength() const
 {
+    if (cellLevel_.size() != mesh_.nCells())
+    {
+        FatalErrorIn("scalar hexRef8::getLevel0EdgeLength() const")
+            << "Number of cells in mesh:" << mesh_.nCells()
+            << " does not equal size of cellLevel:" << cellLevel_.size()
+            << endl
+            << "This might be because of a restart with inconsistent cellLevel."
+            << abort(FatalError);
+    }
+
     // Determine minimum edge length per refinement level
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 

From 57af002275b27e2b1dceb93abfc8e180941781f4 Mon Sep 17 00:00:00 2001
From: Vuko Vukcevic <v.vukcevic@wikki.co.uk>
Date: Fri, 28 Apr 2017 07:13:40 +0200
Subject: [PATCH 2/2] No slip wall boundary condition

---
 src/finiteVolume/Make/files                   |   1 +
 .../noSlipWall/noSlipWallFvPatchVectorField.C | 145 ++++++++++++++++++
 .../noSlipWall/noSlipWallFvPatchVectorField.H | 143 +++++++++++++++++
 3 files changed, 289 insertions(+)
 create mode 100644 src/finiteVolume/fields/fvPatchFields/derived/noSlipWall/noSlipWallFvPatchVectorField.C
 create mode 100644 src/finiteVolume/fields/fvPatchFields/derived/noSlipWall/noSlipWallFvPatchVectorField.H

diff --git a/src/finiteVolume/Make/files b/src/finiteVolume/Make/files
index 6503041fd..a5cd573b8 100644
--- a/src/finiteVolume/Make/files
+++ b/src/finiteVolume/Make/files
@@ -185,6 +185,7 @@ $(derivedFvPatchFields)/waveTransmissive/waveTransmissiveFvPatchFields.C
 $(derivedFvPatchFields)/uniformDensityHydrostaticPressure/uniformDensityHydrostaticPressureFvPatchScalarField.C
 $(derivedFvPatchFields)/pulseFixedValue/pulseFixedValueFvPatchFields.C
 $(derivedFvPatchFields)/waveTransmissiveInlet/waveTransmissiveInletFvPatchFields.C
+$(derivedFvPatchFields)/noSlipWall/noSlipWallFvPatchVectorField.C
 
 fvsPatchFields = fields/fvsPatchFields
 $(fvsPatchFields)/fvsPatchField/fvsPatchFields.C
diff --git a/src/finiteVolume/fields/fvPatchFields/derived/noSlipWall/noSlipWallFvPatchVectorField.C b/src/finiteVolume/fields/fvPatchFields/derived/noSlipWall/noSlipWallFvPatchVectorField.C
new file mode 100644
index 000000000..1bf1f0ea1
--- /dev/null
+++ b/src/finiteVolume/fields/fvPatchFields/derived/noSlipWall/noSlipWallFvPatchVectorField.C
@@ -0,0 +1,145 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | foam-extend: Open Source CFD
+   \\    /   O peration     | Version:     4.0
+    \\  /    A nd           | Web:         http://www.foam-extend.org
+     \\/     M anipulation  | For copyright notice see file Copyright
+-------------------------------------------------------------------------------
+License
+    This file is part of foam-extend.
+
+    foam-extend is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation, either version 3 of the License, or (at your
+    option) any later version.
+
+    foam-extend is distributed in the hope that it will be useful, but
+    WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+    General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "noSlipWallFvPatchVectorField.H"
+#include "addToRunTimeSelectionTable.H"
+#include "fvPatchFieldMapper.H"
+#include "volFields.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+noSlipWallFvPatchVectorField::noSlipWallFvPatchVectorField
+(
+    const fvPatch& p,
+    const DimensionedField<vector, volMesh>& iF
+)
+:
+    fixedValueFvPatchVectorField(p, iF)
+{}
+
+
+noSlipWallFvPatchVectorField::noSlipWallFvPatchVectorField
+(
+    const fvPatch& p,
+    const DimensionedField<vector, volMesh>& iF,
+    const dictionary& dict
+)
+:
+    fixedValueFvPatchVectorField(p, iF, dict)
+{}
+
+
+noSlipWallFvPatchVectorField::noSlipWallFvPatchVectorField
+(
+    const noSlipWallFvPatchVectorField& ptf,
+    const fvPatch& p,
+    const DimensionedField<vector, volMesh>& iF,
+    const fvPatchFieldMapper& mapper
+)
+:
+    fixedValueFvPatchVectorField(ptf, p, iF, mapper)
+{}
+
+
+noSlipWallFvPatchVectorField::noSlipWallFvPatchVectorField
+(
+    const noSlipWallFvPatchVectorField& ptf
+)
+:
+    fixedValueFvPatchVectorField(ptf)
+{}
+
+
+noSlipWallFvPatchVectorField::noSlipWallFvPatchVectorField
+(
+    const noSlipWallFvPatchVectorField& ptf,
+    const DimensionedField<vector, volMesh>& iF
+)
+:
+    fixedValueFvPatchVectorField(ptf, iF)
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+tmp<vectorField> noSlipWallFvPatchVectorField::gradientInternalCoeffs() const
+{
+    const vectorField nPatch = this->patch().nf();
+
+    // Calculate the diagonal part of the transformation tensor
+    vectorField implicitCoeffs(nPatch.size(), vector::zero);
+    contractLinear(implicitCoeffs, I - nPatch*nPatch);
+
+    return -this->patch().deltaCoeffs()*implicitCoeffs;
+}
+
+
+tmp<vectorField> noSlipWallFvPatchVectorField::gradientBoundaryCoeffs() const
+{
+    // Get necessary field
+    const vectorField& UPatch = *this;
+    const vectorField nPatch = this->patch().nf();
+    const tensorField T = I - nPatch*nPatch;
+    const vectorField UPatchInternal = this->patchInternalField();
+
+    // Calculate the explicit contribution in the loop for performance
+    vectorField explicitCoeffs(nPatch.size(), vector::zero);
+
+    forAll(explicitCoeffs, faceI)
+    {
+        const tensor& curT = T[faceI];
+        tensor curDiagT(tensor::zero);
+        expandLinear(curDiagT, contractLinear(curT));
+
+        explicitCoeffs[faceI] =
+            (curT & UPatch[faceI])
+          - (
+                (curT - curDiagT)
+              & UPatchInternal[faceI]
+            );
+    }
+
+    return this->patch().deltaCoeffs()*explicitCoeffs;
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+makePatchTypeField
+(
+    fvPatchVectorField,
+    noSlipWallFvPatchVectorField
+);
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fields/fvPatchFields/derived/noSlipWall/noSlipWallFvPatchVectorField.H b/src/finiteVolume/fields/fvPatchFields/derived/noSlipWall/noSlipWallFvPatchVectorField.H
new file mode 100644
index 000000000..acbd20e86
--- /dev/null
+++ b/src/finiteVolume/fields/fvPatchFields/derived/noSlipWall/noSlipWallFvPatchVectorField.H
@@ -0,0 +1,143 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | foam-extend: Open Source CFD
+   \\    /   O peration     | Version:     4.0
+    \\  /    A nd           | Web:         http://www.foam-extend.org
+     \\/     M anipulation  | For copyright notice see file Copyright
+-------------------------------------------------------------------------------
+License
+    This file is part of foam-extend.
+
+    foam-extend is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation, either version 3 of the License, or (at your
+    option) any later version.
+
+    foam-extend is distributed in the hope that it will be useful, but
+    WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+    General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::noSlipWallFvPatchVectorField
+
+Description
+    Foam::noSlipWallFvPatchVectorField
+
+SourceFiles
+    noSlipWallFvPatchVectorField.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef noSlipWallFvPatchVectorField_H
+#define noSlipWallFvPatchVectorField_H
+
+#include "fvPatchFields.H"
+#include "fixedValueFvPatchFields.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                  Class noSlipWallFvPatchVectorField Declaration
+\*---------------------------------------------------------------------------*/
+
+class noSlipWallFvPatchVectorField
+:
+    public fixedValueFvPatchVectorField
+{
+
+public:
+
+    //- Runtime type information
+    TypeName("noSlipWall");
+
+
+    // Constructors
+
+        //- Construct from patch and internal field
+        noSlipWallFvPatchVectorField
+        (
+            const fvPatch&,
+            const DimensionedField<vector, volMesh>&
+        );
+
+        //- Construct from patch, internal field and dictionary
+        noSlipWallFvPatchVectorField
+        (
+            const fvPatch&,
+            const DimensionedField<vector, volMesh>&,
+            const dictionary&
+        );
+
+        //- Construct by mapping the given noSlipWallFvPatchVectorField
+        //  onto a new patch
+        noSlipWallFvPatchVectorField
+        (
+            const noSlipWallFvPatchVectorField&,
+            const fvPatch&,
+            const DimensionedField<vector, volMesh>&,
+            const fvPatchFieldMapper&
+        );
+
+        //- Construct as copy
+        noSlipWallFvPatchVectorField
+        (
+            const noSlipWallFvPatchVectorField&
+        );
+
+        //- Construct and return a clone
+        virtual tmp<fvPatchVectorField> clone() const
+        {
+            return tmp<fvPatchVectorField>
+            (
+                new noSlipWallFvPatchVectorField(*this)
+            );
+        }
+
+        //- Construct as copy setting internal field reference
+        noSlipWallFvPatchVectorField
+        (
+            const noSlipWallFvPatchVectorField&,
+            const DimensionedField<vector, volMesh>&
+        );
+
+        //- Construct and return a clone setting internal field reference
+        virtual tmp<fvPatchVectorField> clone
+        (
+            const DimensionedField<vector, volMesh>& iF
+        ) const
+        {
+            return tmp<fvPatchVectorField>
+            (
+                new noSlipWallFvPatchVectorField(*this, iF)
+            );
+        }
+
+
+    // Member functions
+
+            //- Return the matrix diagonal coefficients corresponding to the
+            //  evaluation of the gradient of this patchField
+            virtual tmp<vectorField> gradientInternalCoeffs() const;
+
+            //- Return the matrix source coefficients corresponding to the
+            //  evaluation of the gradient of this patchField
+            virtual tmp<vectorField> gradientBoundaryCoeffs() const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //