diff --git a/src/OpenFOAM/meshes/MeshObject/MeshObject.C b/src/OpenFOAM/meshes/MeshObject/MeshObject.C
index c089e0e531ead8fc903bd0641b73a97b1297d7af..fddfadfa59e0c50e6037488a8c8398db5fc2d26e 100644
--- a/src/OpenFOAM/meshes/MeshObject/MeshObject.C
+++ b/src/OpenFOAM/meshes/MeshObject/MeshObject.C
@@ -83,6 +83,26 @@ const Type& Foam::MeshObject<Mesh, Type>::New
}
+template<class Mesh, class Type>
+template<class Data1, class Data2>
+const Type& Foam::MeshObject<Mesh, Type>::New
+(
+ const Mesh& mesh,
+ const Data1& d1,
+ const Data2& d2
+)
+{
+ if (!mesh.db().objectRegistry::foundObject<Type>(Type::typeName))
+ {
+ return store(new Type(mesh, d1, d2));
+ }
+ else
+ {
+ return mesh.db().objectRegistry::lookupObject<Type>(Type::typeName);
+ }
+}
+
+
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * //
template<class Mesh, class Type>
diff --git a/src/OpenFOAM/meshes/MeshObject/MeshObject.H b/src/OpenFOAM/meshes/MeshObject/MeshObject.H
index 595f2f1df4edf577495bb277686e2632c2c0c46c..b142b47f968a7d165fec3adf0025fede3687e556 100644
--- a/src/OpenFOAM/meshes/MeshObject/MeshObject.H
+++ b/src/OpenFOAM/meshes/MeshObject/MeshObject.H
@@ -72,6 +72,9 @@ public:
template<class Data>
static const Type& New(const Mesh& mesh, const Data& d);
+ template<class Data1, class Data2>
+ static const Type& New(const Mesh& mesh, const Data1&, const Data2&);
+
// Destructor
diff --git a/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistribute.C b/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistribute.C
index 5421b328d8e5d585bd0d0d55d3ef13977a86ba3c..5ef7b693a560a736933fbcd70d0e213890184377 100644
--- a/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistribute.C
+++ b/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistribute.C
@@ -271,4 +271,144 @@ Foam::mapDistribute::mapDistribute
}
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+void Foam::mapDistribute::compact(const boolList& elemIsUsed)
+{
+ // 1. send back to sender. Have him delete the corresponding element
+ // from the submap and do the same to the constructMap locally
+ // (and in same order).
+
+ // Send elemIsUsed field to neighbour. Use nonblocking code from
+ // mapDistribute but in reverse order.
+ {
+ List<boolList> sendFields(Pstream::nProcs());
+
+ for (label domain = 0; domain < Pstream::nProcs(); domain++)
+ {
+ const labelList& map = constructMap_[domain];
+
+ if (domain != Pstream::myProcNo() && map.size() > 0)
+ {
+ boolList& subField = sendFields[domain];
+ subField.setSize(map.size());
+ forAll(map, i)
+ {
+ subField[i] = elemIsUsed[map[i]];
+ }
+
+ OPstream::write
+ (
+ Pstream::nonBlocking,
+ domain,
+ reinterpret_cast<const char*>(subField.begin()),
+ subField.size()*sizeof(bool)
+ );
+ }
+ }
+
+ // Set up receives from neighbours
+
+ List<boolList> recvFields(Pstream::nProcs());
+
+ for (label domain = 0; domain < Pstream::nProcs(); domain++)
+ {
+ const labelList& map = subMap_[domain];
+
+ if (domain != Pstream::myProcNo() && map.size() > 0)
+ {
+ recvFields[domain].setSize(map.size());
+ IPstream::read
+ (
+ Pstream::nonBlocking,
+ domain,
+ reinterpret_cast<char*>(recvFields[domain].begin()),
+ recvFields[domain].size()*sizeof(bool)
+ );
+ }
+ }
+
+
+ // Set up 'send' to myself - write directly into recvFields
+
+ {
+ const labelList& map = constructMap_[Pstream::myProcNo()];
+
+ recvFields[Pstream::myProcNo()].setSize(map.size());
+ forAll(map, i)
+ {
+ recvFields[Pstream::myProcNo()][i] = elemIsUsed[map[i]];
+ }
+ }
+
+
+ // Wait for all to finish
+
+ OPstream::waitRequests();
+ IPstream::waitRequests();
+
+
+ // Compact out all submap entries that are referring to unused elements
+ for (label domain = 0; domain < Pstream::nProcs(); domain++)
+ {
+ const labelList& map = subMap_[domain];
+
+ labelList newMap(map.size());
+ label newI = 0;
+
+ forAll(map, i)
+ {
+ if (recvFields[domain][i])
+ {
+ // So element is used on destination side
+ newMap[newI++] = map[i];
+ }
+ }
+ if (newI < map.size())
+ {
+ newMap.setSize(newI);
+ subMap_[domain].transfer(newMap);
+ }
+ }
+ }
+
+
+ // 2. remove from construct map - since end-result (element in elemIsUsed)
+ // not used.
+
+ label maxConstructIndex = -1;
+
+ for (label domain = 0; domain < Pstream::nProcs(); domain++)
+ {
+ const labelList& map = constructMap_[domain];
+
+ labelList newMap(map.size());
+ label newI = 0;
+
+ forAll(map, i)
+ {
+ label destinationI = map[i];
+
+ // Is element is used on destination side
+ if (elemIsUsed[destinationI])
+ {
+ maxConstructIndex = max(maxConstructIndex, destinationI);
+
+ newMap[newI++] = destinationI;
+ }
+ }
+ if (newI < map.size())
+ {
+ newMap.setSize(newI);
+ constructMap_[domain].transfer(newMap);
+ }
+ }
+
+ constructSize_ = maxConstructIndex+1;
+
+ // Clear the schedule (note:not necessary if nothing changed)
+ schedulePtr_.clear();
+}
+
+
// ************************************************************************* //
diff --git a/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistribute.H b/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistribute.H
index 748834ca4d28ef43dcab1816d451f43d4bfdacb7..8427670a25c77705dfadd8b32da2c0e2acb976d3 100644
--- a/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistribute.H
+++ b/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistribute.H
@@ -51,6 +51,7 @@ SourceFiles
#include "labelList.H"
#include "labelPair.H"
#include "Pstream.H"
+#include "boolList.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -155,6 +156,12 @@ public:
// Other
+ //- Compact maps. Gets per field a bool whether it is used (locally)
+ // and works out itself what this side and sender side can remove
+ // from maps.
+ void compact(const boolList& elemIsUsed);
+
+
//- Distribute data. Note:schedule only used for Pstream::scheduled
// for now, all others just use send-to-all, receive-from-all.
template<class T>
@@ -168,6 +175,21 @@ public:
List<T>&
);
+ //- Distribute data. If multiple processors writing to same
+ // position adds contributions using cop.
+ template<class T, class CombineOp>
+ static void distribute
+ (
+ const Pstream::commsTypes commsType,
+ const List<labelPair>& schedule,
+ const label constructSize,
+ const labelListList& subMap,
+ const labelListList& constructMap,
+ List<T>&,
+ const CombineOp& cop,
+ const T& nullValue
+ );
+
//- Distribute data using default commsType.
template<class T>
void distribute(List<T>& fld) const
diff --git a/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistributeTemplates.C b/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistributeTemplates.C
index 5f16ad1335fb2eb58c07f2c38b5f05214dfaee4d..8b4598bfb8ad0ffec1661a184ddc21290e20e842 100644
--- a/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistributeTemplates.C
+++ b/src/OpenFOAM/meshes/polyMesh/mapPolyMesh/mapDistribute/mapDistributeTemplates.C
@@ -268,8 +268,122 @@ void Foam::mapDistribute::distribute
}
+ // Set up 'send' to myself
+
+ {
+ const labelList& map = subMap[Pstream::myProcNo()];
+
+ List<T>& subField = sendFields[Pstream::myProcNo()];
+ subField.setSize(map.size());
+ forAll(map, i)
+ {
+ subField[i] = field[map[i]];
+ }
+ }
+
+
// Combine bits. Note that can reuse field storage
+ field.setSize(constructSize);
+
+
+ // Receive sub field from myself (sendFields[Pstream::myProcNo()])
+ {
+ const labelList& map = constructMap[Pstream::myProcNo()];
+ const List<T>& subField = sendFields[Pstream::myProcNo()];
+
+ forAll(map, i)
+ {
+ field[map[i]] = subField[i];
+ }
+ }
+
+
+ // Wait for all to finish
+
+ OPstream::waitRequests();
+ IPstream::waitRequests();
+
+ // Collect neighbour fields
+
+ for (label domain = 0; domain < Pstream::nProcs(); domain++)
+ {
+ const labelList& map = constructMap[domain];
+
+ if (domain != Pstream::myProcNo() && map.size() > 0)
+ {
+ if (recvFields[domain].size() != map.size())
+ {
+ FatalErrorIn
+ (
+ "template<class T>\n"
+ "void mapDistribute::distribute\n"
+ "(\n"
+ " const Pstream::commsTypes commsType,\n"
+ " const List<labelPair>& schedule,\n"
+ " const label constructSize,\n"
+ " const labelListList& subMap,\n"
+ " const labelListList& constructMap,\n"
+ " List<T>& field\n"
+ ")\n"
+ ) << "Expected from processor " << domain
+ << " " << map.size() << " but received "
+ << recvFields[domain].size() << " elements."
+ << abort(FatalError);
+ }
+
+ forAll(map, i)
+ {
+ field[map[i]] = recvFields[domain][i];
+ }
+ }
+ }
+ }
+ else
+ {
+ FatalErrorIn("mapDistribute::distribute(..)")
+ << "Unknown communication schedule " << commsType
+ << abort(FatalError);
+ }
+}
+
+
+// Distribute list.
+template<class T, class CombineOp>
+void Foam::mapDistribute::distribute
+(
+ const Pstream::commsTypes commsType,
+ const List<labelPair>& schedule,
+ const label constructSize,
+ const labelListList& subMap,
+ const labelListList& constructMap,
+ List<T>& field,
+ const CombineOp& cop,
+ const T& nullValue
+)
+{
+ if (commsType == Pstream::blocking)
+ {
+ // Since buffered sending can reuse the field to collect the
+ // received data.
+
+ // Send sub field to neighbour
+ for (label domain = 0; domain < Pstream::nProcs(); domain++)
+ {
+ const labelList& map = subMap[domain];
+
+ if (domain != Pstream::myProcNo() && map.size() > 0)
+ {
+ List<T> subField(map.size());
+ forAll(map, i)
+ {
+ subField[i] = field[map[i]];
+ }
+ OPstream toNbr(Pstream::blocking, domain);
+ toNbr << subField;
+ }
+ }
+
// Subset myself
const labelList& mySubMap = subMap[Pstream::myProcNo()];
@@ -279,14 +393,231 @@ void Foam::mapDistribute::distribute
subField[i] = field[mySubMap[i]];
}
+ // Receive sub field from myself (subField)
+ const labelList& map = constructMap[Pstream::myProcNo()];
+
field.setSize(constructSize);
+ field = nullValue;
+
+ forAll(map, i)
+ {
+ cop(field[map[i]], subField[i]);
+ }
+
+ // Receive sub field from neighbour
+ for (label domain = 0; domain < Pstream::nProcs(); domain++)
+ {
+ const labelList& map = constructMap[domain];
+
+ if (domain != Pstream::myProcNo() && map.size() > 0)
+ {
+ IPstream fromNbr(Pstream::blocking, domain);
+ List<T> subField(fromNbr);
+
+ if (subField.size() != map.size())
+ {
+ FatalErrorIn
+ (
+ "template<class T>\n"
+ "void mapDistribute::distribute\n"
+ "(\n"
+ " const Pstream::commsTypes commsType,\n"
+ " const List<labelPair>& schedule,\n"
+ " const label constructSize,\n"
+ " const labelListList& subMap,\n"
+ " const labelListList& constructMap,\n"
+ " List<T>& field\n"
+ ")\n"
+ ) << "Expected from processor " << domain
+ << " " << map.size() << " but received "
+ << subField.size() << " elements."
+ << abort(FatalError);
+ }
+
+ forAll(map, i)
+ {
+ cop(field[map[i]], subField[i]);
+ }
+ }
+ }
+ }
+ else if (commsType == Pstream::scheduled)
+ {
+ // Need to make sure I don't overwrite field with received data
+ // since the data might need to be sent to another processor. So
+ // allocate a new field for the results.
+ List<T> newField(constructSize, nullValue);
+
+ // Subset myself
+ const labelList& mySubMap = subMap[Pstream::myProcNo()];
+
+ List<T> subField(mySubMap.size());
+ forAll(mySubMap, i)
+ {
+ subField[i] = field[mySubMap[i]];
+ }
// Receive sub field from myself (subField)
const labelList& map = constructMap[Pstream::myProcNo()];
forAll(map, i)
{
- field[map[i]] = subField[i];
+ cop(newField[map[i]], subField[i]);
+ }
+
+ // Schedule will already have pruned 0-sized comms
+ forAll(schedule, i)
+ {
+ const labelPair& twoProcs = schedule[i];
+ label sendProc = twoProcs[0];
+ label recvProc = twoProcs[1];
+
+ if (Pstream::myProcNo() == sendProc)
+ {
+ // I am sender. Send to recvProc.
+ const labelList& map = subMap[recvProc];
+
+ List<T> subField(map.size());
+ forAll(map, i)
+ {
+ subField[i] = field[map[i]];
+ }
+
+ OPstream toNbr(Pstream::scheduled, recvProc);
+ toNbr << subField;
+ }
+ else
+ {
+ // I am receiver. Receive from sendProc.
+ IPstream fromNbr(Pstream::scheduled, sendProc);
+ List<T> subField(fromNbr);
+
+ const labelList& map = constructMap[sendProc];
+
+ if (subField.size() != map.size())
+ {
+ FatalErrorIn
+ (
+ "template<class T>\n"
+ "void mapDistribute::distribute\n"
+ "(\n"
+ " const Pstream::commsTypes commsType,\n"
+ " const List<labelPair>& schedule,\n"
+ " const label constructSize,\n"
+ " const labelListList& subMap,\n"
+ " const labelListList& constructMap,\n"
+ " List<T>& field\n"
+ ")\n"
+ ) << "Expected from processor " << sendProc
+ << " " << map.size() << " but received "
+ << subField.size() << " elements."
+ << abort(FatalError);
+ }
+
+ forAll(map, i)
+ {
+ cop(newField[map[i]], subField[i]);
+ }
+ }
+ }
+ field.transfer(newField);
+ }
+ else if (commsType == Pstream::nonBlocking)
+ {
+ if (!contiguous<T>())
+ {
+ FatalErrorIn
+ (
+ "template<class T>\n"
+ "void mapDistribute::distribute\n"
+ "(\n"
+ " const Pstream::commsTypes commsType,\n"
+ " const List<labelPair>& schedule,\n"
+ " const label constructSize,\n"
+ " const labelListList& subMap,\n"
+ " const labelListList& constructMap,\n"
+ " List<T>& field\n"
+ ")\n"
+ ) << "Non-blocking only supported for contiguous data."
+ << exit(FatalError);
+ }
+
+ // Set up sends to neighbours
+
+ List<List<T > > sendFields(Pstream::nProcs());
+
+ for (label domain = 0; domain < Pstream::nProcs(); domain++)
+ {
+ const labelList& map = subMap[domain];
+
+ if (domain != Pstream::myProcNo() && map.size() > 0)
+ {
+ List<T>& subField = sendFields[domain];
+ subField.setSize(map.size());
+ forAll(map, i)
+ {
+ subField[i] = field[map[i]];
+ }
+
+ OPstream::write
+ (
+ Pstream::nonBlocking,
+ domain,
+ reinterpret_cast<const char*>(subField.begin()),
+ subField.size()*sizeof(T)
+ );
+ }
+ }
+
+ // Set up receives from neighbours
+
+ List<List<T > > recvFields(Pstream::nProcs());
+
+ for (label domain = 0; domain < Pstream::nProcs(); domain++)
+ {
+ const labelList& map = constructMap[domain];
+
+ if (domain != Pstream::myProcNo() && map.size() > 0)
+ {
+ recvFields[domain].setSize(map.size());
+ IPstream::read
+ (
+ Pstream::nonBlocking,
+ domain,
+ reinterpret_cast<char*>(recvFields[domain].begin()),
+ recvFields[domain].size()*sizeof(T)
+ );
+ }
+ }
+
+ // Set up 'send' to myself
+
+ {
+ const labelList& map = subMap[Pstream::myProcNo()];
+
+ List<T>& subField = sendFields[Pstream::myProcNo()];
+ subField.setSize(map.size());
+ forAll(map, i)
+ {
+ subField[i] = field[map[i]];
+ }
+ }
+
+
+ // Combine bits. Note that can reuse field storage
+
+ field.setSize(constructSize);
+ field = nullValue;
+
+ // Receive sub field from myself (subField)
+ {
+ const labelList& map = constructMap[Pstream::myProcNo()];
+ const List<T>& subField = sendFields[Pstream::myProcNo()];
+
+ forAll(map, i)
+ {
+ cop(field[map[i]], subField[i]);
+ }
}
@@ -325,7 +656,7 @@ void Foam::mapDistribute::distribute
forAll(map, i)
{
- field[map[i]] = recvFields[domain][i];
+ cop(field[map[i]], recvFields[domain][i]);
}
}
}
diff --git a/src/finiteVolume/Make/files b/src/finiteVolume/Make/files
index 923e4c50e60b14f7b95c7daad7b1ab04a90e98aa..7bf27115a57beb33d3ce99c1b17b5204e9b60038 100644
--- a/src/finiteVolume/Make/files
+++ b/src/finiteVolume/Make/files
@@ -40,6 +40,24 @@ $(fvMeshMapper)/fvSurfaceMapper.C
extendedStencil = fvMesh/extendedStencil
$(extendedStencil)/extendedStencil.C
+$(extendedStencil)/extendedUpwindStencil.C
+$(extendedStencil)/extendedCentredStencil.C
+
+$(extendedStencil)/faceStencil/faceStencil.C
+$(extendedStencil)/faceStencil/faceEdgeCellStencil.C
+$(extendedStencil)/faceStencil/cellFaceCellStencil.C
+$(extendedStencil)/faceStencil/cellPointCellStencil.C
+$(extendedStencil)/faceStencil/cellEdgeCellStencil.C
+
+$(extendedStencil)/extendedStencilMeshObjects/centredCECStencilObject.C
+$(extendedStencil)/extendedStencilMeshObjects/centredCFCStencilObject.C
+$(extendedStencil)/extendedStencilMeshObjects/centredCPCStencilObject.C
+$(extendedStencil)/extendedStencilMeshObjects/centredFECStencilObject.C
+
+$(extendedStencil)/extendedStencilMeshObjects/upwindCECStencilObject.C
+$(extendedStencil)/extendedStencilMeshObjects/upwindCFCStencilObject.C
+$(extendedStencil)/extendedStencilMeshObjects/upwindCPCStencilObject.C
+$(extendedStencil)/extendedStencilMeshObjects/upwindFECStencilObject.C
fvPatchFields = fields/fvPatchFields
$(fvPatchFields)/fvPatchField/fvPatchFields.C
@@ -168,6 +186,8 @@ $(schemes)/harmonic/harmonic.C
$(schemes)/localBlended/localBlended.C
$(schemes)/localMax/localMax.C
$(schemes)/localMin/localMin.C
+//$(schemes)/linearFit/linearFit.C
+//$(schemes)/linearFit/linearFitData.C
$(schemes)/quadraticFit/quadraticFit.C
$(schemes)/quadraticFit/quadraticFitData.C
@@ -251,8 +271,8 @@ $(snGradSchemes)/snGradScheme/snGradSchemes.C
$(snGradSchemes)/correctedSnGrad/correctedSnGrads.C
$(snGradSchemes)/limitedSnGrad/limitedSnGrads.C
$(snGradSchemes)/uncorrectedSnGrad/uncorrectedSnGrads.C
-$(snGradSchemes)/quadraticFitSnGrad/quadraticFitSnGradData.C
-$(snGradSchemes)/quadraticFitSnGrad/quadraticFitSnGrads.C
+//$(snGradSchemes)/quadraticFitSnGrad/quadraticFitSnGradData.C
+//$(snGradSchemes)/quadraticFitSnGrad/quadraticFitSnGrads.C
convectionSchemes = finiteVolume/convectionSchemes
$(convectionSchemes)/convectionScheme/convectionSchemes.C
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedCentredStencil.C b/src/finiteVolume/fvMesh/extendedStencil/extendedCentredStencil.C
new file mode 100644
index 0000000000000000000000000000000000000000..956e796e18245219970be90dee4f7a06dce13118
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedCentredStencil.C
@@ -0,0 +1,65 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "mapDistribute.H"
+#include "extendedCentredStencil.H"
+#include "faceStencil.H"
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::extendedCentredStencil::extendedCentredStencil(const faceStencil& stencil)
+:
+ extendedStencil(stencil.mesh())
+{
+ stencil_ = stencil;
+
+ // Calculate distribute map (also renumbers stencil)
+ mapPtr_ = calcDistributeMap(stencil.globalNumbering(), stencil_);
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+// Per face which elements of the stencil to keep.
+void Foam::extendedCentredStencil::compact()
+{
+ boolList isInStencil(map().constructSize(), false);
+
+ forAll(stencil_, faceI)
+ {
+ const labelList& stencilCells = stencil_[faceI];
+
+ forAll(stencilCells, i)
+ {
+ isInStencil[stencilCells[i]] = true;
+ }
+ }
+
+ mapPtr_().compact(isInStencil);
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedCentredStencil.H b/src/finiteVolume/fvMesh/extendedStencil/extendedCentredStencil.H
new file mode 100644
index 0000000000000000000000000000000000000000..0ae37d5eb6c8e83ea015117fe2d1a429d35badc5
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedCentredStencil.H
@@ -0,0 +1,137 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::extendedCentredStencil
+
+Description
+
+SourceFiles
+ extendedCentredStencil.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef extendedCentredStencil_H
+#define extendedCentredStencil_H
+
+#include "extendedStencil.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class faceStencil;
+
+/*---------------------------------------------------------------------------*\
+ Class extendedCentredStencil Declaration
+\*---------------------------------------------------------------------------*/
+
+class extendedCentredStencil
+:
+ public extendedStencil
+{
+ // Private data
+
+ //- Swap map for getting neigbouring data
+ autoPtr<mapDistribute> mapPtr_;
+
+ //- Per face the stencil.
+ labelListList stencil_;
+
+
+ // Private Member Functions
+
+ //- Disallow default bitwise copy construct
+ extendedCentredStencil(const extendedCentredStencil&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const extendedCentredStencil&);
+
+
+public:
+
+ // Constructors
+
+ //- Construct from uncompacted face stencil
+ explicit extendedCentredStencil(const faceStencil&);
+
+
+ // Member Functions
+
+ //- Return reference to the parallel distribution map
+ const mapDistribute& map() const
+ {
+ return mapPtr_();
+ }
+
+ //- Return reference to the stencil
+ const labelListList& stencil() const
+ {
+ return stencil_;
+ }
+
+ //- After removing elements from the stencil adapt the schedule (map).
+ void compact();
+
+ //- Use map to get the data into stencil order
+ template<class T>
+ void collectData
+ (
+ const GeometricField<T, fvPatchField, volMesh>& fld,
+ List<List<T> >& stencilFld
+ ) const
+ {
+ extendedStencil::collectData(map(), stencil(), fld, stencilFld);
+ }
+
+ //- Sum vol field contributions to create face values
+ template<class Type>
+ tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > weightedSum
+ (
+ const GeometricField<Type, fvPatchField, volMesh>& fld,
+ const List<List<scalar> >& stencilWeights
+ ) const
+ {
+ return extendedStencil::weightedSum
+ (
+ map(),
+ stencil(),
+ fld,
+ stencilWeights
+ );
+ }
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencil.C b/src/finiteVolume/fvMesh/extendedStencil/extendedStencil.C
index bb4c6c1534c66c5f32911fd57bd62cd2acd040e5..c2273d4a52f44ccbd104f5e526ac59729f87f3cf 100644
--- a/src/finiteVolume/fvMesh/extendedStencil/extendedStencil.C
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencil.C
@@ -32,154 +32,83 @@ License
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
// Calculates per face a list of global cell/face indices.
-void Foam::extendedStencil::calcFaceStencils
+void Foam::extendedStencil::calcFaceStencil
(
- const polyMesh& mesh,
- const globalIndex& globalNumbering
+ const labelListList& globalCellCells,
+ labelListList& faceStencil
)
{
- const polyBoundaryMesh& patches = mesh.boundaryMesh();
- const label nBnd = mesh.nFaces()-mesh.nInternalFaces();
- const labelList& own = mesh.faceOwner();
- const labelList& nei = mesh.faceNeighbour();
+ const polyBoundaryMesh& patches = mesh_.boundaryMesh();
+ const label nBnd = mesh_.nFaces()-mesh_.nInternalFaces();
+ const labelList& own = mesh_.faceOwner();
+ const labelList& nei = mesh_.faceNeighbour();
- // Determine neighbouring global cell or boundary face
- // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+ // Determine neighbouring global cell Cells
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
- labelList neiGlobal(nBnd);
+ labelListList neiGlobalCellCells(nBnd);
forAll(patches, patchI)
{
const polyPatch& pp = patches[patchI];
- label faceI = pp.start();
if (pp.coupled())
{
- // For coupled faces get the cell on the other side
- forAll(pp, i)
- {
- label bFaceI = faceI-mesh.nInternalFaces();
- neiGlobal[bFaceI] = globalNumbering.toGlobal(own[faceI]);
- faceI++;
- }
- }
- else if (isA<emptyPolyPatch>(pp))
- {
- forAll(pp, i)
- {
- label bFaceI = faceI-mesh.nInternalFaces();
- neiGlobal[bFaceI] = -1;
- faceI++;
- }
- }
- else
- {
- // For noncoupled faces get the boundary face.
+ label faceI = pp.start();
+
forAll(pp, i)
{
- label bFaceI = faceI-mesh.nInternalFaces();
- neiGlobal[bFaceI] =
- globalNumbering.toGlobal(mesh.nCells()+bFaceI);
+ neiGlobalCellCells[faceI-mesh_.nInternalFaces()] =
+ globalCellCells[own[faceI]];
faceI++;
}
}
}
- syncTools::swapBoundaryFaceList(mesh, neiGlobal, false);
-
-
- // Determine cellCells in global numbering
- // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
- labelListList globalCellCells(mesh.nCells());
- forAll(globalCellCells, cellI)
- {
- const cell& cFaces = mesh.cells()[cellI];
-
- labelList& cCells = globalCellCells[cellI];
-
- cCells.setSize(cFaces.size());
-
- // Collect neighbouring cells/faces
- label nNbr = 0;
- forAll(cFaces, i)
- {
- label faceI = cFaces[i];
-
- if (mesh.isInternalFace(faceI))
- {
- label nbrCellI = own[faceI];
- if (nbrCellI == cellI)
- {
- nbrCellI = nei[faceI];
- }
- cCells[nNbr++] = globalNumbering.toGlobal(nbrCellI);
- }
- else
- {
- label nbrCellI = neiGlobal[faceI-mesh.nInternalFaces()];
- if (nbrCellI != -1)
- {
- cCells[nNbr++] = nbrCellI;
- }
- }
- }
- cCells.setSize(nNbr);
- }
-
-
- // Determine neighbouring global cell Cells
- // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
- labelListList neiGlobalCellCells(nBnd);
- for (label faceI = mesh.nInternalFaces(); faceI < mesh.nFaces(); faceI++)
- {
- neiGlobalCellCells[faceI-mesh.nInternalFaces()] =
- globalCellCells[own[faceI]];
- }
- syncTools::swapBoundaryFaceList(mesh, neiGlobalCellCells, false);
+ syncTools::swapBoundaryFaceList(mesh_, neiGlobalCellCells, false);
// Construct stencil in global numbering
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
- stencil_.setSize(mesh.nFaces());
+ faceStencil.setSize(mesh_.nFaces());
- labelHashSet faceStencil;
+ labelHashSet faceStencilSet;
- for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
+ for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
{
- faceStencil.clear();
- label globalOwn = globalNumbering.toGlobal(own[faceI]);
- faceStencil.insert(globalOwn);
+ faceStencilSet.clear();
+
const labelList& ownCCells = globalCellCells[own[faceI]];
+ label globalOwn = ownCCells[0];
+ // Insert cellCells
forAll(ownCCells, i)
{
- faceStencil.insert(ownCCells[i]);
+ faceStencilSet.insert(ownCCells[i]);
}
- label globalNei = globalNumbering.toGlobal(nei[faceI]);
- faceStencil.insert(globalNei);
const labelList& neiCCells = globalCellCells[nei[faceI]];
+ label globalNei = neiCCells[0];
+ // Insert cellCells
forAll(neiCCells, i)
{
- faceStencil.insert(neiCCells[i]);
+ faceStencilSet.insert(neiCCells[i]);
}
// Guarantee owner first, neighbour second.
- stencil_[faceI].setSize(faceStencil.size());
+ faceStencil[faceI].setSize(faceStencilSet.size());
label n = 0;
- stencil_[faceI][n++] = globalOwn;
- stencil_[faceI][n++] = globalNei;
- forAllConstIter(labelHashSet, faceStencil, iter)
+ faceStencil[faceI][n++] = globalOwn;
+ faceStencil[faceI][n++] = globalNei;
+ forAllConstIter(labelHashSet, faceStencilSet, iter)
{
if (iter.key() != globalOwn && iter.key() != globalNei)
{
- stencil_[faceI][n++] = iter.key();
+ faceStencil[faceI][n++] = iter.key();
}
}
- //Pout<< "internalface:" << faceI << " toc:" << faceStencil.toc()
- // << " stencil:" << stencil_[faceI] << endl;
+ //Pout<< "internalface:" << faceI << " toc:" << faceStencilSet.toc()
+ // << " faceStencil:" << faceStencil[faceI] << endl;
}
forAll(patches, patchI)
{
@@ -190,41 +119,40 @@ void Foam::extendedStencil::calcFaceStencils
{
forAll(pp, i)
{
- faceStencil.clear();
- label globalOwn = globalNumbering.toGlobal(own[faceI]);
- faceStencil.insert(globalOwn);
+ faceStencilSet.clear();
+
const labelList& ownCCells = globalCellCells[own[faceI]];
+ label globalOwn = ownCCells[0];
forAll(ownCCells, i)
{
- faceStencil.insert(ownCCells[i]);
+ faceStencilSet.insert(ownCCells[i]);
}
- // Get the coupled cell
- label globalNei = neiGlobal[faceI-mesh.nInternalFaces()];
- faceStencil.insert(globalNei);
+
// And the neighbours of the coupled cell
const labelList& neiCCells =
- neiGlobalCellCells[faceI-mesh.nInternalFaces()];
+ neiGlobalCellCells[faceI-mesh_.nInternalFaces()];
+ label globalNei = neiCCells[0];
forAll(neiCCells, i)
{
- faceStencil.insert(neiCCells[i]);
+ faceStencilSet.insert(neiCCells[i]);
}
// Guarantee owner first, neighbour second.
- stencil_[faceI].setSize(faceStencil.size());
+ faceStencil[faceI].setSize(faceStencilSet.size());
label n = 0;
- stencil_[faceI][n++] = globalOwn;
- stencil_[faceI][n++] = globalNei;
- forAllConstIter(labelHashSet, faceStencil, iter)
+ faceStencil[faceI][n++] = globalOwn;
+ faceStencil[faceI][n++] = globalNei;
+ forAllConstIter(labelHashSet, faceStencilSet, iter)
{
if (iter.key() != globalOwn && iter.key() != globalNei)
{
- stencil_[faceI][n++] = iter.key();
+ faceStencil[faceI][n++] = iter.key();
}
}
//Pout<< "coupledface:" << faceI
- // << " toc:" << faceStencil.toc()
- // << " stencil:" << stencil_[faceI] << endl;
+ // << " toc:" << faceStencilSet.toc()
+ // << " faceStencil:" << faceStencil[faceI] << endl;
faceI++;
}
@@ -233,31 +161,30 @@ void Foam::extendedStencil::calcFaceStencils
{
forAll(pp, i)
{
- faceStencil.clear();
- label globalOwn = globalNumbering.toGlobal(own[faceI]);
- faceStencil.insert(globalOwn);
+ faceStencilSet.clear();
+
const labelList& ownCCells = globalCellCells[own[faceI]];
+ label globalOwn = ownCCells[0];
forAll(ownCCells, i)
{
- faceStencil.insert(ownCCells[i]);
+ faceStencilSet.insert(ownCCells[i]);
}
-
// Guarantee owner first, neighbour second.
- stencil_[faceI].setSize(faceStencil.size());
+ faceStencil[faceI].setSize(faceStencilSet.size());
label n = 0;
- stencil_[faceI][n++] = globalOwn;
- forAllConstIter(labelHashSet, faceStencil, iter)
+ faceStencil[faceI][n++] = globalOwn;
+ forAllConstIter(labelHashSet, faceStencilSet, iter)
{
if (iter.key() != globalOwn)
{
- stencil_[faceI][n++] = iter.key();
+ faceStencil[faceI][n++] = iter.key();
}
}
//Pout<< "boundaryface:" << faceI
- // << " toc:" << faceStencil.toc()
- // << " stencil:" << stencil_[faceI] << endl;
+ // << " toc:" << faceStencilSet.toc()
+ // << " faceStencil:" << faceStencil[faceI] << endl;
faceI++;
}
@@ -266,34 +193,13 @@ void Foam::extendedStencil::calcFaceStencils
}
-// Calculates extended stencil. This is per face
-// - owner
-// - cellCells of owner
-// - neighbour
-// - cellCells of neighbour
-// It comes in two parts:
-// - a map which collects/distributes all necessary data in a compact array
-// - the stencil (a labelList per face) which is a set of indices into this
-// compact array.
-// The compact array is laid out as follows:
-// - first data for current processor (Pstream::myProcNo())
-// - all cells
-// - all boundary faces
-// - then per processor
-// - all used cells and boundary faces
-void Foam::extendedStencil::calcExtendedFaceStencil(const polyMesh& mesh)
+Foam::autoPtr<Foam::mapDistribute> Foam::extendedStencil::calcDistributeMap
+(
+ const globalIndex& globalNumbering,
+ labelListList& faceStencil
+)
{
- const label nBnd = mesh.nFaces()-mesh.nInternalFaces();
-
- // Global numbering for cells and boundary faces
- // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
- globalIndex globalNumbering(mesh.nCells()+nBnd);
-
-
- // Calculate stencil in global cell indices
- // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
- calcFaceStencils(mesh, globalNumbering);
+ const label nBnd = mesh_.nFaces()-mesh_.nInternalFaces();
// Convert stencil to schedule
@@ -309,8 +215,8 @@ void Foam::extendedStencil::calcExtendedFaceStencil(const polyMesh& mesh)
// these are always all needed.
List<Map<label> > globalToProc(Pstream::nProcs());
{
- const labelList& procPatchMap = mesh.globalData().procPatchMap();
- const polyBoundaryMesh& patches = mesh.boundaryMesh();
+ const labelList& procPatchMap = mesh_.globalData().procPatchMap();
+ const polyBoundaryMesh& patches = mesh_.boundaryMesh();
// Presize with (as estimate) size of patch to neighbour.
forAll(procPatchMap, procI)
@@ -325,9 +231,9 @@ void Foam::extendedStencil::calcExtendedFaceStencil(const polyMesh& mesh)
}
// Collect all (non-local) globalcells/faces needed.
- forAll(stencil_, faceI)
+ forAll(faceStencil, faceI)
{
- const labelList& stencilCells = stencil_[faceI];
+ const labelList& stencilCells = faceStencil[faceI];
forAll(stencilCells, i)
{
@@ -359,17 +265,17 @@ void Foam::extendedStencil::calcExtendedFaceStencil(const polyMesh& mesh)
}
- // forAll(globalToProc, procI)
- // {
- // Pout<< "From processor:" << procI << " want cells/faces:" << endl;
- // forAllConstIter(Map<label>, globalToProc[procI], iter)
- // {
- // Pout<< " global:" << iter.key()
- // << " local:" << globalNumbering.toLocal(procI, iter.key())
- // << endl;
- // }
- // Pout<< endl;
- // }
+ //forAll(globalToProc, procI)
+ //{
+ // Pout<< "From processor:" << procI << " want cells/faces:" << endl;
+ // forAllConstIter(Map<label>, globalToProc[procI], iter)
+ // {
+ // Pout<< " global:" << iter.key()
+ // << " local:" << globalNumbering.toLocal(procI, iter.key())
+ // << endl;
+ // }
+ // Pout<< endl;
+ //}
}
@@ -379,19 +285,15 @@ void Foam::extendedStencil::calcExtendedFaceStencil(const polyMesh& mesh)
labelList compactStart(Pstream::nProcs());
compactStart[Pstream::myProcNo()] = 0;
- label nCompact = mesh.nCells()+nBnd;
+ label nCompact = mesh_.nCells()+nBnd;
forAll(compactStart, procI)
{
if (procI != Pstream::myProcNo())
{
compactStart[procI] = nCompact;
nCompact += globalToProc[procI].size();
-
- // Pout<< "Data wanted from " << procI << " starts at "
- // << compactStart[procI] << endl;
}
}
- // Pout<< "Overall cells needed:" << nCompact << endl;
// 3. Find out what to receive/send in compact addressing.
@@ -411,11 +313,6 @@ void Foam::extendedStencil::calcExtendedFaceStencil(const polyMesh& mesh)
i++;
}
- // Pout<< "From proc:" << procI
- // << " I need (globalcells):" << wantedGlobals
- // << " which are my compact:" << recvCompact[procI]
- // << endl;
-
// Send the global cell numbers I need from procI
OPstream str(Pstream::blocking, procI);
str << wantedGlobals;
@@ -424,7 +321,7 @@ void Foam::extendedStencil::calcExtendedFaceStencil(const polyMesh& mesh)
{
recvCompact[procI] =
compactStart[procI]
- + identity(mesh.nCells()+nBnd);
+ + identity(mesh_.nCells()+nBnd);
}
}
labelListList sendCompact(Pstream::nProcs());
@@ -455,9 +352,9 @@ void Foam::extendedStencil::calcExtendedFaceStencil(const polyMesh& mesh)
}
// Convert stencil to compact numbering
- forAll(stencil_, faceI)
+ forAll(faceStencil, faceI)
{
- labelList& stencilCells = stencil_[faceI];
+ labelList& stencilCells = faceStencil[faceI];
forAll(stencilCells, i)
{
@@ -476,10 +373,9 @@ void Foam::extendedStencil::calcExtendedFaceStencil(const polyMesh& mesh)
}
}
- // Pout<< "***stencil_:" << stencil_ << endl;
// Constuct map for distribution of compact data.
- mapPtr_.reset
+ return autoPtr<mapDistribute>
(
new mapDistribute
(
@@ -494,32 +390,10 @@ void Foam::extendedStencil::calcExtendedFaceStencil(const polyMesh& mesh)
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-Foam::extendedStencil::extendedStencil
-(
- const mapDistribute& map,
- const labelListList& stencil
-)
+Foam::extendedStencil::extendedStencil(const polyMesh& mesh)
:
- mapPtr_
- (
- autoPtr<mapDistribute>
- (
- new mapDistribute
- (
- map.constructSize(),
- map.subMap(),
- map.constructMap()
- )
- )
- ),
- stencil_(stencil)
+ mesh_(mesh)
{}
-Foam::extendedStencil::extendedStencil(const polyMesh& mesh)
-{
- calcExtendedFaceStencil(mesh);
-}
-
-
// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencil.H b/src/finiteVolume/fvMesh/extendedStencil/extendedStencil.H
index ad2f4e1b0b0c98ff746733e2b247bbe8fc907ac3..12d3179e62fc1f2a7ffdca439e0164887a952816 100644
--- a/src/finiteVolume/fvMesh/extendedStencil/extendedStencil.H
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencil.H
@@ -42,6 +42,7 @@ Description
SourceFiles
extendedStencil.C
+ extendedStencilTemplates.C
\*---------------------------------------------------------------------------*/
@@ -65,22 +66,33 @@ class globalIndex;
class extendedStencil
{
- // Private data
+protected:
- //- Swap map for getting neigbouring data
- autoPtr<mapDistribute> mapPtr_;
+ // Protected data
- //- Per face the stencil.
- labelListList stencil_;
+ const polyMesh& mesh_;
- // Private Member Functions
+ // Protected Member Functions
+
+ //- Collect cell neighbours into extended stencil
+ void calcFaceStencil
+ (
+ const labelListList& globalCellCells,
+ labelListList& faceStencil
+ );
+
+ //- Calculate distribute map
+ autoPtr<mapDistribute> calcDistributeMap
+ (
+ const globalIndex& globalNumbering,
+ labelListList& faceStencil
+ );
- void calcFaceStencils(const polyMesh&, const globalIndex&);
- //- Calculate the stencil (but not weights)
- void calcExtendedFaceStencil(const polyMesh&);
+private:
+ // Private Member Functions
//- Disallow default bitwise copy construct
extendedStencil(const extendedStencil&);
@@ -93,43 +105,32 @@ public:
// Constructors
- //- Construct from components
- extendedStencil(const mapDistribute& map, const labelListList&);
-
- //- Construct from all cells and boundary faces
- extendedStencil(const polyMesh&);
-
+ //- Construct from mesh
+ explicit extendedStencil(const polyMesh&);
// Member Functions
- //- Return reference to the parallel distribution map
- const mapDistribute& map() const
- {
- return mapPtr_();
- }
-
- //- Return reference to the stencil
- const labelListList& stencil() const
- {
- return stencil_;
- }
-
//- Use map to get the data into stencil order
template<class T>
- void collectData
+ static void collectData
(
+ const mapDistribute& map,
+ const labelListList& stencil,
const GeometricField<T, fvPatchField, volMesh>& fld,
List<List<T> >& stencilFld
- ) const;
+ );
- //- Calculate weighted sum of vol field
+ //- Sum vol field contributions to create face values
template<class Type>
- tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > weightedSum
+ static tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >
+ weightedSum
(
+ const mapDistribute& map,
+ const labelListList& stencil,
const GeometricField<Type, fvPatchField, volMesh>& fld,
const List<List<scalar> >& stencilWeights
- ) const;
+ );
};
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCECStencilObject.C b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCECStencilObject.C
new file mode 100644
index 0000000000000000000000000000000000000000..2c5333bd92fc1082ef4f0ce2f46289e064055039
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCECStencilObject.C
@@ -0,0 +1,38 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "centredCECStencilObject.H"
+
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(centredCECStencilObject, 0);
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCECStencilObject.H b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCECStencilObject.H
new file mode 100644
index 0000000000000000000000000000000000000000..adbc9ec5b3b750ef44527821591d31821688ab8b
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCECStencilObject.H
@@ -0,0 +1,88 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::centredCECStencilObject
+
+Description
+
+SourceFiles
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef centredCECStencilObject_H
+#define centredCECStencilObject_H
+
+#include "extendedCentredStencil.H"
+#include "cellEdgeCellStencil.H"
+#include "MeshObject.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class centredCECStencilObject Declaration
+\*---------------------------------------------------------------------------*/
+
+class centredCECStencilObject
+:
+ public MeshObject<fvMesh, centredCECStencilObject>,
+ public extendedCentredStencil
+{
+
+public:
+
+ TypeName("centredCECStencil");
+
+ // Constructors
+
+ //- Construct from uncompacted face stencil
+ explicit centredCECStencilObject
+ (
+ const fvMesh& mesh
+ )
+ :
+ MeshObject<fvMesh, centredCECStencilObject>(mesh),
+ extendedCentredStencil(cellEdgeCellStencil(mesh))
+ {}
+
+
+ // Destructor
+
+ virtual ~centredCECStencilObject()
+ {}
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCFCStencilObject.C b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCFCStencilObject.C
new file mode 100644
index 0000000000000000000000000000000000000000..e55e3ed2c165ede70fda44abc1a8720ddbda230c
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCFCStencilObject.C
@@ -0,0 +1,38 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "centredCFCStencilObject.H"
+
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(centredCFCStencilObject, 0);
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCFCStencilObject.H b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCFCStencilObject.H
new file mode 100644
index 0000000000000000000000000000000000000000..704f9a92d459394366abf605d3f301276495c0b0
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCFCStencilObject.H
@@ -0,0 +1,88 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::centredCFCStencilObject
+
+Description
+
+SourceFiles
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef centredCFCStencilObject_H
+#define centredCFCStencilObject_H
+
+#include "extendedCentredStencil.H"
+#include "cellFaceCellStencil.H"
+#include "MeshObject.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class centredCFCStencilObject Declaration
+\*---------------------------------------------------------------------------*/
+
+class centredCFCStencilObject
+:
+ public MeshObject<fvMesh, centredCFCStencilObject>,
+ public extendedCentredStencil
+{
+
+public:
+
+ TypeName("centredCFCStencil");
+
+ // Constructors
+
+ //- Construct from uncompacted face stencil
+ explicit centredCFCStencilObject
+ (
+ const fvMesh& mesh
+ )
+ :
+ MeshObject<fvMesh, centredCFCStencilObject>(mesh),
+ extendedCentredStencil(cellFaceCellStencil(mesh))
+ {}
+
+
+ // Destructor
+
+ virtual ~centredCFCStencilObject()
+ {}
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCPCStencilObject.C b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCPCStencilObject.C
new file mode 100644
index 0000000000000000000000000000000000000000..09af9fb3ef8b3bdd7c8f6e23bd6f08698c89e8d4
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCPCStencilObject.C
@@ -0,0 +1,38 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "centredCPCStencilObject.H"
+
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(centredCPCStencilObject, 0);
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCPCStencilObject.H b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCPCStencilObject.H
new file mode 100644
index 0000000000000000000000000000000000000000..049d214437ae2441d77d30b95ff22b3231a9c57c
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredCPCStencilObject.H
@@ -0,0 +1,88 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::centredCPCStencilObject
+
+Description
+
+SourceFiles
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef centredCPCStencilObject_H
+#define centredCPCStencilObject_H
+
+#include "extendedCentredStencil.H"
+#include "cellPointCellStencil.H"
+#include "MeshObject.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class centredCPCStencilObject Declaration
+\*---------------------------------------------------------------------------*/
+
+class centredCPCStencilObject
+:
+ public MeshObject<fvMesh, centredCPCStencilObject>,
+ public extendedCentredStencil
+{
+
+public:
+
+ TypeName("centredCFCStencil");
+
+ // Constructors
+
+ //- Construct from uncompacted face stencil
+ explicit centredCPCStencilObject
+ (
+ const fvMesh& mesh
+ )
+ :
+ MeshObject<fvMesh, centredCPCStencilObject>(mesh),
+ extendedCentredStencil(cellPointCellStencil(mesh))
+ {}
+
+
+ // Destructor
+
+ virtual ~centredCPCStencilObject()
+ {}
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredFECStencilObject.C b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredFECStencilObject.C
new file mode 100644
index 0000000000000000000000000000000000000000..f589b90ff101c640164e3de2bc11a19bd7700425
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredFECStencilObject.C
@@ -0,0 +1,38 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "centredFECStencilObject.H"
+
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(centredFECStencilObject, 0);
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredFECStencilObject.H b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredFECStencilObject.H
new file mode 100644
index 0000000000000000000000000000000000000000..6f7acf78adfb15d03e6f5b93019c40cb65467fe4
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/centredFECStencilObject.H
@@ -0,0 +1,88 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::centredFECStencilObject
+
+Description
+
+SourceFiles
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef centredFECStencilObject_H
+#define centredFECStencilObject_H
+
+#include "extendedCentredStencil.H"
+#include "faceEdgeCellStencil.H"
+#include "MeshObject.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class centredFECStencilObject Declaration
+\*---------------------------------------------------------------------------*/
+
+class centredFECStencilObject
+:
+ public MeshObject<fvMesh, centredFECStencilObject>,
+ public extendedCentredStencil
+{
+
+public:
+
+ TypeName("centredCFCStencil");
+
+ // Constructors
+
+ //- Construct from uncompacted face stencil
+ explicit centredFECStencilObject
+ (
+ const fvMesh& mesh
+ )
+ :
+ MeshObject<fvMesh, centredFECStencilObject>(mesh),
+ extendedCentredStencil(faceEdgeCellStencil(mesh))
+ {}
+
+
+ // Destructor
+
+ virtual ~centredFECStencilObject()
+ {}
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCECStencilObject.C b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCECStencilObject.C
new file mode 100644
index 0000000000000000000000000000000000000000..a869b71d7e133b2c46a466f975a2eac00b5e3782
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCECStencilObject.C
@@ -0,0 +1,38 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "upwindCECStencilObject.H"
+
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(upwindCECStencilObject, 0);
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCECStencilObject.H b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCECStencilObject.H
new file mode 100644
index 0000000000000000000000000000000000000000..f964faf9d0019fbd37ba5c070e697f6663ea3d8a
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCECStencilObject.H
@@ -0,0 +1,89 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::upwindCECStencilObject
+
+Description
+
+SourceFiles
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef upwindCECStencilObject_H
+#define upwindCECStencilObject_H
+
+#include "extendedUpwindStencil.H"
+#include "cellEdgeCellStencil.H"
+#include "MeshObject.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class upwindCECStencilObject Declaration
+\*---------------------------------------------------------------------------*/
+
+class upwindCECStencilObject
+:
+ public MeshObject<fvMesh, upwindCECStencilObject>,
+ public extendedUpwindStencil
+{
+
+public:
+
+ TypeName("upwindCFCStencil");
+
+ // Constructors
+
+ //- Construct from uncompacted face stencil
+ explicit upwindCECStencilObject
+ (
+ const fvMesh& mesh,
+ const scalar minOpposedness
+ )
+ :
+ MeshObject<fvMesh, upwindCECStencilObject>(mesh),
+ extendedUpwindStencil(cellEdgeCellStencil(mesh), minOpposedness)
+ {}
+
+
+ // Destructor
+
+ virtual ~upwindCECStencilObject()
+ {}
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCFCStencilObject.C b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCFCStencilObject.C
new file mode 100644
index 0000000000000000000000000000000000000000..5637b4f9ea10b70e2bd5cdf643d7fc9126652d5b
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCFCStencilObject.C
@@ -0,0 +1,38 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "upwindCFCStencilObject.H"
+
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(upwindCFCStencilObject, 0);
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCFCStencilObject.H b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCFCStencilObject.H
new file mode 100644
index 0000000000000000000000000000000000000000..a38d274b487b466af2b2190e014c984f337176c0
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCFCStencilObject.H
@@ -0,0 +1,89 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::upwindCFCStencilObject
+
+Description
+
+SourceFiles
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef upwindCFCStencilObject_H
+#define upwindCFCStencilObject_H
+
+#include "extendedUpwindStencil.H"
+#include "cellFaceCellStencil.H"
+#include "MeshObject.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class upwindCFCStencilObject Declaration
+\*---------------------------------------------------------------------------*/
+
+class upwindCFCStencilObject
+:
+ public MeshObject<fvMesh, upwindCFCStencilObject>,
+ public extendedUpwindStencil
+{
+
+public:
+
+ TypeName("upwindCFCStencil");
+
+ // Constructors
+
+ //- Construct from uncompacted face stencil
+ explicit upwindCFCStencilObject
+ (
+ const fvMesh& mesh,
+ const scalar minOpposedness
+ )
+ :
+ MeshObject<fvMesh, upwindCFCStencilObject>(mesh),
+ extendedUpwindStencil(cellFaceCellStencil(mesh), minOpposedness)
+ {}
+
+
+ // Destructor
+
+ virtual ~upwindCFCStencilObject()
+ {}
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCPCStencilObject.C b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCPCStencilObject.C
new file mode 100644
index 0000000000000000000000000000000000000000..f75571f3e768bd3bbe8f34f9ac0a58b753ef01de
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCPCStencilObject.C
@@ -0,0 +1,38 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "upwindCPCStencilObject.H"
+
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(upwindCPCStencilObject, 0);
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCPCStencilObject.H b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCPCStencilObject.H
new file mode 100644
index 0000000000000000000000000000000000000000..0d86ef6422b7598447689abbc9befaf5d5a782d1
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindCPCStencilObject.H
@@ -0,0 +1,89 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::upwindCPCStencilObject
+
+Description
+
+SourceFiles
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef upwindCPCStencilObject_H
+#define upwindCPCStencilObject_H
+
+#include "extendedUpwindStencil.H"
+#include "cellPointCellStencil.H"
+#include "MeshObject.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class upwindCPCStencilObject Declaration
+\*---------------------------------------------------------------------------*/
+
+class upwindCPCStencilObject
+:
+ public MeshObject<fvMesh, upwindCPCStencilObject>,
+ public extendedUpwindStencil
+{
+
+public:
+
+ TypeName("upwindCFCStencil");
+
+ // Constructors
+
+ //- Construct from uncompacted face stencil
+ explicit upwindCPCStencilObject
+ (
+ const fvMesh& mesh,
+ const scalar minOpposedness
+ )
+ :
+ MeshObject<fvMesh, upwindCPCStencilObject>(mesh),
+ extendedUpwindStencil(cellPointCellStencil(mesh), minOpposedness)
+ {}
+
+
+ // Destructor
+
+ virtual ~upwindCPCStencilObject()
+ {}
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindFECStencilObject.C b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindFECStencilObject.C
new file mode 100644
index 0000000000000000000000000000000000000000..43b71817c2103ca7ea8809f4bd76adbd9aea861e
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindFECStencilObject.C
@@ -0,0 +1,38 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "upwindFECStencilObject.H"
+
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(upwindFECStencilObject, 0);
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindFECStencilObject.H b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindFECStencilObject.H
new file mode 100644
index 0000000000000000000000000000000000000000..5b2627ebbe2acb4083fddf5274e481e84f8dac73
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilMeshObjects/upwindFECStencilObject.H
@@ -0,0 +1,89 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::upwindFECStencilObject
+
+Description
+
+SourceFiles
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef upwindFECStencilObject_H
+#define upwindFECStencilObject_H
+
+#include "extendedUpwindStencil.H"
+#include "faceEdgeCellStencil.H"
+#include "MeshObject.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class upwindFECStencilObject Declaration
+\*---------------------------------------------------------------------------*/
+
+class upwindFECStencilObject
+:
+ public MeshObject<fvMesh, upwindFECStencilObject>,
+ public extendedUpwindStencil
+{
+
+public:
+
+ TypeName("upwindCFCStencil");
+
+ // Constructors
+
+ //- Construct from uncompacted face stencil
+ explicit upwindFECStencilObject
+ (
+ const fvMesh& mesh,
+ const scalar minOpposedness
+ )
+ :
+ MeshObject<fvMesh, upwindFECStencilObject>(mesh),
+ extendedUpwindStencil(faceEdgeCellStencil(mesh), minOpposedness)
+ {}
+
+
+ // Destructor
+
+ virtual ~upwindFECStencilObject()
+ {}
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilTemplates.C b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilTemplates.C
index 31adb48392ac34d125b7592394ffc1b1d7f57ab6..480cf15a2a9412224c5e465dafde0fb9fc766b31 100644
--- a/src/finiteVolume/fvMesh/extendedStencil/extendedStencilTemplates.C
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedStencilTemplates.C
@@ -31,12 +31,14 @@ License
template<class Type>
void Foam::extendedStencil::collectData
(
+ const mapDistribute& map,
+ const labelListList& stencil,
const GeometricField<Type, fvPatchField, volMesh>& fld,
List<List<Type> >& stencilFld
-) const
+)
{
// 1. Construct cell data in compact addressing
- List<Type> compactFld(map().constructSize(), pTraits<Type>::zero);
+ List<Type> compactFld(map.constructSize(), pTraits<Type>::zero);
// Insert my internal values
forAll(fld, cellI)
@@ -56,14 +58,14 @@ void Foam::extendedStencil::collectData
}
// Do all swapping
- map().distribute(compactFld);
+ map.distribute(compactFld);
// 2. Pull to stencil
- stencilFld.setSize(stencil_.size());
+ stencilFld.setSize(stencil.size());
- forAll(stencil_, faceI)
+ forAll(stencil, faceI)
{
- const labelList& compactCells = stencil_[faceI];
+ const labelList& compactCells = stencil[faceI];
stencilFld[faceI].setSize(compactCells.size());
@@ -79,15 +81,17 @@ template<class Type>
Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh> >
Foam::extendedStencil::weightedSum
(
+ const mapDistribute& map,
+ const labelListList& stencil,
const GeometricField<Type, fvPatchField, volMesh>& fld,
const List<List<scalar> >& stencilWeights
-) const
+)
{
const fvMesh& mesh = fld.mesh();
// Collect internal and boundary values
List<List<Type> > stencilFld;
- collectData(fld, stencilFld);
+ collectData(map, stencil, fld, stencilFld);
tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > tsfCorr
(
@@ -110,6 +114,7 @@ Foam::extendedStencil::weightedSum
);
GeometricField<Type, fvsPatchField, surfaceMesh>& sf = tsfCorr();
+ // Internal faces
for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
{
const List<Type>& stField = stencilFld[faceI];
@@ -121,33 +126,30 @@ Foam::extendedStencil::weightedSum
}
}
- // Coupled boundaries
- /*
+ // Boundaries. Either constrained or calculated so assign value
+ // directly (instead of nicely using operator==)
typename GeometricField<Type, fvsPatchField, surfaceMesh>::
GeometricBoundaryField& bSfCorr = sf.boundaryField();
+
forAll(bSfCorr, patchi)
{
fvsPatchField<Type>& pSfCorr = bSfCorr[patchi];
- if (pSfCorr.coupled())
+ label faceI = pSfCorr.patch().patch().start();
+
+ forAll(pSfCorr, i)
{
- label faceI = pSfCorr.patch().patch().start();
+ const List<Type>& stField = stencilFld[faceI];
+ const List<scalar>& stWeight = stencilWeights[faceI];
- forAll(pSfCorr, i)
+ forAll(stField, j)
{
- const List<Type>& stField = stencilFld[faceI];
- const List<scalar>& stWeight = stencilWeights[faceI];
-
- forAll(stField, j)
- {
- pSfCorr[i] += stField[j]*stWeight[j];
- }
-
- faceI++;
+ pSfCorr[i] += stField[j]*stWeight[j];
}
+
+ faceI++;
}
}
- */
return tsfCorr;
}
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedUpwindStencil.C b/src/finiteVolume/fvMesh/extendedStencil/extendedUpwindStencil.C
new file mode 100644
index 0000000000000000000000000000000000000000..64991775ba354cbe5e16ea74c3e6aa5aee978b39
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedUpwindStencil.C
@@ -0,0 +1,429 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "extendedUpwindStencil.H"
+#include "faceStencil.H"
+#include "syncTools.H"
+#include "SortableList.H"
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+void Foam::extendedUpwindStencil::selectOppositeFaces
+(
+ const boolList& nonEmptyFace,
+ const scalar minOpposedness,
+ const label faceI,
+ const label cellI,
+ DynamicList<label>& oppositeFaces
+) const
+{
+ const vectorField& areas = mesh_.faceAreas();
+ const labelList& own = mesh_.faceOwner();
+ const cell& cFaces = mesh_.cells()[cellI];
+
+ SortableList<scalar> opposedness(cFaces.size(), -GREAT);
+
+ // Pick up all the faces that oppose this one.
+ forAll(cFaces, i)
+ {
+ label otherFaceI = cFaces[i];
+
+ if (otherFaceI != faceI && nonEmptyFace[otherFaceI])
+ {
+ if ((own[otherFaceI] == cellI) == (own[faceI] == cellI))
+ {
+ opposedness[i] = -(areas[otherFaceI] & areas[faceI]);
+ }
+ else
+ {
+ opposedness[i] = (areas[otherFaceI] & areas[faceI]);
+ }
+ }
+ }
+
+ label sz = opposedness.size();
+
+ oppositeFaces.clear();
+
+ scalar myAreaSqr = magSqr(areas[faceI]);
+
+ if (myAreaSqr > VSMALL)
+ {
+ forAll(opposedness, i)
+ {
+ opposedness[i] /= myAreaSqr;
+ }
+ // Sort in incrementing order
+ opposedness.sort();
+
+ // Pick largest no matter what
+ oppositeFaces.append(cFaces[opposedness.indices()[sz-1]]);
+
+ for (label i = sz-2; i >= 0; --i)
+ {
+ if (opposedness[i] < minOpposedness)
+ {
+ break;
+ }
+ oppositeFaces.append(cFaces[opposedness.indices()[i]]);
+ }
+ }
+ else
+ {
+ // Sort in incrementing order
+ opposedness.sort();
+
+ // Tiny face. Do what?
+ // Pick largest no matter what
+ oppositeFaces.append(cFaces[opposedness.indices()[sz-1]]);
+ }
+}
+
+
+void Foam::extendedUpwindStencil::transportStencil
+(
+ const boolList& nonEmptyFace,
+ const labelListList& faceStencil,
+ const scalar minOpposedness,
+ const label faceI,
+ const label cellI,
+ const bool stencilHasNeighbour,
+
+ DynamicList<label>& oppositeFaces,
+ labelHashSet& faceStencilSet,
+ labelList& transportedStencil
+) const
+{
+ label globalOwn = faceStencil[faceI][0];
+ label globalNei = -1;
+ if (stencilHasNeighbour && faceStencil[faceI].size() >= 2)
+ {
+ globalNei = faceStencil[faceI][1];
+ }
+
+
+ selectOppositeFaces
+ (
+ nonEmptyFace,
+ minOpposedness,
+ faceI,
+ cellI,
+ oppositeFaces
+ );
+
+ // Collect all stencils of oppositefaces
+ faceStencilSet.clear();
+ forAll(oppositeFaces, i)
+ {
+ const labelList& fStencil = faceStencil[oppositeFaces[i]];
+
+ forAll(fStencil, j)
+ {
+ label globalI = fStencil[j];
+
+ if (globalI != globalOwn && globalI != globalNei)
+ {
+ faceStencilSet.insert(globalI);
+ }
+ }
+ }
+
+ // Add my owner and neighbour first.
+ if (stencilHasNeighbour)
+ {
+ transportedStencil.setSize(faceStencilSet.size()+2);
+ label n = 0;
+ transportedStencil[n++] = globalOwn;
+ transportedStencil[n++] = globalNei;
+
+ forAllConstIter(labelHashSet, faceStencilSet, iter)
+ {
+ if (iter.key() != globalOwn && iter.key() != globalNei)
+ {
+ transportedStencil[n++] = iter.key();
+ }
+ }
+ if (n != transportedStencil.size())
+ {
+ FatalErrorIn("extendedUpwindStencil::transportStencil(..)")
+ << "problem:" << faceStencilSet
+ << abort(FatalError);
+ }
+ }
+ else
+ {
+ transportedStencil.setSize(faceStencilSet.size()+1);
+ label n = 0;
+ transportedStencil[n++] = globalOwn;
+
+ forAllConstIter(labelHashSet, faceStencilSet, iter)
+ {
+ if (iter.key() != globalOwn)
+ {
+ transportedStencil[n++] = iter.key();
+ }
+ }
+ if (n != transportedStencil.size())
+ {
+ FatalErrorIn("extendedUpwindStencil::transportStencil(..)")
+ << "problem:" << faceStencilSet
+ << abort(FatalError);
+ }
+ }
+}
+
+
+void Foam::extendedUpwindStencil::transportStencils
+(
+ const labelListList& faceStencil,
+ const scalar minOpposedness,
+ labelListList& ownStencil,
+ labelListList& neiStencil
+)
+{
+ const polyBoundaryMesh& patches = mesh_.boundaryMesh();
+ const label nBnd = mesh_.nFaces()-mesh_.nInternalFaces();
+ const labelList& own = mesh_.faceOwner();
+ const labelList& nei = mesh_.faceNeighbour();
+
+ // Work arrays
+ DynamicList<label> oppositeFaces;
+ labelHashSet faceStencilSet;
+
+
+ // For quick detection of empty faces
+ boolList nonEmptyFace(mesh_.nFaces(), true);
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
+
+ if (isA<emptyPolyPatch>(pp))
+ {
+ label faceI = pp.start();
+ forAll(pp, i)
+ {
+ nonEmptyFace[faceI++] = false;
+ }
+ }
+ }
+
+
+ // Do the owner side
+ // ~~~~~~~~~~~~~~~~~
+ // stencil is synchronised at entry so no need to swap.
+
+ ownStencil.setSize(mesh_.nFaces());
+
+ // Internal faces
+ for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
+ {
+ transportStencil
+ (
+ nonEmptyFace,
+ faceStencil,
+ minOpposedness,
+ faceI,
+ own[faceI],
+ true, //stencilHasNeighbour
+ oppositeFaces,
+ faceStencilSet,
+ ownStencil[faceI]
+ );
+ }
+ // Boundary faces
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
+ label faceI = pp.start();
+
+ if (pp.coupled())
+ {
+ forAll(pp, i)
+ {
+ transportStencil
+ (
+ nonEmptyFace,
+ faceStencil,
+ minOpposedness,
+ faceI,
+ own[faceI],
+ true, //stencilHasNeighbour
+
+ oppositeFaces,
+ faceStencilSet,
+ ownStencil[faceI]
+ );
+ faceI++;
+ }
+ }
+ else if (!isA<emptyPolyPatch>(pp))
+ {
+ forAll(pp, i)
+ {
+ // faceStencil does not contain neighbour
+ transportStencil
+ (
+ nonEmptyFace,
+ faceStencil,
+ minOpposedness,
+ faceI,
+ own[faceI],
+ false, //stencilHasNeighbour
+
+ oppositeFaces,
+ faceStencilSet,
+ ownStencil[faceI]
+ );
+ faceI++;
+ }
+ }
+ }
+
+
+ // Swap coupled boundary stencil
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ labelListList neiBndStencil(nBnd);
+ for (label faceI = mesh_.nInternalFaces(); faceI < mesh_.nFaces(); faceI++)
+ {
+ neiBndStencil[faceI-mesh_.nInternalFaces()] = ownStencil[faceI];
+ }
+ syncTools::swapBoundaryFaceList(mesh_, neiBndStencil, false);
+
+
+ // Do the neighbour side
+ // ~~~~~~~~~~~~~~~~~~~~~
+ // - internal faces : get opposite faces on neighbour side
+ // - boundary faces : empty
+ // - coupled faces : in neiBndStencil
+
+ neiStencil.setSize(mesh_.nFaces());
+
+ // Internal faces
+ for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
+ {
+ transportStencil
+ (
+ nonEmptyFace,
+ faceStencil,
+ minOpposedness,
+ faceI,
+ nei[faceI],
+ true, //stencilHasNeighbour
+
+ oppositeFaces,
+ faceStencilSet,
+ neiStencil[faceI]
+ );
+ }
+
+ // Boundary faces
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
+ label faceI = pp.start();
+
+ if (pp.coupled())
+ {
+ forAll(pp, i)
+ {
+ neiStencil[faceI].transfer
+ (
+ neiBndStencil[faceI-mesh_.nInternalFaces()]
+ );
+ faceI++;
+ }
+ }
+ else
+ {
+ // Boundary has empty neighbour stencil
+ }
+ }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::extendedUpwindStencil::extendedUpwindStencil
+(
+ const faceStencil& stencil,
+ const scalar minOpposedness
+)
+:
+ extendedStencil(stencil.mesh())
+{
+ //forAll(stencil, faceI)
+ //{
+ // const labelList& fCells = stencil[faceI];
+ //
+ // Pout<< "Face:" << faceI << " at:" << mesh_.faceCentres()[faceI]
+ // << endl;
+ //
+ // forAll(fCells, i)
+ // {
+ // label globalI = fCells[i];
+ //
+ // if (globalI < mesh_.nCells())
+ // {
+ // Pout<< " cell:" << globalI
+ // << " at:" << mesh_.cellCentres()[globalI] << endl;
+ // }
+ // else
+ // {
+ // label faceI = globalI-mesh_.nCells() + mesh_.nInternalFaces();
+ //
+ // Pout<< " boundary:" << faceI
+ // << " at:" << mesh_.faceCentres()[faceI] << endl;
+ // }
+ // }
+ //}
+ //Pout<< endl << endl;
+
+
+ // Transport centred stencil to upwind/downwind face
+ transportStencils
+ (
+ stencil,
+ minOpposedness,
+ ownStencil_,
+ neiStencil_
+ );
+
+ ownMapPtr_ = calcDistributeMap
+ (
+ stencil.globalNumbering(),
+ ownStencil_
+ );
+
+ neiMapPtr_ = calcDistributeMap
+ (
+ stencil.globalNumbering(),
+ neiStencil_
+ );
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedUpwindStencil.H b/src/finiteVolume/fvMesh/extendedStencil/extendedUpwindStencil.H
new file mode 100644
index 0000000000000000000000000000000000000000..62976f893da9f95fa39866d5a7b00af7eefc6019
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedUpwindStencil.H
@@ -0,0 +1,177 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::extendedUpwindStencil
+
+Description
+
+SourceFiles
+ extendedUpwindStencil.C
+ extendedUpwindStencilTemplates.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef extendedUpwindStencil_H
+#define extendedUpwindStencil_H
+
+#include "extendedStencil.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class faceStencil;
+
+/*---------------------------------------------------------------------------*\
+ Class extendedUpwindStencil Declaration
+\*---------------------------------------------------------------------------*/
+
+class extendedUpwindStencil
+:
+ public extendedStencil
+{
+ // Private data
+
+ //- Swap map for getting neigbouring data
+ autoPtr<mapDistribute> ownMapPtr_;
+ autoPtr<mapDistribute> neiMapPtr_;
+
+ //- Per face the stencil.
+ labelListList ownStencil_;
+ labelListList neiStencil_;
+
+
+
+ // Private Member Functions
+
+ //- Find most 'opposite' faces of cell
+ void selectOppositeFaces
+ (
+ const boolList& nonEmptyFace,
+ const scalar minOpposedness,
+ const label faceI,
+ const label cellI,
+ DynamicList<label>& oppositeFaces
+ ) const;
+
+ //- Transport (centred) face stencil to 'opposite' face.
+ void transportStencil
+ (
+ const boolList& nonEmptyFace,
+ const labelListList& faceStencil,
+ const scalar minOpposedness,
+ const label faceI,
+ const label cellI,
+ const bool stencilHasNeighbour,
+
+ DynamicList<label>& oppositeFaces,
+ labelHashSet& faceStencilSet,
+ labelList& transportedStencil
+ ) const;
+
+ //- Transport (centred) face stencil to 'opposite' faces.
+ void transportStencils
+ (
+ const labelListList& faceStencil,
+ const scalar minOpposedness,
+ labelListList& ownStencil,
+ labelListList& neiStencil
+ );
+
+
+ //- Disallow default bitwise copy construct
+ extendedUpwindStencil(const extendedUpwindStencil&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const extendedUpwindStencil&);
+
+
+public:
+
+ // Constructors
+
+ //- Construct from mesh and uncompacted face stencil
+ extendedUpwindStencil
+ (
+ const faceStencil&,
+ const scalar minOpposedness
+ );
+
+
+ // Member Functions
+
+ //- Return reference to the parallel distribution map
+ const mapDistribute& ownMap() const
+ {
+ return ownMapPtr_();
+ }
+
+ //- Return reference to the parallel distribution map
+ const mapDistribute& neiMap() const
+ {
+ return neiMapPtr_();
+ }
+
+ //- Return reference to the stencil
+ const labelListList& ownStencil() const
+ {
+ return ownStencil_;
+ }
+
+ //- Return reference to the stencil
+ const labelListList& neiStencil() const
+ {
+ return neiStencil_;
+ }
+
+ //- Sum vol field contributions to create face values
+ template<class Type>
+ tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > weightedSum
+ (
+ const surfaceScalarField& phi,
+ const GeometricField<Type, fvPatchField, volMesh>& fld,
+ const List<List<scalar> >& ownWeights,
+ const List<List<scalar> >& neiWeights
+ ) const;
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+# include "extendedUpwindStencilTemplates.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/extendedUpwindStencilTemplates.C b/src/finiteVolume/fvMesh/extendedStencil/extendedUpwindStencilTemplates.C
new file mode 100644
index 0000000000000000000000000000000000000000..1cfa79eaa55ff6154579a2db7cbc7f200d6ae20d
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/extendedUpwindStencilTemplates.C
@@ -0,0 +1,138 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "extendedStencil.H"
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+template<class Type>
+Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh> >
+Foam::extendedUpwindStencil::weightedSum
+(
+ const surfaceScalarField& phi,
+ const GeometricField<Type, fvPatchField, volMesh>& fld,
+ const List<List<scalar> >& ownWeights,
+ const List<List<scalar> >& neiWeights
+) const
+{
+ const fvMesh& mesh = fld.mesh();
+
+ // Collect internal and boundary values
+ List<List<Type> > ownFld;
+ collectData(ownMap(), ownStencil(), fld, ownFld);
+ List<List<Type> > neiFld;
+ collectData(neiMap(), neiStencil(), fld, neiFld);
+
+ tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > tsfCorr
+ (
+ new GeometricField<Type, fvsPatchField, surfaceMesh>
+ (
+ IOobject
+ (
+ fld.name(),
+ mesh.time().timeName(),
+ mesh
+ ),
+ mesh,
+ dimensioned<Type>
+ (
+ fld.name(),
+ fld.dimensions(),
+ pTraits<Type>::zero
+ )
+ )
+ );
+ GeometricField<Type, fvsPatchField, surfaceMesh>& sf = tsfCorr();
+
+ // Internal faces
+ for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
+ {
+ if (phi[faceI] > 0)
+ {
+ // Flux out of owner. Use upwind (= owner side) stencil.
+ const List<Type>& stField = ownFld[faceI];
+ const List<scalar>& stWeight = ownWeights[faceI];
+
+ forAll(stField, i)
+ {
+ sf[faceI] += stField[i]*stWeight[i];
+ }
+ }
+ else
+ {
+ const List<Type>& stField = neiFld[faceI];
+ const List<scalar>& stWeight = neiWeights[faceI];
+
+ forAll(stField, i)
+ {
+ sf[faceI] += stField[i]*stWeight[i];
+ }
+ }
+ }
+
+ // Boundaries. Either constrained or calculated so assign value
+ // directly (instead of nicely using operator==)
+ typename GeometricField<Type, fvsPatchField, surfaceMesh>::
+ GeometricBoundaryField& bSfCorr = sf.boundaryField();
+
+ forAll(bSfCorr, patchi)
+ {
+ fvsPatchField<Type>& pSfCorr = bSfCorr[patchi];
+
+ label faceI = pSfCorr.patch().patch().start();
+
+ forAll(pSfCorr, i)
+ {
+ if (phi[faceI] > 0)
+ {
+ // Flux out of owner. Use upwind (= owner side) stencil.
+ const List<Type>& stField = ownFld[faceI];
+ const List<scalar>& stWeight = ownWeights[faceI];
+
+ forAll(stField, j)
+ {
+ pSfCorr[i] += stField[j]*stWeight[j];
+ }
+ }
+ else
+ {
+ const List<Type>& stField = neiFld[faceI];
+ const List<scalar>& stWeight = neiWeights[faceI];
+
+ forAll(stField, j)
+ {
+ pSfCorr[i] += stField[j]*stWeight[j];
+ }
+ }
+ faceI++;
+ }
+ }
+
+ return tsfCorr;
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellEdgeCellStencil.C b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellEdgeCellStencil.C
new file mode 100644
index 0000000000000000000000000000000000000000..a7189c44ed11bc9e278385d7326187cc6686daa9
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellEdgeCellStencil.C
@@ -0,0 +1,209 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "cellEdgeCellStencil.H"
+#include "syncTools.H"
+//#include "meshTools.H"
+//#include "OFstream.H"
+//#include "Time.H"
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+// Calculates per edge the neighbour data (= edgeCells)
+void Foam::cellEdgeCellStencil::calcEdgeBoundaryData
+(
+ const boolList& isValidBFace,
+ const labelList& boundaryEdges,
+ EdgeMap<labelList>& neiGlobal
+) const
+{
+ neiGlobal.resize(2*boundaryEdges.size());
+
+ labelHashSet edgeGlobals;
+
+ forAll(boundaryEdges, i)
+ {
+ label edgeI = boundaryEdges[i];
+
+ neiGlobal.insert
+ (
+ mesh().edges()[edgeI],
+ calcFaceCells
+ (
+ isValidBFace,
+ mesh().edgeFaces(edgeI),
+ edgeGlobals
+ )
+ );
+ }
+
+ syncTools::syncEdgeMap
+ (
+ mesh(),
+ neiGlobal,
+ unionEqOp(),
+ false // apply separation
+ );
+}
+
+
+// Calculates per cell the neighbour data (= cell or boundary in global
+// numbering). First element is always cell itself!
+void Foam::cellEdgeCellStencil::calcCellStencil
+(
+ labelListList& globalCellCells
+) const
+{
+ // Calculate edges on coupled patches
+ labelList boundaryEdges
+ (
+ allCoupledFacesPatch()().meshEdges
+ (
+ mesh().edges(),
+ mesh().pointEdges()
+ )
+ );
+
+ //{
+ // OFstream str(mesh().time().path()/"boundaryEdges.obj");
+ // Pout<< "DUmping boundary edges to " << str.name() << endl;
+ //
+ // label vertI = 0;
+ // forAll(boundaryEdges, i)
+ // {
+ // label edgeI = boundaryEdges[i];
+ // const edge& e = mesh().edges()[edgeI];
+ // const point& p0 = mesh().points()[e[0]];
+ // const point& p1 = mesh().points()[e[1]];
+ //
+ // Pout<< "boundary edge " << edgeI << " between " << p0 << p1
+ // << endl;
+ //
+ // meshTools::writeOBJ(str, p0);
+ // vertI++;
+ // meshTools::writeOBJ(str, p1);
+ // vertI++;
+ // str << "l " << vertI-1 << ' ' << vertI << nl;
+ // }
+ //}
+
+
+ // Mark boundary faces to be included in stencil (i.e. not coupled or empty)
+ boolList isValidBFace;
+ validBoundaryFaces(isValidBFace);
+
+
+ // Swap edgeCells for coupled edges. Note: use EdgeMap for now since we've
+ // got syncTools::syncEdgeMap for those. Should be replaced with Map and
+ // syncTools functionality to handle those.
+ EdgeMap<labelList> neiGlobal;
+ calcEdgeBoundaryData
+ (
+ isValidBFace,
+ boundaryEdges,
+ neiGlobal
+ );
+
+ globalCellCells.setSize(mesh().nCells());
+
+ // Do coupled edges first
+
+ forAll(boundaryEdges, i)
+ {
+ label edgeI = boundaryEdges[i];
+
+ const labelList& eGlobals = neiGlobal[mesh().edges()[edgeI]];
+
+ // Distribute to all edgeCells
+ const labelList& eCells = mesh().edgeCells(edgeI);
+
+ forAll(eCells, j)
+ {
+ label cellI = eCells[j];
+
+ // Insert pGlobals into globalCellCells
+ merge
+ (
+ globalNumbering().toGlobal(cellI),
+ eGlobals,
+ globalCellCells[cellI]
+ );
+ }
+ }
+ neiGlobal.clear();
+
+ // Do remaining edges cells
+ labelHashSet edgeGlobals;
+
+ for (label edgeI = 0; edgeI < mesh().nEdges(); edgeI++)
+ {
+ labelList eGlobals
+ (
+ calcFaceCells
+ (
+ isValidBFace,
+ mesh().edgeFaces(edgeI),
+ edgeGlobals
+ )
+ );
+
+ const labelList& eCells = mesh().edgeCells(edgeI);
+
+ forAll(eCells, j)
+ {
+ label cellI = eCells[j];
+
+ merge
+ (
+ globalNumbering().toGlobal(cellI),
+ eGlobals,
+ globalCellCells[cellI]
+ );
+ }
+ }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::cellEdgeCellStencil::cellEdgeCellStencil(const polyMesh& mesh)
+:
+ faceStencil(mesh)
+{
+ // Calculate per cell the (edge) connected cells (in global numbering)
+ labelListList globalCellCells;
+ calcCellStencil(globalCellCells);
+
+ // Add stencils of neighbouring cells to create faceStencil
+ labelListList faceStencil;
+ calcFaceStencil(globalCellCells, faceStencil);
+
+ // Transfer to *this
+ transfer(faceStencil);
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellEdgeCellStencil.H b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellEdgeCellStencil.H
new file mode 100644
index 0000000000000000000000000000000000000000..4091ffb40f6ab16c9d5a2a2a48624982268f77cd
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellEdgeCellStencil.H
@@ -0,0 +1,95 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::cellEdgeCellStencil
+
+Description
+
+SourceFiles
+ cellEdgeCellStencil.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef cellEdgeCellStencil_H
+#define cellEdgeCellStencil_H
+
+#include "faceStencil.H"
+#include "boolList.H"
+#include "labelHashSet.H"
+#include "Map.H"
+#include "EdgeMap.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class cellEdgeCellStencil Declaration
+\*---------------------------------------------------------------------------*/
+
+class cellEdgeCellStencil
+:
+ public faceStencil
+{
+ // Private Member Functions
+
+ //- Calculates per edge the neighbour data (= edgeCells)
+ void calcEdgeBoundaryData
+ (
+ const boolList& isValidBFace,
+ const labelList& boundaryEdges,
+ EdgeMap<labelList>& neiGlobal
+ ) const;
+
+ void calcCellStencil(labelListList& globalCellCells) const;
+
+
+ //- Disallow default bitwise copy construct
+ cellEdgeCellStencil(const cellEdgeCellStencil&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const cellEdgeCellStencil&);
+
+
+public:
+
+ // Constructors
+
+ //- Construct from all cells and boundary faces
+ explicit cellEdgeCellStencil(const polyMesh&);
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellFaceCellStencil.C b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellFaceCellStencil.C
new file mode 100644
index 0000000000000000000000000000000000000000..cbcc8f1696dc12bdeb6e5b3bc714bf88ada9f77a
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellFaceCellStencil.C
@@ -0,0 +1,167 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "cellFaceCellStencil.H"
+#include "syncTools.H"
+#include "SortableList.H"
+#include "emptyPolyPatch.H"
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+// Calculates per face the neighbour data (= cell or boundary face)
+void Foam::cellFaceCellStencil::calcFaceBoundaryData
+(
+ labelList& neiGlobal
+) const
+{
+ const polyBoundaryMesh& patches = mesh().boundaryMesh();
+ const label nBnd = mesh().nFaces()-mesh().nInternalFaces();
+ const labelList& own = mesh().faceOwner();
+
+ neiGlobal.setSize(nBnd);
+
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
+ label faceI = pp.start();
+
+ if (pp.coupled())
+ {
+ // For coupled faces get the cell on the other side
+ forAll(pp, i)
+ {
+ label bFaceI = faceI-mesh().nInternalFaces();
+ neiGlobal[bFaceI] = globalNumbering().toGlobal(own[faceI]);
+ faceI++;
+ }
+ }
+ else if (isA<emptyPolyPatch>(pp))
+ {
+ forAll(pp, i)
+ {
+ label bFaceI = faceI-mesh().nInternalFaces();
+ neiGlobal[bFaceI] = -1;
+ faceI++;
+ }
+ }
+ else
+ {
+ // For noncoupled faces get the boundary face.
+ forAll(pp, i)
+ {
+ label bFaceI = faceI-mesh().nInternalFaces();
+ neiGlobal[bFaceI] =
+ globalNumbering().toGlobal(mesh().nCells()+bFaceI);
+ faceI++;
+ }
+ }
+ }
+ syncTools::swapBoundaryFaceList(mesh(), neiGlobal, false);
+}
+
+
+// Calculates per cell the neighbour data (= cell or boundary in global
+// numbering). First element is always cell itself!
+void Foam::cellFaceCellStencil::calcCellStencil(labelListList& globalCellCells)
+ const
+{
+ const label nBnd = mesh().nFaces()-mesh().nInternalFaces();
+ const labelList& own = mesh().faceOwner();
+ const labelList& nei = mesh().faceNeighbour();
+
+
+ // Calculate coupled neighbour (in global numbering)
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ labelList neiGlobal(nBnd);
+ calcFaceBoundaryData(neiGlobal);
+
+
+ // Determine cellCells in global numbering
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ globalCellCells.setSize(mesh().nCells());
+ forAll(globalCellCells, cellI)
+ {
+ const cell& cFaces = mesh().cells()[cellI];
+
+ labelList& cCells = globalCellCells[cellI];
+
+ cCells.setSize(cFaces.size()+1);
+
+ label nNbr = 0;
+
+ // Myself
+ cCells[nNbr++] = globalNumbering().toGlobal(cellI);
+
+ // Collect neighbouring cells/faces
+ forAll(cFaces, i)
+ {
+ label faceI = cFaces[i];
+
+ if (mesh().isInternalFace(faceI))
+ {
+ label nbrCellI = own[faceI];
+ if (nbrCellI == cellI)
+ {
+ nbrCellI = nei[faceI];
+ }
+ cCells[nNbr++] = globalNumbering().toGlobal(nbrCellI);
+ }
+ else
+ {
+ label nbrCellI = neiGlobal[faceI-mesh().nInternalFaces()];
+ if (nbrCellI != -1)
+ {
+ cCells[nNbr++] = nbrCellI;
+ }
+ }
+ }
+ cCells.setSize(nNbr);
+ }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::cellFaceCellStencil::cellFaceCellStencil(const polyMesh& mesh)
+:
+ faceStencil(mesh)
+{
+ // Calculate per cell the (face) connected cells (in global numbering)
+ labelListList globalCellCells;
+ calcCellStencil(globalCellCells);
+
+ // Add stencils of neighbouring cells to create faceStencil
+ labelListList faceStencil;
+ calcFaceStencil(globalCellCells, faceStencil);
+
+ // Transfer to *this
+ transfer(faceStencil);
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellFaceCellStencil.H b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellFaceCellStencil.H
new file mode 100644
index 0000000000000000000000000000000000000000..93e9d9117aba8dd73803079c0683c05751290d1f
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellFaceCellStencil.H
@@ -0,0 +1,83 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::cellFaceCellStencil
+
+Description
+
+SourceFiles
+ cellFaceCellStencil.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef cellFaceCellStencil_H
+#define cellFaceCellStencil_H
+
+#include "faceStencil.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class cellFaceCellStencil Declaration
+\*---------------------------------------------------------------------------*/
+
+class cellFaceCellStencil
+:
+ public faceStencil
+{
+ // Private Member Functions
+
+ void calcFaceBoundaryData(labelList& neiGlobal) const;
+
+ void calcCellStencil(labelListList& globalCellCells) const;
+
+ //- Disallow default bitwise copy construct
+ cellFaceCellStencil(const cellFaceCellStencil&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const cellFaceCellStencil&);
+
+public:
+
+ // Constructors
+
+ //- Construct from mesh
+ explicit cellFaceCellStencil(const polyMesh& mesh);
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellPointCellStencil.C b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellPointCellStencil.C
new file mode 100644
index 0000000000000000000000000000000000000000..d284c2634104cfa1e4fcfc69edaedbd6c625910e
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellPointCellStencil.C
@@ -0,0 +1,174 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "cellPointCellStencil.H"
+#include "syncTools.H"
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+// Calculates per point the neighbour data (= pointCells)
+void Foam::cellPointCellStencil::calcPointBoundaryData
+(
+ const boolList& isValidBFace,
+ const labelList& boundaryPoints,
+ Map<labelList>& neiGlobal
+) const
+{
+ neiGlobal.resize(2*boundaryPoints.size());
+
+ labelHashSet pointGlobals;
+
+ forAll(boundaryPoints, i)
+ {
+ label pointI = boundaryPoints[i];
+
+ neiGlobal.insert
+ (
+ pointI,
+ calcFaceCells
+ (
+ isValidBFace,
+ mesh().pointFaces()[pointI],
+ pointGlobals
+ )
+ );
+ }
+
+ syncTools::syncPointMap
+ (
+ mesh(),
+ neiGlobal,
+ unionEqOp(),
+ false // apply separation
+ );
+}
+
+
+// Calculates per cell the neighbour data (= cell or boundary in global
+// numbering). First element is always cell itself!
+void Foam::cellPointCellStencil::calcCellStencil
+(
+ labelListList& globalCellCells
+) const
+{
+ // Calculate points on coupled patches
+ labelList boundaryPoints(allCoupledFacesPatch()().meshPoints());
+
+
+ // Mark boundary faces to be included in stencil (i.e. not coupled or empty)
+ boolList isValidBFace;
+ validBoundaryFaces(isValidBFace);
+
+
+ // Swap pointCells for coupled points
+ Map<labelList> neiGlobal;
+ calcPointBoundaryData
+ (
+ isValidBFace,
+ boundaryPoints,
+ neiGlobal
+ );
+
+ globalCellCells.setSize(mesh().nCells());
+
+ // Do coupled points first
+
+ forAll(boundaryPoints, i)
+ {
+ label pointI = boundaryPoints[i];
+
+ const labelList& pGlobals = neiGlobal[pointI];
+
+ // Distribute to all pointCells
+ const labelList& pCells = mesh().pointCells(pointI);
+
+ forAll(pCells, j)
+ {
+ label cellI = pCells[j];
+
+ // Insert pGlobals into globalCellCells
+ merge
+ (
+ globalNumbering().toGlobal(cellI),
+ pGlobals,
+ globalCellCells[cellI]
+ );
+ }
+ }
+ neiGlobal.clear();
+
+ // Do remaining points cells
+ labelHashSet pointGlobals;
+
+ for (label pointI = 0; pointI < mesh().nPoints(); pointI++)
+ {
+ labelList pGlobals
+ (
+ calcFaceCells
+ (
+ isValidBFace,
+ mesh().pointFaces()[pointI],
+ pointGlobals
+ )
+ );
+
+ const labelList& pCells = mesh().pointCells(pointI);
+
+ forAll(pCells, j)
+ {
+ label cellI = pCells[j];
+
+ merge
+ (
+ globalNumbering().toGlobal(cellI),
+ pGlobals,
+ globalCellCells[cellI]
+ );
+ }
+ }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::cellPointCellStencil::cellPointCellStencil(const polyMesh& mesh)
+:
+ faceStencil(mesh)
+{
+ // Calculate per cell the (point) connected cells (in global numbering)
+ labelListList globalCellCells;
+ calcCellStencil(globalCellCells);
+
+ // Add stencils of neighbouring cells to create faceStencil
+ labelListList faceStencil;
+ calcFaceStencil(globalCellCells, faceStencil);
+
+ // Transfer to *this
+ transfer(faceStencil);
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellPointCellStencil.H b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellPointCellStencil.H
new file mode 100644
index 0000000000000000000000000000000000000000..301a430899847a0843c28dadf4b0b923e591a501
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/cellPointCellStencil.H
@@ -0,0 +1,94 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::cellPointCellStencil
+
+Description
+
+SourceFiles
+ cellPointCellStencil.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef cellPointCellStencil_H
+#define cellPointCellStencil_H
+
+#include "faceStencil.H"
+#include "boolList.H"
+#include "labelHashSet.H"
+#include "Map.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class cellPointCellStencil Declaration
+\*---------------------------------------------------------------------------*/
+
+class cellPointCellStencil
+:
+ public faceStencil
+{
+ // Private Member Functions
+
+ //- Calculates per point the neighbour data (= pointCells)
+ void calcPointBoundaryData
+ (
+ const boolList& isValidBFace,
+ const labelList& boundaryPoints,
+ Map<labelList>& neiGlobal
+ ) const;
+
+ void calcCellStencil(labelListList& globalCellCells) const;
+
+
+ //- Disallow default bitwise copy construct
+ cellPointCellStencil(const cellPointCellStencil&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const cellPointCellStencil&);
+
+
+public:
+
+ // Constructors
+
+ //- Construct from all cells and boundary faces
+ explicit cellPointCellStencil(const polyMesh&);
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/faceStencil/faceEdgeCellStencil.C b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/faceEdgeCellStencil.C
new file mode 100644
index 0000000000000000000000000000000000000000..ac78ba76a41409f0e7bcc9bca06271e963eefcb3
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/faceEdgeCellStencil.C
@@ -0,0 +1,445 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "faceEdgeCellStencil.H"
+#include "syncTools.H"
+#include "emptyPolyPatch.H"
+//#include "meshTools.H"
+//#include "OFstream.H"
+//#include "Time.H"
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+// Calculates per edge the neighbour data (= edgeCells)
+void Foam::faceEdgeCellStencil::calcEdgeBoundaryData
+(
+ const boolList& isValidBFace,
+ const labelList& boundaryEdges,
+ EdgeMap<labelList>& neiGlobal
+) const
+{
+ neiGlobal.resize(2*boundaryEdges.size());
+
+ labelHashSet edgeGlobals;
+
+ forAll(boundaryEdges, i)
+ {
+ label edgeI = boundaryEdges[i];
+
+ neiGlobal.insert
+ (
+ mesh().edges()[edgeI],
+ calcFaceCells
+ (
+ isValidBFace,
+ mesh().edgeFaces(edgeI),
+ edgeGlobals
+ )
+ );
+ }
+
+ syncTools::syncEdgeMap
+ (
+ mesh(),
+ neiGlobal,
+ unionEqOp(),
+ false // apply separation
+ );
+}
+
+
+// Calculates per face the edge connected data (= cell or boundary in global
+// numbering).
+void Foam::faceEdgeCellStencil::calcFaceStencil
+(
+ labelListList& faceStencil
+) const
+{
+ const polyBoundaryMesh& patches = mesh().boundaryMesh();
+ const label nBnd = mesh().nFaces()-mesh().nInternalFaces();
+ const labelList& own = mesh().faceOwner();
+ const labelList& nei = mesh().faceNeighbour();
+
+
+
+ // Determine neighbouring global cell
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ labelList neiGlobalCell(nBnd);
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
+
+ if (pp.coupled())
+ {
+ label faceI = pp.start();
+
+ forAll(pp, i)
+ {
+ neiGlobalCell[faceI-mesh().nInternalFaces()] =
+ globalNumbering().toGlobal(own[faceI]);
+ faceI++;
+ }
+ }
+ }
+ syncTools::swapBoundaryFaceList(mesh(), neiGlobalCell, false);
+
+
+
+ // Determine on coupled edges the edge cells on the other side
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ // Calculate edges on coupled patches
+ labelList boundaryEdges
+ (
+ allCoupledFacesPatch()().meshEdges
+ (
+ mesh().edges(),
+ mesh().pointEdges()
+ )
+ );
+
+ // Mark boundary faces to be included in stencil (i.e. not coupled or empty)
+ boolList isValidBFace;
+ validBoundaryFaces(isValidBFace);
+
+ // Swap edgeCells for coupled edges. Note: use EdgeMap for now since we've
+ // got syncTools::syncEdgeMap for those. Should be replaced with Map and
+ // syncTools functionality to handle those.
+ EdgeMap<labelList> neiGlobal;
+ calcEdgeBoundaryData
+ (
+ isValidBFace,
+ boundaryEdges,
+ neiGlobal
+ );
+
+
+
+ // Construct stencil in global numbering
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ faceStencil.setSize(mesh().nFaces());
+
+ // Do coupled edges first
+
+ forAll(boundaryEdges, i)
+ {
+ label edgeI = boundaryEdges[i];
+
+ const labelList& eGlobals = neiGlobal[mesh().edges()[edgeI]];
+
+ // Distribute to all edgeFaces
+ const labelList& eFaces = mesh().edgeFaces(edgeI);
+
+ forAll(eFaces, j)
+ {
+ label faceI = eFaces[j];
+
+ // Insert eGlobals into faceStencil.
+ merge(-1, -1, eGlobals, faceStencil[faceI]);
+ }
+ }
+ neiGlobal.clear();
+
+
+ // Do remaining edges by looping over all faces
+
+ // Work arrays
+ DynamicList<label> fEdgesSet;
+ DynamicList<label> eFacesSet;
+ labelHashSet faceStencilSet;
+
+ for (label faceI = 0; faceI < mesh().nInternalFaces(); faceI++)
+ {
+ label globalOwn = globalNumbering().toGlobal(own[faceI]);
+ label globalNei = globalNumbering().toGlobal(nei[faceI]);
+
+ // Convert any existing faceStencil (from coupled edges) into
+ // set and operate on this.
+
+ faceStencilSet.clear();
+
+ // Insert all but global owner and neighbour
+ forAll(faceStencil[faceI], i)
+ {
+ label globalI = faceStencil[faceI][i];
+ if (globalI != globalOwn && globalI != globalNei)
+ {
+ faceStencilSet.insert(globalI);
+ }
+ }
+ faceStencil[faceI].clear();
+
+ // Collect all edge connected (internal) cells
+ const labelList& fEdges = mesh().faceEdges(faceI, fEdgesSet);
+
+ forAll(fEdges, i)
+ {
+ label edgeI = fEdges[i];
+
+ insertFaceCells
+ (
+ globalOwn,
+ globalNei,
+ isValidBFace,
+ mesh().edgeFaces(edgeI, eFacesSet),
+ faceStencilSet
+ );
+ }
+
+ // Extract, guarantee owner first, neighbour second.
+ faceStencil[faceI].setSize(faceStencilSet.size()+2);
+ label n = 0;
+ faceStencil[faceI][n++] = globalOwn;
+ faceStencil[faceI][n++] = globalNei;
+ forAllConstIter(labelHashSet, faceStencilSet, iter)
+ {
+ if (iter.key() == globalOwn || iter.key() == globalNei)
+ {
+ FatalErrorIn("faceEdgeCellStencil::calcFaceStencil(..)")
+ << "problem:" << faceStencilSet
+ << abort(FatalError);
+ }
+ faceStencil[faceI][n++] = iter.key();
+ }
+ }
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
+ label faceI = pp.start();
+
+ if (pp.coupled())
+ {
+ forAll(pp, i)
+ {
+ label globalOwn = globalNumbering().toGlobal(own[faceI]);
+ label globalNei = neiGlobalCell[faceI-mesh().nInternalFaces()];
+
+ // Convert any existing faceStencil (from coupled edges) into
+ // set and operate on this.
+
+ faceStencilSet.clear();
+
+ // Insert all but global owner and neighbour
+ forAll(faceStencil[faceI], i)
+ {
+ label globalI = faceStencil[faceI][i];
+ if (globalI != globalOwn && globalI != globalNei)
+ {
+ faceStencilSet.insert(globalI);
+ }
+ }
+ faceStencil[faceI].clear();
+
+ // Collect all edge connected (internal) cells
+ const labelList& fEdges = mesh().faceEdges(faceI, fEdgesSet);
+
+ forAll(fEdges, i)
+ {
+ label edgeI = fEdges[i];
+
+ insertFaceCells
+ (
+ globalOwn,
+ globalNei,
+ isValidBFace,
+ mesh().edgeFaces(edgeI, eFacesSet),
+ faceStencilSet
+ );
+ }
+
+ // Extract, guarantee owner first, neighbour second.
+ faceStencil[faceI].setSize(faceStencilSet.size()+2);
+ label n = 0;
+ faceStencil[faceI][n++] = globalOwn;
+ faceStencil[faceI][n++] = globalNei;
+ forAllConstIter(labelHashSet, faceStencilSet, iter)
+ {
+ if (iter.key() == globalOwn || iter.key() == globalNei)
+ {
+ FatalErrorIn("faceEdgeCellStencil::calcFaceStencil(..)")
+ << "problem:" << faceStencilSet
+ << abort(FatalError);
+ }
+ faceStencil[faceI][n++] = iter.key();
+ }
+
+ if (n != faceStencil[faceI].size())
+ {
+ FatalErrorIn("problem") << "n:" << n
+ << " size:" << faceStencil[faceI].size()
+ << abort(FatalError);
+ }
+
+ faceI++;
+ }
+ }
+ else if (!isA<emptyPolyPatch>(pp))
+ {
+ forAll(pp, i)
+ {
+ label globalOwn = globalNumbering().toGlobal(own[faceI]);
+
+ // Convert any existing faceStencil (from coupled edges) into
+ // set and operate on this.
+
+ faceStencilSet.clear();
+
+ // Insert all but global owner and neighbour
+ forAll(faceStencil[faceI], i)
+ {
+ label globalI = faceStencil[faceI][i];
+ if (globalI != globalOwn)
+ {
+ faceStencilSet.insert(globalI);
+ }
+ }
+ faceStencil[faceI].clear();
+
+ // Collect all edge connected (internal) cells
+ const labelList& fEdges = mesh().faceEdges(faceI, fEdgesSet);
+
+ forAll(fEdges, i)
+ {
+ label edgeI = fEdges[i];
+
+ insertFaceCells
+ (
+ globalOwn,
+ -1,
+ isValidBFace,
+ mesh().edgeFaces(edgeI, eFacesSet),
+ faceStencilSet
+ );
+ }
+
+ // Extract, guarantee owner first, neighbour second.
+ faceStencil[faceI].setSize(faceStencilSet.size()+1);
+ label n = 0;
+ faceStencil[faceI][n++] = globalOwn;
+ forAllConstIter(labelHashSet, faceStencilSet, iter)
+ {
+ if (iter.key() == globalOwn)
+ {
+ FatalErrorIn("faceEdgeCellStencil::calcFaceStencil(..)")
+ << "problem:" << faceStencilSet
+ << abort(FatalError);
+ }
+ faceStencil[faceI][n++] = iter.key();
+ }
+
+ faceI++;
+ }
+ }
+ }
+
+
+ for (label faceI = 0; faceI < mesh().nInternalFaces(); faceI++)
+ {
+ label globalOwn = globalNumbering().toGlobal(own[faceI]);
+ if (faceStencil[faceI][0] != globalOwn)
+ {
+ FatalErrorIn("faceEdgeCellStencil::calcFaceStencil(..)")
+ << "problem:" << faceStencil[faceI]
+ << " globalOwn:" << globalOwn
+ << abort(FatalError);
+ }
+ label globalNei = globalNumbering().toGlobal(nei[faceI]);
+ if (faceStencil[faceI][1] != globalNei)
+ {
+ FatalErrorIn("faceEdgeCellStencil::calcFaceStencil(..)")
+ << "problem:" << faceStencil[faceI]
+ << " globalNei:" << globalNei
+ << abort(FatalError);
+ }
+ }
+
+
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
+
+ if (pp.coupled())
+ {
+ forAll(pp, i)
+ {
+ label faceI = pp.start()+i;
+
+ label globalOwn = globalNumbering().toGlobal(own[faceI]);
+ if (faceStencil[faceI][0] != globalOwn)
+ {
+ FatalErrorIn("faceEdgeCellStencil::calcFaceStencil(..)")
+ << "problem:" << faceStencil[faceI]
+ << " globalOwn:" << globalOwn
+ << abort(FatalError);
+ }
+ label globalNei = neiGlobalCell[faceI-mesh().nInternalFaces()];
+ if (faceStencil[faceI][1] != globalNei)
+ {
+ FatalErrorIn("faceEdgeCellStencil::calcFaceStencil(..)")
+ << "problem:" << faceStencil[faceI]
+ << " globalNei:" << globalNei
+ << abort(FatalError);
+ }
+ }
+ }
+ else if (!isA<emptyPolyPatch>(pp))
+ {
+ forAll(pp, i)
+ {
+ label faceI = pp.start()+i;
+
+ label globalOwn = globalNumbering().toGlobal(own[faceI]);
+ if (faceStencil[faceI][0] != globalOwn)
+ {
+ FatalErrorIn("faceEdgeCellStencil::calcFaceStencil(..)")
+ << "problem:" << faceStencil[faceI]
+ << " globalOwn:" << globalOwn
+ << abort(FatalError);
+ }
+ }
+ }
+ }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::faceEdgeCellStencil::faceEdgeCellStencil(const polyMesh& mesh)
+:
+ faceStencil(mesh)
+{
+ // Calculate per face the (edge) connected cells (in global numbering)
+ labelListList faceStencil;
+ calcFaceStencil(faceStencil);
+
+ // Transfer to *this
+ transfer(faceStencil);
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/faceStencil/faceEdgeCellStencil.H b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/faceEdgeCellStencil.H
new file mode 100644
index 0000000000000000000000000000000000000000..be5ab855a041ad41a9a19d9b2d4d6824635f073f
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/faceEdgeCellStencil.H
@@ -0,0 +1,96 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::faceEdgeCellStencil
+
+Description
+ All cells connected via edge to face.
+
+SourceFiles
+ faceEdgeCellStencil.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef faceEdgeCellStencil_H
+#define faceEdgeCellStencil_H
+
+#include "faceStencil.H"
+#include "boolList.H"
+#include "labelHashSet.H"
+#include "Map.H"
+#include "EdgeMap.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class faceEdgeCellStencil Declaration
+\*---------------------------------------------------------------------------*/
+
+class faceEdgeCellStencil
+:
+ public faceStencil
+{
+ // Private Member Functions
+
+ //- Calculates per edge the neighbour data (= edgeCells)
+ void calcEdgeBoundaryData
+ (
+ const boolList& isValidBFace,
+ const labelList& boundaryEdges,
+ EdgeMap<labelList>& neiGlobal
+ ) const;
+
+ void calcFaceStencil(labelListList& faceStencil) const;
+
+
+ //- Disallow default bitwise copy construct
+ faceEdgeCellStencil(const faceEdgeCellStencil&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const faceEdgeCellStencil&);
+
+
+public:
+
+ // Constructors
+
+ //- Construct from all cells and boundary faces
+ explicit faceEdgeCellStencil(const polyMesh&);
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/faceStencil/faceStencil.C b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/faceStencil.C
new file mode 100644
index 0000000000000000000000000000000000000000..b3d72c3b5a1d258a49d74cef6f995c05ef6440ec
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/faceStencil.C
@@ -0,0 +1,512 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "faceStencil.H"
+#include "syncTools.H"
+#include "SortableList.H"
+#include "emptyPolyPatch.H"
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+// Merge two list and guarantee global0,global1 are first.
+void Foam::faceStencil::merge
+(
+ const label global0,
+ const label global1,
+ const labelList& listA,
+ labelList& listB
+)
+{
+ sort(listB);
+
+ // See if global0, global1 already present in listB
+ label nGlobalInsert = 0;
+
+ if (global0 != -1)
+ {
+ label index0 = findSortedIndex(listB, global0);
+ if (index0 == -1)
+ {
+ nGlobalInsert++;
+ }
+ }
+
+ if (global1 != -1)
+ {
+ label index1 = findSortedIndex(listB, global1);
+ if (index1 == -1)
+ {
+ nGlobalInsert++;
+ }
+ }
+
+
+ // For all in listA see if they are present
+ label nInsert = 0;
+
+ forAll(listA, i)
+ {
+ label elem = listA[i];
+
+ if (elem != global0 && elem != global1)
+ {
+ if (findSortedIndex(listB, elem) == -1)
+ {
+ nInsert++;
+ }
+ }
+ }
+
+ // Extend B with nInsert and whether global0,global1 need to be inserted.
+ labelList result(listB.size() + nGlobalInsert + nInsert);
+
+ label resultI = 0;
+
+ // Insert global0,1 first
+ if (global0 != -1)
+ {
+ result[resultI++] = global0;
+ }
+ if (global1 != -1)
+ {
+ result[resultI++] = global1;
+ }
+
+
+ // Insert listB
+ forAll(listB, i)
+ {
+ label elem = listB[i];
+
+ if (elem != global0 && elem != global1)
+ {
+ result[resultI++] = elem;
+ }
+ }
+
+
+ // Insert listA
+ forAll(listA, i)
+ {
+ label elem = listA[i];
+
+ if (elem != global0 && elem != global1)
+ {
+ if (findSortedIndex(listB, elem) == -1)
+ {
+ result[resultI++] = elem;
+ }
+ }
+ }
+
+ if (resultI != result.size())
+ {
+ FatalErrorIn("faceStencil::merge(..)")
+ << "problem" << abort(FatalError);
+ }
+
+ listB.transfer(result);
+}
+
+
+// Merge two list and guarantee globalI is first.
+void Foam::faceStencil::merge
+(
+ const label globalI,
+ const labelList& pGlobals,
+ labelList& cCells
+)
+{
+ labelHashSet set;
+ forAll(cCells, i)
+ {
+ if (cCells[i] != globalI)
+ {
+ set.insert(cCells[i]);
+ }
+ }
+
+ forAll(pGlobals, i)
+ {
+ if (pGlobals[i] != globalI)
+ {
+ set.insert(pGlobals[i]);
+ }
+ }
+
+ cCells.setSize(set.size()+1);
+ label n = 0;
+ cCells[n++] = globalI;
+
+ forAllConstIter(labelHashSet, set, iter)
+ {
+ cCells[n++] = iter.key();
+ }
+}
+
+
+void Foam::faceStencil::validBoundaryFaces(boolList& isValidBFace) const
+{
+ const polyBoundaryMesh& patches = mesh().boundaryMesh();
+
+ isValidBFace.setSize(mesh().nFaces()-mesh().nInternalFaces(), true);
+
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
+
+ if (pp.coupled() || isA<emptyPolyPatch>(pp))
+ {
+ label bFaceI = pp.start()-mesh().nInternalFaces();
+ forAll(pp, i)
+ {
+ isValidBFace[bFaceI++] = false;
+ }
+ }
+ }
+}
+
+
+Foam::autoPtr<Foam::indirectPrimitivePatch>
+Foam::faceStencil::allCoupledFacesPatch() const
+{
+ const polyBoundaryMesh& patches = mesh().boundaryMesh();
+
+ label nCoupled = 0;
+
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
+
+ if (pp.coupled())
+ {
+ nCoupled += pp.size();
+ }
+ }
+ labelList coupledFaces(nCoupled);
+ nCoupled = 0;
+
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
+
+ if (pp.coupled())
+ {
+ label faceI = pp.start();
+
+ forAll(pp, i)
+ {
+ coupledFaces[nCoupled++] = faceI++;
+ }
+ }
+ }
+
+ return autoPtr<indirectPrimitivePatch>
+ (
+ new indirectPrimitivePatch
+ (
+ IndirectList<face>
+ (
+ mesh().faces(),
+ coupledFaces
+ ),
+ mesh().points()
+ )
+ );
+}
+
+
+void Foam::faceStencil::unionEqOp::operator()
+(
+ labelList& x,
+ const labelList& y
+) const
+{
+ if (y.size() > 0)
+ {
+ if (x.size() == 0)
+ {
+ x = y;
+ }
+ else
+ {
+ labelHashSet set(x);
+ forAll(y, i)
+ {
+ set.insert(y[i]);
+ }
+ x = set.toc();
+ }
+ }
+}
+
+
+void Foam::faceStencil::insertFaceCells
+(
+ const label exclude0,
+ const label exclude1,
+ const boolList& isValidBFace,
+ const labelList& faceLabels,
+ labelHashSet& globals
+) const
+{
+ const labelList& own = mesh().faceOwner();
+ const labelList& nei = mesh().faceNeighbour();
+
+ forAll(faceLabels, i)
+ {
+ label faceI = faceLabels[i];
+
+ label globalOwn = globalNumbering().toGlobal(own[faceI]);
+ if (globalOwn != exclude0 && globalOwn != exclude1)
+ {
+ globals.insert(globalOwn);
+ }
+
+ if (mesh().isInternalFace(faceI))
+ {
+ label globalNei = globalNumbering().toGlobal(nei[faceI]);
+ if (globalNei != exclude0 && globalNei != exclude1)
+ {
+ globals.insert(globalNei);
+ }
+ }
+ else
+ {
+ label bFaceI = faceI-mesh().nInternalFaces();
+
+ if (isValidBFace[bFaceI])
+ {
+ label globalI = globalNumbering().toGlobal
+ (
+ mesh().nCells()
+ + bFaceI
+ );
+
+ if (globalI != exclude0 && globalI != exclude1)
+ {
+ globals.insert(globalI);
+ }
+ }
+ }
+ }
+}
+
+
+Foam::labelList Foam::faceStencil::calcFaceCells
+(
+ const boolList& isValidBFace,
+ const labelList& faceLabels,
+ labelHashSet& globals
+) const
+{
+ globals.clear();
+
+ insertFaceCells
+ (
+ -1,
+ -1,
+ isValidBFace,
+ faceLabels,
+ globals
+ );
+
+ return globals.toc();
+}
+
+
+// Calculates per face a list of global cell/face indices.
+void Foam::faceStencil::calcFaceStencil
+(
+ const labelListList& globalCellCells,
+ labelListList& faceStencil
+) const
+{
+ const polyBoundaryMesh& patches = mesh_.boundaryMesh();
+ const label nBnd = mesh_.nFaces()-mesh_.nInternalFaces();
+ const labelList& own = mesh_.faceOwner();
+ const labelList& nei = mesh_.faceNeighbour();
+
+
+ // Determine neighbouring global cell Cells
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ labelListList neiGlobalCellCells(nBnd);
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
+
+ if (pp.coupled())
+ {
+ label faceI = pp.start();
+
+ forAll(pp, i)
+ {
+ neiGlobalCellCells[faceI-mesh_.nInternalFaces()] =
+ globalCellCells[own[faceI]];
+ faceI++;
+ }
+ }
+ }
+ syncTools::swapBoundaryFaceList(mesh_, neiGlobalCellCells, false);
+
+
+
+ // Construct stencil in global numbering
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ faceStencil.setSize(mesh_.nFaces());
+
+ labelHashSet faceStencilSet;
+
+ for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
+ {
+ faceStencilSet.clear();
+
+ const labelList& ownCCells = globalCellCells[own[faceI]];
+ label globalOwn = ownCCells[0];
+ // Insert cellCells
+ forAll(ownCCells, i)
+ {
+ faceStencilSet.insert(ownCCells[i]);
+ }
+
+ const labelList& neiCCells = globalCellCells[nei[faceI]];
+ label globalNei = neiCCells[0];
+ // Insert cellCells
+ forAll(neiCCells, i)
+ {
+ faceStencilSet.insert(neiCCells[i]);
+ }
+
+ // Guarantee owner first, neighbour second.
+ faceStencil[faceI].setSize(faceStencilSet.size());
+ label n = 0;
+ faceStencil[faceI][n++] = globalOwn;
+ faceStencil[faceI][n++] = globalNei;
+ forAllConstIter(labelHashSet, faceStencilSet, iter)
+ {
+ if (iter.key() != globalOwn && iter.key() != globalNei)
+ {
+ faceStencil[faceI][n++] = iter.key();
+ }
+ }
+ //Pout<< "internalface:" << faceI << " toc:" << faceStencilSet.toc()
+ // << " faceStencil:" << faceStencil[faceI] << endl;
+ }
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
+ label faceI = pp.start();
+
+ if (pp.coupled())
+ {
+ forAll(pp, i)
+ {
+ faceStencilSet.clear();
+
+ const labelList& ownCCells = globalCellCells[own[faceI]];
+ label globalOwn = ownCCells[0];
+ forAll(ownCCells, i)
+ {
+ faceStencilSet.insert(ownCCells[i]);
+ }
+
+ // And the neighbours of the coupled cell
+ const labelList& neiCCells =
+ neiGlobalCellCells[faceI-mesh_.nInternalFaces()];
+ label globalNei = neiCCells[0];
+ forAll(neiCCells, i)
+ {
+ faceStencilSet.insert(neiCCells[i]);
+ }
+
+ // Guarantee owner first, neighbour second.
+ faceStencil[faceI].setSize(faceStencilSet.size());
+ label n = 0;
+ faceStencil[faceI][n++] = globalOwn;
+ faceStencil[faceI][n++] = globalNei;
+ forAllConstIter(labelHashSet, faceStencilSet, iter)
+ {
+ if (iter.key() != globalOwn && iter.key() != globalNei)
+ {
+ faceStencil[faceI][n++] = iter.key();
+ }
+ }
+
+ //Pout<< "coupledface:" << faceI
+ // << " toc:" << faceStencilSet.toc()
+ // << " faceStencil:" << faceStencil[faceI] << endl;
+
+ faceI++;
+ }
+ }
+ else if (!isA<emptyPolyPatch>(pp))
+ {
+ forAll(pp, i)
+ {
+ faceStencilSet.clear();
+
+ const labelList& ownCCells = globalCellCells[own[faceI]];
+ label globalOwn = ownCCells[0];
+ forAll(ownCCells, i)
+ {
+ faceStencilSet.insert(ownCCells[i]);
+ }
+
+ // Guarantee owner first
+ faceStencil[faceI].setSize(faceStencilSet.size());
+ label n = 0;
+ faceStencil[faceI][n++] = globalOwn;
+ forAllConstIter(labelHashSet, faceStencilSet, iter)
+ {
+ if (iter.key() != globalOwn)
+ {
+ faceStencil[faceI][n++] = iter.key();
+ }
+ }
+
+ //Pout<< "boundaryface:" << faceI
+ // << " toc:" << faceStencilSet.toc()
+ // << " faceStencil:" << faceStencil[faceI] << endl;
+
+ faceI++;
+ }
+ }
+ }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::faceStencil::faceStencil(const polyMesh& mesh)
+:
+ mesh_(mesh),
+ globalNumbering_(mesh_.nCells()+mesh_.nFaces()-mesh_.nInternalFaces())
+{}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/fvMesh/extendedStencil/faceStencil/faceStencil.H b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/faceStencil.H
new file mode 100644
index 0000000000000000000000000000000000000000..058d1e61ce56cec2ef8e159e1804d6b5dbd8583b
--- /dev/null
+++ b/src/finiteVolume/fvMesh/extendedStencil/faceStencil/faceStencil.H
@@ -0,0 +1,160 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::faceStencil
+
+Description
+ baseclass for extended face-cell addressing.
+
+SourceFiles
+ faceStencil.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef faceStencil_H
+#define faceStencil_H
+
+#include "globalIndex.H"
+#include "boolList.H"
+#include "labelHashSet.H"
+#include "indirectPrimitivePatch.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class polyMesh;
+
+/*---------------------------------------------------------------------------*\
+ Class faceStencil Declaration
+\*---------------------------------------------------------------------------*/
+
+class faceStencil
+:
+ public labelListList
+{
+ // Private data
+
+ const polyMesh& mesh_;
+
+ //- Global numbering for cells and boundary faces
+ const globalIndex globalNumbering_;
+
+
+ // Private Member Functions
+
+ //- Disallow default bitwise assignment
+ void operator=(const faceStencil&);
+
+
+
+protected:
+
+ //- Merge two lists.
+ static void merge
+ (
+ const label,
+ const label,
+ const labelList&,
+ labelList&
+ );
+
+ //- Merge two lists.
+ static void merge(const label, const labelList&, labelList&);
+
+ //- Valid boundary faces (not empty and not coupled)
+ void validBoundaryFaces(boolList& isValidBFace) const;
+
+ //- Return patch of all coupled faces.
+ autoPtr<indirectPrimitivePatch> allCoupledFacesPatch() const;
+
+ //- Combine operator for labelLists
+ class unionEqOp
+ {
+ public:
+ void operator()( labelList& x, const labelList& y ) const;
+ };
+
+ //- Collect cell neighbours of faces in global numbering
+ void insertFaceCells
+ (
+ const label exclude0,
+ const label exclude1,
+ const boolList& nonEmptyFace,
+ const labelList& faceLabels,
+ labelHashSet& globals
+ ) const;
+
+ //- Collect cell neighbours of faces in global numbering
+ labelList calcFaceCells
+ (
+ const boolList& nonEmptyFace,
+ const labelList& faceLabels,
+ labelHashSet& globals
+ ) const;
+
+
+ //- Collect cell neighbours into extended stencil
+ void calcFaceStencil
+ (
+ const labelListList& globalCellCells,
+ labelListList& faceStencil
+ ) const;
+
+
+public:
+
+ // Constructors
+
+ //- Construct from mesh
+ explicit faceStencil(const polyMesh&);
+
+
+ // Member Functions
+
+ const polyMesh& mesh() const
+ {
+ return mesh_;
+ }
+
+ //- Global numbering for cells and boundary faces
+ const globalIndex& globalNumbering() const
+ {
+ return globalNumbering_;
+ }
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearFit/linearFit.C b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearFit/linearFit.C
new file mode 100644
index 0000000000000000000000000000000000000000..a7dead621f9dfa472a2ad370487a990511a72e2c
--- /dev/null
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearFit/linearFit.C
@@ -0,0 +1,36 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "linearFit.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ makeSurfaceInterpolationScheme(linearFit);
+}
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearFit/linearFit.H b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearFit/linearFit.H
new file mode 100644
index 0000000000000000000000000000000000000000..882f3e94bcc2dc046c226062820b1144a1a4e5ff
--- /dev/null
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearFit/linearFit.H
@@ -0,0 +1,138 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ linearFit
+
+Description
+ Linear fit interpolation scheme which applies an explicit correction to
+ linear.
+
+SourceFiles
+ linearFit.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef linearFit_H
+#define linearFit_H
+
+#include "linear.H"
+#include "linearFitData.H"
+#include "extendedStencil.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class linearFit Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class Type>
+class linearFit
+:
+ public linear<Type>
+{
+ // Private Data
+ const scalar centralWeight_;
+
+ // Private Member Functions
+
+ //- Disallow default bitwise copy construct
+ linearFit(const linearFit&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const linearFit&);
+
+
+public:
+
+ //- Runtime type information
+ TypeName("linearFit");
+
+
+ // Constructors
+
+ //- Construct from mesh and Istream
+ linearFit(const fvMesh& mesh, Istream& is)
+ :
+ linear<Type>(mesh),
+ centralWeight_(readScalar(is))
+ {}
+
+
+ //- Construct from mesh, faceFlux and Istream
+ linearFit
+ (
+ const fvMesh& mesh,
+ const surfaceScalarField& faceFlux,
+ Istream& is
+ )
+ :
+ linear<Type>(mesh),
+ centralWeight_(readScalar(is))
+ {}
+
+
+ // Member Functions
+
+ //- Return true if this scheme uses an explicit correction
+ virtual bool corrected() const
+ {
+ return true;
+ }
+
+ //- Return the explicit correction to the face-interpolate
+ virtual tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >
+ correction
+ (
+ const GeometricField<Type, fvPatchField, volMesh>& vf
+ ) const
+ {
+ const fvMesh& mesh = this->mesh();
+
+ const linearFitData& cfd = linearFitData::New
+ (
+ mesh,
+ centralWeight_
+ );
+
+ const extendedStencil& stencil = cfd.stencil();
+ const List<scalarList>& f = cfd.fit();
+
+ return stencil.weightedSum(vf, f);
+ }
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearFit/linearFitData.C b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearFit/linearFitData.C
new file mode 100644
index 0000000000000000000000000000000000000000..996d915e1cc56aad3734772bbf7be4e5667f0e90
--- /dev/null
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearFit/linearFitData.C
@@ -0,0 +1,371 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "linearFitData.H"
+#include "surfaceFields.H"
+#include "volFields.H"
+#include "SVD.H"
+#include "syncTools.H"
+
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(linearFitData, 0);
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+static int count = 0;
+
+Foam::linearFitData::linearFitData
+(
+ const fvMesh& mesh,
+ const scalar cWeight
+)
+:
+ MeshObject<fvMesh, linearFitData>(mesh),
+ centralWeight_(cWeight),
+# ifdef SPHERICAL_GEOMETRY
+ dim_(2),
+# else
+ dim_(mesh.nGeometricD()),
+# endif
+ minSize_
+ (
+ dim_ == 1 ? 2 :
+ dim_ == 2 ? 3 :
+ dim_ == 3 ? 4 : 0
+ ),
+ stencil_(mesh),
+ fit_(mesh.nInternalFaces())
+{
+ if (debug)
+ {
+ Info << "Contructing linearFitData" << endl;
+ }
+
+ // check input
+ if (centralWeight_ < 1 - SMALL)
+ {
+ FatalErrorIn("linearFitData::linearFitData")
+ << "centralWeight requested = " << centralWeight_
+ << " should not be less than one"
+ << exit(FatalError);
+ }
+
+ if (minSize_ == 0)
+ {
+ FatalErrorIn("linearFitSnGradData")
+ << " dimension must be 1,2 or 3, not" << dim_ << exit(FatalError);
+ }
+
+ // store the polynomial size for each cell to write out
+ surfaceScalarField interpPolySize
+ (
+ IOobject
+ (
+ "linearFitInterpPolySize",
+ "constant",
+ mesh,
+ IOobject::NO_READ,
+ IOobject::NO_WRITE
+ ),
+ mesh,
+ dimensionedScalar("linearFitInterpPolySize", dimless, scalar(0))
+ );
+
+ // Get the cell/face centres in stencil order.
+ // Centred face stencils no good for triangles of tets. Need bigger stencils
+ List<List<point> > stencilPoints(stencil_.stencil().size());
+ stencil_.collectData
+ (
+ mesh.C(),
+ stencilPoints
+ );
+
+ // find the fit coefficients for every face in the mesh
+
+ for(label faci = 0; faci < mesh.nInternalFaces(); faci++)
+ {
+ interpPolySize[faci] = calcFit(stencilPoints[faci], faci);
+ }
+
+ Pout<< "count = " << count << endl;
+
+ if (debug)
+ {
+ Info<< "linearFitData::linearFitData() :"
+ << "Finished constructing polynomialFit data"
+ << endl;
+
+ interpPolySize.write();
+ }
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+void Foam::linearFitData::findFaceDirs
+(
+ vector& idir, // value changed in return
+ vector& jdir, // value changed in return
+ vector& kdir, // value changed in return
+ const fvMesh& mesh,
+ const label faci
+)
+{
+ idir = mesh.Sf()[faci];
+ //idir = mesh.C()[mesh.neighbour()[faci]] - mesh.C()[mesh.owner()[faci]];
+ idir /= mag(idir);
+
+# ifndef SPHERICAL_GEOMETRY
+ if (mesh.nGeometricD() <= 2) // find the normal direcion
+ {
+ if (mesh.directions()[0] == -1)
+ {
+ kdir = vector(1, 0, 0);
+ }
+ else if (mesh.directions()[1] == -1)
+ {
+ kdir = vector(0, 1, 0);
+ }
+ else
+ {
+ kdir = vector(0, 0, 1);
+ }
+ }
+ else // 3D so find a direction in the place of the face
+ {
+ const face& f = mesh.faces()[faci];
+ kdir = mesh.points()[f[0]] - mesh.points()[f[1]];
+ }
+# else
+ // Spherical geometry so kdir is the radial direction
+ kdir = mesh.Cf()[faci];
+# endif
+
+ if (mesh.nGeometricD() == 3)
+ {
+ // Remove the idir component from kdir and normalise
+ kdir -= (idir & kdir)*idir;
+
+ scalar magk = mag(kdir);
+
+ if (magk < SMALL)
+ {
+ FatalErrorIn("findFaceDirs") << " calculated kdir = zero"
+ << exit(FatalError);
+ }
+ else
+ {
+ kdir /= magk;
+ }
+ }
+
+ jdir = kdir ^ idir;
+}
+
+
+Foam::label Foam::linearFitData::calcFit
+(
+ const List<point>& C,
+ const label faci
+)
+{
+ vector idir(1,0,0);
+ vector jdir(0,1,0);
+ vector kdir(0,0,1);
+ findFaceDirs(idir, jdir, kdir, mesh(), faci);
+
+ scalarList wts(C.size(), scalar(1));
+ wts[0] = centralWeight_;
+ wts[1] = centralWeight_;
+
+ point p0 = mesh().faceCentres()[faci];
+ scalar scale = 0;
+
+ // calculate the matrix of the polynomial components
+ scalarRectangularMatrix B(C.size(), minSize_, scalar(0));
+
+ for(label ip = 0; ip < C.size(); ip++)
+ {
+ const point& p = C[ip];
+
+ scalar px = (p - p0)&idir;
+ scalar py = (p - p0)&jdir;
+# ifndef SPHERICAL_GEOMETRY
+ scalar pz = (p - p0)&kdir;
+# else
+ scalar pz = mag(p) - mag(p0);
+# endif
+
+ if (ip == 0)
+ {
+ scale = max(max(mag(px), mag(py)), mag(pz));
+ }
+
+ px /= scale;
+ py /= scale;
+ pz /= scale;
+
+ label is = 0;
+
+ B[ip][is++] = wts[0]*wts[ip];
+ B[ip][is++] = wts[0]*wts[ip]*px;
+
+ if (dim_ >= 2)
+ {
+ B[ip][is++] = wts[ip]*py;
+ }
+ if (dim_ == 3)
+ {
+ B[ip][is++] = wts[ip]*pz;
+ }
+ }
+
+ // Set the fit
+ label stencilSize = C.size();
+ fit_[faci].setSize(stencilSize);
+ scalarList singVals(minSize_);
+ label nSVDzeros = 0;
+
+ const GeometricField<scalar, fvsPatchField, surfaceMesh>& w =
+ mesh().surfaceInterpolation::weights();
+
+ bool goodFit = false;
+ for(int iIt = 0; iIt < 10 && !goodFit; iIt++)
+ {
+ SVD svd(B, SMALL);
+
+ scalar fit0 = wts[0]*wts[0]*svd.VSinvUt()[0][0];
+ scalar fit1 = wts[0]*wts[1]*svd.VSinvUt()[0][1];
+
+ //goodFit = (fit0 > 0 && fit1 > 0);
+
+ goodFit =
+ (mag(fit0 - w[faci])/w[faci] < 0.5)
+ && (mag(fit1 - (1 - w[faci]))/(1 - w[faci]) < 0.5);
+
+ //scalar w0Err = fit0/w[faci];
+ //scalar w1Err = fit1/(1 - w[faci]);
+
+ //goodFit =
+ // (w0Err > 0.5 && w0Err < 1.5)
+ // && (w1Err > 0.5 && w1Err < 1.5);
+
+ if (goodFit)
+ {
+ fit_[faci][0] = fit0;
+ fit_[faci][1] = fit1;
+
+ for(label i=2; i<stencilSize; i++)
+ {
+ fit_[faci][i] = wts[0]*wts[i]*svd.VSinvUt()[0][i];
+ }
+
+ singVals = svd.S();
+ nSVDzeros = svd.nZeros();
+ }
+ else // (not good fit so increase weight in the centre and for linear)
+ {
+ wts[0] *= 10;
+ wts[1] *= 10;
+
+ for(label i = 0; i < B.n(); i++)
+ {
+ B[i][0] *= 10;
+ B[i][1] *= 10;
+ }
+
+ for(label j = 0; j < B.m(); j++)
+ {
+ B[0][j] *= 10;
+ B[1][j] *= 10;
+ }
+ }
+ }
+
+ // static const scalar L = 0.1;
+ // static const scalar R = 0.2;
+
+ // static const scalar beta = 1.0/(R - L);
+ // static const scalar alpha = R*beta;
+
+ if (goodFit)
+ {
+ if ((mag(fit_[faci][0] - w[faci])/w[faci] < 0.15)
+ && (mag(fit_[faci][1] - (1 - w[faci]))/(1 - w[faci]) < 0.15))
+ {
+ count++;
+ //Pout<< "fit " << mag(fit_[faci][0] - w[faci])/w[faci] << " " << mag(fit_[faci][1] - (1 - w[faci]))/(1 - w[faci]) << endl;
+ }
+
+ // scalar limiter =
+ // max
+ // (
+ // min
+ // (
+ // min(alpha - beta*mag(fit_[faci][0] - w[faci])/w[faci], 1),
+ // min(alpha - beta*mag(fit_[faci][1] - (1 - w[faci]))/(1 - w[faci]), 1)
+ // ), 0
+ // );
+
+ // Remove the uncorrected linear coefficients
+ fit_[faci][0] -= w[faci];
+ fit_[faci][1] -= 1 - w[faci];
+
+ // if (limiter < 0.99)
+ // {
+ // for(label i = 0; i < stencilSize; i++)
+ // {
+ // fit_[faci][i] *= limiter;
+ // }
+ // }
+ }
+ else
+ {
+ Pout<< "Could not fit face " << faci
+ << " " << fit_[faci][0] << " " << w[faci]
+ << " " << fit_[faci][1] << " " << 1 - w[faci]<< endl;
+ fit_[faci] = 0;
+ }
+
+ return minSize_ - nSVDzeros;
+}
+
+
+bool Foam::linearFitData::movePoints()
+{
+ notImplemented("linearFitData::movePoints()");
+
+ return true;
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearFit/linearFitData.H b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearFit/linearFitData.H
new file mode 100644
index 0000000000000000000000000000000000000000..7a2bc7d28e31888f059dec2e35e862daf72fc95c
--- /dev/null
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/linearFit/linearFitData.H
@@ -0,0 +1,140 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM 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 2 of the License, or (at your
+ option) any later version.
+
+ OpenFOAM 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 OpenFOAM; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ linearFitData
+
+Description
+ Data for the linear fit correction interpolation scheme
+
+SourceFiles
+ linearFitData.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef linearFitData_H
+#define linearFitData_H
+
+#include "MeshObject.H"
+#include "fvMesh.H"
+#include "extendedStencil.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class globalIndex;
+
+/*---------------------------------------------------------------------------*\
+ Class linearFitData Declaration
+\*---------------------------------------------------------------------------*/
+
+class linearFitData
+:
+ public MeshObject<fvMesh, linearFitData>
+{
+ // Private data
+
+ //- weights for central stencil
+ const scalar centralWeight_;
+
+ //- dimensionality of the geometry
+ const label dim_;
+
+ //- minimum stencil size
+ const label minSize_;
+
+ //- Extended stencil addressing
+ extendedStencil stencil_;
+
+ //- For each cell in the mesh store the values which multiply the
+ // values of the stencil to obtain the gradient for each direction
+ List<scalarList> fit_;
+
+
+ // Private member functions
+
+ //- Find the normal direction and i, j and k directions for face faci
+ static void findFaceDirs
+ (
+ vector& idir, // value changed in return
+ vector& jdir, // value changed in return
+ vector& kdir, // value changed in return
+ const fvMesh& mesh,
+ const label faci
+ );
+
+ label calcFit(const List<point>&, const label faci);
+
+
+public:
+
+ TypeName("linearFitData");
+
+
+ // Constructors
+
+ explicit linearFitData
+ (
+ const fvMesh& mesh,
+ scalar cWeightDim
+ );
+
+
+ // Destructor
+
+ virtual ~linearFitData()
+ {}
+
+
+ // Member functions
+
+
+ //- Return reference to the stencil
+ const extendedStencil& stencil() const
+ {
+ return stencil_;
+ }
+
+ //- Return reference to fit coefficients
+ const List<scalarList>& fit() const
+ {
+ return fit_;
+ }
+
+ //- Delete the data when the mesh moves not implemented
+ virtual bool movePoints();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/quadraticFit/quadraticFit.H b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/quadraticFit/quadraticFit.H
index 9db4996862881c924150ceaf3bbc762087a955e4..ff9a668f3d30a1264fe2cfb2da240595a3b79c2f 100644
--- a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/quadraticFit/quadraticFit.H
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/quadraticFit/quadraticFit.H
@@ -37,9 +37,9 @@ SourceFiles
#ifndef quadraticFit_H
#define quadraticFit_H
-#include "linear.H"
#include "quadraticFitData.H"
-#include "extendedStencil.H"
+#include "linear.H"
+#include "centredCFCStencilObject.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -113,13 +113,19 @@ public:
{
const fvMesh& mesh = this->mesh();
+ const extendedCentredStencil& stencil =
+ centredCFCStencilObject::New
+ (
+ mesh
+ );
+
const quadraticFitData& cfd = quadraticFitData::New
(
mesh,
+ stencil,
centralWeight_
);
- const extendedStencil& stencil = cfd.stencil();
const List<scalarList>& f = cfd.fit();
return stencil.weightedSum(vf, f);
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/quadraticFit/quadraticFitData.C b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/quadraticFit/quadraticFitData.C
index 52cf371c70caf34814d99c9e32d2d13ff22df2e2..c5269877d8d9741295e80c2dba77703f9ed5944c 100644
--- a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/quadraticFit/quadraticFitData.C
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/quadraticFit/quadraticFitData.C
@@ -29,7 +29,7 @@ License
#include "volFields.H"
#include "SVD.H"
#include "syncTools.H"
-
+#include "extendedCentredStencil.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -46,6 +46,7 @@ static int count = 0;
Foam::quadraticFitData::quadraticFitData
(
const fvMesh& mesh,
+ const extendedCentredStencil& stencil,
const scalar cWeight
)
:
@@ -62,7 +63,6 @@ Foam::quadraticFitData::quadraticFitData
dim_ == 2 ? 6 :
dim_ == 3 ? 9 : 0
),
- stencil_(mesh),
fit_(mesh.nInternalFaces())
{
if (debug)
@@ -102,8 +102,8 @@ Foam::quadraticFitData::quadraticFitData
// Get the cell/face centres in stencil order.
// Centred face stencils no good for triangles of tets. Need bigger stencils
- List<List<point> > stencilPoints(stencil_.stencil().size());
- stencil_.collectData
+ List<List<point> > stencilPoints(mesh.nFaces());
+ stencil.collectData
(
mesh.C(),
stencilPoints
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/quadraticFit/quadraticFitData.H b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/quadraticFit/quadraticFitData.H
index 19859ea621cde6d18e6461f60aae6a7c49387d04..b44a2e4d90f4ff84d36dc360f97e34bd7ec000da 100644
--- a/src/finiteVolume/interpolation/surfaceInterpolation/schemes/quadraticFit/quadraticFitData.H
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/schemes/quadraticFit/quadraticFitData.H
@@ -38,14 +38,13 @@ SourceFiles
#include "MeshObject.H"
#include "fvMesh.H"
-#include "extendedStencil.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
-class globalIndex;
+class extendedCentredStencil;
/*---------------------------------------------------------------------------*\
Class quadraticFitData Declaration
@@ -66,9 +65,6 @@ class quadraticFitData
//- minimum stencil size
const label minSize_;
- //- Extended stencil addressing
- extendedStencil stencil_;
-
//- For each cell in the mesh store the values which multiply the
// values of the stencil to obtain the gradient for each direction
List<scalarList> fit_;
@@ -99,6 +95,7 @@ public:
explicit quadraticFitData
(
const fvMesh& mesh,
+ const extendedCentredStencil& stencil,
scalar cWeightDim
);
@@ -112,12 +109,6 @@ public:
// Member functions
- //- Return reference to the stencil
- const extendedStencil& stencil() const
- {
- return stencil_;
- }
-
//- Return reference to fit coefficients
const List<scalarList>& fit() const
{
diff --git a/src/sampling/probes/probes.C b/src/sampling/probes/probes.C
index e0576b60348c4c4c2bb70df5e98bdbee1d24103b..26b07afc595ac313e425116ff082dafd22fd5504 100644
--- a/src/sampling/probes/probes.C
+++ b/src/sampling/probes/probes.C
@@ -240,7 +240,7 @@ bool Foam::probes::checkFieldTypes()
forAll(probeLocations_, probeI)
{
- *sPtr<< setw(w) << probeLocations_[probeI][cmpt];
+ *sPtr<< ' ' << setw(w) << probeLocations_[probeI][cmpt];
}
*sPtr << endl;
}
diff --git a/src/sampling/probes/probesTemplates.C b/src/sampling/probes/probesTemplates.C
index 109f290ce1fc2c9d9e856e55cc3a5828a046f05d..fc90ed1dcab5bf86fb1abb1723e6f287157d5ba4 100644
--- a/src/sampling/probes/probesTemplates.C
+++ b/src/sampling/probes/probesTemplates.C
@@ -109,7 +109,7 @@ void Foam::probes::sampleAndWrite
forAll(values, probeI)
{
- probeStream << setw(w) << values[probeI];
+ probeStream << ' ' << setw(w) << values[probeI];
}
probeStream << endl;
}