ENH: Added and abstracted out new AMI methods - direct, mapNearest and faceAreaWeight

This commit is contained in:
andy
2013-04-11 14:52:53 +01:00
parent 652e318422
commit 9817cb318d
10 changed files with 2308 additions and 0 deletions

View File

@ -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()
{}
// ************************************************************************* //

View File

@ -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
// ************************************************************************* //

View File

@ -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_;
}
// ************************************************************************* //

View File

@ -0,0 +1,79 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / 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>
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
)
{
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
)
);
}
// ************************************************************************* //

View File

@ -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(srcAddr[i].size(), magSf);
}
forAll(tgtAddr, i)
{
scalar magSf = this->tgtMagSf_[i];
tgtAddress[i].transfer(tgtAddr[i]);
tgtWeights[i] = scalarList(tgtAddr[i].size(), magSf);
}
}
// ************************************************************************* //

View File

@ -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
// ************************************************************************* //

View File

@ -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]);
}
}
// ************************************************************************* //

View File

@ -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("faceAreaWeight");
// 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
// ************************************************************************* //

View File

@ -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(srcAddr[i].size(), magSf);
}
forAll(tgtAddr, i)
{
scalar magSf = this->tgtMagSf_[i];
tgtAddress[i].transfer(tgtAddr[i]);
tgtWeights[i] = scalarList(tgtAddr[i].size(), magSf);
}
}
// ************************************************************************* //

View File

@ -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
// ************************************************************************* //