diff --git a/applications/solvers/incompressible/pimpleFoam/pimpleDyMFoam/UEqn.H b/applications/solvers/incompressible/pimpleFoam/pimpleDyMFoam/UEqn.H
index 430e60c8482ea288a734f51f5a35dcf434dde0bf..eaaa3e60c9015afad18bdc8ff96bb3791639c951 100644
--- a/applications/solvers/incompressible/pimpleFoam/pimpleDyMFoam/UEqn.H
+++ b/applications/solvers/incompressible/pimpleFoam/pimpleDyMFoam/UEqn.H
@@ -5,7 +5,7 @@ tmp<fvVectorMatrix> UEqn
fvm::ddt(U)
+ fvm::div(phi, U)
+ turbulence->divDevReff(U)
- ==
+ ==
fvOptions(U)
);
@@ -17,7 +17,7 @@ rAU = 1.0/UEqn().A();
if (pimple.momentumPredictor())
{
- solve(UEqn() == -fvc::grad(p) + fvOptions(U));
+ solve(UEqn() == -fvc::grad(p));
fvOptions.correct(U);
}
diff --git a/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C
index 2406cb667c45098ebd0ae54dfbf61f325e48d68c..14bf0ee9eb9d5eec0442060b3ca8a308f2bd1e74 100644
--- a/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C
+++ b/applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -47,6 +47,8 @@ Description
#include "snapParameters.H"
#include "layerParameters.H"
#include "vtkSetWriter.H"
+#include "faceSet.H"
+#include "motionSmoother.H"
using namespace Foam;
@@ -675,6 +677,31 @@ int main(int argc, char *argv[])
}
+ {
+ // Check final mesh
+ Info<< "Checking final mesh ..." << endl;
+ faceSet wrongFaces(mesh, "wrongFaces", mesh.nFaces()/100);
+ motionSmoother::checkMesh(false, mesh, motionDict, wrongFaces);
+ const label nErrors = returnReduce
+ (
+ wrongFaces.size(),
+ sumOp<label>()
+ );
+
+ if (nErrors > 0)
+ {
+ Info<< "Finished meshing with " << nErrors << " illegal faces"
+ << " (concave, zero area or negative cell pyramid volume)"
+ << endl;
+ wrongFaces.write();
+ }
+ else
+ {
+ Info<< "Finished meshing without any errors" << endl;
+ }
+ }
+
+
Info<< "Finished meshing in = "
<< runTime.elapsedCpuTime() << " s." << endl;
diff --git a/applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C b/applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C
index ef521a257d2238f164acab87bf6e4307d8c8218a..cf503b057483308a9901df7500dce05572a5d308 100644
--- a/applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C
+++ b/applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C
@@ -832,5 +832,39 @@ Foam::label Foam::checkGeometry(const polyMesh& mesh, const bool allGeometry)
}
}
+ if (allGeometry)
+ {
+ faceSet faces(mesh, "lowWeightFaces", mesh.nFaces()/100);
+ if (mesh.checkFaceWeight(true, 0.05, &faces))
+ {
+ noFailedChecks++;
+
+ label nFaces = returnReduce(faces.size(), sumOp<label>());
+
+ Info<< " <<Writing " << nFaces
+ << " faces with low interpolation weights to set "
+ << faces.name() << endl;
+ faces.instance() = mesh.pointsInstance();
+ faces.write();
+ }
+ }
+
+ if (allGeometry)
+ {
+ faceSet faces(mesh, "lowVolRatioFaces", mesh.nFaces()/100);
+ if (mesh.checkVolRatio(true, 0.05, &faces))
+ {
+ noFailedChecks++;
+
+ label nFaces = returnReduce(faces.size(), sumOp<label>());
+
+ Info<< " <<Writing " << nFaces
+ << " faces with low volume ratio cells to set "
+ << faces.name() << endl;
+ faces.instance() = mesh.pointsInstance();
+ faces.write();
+ }
+ }
+
return noFailedChecks;
}
diff --git a/applications/utilities/mesh/manipulation/createBaffles/createBafflesDict b/applications/utilities/mesh/manipulation/createBaffles/createBafflesDict
index 99c03b4cfc7429af81364a42f2613eec1ee21bb9..b392f0764278002f565c875be900379d73f6e039 100644
--- a/applications/utilities/mesh/manipulation/createBaffles/createBafflesDict
+++ b/applications/utilities/mesh/manipulation/createBaffles/createBafflesDict
@@ -138,6 +138,8 @@ baffles
//- Select faces and orientation through a searchableSurface
type searchableSurface;
surface searchablePlate;
+ //name sphere.stl; // name if surface=triSurfaceMesh
+
origin (0.099 -0.006 0.004);
span (0 0.012 0.012);
diff --git a/applications/utilities/mesh/manipulation/createBaffles/faceSelection/searchableSurfaceSelection.C b/applications/utilities/mesh/manipulation/createBaffles/faceSelection/searchableSurfaceSelection.C
index 04dbe07929ed96a0474109c32fdd149baf0a88ee..e7439d52f68ed987b84f1bddbb799ca3706eb58c 100644
--- a/applications/utilities/mesh/manipulation/createBaffles/faceSelection/searchableSurfaceSelection.C
+++ b/applications/utilities/mesh/manipulation/createBaffles/faceSelection/searchableSurfaceSelection.C
@@ -28,6 +28,7 @@ License
#include "syncTools.H"
#include "searchableSurface.H"
#include "fvMesh.H"
+#include "Time.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -61,7 +62,15 @@ Foam::faceSelections::searchableSurfaceSelection::searchableSurfaceSelection
searchableSurface::New
(
word(dict.lookup("surface")),
- mesh.objectRegistry::db(),
+ IOobject
+ (
+ dict.lookupOrDefault("name", mesh.objectRegistry::db().name()),
+ mesh.time().constant(),
+ "triSurface",
+ mesh.objectRegistry::db(),
+ IOobject::MUST_READ,
+ IOobject::NO_WRITE
+ ),
dict
)
)
diff --git a/applications/utilities/mesh/manipulation/topoSet/topoSetDict b/applications/utilities/mesh/manipulation/topoSet/topoSetDict
index 57514a205352c76e8b6f0e167a8a8739c436c808..11dccf942b0f35486f18d55f8fca470b437c465f 100644
--- a/applications/utilities/mesh/manipulation/topoSet/topoSetDict
+++ b/applications/utilities/mesh/manipulation/topoSet/topoSetDict
@@ -378,6 +378,7 @@ FoamFile
// surface searchableSphere;
// centre (0.05 0.05 0.005);
// radius 0.025;
+// //name sphere.stl; // Optional name if surface triSurfaceMesh
// }
// }
//
diff --git a/applications/utilities/postProcessing/graphics/PV398Readers/Allwmake b/applications/utilities/postProcessing/graphics/PV398Readers/Allwmake
index 33d6bd18dd866fc70eae29df18d7075cdc307206..1469999cb6fcaac52d16d5c5a0be15c74fb6a483 100755
--- a/applications/utilities/postProcessing/graphics/PV398Readers/Allwmake
+++ b/applications/utilities/postProcessing/graphics/PV398Readers/Allwmake
@@ -2,21 +2,26 @@
cd ${0%/*} || exit 1 # run from this directory
#set -x
-if [ -d "$ParaView_DIR" -a -r "$ParaView_DIR" ]
-then
- [ -n "$PV_PLUGIN_PATH" ] || {
- echo "$0 : PV_PLUGIN_PATH not valid - it is unset"
- exit 1
- }
+if [ "$ParaView_VERSION" == "3.98.1" ]
+ then
+ if [ -d "$ParaView_DIR" -a -r "$ParaView_DIR" ]
+ then
+ [ -n "$PV_PLUGIN_PATH" ] || {
+ echo "$0 : PV_PLUGIN_PATH not valid - it is unset"
+ exit 1
+ }
- # ensure CMake gets the correct C++ compiler
- [ -n "$WM_CXX" ] && export CXX="$WM_CXX"
+ # ensure CMake gets the correct C++ compiler
+ [ -n "$WM_CXX" ] && export CXX="$WM_CXX"
- wmake libso vtkPV398Readers
- PV398blockMeshReader/Allwmake
- PV398FoamReader/Allwmake
+ wmake libso vtkPV398Readers
+ PV398blockMeshReader/Allwmake
+ PV398FoamReader/Allwmake
+ else
+ echo "ERROR: ParaView not found in $ParaView_DIR"
+ fi
else
- echo "ERROR: ParaView not found in $ParaView_DIR"
+ echo "WARN: PV398 readers not building: ParaView_VERSION=$ParaView_VERSION"
fi
# ----------------------------------------------------------------- end-of-file
diff --git a/applications/utilities/postProcessing/graphics/PV3Readers/Allwmake b/applications/utilities/postProcessing/graphics/PV3Readers/Allwmake
index 1568f8afb808e9df76a3c9e244ccce0cb0aba9f5..03b352116e7dd678042468da06d6640c637f6cfe 100755
--- a/applications/utilities/postProcessing/graphics/PV3Readers/Allwmake
+++ b/applications/utilities/postProcessing/graphics/PV3Readers/Allwmake
@@ -2,21 +2,26 @@
cd ${0%/*} || exit 1 # run from this directory
#set -x
-if [ -d "$ParaView_DIR" -a -r "$ParaView_DIR" ]
-then
- [ -n "$PV_PLUGIN_PATH" ] || {
- echo "$0 : PV_PLUGIN_PATH not valid - it is unset"
- exit 1
- }
+if [ "$ParaView_VERSION" != "3.98.1" ]
+ then
+ if [ -d "$ParaView_DIR" -a -r "$ParaView_DIR" ]
+ then
+ [ -n "$PV_PLUGIN_PATH" ] || {
+ echo "$0 : PV_PLUGIN_PATH not valid - it is unset"
+ exit 1
+ }
- # ensure CMake gets the correct C++ compiler
- [ -n "$WM_CXX" ] && export CXX="$WM_CXX"
+ # ensure CMake gets the correct C++ compiler
+ [ -n "$WM_CXX" ] && export CXX="$WM_CXX"
- wmake libso vtkPV3Readers
- PV3blockMeshReader/Allwmake
- PV3FoamReader/Allwmake
+ wmake libso vtkPV3Readers
+ PV3blockMeshReader/Allwmake
+ PV3FoamReader/Allwmake
+ else
+ echo "ERROR: ParaView not found in $ParaView_DIR"
+ fi
else
- echo "ERROR: ParaView not found in $ParaView_DIR"
+ echo "WARN: PV3 readers not building: ParaView_VERSION=$ParaView_VERSION"
fi
# ----------------------------------------------------------------- end-of-file
diff --git a/applications/utilities/postProcessing/graphics/PV3Readers/PV3FoamReader/vtkPV3Foam/vtkPV3FoamMesh.C b/applications/utilities/postProcessing/graphics/PV3Readers/PV3FoamReader/vtkPV3Foam/vtkPV3FoamMesh.C
index b80acb4949984b50d1145d06d96c20d65f4eb012..bbb7072b22f24690bfe262350b22852c791b2cf0 100644
--- a/applications/utilities/postProcessing/graphics/PV3Readers/PV3FoamReader/vtkPV3Foam/vtkPV3FoamMesh.C
+++ b/applications/utilities/postProcessing/graphics/PV3Readers/PV3FoamReader/vtkPV3Foam/vtkPV3FoamMesh.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -179,7 +179,10 @@ void Foam::vtkPV3Foam::convertMeshPatches
const word patchName = getPartName(partId);
- labelHashSet patchIds(patches.patchSet(List<wordRe>(1, patchName)));
+ labelHashSet patchIds
+ (
+ patches.patchSet(List<wordRe>(1, wordRe(patchName)))
+ );
if (debug)
{
diff --git a/src/OpenFOAM/interpolations/interpolationTable/tableReaders/csv/csvTableReader.C b/src/OpenFOAM/interpolations/interpolationTable/tableReaders/csv/csvTableReader.C
index b71e274842ad036d6d84a72fd685fc80ba4c965d..09faf41fd7b5e1f4349bb7b9a6cc7cf5858c22c9 100644
--- a/src/OpenFOAM/interpolations/interpolationTable/tableReaders/csv/csvTableReader.C
+++ b/src/OpenFOAM/interpolations/interpolationTable/tableReaders/csv/csvTableReader.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -182,8 +182,17 @@ void Foam::csvTableReader<Type>::write(Ostream& os) const
<< headerLine_ << token::END_STATEMENT << nl;
os.writeKeyword("timeColumn")
<< timeColumn_ << token::END_STATEMENT << nl;
- os.writeKeyword("valueColumns")
- << componentColumns_ << token::END_STATEMENT << nl;
+
+ // Force writing labelList in ascii
+ os.writeKeyword("valueColumns");
+ if (os.format() == IOstream::BINARY)
+ {
+ os.format(IOstream::ASCII);
+ os << componentColumns_;
+ os.format(IOstream::BINARY);
+ }
+ os << token::END_STATEMENT << nl;
+
os.writeKeyword("separator")
<< string(separator_) << token::END_STATEMENT << nl;
}
diff --git a/src/OpenFOAM/meshes/polyMesh/polyMesh.H b/src/OpenFOAM/meshes/polyMesh/polyMesh.H
index 76c5f1c8a3777cc65b520d2a720fc1fc66ddd895..a416fafc301d52778034467bdf9f83e649956d08 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyMesh.H
+++ b/src/OpenFOAM/meshes/polyMesh/polyMesh.H
@@ -265,6 +265,24 @@ private:
const Vector<label>& meshD
) const;
+ bool checkFaceWeight
+ (
+ const vectorField& fCtrs,
+ const vectorField& fAreas,
+ const vectorField& cellCtrs,
+ const bool report,
+ const scalar minWeight,
+ labelHashSet* setPtr
+ ) const;
+
+ bool checkVolRatio
+ (
+ const scalarField& cellVols,
+ const bool report,
+ const scalar minRatio,
+ labelHashSet* setPtr
+ ) const;
+
public:
// Public typedefs
@@ -612,6 +630,22 @@ public:
const bool detailedReport = false
) const;
+ //- Check for face weights
+ virtual bool checkFaceWeight
+ (
+ const bool report,
+ const scalar minWeight = 0.05,
+ labelHashSet* setPtr = NULL
+ ) const;
+
+ //- Check for neighbouring cell volumes
+ virtual bool checkVolRatio
+ (
+ const bool report,
+ const scalar minRatio = 0.01,
+ labelHashSet* setPtr = NULL
+ ) const;
+
// Helper functions
diff --git a/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshCheck.C b/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshCheck.C
index f7daba04955b96710f570267d344519909a3d9f9..2b241bf6cddc84d2420a2cfe0361abee466dc4aa 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshCheck.C
+++ b/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshCheck.C
@@ -142,7 +142,8 @@ bool Foam::polyMesh::checkFaceOrthogonality
if (severeNonOrth > 0)
{
- Info<< " *Number of severely non-orthogonal faces: "
+ Info<< " *Number of severely non-orthogonal (> "
+ << primitiveMesh::nonOrthThreshold_ << " degrees) faces: "
<< severeNonOrth << "." << endl;
}
}
@@ -427,6 +428,8 @@ bool Foam::polyMesh::checkCellDeterminant
const Vector<label>& meshD
) const
{
+ const scalar warnDet = 1e-3;
+
if (debug)
{
Info<< "bool polyMesh::checkCellDeterminant(const bool"
@@ -450,7 +453,7 @@ bool Foam::polyMesh::checkCellDeterminant
forAll (cellDeterminant, cellI)
{
- if (cellDeterminant[cellI] < 1e-3)
+ if (cellDeterminant[cellI] < warnDet)
{
if (setPtr)
{
@@ -480,7 +483,8 @@ bool Foam::polyMesh::checkCellDeterminant
{
if (debug || report)
{
- Info<< " ***Cells with small determinant found, number of cells: "
+ Info<< " ***Cells with small determinant (< "
+ << warnDet << ") found, number of cells: "
<< nErrorCells << endl;
}
@@ -500,6 +504,192 @@ bool Foam::polyMesh::checkCellDeterminant
}
+bool Foam::polyMesh::checkFaceWeight
+(
+ const vectorField& fCtrs,
+ const vectorField& fAreas,
+ const vectorField& cellCtrs,
+ const bool report,
+ const scalar minWeight,
+ labelHashSet* setPtr
+) const
+{
+ if (debug)
+ {
+ Info<< "bool polyMesh::checkFaceWeight(const bool"
+ << ", labelHashSet*) const: "
+ << "checking for low face interpolation weights" << endl;
+ }
+
+ tmp<scalarField> tfaceWght = polyMeshTools::faceWeights
+ (
+ *this,
+ fCtrs,
+ fAreas,
+ cellCtrs
+ );
+ scalarField& faceWght = tfaceWght();
+
+
+ label nErrorFaces = 0;
+ scalar minDet = GREAT;
+ scalar sumDet = 0.0;
+ label nSummed = 0;
+
+ // Statistics only for internal and masters of coupled faces
+ PackedBoolList isMasterFace(syncTools::getInternalOrMasterFaces(*this));
+
+ forAll (faceWght, faceI)
+ {
+ if (faceWght[faceI] < minWeight)
+ {
+ // Note: insert both sides of coupled faces
+ if (setPtr)
+ {
+ setPtr->insert(faceI);
+ }
+
+ nErrorFaces++;
+ }
+
+ // Note: statistics only on master of coupled faces
+ if (isMasterFace[faceI])
+ {
+ minDet = min(minDet, faceWght[faceI]);
+ sumDet += faceWght[faceI];
+ nSummed++;
+ }
+ }
+
+ reduce(nErrorFaces, sumOp<label>());
+ reduce(minDet, minOp<scalar>());
+ reduce(sumDet, sumOp<scalar>());
+ reduce(nSummed, sumOp<label>());
+
+ if (debug || report)
+ {
+ if (nSummed > 0)
+ {
+ Info<< " Face interpolation weight : minimum: " << minDet
+ << " average: " << sumDet/nSummed
+ << endl;
+ }
+ }
+
+ if (nErrorFaces > 0)
+ {
+ if (debug || report)
+ {
+ Info<< " ***Faces with small interpolation weight (< " << minWeight
+ << ") found, number of faces: "
+ << nErrorFaces << endl;
+ }
+
+ return true;
+ }
+ else
+ {
+ if (debug || report)
+ {
+ Info<< " Face interpolation weight check OK." << endl;
+ }
+
+ return false;
+ }
+
+ return false;
+}
+
+
+bool Foam::polyMesh::checkVolRatio
+(
+ const scalarField& cellVols,
+ const bool report,
+ const scalar minRatio,
+ labelHashSet* setPtr
+) const
+{
+ if (debug)
+ {
+ Info<< "bool polyMesh::checkVolRatio(const bool"
+ << ", labelHashSet*) const: "
+ << "checking for volume ratio < " << minRatio << endl;
+ }
+
+ tmp<scalarField> tvolRatio = polyMeshTools::volRatio(*this, cellVols);
+ scalarField& volRatio = tvolRatio();
+
+
+ label nErrorFaces = 0;
+ scalar minDet = GREAT;
+ scalar sumDet = 0.0;
+ label nSummed = 0;
+
+ // Statistics only for internal and masters of coupled faces
+ PackedBoolList isMasterFace(syncTools::getInternalOrMasterFaces(*this));
+
+ forAll (volRatio, faceI)
+ {
+ if (volRatio[faceI] < minRatio)
+ {
+ // Note: insert both sides of coupled faces
+ if (setPtr)
+ {
+ setPtr->insert(faceI);
+ }
+
+ nErrorFaces++;
+ }
+
+ // Note: statistics only on master of coupled faces
+ if (isMasterFace[faceI])
+ {
+ minDet = min(minDet, volRatio[faceI]);
+ sumDet += volRatio[faceI];
+ nSummed++;
+ }
+ }
+
+ reduce(nErrorFaces, sumOp<label>());
+ reduce(minDet, minOp<scalar>());
+ reduce(sumDet, sumOp<scalar>());
+ reduce(nSummed, sumOp<label>());
+
+ if (debug || report)
+ {
+ if (nSummed > 0)
+ {
+ Info<< " Face volume ratio : minimum: " << minDet
+ << " average: " << sumDet/nSummed
+ << endl;
+ }
+ }
+
+ if (nErrorFaces > 0)
+ {
+ if (debug || report)
+ {
+ Info<< " ***Faces with small volume ratio (< " << minRatio
+ << ") found, number of faces: "
+ << nErrorFaces << endl;
+ }
+
+ return true;
+ }
+ else
+ {
+ if (debug || report)
+ {
+ Info<< " Face volume ratio check OK." << endl;
+ }
+
+ return false;
+ }
+
+ return false;
+}
+
+
//- Could override checkClosedBoundary to not look at (collocated!) coupled
// faces
//bool Foam::polyMesh::checkClosedBoundary(const bool report) const
@@ -582,6 +772,36 @@ bool Foam::polyMesh::checkCellDeterminant
}
+bool Foam::polyMesh::checkFaceWeight
+(
+ const bool report,
+ const scalar minWeight,
+ labelHashSet* setPtr
+) const
+{
+ return checkFaceWeight
+ (
+ faceCentres(),
+ faceAreas(),
+ cellCentres(),
+ report,
+ minWeight,
+ setPtr
+ );
+}
+
+
+bool Foam::polyMesh::checkVolRatio
+(
+ const bool report,
+ const scalar minRatio,
+ labelHashSet* setPtr
+) const
+{
+ return checkVolRatio(cellVolumes(), report, minRatio, setPtr);
+}
+
+
bool Foam::polyMesh::checkMeshMotion
(
const pointField& newPoints,
diff --git a/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshTools.C b/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshTools.C
index 55f6a907676e7f241fd34c11befaa23a5de4332f..59947b0fbe2528a6045ce36daf90570f8cdbbd55 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshTools.C
+++ b/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshTools.C
@@ -187,4 +187,112 @@ Foam::tmp<Foam::scalarField> Foam::polyMeshTools::faceSkewness
}
+Foam::tmp<Foam::scalarField> Foam::polyMeshTools::faceWeights
+(
+ const polyMesh& mesh,
+ const vectorField& fCtrs,
+ const vectorField& fAreas,
+ const vectorField& cellCtrs
+)
+{
+ const labelList& own = mesh.faceOwner();
+ const labelList& nei = mesh.faceNeighbour();
+ const polyBoundaryMesh& pbm = mesh.boundaryMesh();
+
+ tmp<scalarField> tweight(new scalarField(mesh.nFaces(), 1.0));
+ scalarField& weight = tweight();
+
+ // Internal faces
+ forAll(nei, faceI)
+ {
+ const point& fc = fCtrs[faceI];
+ const vector& fa = fAreas[faceI];
+
+ scalar dOwn = mag(fa & (fc-cellCtrs[own[faceI]]));
+ scalar dNei = mag(fa & (cellCtrs[nei[faceI]]-fc));
+
+ weight[faceI] = min(dNei,dOwn)/(dNei+dOwn+VSMALL);
+ }
+
+
+ // Coupled faces
+
+ pointField neiCc;
+ syncTools::swapBoundaryCellPositions(mesh, cellCtrs, neiCc);
+
+ forAll(pbm, patchI)
+ {
+ const polyPatch& pp = pbm[patchI];
+ if (pp.coupled())
+ {
+ forAll(pp, i)
+ {
+ label faceI = pp.start() + i;
+ label bFaceI = faceI - mesh.nInternalFaces();
+
+ const point& fc = fCtrs[faceI];
+ const vector& fa = fAreas[faceI];
+
+ scalar dOwn = mag(fa & (fc-cellCtrs[own[faceI]]));
+ scalar dNei = mag(fa & (neiCc[bFaceI]-fc));
+
+ weight[faceI] = min(dNei,dOwn)/(dNei+dOwn+VSMALL);
+ }
+ }
+ }
+
+ return tweight;
+}
+
+
+Foam::tmp<Foam::scalarField> Foam::polyMeshTools::volRatio
+(
+ const polyMesh& mesh,
+ const scalarField& vol
+)
+{
+ const labelList& own = mesh.faceOwner();
+ const labelList& nei = mesh.faceNeighbour();
+ const polyBoundaryMesh& pbm = mesh.boundaryMesh();
+
+ tmp<scalarField> tratio(new scalarField(mesh.nFaces(), 1.0));
+ scalarField& ratio = tratio();
+
+ // Internal faces
+ forAll(nei, faceI)
+ {
+ scalar volOwn = vol[own[faceI]];
+ scalar volNei = vol[nei[faceI]];
+
+ ratio[faceI] = min(volOwn,volNei)/(max(volOwn, volNei)+VSMALL);
+ }
+
+
+ // Coupled faces
+
+ scalarField neiVol;
+ syncTools::swapBoundaryCellList(mesh, vol, neiVol);
+
+ forAll(pbm, patchI)
+ {
+ const polyPatch& pp = pbm[patchI];
+ if (pp.coupled())
+ {
+ forAll(pp, i)
+ {
+ label faceI = pp.start() + i;
+ label bFaceI = faceI - mesh.nInternalFaces();
+
+ scalar volOwn = vol[own[faceI]];
+ scalar volNei = neiVol[bFaceI];
+
+ ratio[faceI] = min(volOwn,volNei)/(max(volOwn, volNei)+VSMALL);
+ }
+ }
+ }
+
+ return tratio;
+}
+
+
// ************************************************************************* //
diff --git a/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshTools.H b/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshTools.H
index 53c758c3be7967f5a914b39e7e0c12a50bbaf86c..7f23facd04ec02eb0b43abb1b64ee2702f5be88e 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshTools.H
+++ b/src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshTools.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2012-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -95,6 +95,22 @@ public:
// );
// }
+ //- Generate interpolation factors field
+ static tmp<scalarField> faceWeights
+ (
+ const polyMesh& mesh,
+ const vectorField& fCtrs,
+ const vectorField& fAreas,
+ const vectorField& cellCtrs
+ );
+
+ //- Generate volume ratio field
+ static tmp<scalarField> volRatio
+ (
+ const polyMesh& mesh,
+ const scalarField& vol
+ );
+
};
diff --git a/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/processorCyclic/processorCyclicPolyPatch.H b/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/processorCyclic/processorCyclicPolyPatch.H
index 0b9d4f70d94439fcf120b74ea404f75105a5aa53..76a833ba612e9adb0c77353c60913195e8fc853d 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/processorCyclic/processorCyclicPolyPatch.H
+++ b/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/processorCyclic/processorCyclicPolyPatch.H
@@ -45,7 +45,7 @@ namespace Foam
{
/*---------------------------------------------------------------------------*\
- Class processorCyclicPolyPatch Declaration
+ Class processorCyclicPolyPatch Declaration
\*---------------------------------------------------------------------------*/
class processorCyclicPolyPatch
@@ -63,8 +63,6 @@ class processorCyclicPolyPatch
//- Index of originating patch
mutable label referPatchID_;
- // Private member functions
-
protected:
@@ -223,8 +221,7 @@ public:
// Destructor
-
- virtual ~processorCyclicPolyPatch();
+ virtual ~processorCyclicPolyPatch();
// Member functions
diff --git a/src/OpenFOAM/meshes/primitiveShapes/triangle/triangle.H b/src/OpenFOAM/meshes/primitiveShapes/triangle/triangle.H
index 39a61c834b90d9ebefe6388a63283bca55ce739e..de339b172e0d00881d2c14e893a021bc1e26e9e3 100644
--- a/src/OpenFOAM/meshes/primitiveShapes/triangle/triangle.H
+++ b/src/OpenFOAM/meshes/primitiveShapes/triangle/triangle.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -42,6 +42,7 @@ SourceFiles
#include "Random.H"
#include "FixedList.H"
#include "UList.H"
+#include "linePointRef.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -222,6 +223,15 @@ public:
label& nearLabel
) const;
+ //- Return nearest point to line on triangle. Returns hit if
+ // point is inside triangle. Sets edgePoint to point on edge
+ // (hit if nearest is inside line)
+ inline pointHit nearestPoint
+ (
+ const linePointRef& edge,
+ pointHit& edgePoint
+ ) const;
+
// IOstream operators
diff --git a/src/OpenFOAM/meshes/primitiveShapes/triangle/triangleI.H b/src/OpenFOAM/meshes/primitiveShapes/triangle/triangleI.H
index 496442012d0e92c0b1b65ff0842d36357269b5dd..2c945fae5450a1b67e859a6481b4e9ed6407ff98 100644
--- a/src/OpenFOAM/meshes/primitiveShapes/triangle/triangleI.H
+++ b/src/OpenFOAM/meshes/primitiveShapes/triangle/triangleI.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -659,6 +659,132 @@ inline bool Foam::triangle<Point, PointRef>::classify
}
+template<class Point, class PointRef>
+inline Foam::pointHit Foam::triangle<Point, PointRef>::nearestPoint
+(
+ const linePointRef& ln,
+ pointHit& lnInfo
+) const
+{
+ vector q = ln.vec();
+ pointHit triInfo
+ (
+ triangle<Point, PointRef>::intersection
+ (
+ ln.start(),
+ q,
+ intersection::FULL_RAY
+ )
+ );
+
+ if (triInfo.hit())
+ {
+ // Line hits triangle. Find point on line.
+ if (triInfo.distance() > 1)
+ {
+ // Hit beyond endpoint
+ lnInfo.setMiss(true);
+ lnInfo.setPoint(ln.end());
+ scalar dist = Foam::mag(triInfo.hitPoint()-lnInfo.missPoint());
+ lnInfo.setDistance(dist);
+ triInfo.setMiss(true);
+ triInfo.setDistance(dist);
+ }
+ else if (triInfo.distance() < 0)
+ {
+ // Hit beyond startpoint
+ lnInfo.setMiss(true);
+ lnInfo.setPoint(ln.start());
+ scalar dist = Foam::mag(triInfo.hitPoint()-lnInfo.missPoint());
+ lnInfo.setDistance(dist);
+ triInfo.setMiss(true);
+ triInfo.setDistance(dist);
+ }
+ else
+ {
+ // Hit on line
+ lnInfo.setHit();
+ lnInfo.setPoint(triInfo.hitPoint());
+ lnInfo.setDistance(0.0);
+ triInfo.setDistance(0.0);
+ }
+ }
+ else
+ {
+ // Line skips triangle. See which triangle edge it gets closest to
+
+ point nearestEdgePoint;
+ point nearestLinePoint;
+ label minEdgeIndex = 0;
+ scalar minDist = ln.nearestDist
+ (
+ linePointRef(a_, b_),
+ nearestLinePoint,
+ nearestEdgePoint
+ );
+
+ {
+ point linePoint;
+ point triEdgePoint;
+ scalar dist = ln.nearestDist
+ (
+ linePointRef(b_, c_),
+ linePoint,
+ triEdgePoint
+ );
+ if (dist < minDist)
+ {
+ minDist = dist;
+ nearestEdgePoint = triEdgePoint;
+ nearestLinePoint = linePoint;
+ minEdgeIndex = 1;
+ }
+ }
+
+ {
+ point linePoint;
+ point triEdgePoint;
+ scalar dist = ln.nearestDist
+ (
+ linePointRef(c_, a_),
+ linePoint,
+ triEdgePoint
+ );
+ if (dist < minDist)
+ {
+ minDist = dist;
+ nearestEdgePoint = triEdgePoint;
+ nearestLinePoint = linePoint;
+ minEdgeIndex = 2;
+ }
+ }
+
+ lnInfo.setDistance(minDist);
+ triInfo.setDistance(minDist);
+ triInfo.setMiss(false);
+ triInfo.setPoint(nearestEdgePoint);
+
+ // Convert point on line to pointHit
+ if (Foam::mag(nearestLinePoint-ln.start()) < SMALL)
+ {
+ lnInfo.setMiss(true);
+ lnInfo.setPoint(ln.start());
+ }
+ else if (Foam::mag(nearestLinePoint-ln.end()) < SMALL)
+ {
+ lnInfo.setMiss(true);
+ lnInfo.setPoint(ln.end());
+ }
+ else
+ {
+ lnInfo.setHit();
+ lnInfo.setPoint(nearestLinePoint);
+ }
+ }
+ return triInfo;
+}
+
+
// * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * //
template<class Point, class PointRef>
diff --git a/src/OpenFOAM/primitives/functions/DataEntry/CSV/CSVIO.C b/src/OpenFOAM/primitives/functions/DataEntry/CSV/CSVIO.C
index 36595e7dd6d625092c5b5ba3cf63ff3fbc2cfd54..1dd395cee0d227186830f57016d4e9997e22aecf 100644
--- a/src/OpenFOAM/primitives/functions/DataEntry/CSV/CSVIO.C
+++ b/src/OpenFOAM/primitives/functions/DataEntry/CSV/CSVIO.C
@@ -74,8 +74,17 @@ void Foam::CSV<Type>::writeData(Ostream& os) const
os.writeKeyword("nHeaderLine") << nHeaderLine_ << token::END_STATEMENT
<< nl;
os.writeKeyword("refColumn") << refColumn_ << token::END_STATEMENT << nl;
- os.writeKeyword("componentColumns") << componentColumns_
- << token::END_STATEMENT << nl;
+
+ // Force writing labelList in ascii
+ os.writeKeyword("componentColumns");
+ if (os.format() == IOstream::BINARY)
+ {
+ os.format(IOstream::ASCII);
+ os << componentColumns_;
+ os.format(IOstream::BINARY);
+ }
+ os << token::END_STATEMENT << nl;
+
os.writeKeyword("separator") << string(separator_)
<< token::END_STATEMENT << nl;
os.writeKeyword("mergeSeparators") << string(mergeSeparators_)
diff --git a/src/OpenFOAM/primitives/strings/wordRe/wordRe.H b/src/OpenFOAM/primitives/strings/wordRe/wordRe.H
index e45044c9975b3cc951a0d65e52add51492f5cdc0..b0cb62e0b43f490c8c62befe5e2d618dfec493ce 100644
--- a/src/OpenFOAM/primitives/strings/wordRe/wordRe.H
+++ b/src/OpenFOAM/primitives/strings/wordRe/wordRe.H
@@ -120,22 +120,25 @@ public:
inline wordRe(const wordRe&);
//- Construct from keyType
- inline wordRe(const keyType&, const compOption=LITERAL);
+ inline explicit wordRe(const keyType&);
+
+ //- Construct from keyType
+ inline wordRe(const keyType&, const compOption);
//- Construct as copy of word
inline explicit wordRe(const word&);
//- Construct as copy of character array
// Optionally specify how it should be treated.
- inline wordRe(const char*, const compOption = LITERAL);
+ inline explicit wordRe(const char*, const compOption = LITERAL);
//- Construct as copy of string.
// Optionally specify how it should be treated.
- inline wordRe(const string&, const compOption = LITERAL);
+ inline explicit wordRe(const string&, const compOption = LITERAL);
//- Construct as copy of std::string
// Optionally specify how it should be treated.
- inline wordRe(const std::string&, const compOption = LITERAL);
+ inline explicit wordRe(const std::string&, const compOption = LITERAL);
//- Construct from Istream
// Words are treated as literals, strings with an auto-test
diff --git a/src/OpenFOAM/primitives/strings/wordRe/wordReI.H b/src/OpenFOAM/primitives/strings/wordRe/wordReI.H
index 431106716544710f262d89bc37946a8372c2582a..01a3114c6c9b2c6da7c389b501831eaaa7dca735 100644
--- a/src/OpenFOAM/primitives/strings/wordRe/wordReI.H
+++ b/src/OpenFOAM/primitives/strings/wordRe/wordReI.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -65,6 +65,18 @@ inline Foam::wordRe::wordRe(const word& str)
{}
+inline Foam::wordRe::wordRe(const keyType& str)
+:
+ word(str, false),
+ re_()
+{
+ if (str.isPattern())
+ {
+ compile();
+ }
+}
+
+
inline Foam::wordRe::wordRe(const keyType& str, const compOption opt)
:
word(str, false),
diff --git a/src/dynamicMesh/motionSmoother/motionSmoother.C b/src/dynamicMesh/motionSmoother/motionSmoother.C
index c7d5d28e765ada2e25f0562519b16b5b5c7c4d7d..dad3ea59774e79359825640844dcd9ce9d83a3b1 100644
--- a/src/dynamicMesh/motionSmoother/motionSmoother.C
+++ b/src/dynamicMesh/motionSmoother/motionSmoother.C
@@ -678,34 +678,8 @@ void Foam::motionSmoother::correct()
}
-void Foam::motionSmoother::setDisplacement(pointField& patchDisp)
+void Foam::motionSmoother::setDisplacementPatchFields()
{
- // See comment in .H file about shared points.
- const polyBoundaryMesh& patches = mesh_.boundaryMesh();
-
- forAll(patches, patchI)
- {
- const polyPatch& pp = patches[patchI];
-
- if (pp.coupled())
- {
- const labelList& meshPoints = pp.meshPoints();
-
- forAll(meshPoints, i)
- {
- displacement_[meshPoints[i]] = vector::zero;
- }
- }
- }
-
- const labelList& ppMeshPoints = pp_.meshPoints();
-
- // Set internal point data from displacement on combined patch points.
- forAll(ppMeshPoints, patchPointI)
- {
- displacement_[ppMeshPoints[patchPointI]] = patchDisp[patchPointI];
- }
-
// Adapt the fixedValue bc's (i.e. copy internal point data to
// boundaryField for all affected patches)
forAll(adaptPatchIDs_, i)
@@ -765,6 +739,42 @@ void Foam::motionSmoother::setDisplacement(pointField& patchDisp)
displacement_.boundaryField()[patchI] ==
displacement_.boundaryField()[patchI].patchInternalField();
}
+}
+
+
+void Foam::motionSmoother::setDisplacement(pointField& patchDisp)
+{
+ // See comment in .H file about shared points.
+ const polyBoundaryMesh& patches = mesh_.boundaryMesh();
+
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
+
+ if (pp.coupled())
+ {
+ const labelList& meshPoints = pp.meshPoints();
+
+ forAll(meshPoints, i)
+ {
+ displacement_[meshPoints[i]] = vector::zero;
+ }
+ }
+ }
+
+ const labelList& ppMeshPoints = pp_.meshPoints();
+
+ // Set internal point data from displacement on combined patch points.
+ forAll(ppMeshPoints, patchPointI)
+ {
+ displacement_[ppMeshPoints[patchPointI]] = patchDisp[patchPointI];
+ }
+
+
+ // Adapt the fixedValue bc's (i.e. copy internal point data to
+ // boundaryField for all affected patches)
+ setDisplacementPatchFields();
+
if (debug)
{
diff --git a/src/dynamicMesh/motionSmoother/motionSmoother.H b/src/dynamicMesh/motionSmoother/motionSmoother.H
index ee5432ad4ec4218768f8226a95770ecd71915897..fff07c4b94a8a73c1490ef084024a93bb73f8a3b 100644
--- a/src/dynamicMesh/motionSmoother/motionSmoother.H
+++ b/src/dynamicMesh/motionSmoother/motionSmoother.H
@@ -385,6 +385,10 @@ public:
//- Take over existing mesh position.
void correct();
+ //- Set patch fields on displacement to be consistent with
+ // internal values.
+ void setDisplacementPatchFields();
+
//- Set displacement field from displacement on patch points.
// Modify provided displacement to be consistent with actual
// boundary conditions on displacement. Note: resets the
diff --git a/src/finiteVolume/fields/fvPatchFields/derived/mappedField/mappedPatchFieldBase.H b/src/finiteVolume/fields/fvPatchFields/derived/mappedField/mappedPatchFieldBase.H
index 50728410949b23ae732b5110718dfc7a3b2642c3..22a3872effe9d5eb7b007c37a350c91ec9d9ad24 100644
--- a/src/finiteVolume/fields/fvPatchFields/derived/mappedField/mappedPatchFieldBase.H
+++ b/src/finiteVolume/fields/fvPatchFields/derived/mappedField/mappedPatchFieldBase.H
@@ -25,7 +25,18 @@ Class
Foam::mappedPatchFieldBase
Description
- Functionality for sampling fields using mappedPatchBase.
+ Functionality for sampling fields using mappedPatchBase. Every call to
+ mappedField() returns a sampled field, optionally scaled to maintain an
+ area-weighted average.
+
+ Example usage:
+
+ {
+ fieldName T; // default is same as fvPatchField
+ setAverage false;
+ average 1.0; // only if setAverage=true
+ interpolationScheme cellPoint; // default is cell
+ }
SourceFiles
mappedPatchFieldBase.C
diff --git a/src/finiteVolume/fields/fvPatchFields/derived/pressureInletOutletVelocity/pressureInletOutletVelocityFvPatchVectorField.H b/src/finiteVolume/fields/fvPatchFields/derived/pressureInletOutletVelocity/pressureInletOutletVelocityFvPatchVectorField.H
index 1ed301b191ea185db943a58de71215682d509370..49ebb5d77e8ca16dfd2ba119b0b49cd3cf3fb367 100644
--- a/src/finiteVolume/fields/fvPatchFields/derived/pressureInletOutletVelocity/pressureInletOutletVelocityFvPatchVectorField.H
+++ b/src/finiteVolume/fields/fvPatchFields/derived/pressureInletOutletVelocity/pressureInletOutletVelocityFvPatchVectorField.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -29,7 +29,7 @@ Group
Description
This velocity inlet/outlet boundary condition is applied to pressure
- boundaries where the pressure is specified. A zero-gradient condtion is
+ boundaries where the pressure is specified. A zero-gradient condition is
applied for outflow (as defined by the flux); for inflow, the velocity is
obtained from the patch-face normal component of the internal-cell value.
diff --git a/src/finiteVolume/finiteVolume/gradSchemes/fourthGrad/fourthGrad.C b/src/finiteVolume/finiteVolume/gradSchemes/fourthGrad/fourthGrad.C
index a5b561f8fd81947991bab22bb47a6b909fc6f3f0..2293c126db80a3f0ea39bb0ba0c3afe47dea0193 100644
--- a/src/finiteVolume/finiteVolume/gradSchemes/fourthGrad/fourthGrad.C
+++ b/src/finiteVolume/finiteVolume/gradSchemes/fourthGrad/fourthGrad.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -62,7 +62,11 @@ Foam::fv::fourthGrad<Type>::calcGrad
// Assemble the second-order least-square gradient
// Calculate the second-order least-square gradient
tmp<GeometricField<GradType, fvPatchField, volMesh> > tsecondfGrad
- = leastSquaresGrad<Type>(mesh).grad(vsf);
+ = leastSquaresGrad<Type>(mesh).grad
+ (
+ vsf,
+ "leastSquaresGrad(" + vsf.name() + ")"
+ );
const GeometricField<GradType, fvPatchField, volMesh>& secondfGrad =
tsecondfGrad();
diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C
index 461ec5b285c6e4eac516534ca1eb7180cdcf7f15..9350e1251ed9208395e9f7c5f4e83a8a55da31ca 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -477,7 +477,7 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc
const scalar T4 = pow4(this->T_);
td.cloud().radAreaP()[cellI] += dt*np0*ap;
td.cloud().radT4()[cellI] += dt*np0*T4;
- td.cloud().radAreaP()[cellI] += dt*np0*ap*T4;
+ td.cloud().radAreaPT4()[cellI] += dt*np0*ap*T4;
}
}
}
diff --git a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C
index ed60ded9596a6b0545dd8968f5f7782fd2150694..4906280cf9f8c3a3ecae5f2a79e34e7edcc4a295 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C
@@ -466,7 +466,7 @@ void Foam::ReactingParcel<ParcelType>::calc
const scalar T4 = pow4(this->T_);
td.cloud().radAreaP()[cellI] += dt*np0*ap;
td.cloud().radT4()[cellI] += dt*np0*T4;
- td.cloud().radAreaP()[cellI] += dt*np0*ap*T4;
+ td.cloud().radAreaPT4()[cellI] += dt*np0*ap*T4;
}
}
}
diff --git a/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C b/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C
index 7baa54c1eb4ce894b8a691a7b4c06958509cc3d5..5c1dd3566a6dd709efd957369c79bdaeb5e45c14 100644
--- a/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C
+++ b/src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -260,7 +260,7 @@ void Foam::ThermoParcel<ParcelType>::calc
const scalar T4 = pow4(this->T_);
td.cloud().radAreaP()[cellI] += dt*np0*ap;
td.cloud().radT4()[cellI] += dt*np0*T4;
- td.cloud().radAreaP()[cellI] += dt*np0*ap*T4;
+ td.cloud().radAreaPT4()[cellI] += dt*np0*ap*T4;
}
}
}
diff --git a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.C b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.C
index 6574ad9cbd88f9ac84769eca3eb87f2eff793750..8444b6632c8237e5ab6316ec9d36d0e5171ee3a0 100644
--- a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.C
+++ b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -138,6 +138,7 @@ Foam::layerParameters::layerParameters
readLabel(dict.lookup("nSmoothSurfaceNormals"))
),
nSmoothNormals_(readLabel(dict.lookup("nSmoothNormals"))),
+ nSmoothDisplacement_(dict.lookupOrDefault("nSmoothDisplacement", 0)),
nSmoothThickness_(readLabel(dict.lookup("nSmoothThickness"))),
maxFaceThicknessRatio_
(
@@ -278,7 +279,7 @@ Foam::layerParameters::layerParameters
const keyType& key = iter().keyword();
const labelHashSet patchIDs
(
- boundaryMesh.patchSet(List<wordRe>(1, key))
+ boundaryMesh.patchSet(List<wordRe>(1, wordRe(key)))
);
if (patchIDs.size() == 0)
diff --git a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.H b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.H
index 2578c44b687459cadb2daf8e11798fae1b069a7d..4005dd5df82e58192e04e7fc05cdc69fba0c939b 100644
--- a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.H
+++ b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -119,6 +119,8 @@ private:
label nSmoothNormals_;
+ label nSmoothDisplacement_;
+
label nSmoothThickness_;
scalar maxFaceThicknessRatio_;
@@ -275,6 +277,12 @@ public:
return layerTerminationCos_;
}
+ //- Smooth internal displacement
+ label nSmoothDisplacement() const
+ {
+ return nSmoothDisplacement_;
+ }
+
//- Smooth layer thickness over surface patches
label nSmoothThickness() const
{
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C
index f62ea404d9ae2090ea452e6f3fbadcf4d4001ba1..f37d7ae2bdc2f6240f9149b2864ccce8c6b815e7 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.C
@@ -24,128 +24,107 @@ License
\*---------------------------------------------------------------------------*/
#include "AMIInterpolation.H"
+#include "AMIMethod.H"
#include "meshTools.H"
#include "mapDistribute.H"
-// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
template<class SourcePatch, class TargetPatch>
-void Foam::AMIInterpolation<SourcePatch, TargetPatch>::writeIntersectionOBJ
+Foam::word
+Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolationMethodToWord
(
- const scalar area,
- const face& f1,
- const face& f2,
- const pointField& f1Points,
- const pointField& f2Points
-) const
+ const interpolationMethod& im
+)
{
- static label count = 1;
-
- const pointField f1pts = f1.points(f1Points);
- const pointField f2pts = f2.points(f2Points);
-
- Pout<< "Face intersection area (" << count << "):" << nl
- << " f1 face = " << f1 << nl
- << " f1 pts = " << f1pts << nl
- << " f2 face = " << f2 << nl
- << " f2 pts = " << f2pts << nl
- << " area = " << area
- << endl;
-
- OFstream os("areas" + name(count) + ".obj");
-
- forAll(f1pts, i)
- {
- meshTools::writeOBJ(os, f1pts[i]);
- }
- os<< "l";
- forAll(f1pts, i)
- {
- os<< " " << i + 1;
- }
- os<< " 1" << endl;
-
+ word method = "unknown-interpolationMethod";
- forAll(f2pts, i)
+ switch (im)
{
- meshTools::writeOBJ(os, f2pts[i]);
- }
- os<< "l";
- forAll(f2pts, i)
- {
- os<< " " << f1pts.size() + i + 1;
+ case imDirect:
+ {
+ method = "directAMI";
+ break;
+ }
+ case imMapNearest:
+ {
+ method = "mapNearestAMI";
+ break;
+ }
+ case imFaceAreaWeight:
+ {
+ method = "faceAreaWeightAMI";
+ break;
+ }
+ default:
+ {
+ FatalErrorIn
+ (
+ "const Foam::word"
+ "Foam::AMIInterpolation<SourcePatch, TargetPatch>::"
+ "interpolationMethodToWord"
+ "("
+ "const interpolationMethod&"
+ ")"
+ )
+ << "Unhandled interpolationMethod enumeration " << method
+ << abort(FatalError);
+ }
}
- os<< " " << f1pts.size() + 1 << endl;
- count++;
+ return method;
}
template<class SourcePatch, class TargetPatch>
-void Foam::AMIInterpolation<SourcePatch, TargetPatch>::checkPatches
+typename Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolationMethod
+Foam::AMIInterpolation<SourcePatch, TargetPatch>::wordTointerpolationMethod
(
- const SourcePatch& srcPatch,
- const TargetPatch& tgtPatch
-) const
+ const word& im
+)
{
- const scalar maxBoundsError = 0.05;
-
- // check bounds of source and target
- boundBox bbSrc(srcPatch.points(), srcPatch.meshPoints(), true);
- boundBox bbTgt(tgtPatch.points(), tgtPatch.meshPoints(), true);
+ interpolationMethod method = imDirect;
- boundBox bbTgtInf(bbTgt);
- bbTgtInf.inflate(maxBoundsError);
+ wordList methods
+ (
+ IStringStream("(directAMI mapNearestAMI faceAreaWeightAMI)")()
+ );
- if (!bbTgtInf.contains(bbSrc))
+ if (im == "directAMI")
+ {
+ method = imDirect;
+ }
+ else if (im == "mapNearestAMI")
{
- WarningIn
+ method = imMapNearest;
+ }
+ else if (im == "faceAreaWeightAMI")
+ {
+ method = imFaceAreaWeight;
+ }
+ else
+ {
+ FatalErrorIn
(
- "AMIInterpolation<SourcePatch, TargetPatch>::checkPatches"
+ "Foam::AMIInterpolation<SourcePatch, TargetPatch>::"
+ "interpolationMethod"
+ "Foam::AMIInterpolation<SourcePatch, TargetPatch>::"
+ "wordTointerpolationMethod"
"("
- "const SourcePatch&, "
- "const TargetPatch&"
+ "const word&"
")"
- ) << "Source and target patch bounding boxes are not similar" << nl
- << " source box span : " << bbSrc.span() << nl
- << " target box span : " << bbTgt.span() << nl
- << " source box : " << bbSrc << nl
- << " target box : " << bbTgt << nl
- << " inflated target box : " << bbTgtInf << endl;
+ )
+ << "Invalid interpolationMethod " << im
+ << ". Valid methods are:" << methods
+ << exit(FatalError);
}
-}
-
-
-template<class SourcePatch, class TargetPatch>
-void Foam::AMIInterpolation<SourcePatch, TargetPatch>::resetTree
-(
- const TargetPatch& tgtPatch
-)
-{
- // Clear the old octree
- treePtr_.clear();
-
-
- treeBoundBox bb(tgtPatch.points());
- bb.inflate(0.01);
- if (!treePtr_.valid())
- {
- treePtr_.reset
- (
- new indexedOctree<treeType>
- (
- treeType(false, tgtPatch),
- bb, // overall search domain
- 8, // maxLevel
- 10, // leaf size
- 3.0 // duplicity
- )
- );
- }
+ return method;
}
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
template<class SourcePatch, class TargetPatch>
void Foam::AMIInterpolation<SourcePatch, TargetPatch>::projectPointsToSurface
(
@@ -182,7 +161,8 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::projectPointsToSurface
{
FatalErrorIn
(
- "void Foam::projectPointsToSurface"
+ "void Foam::AMIInterpolation<SourcePatch, TargetPatch>::"
+ "projectPointsToSurface"
"("
"const searchableSurface&, "
"pointField&"
@@ -195,635 +175,6 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::projectPointsToSurface
}
-template<class SourcePatch, class TargetPatch>
-Foam::label Foam::AMIInterpolation<SourcePatch, TargetPatch>::findTargetFace
-(
- const label srcFaceI,
- const SourcePatch& srcPatch
-) const
-{
- label targetFaceI = -1;
-
- const pointField& srcPts = srcPatch.points();
- const face& srcFace = srcPatch[srcFaceI];
- const point srcPt = srcFace.centre(srcPts);
- const scalar srcFaceArea = srcMagSf_[srcFaceI];
-
- pointIndexHit sample = treePtr_->findNearest(srcPt, 10.0*srcFaceArea);
-
-
- if (debug)
- {
- Pout<< "Source point = " << srcPt << ", Sample point = "
- << sample.hitPoint() << ", Sample index = " << sample.index()
- << endl;
- }
-
- if (sample.hit())
- {
- targetFaceI = sample.index();
- }
-
- return targetFaceI;
-}
-
-
-template<class SourcePatch, class TargetPatch>
-void Foam::AMIInterpolation<SourcePatch, TargetPatch>::appendNbrFaces
-(
- const label faceI,
- const TargetPatch& patch,
- const DynamicList<label>& visitedFaces,
- DynamicList<label>& faceIDs
-) const
-{
- const labelList& nbrFaces = patch.faceFaces()[faceI];
-
- // filter out faces already visited from src face neighbours
- forAll(nbrFaces, i)
- {
- label nbrFaceI = nbrFaces[i];
- bool valid = true;
- forAll(visitedFaces, j)
- {
- if (nbrFaceI == visitedFaces[j])
- {
- valid = false;
- break;
- }
- }
-
- if (valid)
- {
- forAll(faceIDs, j)
- {
- if (nbrFaceI == faceIDs[j])
- {
- valid = false;
- break;
- }
- }
- }
-
- if (valid)
- {
- faceIDs.append(nbrFaceI);
- }
- }
-}
-
-
-template<class SourcePatch, class TargetPatch>
-bool Foam::AMIInterpolation<SourcePatch, TargetPatch>::processSourceFace
-(
- const SourcePatch& srcPatch,
- const TargetPatch& tgtPatch,
- const label srcFaceI,
- const label tgtStartFaceI,
-
- // list of tgt face neighbour faces
- DynamicList<label>& nbrFaces,
- // list of faces currently visited for srcFaceI to avoid multiple hits
- DynamicList<label>& visitedFaces,
-
- // temporary storage for addressing and weights
- List<DynamicList<label> >& srcAddr,
- List<DynamicList<scalar> >& srcWght,
- List<DynamicList<label> >& tgtAddr,
- List<DynamicList<scalar> >& tgtWght
-)
-{
- nbrFaces.clear();
- visitedFaces.clear();
-
- // append initial target face and neighbours
- nbrFaces.append(tgtStartFaceI);
- appendNbrFaces(tgtStartFaceI, tgtPatch, visitedFaces, nbrFaces);
-
- bool faceProcessed = false;
-
- do
- {
- // process new target face
- label tgtFaceI = nbrFaces.remove();
- visitedFaces.append(tgtFaceI);
- scalar area = interArea(srcFaceI, tgtFaceI, srcPatch, tgtPatch);
-
- // store when intersection area > 0
- if (area > 0)
- {
- srcAddr[srcFaceI].append(tgtFaceI);
- srcWght[srcFaceI].append(area);
-
- tgtAddr[tgtFaceI].append(srcFaceI);
- tgtWght[tgtFaceI].append(area);
-
- appendNbrFaces(tgtFaceI, tgtPatch, visitedFaces, nbrFaces);
-
- faceProcessed = true;
- }
-
- } while (nbrFaces.size() > 0);
-
- return faceProcessed;
-}
-
-
-template<class SourcePatch, class TargetPatch>
-void Foam::AMIInterpolation<SourcePatch, TargetPatch>::setNextFaces
-(
- label& startSeedI,
- label& srcFaceI,
- label& tgtFaceI,
- const SourcePatch& srcPatch0,
- const TargetPatch& tgtPatch0,
- const boolList& mapFlag,
- labelList& seedFaces,
- const DynamicList<label>& visitedFaces
-) const
-{
- const labelList& srcNbrFaces = srcPatch0.faceFaces()[srcFaceI];
-
- // set possible seeds for later use
- bool valuesSet = false;
- forAll(srcNbrFaces, i)
- {
- label faceS = srcNbrFaces[i];
-
- if (mapFlag[faceS] && seedFaces[faceS] == -1)
- {
- forAll(visitedFaces, j)
- {
- label faceT = visitedFaces[j];
- scalar area = interArea(faceS, faceT, srcPatch0, tgtPatch0);
-
- // Check that faces have enough overlap for robust walking
- if (area/srcMagSf_[srcFaceI] > faceAreaIntersect::tolerance())
- {
- // TODO - throwing area away - re-use in next iteration?
-
- seedFaces[faceS] = faceT;
-
- if (!valuesSet)
- {
- srcFaceI = faceS;
- tgtFaceI = faceT;
- valuesSet = true;
- }
- }
- }
- }
- }
-
- // set next src and tgt faces if not set above
- if (valuesSet)
- {
- return;
- }
- else
- {
- // try to use existing seed
- bool foundNextSeed = false;
- for (label faceI = startSeedI; faceI < mapFlag.size(); faceI++)
- {
- if (mapFlag[faceI])
- {
- if (!foundNextSeed)
- {
- startSeedI = faceI;
- foundNextSeed = true;
- }
-
- if (seedFaces[faceI] != -1)
- {
- srcFaceI = faceI;
- tgtFaceI = seedFaces[faceI];
-
- return;
- }
- }
- }
-
- // perform new search to find match
- if (debug)
- {
- Pout<< "Advancing front stalled: searching for new "
- << "target face" << endl;
- }
-
- foundNextSeed = false;
- for (label faceI = startSeedI; faceI < mapFlag.size(); faceI++)
- {
- if (mapFlag[faceI])
- {
- if (!foundNextSeed)
- {
- startSeedI = faceI + 1;
- foundNextSeed = true;
- }
-
- srcFaceI = faceI;
- tgtFaceI = findTargetFace(srcFaceI, srcPatch0);
-
- if (tgtFaceI >= 0)
- {
- return;
- }
- }
- }
-
- FatalErrorIn
- (
- "void Foam::AMIInterpolation<SourcePatch, TargetPatch>::"
- "setNextFaces"
- "("
- "label&, "
- "label&, "
- "label&, "
- "const SourcePatch&, "
- "const TargetPatch&, "
- "const boolList&, "
- "labelList&, "
- "const DynamicList<label>&"
- ") const"
- ) << "Unable to set source and target faces" << abort(FatalError);
- }
-}
-
-
-template<class SourcePatch, class TargetPatch>
-Foam::scalar Foam::AMIInterpolation<SourcePatch, TargetPatch>::interArea
-(
- const label srcFaceI,
- const label tgtFaceI,
- const SourcePatch& srcPatch,
- const TargetPatch& tgtPatch
-) const
-{
- const pointField& srcPoints = srcPatch.points();
- const pointField& tgtPoints = tgtPatch.points();
-
- // references to candidate faces
- const face& src = srcPatch[srcFaceI];
- const face& tgt = tgtPatch[tgtFaceI];
-
- // quick reject if either face has zero area
- // Note: do not used stored face areas for target patch
- if ((srcMagSf_[srcFaceI] < ROOTVSMALL) || (tgt.mag(tgtPoints) < ROOTVSMALL))
- {
- return 0.0;
- }
-
- // create intersection object
- faceAreaIntersect inter(srcPoints, tgtPoints, reverseTarget_);
-
- // crude resultant norm
- vector n(-src.normal(srcPoints));
- if (reverseTarget_)
- {
- n -= tgt.normal(tgtPoints);
- }
- else
- {
- n += tgt.normal(tgtPoints);
- }
- n *= 0.5;
-
- scalar area = 0;
- if (mag(n) > ROOTVSMALL)
- {
- area = inter.calc(src, tgt, n, triMode_);
- }
- else
- {
- WarningIn
- (
- "void Foam::AMIInterpolation<SourcePatch, TargetPatch>::"
- "interArea"
- "("
- "const label, "
- "const label, "
- "const SourcePatch&, "
- "const TargetPatch&"
- ") const"
- ) << "Invalid normal for source face " << srcFaceI
- << " points " << UIndirectList<point>(srcPoints, src)
- << " target face " << tgtFaceI
- << " points " << UIndirectList<point>(tgtPoints, tgt)
- << endl;
- }
-
-
- if ((debug > 1) && (area > 0))
- {
- writeIntersectionOBJ(area, src, tgt, srcPoints, tgtPoints);
- }
-
- return area;
-}
-
-
-template<class SourcePatch, class TargetPatch>
-void Foam::AMIInterpolation<SourcePatch, TargetPatch>::
-restartUncoveredSourceFace
-(
- const SourcePatch& srcPatch,
- const TargetPatch& tgtPatch,
- List<DynamicList<label> >& srcAddr,
- List<DynamicList<scalar> >& srcWght,
- List<DynamicList<label> >& tgtAddr,
- List<DynamicList<scalar> >& tgtWght
-)
-{
- // Collect all src faces with a low weight
- // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
- labelHashSet lowWeightFaces(100);
- forAll(srcWght, srcFaceI)
- {
- scalar s = sum(srcWght[srcFaceI]);
- scalar t = s/srcMagSf_[srcFaceI];
-
- if (t < 0.5)
- {
- lowWeightFaces.insert(srcFaceI);
- }
- }
-
- if (debug)
- {
- Pout<< "AMIInterpolation: restarting search on "
- << lowWeightFaces.size() << " faces since sum of weights < 0.5"
- << endl;
- }
-
- if (lowWeightFaces.size() > 0)
- {
- // Erase all the lowWeight source faces from the target
- // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
- DynamicList<label> okSrcFaces(10);
- DynamicList<scalar> okSrcWeights(10);
- forAll(tgtAddr, tgtFaceI)
- {
- okSrcFaces.clear();
- okSrcWeights.clear();
- DynamicList<label>& srcFaces = tgtAddr[tgtFaceI];
- DynamicList<scalar>& srcWeights = tgtWght[tgtFaceI];
- forAll(srcFaces, i)
- {
- if (!lowWeightFaces.found(srcFaces[i]))
- {
- okSrcFaces.append(srcFaces[i]);
- okSrcWeights.append(srcWeights[i]);
- }
- }
- if (okSrcFaces.size() < srcFaces.size())
- {
- srcFaces.transfer(okSrcFaces);
- srcWeights.transfer(okSrcWeights);
- }
- }
-
-
-
- // Restart search from best hit
- // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-
- // list of tgt face neighbour faces
- DynamicList<label> nbrFaces(10);
-
- // list of faces currently visited for srcFaceI to avoid multiple hits
- DynamicList<label> visitedFaces(10);
-
- forAllConstIter(labelHashSet, lowWeightFaces, iter)
- {
- label srcFaceI = iter.key();
- label tgtFaceI = findTargetFace(srcFaceI, srcPatch);
- if (tgtFaceI != -1)
- {
- //bool faceProcessed =
- processSourceFace
- (
- srcPatch,
- tgtPatch,
- srcFaceI,
- tgtFaceI,
-
- nbrFaces,
- visitedFaces,
-
- srcAddr,
- srcWght,
- tgtAddr,
- tgtWght
- );
- // ? Check faceProcessed to see if restarting has worked.
- }
- }
- }
-}
-
-
-template<class SourcePatch, class TargetPatch>
-void Foam::AMIInterpolation<SourcePatch, TargetPatch>::calcAddressing
-(
- const SourcePatch& srcPatch,
- const TargetPatch& tgtPatch,
- label srcFaceI,
- label tgtFaceI
-)
-{
- // Pre-size to handle early exit
- srcAddress_.setSize(srcPatch.size());
- srcWeights_.setSize(srcPatch.size());
- tgtAddress_.setSize(tgtPatch.size());
- tgtWeights_.setSize(tgtPatch.size());
-
-
- if (debug && (!srcPatch.size() || !tgtPatch.size()))
- {
- Pout<< "AMI: Patches not on processor: Source faces = "
- << srcPatch.size() << ", target faces = " << tgtPatch.size()
- << endl;
- }
-
- if (!srcPatch.size())
- {
- return;
- }
- else if (!tgtPatch.size())
- {
- WarningIn
- (
- "void Foam::AMIInterpolation<SourcePatch, TargetPatch>::"
- "calcAddressing"
- "("
- "const SourcePatch&, "
- "const TargetPatch&, "
- "label, "
- "label"
- ")"
- ) << "Have " << srcPatch.size() << " source faces but no target faces."
- << endl;
-
- return;
- }
-
- resetTree(tgtPatch);
-
- // temporary storage for addressing and weights
- List<DynamicList<label> > srcAddr(srcPatch.size());
- List<DynamicList<scalar> > srcWght(srcPatch.size());
- List<DynamicList<label> > tgtAddr(tgtPatch.size());
- List<DynamicList<scalar> > tgtWght(tgtPatch.size());
-
-
- // find initial face match using brute force/octree search
- // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
- if ((srcFaceI == -1) || (tgtFaceI == -1))
- {
- srcFaceI = 0;
- tgtFaceI = 0;
- bool foundFace = false;
- forAll(srcPatch, faceI)
- {
- tgtFaceI = findTargetFace(faceI, srcPatch);
- if (tgtFaceI >= 0)
- {
- srcFaceI = faceI;
- foundFace = true;
- break;
- }
- }
-
- if (!foundFace)
- {
- FatalErrorIn
- (
- "void Foam::AMIInterpolation<SourcePatch, TargetPatch>::"
- "calcAddressing"
- "("
- "const SourcePatch&, "
- "const TargetPatch&, "
- "label, "
- "label"
- ")"
- ) << "Unable to find initial target face" << abort(FatalError);
- }
- }
-
- if (debug)
- {
- Pout<< "AMI: initial target face = " << tgtFaceI << endl;
- }
-
-
- // construct weights and addressing
- // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
- label nFacesRemaining = srcPatch.size();
-
- // list of tgt face neighbour faces
- DynamicList<label> nbrFaces(10);
-
- // list of faces currently visited for srcFaceI to avoid multiple hits
- DynamicList<label> visitedFaces(10);
-
- // list to keep track of tgt faces used to seed src faces
- labelList seedFaces(nFacesRemaining, -1);
- seedFaces[srcFaceI] = tgtFaceI;
-
- // list to keep track of whether src face can be mapped
- boolList mapFlag(nFacesRemaining, true);
-
- // reset starting seed
- label startSeedI = 0;
-
- DynamicList<label> nonOverlapFaces;
- do
- {
- // Do advancing front starting from srcFaceI,tgtFaceI
- bool faceProcessed = processSourceFace
- (
- srcPatch,
- tgtPatch,
- srcFaceI,
- tgtFaceI,
-
- nbrFaces,
- visitedFaces,
-
- srcAddr,
- srcWght,
- tgtAddr,
- tgtWght
- );
-
- mapFlag[srcFaceI] = false;
-
- nFacesRemaining--;
-
- if (!faceProcessed)
- {
- nonOverlapFaces.append(srcFaceI);
- }
-
- // choose new src face from current src face neighbour
- if (nFacesRemaining > 0)
- {
- setNextFaces
- (
- startSeedI,
- srcFaceI,
- tgtFaceI,
- srcPatch,
- tgtPatch,
- mapFlag,
- seedFaces,
- visitedFaces
- );
- }
- } while (nFacesRemaining > 0);
-
- if (nonOverlapFaces.size() != 0)
- {
- Pout<< "AMI: " << nonOverlapFaces.size()
- << " non-overlap faces identified"
- << endl;
-
- srcNonOverlap_.transfer(nonOverlapFaces);
- }
-
-
- // Check for badly covered faces
- if (debug)
- {
- restartUncoveredSourceFace
- (
- srcPatch,
- tgtPatch,
- srcAddr,
- srcWght,
- tgtAddr,
- tgtWght
- );
- }
-
-
- // transfer data to persistent storage
- forAll(srcAddr, i)
- {
- srcAddress_[i].transfer(srcAddr[i]);
- srcWeights_[i].transfer(srcWght[i]);
- }
- forAll(tgtAddr, i)
- {
- tgtAddress_[i].transfer(tgtAddr[i]);
- tgtWeights_[i].transfer(tgtWght[i]);
- }
-}
-
-
template<class SourcePatch, class TargetPatch>
void Foam::AMIInterpolation<SourcePatch, TargetPatch>::normaliseWeights
(
@@ -857,7 +208,8 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::normaliseWeights
if (nFace)
{
- IInfo<< "AMI: Patch " << patchName << " weights min/max/average = "
+ IInfo<< "AMI: Patch " << patchName
+ << " sum(weights) min/max/average = "
<< gMin(wghtSum) << ", "
<< gMax(wghtSum) << ", "
<< gAverage(wghtSum) << endl;
@@ -1136,30 +488,23 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::AMIInterpolation
const SourcePatch& srcPatch,
const TargetPatch& tgtPatch,
const faceAreaIntersect::triangulationMode& triMode,
+ const interpolationMethod& method,
const bool reverseTarget
)
:
+ method_(method),
reverseTarget_(reverseTarget),
singlePatchProc_(-999),
srcAddress_(),
srcWeights_(),
srcWeightsSum_(),
- srcNonOverlap_(),
tgtAddress_(),
tgtWeights_(),
tgtWeightsSum_(),
- treePtr_(NULL),
triMode_(triMode),
srcMapPtr_(NULL),
tgtMapPtr_(NULL)
{
- label srcSize = returnReduce(srcPatch.size(), sumOp<label>());
- label tgtSize = returnReduce(tgtPatch.size(), sumOp<label>());
-
- IInfo<< "AMI: Creating addressing and weights between "
- << srcSize << " source faces and " << tgtSize << " target faces"
- << endl;
-
update(srcPatch, tgtPatch);
}
@@ -1171,30 +516,23 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::AMIInterpolation
const TargetPatch& tgtPatch,
const autoPtr<searchableSurface>& surfPtr,
const faceAreaIntersect::triangulationMode& triMode,
+ const interpolationMethod& method,
const bool reverseTarget
)
:
+ method_(method),
reverseTarget_(reverseTarget),
singlePatchProc_(-999),
srcAddress_(),
srcWeights_(),
srcWeightsSum_(),
- srcNonOverlap_(),
tgtAddress_(),
tgtWeights_(),
tgtWeightsSum_(),
- treePtr_(NULL),
triMode_(triMode),
srcMapPtr_(NULL),
tgtMapPtr_(NULL)
{
- label srcSize = returnReduce(srcPatch.size(), sumOp<label>());
- label tgtSize = returnReduce(tgtPatch.size(), sumOp<label>());
-
- IInfo<< "AMI: Creating addressing and weights between "
- << srcSize << " source faces and " << tgtSize << " target faces"
- << endl;
-
if (surfPtr.valid())
{
// create new patches for source and target
@@ -1264,16 +602,15 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::AMIInterpolation
const labelList& targetRestrictAddressing
)
:
+ method_(fineAMI.method_),
reverseTarget_(fineAMI.reverseTarget_),
singlePatchProc_(fineAMI.singlePatchProc_),
srcAddress_(),
srcWeights_(),
srcWeightsSum_(),
- srcNonOverlap_(),
tgtAddress_(),
tgtWeights_(),
tgtWeightsSum_(),
- treePtr_(NULL),
triMode_(fineAMI.triMode_),
srcMapPtr_(NULL),
tgtMapPtr_(NULL)
@@ -1456,11 +793,32 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::update
newTgtPoints
);
- checkPatches(srcPatch, newTgtPatch);
-
// calculate AMI interpolation
- calcAddressing(srcPatch, newTgtPatch);
+ {
+ autoPtr<AMIMethod<SourcePatch, TargetPatch> > AMIPtr
+ (
+ AMIMethod<SourcePatch, TargetPatch>::New
+ (
+ interpolationMethodToWord(method_),
+ srcPatch,
+ newTgtPatch,
+ srcMagSf_,
+ tgtMagSf_,
+ triMode_,
+ reverseTarget_
+ )
+ );
+
+ AMIPtr->calculate
+ (
+ srcAddress_,
+ srcWeights_,
+ tgtAddress_,
+ tgtWeights_
+ );
+ }
+
// Now
// ~~~
@@ -1552,9 +910,31 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::update
}
else
{
- checkPatches(srcPatch, tgtPatch);
- calcAddressing(srcPatch, tgtPatch);
+ // calculate AMI interpolation
+ {
+ autoPtr<AMIMethod<SourcePatch, TargetPatch> > AMIPtr
+ (
+ AMIMethod<SourcePatch, TargetPatch>::New
+ (
+ interpolationMethodToWord(method_),
+ srcPatch,
+ tgtPatch,
+ srcMagSf_,
+ tgtMagSf_,
+ triMode_,
+ reverseTarget_
+ )
+ );
+
+ AMIPtr->calculate
+ (
+ srcAddress_,
+ srcWeights_,
+ tgtAddress_,
+ tgtWeights_
+ );
+ }
normaliseWeights
(
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H
index 9f94fa41169965fbd12bece06b9c0afc0edb6563..e31eb91e1c6105deaaee94a117ce08f38d769422 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H
@@ -50,14 +50,11 @@ SourceFiles
#define AMIInterpolation_H
#include "className.H"
-#include "DynamicList.H"
#include "searchableSurface.H"
+#include "treeBoundBoxList.H"
#include "boolList.H"
#include "primitivePatch.H"
#include "faceAreaIntersect.H"
-#include "indexedOctree.H"
-#include "treeDataPrimitivePatch.H"
-#include "treeBoundBoxList.H"
#include "globalIndex.H"
#include "ops.H"
@@ -82,12 +79,38 @@ class AMIInterpolation
:
public AMIInterpolationName
{
- //- local typedef to octree tree-type
- typedef treeDataPrimitivePatch<TargetPatch> treeType;
+public:
+
+ // Public data types
+
+ //- Enumeration specifying interpolation method
+ enum interpolationMethod
+ {
+ imDirect,
+ imMapNearest,
+ imFaceAreaWeight
+ };
+
+ //- Convert interpolationMethod to word representation
+ static word interpolationMethodToWord
+ (
+ const interpolationMethod& method
+ );
+ //- Convert word to interpolationMethod
+ static interpolationMethod wordTointerpolationMethod
+ (
+ const word& method
+ );
+
+
+private:
// Private data
+ //- Interpolation method
+ interpolationMethod method_;
+
//- Flag to indicate that the two patches are co-directional and
// that the orientation of the target patch should be reversed
const bool reverseTarget_;
@@ -96,6 +119,7 @@ class AMIInterpolation
// cases
label singlePatchProc_;
+
// Source patch
//- Source face areas
@@ -110,10 +134,6 @@ class AMIInterpolation
//- Sum of weights of target faces per source face
scalarField srcWeightsSum_;
- //- Labels of faces that are not overlapped by any target faces
- // (should be empty for correct functioning)
- labelList srcNonOverlap_;
-
// Target patch
@@ -130,9 +150,6 @@ class AMIInterpolation
scalarField tgtWeightsSum_;
- //- Octree used to find face seeds
- autoPtr<indexedOctree<treeType> > treePtr_;
-
//- Face triangulation mode
const faceAreaIntersect::triangulationMode triMode_;
@@ -152,30 +169,6 @@ class AMIInterpolation
void operator=(const AMIInterpolation&);
- // Helper functions
-
- //- Write triangle intersection to OBJ file
- void writeIntersectionOBJ
- (
- const scalar area,
- const face& f1,
- const face& f2,
- const pointField& f1Points,
- const pointField& f2Points
- ) const;
-
- //- Check that patches are valid
- void checkPatches
- (
- const SourcePatch& srcPatch,
- const TargetPatch& tgtPatch
- ) const;
-
- //- Reset the octree for the traget patch face search
- void resetTree(const TargetPatch& tgtPatch);
-
-
-
// Parallel functionality
//- Calculate if patches are on multiple processors
@@ -229,83 +222,8 @@ class AMIInterpolation
) const;
- // Marching front
-
- //- Find face on target patch that overlaps source face
- label findTargetFace
- (
- const label srcFaceI,
- const SourcePatch& srcPatch
- ) const;
-
- //- Add faces neighbouring faceI to the ID list
- void appendNbrFaces
- (
- const label faceI,
- const TargetPatch& patch,
- const DynamicList<label>& visitedFaces,
- DynamicList<label>& faceIDs
- ) const;
-
- bool processSourceFace
- (
- const SourcePatch& srcPatch,
- const TargetPatch& tgtPatch,
- const label srcFaceI,
- const label tgtStartFaceI,
-
- DynamicList<label>& nbrFaces,
- DynamicList<label>& visitedFaces,
- List<DynamicList<label> >& srcAddr,
- List<DynamicList<scalar> >& srcWght,
- List<DynamicList<label> >& tgtAddr,
- List<DynamicList<scalar> >& tgtWght
- );
-
- void restartUncoveredSourceFace
- (
- const SourcePatch& srcPatch,
- const TargetPatch& tgtPatch,
- List<DynamicList<label> >& srcAddr,
- List<DynamicList<scalar> >& srcWght,
- List<DynamicList<label> >& tgtAddr,
- List<DynamicList<scalar> >& tgtWght
- );
-
- //- Set the source and target seed faces
- void setNextFaces
- (
- label& startSeedI,
- label& srcFaceI,
- label& tgtFaceI,
- const SourcePatch& srcPatch0,
- const TargetPatch& tgtPatch0,
- const boolList& mapFlag,
- labelList& seedFaces,
- const DynamicList<label>& visitedFaces
- ) const;
-
-
// Evaluation
- //- Area of intersection between source and target faces
- scalar interArea
- (
- const label srcFaceI,
- const label tgtFaceI,
- const SourcePatch& srcPatch,
- const TargetPatch& tgtPatch
- ) const;
-
- //- Calculate addressing
- void calcAddressing
- (
- const SourcePatch& srcPatch,
- const TargetPatch& tgtPatch,
- label srcFaceI = -1,
- label tgtFaceI = -1
- );
-
//- Normalise the (area) weights - suppresses numerical error in
// weights calculation
// NOTE: if area weights are incorrect by 'a significant amount'
@@ -352,6 +270,7 @@ public:
const SourcePatch& srcPatch,
const TargetPatch& tgtPatch,
const faceAreaIntersect::triangulationMode& triMode,
+ const interpolationMethod& method = imFaceAreaWeight,
const bool reverseTarget = false
);
@@ -362,6 +281,7 @@ public:
const TargetPatch& tgtPatch,
const autoPtr<searchableSurface>& surf,
const faceAreaIntersect::triangulationMode& triMode,
+ const interpolationMethod& method = imFaceAreaWeight,
const bool reverseTarget = false
);
@@ -378,6 +298,12 @@ public:
//- Destructor
~AMIInterpolation();
+ // Typedef to SourcePatch type this AMIInterplation is instantiated on
+ typedef SourcePatch sourcePatchType;
+
+ // Typedef to TargetPatch type this AMIInterplation is instantiated on
+ typedef TargetPatch targetPatchType;
+
// Member Functions
@@ -403,10 +329,6 @@ public:
// patch weights (i.e. the sum before normalisation)
inline const scalarField& srcWeightsSum() const;
- //- Labels of faces that are not overlapped by any target faces
- // (should be empty for correct functioning)
- inline const labelList& srcNonOverlap() const;
-
//- Source map pointer - valid only if singlePatchProc = -1
// This gets source data into a form to be consumed by
// tgtAddress, tgtWeights
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationI.H
index f2782bad886aef725f81f9b1fb74c6694b085092..84291ce326c799bb781bd44ed9f6fd1871614bc1 100644
--- a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationI.H
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationI.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -55,14 +55,6 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcWeights() const
}
-template<class SourcePatch, class TargetPatch>
-inline const Foam::labelList&
-Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcNonOverlap() const
-{
- return srcNonOverlap_;
-}
-
-
template<class SourcePatch, class TargetPatch>
inline const Foam::scalarField&
Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcWeightsSum() const
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.C
new file mode 100644
index 0000000000000000000000000000000000000000..03f6355870b2b89e37940ad64b89bb0f48d1e1d7
--- /dev/null
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.C
@@ -0,0 +1,339 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "AMIMethod.H"
+#include "meshTools.H"
+#include "mapDistribute.H"
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+void Foam::AMIMethod<SourcePatch, TargetPatch>::checkPatches() const
+{
+ if (debug && (!srcPatch_.size() || !tgtPatch_.size()))
+ {
+ Pout<< "AMI: Patches not on processor: Source faces = "
+ << srcPatch_.size() << ", target faces = " << tgtPatch_.size()
+ << endl;
+ }
+
+
+ const scalar maxBoundsError = 0.05;
+
+ // check bounds of source and target
+ boundBox bbSrc(srcPatch_.points(), srcPatch_.meshPoints(), true);
+ boundBox bbTgt(tgtPatch_.points(), tgtPatch_.meshPoints(), true);
+
+ boundBox bbTgtInf(bbTgt);
+ bbTgtInf.inflate(maxBoundsError);
+
+ if (!bbTgtInf.contains(bbSrc))
+ {
+ WarningIn("AMIMethod<SourcePatch, TargetPatch>::checkPatches()")
+ << "Source and target patch bounding boxes are not similar" << nl
+ << " source box span : " << bbSrc.span() << nl
+ << " target box span : " << bbTgt.span() << nl
+ << " source box : " << bbSrc << nl
+ << " target box : " << bbTgt << nl
+ << " inflated target box : " << bbTgtInf << endl;
+ }
+}
+
+
+template<class SourcePatch, class TargetPatch>
+bool Foam::AMIMethod<SourcePatch, TargetPatch>::initialise
+(
+ labelListList& srcAddress,
+ scalarListList& srcWeights,
+ labelListList& tgtAddress,
+ scalarListList& tgtWeights,
+ label& srcFaceI,
+ label& tgtFaceI
+)
+{
+ // set initial sizes for weights and addressing - must be done even if
+ // returns false below
+ srcAddress.setSize(srcPatch_.size());
+ srcWeights.setSize(srcPatch_.size());
+ tgtAddress.setSize(tgtPatch_.size());
+ tgtWeights.setSize(tgtPatch_.size());
+
+ // check that patch sizes are valid
+ if (!srcPatch_.size())
+ {
+ return false;
+ }
+ else if (!tgtPatch_.size())
+ {
+ WarningIn
+ (
+ "void Foam::AMIMethod<SourcePatch, TargetPatch>::initialise"
+ "("
+ "label&, "
+ "label&"
+ ")"
+ )
+ << srcPatch_.size() << " source faces but no target faces" << endl;
+
+ return false;
+ }
+
+ // reset the octree
+ resetTree();
+
+ // find initial face match using brute force/octree search
+ if ((srcFaceI == -1) || (tgtFaceI == -1))
+ {
+ srcFaceI = 0;
+ tgtFaceI = 0;
+ bool foundFace = false;
+ forAll(srcPatch_, faceI)
+ {
+ tgtFaceI = findTargetFace(faceI);
+ if (tgtFaceI >= 0)
+ {
+ srcFaceI = faceI;
+ foundFace = true;
+ break;
+ }
+ }
+
+ if (!foundFace)
+ {
+ FatalErrorIn
+ (
+ "void Foam::AMIMethod<SourcePatch, TargetPatch>::initialise"
+ "("
+ "label&, "
+ "label&"
+ ")"
+ ) << "Unable to find initial target face" << abort(FatalError);
+ }
+ }
+
+ if (debug)
+ {
+ Pout<< "AMI: initial target face = " << tgtFaceI << endl;
+ }
+
+ return true;
+}
+
+
+template<class SourcePatch, class TargetPatch>
+void Foam::AMIMethod<SourcePatch, TargetPatch>::writeIntersectionOBJ
+(
+ const scalar area,
+ const face& f1,
+ const face& f2,
+ const pointField& f1Points,
+ const pointField& f2Points
+) const
+{
+ static label count = 1;
+
+ const pointField f1pts = f1.points(f1Points);
+ const pointField f2pts = f2.points(f2Points);
+
+ Pout<< "Face intersection area (" << count << "):" << nl
+ << " f1 face = " << f1 << nl
+ << " f1 pts = " << f1pts << nl
+ << " f2 face = " << f2 << nl
+ << " f2 pts = " << f2pts << nl
+ << " area = " << area
+ << endl;
+
+ OFstream os("areas" + name(count) + ".obj");
+
+ forAll(f1pts, i)
+ {
+ meshTools::writeOBJ(os, f1pts[i]);
+ }
+ os<< "l";
+ forAll(f1pts, i)
+ {
+ os<< " " << i + 1;
+ }
+ os<< " 1" << endl;
+
+
+ forAll(f2pts, i)
+ {
+ meshTools::writeOBJ(os, f2pts[i]);
+ }
+ os<< "l";
+ forAll(f2pts, i)
+ {
+ os<< " " << f1pts.size() + i + 1;
+ }
+ os<< " " << f1pts.size() + 1 << endl;
+
+ count++;
+}
+
+
+template<class SourcePatch, class TargetPatch>
+void Foam::AMIMethod<SourcePatch, TargetPatch>::resetTree()
+{
+ // Clear the old octree
+ treePtr_.clear();
+
+ treeBoundBox bb(tgtPatch_.points());
+ bb.inflate(0.01);
+
+ if (!treePtr_.valid())
+ {
+ treePtr_.reset
+ (
+ new indexedOctree<treeType>
+ (
+ treeType(false, tgtPatch_),
+ bb, // overall search domain
+ 8, // maxLevel
+ 10, // leaf size
+ 3.0 // duplicity
+ )
+ );
+ }
+}
+
+
+template<class SourcePatch, class TargetPatch>
+Foam::label Foam::AMIMethod<SourcePatch, TargetPatch>::findTargetFace
+(
+ const label srcFaceI
+) const
+{
+ label targetFaceI = -1;
+
+ const pointField& srcPts = srcPatch_.points();
+ const face& srcFace = srcPatch_[srcFaceI];
+ const point srcPt = srcFace.centre(srcPts);
+ const scalar srcFaceArea = srcMagSf_[srcFaceI];
+
+ pointIndexHit sample = treePtr_->findNearest(srcPt, 10.0*srcFaceArea);
+
+
+ if (debug)
+ {
+ Pout<< "Source point = " << srcPt << ", Sample point = "
+ << sample.hitPoint() << ", Sample index = " << sample.index()
+ << endl;
+ }
+
+ if (sample.hit())
+ {
+ targetFaceI = sample.index();
+ }
+
+ return targetFaceI;
+}
+
+
+template<class SourcePatch, class TargetPatch>
+void Foam::AMIMethod<SourcePatch, TargetPatch>::appendNbrFaces
+(
+ const label faceI,
+ const TargetPatch& patch,
+ const DynamicList<label>& visitedFaces,
+ DynamicList<label>& faceIDs
+) const
+{
+ const labelList& nbrFaces = patch.faceFaces()[faceI];
+
+ // filter out faces already visited from src face neighbours
+ forAll(nbrFaces, i)
+ {
+ label nbrFaceI = nbrFaces[i];
+ bool valid = true;
+ forAll(visitedFaces, j)
+ {
+ if (nbrFaceI == visitedFaces[j])
+ {
+ valid = false;
+ break;
+ }
+ }
+
+ if (valid)
+ {
+ forAll(faceIDs, j)
+ {
+ if (nbrFaceI == faceIDs[j])
+ {
+ valid = false;
+ break;
+ }
+ }
+ }
+
+ if (valid)
+ {
+ faceIDs.append(nbrFaceI);
+ }
+ }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::AMIMethod<SourcePatch, TargetPatch>::AMIMethod
+(
+ const SourcePatch& srcPatch,
+ const TargetPatch& tgtPatch,
+ const scalarField& srcMagSf,
+ const scalarField& tgtMagSf,
+ const faceAreaIntersect::triangulationMode& triMode,
+ const bool reverseTarget
+)
+:
+ srcPatch_(srcPatch),
+ tgtPatch_(tgtPatch),
+ reverseTarget_(reverseTarget),
+ srcMagSf_(srcMagSf),
+ tgtMagSf_(tgtMagSf),
+ srcNonOverlap_(),
+ triMode_(triMode)
+{
+ checkPatches();
+
+ label srcSize = returnReduce(srcPatch_.size(), sumOp<label>());
+ label tgtSize = returnReduce(tgtPatch_.size(), sumOp<label>());
+
+ IInfo<< "AMI: Creating addressing and weights between "
+ << srcSize << " source faces and " << tgtSize << " target faces"
+ << endl;
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::AMIMethod<SourcePatch, TargetPatch>::~AMIMethod()
+{}
+
+
+// ************************************************************************* //
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.H
new file mode 100644
index 0000000000000000000000000000000000000000..7b6482b7c6cfe1850461b39895f088000f4e3119
--- /dev/null
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.H
@@ -0,0 +1,268 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::AMIMethod
+
+Description
+ Base class for Arbitrary Mesh Interface (AMI) methods
+
+SourceFiles
+ AMIMethod.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef AMIMethod_H
+#define AMIMethod_H
+
+#include "className.H"
+#include "DynamicList.H"
+#include "faceAreaIntersect.H"
+#include "indexedOctree.H"
+#include "treeDataPrimitivePatch.H"
+#include "treeBoundBoxList.H"
+#include "runTimeSelectionTables.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class AMIMethod Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class SourcePatch, class TargetPatch>
+class AMIMethod
+{
+
+private:
+
+ // Private Member Functions
+
+ //- Disallow default bitwise copy construct
+ AMIMethod(const AMIMethod&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const AMIMethod&);
+
+
+protected:
+
+ //- local typedef to octree tree-type
+ typedef treeDataPrimitivePatch<TargetPatch> treeType;
+
+
+ // Protected data
+
+ //- Reference to source patch
+ const SourcePatch& srcPatch_;
+
+ //- Reference to target patch
+ const SourcePatch& tgtPatch_;
+
+ //- Flag to indicate that the two patches are co-directional and
+ // that the orientation of the target patch should be reversed
+ const bool reverseTarget_;
+
+ //- Source face areas
+ const scalarField& srcMagSf_;
+
+ //- Target face areas
+ const scalarField& tgtMagSf_;
+
+ //- Labels of faces that are not overlapped by any target faces
+ // (should be empty for correct functioning)
+ labelList srcNonOverlap_;
+
+ //- Octree used to find face seeds
+ autoPtr<indexedOctree<treeType> > treePtr_;
+
+ //- Face triangulation mode
+ const faceAreaIntersect::triangulationMode triMode_;
+
+
+ // Protected Member Functions
+
+ // Helper functions
+
+ //- Check AMI patch coupling
+ void checkPatches() const;
+
+ //- Initialise and return true if all ok
+ bool initialise
+ (
+ labelListList& srcAddress,
+ scalarListList& srcWeights,
+ labelListList& tgtAddress,
+ scalarListList& tgtWeights,
+ label& srcFaceI,
+ label& tgtFaceI
+ );
+
+ //- Write triangle intersection to OBJ file
+ void writeIntersectionOBJ
+ (
+ const scalar area,
+ const face& f1,
+ const face& f2,
+ const pointField& f1Points,
+ const pointField& f2Points
+ ) const;
+
+
+ // Common AMI method functions
+
+ //- Reset the octree for the target patch face search
+ void resetTree();
+
+ //- Find face on target patch that overlaps source face
+ label findTargetFace(const label srcFaceI) const;
+
+ //- Add faces neighbouring faceI to the ID list
+ void appendNbrFaces
+ (
+ const label faceI,
+ const TargetPatch& patch,
+ const DynamicList<label>& visitedFaces,
+ DynamicList<label>& faceIDs
+ ) const;
+
+
+public:
+
+ //- Runtime type information
+ TypeName("AMIMethod");
+
+ //- Declare runtime constructor selection table
+ declareRunTimeSelectionTable
+ (
+ autoPtr,
+ AMIMethod,
+ components,
+ (
+ const SourcePatch& srcPatch,
+ const TargetPatch& tgtPatch,
+ const scalarField& srcMagSf,
+ const scalarField& tgtMagSf,
+ const faceAreaIntersect::triangulationMode& triMode,
+ const bool reverseTarget
+ ),
+ (srcPatch, tgtPatch, srcMagSf, tgtMagSf, triMode, reverseTarget)
+ );
+
+ //- Construct from components
+ AMIMethod
+ (
+ const SourcePatch& srcPatch,
+ const TargetPatch& tgtPatch,
+ const scalarField& srcMagSf,
+ const scalarField& tgtMagSf,
+ const faceAreaIntersect::triangulationMode& triMode,
+ const bool reverseTarget
+ );
+
+ //- Selector
+ static autoPtr<AMIMethod> New
+ (
+ const word& methodName,
+ const SourcePatch& srcPatch,
+ const TargetPatch& tgtPatch,
+ const scalarField& srcMagSf,
+ const scalarField& tgtMagSf,
+ const faceAreaIntersect::triangulationMode& triMode,
+ const bool reverseTarget
+ );
+
+
+ //- Destructor
+ virtual ~AMIMethod();
+
+
+ // Member Functions
+
+ // Access
+
+ //- Labels of faces that are not overlapped by any target faces
+ // Note: this should be empty for correct functioning
+ inline const labelList& srcNonOverlap() const;
+
+
+ // Manipulation
+
+ //- Update addressing and weights
+ virtual void calculate
+ (
+ labelListList& srcAddress,
+ scalarListList& srcWeights,
+ labelListList& tgtAddress,
+ scalarListList& tgtWeights,
+ label srcFaceI = -1,
+ label tgtFaceI = -1
+ ) = 0;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#define makeAMIMethod(AMIType) \
+ \
+ typedef AMIMethod<AMIType::sourcePatchType,AMIType::targetPatchType> \
+ AMIMethod##AMIType; \
+ \
+ defineNamedTemplateTypeNameAndDebug(AMIMethod##AMIType, 0); \
+ defineTemplateRunTimeSelectionTable(AMIMethod##AMIType, components);
+
+
+#define makeAMIMethodType(AMIType, Method) \
+ \
+ typedef Method<AMIType::sourcePatchType,AMIType::targetPatchType> \
+ Method##AMIType; \
+ \
+ defineNamedTemplateTypeNameAndDebug(Method##AMIType, 0); \
+ \
+ AMIMethod<AMIType::sourcePatchType,AMIType::targetPatchType>:: \
+ addcomponentsConstructorToTable<Method##AMIType> \
+ add##Method##AMIType##ConstructorToTable_;
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "AMIMethodI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+# include "AMIMethod.C"
+# include "AMIMethodNew.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodI.H
new file mode 100644
index 0000000000000000000000000000000000000000..6a86eccf9c6f383fbef90ccf0ca857297098d3f7
--- /dev/null
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodI.H
@@ -0,0 +1,34 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+template<class SourcePatch, class TargetPatch>
+inline const Foam::labelList&
+Foam::AMIMethod<SourcePatch, TargetPatch>::srcNonOverlap() const
+{
+ return srcNonOverlap_;
+}
+
+
+// ************************************************************************* //
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodNew.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodNew.C
new file mode 100644
index 0000000000000000000000000000000000000000..39cc1d3f8f1eff1039c5ff8169e7cfad38e5004b
--- /dev/null
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodNew.C
@@ -0,0 +1,84 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::autoPtr<Foam::AMIMethod<SourcePatch, TargetPatch> >
+Foam::AMIMethod<SourcePatch, TargetPatch>::New
+(
+ const word& methodName,
+ const SourcePatch& srcPatch,
+ const TargetPatch& tgtPatch,
+ const scalarField& srcMagSf,
+ const scalarField& tgtMagSf,
+ const faceAreaIntersect::triangulationMode& triMode,
+ const bool reverseTarget
+)
+{
+ if (debug)
+ {
+ Info<< "Selecting AMIMethod " << methodName << endl;
+ }
+
+ typename componentsConstructorTable::iterator cstrIter =
+ componentsConstructorTablePtr_->find(methodName);
+
+ if (cstrIter == componentsConstructorTablePtr_->end())
+ {
+ FatalErrorIn
+ (
+ "AMIMethod<SourcePatch, TargetPatch>::New"
+ "("
+ "const word&, "
+ "const SourcePatch&, "
+ "const TargetPatch&, "
+ "const scalarField&, "
+ "const scalarField&, "
+ "const faceAreaIntersect::triangulationMode&, "
+ "const bool"
+ ")"
+ ) << "Unknown AMIMethod type "
+ << methodName << nl << nl
+ << "Valid AMIMethod types are:" << nl
+ << componentsConstructorTablePtr_->sortedToc() << exit(FatalError);
+ }
+
+ return autoPtr<AMIMethod<SourcePatch, TargetPatch> >
+ (
+ cstrIter()
+ (
+ srcPatch,
+ tgtPatch,
+ srcMagSf,
+ tgtMagSf,
+ triMode,
+ reverseTarget
+ )
+ );
+}
+
+
+// ************************************************************************* //
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.C
new file mode 100644
index 0000000000000000000000000000000000000000..fa5d0b73e12ab3d23e742fbada71dd7436d3807d
--- /dev/null
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.C
@@ -0,0 +1,224 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "directAMI.H"
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+void Foam::directAMI<SourcePatch, TargetPatch>::appendToDirectSeeds
+(
+ boolList& mapFlag,
+ labelList& srcTgtSeed,
+ DynamicList<label>& srcSeeds,
+ DynamicList<label>& nonOverlapFaces,
+ label& srcFaceI,
+ label& tgtFaceI
+) const
+{
+ const labelList& srcNbr = this->srcPatch_.faceFaces()[srcFaceI];
+ const labelList& tgtNbr = this->tgtPatch_.faceFaces()[tgtFaceI];
+
+ const pointField& srcPoints = this->srcPatch_.points();
+ const pointField& tgtPoints = this->tgtPatch_.points();
+
+ const vectorField& srcCf = this->srcPatch_.faceCentres();
+
+ forAll(srcNbr, i)
+ {
+ label srcI = srcNbr[i];
+
+ if (mapFlag[srcI] && srcTgtSeed[srcI] == -1)
+ {
+ const face& srcF = this->srcPatch_[srcI];
+ const vector srcN = srcF.normal(srcPoints);
+
+ bool found = false;
+ forAll(tgtNbr, j)
+ {
+ label tgtI = tgtNbr[j];
+ const face& tgtF = this->tgtPatch_[tgtI];
+ pointHit ray = tgtF.ray(srcCf[srcI], srcN, tgtPoints);
+
+ if (ray.hit())
+ {
+ // new match - append to lists
+ found = true;
+
+ srcTgtSeed[srcI] = tgtI;
+ srcSeeds.append(srcI);
+
+ break;
+ }
+ }
+
+ if (!found)
+ {
+ // no match available for source face srcI
+ mapFlag[srcI] = false;
+ nonOverlapFaces.append(srcI);
+ }
+ }
+ }
+
+ if (srcSeeds.size())
+ {
+ srcFaceI = srcSeeds.remove();
+ tgtFaceI = srcTgtSeed[srcFaceI];
+ }
+ else
+ {
+ srcFaceI = -1;
+ tgtFaceI = -1;
+ }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::directAMI<SourcePatch, TargetPatch>::directAMI
+(
+ const SourcePatch& srcPatch,
+ const TargetPatch& tgtPatch,
+ const scalarField& srcMagSf,
+ const scalarField& tgtMagSf,
+ const faceAreaIntersect::triangulationMode& triMode,
+ const bool reverseTarget
+)
+:
+ AMIMethod<SourcePatch, TargetPatch>
+ (
+ srcPatch,
+ tgtPatch,
+ srcMagSf,
+ tgtMagSf,
+ triMode,
+ reverseTarget
+ )
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::directAMI<SourcePatch, TargetPatch>::~directAMI()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+void Foam::directAMI<SourcePatch, TargetPatch>::calculate
+(
+ labelListList& srcAddress,
+ scalarListList& srcWeights,
+ labelListList& tgtAddress,
+ scalarListList& tgtWeights,
+ label srcFaceI,
+ label tgtFaceI
+)
+{
+ bool ok =
+ this->initialise
+ (
+ srcAddress,
+ srcWeights,
+ tgtAddress,
+ tgtWeights,
+ srcFaceI,
+ tgtFaceI
+ );
+
+ if (!ok)
+ {
+ return;
+ }
+
+
+ // temporary storage for addressing and weights
+ List<DynamicList<label> > srcAddr(this->srcPatch_.size());
+ List<DynamicList<label> > tgtAddr(this->tgtPatch_.size());
+
+
+ // construct weights and addressing
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ // list of faces currently visited for srcFaceI to avoid multiple hits
+ DynamicList<label> srcSeeds(10);
+
+ // list to keep track of tgt faces used to seed src faces
+ labelList srcTgtSeed(srcAddr.size(), -1);
+ srcTgtSeed[srcFaceI] = tgtFaceI;
+
+ // list to keep track of whether src face can be mapped
+ boolList mapFlag(srcAddr.size(), true);
+
+ DynamicList<label> nonOverlapFaces;
+ do
+ {
+ srcAddr[srcFaceI].append(tgtFaceI);
+ tgtAddr[tgtFaceI].append(srcFaceI);
+
+ mapFlag[srcFaceI] = false;
+
+ // Do advancing front starting from srcFaceI, tgtFaceI
+ appendToDirectSeeds
+ (
+ mapFlag,
+ srcTgtSeed,
+ srcSeeds,
+ nonOverlapFaces,
+ srcFaceI,
+ tgtFaceI
+ );
+ } while (srcFaceI >= 0);
+
+ if (nonOverlapFaces.size() != 0)
+ {
+ Pout<< " AMI: " << nonOverlapFaces.size()
+ << " non-overlap faces identified"
+ << endl;
+
+ this->srcNonOverlap_.transfer(nonOverlapFaces);
+ }
+
+ // transfer data to persistent storage
+ forAll(srcAddr, i)
+ {
+ scalar magSf = this->srcMagSf_[i];
+ srcAddress[i].transfer(srcAddr[i]);
+ srcWeights[i] = scalarList(1, magSf);
+ }
+ forAll(tgtAddr, i)
+ {
+ scalar magSf = this->tgtMagSf_[i];
+ tgtAddress[i].transfer(tgtAddr[i]);
+ tgtWeights[i] = scalarList(1, magSf);
+ }
+}
+
+
+// ************************************************************************* //
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.H
new file mode 100644
index 0000000000000000000000000000000000000000..e7842ff9f8957b2557bf9245a6166b608655f38d
--- /dev/null
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.H
@@ -0,0 +1,144 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::directAMI
+
+Description
+ Direct mapped Arbitrary Mesh Interface (AMI) method
+
+SourceFiles
+ directAMI.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef directAMI_H
+#define directAMI_H
+
+#include "AMIMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class directAMI Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class SourcePatch, class TargetPatch>
+class directAMI
+:
+ public AMIMethod<SourcePatch, TargetPatch>
+{
+
+private:
+
+ // Private Member Functions
+
+ //- Disallow default bitwise copy construct
+ directAMI(const directAMI&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const directAMI&);
+
+ // Marching front
+
+ //- Append to list of src face seed indices
+ void appendToDirectSeeds
+ (
+ boolList& mapFlag,
+ labelList& srcTgtSeed,
+ DynamicList<label>& srcSeeds,
+ DynamicList<label>& nonOverlapFaces,
+ label& srcFaceI,
+ label& tgtFaceI
+ ) const;
+
+
+ // Evaluation
+
+ //- Area of intersection between source and target faces
+ scalar interArea
+ (
+ const label srcFaceI,
+ const label tgtFaceI
+ ) const;
+
+
+public:
+
+ //- Runtime type information
+ TypeName("directAMI");
+
+
+ // Constructors
+
+ //- Construct from components
+ directAMI
+ (
+ const SourcePatch& srcPatch,
+ const TargetPatch& tgtPatch,
+ const scalarField& srcMagSf,
+ const scalarField& tgtMagSf,
+ const faceAreaIntersect::triangulationMode& triMode,
+ const bool reverseTarget = false
+ );
+
+
+ //- Destructor
+ virtual ~directAMI();
+
+
+ // Member Functions
+
+ // Manipulation
+
+ //- Update addressing and weights
+ virtual void calculate
+ (
+ labelListList& srcAddress,
+ scalarListList& srcWeights,
+ labelListList& tgtAddress,
+ scalarListList& tgtWeights,
+ label srcFaceI = -1,
+ label tgtFaceI = -1
+ );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+# include "directAMI.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.C
new file mode 100644
index 0000000000000000000000000000000000000000..6aea34e945827dbce24cc272e9dba3c811cac77f
--- /dev/null
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.C
@@ -0,0 +1,553 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "faceAreaWeightAMI.H"
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::processSourceFace
+(
+ const label srcFaceI,
+ const label tgtStartFaceI,
+
+ // list of tgt face neighbour faces
+ DynamicList<label>& nbrFaces,
+ // list of faces currently visited for srcFaceI to avoid multiple hits
+ DynamicList<label>& visitedFaces,
+
+ // temporary storage for addressing and weights
+ List<DynamicList<label> >& srcAddr,
+ List<DynamicList<scalar> >& srcWght,
+ List<DynamicList<label> >& tgtAddr,
+ List<DynamicList<scalar> >& tgtWght
+)
+{
+ nbrFaces.clear();
+ visitedFaces.clear();
+
+ // append initial target face and neighbours
+ nbrFaces.append(tgtStartFaceI);
+ this->appendNbrFaces
+ (
+ tgtStartFaceI,
+ this->tgtPatch_,
+ visitedFaces,
+ nbrFaces
+ );
+
+ bool faceProcessed = false;
+
+ do
+ {
+ // process new target face
+ label tgtFaceI = nbrFaces.remove();
+ visitedFaces.append(tgtFaceI);
+ scalar area = interArea(srcFaceI, tgtFaceI);
+
+ // store when intersection area > 0
+ if (area > 0)
+ {
+ srcAddr[srcFaceI].append(tgtFaceI);
+ srcWght[srcFaceI].append(area);
+
+ tgtAddr[tgtFaceI].append(srcFaceI);
+ tgtWght[tgtFaceI].append(area);
+
+ this->appendNbrFaces
+ (
+ tgtFaceI,
+ this->tgtPatch_,
+ visitedFaces,
+ nbrFaces
+ );
+
+ faceProcessed = true;
+ }
+
+ } while (nbrFaces.size() > 0);
+
+ return faceProcessed;
+}
+
+
+template<class SourcePatch, class TargetPatch>
+void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces
+(
+ label& startSeedI,
+ label& srcFaceI,
+ label& tgtFaceI,
+ const boolList& mapFlag,
+ labelList& seedFaces,
+ const DynamicList<label>& visitedFaces
+) const
+{
+ const labelList& srcNbrFaces = this->srcPatch_.faceFaces()[srcFaceI];
+
+ // set possible seeds for later use
+ bool valuesSet = false;
+ forAll(srcNbrFaces, i)
+ {
+ label faceS = srcNbrFaces[i];
+
+ if (mapFlag[faceS] && seedFaces[faceS] == -1)
+ {
+ forAll(visitedFaces, j)
+ {
+ label faceT = visitedFaces[j];
+ scalar area = interArea(faceS, faceT);
+ scalar areaTotal = this->srcMagSf_[srcFaceI];
+
+ // Check that faces have enough overlap for robust walking
+ if (area/areaTotal > faceAreaIntersect::tolerance())
+ {
+ // TODO - throwing area away - re-use in next iteration?
+
+ seedFaces[faceS] = faceT;
+
+ if (!valuesSet)
+ {
+ srcFaceI = faceS;
+ tgtFaceI = faceT;
+ valuesSet = true;
+ }
+ }
+ }
+ }
+ }
+
+ // set next src and tgt faces if not set above
+ if (valuesSet)
+ {
+ return;
+ }
+ else
+ {
+ // try to use existing seed
+ bool foundNextSeed = false;
+ for (label faceI = startSeedI; faceI < mapFlag.size(); faceI++)
+ {
+ if (mapFlag[faceI])
+ {
+ if (!foundNextSeed)
+ {
+ startSeedI = faceI;
+ foundNextSeed = true;
+ }
+
+ if (seedFaces[faceI] != -1)
+ {
+ srcFaceI = faceI;
+ tgtFaceI = seedFaces[faceI];
+
+ return;
+ }
+ }
+ }
+
+ // perform new search to find match
+ if (debug)
+ {
+ Pout<< "Advancing front stalled: searching for new "
+ << "target face" << endl;
+ }
+
+ foundNextSeed = false;
+ for (label faceI = startSeedI; faceI < mapFlag.size(); faceI++)
+ {
+ if (mapFlag[faceI])
+ {
+ if (!foundNextSeed)
+ {
+ startSeedI = faceI + 1;
+ foundNextSeed = true;
+ }
+
+ srcFaceI = faceI;
+ tgtFaceI = this->findTargetFace(srcFaceI);
+
+ if (tgtFaceI >= 0)
+ {
+ return;
+ }
+ }
+ }
+
+ FatalErrorIn
+ (
+ "void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::"
+ "setNextFaces"
+ "("
+ "label&, "
+ "label&, "
+ "label&, "
+ "const boolList&, "
+ "labelList&, "
+ "const DynamicList<label>&"
+ ") const"
+ ) << "Unable to set source and target faces" << abort(FatalError);
+ }
+}
+
+
+template<class SourcePatch, class TargetPatch>
+Foam::scalar Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::interArea
+(
+ const label srcFaceI,
+ const label tgtFaceI
+) const
+{
+ const pointField& srcPoints = this->srcPatch_.points();
+ const pointField& tgtPoints = this->tgtPatch_.points();
+
+ // references to candidate faces
+ const face& src = this->srcPatch_[srcFaceI];
+ const face& tgt = this->tgtPatch_[tgtFaceI];
+
+ // quick reject if either face has zero area
+ // Note: do not used stored face areas for target patch
+ if
+ (
+ (this->srcMagSf_[srcFaceI] < ROOTVSMALL)
+ || (tgt.mag(tgtPoints) < ROOTVSMALL)
+ )
+ {
+ return 0.0;
+ }
+
+ // create intersection object
+ faceAreaIntersect inter(srcPoints, tgtPoints, this->reverseTarget_);
+
+ // crude resultant norm
+ vector n(-src.normal(srcPoints));
+ if (this->reverseTarget_)
+ {
+ n -= tgt.normal(tgtPoints);
+ }
+ else
+ {
+ n += tgt.normal(tgtPoints);
+ }
+ n *= 0.5;
+
+ scalar area = 0;
+ if (mag(n) > ROOTVSMALL)
+ {
+ area = inter.calc(src, tgt, n, this->triMode_);
+ }
+ else
+ {
+ WarningIn
+ (
+ "void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::"
+ "interArea"
+ "("
+ "const label, "
+ "const label"
+ ") const"
+ ) << "Invalid normal for source face " << srcFaceI
+ << " points " << UIndirectList<point>(srcPoints, src)
+ << " target face " << tgtFaceI
+ << " points " << UIndirectList<point>(tgtPoints, tgt)
+ << endl;
+ }
+
+
+ if ((debug > 1) && (area > 0))
+ {
+ this->writeIntersectionOBJ(area, src, tgt, srcPoints, tgtPoints);
+ }
+
+ return area;
+}
+
+
+template<class SourcePatch, class TargetPatch>
+void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::
+restartUncoveredSourceFace
+(
+ List<DynamicList<label> >& srcAddr,
+ List<DynamicList<scalar> >& srcWght,
+ List<DynamicList<label> >& tgtAddr,
+ List<DynamicList<scalar> >& tgtWght
+)
+{
+ // Collect all src faces with a low weight
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ labelHashSet lowWeightFaces(100);
+ forAll(srcWght, srcFaceI)
+ {
+ scalar s = sum(srcWght[srcFaceI]);
+ scalar t = s/this->srcMagSf_[srcFaceI];
+
+ if (t < 0.5)
+ {
+ lowWeightFaces.insert(srcFaceI);
+ }
+ }
+
+ if (debug)
+ {
+ Pout<< "faceAreaWeightAMI: restarting search on "
+ << lowWeightFaces.size() << " faces since sum of weights < 0.5"
+ << endl;
+ }
+
+ if (lowWeightFaces.size() > 0)
+ {
+ // Erase all the lowWeight source faces from the target
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ DynamicList<label> okSrcFaces(10);
+ DynamicList<scalar> okSrcWeights(10);
+ forAll(tgtAddr, tgtFaceI)
+ {
+ okSrcFaces.clear();
+ okSrcWeights.clear();
+ DynamicList<label>& srcFaces = tgtAddr[tgtFaceI];
+ DynamicList<scalar>& srcWeights = tgtWght[tgtFaceI];
+ forAll(srcFaces, i)
+ {
+ if (!lowWeightFaces.found(srcFaces[i]))
+ {
+ okSrcFaces.append(srcFaces[i]);
+ okSrcWeights.append(srcWeights[i]);
+ }
+ }
+ if (okSrcFaces.size() < srcFaces.size())
+ {
+ srcFaces.transfer(okSrcFaces);
+ srcWeights.transfer(okSrcWeights);
+ }
+ }
+
+
+
+ // Restart search from best hit
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ // list of tgt face neighbour faces
+ DynamicList<label> nbrFaces(10);
+
+ // list of faces currently visited for srcFaceI to avoid multiple hits
+ DynamicList<label> visitedFaces(10);
+
+ forAllConstIter(labelHashSet, lowWeightFaces, iter)
+ {
+ label srcFaceI = iter.key();
+ label tgtFaceI = this->findTargetFace(srcFaceI);
+ if (tgtFaceI != -1)
+ {
+ //bool faceProcessed =
+ processSourceFace
+ (
+ srcFaceI,
+ tgtFaceI,
+
+ nbrFaces,
+ visitedFaces,
+
+ srcAddr,
+ srcWght,
+ tgtAddr,
+ tgtWght
+ );
+ // ? Check faceProcessed to see if restarting has worked.
+ }
+ }
+ }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::faceAreaWeightAMI
+(
+ const SourcePatch& srcPatch,
+ const TargetPatch& tgtPatch,
+ const scalarField& srcMagSf,
+ const scalarField& tgtMagSf,
+ const faceAreaIntersect::triangulationMode& triMode,
+ const bool reverseTarget
+)
+:
+ AMIMethod<SourcePatch, TargetPatch>
+ (
+ srcPatch,
+ tgtPatch,
+ srcMagSf,
+ tgtMagSf,
+ triMode,
+ reverseTarget
+ )
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::~faceAreaWeightAMI()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calculate
+(
+ labelListList& srcAddress,
+ scalarListList& srcWeights,
+ labelListList& tgtAddress,
+ scalarListList& tgtWeights,
+ label srcFaceI,
+ label tgtFaceI
+)
+{
+ bool ok =
+ this->initialise
+ (
+ srcAddress,
+ srcWeights,
+ tgtAddress,
+ tgtWeights,
+ srcFaceI,
+ tgtFaceI
+ );
+
+ if (!ok)
+ {
+ return;
+ }
+
+ // temporary storage for addressing and weights
+ List<DynamicList<label> > srcAddr(this->srcPatch_.size());
+ List<DynamicList<scalar> > srcWght(srcAddr.size());
+ List<DynamicList<label> > tgtAddr(this->tgtPatch_.size());
+ List<DynamicList<scalar> > tgtWght(tgtAddr.size());
+
+
+ // construct weights and addressing
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ label nFacesRemaining = srcAddr.size();
+
+ // list of tgt face neighbour faces
+ DynamicList<label> nbrFaces(10);
+
+ // list of faces currently visited for srcFaceI to avoid multiple hits
+ DynamicList<label> visitedFaces(10);
+
+ // list to keep track of tgt faces used to seed src faces
+ labelList seedFaces(nFacesRemaining, -1);
+ seedFaces[srcFaceI] = tgtFaceI;
+
+ // list to keep track of whether src face can be mapped
+ boolList mapFlag(nFacesRemaining, true);
+
+ // reset starting seed
+ label startSeedI = 0;
+
+ DynamicList<label> nonOverlapFaces;
+ do
+ {
+ // Do advancing front starting from srcFaceI,tgtFaceI
+ bool faceProcessed = processSourceFace
+ (
+ srcFaceI,
+ tgtFaceI,
+
+ nbrFaces,
+ visitedFaces,
+
+ srcAddr,
+ srcWght,
+ tgtAddr,
+ tgtWght
+ );
+
+ mapFlag[srcFaceI] = false;
+
+ nFacesRemaining--;
+
+ if (!faceProcessed)
+ {
+ nonOverlapFaces.append(srcFaceI);
+ }
+
+ // choose new src face from current src face neighbour
+ if (nFacesRemaining > 0)
+ {
+ setNextFaces
+ (
+ startSeedI,
+ srcFaceI,
+ tgtFaceI,
+ mapFlag,
+ seedFaces,
+ visitedFaces
+ );
+ }
+ } while (nFacesRemaining > 0);
+
+ if (nonOverlapFaces.size() != 0)
+ {
+ Pout<< " AMI: " << nonOverlapFaces.size()
+ << " non-overlap faces identified"
+ << endl;
+
+ this->srcNonOverlap_.transfer(nonOverlapFaces);
+ }
+
+
+ // Check for badly covered faces
+ if (debug)
+ {
+ restartUncoveredSourceFace
+ (
+ srcAddr,
+ srcWght,
+ tgtAddr,
+ tgtWght
+ );
+ }
+
+
+ // transfer data to persistent storage
+ forAll(srcAddr, i)
+ {
+ srcAddress[i].transfer(srcAddr[i]);
+ srcWeights[i].transfer(srcWght[i]);
+ }
+ forAll(tgtAddr, i)
+ {
+ tgtAddress[i].transfer(tgtAddr[i]);
+ tgtWeights[i].transfer(tgtWght[i]);
+ }
+}
+
+
+// ************************************************************************* //
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.H
new file mode 100644
index 0000000000000000000000000000000000000000..82061d2727ae239bbcc238a6ddf009f72940108f
--- /dev/null
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.H
@@ -0,0 +1,166 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::faceAreaWeightAMI
+
+Description
+ Face area weighted Arbitrary Mesh Interface (AMI) method
+
+SourceFiles
+ faceAreaWeightAMI.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef faceAreaWeightAMI_H
+#define faceAreaWeightAMI_H
+
+#include "AMIMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class faceAreaWeightAMI Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class SourcePatch, class TargetPatch>
+class faceAreaWeightAMI
+:
+ public AMIMethod<SourcePatch, TargetPatch>
+{
+
+private:
+
+ // Private Member Functions
+
+ //- Disallow default bitwise copy construct
+ faceAreaWeightAMI(const faceAreaWeightAMI&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const faceAreaWeightAMI&);
+
+ // Marching front
+
+ //- Determine overlap contributions for source face srcFaceI
+ bool processSourceFace
+ (
+ const label srcFaceI,
+ const label tgtStartFaceI,
+ DynamicList<label>& nbrFaces,
+ DynamicList<label>& visitedFaces,
+ List<DynamicList<label> >& srcAddr,
+ List<DynamicList<scalar> >& srcWght,
+ List<DynamicList<label> >& tgtAddr,
+ List<DynamicList<scalar> >& tgtWght
+ );
+
+ //- Attempt to re-evaluate source faces that have not been included
+ void restartUncoveredSourceFace
+ (
+ List<DynamicList<label> >& srcAddr,
+ List<DynamicList<scalar> >& srcWght,
+ List<DynamicList<label> >& tgtAddr,
+ List<DynamicList<scalar> >& tgtWght
+ );
+
+ //- Set the source and target seed faces
+ void setNextFaces
+ (
+ label& startSeedI,
+ label& srcFaceI,
+ label& tgtFaceI,
+ const boolList& mapFlag,
+ labelList& seedFaces,
+ const DynamicList<label>& visitedFaces
+ ) const;
+
+
+ // Evaluation
+
+ //- Area of intersection between source and target faces
+ scalar interArea
+ (
+ const label srcFaceI,
+ const label tgtFaceI
+ ) const;
+
+
+public:
+
+ //- Runtime type information
+ TypeName("faceAreaWeightAMI");
+
+
+ // Constructors
+
+ //- Construct from components
+ faceAreaWeightAMI
+ (
+ const SourcePatch& srcPatch,
+ const TargetPatch& tgtPatch,
+ const scalarField& srcMagSf,
+ const scalarField& tgtMagSf,
+ const faceAreaIntersect::triangulationMode& triMode,
+ const bool reverseTarget = false
+ );
+
+
+ //- Destructor
+ virtual ~faceAreaWeightAMI();
+
+
+ // Member Functions
+
+ // Manipulation
+
+ //- Update addressing and weights
+ virtual void calculate
+ (
+ labelListList& srcAddress,
+ scalarListList& srcWeights,
+ labelListList& tgtAddress,
+ scalarListList& tgtWeights,
+ label srcFaceI = -1,
+ label tgtFaceI = -1
+ );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+# include "faceAreaWeightAMI.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.C
new file mode 100644
index 0000000000000000000000000000000000000000..61ed01727e2a42e8109661aae7f76befbcd8dd32
--- /dev/null
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.C
@@ -0,0 +1,343 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "mapNearestAMI.H"
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+void Foam::mapNearestAMI<SourcePatch, TargetPatch>::findNearestFace
+(
+ const SourcePatch& srcPatch,
+ const TargetPatch& tgtPatch,
+ const label& srcFaceI,
+ label& tgtFaceI
+) const
+{
+ const vectorField& srcCf = srcPatch.faceCentres();
+ const vectorField& tgtCf = tgtPatch.faceCentres();
+
+ const vector srcP = srcCf[srcFaceI];
+
+ DynamicList<label> tgtFaces(10);
+ tgtFaces.append(tgtFaceI);
+
+ DynamicList<label> visitedFaces(10);
+
+ scalar d = GREAT;
+
+ do
+ {
+ label tgtI = tgtFaces.remove();
+ visitedFaces.append(tgtI);
+
+ scalar dTest = magSqr(tgtCf[tgtI] - srcP);
+ if (dTest < d)
+ {
+ tgtFaceI = tgtI;
+ d = dTest;
+
+ this->appendNbrFaces
+ (
+ tgtFaceI,
+ tgtPatch,
+ visitedFaces,
+ tgtFaces
+ );
+ }
+
+ } while (tgtFaces.size() > 0);
+}
+
+
+template<class SourcePatch, class TargetPatch>
+void Foam::mapNearestAMI<SourcePatch, TargetPatch>::setNextNearestFaces
+(
+ boolList& mapFlag,
+ label& startSeedI,
+ label& srcFaceI,
+ label& tgtFaceI
+) const
+{
+ const labelList& srcNbr = this->srcPatch_.faceFaces()[srcFaceI];
+
+ srcFaceI = -1;
+
+ forAll(srcNbr, i)
+ {
+ label faceI = srcNbr[i];
+ if (mapFlag[faceI])
+ {
+ srcFaceI = faceI;
+ startSeedI = faceI + 1;
+
+ return;
+ }
+ }
+
+ forAll(mapFlag, srcFaceI)
+ {
+ if (!mapFlag[srcFaceI])
+ {
+ tgtFaceI = this->findTargetFace(srcFaceI);
+ if (tgtFaceI == -1)
+ {
+ const vectorField& srcCf = this->srcPatch_.faceCentres();
+
+ FatalErrorIn
+ (
+ "void Foam::mapNearestAMI<SourcePatch, TargetPatch>::"
+ "setNextNearestFaces"
+ "("
+ "boolList&, "
+ "label&, "
+ "label&, "
+ "label&"
+ ") const"
+ )
+ << "Unable to find target face for source face "
+ << srcFaceI << " with face centre " << srcCf[srcFaceI]
+ << abort(FatalError);
+ }
+
+ break;
+ }
+ }
+}
+
+
+template<class SourcePatch, class TargetPatch>
+Foam::label Foam::mapNearestAMI<SourcePatch, TargetPatch>::findMappedSrcFace
+(
+ const label tgtFaceI,
+ const List<DynamicList<label> >& tgtToSrc
+) const
+{
+ DynamicList<label> testFaces(10);
+ DynamicList<label> visitedFaces(10);
+
+ testFaces.append(tgtFaceI);
+
+ do
+ {
+ // search target tgtFaceI neighbours for match with source face
+ label tgtI = testFaces.remove();
+
+ if (findIndex(visitedFaces, tgtI) == -1)
+ {
+ visitedFaces.append(tgtI);
+
+ if (tgtToSrc[tgtI].size())
+ {
+ return tgtToSrc[tgtI][0];
+ }
+ else
+ {
+ const labelList& nbrFaces = this->tgtPatch_.faceFaces()[tgtI];
+
+ forAll(nbrFaces, i)
+ {
+ if (findIndex(visitedFaces, nbrFaces[i]) == -1)
+ {
+ testFaces.append(nbrFaces[i]);
+ }
+ }
+ }
+ }
+ } while (testFaces.size());
+
+ // did not find any match - should not be possible to get here!
+ return -1;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::mapNearestAMI<SourcePatch, TargetPatch>::mapNearestAMI
+(
+ const SourcePatch& srcPatch,
+ const TargetPatch& tgtPatch,
+ const scalarField& srcMagSf,
+ const scalarField& tgtMagSf,
+ const faceAreaIntersect::triangulationMode& triMode,
+ const bool reverseTarget
+)
+:
+ AMIMethod<SourcePatch, TargetPatch>
+ (
+ srcPatch,
+ tgtPatch,
+ srcMagSf,
+ tgtMagSf,
+ triMode,
+ reverseTarget
+ )
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::mapNearestAMI<SourcePatch, TargetPatch>::~mapNearestAMI()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+void Foam::mapNearestAMI<SourcePatch, TargetPatch>::calculate
+(
+ labelListList& srcAddress,
+ scalarListList& srcWeights,
+ labelListList& tgtAddress,
+ scalarListList& tgtWeights,
+ label srcFaceI,
+ label tgtFaceI
+)
+{
+ bool ok =
+ this->initialise
+ (
+ srcAddress,
+ srcWeights,
+ tgtAddress,
+ tgtWeights,
+ srcFaceI,
+ tgtFaceI
+ );
+
+ if (!ok)
+ {
+ return;
+ }
+
+
+ // temporary storage for addressing and weights
+ List<DynamicList<label> > srcAddr(this->srcPatch_.size());
+ List<DynamicList<label> > tgtAddr(this->tgtPatch_.size());
+
+
+ // construct weights and addressing
+ // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+ // list to keep track of whether src face can be mapped
+ boolList mapFlag(srcAddr.size(), true);
+
+ // reset starting seed
+ label startSeedI = 0;
+
+ DynamicList<label> nonOverlapFaces;
+ do
+ {
+ findNearestFace(this->srcPatch_, this->tgtPatch_, srcFaceI, tgtFaceI);
+
+ srcAddr[srcFaceI].append(tgtFaceI);
+ tgtAddr[tgtFaceI].append(srcFaceI);
+
+ mapFlag[srcFaceI] = false;
+
+ // Do advancing front starting from srcFaceI, tgtFaceI
+ setNextNearestFaces
+ (
+ mapFlag,
+ startSeedI,
+ srcFaceI,
+ tgtFaceI
+ );
+ } while (srcFaceI >= 0);
+
+
+ // for the case of multiple source faces per target face, select the
+ // nearest source face only and discard the others
+ const vectorField& srcCf = this->srcPatch_.faceCentres();
+ const vectorField& tgtCf = this->tgtPatch_.faceCentres();
+
+ forAll(tgtAddr, targetFaceI)
+ {
+ if (tgtAddr[targetFaceI].size() > 1)
+ {
+ const vector& tgtC = tgtCf[tgtFaceI];
+
+ DynamicList<label>& srcFaces = tgtAddr[targetFaceI];
+
+ label srcFaceI = srcFaces[0];
+ scalar d = magSqr(tgtC - srcCf[srcFaceI]);
+
+ for (label i = 1; i < srcFaces.size(); i++)
+ {
+ label srcI = srcFaces[i];
+ scalar dNew = magSqr(tgtC - srcCf[srcI]);
+ if (dNew < d)
+ {
+ d = dNew;
+ srcFaceI = srcI;
+ }
+ }
+
+ srcFaces.clear();
+ srcFaces.append(srcFaceI);
+ }
+ }
+
+ // If there are more target faces than source faces, some target faces
+ // might not yet be mapped
+ forAll(tgtAddr, tgtFaceI)
+ {
+ if (tgtAddr[tgtFaceI].empty())
+ {
+ label srcFaceI = findMappedSrcFace(tgtFaceI, tgtAddr);
+
+ // note - reversed search from src->tgt to tgt->src
+ findNearestFace
+ (
+ this->tgtPatch_,
+ this->srcPatch_,
+ tgtFaceI,
+ srcFaceI
+ );
+
+ tgtAddr[tgtFaceI].append(srcFaceI);
+ }
+ }
+
+
+ // transfer data to persistent storage
+ forAll(srcAddr, i)
+ {
+ scalar magSf = this->srcMagSf_[i];
+ srcAddress[i].transfer(srcAddr[i]);
+ srcWeights[i] = scalarList(1, magSf);
+ }
+ forAll(tgtAddr, i)
+ {
+ scalar magSf = this->tgtMagSf_[i];
+ tgtAddress[i].transfer(tgtAddr[i]);
+ tgtWeights[i] = scalarList(1, magSf);
+ }
+}
+
+
+// ************************************************************************* //
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.H b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.H
new file mode 100644
index 0000000000000000000000000000000000000000..0daf26662d27a98473c4b2e233ed405f65fa2aba
--- /dev/null
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.H
@@ -0,0 +1,158 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::mapNearestAMI
+
+Description
+ Nearest-mapping Arbitrary Mesh Interface (AMI) method
+
+SourceFiles
+ mapNearestAMI.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef mapNearestAMI_H
+#define mapNearestAMI_H
+
+#include "AMIMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class mapNearestAMI Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class SourcePatch, class TargetPatch>
+class mapNearestAMI
+:
+ public AMIMethod<SourcePatch, TargetPatch>
+{
+
+private:
+
+ // Private Member Functions
+
+ //- Disallow default bitwise copy construct
+ mapNearestAMI(const mapNearestAMI&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const mapNearestAMI&);
+
+ // Marching front
+
+ //- Find nearest target face for source face srcFaceI
+ void findNearestFace
+ (
+ const SourcePatch& srcPatch,
+ const TargetPatch& tgtPatch,
+ const label& srcFaceI,
+ label& tgtFaceI
+ ) const;
+
+ //- Determine next source-target face pair
+ void setNextNearestFaces
+ (
+ boolList& mapFlag,
+ label& startSeedI,
+ label& srcFaceI,
+ label& tgtFaceI
+ ) const;
+
+ //- Find mapped source face
+ label findMappedSrcFace
+ (
+ const label tgtFaceI,
+ const List<DynamicList<label> >& tgtToSrc
+ ) const;
+
+
+ // Evaluation
+
+ //- Area of intersection between source and target faces
+ scalar interArea
+ (
+ const label srcFaceI,
+ const label tgtFaceI
+ ) const;
+
+
+public:
+
+ //- Runtime type information
+ TypeName("mapNearestAMI");
+
+
+ // Constructors
+
+ //- Construct from components
+ mapNearestAMI
+ (
+ const SourcePatch& srcPatch,
+ const TargetPatch& tgtPatch,
+ const scalarField& srcMagSf,
+ const scalarField& tgtMagSf,
+ const faceAreaIntersect::triangulationMode& triMode,
+ const bool reverseTarget = false
+ );
+
+
+ //- Destructor
+ virtual ~mapNearestAMI();
+
+
+ // Member Functions
+
+ // Manipulation
+
+ //- Update addressing and weights
+ virtual void calculate
+ (
+ labelListList& srcAddress,
+ scalarListList& srcWeights,
+ labelListList& tgtAddress,
+ scalarListList& tgtWeights,
+ label srcFaceI = -1,
+ label tgtFaceI = -1
+ );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+# include "mapNearestAMI.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/meshTools/AMIInterpolation/AMIInterpolation/AMIPatchToPatchInterpolation.C b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIPatchToPatchInterpolation.C
new file mode 100644
index 0000000000000000000000000000000000000000..6cb5678aafd2bdbbdb77e97bf9c20938e3a0dd86
--- /dev/null
+++ b/src/meshTools/AMIInterpolation/AMIInterpolation/AMIPatchToPatchInterpolation.C
@@ -0,0 +1,44 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "AMIPatchToPatchInterpolation.H"
+#include "AMIMethod.H"
+#include "directAMI.H"
+#include "mapNearestAMI.H"
+#include "faceAreaWeightAMI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ makeAMIMethod(AMIPatchToPatchInterpolation);
+
+ makeAMIMethodType(AMIPatchToPatchInterpolation, directAMI);
+ makeAMIMethodType(AMIPatchToPatchInterpolation, mapNearestAMI);
+ makeAMIMethodType(AMIPatchToPatchInterpolation, faceAreaWeightAMI);
+}
+
+
+// ************************************************************************* //
diff --git a/src/meshTools/AMIInterpolation/patches/cyclic/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C b/src/meshTools/AMIInterpolation/patches/cyclic/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C
index 1d37199cd3b6087ba6e97f08d3ff3b45e7e99b56..0458dad9dee1edb6e0c8180dba59b4f1162a618c 100644
--- a/src/meshTools/AMIInterpolation/patches/cyclic/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C
+++ b/src/meshTools/AMIInterpolation/patches/cyclic/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C
@@ -283,6 +283,7 @@ void Foam::cyclicAMIPolyPatch::resetAMI() const
nbrPatch0,
surfPtr(),
faceAreaIntersect::tmMesh,
+ AMIPatchToPatchInterpolation::imFaceAreaWeight,
AMIReverse_
)
);
diff --git a/src/meshTools/Make/files b/src/meshTools/Make/files
index d46d774d887cb783668c650821021693af74a76b..f11e6fc4712a237cb2bf553d83850e0d9b79d10e 100644
--- a/src/meshTools/Make/files
+++ b/src/meshTools/Make/files
@@ -170,6 +170,7 @@ twoDPointCorrector/twoDPointCorrector.C
AMI=AMIInterpolation
$(AMI)/AMIInterpolation/AMIInterpolationName.C
+$(AMI)/AMIInterpolation/AMIPatchToPatchInterpolation.C
$(AMI)/faceAreaIntersect/faceAreaIntersect.C
$(AMI)/GAMG/interfaces/cyclicAMIGAMGInterface/cyclicAMIGAMGInterface.C
$(AMI)/GAMG/interfaceFields/cyclicAMIGAMGInterfaceField/cyclicAMIGAMGInterfaceField.C
diff --git a/src/meshTools/indexedOctree/treeDataPoint.C b/src/meshTools/indexedOctree/treeDataPoint.C
index 7c3a831159c86e18c112de53c8db1f92fbefa0ea..1d7f3d133ec2710e5ad64f1adb12723456db8cd0 100644
--- a/src/meshTools/indexedOctree/treeDataPoint.C
+++ b/src/meshTools/indexedOctree/treeDataPoint.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -160,7 +160,11 @@ void Foam::treeDataPoint::findNearest
) const
{
// Best so far
- scalar nearestDistSqr = magSqr(linePoint - nearestPoint);
+ scalar nearestDistSqr = GREAT;
+ if (minIndex >= 0)
+ {
+ nearestDistSqr = magSqr(linePoint - nearestPoint);
+ }
forAll(indices, i)
{
diff --git a/src/meshTools/indexedOctree/treeDataTriSurface.C b/src/meshTools/indexedOctree/treeDataTriSurface.C
index fb1a4645f690e96b62256556b420e69800697f40..0aa567f68b0ce38a80b6bf4576dc25f25accaa6d 100644
--- a/src/meshTools/indexedOctree/treeDataTriSurface.C
+++ b/src/meshTools/indexedOctree/treeDataTriSurface.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -404,11 +404,50 @@ void Foam::treeDataTriSurface::findNearest
point& nearestPoint
) const
{
- notImplemented
- (
- "treeDataTriSurface::findNearest(const labelUList&"
- ", const linePointRef&, treeBoundBox&, label&, point&, point&) const"
- );
+ // Best so far
+ scalar nearestDistSqr = VGREAT;
+ if (minIndex >= 0)
+ {
+ nearestDistSqr = magSqr(linePoint - nearestPoint);
+ }
+
+ const pointField& points = surface_.points();
+
+ forAll(indices, i)
+ {
+ label index = indices[i];
+ const triSurface::FaceType& f = surface_[index];
+
+ triPointRef tri(f.tri(points));
+
+ pointHit lineInfo(point::max);
+ pointHit pHit = tri.nearestPoint(ln, lineInfo);
+
+ scalar distSqr = sqr(pHit.distance());
+
+ if (distSqr < nearestDistSqr)
+ {
+ nearestDistSqr = distSqr;
+ minIndex = index;
+ linePoint = lineInfo.rawPoint();
+ nearestPoint = pHit.rawPoint();
+
+ {
+ point& minPt = tightest.min();
+ minPt = min(ln.start(), ln.end());
+ minPt.x() -= pHit.distance();
+ minPt.y() -= pHit.distance();
+ minPt.z() -= pHit.distance();
+ }
+ {
+ point& maxPt = tightest.max();
+ maxPt = max(ln.start(), ln.end());
+ maxPt.x() += pHit.distance();
+ maxPt.y() += pHit.distance();
+ maxPt.z() += pHit.distance();
+ }
+ }
+ }
}
diff --git a/src/meshTools/mappedPatches/mappedPolyPatch/mappedPatchBase.C b/src/meshTools/mappedPatches/mappedPolyPatch/mappedPatchBase.C
index 304501870181cf88709362df778bd3f7d7c7f714..4834c1a81a0e12e5da3f9ef7209139fc53548a48 100644
--- a/src/meshTools/mappedPatches/mappedPolyPatch/mappedPatchBase.C
+++ b/src/meshTools/mappedPatches/mappedPolyPatch/mappedPatchBase.C
@@ -849,6 +849,7 @@ void Foam::mappedPatchBase::calcAMI() const
samplePolyPatch(), // nbrPatch0,
surfPtr(),
faceAreaIntersect::tmMesh,
+ AMIPatchToPatchInterpolation::imFaceAreaWeight,
AMIReverse_
)
);
diff --git a/src/meshTools/regionCoupled/patches/regionCoupledPolyPatch/regionCoupledBase.C b/src/meshTools/regionCoupled/patches/regionCoupledPolyPatch/regionCoupledBase.C
index de7215e3da639686cda4088b2a538f12e00d58a7..5ee79489341cb939fcc59e83794a38b8fd7dc0c8 100644
--- a/src/meshTools/regionCoupled/patches/regionCoupledPolyPatch/regionCoupledBase.C
+++ b/src/meshTools/regionCoupled/patches/regionCoupledPolyPatch/regionCoupledBase.C
@@ -91,6 +91,7 @@ void Foam::regionCoupledBase::resetAMI() const
nbrPatch0,
surfPtr(),
faceAreaIntersect::tmMesh,
+ AMIPatchToPatchInterpolation::imFaceAreaWeight,
AMIReverse_
)
);
diff --git a/src/meshTools/sets/faceZoneSources/searchableSurfaceToFaceZone/searchableSurfaceToFaceZone.C b/src/meshTools/sets/faceZoneSources/searchableSurfaceToFaceZone/searchableSurfaceToFaceZone.C
index 039e14555f8515fa01346f9cc2a643d563606bab..c9ce92f25b040bb31023d61714a77715b7677577 100644
--- a/src/meshTools/sets/faceZoneSources/searchableSurfaceToFaceZone/searchableSurfaceToFaceZone.C
+++ b/src/meshTools/sets/faceZoneSources/searchableSurfaceToFaceZone/searchableSurfaceToFaceZone.C
@@ -28,6 +28,7 @@ License
#include "faceZoneSet.H"
#include "searchableSurface.H"
#include "syncTools.H"
+#include "Time.H"
#include "addToRunTimeSelectionTable.H"
@@ -69,7 +70,15 @@ Foam::searchableSurfaceToFaceZone::searchableSurfaceToFaceZone
searchableSurface::New
(
word(dict.lookup("surface")),
- mesh.objectRegistry::db(),
+ IOobject
+ (
+ dict.lookupOrDefault("name", mesh.objectRegistry::db().name()),
+ mesh.time().constant(),
+ "triSurface",
+ mesh.objectRegistry::db(),
+ IOobject::MUST_READ,
+ IOobject::NO_WRITE
+ ),
dict
)
)
diff --git a/src/regionModels/regionModel/regionModel/regionModel.C b/src/regionModels/regionModel/regionModel/regionModel.C
index 0590e4e418ce0b2fe9784e0b867ff8088b127236..125f48c63495a0b37538966f65182f5def1e28e5 100644
--- a/src/regionModels/regionModel/regionModel/regionModel.C
+++ b/src/regionModels/regionModel/regionModel/regionModel.C
@@ -221,6 +221,7 @@ Foam::regionModels::regionModel::interRegionAMI
p,
nbrP,
faceAreaIntersect::tmMesh,
+ AMIPatchToPatchInterpolation::imFaceAreaWeight,
flip
)
);
@@ -261,6 +262,7 @@ Foam::regionModels::regionModel::interRegionAMI
p,
nbrP,
faceAreaIntersect::tmMesh,
+ AMIPatchToPatchInterpolation::imFaceAreaWeight,
flip
)
);
@@ -318,7 +320,7 @@ Foam::label Foam::regionModels::regionModel::nbrCoupledPatchID
(
"Foam::label Foam::regionModels::regionModel::nbrCoupledPatchID"
"("
- "const regionModel& , "
+ "const regionModel&, "
"const label"
") const"
)
diff --git a/src/sampling/Make/files b/src/sampling/Make/files
index 309cf09d7ab39e79b984b4c0e7e8cce10f873b09..6d5b7ad19ead95add8c5b765617b1d571e7709d0 100644
--- a/src/sampling/Make/files
+++ b/src/sampling/Make/files
@@ -60,10 +60,14 @@ $(meshToMesh)/calculateMeshToMeshAddressing.C
$(meshToMesh)/calculateMeshToMeshWeights.C
meshToMeshNew = meshToMeshInterpolation/meshToMeshNew
-$(meshToMeshNew)/calcDirect.C
-$(meshToMeshNew)/calcMapNearest.C
-$(meshToMeshNew)/calcCellVolumeWeight.C
$(meshToMeshNew)/meshToMeshNew.C
$(meshToMeshNew)/meshToMeshNewParallelOps.C
+meshToMeshNewMethods = meshToMeshInterpolation/meshToMeshNew/calcMethod
+$(meshToMeshNewMethods)/meshToMeshMethod/meshToMeshMethod.C
+$(meshToMeshNewMethods)/meshToMeshMethod/meshToMeshMethodNew.C
+$(meshToMeshNewMethods)/cellVolumeWeight/cellVolumeWeightMethod.C
+$(meshToMeshNewMethods)/direct/directMethod.C
+$(meshToMeshNewMethods)/mapNearest/mapNearestMethod.C
+
LIB = $(FOAM_LIBBIN)/libsampling
diff --git a/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcDirect.C b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcDirect.C
deleted file mode 100644
index 8e6a17511ef590af5f8eddfa339a73994a146109..0000000000000000000000000000000000000000
--- a/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcDirect.C
+++ /dev/null
@@ -1,159 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2012-2013 OpenFOAM Foundation
- \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "meshToMeshNew.H"
-
-// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
-
-void Foam::meshToMeshNew::calcDirect
-(
- const polyMesh& src,
- const polyMesh& tgt,
- const label srcSeedI,
- const label tgtSeedI
-)
-{
- // store a list of src cells already mapped
- boolList srcSeedFlag(src.nCells(), true);
- labelList srcTgtSeed(src.nCells(), -1);
-
- List<DynamicList<label> > srcToTgt(src.nCells());
- List<DynamicList<label> > tgtToSrc(tgt.nCells());
-
- DynamicList<label> srcSeeds;
-
- const scalarField& srcVc = src.cellVolumes();
- const scalarField& tgtVc = tgt.cellVolumes();
-
- label srcCellI = srcSeedI;
- label tgtCellI = tgtSeedI;
-
- do
- {
- // store src/tgt cell pair
- srcToTgt[srcCellI].append(tgtCellI);
- tgtToSrc[tgtCellI].append(srcCellI);
-
- // mark source cell srcSeedI as matched
- srcSeedFlag[srcCellI] = false;
-
- // accumulate intersection volume
- V_ += srcVc[srcCellI];
-
- // find new source seed cell
- appendToDirectSeeds
- (
- src,
- tgt,
- srcSeedFlag,
- srcTgtSeed,
- srcSeeds,
- srcCellI,
- tgtCellI
- );
- }
- while (srcCellI >= 0);
-
- // transfer addressing into persistent storage
- forAll(srcToTgtCellAddr_, i)
- {
- scalar v = srcVc[i];
- srcToTgtCellAddr_[i].transfer(srcToTgt[i]);
- srcToTgtCellWght_[i] = scalarList(srcToTgtCellAddr_[i].size(), v);
- }
-
- forAll(tgtToSrcCellAddr_, i)
- {
- scalar v = tgtVc[i];
- tgtToSrcCellAddr_[i].transfer(tgtToSrc[i]);
- tgtToSrcCellWght_[i] = scalarList(tgtToSrcCellAddr_[i].size(), v);
- }
-}
-
-
-void Foam::meshToMeshNew::appendToDirectSeeds
-(
- const polyMesh& src,
- const polyMesh& tgt,
- boolList& mapFlag,
- labelList& srcTgtSeed,
- DynamicList<label>& srcSeeds,
- label& srcSeedI,
- label& tgtSeedI
-) const
-{
- const labelList& srcNbr = src.cellCells()[srcSeedI];
- const labelList& tgtNbr = tgt.cellCells()[tgtSeedI];
-
- const vectorField& srcCentre = src.cellCentres();
-
- forAll(srcNbr, i)
- {
- label srcI = srcNbr[i];
-
- if (mapFlag[srcI] && (srcTgtSeed[srcI] == -1))
- {
- // source cell srcI not yet mapped
-
- // identfy if target cell exists for source cell srcI
- bool found = false;
- forAll(tgtNbr, j)
- {
- label tgtI = tgtNbr[j];
-
- if (tgt.pointInCell(srcCentre[srcI], tgtI))
- {
- // new match - append to lists
- found = true;
-
- srcTgtSeed[srcI] = tgtI;
- srcSeeds.append(srcI);
-
- break;
- }
- }
-
- if (!found)
- {
- // no match available for source cell srcI
- mapFlag[srcI] = false;
- }
- }
- }
-
- if (srcSeeds.size())
- {
- srcSeedI = srcSeeds.remove();
- tgtSeedI = srcTgtSeed[srcSeedI];
- }
- else
- {
- srcSeedI = -1;
- tgtSeedI = -1;
- }
-}
-
-
-// ************************************************************************* //
diff --git a/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMapNearest.C b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMapNearest.C
deleted file mode 100644
index 0276ee560495f8679eda7c7aa7b7710a41030aa3..0000000000000000000000000000000000000000
--- a/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMapNearest.C
+++ /dev/null
@@ -1,265 +0,0 @@
-/*---------------------------------------------------------------------------*\
- ========= |
- \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
- \\ / O peration |
- \\ / A nd | Copyright (C) 2012-2013 OpenFOAM Foundation
- \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "meshToMeshNew.H"
-#include "ListOps.H"
-
-// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
-
-void Foam::meshToMeshNew::calcMapNearest
-(
- const polyMesh& src,
- const polyMesh& tgt,
- const label srcSeedI,
- const label tgtSeedI,
- const labelList& srcCellIDs,
- boolList& mapFlag,
- label& startSeedI
-)
-{
- List<DynamicList<label> > srcToTgt(src.nCells());
- List<DynamicList<label> > tgtToSrc(tgt.nCells());
-
- const scalarField& srcVc = src.cellVolumes();
- const scalarField& tgtVc = tgt.cellVolumes();
-
- label srcCellI = srcSeedI;
- label tgtCellI = tgtSeedI;
-
- do
- {
- // find nearest tgt cell
- findNearestCell(src, tgt, srcCellI, tgtCellI);
-
- // store src/tgt cell pair
- srcToTgt[srcCellI].append(tgtCellI);
- tgtToSrc[tgtCellI].append(srcCellI);
-
- // mark source cell srcCellI and tgtCellI as matched
- mapFlag[srcCellI] = false;
-
- // accumulate intersection volume
- V_ += srcVc[srcCellI];
-
- // find new source cell
- setNextNearestCells
- (
- startSeedI,
- srcCellI,
- tgtCellI,
- mapFlag,
- src,
- tgt,
- srcCellIDs
- );
- }
- while (srcCellI >= 0);
-
-
- // for the case of multiple source cells per target cell, select the
- // nearest source cell only and discard the others
- const vectorField& srcCc = src.cellCentres();
- const vectorField& tgtCc = tgt.cellCentres();
-
- forAll(tgtToSrc, targetCellI)
- {
- if (tgtToSrc[targetCellI].size() > 1)
- {
- const vector& tgtC = tgtCc[tgtCellI];
-
- DynamicList<label>& srcCells = tgtToSrc[targetCellI];
-
- label srcCellI = srcCells[0];
- scalar d = magSqr(tgtC - srcCc[srcCellI]);
-
- for (label i = 1; i < srcCells.size(); i++)
- {
- label srcI = srcCells[i];
- scalar dNew = magSqr(tgtC - srcCc[srcI]);
- if (dNew < d)
- {
- d = dNew;
- srcCellI = srcI;
- }
- }
-
- srcCells.clear();
- srcCells.append(srcCellI);
- }
- }
-
- // If there are more target cells than source cells, some target cells
- // might not yet be mapped
- forAll(tgtToSrc, tgtCellI)
- {
- if (tgtToSrc[tgtCellI].empty())
- {
- label srcCellI = findMappedSrcCell(tgt, tgtCellI, tgtToSrc);
-
- findNearestCell(tgt, src, tgtCellI, srcCellI);
-
- tgtToSrc[tgtCellI].append(srcCellI);
- }
- }
-
- // transfer addressing into persistent storage
- forAll(srcToTgtCellAddr_, i)
- {
- scalar v = srcVc[i];
- srcToTgtCellWght_[i] = scalarList(srcToTgt[i].size(), v);
- srcToTgtCellAddr_[i].transfer(srcToTgt[i]);
- }
-
- forAll(tgtToSrcCellAddr_, i)
- {
- scalar v = tgtVc[i];
- tgtToSrcCellWght_[i] = scalarList(tgtToSrc[i].size(), v);
- tgtToSrcCellAddr_[i].transfer(tgtToSrc[i]);
- }
-}
-
-
-void Foam::meshToMeshNew::findNearestCell
-(
- const polyMesh& src,
- const polyMesh& tgt,
- const label srcCellI,
- label& tgtCellI
-)
-{
- const vectorField& srcC = src.cellCentres();
- const vectorField& tgtC = tgt.cellCentres();
-
- const vector& srcP = srcC[srcCellI];
-
- DynamicList<label> tgtCells(10);
- tgtCells.append(tgtCellI);
-
- DynamicList<label> visitedCells(10);
-
- scalar d = GREAT;
-
- do
- {
- label tgtI = tgtCells.remove();
- visitedCells.append(tgtI);
-
- scalar dTest = magSqr(tgtC[tgtI] - srcP);
- if (dTest < d)
- {
- tgtCellI = tgtI;
- d = dTest;
- appendNbrCells(tgtCellI, tgt, visitedCells, tgtCells);
- }
-
- } while (tgtCells.size() > 0);
-}
-
-
-void Foam::meshToMeshNew::setNextNearestCells
-(
- label& startSeedI,
- label& srcCellI,
- label& tgtCellI,
- boolList& mapFlag,
- const polyMesh& src,
- const polyMesh& tgt,
- const labelList& srcCellIDs
-)
-{
- const labelList& srcNbr = src.cellCells()[srcCellI];
-
- srcCellI = -1;
- forAll(srcNbr, i)
- {
- label cellI = srcNbr[i];
- if (mapFlag[cellI])
- {
- srcCellI = cellI;
- startSeedI = cellI + 1;
-
- return;
- }
- }
-
- (void)findInitialSeeds
- (
- src,
- tgt,
- srcCellIDs,
- mapFlag,
- startSeedI,
- srcCellI,
- tgtCellI
- );
-}
-
-
-Foam::label Foam::meshToMeshNew::findMappedSrcCell
-(
- const polyMesh& tgt,
- const label tgtCellI,
- const List<DynamicList<label> >& tgtToSrc
-) const
-{
- DynamicList<label> testCells(10);
- DynamicList<label> visitedCells(10);
-
- testCells.append(tgtCellI);
-
- do
- {
- // search target tgtCellI neighbours for match with source cell
- label tgtI = testCells.remove();
-
- if (findIndex(visitedCells, tgtI) == -1)
- {
- visitedCells.append(tgtI);
-
- if (tgtToSrc[tgtI].size())
- {
- return tgtToSrc[tgtI][0];
- }
- else
- {
- const labelList& nbrCells = tgt.cellCells()[tgtI];
-
- forAll(nbrCells, i)
- {
- if (findIndex(visitedCells, nbrCells[i]) == -1)
- {
- testCells.append(nbrCells[i]);
- }
- }
- }
- }
- } while (testCells.size());
-
- // did not find any match - should not be possible to get here!
- return -1;
-}
-
-
-// ************************************************************************* //
diff --git a/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcCellVolumeWeight.C b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/cellVolumeWeight/cellVolumeWeightMethod.C
similarity index 53%
rename from src/sampling/meshToMeshInterpolation/meshToMeshNew/calcCellVolumeWeight.C
rename to src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/cellVolumeWeight/cellVolumeWeightMethod.C
index 9cd8e417f1c497ee21f638702c6d7b56edcd3553..67694afd9d02bc8d46da253b01650cfc706834fd 100644
--- a/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcCellVolumeWeight.C
+++ b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/cellVolumeWeight/cellVolumeWeightMethod.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2012-2013 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -23,15 +23,79 @@ License
\*---------------------------------------------------------------------------*/
-#include "meshToMeshNew.H"
-#include "tetOverlapVolume.H"
+#include "cellVolumeWeightMethod.H"
+#include "indexedOctree.H"
+#include "treeDataCell.H"
+#include "addToRunTimeSelectionTable.H"
-// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-void Foam::meshToMeshNew::calcCellVolumeWeight
+namespace Foam
+{
+ defineTypeNameAndDebug(cellVolumeWeightMethod, 0);
+ addToRunTimeSelectionTable
+ (
+ meshToMeshMethod,
+ cellVolumeWeightMethod,
+ components
+ );
+}
+
+// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
+
+bool Foam::cellVolumeWeightMethod::findInitialSeeds
(
- const polyMesh& src,
- const polyMesh& tgt,
+ const labelList& srcCellIDs,
+ const boolList& mapFlag,
+ const label startSeedI,
+ label& srcSeedI,
+ label& tgtSeedI
+) const
+{
+ const cellList& srcCells = src_.cells();
+ const faceList& srcFaces = src_.faces();
+ const pointField& srcPts = src_.points();
+
+ for (label i = startSeedI; i < srcCellIDs.size(); i++)
+ {
+ label srcI = srcCellIDs[i];
+
+ if (mapFlag[srcI])
+ {
+ const pointField
+ pts(srcCells[srcI].points(srcFaces, srcPts).xfer());
+
+ forAll(pts, ptI)
+ {
+ const point& pt = pts[ptI];
+ label tgtI = tgt_.cellTree().findInside(pt);
+
+ if (tgtI != -1 && intersect(srcI, tgtI))
+ {
+ srcSeedI = srcI;
+ tgtSeedI = tgtI;
+
+ return true;
+ }
+ }
+ }
+ }
+
+ if (debug)
+ {
+ Pout<< "could not find starting seed" << endl;
+ }
+
+ return false;
+}
+
+
+void Foam::cellVolumeWeightMethod::calculateAddressing
+(
+ labelListList& srcToTgtCellAddr,
+ scalarListList& srcToTgtCellWght,
+ labelListList& tgtToSrcCellAddr,
+ scalarListList& tgtToSrcCellWght,
const label srcSeedI,
const label tgtSeedI,
const labelList& srcCellIDs,
@@ -42,11 +106,11 @@ void Foam::meshToMeshNew::calcCellVolumeWeight
label srcCellI = srcSeedI;
label tgtCellI = tgtSeedI;
- List<DynamicList<label> > srcToTgtAddr(src.nCells());
- List<DynamicList<scalar> > srcToTgtWght(src.nCells());
+ List<DynamicList<label> > srcToTgtAddr(src_.nCells());
+ List<DynamicList<scalar> > srcToTgtWght(src_.nCells());
- List<DynamicList<label> > tgtToSrcAddr(tgt.nCells());
- List<DynamicList<scalar> > tgtToSrcWght(tgt.nCells());
+ List<DynamicList<label> > tgtToSrcAddr(tgt_.nCells());
+ List<DynamicList<scalar> > tgtToSrcWght(tgt_.nCells());
// list of tgt cell neighbour cells
DynamicList<label> nbrTgtCells(10);
@@ -55,10 +119,10 @@ void Foam::meshToMeshNew::calcCellVolumeWeight
DynamicList<label> visitedTgtCells(10);
// list to keep track of tgt cells used to seed src cells
- labelList seedCells(src.nCells(), -1);
+ labelList seedCells(src_.nCells(), -1);
seedCells[srcCellI] = tgtCellI;
- const scalarField& srcVol = src.cellVolumes();
+ const scalarField& srcVol = src_.cellVolumes();
do
{
@@ -67,14 +131,14 @@ void Foam::meshToMeshNew::calcCellVolumeWeight
// append initial target cell and neighbours
nbrTgtCells.append(tgtCellI);
- appendNbrCells(tgtCellI, tgt, visitedTgtCells, nbrTgtCells);
+ appendNbrCells(tgtCellI, tgt_, visitedTgtCells, nbrTgtCells);
do
{
tgtCellI = nbrTgtCells.remove();
visitedTgtCells.append(tgtCellI);
- scalar vol = interVol(src, tgt, srcCellI, tgtCellI);
+ scalar vol = interVol(srcCellI, tgtCellI);
// accumulate addressing and weights for valid intersection
if (vol/srcVol[srcCellI] > tolerance_)
@@ -86,7 +150,7 @@ void Foam::meshToMeshNew::calcCellVolumeWeight
tgtToSrcAddr[tgtCellI].append(srcCellI);
tgtToSrcWght[tgtCellI].append(vol);
- appendNbrCells(tgtCellI, tgt, visitedTgtCells, nbrTgtCells);
+ appendNbrCells(tgtCellI, tgt_, visitedTgtCells, nbrTgtCells);
// accumulate intersection volume
V_ += vol;
@@ -102,8 +166,6 @@ void Foam::meshToMeshNew::calcCellVolumeWeight
startSeedI,
srcCellI,
tgtCellI,
- src,
- tgt,
srcCellIDs,
mapFlag,
visitedTgtCells,
@@ -113,34 +175,32 @@ void Foam::meshToMeshNew::calcCellVolumeWeight
while (srcCellI != -1);
// transfer addressing into persistent storage
- forAll(srcToTgtCellAddr_, i)
+ forAll(srcToTgtCellAddr, i)
{
- srcToTgtCellAddr_[i].transfer(srcToTgtAddr[i]);
- srcToTgtCellWght_[i].transfer(srcToTgtWght[i]);
+ srcToTgtCellAddr[i].transfer(srcToTgtAddr[i]);
+ srcToTgtCellWght[i].transfer(srcToTgtWght[i]);
}
- forAll(tgtToSrcCellAddr_, i)
+ forAll(tgtToSrcCellAddr, i)
{
- tgtToSrcCellAddr_[i].transfer(tgtToSrcAddr[i]);
- tgtToSrcCellWght_[i].transfer(tgtToSrcWght[i]);
+ tgtToSrcCellAddr[i].transfer(tgtToSrcAddr[i]);
+ tgtToSrcCellWght[i].transfer(tgtToSrcWght[i]);
}
}
-void Foam::meshToMeshNew::setNextCells
+void Foam::cellVolumeWeightMethod::setNextCells
(
label& startSeedI,
label& srcCellI,
label& tgtCellI,
- const polyMesh& src,
- const polyMesh& tgt,
const labelList& srcCellIDs,
const boolList& mapFlag,
const DynamicList<label>& visitedCells,
labelList& seedCells
) const
{
- const labelList& srcNbrCells = src.cellCells()[srcCellI];
+ const labelList& srcNbrCells = src_.cellCells()[srcCellI];
// set possible seeds for later use by querying all src cell neighbours
// with all visited target cells
@@ -155,7 +215,7 @@ void Foam::meshToMeshNew::setNextCells
{
label cellT = visitedCells[j];
- if (intersect(src, tgt, cellS, cellT))
+ if (intersect(cellS, cellT))
{
seedCells[cellS] = cellT;
@@ -211,8 +271,6 @@ void Foam::meshToMeshNew::setNextCells
bool restart =
findInitialSeeds
(
- src,
- tgt,
srcCellIDs,
mapFlag,
startSeedI,
@@ -233,68 +291,91 @@ void Foam::meshToMeshNew::setNextCells
}
-bool Foam::meshToMeshNew::intersect
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::cellVolumeWeightMethod::cellVolumeWeightMethod
(
const polyMesh& src,
- const polyMesh& tgt,
- const label srcCellI,
- const label tgtCellI
-) const
-{
- scalar threshold = tolerance_*src.cellVolumes()[srcCellI];
+ const polyMesh& tgt
+)
+:
+ meshToMeshMethod(src, tgt)
+{}
- tetOverlapVolume overlapEngine;
- treeBoundBox bbTgtCell
- (
- pointField
- (
- tgt.points(),
- tgt.cellPoints()[tgtCellI]
- )
- );
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::cellVolumeWeightMethod::~cellVolumeWeightMethod()
+{}
- return overlapEngine.cellCellOverlapMinDecomp
- (
- src,
- srcCellI,
- tgt,
- tgtCellI,
- bbTgtCell,
- threshold
- );
-}
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
-Foam::scalar Foam::meshToMeshNew::interVol
+void Foam::cellVolumeWeightMethod::calculate
(
- const polyMesh& src,
- const polyMesh& tgt,
- const label srcCellI,
- const label tgtCellI
-) const
+ labelListList& srcToTgtAddr,
+ scalarListList& srcToTgtWght,
+ labelListList& tgtToSrcAddr,
+ scalarListList& tgtToSrcWght
+)
{
- tetOverlapVolume overlapEngine;
-
- treeBoundBox bbTgtCell
+ bool ok = initialise
(
- pointField
- (
- tgt.points(),
- tgt.cellPoints()[tgtCellI]
- )
+ srcToTgtAddr,
+ srcToTgtWght,
+ tgtToSrcAddr,
+ tgtToSrcWght
);
- scalar vol = overlapEngine.cellCellOverlapVolumeMinDecomp
- (
- src,
- srcCellI,
- tgt,
- tgtCellI,
- bbTgtCell
- );
+ if (!ok)
+ {
+ return;
+ }
+
+ // (potentially) participating source mesh cells
+ const labelList srcCellIDs(maskCells());
- return vol;
+ // list to keep track of whether src cell can be mapped
+ boolList mapFlag(src_.nCells(), false);
+ UIndirectList<bool>(mapFlag, srcCellIDs) = true;
+
+ // find initial point in tgt mesh
+ label srcSeedI = -1;
+ label tgtSeedI = -1;
+ label startSeedI = 0;
+
+ bool startWalk =
+ findInitialSeeds
+ (
+ srcCellIDs,
+ mapFlag,
+ startSeedI,
+ srcSeedI,
+ tgtSeedI
+ );
+
+ if (startWalk)
+ {
+ calculateAddressing
+ (
+ srcToTgtAddr,
+ srcToTgtWght,
+ tgtToSrcAddr,
+ tgtToSrcWght,
+ srcSeedI,
+ tgtSeedI,
+ srcCellIDs,
+ mapFlag,
+ startSeedI
+ );
+ }
+ else
+ {
+ // if meshes are collocated, after inflating the source mesh bounding
+ // box tgt mesh cells may be transferred, but may still not overlap
+ // with the source mesh
+ return;
+ }
}
diff --git a/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/cellVolumeWeight/cellVolumeWeightMethod.H b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/cellVolumeWeight/cellVolumeWeightMethod.H
new file mode 100644
index 0000000000000000000000000000000000000000..ad38d76246579cfb311ed4a5e4afdf368fa20e77
--- /dev/null
+++ b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/cellVolumeWeight/cellVolumeWeightMethod.H
@@ -0,0 +1,138 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::cellVolumeWeightMethod
+
+Description
+ Cell-volume-weighted mesh-to-mesh interpolation class
+
+ Volume conservative.
+
+SourceFiles
+ cellVolumeWeightMethod.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef cellVolumeWeightMethod_H
+#define cellVolumeWeightMethod_H
+
+#include "meshToMeshMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class cellVolumeWeightMethod Declaration
+\*---------------------------------------------------------------------------*/
+
+class cellVolumeWeightMethod
+:
+ public meshToMeshMethod
+{
+protected:
+
+ // Protected Member Functions
+
+ //- Find indices of overlapping cells in src and tgt meshes - returns
+ // true if found a matching pair
+ bool findInitialSeeds
+ (
+ const labelList& srcCellIDs,
+ const boolList& mapFlag,
+ const label startSeedI,
+ label& srcSeedI,
+ label& tgtSeedI
+ ) const;
+
+ //- Calculate the mesh-to-mesh addressing and weights
+ void calculateAddressing
+ (
+ labelListList& srcToTgtCellAddr,
+ scalarListList& srcToTgtCellWght,
+ labelListList& tgtToSrcCellAddr,
+ scalarListList& tgtToSrcCellWght,
+ const label srcSeedI,
+ const label tgtSeedI,
+ const labelList& srcCellIDs,
+ boolList& mapFlag,
+ label& startSeedI
+ );
+
+ //- Set the next cells in the advancing front algorithm
+ void setNextCells
+ (
+ label& startSeedI,
+ label& srcCellI,
+ label& tgtCellI,
+ const labelList& srcCellIDs,
+ const boolList& mapFlag,
+ const DynamicList<label>& visitedCells,
+ labelList& seedCells
+ ) const;
+
+ //- Disallow default bitwise copy construct
+ cellVolumeWeightMethod(const cellVolumeWeightMethod&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const cellVolumeWeightMethod&);
+
+
+public:
+
+ //- Run-time type information
+ TypeName("cellVolumeWeight");
+
+ //- Construct from source and target meshes
+ cellVolumeWeightMethod(const polyMesh& src, const polyMesh& tgt);
+
+ //- Destructor
+ virtual ~cellVolumeWeightMethod();
+
+
+ // Member Functions
+
+ // Evaluate
+
+ //- Calculate addressing and weights
+ virtual void calculate
+ (
+ labelListList& srcToTgtAddr,
+ scalarListList& srcToTgtWght,
+ labelListList& tgtToTgtAddr,
+ scalarListList& tgtToTgtWght
+ );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/direct/directMethod.C b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/direct/directMethod.C
new file mode 100644
index 0000000000000000000000000000000000000000..9f3674700cd958090647ac3be3c7cfe7847304c0
--- /dev/null
+++ b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/direct/directMethod.C
@@ -0,0 +1,303 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "directMethod.H"
+#include "indexedOctree.H"
+#include "treeDataCell.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(directMethod, 0);
+ addToRunTimeSelectionTable(meshToMeshMethod, directMethod, components);
+}
+
+// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
+
+bool Foam::directMethod::findInitialSeeds
+(
+ const labelList& srcCellIDs,
+ const boolList& mapFlag,
+ const label startSeedI,
+ label& srcSeedI,
+ label& tgtSeedI
+) const
+{
+ const cellList& srcCells = src_.cells();
+ const faceList& srcFaces = src_.faces();
+ const pointField& srcPts = src_.points();
+
+ for (label i = startSeedI; i < srcCellIDs.size(); i++)
+ {
+ label srcI = srcCellIDs[i];
+
+ if (mapFlag[srcI])
+ {
+ const pointField
+ pts(srcCells[srcI].points(srcFaces, srcPts).xfer());
+
+ forAll(pts, ptI)
+ {
+ const point& pt = pts[ptI];
+ label tgtI = tgt_.cellTree().findInside(pt);
+
+ if (tgtI != -1 && intersect(srcI, tgtI))
+ {
+ srcSeedI = srcI;
+ tgtSeedI = tgtI;
+
+ return true;
+ }
+ }
+ }
+ }
+
+ if (debug)
+ {
+ Pout<< "could not find starting seed" << endl;
+ }
+
+ return false;
+}
+
+
+void Foam::directMethod::calculateAddressing
+(
+ labelListList& srcToTgtCellAddr,
+ scalarListList& srcToTgtCellWght,
+ labelListList& tgtToSrcCellAddr,
+ scalarListList& tgtToSrcCellWght,
+ const label srcSeedI,
+ const label tgtSeedI,
+ const labelList& srcCellIDs, // not used
+ boolList& mapFlag,
+ label& startSeedI
+)
+{
+ // store a list of src cells already mapped
+ labelList srcTgtSeed(src_.nCells(), -1);
+
+ List<DynamicList<label> > srcToTgt(src_.nCells());
+ List<DynamicList<label> > tgtToSrc(tgt_.nCells());
+
+ DynamicList<label> srcSeeds(10);
+
+ const scalarField& srcVc = src_.cellVolumes();
+ const scalarField& tgtVc = tgt_.cellVolumes();
+
+ label srcCellI = srcSeedI;
+ label tgtCellI = tgtSeedI;
+
+ do
+ {
+ // store src/tgt cell pair
+ srcToTgt[srcCellI].append(tgtCellI);
+ tgtToSrc[tgtCellI].append(srcCellI);
+
+ // mark source cell srcSeedI as matched
+ mapFlag[srcCellI] = false;
+
+ // accumulate intersection volume
+ V_ += srcVc[srcCellI];
+
+ // find new source seed cell
+ appendToDirectSeeds
+ (
+ mapFlag,
+ srcTgtSeed,
+ srcSeeds,
+ srcCellI,
+ tgtCellI
+ );
+ }
+ while (srcCellI >= 0);
+
+ // transfer addressing into persistent storage
+ forAll(srcToTgtCellAddr, i)
+ {
+ scalar v = srcVc[i];
+ srcToTgtCellAddr[i].transfer(srcToTgt[i]);
+ srcToTgtCellWght[i] = scalarList(1, v);
+ }
+
+ forAll(tgtToSrcCellAddr, i)
+ {
+ scalar v = tgtVc[i];
+ tgtToSrcCellAddr[i].transfer(tgtToSrc[i]);
+ tgtToSrcCellWght[i] = scalarList(1, v);
+ }
+}
+
+
+void Foam::directMethod::appendToDirectSeeds
+(
+ boolList& mapFlag,
+ labelList& srcTgtSeed,
+ DynamicList<label>& srcSeeds,
+ label& srcSeedI,
+ label& tgtSeedI
+) const
+{
+ const labelList& srcNbr = src_.cellCells()[srcSeedI];
+ const labelList& tgtNbr = tgt_.cellCells()[tgtSeedI];
+
+ const vectorField& srcCentre = src_.cellCentres();
+
+ forAll(srcNbr, i)
+ {
+ label srcI = srcNbr[i];
+
+ if (mapFlag[srcI] && (srcTgtSeed[srcI] == -1))
+ {
+ // source cell srcI not yet mapped
+
+ // identfy if target cell exists for source cell srcI
+ bool found = false;
+ forAll(tgtNbr, j)
+ {
+ label tgtI = tgtNbr[j];
+
+ if (tgt_.pointInCell(srcCentre[srcI], tgtI))
+ {
+ // new match - append to lists
+ found = true;
+
+ srcTgtSeed[srcI] = tgtI;
+ srcSeeds.append(srcI);
+
+ break;
+ }
+ }
+
+ if (!found)
+ {
+ // no match available for source cell srcI
+ mapFlag[srcI] = false;
+ }
+ }
+ }
+
+ if (srcSeeds.size())
+ {
+ srcSeedI = srcSeeds.remove();
+ tgtSeedI = srcTgtSeed[srcSeedI];
+ }
+ else
+ {
+ srcSeedI = -1;
+ tgtSeedI = -1;
+ }
+}
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::directMethod::directMethod
+(
+ const polyMesh& src,
+ const polyMesh& tgt
+)
+:
+ meshToMeshMethod(src, tgt)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::directMethod::~directMethod()
+{}
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+void Foam::directMethod::calculate
+(
+ labelListList& srcToTgtAddr,
+ scalarListList& srcToTgtWght,
+ labelListList& tgtToSrcAddr,
+ scalarListList& tgtToSrcWght
+)
+{
+ bool ok = initialise
+ (
+ srcToTgtAddr,
+ srcToTgtWght,
+ tgtToSrcAddr,
+ tgtToSrcWght
+ );
+
+ if (!ok)
+ {
+ return;
+ }
+
+ // (potentially) participating source mesh cells
+ const labelList srcCellIDs(maskCells());
+
+ // list to keep track of whether src cell can be mapped
+ boolList mapFlag(src_.nCells(), false);
+ UIndirectList<bool>(mapFlag, srcCellIDs) = true;
+
+ // find initial point in tgt mesh
+ label srcSeedI = -1;
+ label tgtSeedI = -1;
+ label startSeedI = 0;
+
+ bool startWalk =
+ findInitialSeeds
+ (
+ srcCellIDs,
+ mapFlag,
+ startSeedI,
+ srcSeedI,
+ tgtSeedI
+ );
+
+ if (startWalk)
+ {
+ calculateAddressing
+ (
+ srcToTgtAddr,
+ srcToTgtWght,
+ tgtToSrcAddr,
+ tgtToSrcWght,
+ srcSeedI,
+ tgtSeedI,
+ srcCellIDs,
+ mapFlag,
+ startSeedI
+ );
+ }
+ else
+ {
+ // if meshes are collocated, after inflating the source mesh bounding
+ // box tgt mesh cells may be transferred, but may still not overlap
+ // with the source mesh
+ return;
+ }
+}
+
+
+// ************************************************************************* //
diff --git a/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/direct/directMethod.H b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/direct/directMethod.H
new file mode 100644
index 0000000000000000000000000000000000000000..1c9f489a83c3fc8f19fb6db0757bb0c00d361c31
--- /dev/null
+++ b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/direct/directMethod.H
@@ -0,0 +1,136 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::directMethod
+
+Description
+ Direct (one-to-one cell correspondence) mesh-to-mesh interpolation class
+
+ Volume conservative.
+
+SourceFiles
+ directMethod.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef directMethod_H
+#define directMethod_H
+
+#include "meshToMeshMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class directMethod Declaration
+\*---------------------------------------------------------------------------*/
+
+class directMethod
+:
+ public meshToMeshMethod
+{
+protected:
+
+ // Protected Member Functions
+
+ //- Find indices of overlapping cells in src and tgt meshes - returns
+ // true if found a matching pair
+ bool findInitialSeeds
+ (
+ const labelList& srcCellIDs,
+ const boolList& mapFlag,
+ const label startSeedI,
+ label& srcSeedI,
+ label& tgtSeedI
+ ) const;
+
+ //- Calculate the mesh-to-mesh addressing and weights
+ void calculateAddressing
+ (
+ labelListList& srcToTgtCellAddr,
+ scalarListList& srcToTgtCellWght,
+ labelListList& tgtToSrcCellAddr,
+ scalarListList& tgtToSrcCellWght,
+ const label srcSeedI,
+ const label tgtSeedI,
+ const labelList& srcCellIDs,
+ boolList& mapFlag,
+ label& startSeedI
+ );
+
+ //- Append to list of src mesh seed indices
+ void appendToDirectSeeds
+ (
+ boolList& mapFlag,
+ labelList& srcTgtSeed,
+ DynamicList<label>& srcSeeds,
+ label& srcSeedI,
+ label& tgtSeedI
+ ) const;
+
+ //- Disallow default bitwise copy construct
+ directMethod(const directMethod&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const directMethod&);
+
+
+public:
+
+ //- Run-time type information
+ TypeName("direct");
+
+ //- Construct from source and target meshes
+ directMethod(const polyMesh& src, const polyMesh& tgt);
+
+ //- Destructor
+ virtual ~directMethod();
+
+
+ // Member Functions
+
+ // Evaluate
+
+ //- Calculate addressing and weights
+ virtual void calculate
+ (
+ labelListList& srcToTgtAddr,
+ scalarListList& srcToTgtWght,
+ labelListList& tgtToTgtAddr,
+ scalarListList& tgtToTgtWght
+ );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/mapNearest/mapNearestMethod.C b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/mapNearest/mapNearestMethod.C
new file mode 100644
index 0000000000000000000000000000000000000000..9327c3abb193b2c50e180056dba313d477f1e494
--- /dev/null
+++ b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/mapNearest/mapNearestMethod.C
@@ -0,0 +1,424 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "mapNearestMethod.H"
+#include "pointIndexHit.H"
+#include "indexedOctree.H"
+#include "treeDataCell.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(mapNearestMethod, 0);
+ addToRunTimeSelectionTable(meshToMeshMethod, mapNearestMethod, components);
+}
+
+// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
+
+bool Foam::mapNearestMethod::findInitialSeeds
+(
+ const labelList& srcCellIDs,
+ const boolList& mapFlag,
+ const label startSeedI,
+ label& srcSeedI,
+ label& tgtSeedI
+) const
+{
+ const cellList& srcCells = src_.cells();
+ const faceList& srcFaces = src_.faces();
+ const pointField& srcPts = src_.points();
+
+ for (label i = startSeedI; i < srcCellIDs.size(); i++)
+ {
+ label srcI = srcCellIDs[i];
+
+ if (mapFlag[srcI])
+ {
+ const pointField
+ pts(srcCells[srcI].points(srcFaces, srcPts).xfer());
+
+ const point& pt = pts[0];
+ pointIndexHit hit = tgt_.cellTree().findNearest(pt, GREAT);
+
+ if (hit.hit())
+ {
+ srcSeedI = srcI;
+ tgtSeedI = hit.index();
+
+ return true;
+ }
+ else
+ {
+ FatalErrorIn
+ (
+ "bool Foam::mapNearestMethod::findInitialSeeds"
+ "("
+ "const labelList&, "
+ "const boolList&, "
+ "const label, "
+ "label&, "
+ "label&"
+ ") const"
+ )
+ << "Unable to find nearest target cell"
+ << " for source cell " << srcI
+ << " with centre "
+ << srcCells[srcI].centre(srcPts, srcFaces)
+ << abort(FatalError);
+ }
+
+ break;
+ }
+ }
+
+ if (debug)
+ {
+ Pout<< "could not find starting seed" << endl;
+ }
+
+ return false;
+}
+
+
+void Foam::mapNearestMethod::calculateAddressing
+(
+ labelListList& srcToTgtCellAddr,
+ scalarListList& srcToTgtCellWght,
+ labelListList& tgtToSrcCellAddr,
+ scalarListList& tgtToSrcCellWght,
+ const label srcSeedI,
+ const label tgtSeedI,
+ const labelList& srcCellIDs,
+ boolList& mapFlag,
+ label& startSeedI
+)
+{
+ List<DynamicList<label> > srcToTgt(src_.nCells());
+ List<DynamicList<label> > tgtToSrc(tgt_.nCells());
+
+ const scalarField& srcVc = src_.cellVolumes();
+ const scalarField& tgtVc = tgt_.cellVolumes();
+
+ label srcCellI = srcSeedI;
+ label tgtCellI = tgtSeedI;
+
+ do
+ {
+ // find nearest tgt cell
+ findNearestCell(src_, tgt_, srcCellI, tgtCellI);
+
+ // store src/tgt cell pair
+ srcToTgt[srcCellI].append(tgtCellI);
+ tgtToSrc[tgtCellI].append(srcCellI);
+
+ // mark source cell srcCellI and tgtCellI as matched
+ mapFlag[srcCellI] = false;
+
+ // accumulate intersection volume
+ V_ += srcVc[srcCellI];
+
+ // find new source cell
+ setNextNearestCells
+ (
+ startSeedI,
+ srcCellI,
+ tgtCellI,
+ mapFlag,
+ srcCellIDs
+ );
+ }
+ while (srcCellI >= 0);
+
+
+ // for the case of multiple source cells per target cell, select the
+ // nearest source cell only and discard the others
+ const vectorField& srcCc = src_.cellCentres();
+ const vectorField& tgtCc = tgt_.cellCentres();
+
+ forAll(tgtToSrc, targetCellI)
+ {
+ if (tgtToSrc[targetCellI].size() > 1)
+ {
+ const vector& tgtC = tgtCc[tgtCellI];
+
+ DynamicList<label>& srcCells = tgtToSrc[targetCellI];
+
+ label srcCellI = srcCells[0];
+ scalar d = magSqr(tgtC - srcCc[srcCellI]);
+
+ for (label i = 1; i < srcCells.size(); i++)
+ {
+ label srcI = srcCells[i];
+ scalar dNew = magSqr(tgtC - srcCc[srcI]);
+ if (dNew < d)
+ {
+ d = dNew;
+ srcCellI = srcI;
+ }
+ }
+
+ srcCells.clear();
+ srcCells.append(srcCellI);
+ }
+ }
+
+ // If there are more target cells than source cells, some target cells
+ // might not yet be mapped
+ forAll(tgtToSrc, tgtCellI)
+ {
+ if (tgtToSrc[tgtCellI].empty())
+ {
+ label srcCellI = findMappedSrcCell(tgtCellI, tgtToSrc);
+
+ findNearestCell(tgt_, src_, tgtCellI, srcCellI);
+
+ tgtToSrc[tgtCellI].append(srcCellI);
+ }
+ }
+
+ // transfer addressing into persistent storage
+ forAll(srcToTgtCellAddr, i)
+ {
+ scalar v = srcVc[i];
+ srcToTgtCellAddr[i].transfer(srcToTgt[i]);
+ srcToTgtCellWght[i] = scalarList(1, v);
+ }
+
+ forAll(tgtToSrcCellAddr, i)
+ {
+ scalar v = tgtVc[i];
+ tgtToSrcCellAddr[i].transfer(tgtToSrc[i]);
+ tgtToSrcCellWght[i] = scalarList(1, v);
+ }
+}
+
+
+void Foam::mapNearestMethod::findNearestCell
+(
+ const polyMesh& mesh1,
+ const polyMesh& mesh2,
+ const label cell1,
+ label& cell2
+) const
+{
+ const vectorField& Cc1 = mesh1.cellCentres();
+ const vectorField& Cc2 = mesh2.cellCentres();
+
+ const vector& p1 = Cc1[cell1];
+
+ DynamicList<label> cells2(10);
+ cells2.append(cell2);
+
+ DynamicList<label> visitedCells(10);
+
+ scalar d = GREAT;
+
+ do
+ {
+ label c2 = cells2.remove();
+ visitedCells.append(c2);
+
+ scalar dTest = magSqr(Cc2[c2] - p1);
+ if (dTest < d)
+ {
+ cell2 = c2;
+ d = dTest;
+ appendNbrCells(cell2, mesh2, visitedCells, cells2);
+ }
+
+ } while (cells2.size() > 0);
+}
+
+
+void Foam::mapNearestMethod::setNextNearestCells
+(
+ label& startSeedI,
+ label& srcCellI,
+ label& tgtCellI,
+ boolList& mapFlag,
+ const labelList& srcCellIDs
+) const
+{
+ const labelList& srcNbr = src_.cellCells()[srcCellI];
+
+ srcCellI = -1;
+ forAll(srcNbr, i)
+ {
+ label cellI = srcNbr[i];
+ if (mapFlag[cellI])
+ {
+ srcCellI = cellI;
+ startSeedI = cellI + 1;
+
+ return;
+ }
+ }
+
+ (void)findInitialSeeds
+ (
+ srcCellIDs,
+ mapFlag,
+ startSeedI,
+ srcCellI,
+ tgtCellI
+ );
+}
+
+
+Foam::label Foam::mapNearestMethod::findMappedSrcCell
+(
+ const label tgtCellI,
+ const List<DynamicList<label> >& tgtToSrc
+) const
+{
+ DynamicList<label> testCells(10);
+ DynamicList<label> visitedCells(10);
+
+ testCells.append(tgtCellI);
+
+ do
+ {
+ // search target tgtCellI neighbours for match with source cell
+ label tgtI = testCells.remove();
+
+ if (findIndex(visitedCells, tgtI) == -1)
+ {
+ visitedCells.append(tgtI);
+
+ if (tgtToSrc[tgtI].size())
+ {
+ return tgtToSrc[tgtI][0];
+ }
+ else
+ {
+ const labelList& nbrCells = tgt_.cellCells()[tgtI];
+
+ forAll(nbrCells, i)
+ {
+ if (findIndex(visitedCells, nbrCells[i]) == -1)
+ {
+ testCells.append(nbrCells[i]);
+ }
+ }
+ }
+ }
+ } while (testCells.size());
+
+ // did not find any match - should not be possible to get here!
+ return -1;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::mapNearestMethod::mapNearestMethod
+(
+ const polyMesh& src,
+ const polyMesh& tgt
+)
+:
+ meshToMeshMethod(src, tgt)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::mapNearestMethod::~mapNearestMethod()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+void Foam::mapNearestMethod::calculate
+(
+ labelListList& srcToTgtAddr,
+ scalarListList& srcToTgtWght,
+ labelListList& tgtToSrcAddr,
+ scalarListList& tgtToSrcWght
+)
+{
+ bool ok = initialise
+ (
+ srcToTgtAddr,
+ srcToTgtWght,
+ tgtToSrcAddr,
+ tgtToSrcWght
+ );
+
+ if (!ok)
+ {
+ return;
+ }
+
+ // (potentially) participating source mesh cells
+ const labelList srcCellIDs(maskCells());
+
+ // list to keep track of whether src cell can be mapped
+ boolList mapFlag(src_.nCells(), false);
+ UIndirectList<bool>(mapFlag, srcCellIDs) = true;
+
+ // find initial point in tgt mesh
+ label srcSeedI = -1;
+ label tgtSeedI = -1;
+ label startSeedI = 0;
+
+ bool startWalk =
+ findInitialSeeds
+ (
+ srcCellIDs,
+ mapFlag,
+ startSeedI,
+ srcSeedI,
+ tgtSeedI
+ );
+
+ if (startWalk)
+ {
+ calculateAddressing
+ (
+ srcToTgtAddr,
+ srcToTgtWght,
+ tgtToSrcAddr,
+ tgtToSrcWght,
+ srcSeedI,
+ tgtSeedI,
+ srcCellIDs,
+ mapFlag,
+ startSeedI
+ );
+ }
+ else
+ {
+ // if meshes are collocated, after inflating the source mesh bounding
+ // box tgt mesh cells may be transferred, but may still not overlap
+ // with the source mesh
+ return;
+ }
+}
+
+
+// ************************************************************************* //
diff --git a/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/mapNearest/mapNearestMethod.H b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/mapNearest/mapNearestMethod.H
new file mode 100644
index 0000000000000000000000000000000000000000..813c46683ff5671e27ac3fd0bd3f6e44914e70a1
--- /dev/null
+++ b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/mapNearest/mapNearestMethod.H
@@ -0,0 +1,152 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::mapNearestMethod
+
+Description
+ Map nearest mesh-to-mesh interpolation class
+
+ Not volume conservative.
+
+SourceFiles
+ mapNearestMethod.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef mapNearestMethod_H
+#define mapNearestMethod_H
+
+#include "meshToMeshMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class mapNearestMethod Declaration
+\*---------------------------------------------------------------------------*/
+
+class mapNearestMethod
+:
+ public meshToMeshMethod
+{
+protected:
+
+ // Protected Member Functions
+
+ //- Find indices of overlapping cells in src and tgt meshes - returns
+ // true if found a matching pair
+ bool findInitialSeeds
+ (
+ const labelList& srcCellIDs,
+ const boolList& mapFlag,
+ const label startSeedI,
+ label& srcSeedI,
+ label& tgtSeedI
+ ) const;
+
+ //- Calculate the mesh-to-mesh addressing and weights
+ void calculateAddressing
+ (
+ labelListList& srcToTgtCellAddr,
+ scalarListList& srcToTgtCellWght,
+ labelListList& tgtToSrcCellAddr,
+ scalarListList& tgtToSrcCellWght,
+ const label srcSeedI,
+ const label tgtSeedI,
+ const labelList& srcCellIDs,
+ boolList& mapFlag,
+ label& startSeedI
+ );
+
+ //- Find the nearest cell on mesh2 for cell1 on mesh1
+ void findNearestCell
+ (
+ const polyMesh& mesh1,
+ const polyMesh& mesh2,
+ const label cell1,
+ label& cell2
+ ) const;
+
+ //- Set the next cells for the marching front algorithm
+ void setNextNearestCells
+ (
+ label& startSeedI,
+ label& srcCellI,
+ label& tgtCellI,
+ boolList& mapFlag,
+ const labelList& srcCellIDs
+ ) const;
+
+ //- Find a source cell mapped to target cell tgtCellI
+ label findMappedSrcCell
+ (
+ const label tgtCellI,
+ const List<DynamicList<label> >& tgtToSrc
+ ) const;
+
+ //- Disallow default bitwise copy construct
+ mapNearestMethod(const mapNearestMethod&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const mapNearestMethod&);
+
+
+public:
+
+ //- Run-time type information
+ TypeName("mapNearest");
+
+ //- Construct from source and target meshes
+ mapNearestMethod(const polyMesh& src, const polyMesh& tgt);
+
+ //- Destructor
+ virtual ~mapNearestMethod();
+
+
+ // Member Functions
+
+ // Evaluate
+
+ //- Calculate addressing and weights
+ virtual void calculate
+ (
+ labelListList& srcToTgtAddr,
+ scalarListList& srcToTgtWght,
+ labelListList& tgtToTgtAddr,
+ scalarListList& tgtToTgtWght
+ );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethod.C b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethod.C
new file mode 100644
index 0000000000000000000000000000000000000000..29e066dd542937b69a6872ca663bd0c9e1622c18
--- /dev/null
+++ b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethod.C
@@ -0,0 +1,274 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "meshToMeshMethod.H"
+#include "tetOverlapVolume.H"
+#include "OFstream.H"
+#include "Time.H"
+#include "treeBoundBox.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(meshToMeshMethod, 0);
+ defineRunTimeSelectionTable(meshToMeshMethod, components);
+}
+
+Foam::scalar Foam::meshToMeshMethod::tolerance_ = 1e-6;
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+Foam::labelList Foam::meshToMeshMethod::maskCells() const
+{
+ boundBox intersectBb
+ (
+ max(src_.bounds().min(), tgt_.bounds().min()),
+ min(src_.bounds().max(), tgt_.bounds().max())
+ );
+
+ intersectBb.inflate(0.01);
+
+ const cellList& srcCells = src_.cells();
+ const faceList& srcFaces = src_.faces();
+ const pointField& srcPts = src_.points();
+
+ DynamicList<label> cells(src_.nCells());
+ forAll(srcCells, srcI)
+ {
+ boundBox cellBb(srcCells[srcI].points(srcFaces, srcPts), false);
+ if (intersectBb.overlaps(cellBb))
+ {
+ cells.append(srcI);
+ }
+ }
+
+ if (debug)
+ {
+ Pout<< "participating source mesh cells: " << cells.size() << endl;
+ }
+
+ return cells;
+}
+
+
+bool Foam::meshToMeshMethod::intersect
+(
+ const label srcCellI,
+ const label tgtCellI
+) const
+{
+ scalar threshold = tolerance_*src_.cellVolumes()[srcCellI];
+
+ tetOverlapVolume overlapEngine;
+
+ treeBoundBox bbTgtCell
+ (
+ pointField
+ (
+ tgt_.points(),
+ tgt_.cellPoints()[tgtCellI]
+ )
+ );
+
+ return overlapEngine.cellCellOverlapMinDecomp
+ (
+ src_,
+ srcCellI,
+ tgt_,
+ tgtCellI,
+ bbTgtCell,
+ threshold
+ );
+}
+
+
+Foam::scalar Foam::meshToMeshMethod::interVol
+(
+ const label srcCellI,
+ const label tgtCellI
+) const
+{
+ tetOverlapVolume overlapEngine;
+
+ treeBoundBox bbTgtCell
+ (
+ pointField
+ (
+ tgt_.points(),
+ tgt_.cellPoints()[tgtCellI]
+ )
+ );
+
+ scalar vol = overlapEngine.cellCellOverlapVolumeMinDecomp
+ (
+ src_,
+ srcCellI,
+ tgt_,
+ tgtCellI,
+ bbTgtCell
+ );
+
+ return vol;
+}
+
+
+void Foam::meshToMeshMethod::appendNbrCells
+(
+ const label cellI,
+ const polyMesh& mesh,
+ const DynamicList<label>& visitedCells,
+ DynamicList<label>& nbrCellIDs
+) const
+{
+ const labelList& nbrCells = mesh.cellCells()[cellI];
+
+ // filter out cells already visited from cell neighbours
+ forAll(nbrCells, i)
+ {
+ label nbrCellI = nbrCells[i];
+
+ if
+ (
+ (findIndex(visitedCells, nbrCellI) == -1)
+ && (findIndex(nbrCellIDs, nbrCellI) == -1)
+ )
+ {
+ nbrCellIDs.append(nbrCellI);
+ }
+ }
+}
+
+
+bool Foam::meshToMeshMethod::initialise
+(
+ labelListList& srcToTgtAddr,
+ scalarListList& srcToTgtWght,
+ labelListList& tgtToSrcAddr,
+ scalarListList& tgtToSrcWght
+) const
+{
+ srcToTgtAddr.setSize(src_.nCells());
+ srcToTgtWght.setSize(src_.nCells());
+ tgtToSrcAddr.setSize(tgt_.nCells());
+ tgtToSrcWght.setSize(tgt_.nCells());
+
+ if (!src_.nCells())
+ {
+ return false;
+ }
+ else if (!tgt_.nCells())
+ {
+ if (debug)
+ {
+ Pout<< "mesh interpolation: hhave " << src_.nCells() << " source "
+ << " cells but no target cells" << endl;
+ }
+
+ return false;
+ }
+
+ return true;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::meshToMeshMethod::meshToMeshMethod
+(
+ const polyMesh& src,
+ const polyMesh& tgt
+)
+:
+ src_(src),
+ tgt_(tgt),
+ V_(0.0)
+{
+ if (!src_.nCells() || !tgt_.nCells())
+ {
+ if (debug)
+ {
+ Pout<< "mesh interpolation: cells not on processor: Source cells = "
+ << src_.nCells() << ", target cells = " << tgt_.nCells()
+ << endl;
+ }
+ }
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::meshToMeshMethod::~meshToMeshMethod()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+void Foam::meshToMeshMethod::writeConnectivity
+(
+ const polyMesh& mesh1,
+ const polyMesh& mesh2,
+ const labelListList& mesh1ToMesh2Addr
+) const
+{
+ Pout<< "Source size = " << mesh1.nCells() << endl;
+ Pout<< "Target size = " << mesh2.nCells() << endl;
+
+ word fName("addressing_" + mesh1.name() + "_to_" + mesh2.name());
+
+ if (Pstream::parRun())
+ {
+ fName = fName + "_proc" + Foam::name(Pstream::myProcNo());
+ }
+
+ OFstream os(src_.time().path()/fName + ".obj");
+
+ label vertI = 0;
+ forAll(mesh1ToMesh2Addr, i)
+ {
+ const labelList& addr = mesh1ToMesh2Addr[i];
+ forAll(addr, j)
+ {
+ label cellI = addr[j];
+ const vector& c0 = mesh1.cellCentres()[i];
+
+ const cell& c = mesh2.cells()[cellI];
+ const pointField pts(c.points(mesh2.faces(), mesh2.points()));
+ forAll(pts, j)
+ {
+ const point& p = pts[j];
+ os << "v " << p.x() << ' ' << p.y() << ' ' << p.z() << nl;
+ vertI++;
+ os << "v " << c0.x() << ' ' << c0.y() << ' ' << c0.z()
+ << nl;
+ vertI++;
+ os << "l " << vertI - 1 << ' ' << vertI << nl;
+ }
+ }
+ }
+}
+
+
+// ************************************************************************* //
diff --git a/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethod.H b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethod.H
new file mode 100644
index 0000000000000000000000000000000000000000..3338fa5349569629ec72684a0304702c7323661e
--- /dev/null
+++ b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethod.H
@@ -0,0 +1,188 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::meshToMeshMethod
+
+Description
+ Base class for mesh-to-mesh calculation methods
+
+SourceFiles
+ meshToMeshMethod.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef meshToMeshMethod_H
+#define meshToMeshMethod_H
+
+#include "polyMesh.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class meshToMeshMethod Declaration
+\*---------------------------------------------------------------------------*/
+
+class meshToMeshMethod
+{
+
+protected:
+
+ // Protected data
+
+ //- Reference to the source mesh
+ const polyMesh& src_;
+
+ //- Reference to the target mesh
+ const polyMesh& tgt_;
+
+ //- Cell total volume in overlap region [m3]
+ scalar V_;
+
+ //- Tolerance used in volume overlap calculations
+ static scalar tolerance_;
+
+
+ // Protected Member Functions
+
+ //- Return src cell IDs for the overlap region
+ labelList maskCells() const;
+
+ //- Return the true if cells intersect
+ bool intersect(const label srcCellI, const label tgtCellI) const;
+
+ //- Return the intersection volume between two cells
+ scalar interVol(const label srcCellI, const label tgtCellI) const;
+
+ //- Append target cell neihgbour cells to cellIDs list
+ void appendNbrCells
+ (
+ const label tgtCellI,
+ const polyMesh& mesh,
+ const DynamicList<label>& visitedTgtCells,
+ DynamicList<label>& nbrTgtCellIDs
+ ) const;
+
+ bool initialise
+ (
+ labelListList& srcToTgtAddr,
+ scalarListList& srcToTgtWght,
+ labelListList& tgtToTgtAddr,
+ scalarListList& tgtToTgtWght
+ ) const;
+
+ //- Disallow default bitwise copy construct
+ meshToMeshMethod(const meshToMeshMethod&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const meshToMeshMethod&);
+
+
+public:
+
+ //- Run-time type information
+ TypeName("meshToMeshMethod");
+
+ //- Declare runtime constructor selection table
+ declareRunTimeSelectionTable
+ (
+ autoPtr,
+ meshToMeshMethod,
+ components,
+ (
+ const polyMesh& src,
+ const polyMesh& tgt
+ ),
+ (src, tgt)
+ );
+
+ //- Construct from source and target meshes
+ meshToMeshMethod(const polyMesh& src, const polyMesh& tgt);
+
+ //- Selector
+ static autoPtr<meshToMeshMethod> New
+ (
+ const word& methodName,
+ const polyMesh& src,
+ const polyMesh& tgt
+ );
+
+
+ //- Destructor
+ virtual ~meshToMeshMethod();
+
+
+ // Member Functions
+
+ // Evaluate
+
+ //- Calculate addressing and weights
+ virtual void calculate
+ (
+ labelListList& srcToTgtAddr,
+ scalarListList& srcToTgtWght,
+ labelListList& tgtToTgtAddr,
+ scalarListList& tgtToTgtWght
+ ) = 0;
+
+
+ // Access
+
+ //- Return const access to the source mesh
+ inline const polyMesh& src() const;
+
+ //- Return const access to the target mesh
+ inline const polyMesh& tgt() const;
+
+ //- Return const access to the overlap volume
+ inline scalar V() const;
+
+
+ // Check
+
+ //- Write the connectivity (debugging)
+ void writeConnectivity
+ (
+ const polyMesh& mesh1,
+ const polyMesh& mesh2,
+ const labelListList& mesh1ToMesh2Addr
+ ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "meshToMeshMethodI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethodI.H b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethodI.H
new file mode 100644
index 0000000000000000000000000000000000000000..a1f98a9123c905768b0b0f6c75b2e6e49268ecf7
--- /dev/null
+++ b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethodI.H
@@ -0,0 +1,44 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+const Foam::polyMesh& Foam::meshToMeshMethod::src() const
+{
+ return src_;
+}
+
+
+const Foam::polyMesh& Foam::meshToMeshMethod::tgt() const
+{
+ return tgt_;
+}
+
+
+Foam::scalar Foam::meshToMeshMethod::V() const
+{
+ return V_;
+}
+
+
+// ************************************************************************* //
diff --git a/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethodNew.C b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethodNew.C
new file mode 100644
index 0000000000000000000000000000000000000000..befb3e6debdd70685e9c1ac6e0698deb3e25a34d
--- /dev/null
+++ b/src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethodNew.C
@@ -0,0 +1,65 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "meshToMeshMethod.H"
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+Foam::autoPtr<Foam::meshToMeshMethod> Foam::meshToMeshMethod::New
+(
+ const word& methodName,
+ const polyMesh& src,
+ const polyMesh& tgt
+)
+{
+ if (debug)
+ {
+ Info<< "Selecting AMIMethod " << methodName << endl;
+ }
+
+ componentsConstructorTable::iterator cstrIter =
+ componentsConstructorTablePtr_->find(methodName);
+
+ if (cstrIter == componentsConstructorTablePtr_->end())
+ {
+ FatalErrorIn
+ (
+ "Foam::autoPtr<Foam::meshToMeshMethod> Foam::meshToMeshMethod::New"
+ "("
+ "const word&, "
+ "const polyMesh&, "
+ "const polyMesh&"
+ ")"
+ ) << "Unknown meshToMesh type "
+ << methodName << nl << nl
+ << "Valid meshToMesh types are:" << nl
+ << componentsConstructorTablePtr_->sortedToc() << exit(FatalError);
+ }
+
+ return autoPtr<meshToMeshMethod>(cstrIter()(src, tgt));
+}
+
+
+// ************************************************************************* //
diff --git a/src/sampling/meshToMeshInterpolation/meshToMeshNew/meshToMeshNew.C b/src/sampling/meshToMeshInterpolation/meshToMeshNew/meshToMeshNew.C
index 4add751668f357bbc37cbd493660df5ba5820500..27f0bd592c34f61382989fa2fd517f7b77f61f08 100644
--- a/src/sampling/meshToMeshInterpolation/meshToMeshNew/meshToMeshNew.C
+++ b/src/sampling/meshToMeshInterpolation/meshToMeshNew/meshToMeshNew.C
@@ -24,14 +24,9 @@ License
\*---------------------------------------------------------------------------*/
#include "meshToMeshNew.H"
-#include "OFstream.H"
#include "Time.H"
#include "globalIndex.H"
-#include "mergePoints.H"
-#include "treeBoundBox.H"
-#include "indexedOctree.H"
-#include "treeDataCell.H"
-#include "ListOps.H"
+#include "meshToMeshMethod.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -55,55 +50,9 @@ namespace Foam
meshToMeshNew::interpolationMethodNames_;
}
-Foam::scalar Foam::meshToMeshNew::tolerance_ = 1e-6;
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
-void Foam::meshToMeshNew::writeConnectivity
-(
- const polyMesh& src,
- const polyMesh& tgt,
- const labelListList& srcToTargetAddr
-) const
-{
- Pout<< "Source size = " << src.nCells() << endl;
- Pout<< "Target size = " << tgt.nCells() << endl;
-
- word fName("addressing_" + src.name() + "_to_" + tgt.name());
-
- if (Pstream::parRun())
- {
- fName = fName + "_proc" + Foam::name(Pstream::myProcNo());
- }
-
- OFstream os(src.time().path()/fName + ".obj");
-
- label vertI = 0;
- forAll(srcToTargetAddr, i)
- {
- const labelList& tgtAddress = srcToTargetAddr[i];
- forAll(tgtAddress, j)
- {
- label tgtI = tgtAddress[j];
- const vector& c0 = src.cellCentres()[i];
-
- const cell& c = tgt.cells()[tgtI];
- const pointField pts(c.points(tgt.faces(), tgt.points()));
- forAll(pts, j)
- {
- const point& p = pts[j];
- os << "v " << p.x() << ' ' << p.y() << ' ' << p.z() << nl;
- vertI++;
- os << "v " << c0.x() << ' ' << c0.y() << ' ' << c0.z()
- << nl;
- vertI++;
- os << "l " << vertI - 1 << ' ' << vertI << nl;
- }
- }
- }
-}
-
-
Foam::labelList Foam::meshToMeshNew::maskCells
(
const polyMesh& src,
@@ -141,82 +90,6 @@ Foam::labelList Foam::meshToMeshNew::maskCells
}
-bool Foam::meshToMeshNew::findInitialSeeds
-(
- const polyMesh& src,
- const polyMesh& tgt,
- const labelList& srcCellIDs,
- const boolList& mapFlag,
- const label startSeedI,
- label& srcSeedI,
- label& tgtSeedI
-) const
-{
- const cellList& srcCells = src.cells();
- const faceList& srcFaces = src.faces();
- const pointField& srcPts = src.points();
-
- for (label i = startSeedI; i < srcCellIDs.size(); i++)
- {
- label srcI = srcCellIDs[i];
-
- if (mapFlag[srcI])
- {
- const pointField
- pts(srcCells[srcI].points(srcFaces, srcPts).xfer());
-
- forAll(pts, ptI)
- {
- const point& pt = pts[ptI];
- label tgtI = tgt.cellTree().findInside(pt);
-
- if (tgtI != -1 && intersect(src, tgt, srcI, tgtI))
- {
- srcSeedI = srcI;
- tgtSeedI = tgtI;
-
- return true;
- }
- }
- }
- }
-
- if (debug)
- {
- Pout<< "could not find starting seed" << endl;
- }
-
- return false;
-}
-
-
-void Foam::meshToMeshNew::appendNbrCells
-(
- const label cellI,
- const polyMesh& mesh,
- const DynamicList<label>& visitedCells,
- DynamicList<label>& nbrCellIDs
-) const
-{
- const labelList& nbrCells = mesh.cellCells()[cellI];
-
- // filter out cells already visited from cell neighbours
- forAll(nbrCells, i)
- {
- label nbrCellI = nbrCells[i];
-
- if
- (
- (findIndex(visitedCells, nbrCellI) == -1)
- && (findIndex(nbrCellIDs, nbrCellI) == -1)
- )
- {
- nbrCellIDs.append(nbrCellI);
- }
- }
-}
-
-
void Foam::meshToMeshNew::normaliseWeights
(
const word& descriptor,
@@ -260,124 +133,29 @@ void Foam::meshToMeshNew::calcAddressing
const polyMesh& tgt
)
{
- srcToTgtCellAddr_.setSize(src.nCells());
- srcToTgtCellWght_.setSize(src.nCells());
-
- tgtToSrcCellAddr_.setSize(tgt.nCells());
- tgtToSrcCellWght_.setSize(tgt.nCells());
-
- if (!src.nCells() || !tgt.nCells())
- {
- if (debug)
- {
- Pout<< "mesh interpolation: cells not on processor: Source cells = "
- << src.nCells() << ", target cells = " << tgt.nCells()
- << endl;
- }
- }
-
- if (!src.nCells())
- {
- return;
- }
- else if (!tgt.nCells())
- {
- if (debug)
- {
- Pout<< "mesh interpolation: hhave " << src.nCells() << " source "
- << " cells but no target cells" << endl;
- }
-
- return;
- }
-
- // (potentially) participating source mesh cells
- const labelList srcCellIDs = maskCells(src, tgt);
-
- // list to keep track of whether src cell can be mapped
- boolList mapFlag(src.nCells(), false);
- UIndirectList<bool>(mapFlag, srcCellIDs) = true;
-
- // find initial point in tgt mesh
- label srcSeedI = -1;
- label tgtSeedI = -1;
- label startSeedI = 0;
-
- bool startWalk =
- findInitialSeeds
+ autoPtr<meshToMeshMethod> methodPtr
+ (
+ meshToMeshMethod::New
(
+ interpolationMethodNames_[method_],
src,
- tgt,
- srcCellIDs,
- mapFlag,
- startSeedI,
- srcSeedI,
- tgtSeedI
- );
-
- if (!startWalk)
- {
- // if meshes are collocated, after inflating the source mesh bounding
- // box tgt mesh cells may be transferred, but may still not overlap
- // with the source mesh
- return;
- }
-
+ tgt
+ )
+ );
- switch (method_)
- {
- case imDirect:
- {
- calcDirect(src, tgt, srcSeedI, tgtSeedI);
- break;
- }
- case imMapNearest:
- {
- calcMapNearest
- (
- src,
- tgt,
- srcSeedI,
- tgtSeedI,
- srcCellIDs,
- mapFlag,
- startSeedI
- );
- break;
- }
- case imCellVolumeWeight:
- {
- calcCellVolumeWeight
- (
- src,
- tgt,
- srcSeedI,
- tgtSeedI,
- srcCellIDs,
- mapFlag,
- startSeedI
- );
- break;
- }
- default:
- {
- FatalErrorIn
- (
- "void Foam::meshToMeshNew::calcAddressing"
- "("
- "const polyMesh&, "
- "const polyMesh&"
- ")"
- )
- << "Unknown interpolation method"
- << abort(FatalError);
- }
- }
+ methodPtr->calculate
+ (
+ srcToTgtCellAddr_,
+ srcToTgtCellWght_,
+ tgtToSrcCellAddr_,
+ tgtToSrcCellWght_
+ );
+ V_ = methodPtr->V();
if (debug)
{
- writeConnectivity(src, tgt, srcToTgtCellAddr_);
+ methodPtr->writeConnectivity(src, tgt, srcToTgtCellAddr_);
}
}
@@ -577,11 +355,59 @@ void Foam::meshToMeshNew::calculate()
}
+Foam::AMIPatchToPatchInterpolation::interpolationMethod
+Foam::meshToMeshNew::interpolationMethodAMI
+(
+ const interpolationMethod method
+) const
+{
+ switch (method_)
+ {
+ case imDirect:
+ {
+ return AMIPatchToPatchInterpolation::imDirect;
+ break;
+ }
+ case imMapNearest:
+ {
+ return AMIPatchToPatchInterpolation::imMapNearest;
+ break;
+ }
+ case imCellVolumeWeight:
+ {
+ return AMIPatchToPatchInterpolation::imFaceAreaWeight;
+ break;
+ }
+ default:
+ {
+ FatalErrorIn
+ (
+ "Foam::AMIPatchToPatchInterpolation::interpolationMethod"
+ "Foam::meshToMeshNew::interpolationMethodAMI"
+ "("
+ "const interpolationMethod method"
+ ") const"
+ )
+ << "Unhandled enumeration " << method_
+ << abort(FatalError);
+ }
+ }
+
+ return AMIPatchToPatchInterpolation::imDirect;
+}
+
+
const Foam::PtrList<Foam::AMIPatchToPatchInterpolation>&
Foam::meshToMeshNew::patchAMIs() const
{
if (patchAMIs_.empty())
{
+ const word amiMethod =
+ AMIPatchToPatchInterpolation::interpolationMethodToWord
+ (
+ interpolationMethodAMI(method_)
+ );
+
patchAMIs_.setSize(srcPatchID_.size());
forAll(srcPatchID_, i)
@@ -593,7 +419,9 @@ Foam::meshToMeshNew::patchAMIs() const
const polyPatch& tgtPP = tgtRegion_.boundaryMesh()[tgtPatchI];
Info<< "Creating AMI between source patch " << srcPP.name()
- << " and target patch " << tgtPP.name() << endl;
+ << " and target patch " << tgtPP.name()
+ << " using " << amiMethod
+ << endl;
Info<< incrIndent;
@@ -605,6 +433,7 @@ Foam::meshToMeshNew::patchAMIs() const
srcPP,
tgtPP,
faceAreaIntersect::tmMesh,
+ interpolationMethodAMI(method_),
true // flip target patch since patch normals are aligned
)
);
diff --git a/src/sampling/meshToMeshInterpolation/meshToMeshNew/meshToMeshNew.H b/src/sampling/meshToMeshInterpolation/meshToMeshNew/meshToMeshNew.H
index 8404ce7b9fa077ecc07ea02671d4673dfeb6981c..ab05a058ae98b7b1efab12e8cf215625bfed3c75 100644
--- a/src/sampling/meshToMeshInterpolation/meshToMeshNew/meshToMeshNew.H
+++ b/src/sampling/meshToMeshInterpolation/meshToMeshNew/meshToMeshNew.H
@@ -27,20 +27,14 @@ Class
Description
Class to calculate the cell-addressing between two overlapping meshes
- Three methods are currently available:
- - direct : 1-to-1 mapping between meshes
- - mapNearest : assign nearest cell values without interpolation
- - cellVolumeWeight : volume consistent mapping
-
- The \c direct and \c cellVolumeWeight options are volume conservative,
- whereas mapNearest is non-conservative.
+ Mapping is performed using a run-time selectable interpolation mothod
+SeeAlso
+ meshToMeshMethod
SourceFiles
- calcDirect.C
- calcMapNearest.C
- calcCellVolumeWeight.C
meshToMeshNew.C
+ meshToMeshNewParallelOps.C
meshToMeshNewTemplates.C
\*---------------------------------------------------------------------------*/
@@ -70,7 +64,7 @@ public:
// Public data types
- //- Enumeration specifying required accuracy
+ //- Enumeration specifying interpolation method
enum interpolationMethod
{
imDirect,
@@ -128,9 +122,6 @@ private:
//- Target map pointer - parallel running only
autoPtr<mapDistribute> tgtMapPtr_;
- //- Tolerance used in volume overlap calculations
- static scalar tolerance_;
-
// Private Member Functions
@@ -138,39 +129,9 @@ private:
template<class Type>
void add(UList<Type>& fld, const label offset) const;
- //- Write the connectivity (debugging)
- void writeConnectivity
- (
- const polyMesh& src,
- const polyMesh& tgt,
- const labelListList& srcToTargetAddr
- ) const;
-
//- Return src cell IDs for the overlap region
labelList maskCells(const polyMesh& src, const polyMesh& tgt) const;
- //- Find indices of overlapping cells in src and tgt meshes - returns
- // true if found a matching pair
- bool findInitialSeeds
- (
- const polyMesh& src,
- const polyMesh& tgt,
- const labelList& srcCellIDs,
- const boolList& mapFlag,
- const label startSeedI,
- label& srcSeedI,
- label& tgtSeedI
- ) const;
-
- //- Append target cell neihgbour cells to cellIDs list
- void appendNbrCells
- (
- const label tgtCellI,
- const polyMesh& tgt,
- const DynamicList<label>& visitedTgtCells,
- DynamicList<label>& nbrTgtCellIDs
- ) const;
-
//- Normalise the interpolation weights
void normaliseWeights
(
@@ -180,121 +141,6 @@ private:
scalarListList& wght
) const;
-
- // Direct (one-to-one) mapping
-
- //- Main driver routine for direct mapping
- void calcDirect
- (
- const polyMesh& src,
- const polyMesh& tgt,
- const label srcSeedI,
- const label tgtSeedI
- );
-
- //- Append to list of src mesh seed indices
- void appendToDirectSeeds
- (
- const polyMesh& src,
- const polyMesh& tgt,
- boolList& mapFlag,
- labelList& srcTgtSeed,
- DynamicList<label>& srcSeeds,
- label& srcSeedI,
- label& tgtSeedI
- ) const;
-
- // Nearest (non-conformal) mapping
-
- //- Main driver routine for nearest-mapping routine
- void calcMapNearest
- (
- const polyMesh& src,
- const polyMesh& tgt,
- const label srcSeedI,
- const label tgtSeedI,
- const labelList& srcCellIDs,
- boolList& mapFlag,
- label& startSeedI
- );
-
- //- Find target cell index of cell closest to source cell
- void findNearestCell
- (
- const polyMesh& src,
- const polyMesh& tgt,
- const label srcCellI,
- label& tgtCellI
- );
-
- //- Set the next pair of cells
- void setNextNearestCells
- (
- label& startSeedI,
- label& srcCellI,
- label& tgtCellI,
- boolList& mapFlag,
- const polyMesh& src,
- const polyMesh& tgt,
- const labelList& srcCellIDs
- );
-
- //- Find source cell for target cell
- label findMappedSrcCell
- (
- const polyMesh& tgt,
- const label tgtCellI,
- const List<DynamicList<label> >& tgtToSrc
- ) const;
-
-
- // Cell volume weighted (non-conformal) interpolation
-
- //- Main driver routine for cell volume weighted interpolation
- void calcCellVolumeWeight
- (
- const polyMesh& src,
- const polyMesh& tgt,
- const label srcSeedI,
- const label tgtSeedI,
- const labelList& srcCellIDs,
- boolList& mapFlag,
- label& startSeedI
- );
-
- //- Set the next cells in the advancing front algorithm
- void setNextCells
- (
- label& startSeedI,
- label& srcCellI,
- label& tgtCellI,
- const polyMesh& src,
- const polyMesh& tgt,
- const labelList& srcCellIDs,
- const boolList& mapFlag,
- const DynamicList<label>& visitedCells,
- labelList& seedCells
- ) const;
-
- //- Return the true if cells intersect
- bool intersect
- (
- const polyMesh& src,
- const polyMesh& tgt,
- const label srcCellI,
- const label tgtCellI
- ) const;
-
- //- Return the intersection volume between two cells
- scalar interVol
- (
- const polyMesh& src,
- const polyMesh& tgt,
- const label srcCellI,
- const label tgtCellI
- ) const;
-
-
//- Calculate the addressing between overalping regions of src and tgt
// meshes
void calcAddressing(const polyMesh& src, const polyMesh& tgt);
@@ -302,6 +148,13 @@ private:
//- Calculate - main driver function
void calculate();
+ //- Conversion between mesh and patch interpolation methods
+ AMIPatchToPatchInterpolation::interpolationMethod
+ interpolationMethodAMI
+ (
+ const interpolationMethod method
+ ) const;
+
//- Return the list of AMIs between source and target patches
const PtrList<AMIPatchToPatchInterpolation>& patchAMIs() const;
diff --git a/src/sampling/sampledSurface/sampledSurfaces/sampledSurfacesTemplates.C b/src/sampling/sampledSurface/sampledSurfaces/sampledSurfacesTemplates.C
index e6ba8d42f198e7b4778bb8468e28e235258ac772..35880afafa1bdb1cb2d64c3d447dd594844f5fcb 100644
--- a/src/sampling/sampledSurface/sampledSurfaces/sampledSurfacesTemplates.C
+++ b/src/sampling/sampledSurface/sampledSurfaces/sampledSurfacesTemplates.C
@@ -168,18 +168,32 @@ void Foam::sampledSurfaces::sampleAndWrite
template<class GeoField>
void Foam::sampledSurfaces::sampleAndWrite(const IOobjectList& objects)
{
+ wordList names;
if (loadFromFiles_)
{
IOobjectList fieldObjects(objects.lookupClass(GeoField::typeName));
+ names = fieldObjects.names();
+ }
+ else
+ {
+ names = mesh_.thisDb().names<GeoField>();
+ }
+
+ labelList nameIDs(findStrings(fieldSelection_, names));
- wordList names(fieldObjects.names());
+ wordHashSet fieldNames(wordList(names, nameIDs));
- labelList fieldNames(findStrings(fieldSelection_, names));
+ forAllConstIter(wordHashSet, fieldNames, iter)
+ {
+ const word& fieldName = iter.key();
- forAll(fieldNames, fieldI)
+ if ((Pstream::master()) && verbose_)
{
- const word& fieldName = names[fieldNames[fieldI]];
+ Pout<< "sampleAndWrite: " << fieldName << endl;
+ }
+ if (loadFromFiles_)
+ {
const GeoField fld
(
IOobject
@@ -192,30 +206,10 @@ void Foam::sampledSurfaces::sampleAndWrite(const IOobjectList& objects)
mesh_
);
- if ((Pstream::master()) && verbose_)
- {
- Pout<< "sampleAndWrite: " << fieldName << endl;
- }
-
sampleAndWrite(fld);
}
- }
- else
- {
- const wordList fieldNames
- (
- mesh_.thisDb().names<GeoField>(fieldSelection_)
- );
-
- forAll(fieldNames, i)
+ else
{
- const word& fieldName = fieldNames[i];
-
- if ((Pstream::master()) && verbose_)
- {
- Pout<< "sampleAndWrite: " << fieldName << endl;
- }
-
sampleAndWrite
(
mesh_.thisDb().lookupObject<GeoField>(fieldName)
diff --git a/src/thermophysicalModels/solidChemistryModel/makeSolidChemistryModel.H b/src/thermophysicalModels/solidChemistryModel/makeSolidChemistryModel.H
index 510897f541db6c18a6adb5f9f4711027cb195416..febcc0ac698adc395da780e6c70a8e8c6567ba64 100644
--- a/src/thermophysicalModels/solidChemistryModel/makeSolidChemistryModel.H
+++ b/src/thermophysicalModels/solidChemistryModel/makeSolidChemistryModel.H
@@ -48,7 +48,8 @@ namespace Foam
defineTemplateTypeNameAndDebugWithName \
( \
sChemistryl##Comp##SThermo, \
- (#sChemistry"<"#Comp"," + SThermo::typeName() + ">").c_str(), \
+ (word(sChemistry::typeName_()) + "<"#Comp"," + SThermo::typeName() + \
+ + ">").c_str(), \
0 \
); \
\
diff --git a/src/thermophysicalModels/solidChemistryModel/solidChemistrySolver/makeSolidChemistrySolverType.H b/src/thermophysicalModels/solidChemistryModel/solidChemistrySolver/makeSolidChemistrySolverType.H
index 3f24e11e4a9cb4f8622dfb687c1e6eee6a50d840..734d51083100def428f7c985368c39bf70854698 100644
--- a/src/thermophysicalModels/solidChemistryModel/solidChemistrySolver/makeSolidChemistrySolverType.H
+++ b/src/thermophysicalModels/solidChemistryModel/solidChemistrySolver/makeSolidChemistrySolverType.H
@@ -47,8 +47,8 @@ namespace Foam
defineTemplateTypeNameAndDebugWithName \
( \
SS##Schem##Comp##SThermo##GThermo, \
- (#SS"<"#Schem"<"#Comp"," + SThermo::typeName() + "," \
- + GThermo::typeName() + ">>").c_str(), \
+ (#SS"<" + word(Schem::typeName_()) +"<"#Comp"," + SThermo::typeName() \
+ + "," + GThermo::typeName() + ">>").c_str(), \
0 \
); \
\
diff --git a/src/triSurface/triSurface/triSurface.H b/src/triSurface/triSurface/triSurface.H
index 62c6cd3d6fbddaf89b8d8a2dcd9860aedb00c10b..b32ca479474e4f51bda07fd8e5a105ed34978aeb 100644
--- a/src/triSurface/triSurface/triSurface.H
+++ b/src/triSurface/triSurface/triSurface.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -175,6 +175,7 @@ class triSurface
void writeDXGeometry(const bool, Ostream&) const;
void writeDXTrailer(Ostream&) const;
+
// Static private functions
//- Convert faces to labelledTri. All get same region.
@@ -203,6 +204,7 @@ class triSurface
//- read non-comment line
static string getLineNoComment(IFstream&);
+
protected:
// Protected Member Functions
@@ -219,6 +221,7 @@ protected:
return static_cast<List<Face>&>(*this);
}
+
public:
// Public typedefs
@@ -287,7 +290,6 @@ public:
triSurface(const triSurface&);
-
//- Destructor
~triSurface();
@@ -324,6 +326,7 @@ public:
// If >2 neighbours: undetermined.
const labelList& edgeOwner() const;
+
// Edit
//- Move points
@@ -381,6 +384,7 @@ public:
labelList& faceMap
) const;
+
// Write
//- Write to Ostream in simple FOAM format
diff --git a/src/turbulenceModels/incompressible/LES/dynLagrangian/dynLagrangian.C b/src/turbulenceModels/incompressible/LES/dynLagrangian/dynLagrangian.C
index 32d23f4eab9caa2b93364f192d30317d829584a9..24ef2ab2316a18b6f14cb3917b35a7cab0bddffd 100644
--- a/src/turbulenceModels/incompressible/LES/dynLagrangian/dynLagrangian.C
+++ b/src/turbulenceModels/incompressible/LES/dynLagrangian/dynLagrangian.C
@@ -47,7 +47,7 @@ void dynLagrangian::updateSubGridScaleFields
const tmp<volTensorField>& gradU
)
{
- nuSgs_ = (flm_/fmm_)*delta()*sqrt(k(gradU));
+ nuSgs_ = (flm_/fmm_)*sqr(delta())*mag(dev(symm(gradU)));
nuSgs_.correctBoundaryConditions();
}
diff --git a/src/turbulenceModels/incompressible/LES/dynLagrangian/dynLagrangian.H b/src/turbulenceModels/incompressible/LES/dynLagrangian/dynLagrangian.H
index 004efeaaed7e9fe92cb047eb008145995c74e7b4..a0c703dd1cc555d600b628b6744749e3c805d010 100644
--- a/src/turbulenceModels/incompressible/LES/dynLagrangian/dynLagrangian.H
+++ b/src/turbulenceModels/incompressible/LES/dynLagrangian/dynLagrangian.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -153,7 +153,10 @@ public:
//- Return SGS kinetic energy
tmp<volScalarField> k(const tmp<volTensorField>& gradU) const
{
- return 2.0*sqr(delta())*magSqr(dev(symm(gradU)));
+ return
+ pow(2.0*flm_/fmm_, 2.0/3.0)
+ * pow(ce_, -2.0/3.0)
+ * sqr(delta())*magSqr(dev(symm(gradU)));
}
//- Return SGS kinetic energy
diff --git a/tutorials/combustion/fireFoam/les/oppositeBurningPanels/constant/panelRegion/chemistryProperties b/tutorials/combustion/fireFoam/les/oppositeBurningPanels/constant/panelRegion/chemistryProperties
index dd30353e6960ee7e6778a66426cf7a55b691694d..9f8ecdd196e68ec48372f8e6325f33d56c43f38f 100644
--- a/tutorials/combustion/fireFoam/les/oppositeBurningPanels/constant/panelRegion/chemistryProperties
+++ b/tutorials/combustion/fireFoam/les/oppositeBurningPanels/constant/panelRegion/chemistryProperties
@@ -18,7 +18,7 @@ FoamFile
chemistryType
{
chemistrySolver ode;
- chemistryThermo pyrolysisChemistryModel;
+ chemistryThermo pyrolysis;
}
chemistry on;