Compare commits

...

1 Commits

Author SHA1 Message Date
8f53fcfdf6 ENH: cyclicAMI: alternative low-weight correction 2022-02-15 16:55:13 +00:00
10 changed files with 1395 additions and 10 deletions

View File

@ -798,6 +798,20 @@ bool Foam::AMIInterpolation::calculate
}
bool Foam::AMIInterpolation::calculate
(
const polyMesh& mesh,
const label srcPatchi,
const primitivePatch& srcPatch,
const label tgtPatchi,
const primitivePatch& tgtPatch,
const autoPtr<searchableSurface>& surfPtr
)
{
return calculate(srcPatch, tgtPatch, surfPtr);
}
void Foam::AMIInterpolation::reset
(
autoPtr<mapDistribute>&& srcToTgtMap,

View File

@ -73,6 +73,8 @@ SourceFiles
namespace Foam
{
class polyMesh;
/*---------------------------------------------------------------------------*\
Class AMIInterpolation Declaration
\*---------------------------------------------------------------------------*/
@ -198,15 +200,6 @@ protected:
) const;
// Access
//- Return the orginal src patch with optionally updated points
inline const primitivePatch& srcPatch0() const;
//- Return the orginal tgt patch with optionally updated points
inline const primitivePatch& tgtPatch0() const;
// Evaluation
//- Normalise the (area) weights - suppresses numerical error in
@ -348,6 +341,15 @@ public:
// Member Functions
// Access
//- Return the orginal src patch with optionally updated points
inline const primitivePatch& srcPatch0() const;
//- Return the orginal tgt patch with optionally updated points
inline const primitivePatch& tgtPatch0() const;
// Access
//- Access to the up-to-date flag
@ -446,6 +448,12 @@ public:
//- patch weights (i.e. the sum before normalisation)
inline const scalarField& tgtWeightsSum() const;
//- Return const access to target patch face centroids
inline const pointListList& tgtCentroids() const;
//- Return access to target patch face centroids
inline pointListList& tgtCentroids();
//- Return access to normalisation factor of target
//- patch weights (i.e. the sum before normalisation)
inline scalarField& tgtWeightsSum();
@ -466,6 +474,17 @@ public:
const autoPtr<searchableSurface>& surfPtr = nullptr
);
//- Update addressing, weights and (optional) centroids
virtual bool calculate
(
const polyMesh& mesh,
const label srcPatchi,
const primitivePatch& srcPatch,
const label tgtPatchi,
const primitivePatch& tgtPatch,
const autoPtr<searchableSurface>& surfPtr = nullptr
);
//- Set the maps, addresses and weights from an external source
void reset
(

View File

@ -222,6 +222,18 @@ inline const Foam::scalarField& Foam::AMIInterpolation::tgtWeightsSum() const
}
inline const Foam::pointListList& Foam::AMIInterpolation::tgtCentroids() const
{
return tgtCentroids_;
}
inline Foam::pointListList& Foam::AMIInterpolation::tgtCentroids()
{
return tgtCentroids_;
}
inline Foam::scalarField& Foam::AMIInterpolation::tgtWeightsSum()
{
return tgtWeightsSum_;

View File

@ -0,0 +1,52 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
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/>.
\*---------------------------------------------------------------------------*/
// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
template<class Type, class PrimitivePatchType>
Foam::Ostream& Foam::operator<<
(
Foam::Ostream& os,
const Foam::edgeTopoDistancesData<Type, PrimitivePatchType>& wDist
)
{
return os << wDist.distance_ << token::SPACE << wDist.data_;
}
template<class Type, class PrimitivePatchType>
Foam::Istream& Foam::operator>>
(
Foam::Istream& is,
Foam::edgeTopoDistancesData<Type, PrimitivePatchType>& wDist
)
{
return is >> wDist.distance_ >> wDist.data_;
}
// ************************************************************************* //

View File

@ -0,0 +1,261 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
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::edgeTopoDistancesData
Description
For use with PatchEdgeFaceWave. Determines topological distance to
starting edges. Templated on passive transported data.
SourceFiles
edgeTopoDistancesDataI.H
edgeTopoDistancesData.C
\*---------------------------------------------------------------------------*/
#ifndef edgeTopoDistancesData_H
#define edgeTopoDistancesData_H
#include "point.H"
#include "tensor.H"
#include "indirectPrimitivePatch.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// Forward Declarations
class polyPatch;
class polyMesh;
template<class Type, class PrimitivePatchType>
class edgeTopoDistancesData;
template<class Type, class PrimitivePatchType>
Istream& operator>>
(
Istream&,
edgeTopoDistancesData<Type, PrimitivePatchType>&
);
template<class Type, class PrimitivePatchType>
Ostream& operator<<
(
Ostream&,
const edgeTopoDistancesData<Type, PrimitivePatchType>&
);
/*---------------------------------------------------------------------------*\
Class edgeTopoDistancesData Declaration
\*---------------------------------------------------------------------------*/
template<class Type, class PrimitivePatchType = indirectPrimitivePatch>
class edgeTopoDistancesData
{
public:
//- Class used to pass additional data in
struct trackData
{
//- Max number of entries stored per element
label n_;
};
protected:
// Protected data
//- Distance (up to size n_)
labelList distance_;
//- Starting data (up to size n_)
List<Type> data_;
public:
typedef Type dataType;
// Constructors
//- Construct null with null constructor for distance and data
inline edgeTopoDistancesData();
//- Construct from distance, data
inline edgeTopoDistancesData
(
const labelList& distance,
const List<Type>& data
);
// Member Functions
// Access
inline const labelList& distance() const
{
return distance_;
}
inline const List<Type>& data() const
{
return data_;
}
// Needed by PatchEdgeFaceWave
//- Check whether origin has been changed at all or
// still contains original (invalid) value.
template<class TrackingData>
inline bool valid(TrackingData& td) const;
//- Apply rotation matrix
template<class TrackingData>
inline void transform
(
const polyMesh& mesh,
const PrimitivePatchType& patch,
const tensor& rotTensor,
const scalar tol,
TrackingData& td
);
//- Influence of face on edge
template<class TrackingData>
inline bool updateEdge
(
const polyMesh& mesh,
const PrimitivePatchType& patch,
const label edgeI,
const label facei,
const edgeTopoDistancesData<Type, PrimitivePatchType>& faceInfo,
const scalar tol,
TrackingData& td
);
//- New information for edge (from e.g. coupled edge)
template<class TrackingData>
inline bool updateEdge
(
const polyMesh& mesh,
const PrimitivePatchType& patch,
const edgeTopoDistancesData<Type, PrimitivePatchType>& edgeInfo,
const bool sameOrientation,
const scalar tol,
TrackingData& td
);
//- Influence of edge on face.
template<class TrackingData>
inline bool updateFace
(
const polyMesh& mesh,
const PrimitivePatchType& patch,
const label facei,
const label edgeI,
const edgeTopoDistancesData<Type, PrimitivePatchType>& edgeInfo,
const scalar tol,
TrackingData& td
);
//- Same (like operator==)
template<class TrackingData>
inline bool equal
(
const edgeTopoDistancesData<Type, PrimitivePatchType>&,
TrackingData&
) const;
// Member Operators
// Needed for List IO
inline bool operator==
(
const edgeTopoDistancesData<Type, PrimitivePatchType>&
) const;
inline bool operator!=
(
const edgeTopoDistancesData<Type, PrimitivePatchType>&
) const;
// IOstream Operators
friend Ostream& operator<< <Type, PrimitivePatchType>
(
Ostream&,
const edgeTopoDistancesData<Type, PrimitivePatchType>&
);
friend Istream& operator>> <Type, PrimitivePatchType>
(
Istream&,
edgeTopoDistancesData<Type, PrimitivePatchType>&
);
};
// * * * * * * * * * * * * * * * * * Traits * * * * * * * * * * * * * * * * //
//- Data are contiguous if data type is contiguous
template<class Type, class PrimitivePatchType>
struct is_contiguous<edgeTopoDistancesData<Type, PrimitivePatchType>> :
std::false_type {};
//- Data are contiguous label if data type is label
template<class Type, class PrimitivePatchType>
struct is_contiguous_label<edgeTopoDistancesData<Type, PrimitivePatchType>> :
std::false_type {};
//- Data are contiguous scalar if data type is scalar
template<class Type, class PrimitivePatchType>
struct is_contiguous_scalar<edgeTopoDistancesData<Type, PrimitivePatchType>> :
std::false_type {};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
#include "edgeTopoDistancesData.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "edgeTopoDistancesDataI.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,224 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
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 "polyMesh.H"
#include "transform.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template<class Type, class PrimitivePatchType>
inline
Foam::edgeTopoDistancesData<Type, PrimitivePatchType>::edgeTopoDistancesData()
:
distance_(),
data_()
{}
template<class Type, class PrimitivePatchType>
inline
Foam::edgeTopoDistancesData<Type, PrimitivePatchType>::edgeTopoDistancesData
(
const labelList& distance,
const List<Type>& data
)
:
distance_(distance),
data_(data)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type, class PrimitivePatchType>
template<class TrackingData>
inline bool Foam::edgeTopoDistancesData<Type, PrimitivePatchType>::valid
(
TrackingData& td
) const
{
return distance_.size();
}
template<class Type, class PrimitivePatchType>
template<class TrackingData>
inline void Foam::edgeTopoDistancesData<Type, PrimitivePatchType>::transform
(
const polyMesh& mesh,
const PrimitivePatchType& patch,
const tensor& rotTensor,
const scalar tol,
TrackingData& td
)
{}
template<class Type, class PrimitivePatchType>
template<class TrackingData>
inline bool Foam::edgeTopoDistancesData<Type, PrimitivePatchType>::updateEdge
(
const polyMesh& mesh,
const PrimitivePatchType& patch,
const label edgeI,
const label facei,
const edgeTopoDistancesData<Type, PrimitivePatchType>& faceInfo,
const scalar tol,
TrackingData& td
)
{
// From face to edge
if (distance_.size() >= td.n_)
{
return false;
}
const label oldSize = data_.size();
forAll(faceInfo.data_, i)
{
const auto& d = faceInfo.data_[i];
if (!data_.found(d))
{
data_.append(d);
distance_.append(faceInfo.distance_[i] + 1);
}
}
return data_.size() > oldSize;
}
template<class Type, class PrimitivePatchType>
template<class TrackingData>
inline bool Foam::edgeTopoDistancesData<Type, PrimitivePatchType>::updateEdge
(
const polyMesh& mesh,
const PrimitivePatchType& patch,
const edgeTopoDistancesData<Type, PrimitivePatchType>& edgeInfo,
const bool sameOrientation,
const scalar tol,
TrackingData& td
)
{
// From edge to edge (e.g. coupled edges)
if (distance_.size() >= td.n_)
{
return false;
}
const label oldSize = data_.size();
forAll(edgeInfo.data_, i)
{
const auto& d = edgeInfo.data_[i];
if (!data_.found(d))
{
data_.append(d);
distance_.append(edgeInfo.distance_[i]);
}
}
return data_.size() > oldSize;
}
template<class Type, class PrimitivePatchType>
template<class TrackingData>
inline bool Foam::edgeTopoDistancesData<Type, PrimitivePatchType>::updateFace
(
const polyMesh& mesh,
const PrimitivePatchType& patch,
const label facei,
const label edgeI,
const edgeTopoDistancesData<Type, PrimitivePatchType>& edgeInfo,
const scalar tol,
TrackingData& td
)
{
// From edge to edge (e.g. coupled edges)
if (distance_.size() >= td.n_)
{
return false;
}
const label oldSize = data_.size();
forAll(edgeInfo.data_, i)
{
const auto& d = edgeInfo.data_[i];
if (!data_.found(d))
{
data_.append(d);
distance_.append(edgeInfo.distance_[i]);
}
}
return data_.size() > oldSize;
}
template<class Type, class PrimitivePatchType>
template<class TrackingData>
inline bool Foam::edgeTopoDistancesData<Type, PrimitivePatchType>::equal
(
const edgeTopoDistancesData<Type, PrimitivePatchType>& rhs,
TrackingData& td
) const
{
return operator==(rhs);
}
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
template<class Type, class PrimitivePatchType>
inline bool Foam::edgeTopoDistancesData<Type, PrimitivePatchType>::operator==
(
const Foam::edgeTopoDistancesData<Type, PrimitivePatchType>& rhs
) const
{
return distance() == rhs.distance() && data() == rhs.data();
}
template<class Type, class PrimitivePatchType>
inline bool Foam::edgeTopoDistancesData<Type, PrimitivePatchType>::operator!=
(
const Foam::edgeTopoDistancesData<Type, PrimitivePatchType>& rhs
) const
{
return !(*this == rhs);
}
// ************************************************************************* //

View File

@ -0,0 +1,578 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
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 "lowWeightCorrection.H"
#include "addToRunTimeSelectionTable.H"
#include "edgeTopoDistancesData.H"
#include "PatchEdgeFaceWave.H"
#include "OBJstream.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(lowWeightCorrection, 0);
addToRunTimeSelectionTable(AMIInterpolation, lowWeightCorrection, dict);
addToRunTimeSelectionTable
(
AMIInterpolation,
lowWeightCorrection,
component
);
}
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::lowWeightCorrection::findNearest
(
const polyMesh& mesh,
const primitivePatch& pp,
const globalIndex& globalFaces,
const labelList& uncoveredFaces,
labelListList& nearestCoveredFaces
) const
{
// Data on all edges and faces
typedef edgeTopoDistancesData<label, primitivePatch> Type;
List<Type> allEdgeInfo(pp.nEdges());
List<Type> allFaceInfo(pp.size());
// Seed all covered faces
bitSet isUncovered(pp.size(), uncoveredFaces);
DynamicList<label> changedEdges(pp.nEdges());
DynamicList<Type> changedInfo(pp.nEdges());
forAll(pp, facei)
{
if (!isUncovered(facei))
{
const labelList& fEdges = pp.faceEdges()[facei];
const Type seed
(
labelList(1, 0),
labelList(1, globalFaces.toGlobal(facei))
);
allFaceInfo[facei] = seed;
for (const label edgei : fEdges)
{
changedEdges.append(edgei);
changedInfo.append(seed);
}
}
}
// Walk
Type::trackData td;
td.n_ = nDonors_;
PatchEdgeFaceWave
<
primitivePatch,
Type,
Type::trackData
> calc
(
mesh,
pp,
changedEdges,
changedInfo,
allEdgeInfo,
allFaceInfo,
0, //returnReduce(pp.nEdges(), sumOp<label>())
td
);
calc.iterate(nIters_); // should be enough iterations?
nearestCoveredFaces.setSize(uncoveredFaces.size());
forAll(uncoveredFaces, i)
{
const label facei = uncoveredFaces[i];
if (allFaceInfo[facei].valid(calc.data()))
{
// Collect donor faces
nearestCoveredFaces[i] = allFaceInfo[facei].data();
}
}
}
void Foam::lowWeightCorrection::calculateStencil
(
const label facei,
const point& sample,
const labelList& covered,
const labelListList& stencilAddressing,
labelListList& addressing,
scalarListList& weights,
scalarField& sumWeights,
const pointField& faceCentres
) const
{
addressing[facei].clear();
weights[facei].clear();
sumWeights[facei] = Zero;
// The (uncovered) facei has local donors
// Add the donors
for (const label coveredSlot : covered)
{
// Get remote slots for the local covered face
const labelList& slots = stencilAddressing[coveredSlot];
for (const label sloti : slots)
{
const point& donorPt = faceCentres[sloti];
const label index = addressing[facei].find(sloti);
if (index == -1)
{
addressing[facei].append(sloti);
weights[facei].append(1.0/mag(sample-donorPt));
}
else
{
weights[facei][index] += 1.0/mag(sample-donorPt);
}
}
}
scalarList& wghts = weights[facei];
const scalar w = sum(wghts);
forAll(wghts, i)
{
wghts[i] /= w;
}
sumWeights[facei] = 1.0;
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::lowWeightCorrection::lowWeightCorrection
(
const dictionary& dict,
const bool reverseTarget
)
:
AMIInterpolation(dict, reverseTarget),
lowWeightCorrection_(dict.get<scalar>("lowWeightCorrection")),
nDonors_(dict.get<label>("nDonors")),
nIters_(dict.get<label>("nIters")),
AMIPtr_
(
AMIInterpolation::New
(
dict.subDict(type() + "Coeffs").get<word>("AMIMethod"),
dict.subDict(type() + "Coeffs"),
reverseTarget
)
)
{
if (nDonors_ <= 0 || nIters_ <= 0)
{
WarningInFunction << "Disabled low-weight correction. Using "
<< AMIPtr_->type() << " AMI interpolation" << endl;
}
}
Foam::lowWeightCorrection::lowWeightCorrection
(
const bool requireMatch,
const bool reverseTarget,
const scalar lowWeightCorrection
)
:
AMIInterpolation(requireMatch, reverseTarget, lowWeightCorrection),
lowWeightCorrection_(-1),
nDonors_(0),
nIters_(0),
AMIPtr_()
{}
Foam::lowWeightCorrection::lowWeightCorrection(const lowWeightCorrection& ami)
:
AMIInterpolation(ami),
lowWeightCorrection_(ami.lowWeightCorrection_),
nDonors_(ami.nDonors_),
nIters_(ami.nIters_),
AMIPtr_(ami.AMIPtr_.clone())
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::lowWeightCorrection::calculate
(
const polyMesh& mesh,
const label srcPatchi,
const primitivePatch& srcPatch,
const label tgtPatchi,
const primitivePatch& tgtPatch,
const autoPtr<searchableSurface>& surfPtr
)
{
if (upToDate_)
{
return false;
}
AMIInterpolation::calculate(srcPatch, tgtPatch, surfPtr);
AMIPtr_->calculate
(
mesh,
srcPatchi,
srcPatch,
tgtPatchi,
tgtPatch,
surfPtr
);
// Take over AMI data
// ~~~~~~~~~~~~~~~~~~
upToDate_ = AMIPtr_->upToDate();
// distributed / singlePatchProc
singlePatchProc_ = AMIPtr_->singlePatchProc();;
// Source patch
srcMagSf_ = AMIPtr_->srcMagSf();
srcAddress_ = AMIPtr_->srcAddress();
srcWeights_ = AMIPtr_->srcWeights();
srcWeightsSum_ = AMIPtr_->srcWeightsSum();
srcCentroids_ = AMIPtr_->srcCentroids();
//TBD: srcPatchPts_ = AMIPtr_->srcPatchPts();
tsrcPatch0_ = AMIPtr_->srcPatch0();
if (AMIPtr_->distributed())
{
srcMapPtr_.reset(new mapDistribute(AMIPtr_->srcMap()));
}
// Target
tgtMagSf_ = AMIPtr_->tgtMagSf();
tgtAddress_ = AMIPtr_->tgtAddress();
tgtWeights_ = AMIPtr_->tgtWeights();
tgtWeightsSum_ = AMIPtr_->tgtWeightsSum();
tgtCentroids_ = AMIPtr_->tgtCentroids();
//TBD: tgtPatchPts_ = AMIPtr_->tgtPatchPts();
ttgtPatch0_ = AMIPtr_->tgtPatch0();
if (AMIPtr_->distributed())
{
tgtMapPtr_.reset(new mapDistribute(AMIPtr_->tgtMap()));
}
// Walk out valid donors
// ~~~~~~~~~~~~~~~~~~~~~
if (nDonors_ > 0 && nIters_ > 0)
{
// Extend source side
{
// Low weight faces
DynamicList<label> uncoveredFaces;
forAll(srcWeightsSum(), facei)
{
if (srcWeightsSum()[facei] < lowWeightCorrection_)
{
uncoveredFaces.append(facei);
}
}
uncoveredFaces.shrink();
// Global indexing for src faces
const globalIndex globalFaces(srcPatch.size());
// Find nearest face with high weight
labelListList nearestCoveredFaces;
findNearest
(
mesh,
mesh.boundaryMesh()[srcPatchi],
globalFaces,
uncoveredFaces,
nearestCoveredFaces
);
// Create map to get remote data into nearestCoveredFaces
// order
List<Map<label>> compactMap;
const mapDistribute uncoveredToPatchMap
(
globalFaces,
nearestCoveredFaces,
compactMap
);
// Get srcPatch faceCentres in stencil ordering
pointField stencilFcs(srcPatch.faceCentres());
uncoveredToPatchMap.distribute(stencilFcs);
// Get addressing (to donors) over in stencil ordering
labelListList stencilAddressing(srcAddress());
uncoveredToPatchMap.distribute(stencilAddressing);
//scalarListList stencilWeights(srcWeights());
//uncoveredToPatchMap.distribute(stencilWeights);
// Target side patch centres
tmp<pointField> totherFcs;
if (AMIPtr_->distributed())
{
totherFcs = new pointField(tgtPatch.faceCentres());
AMIPtr_->tgtMap().distribute(totherFcs.ref());
}
else
{
totherFcs = tgtPatch.faceCentres();
}
const pointField& otherFcs = totherFcs();
//forAll(uncoveredFaces, i)
//{
// const label facei = uncoveredFaces[i];
// const labelList& covered = nearestCoveredFaces[i];
// Pout<< "Uncovered face:" << facei
// << " low weihgt:" << srcWeightsSum()[facei]
// << " at:" << srcPatch.faceCentres()[facei]
// << " has local donors:" << endl;
// for (const label coveredSlot : covered)
// {
// Pout<< " fc:" << stencilFcs[coveredSlot] << nl
// << " which has remotes:"
// << stencilAddressing[coveredSlot] << nl
// << " at:"
// << UIndirectList<point>
// (
// otherFcs,
// stencilAddressing[coveredSlot]
// ) << nl
// //<< " with weights:"
// //<< stencilWeights[coveredSlot] << nl
// << endl;
// }
//}
// Re-do interpolation on uncovered faces
forAll(uncoveredFaces, i)
{
const label facei = uncoveredFaces[i];
calculateStencil
(
facei,
srcPatch.faceCentres()[facei],
nearestCoveredFaces[i],
stencilAddressing,
srcAddress_,
srcWeights_,
srcWeightsSum_,
otherFcs
);
}
//OBJstream os
//(
// mesh.time().path()
// /mesh.boundaryMesh()[srcPatchi].name()+"_src.obj"
//);
//forAll(uncoveredFaces, i)
//{
// const label facei = uncoveredFaces[i];
// const point& fc = srcPatch.faceCentres()[facei];
// const labelList& slots = srcAddress_[facei];
// for (const label sloti : slots)
// {
// os.write(linePointRef(fc, otherFcs[sloti]));
// }
//}
}
// Extend target side
{
// Low weight faces
DynamicList<label> uncoveredFaces;
forAll(tgtWeightsSum(), facei)
{
if (tgtWeightsSum()[facei] < lowWeightCorrection_)
{
uncoveredFaces.append(facei);
}
}
uncoveredFaces.shrink();
// Global indexing for tgt faces
const globalIndex globalFaces(tgtPatch.size());
// Find nearest face with high weight
labelListList nearestCoveredFaces;
findNearest
(
mesh,
mesh.boundaryMesh()[tgtPatchi],
globalFaces,
uncoveredFaces,
nearestCoveredFaces
);
// Create map to get remote data into nearestCoveredFaces
// order
List<Map<label>> compactMap;
const mapDistribute uncoveredToPatchMap
(
globalFaces,
nearestCoveredFaces,
compactMap
);
// Get tgtPatch faceCentres in stencil ordering
pointField stencilFcs(tgtPatch.faceCentres());
uncoveredToPatchMap.distribute(stencilFcs);
// Get addressing (to donors) over in stencil ordering
labelListList stencilAddressing(tgtAddress());
uncoveredToPatchMap.distribute(stencilAddressing);
scalarListList stencilWeights(tgtWeights());
uncoveredToPatchMap.distribute(stencilWeights);
// Target side patch centres
tmp<pointField> totherFcs;
if (AMIPtr_->distributed())
{
totherFcs = new pointField(srcPatch.faceCentres());
AMIPtr_->srcMap().distribute(totherFcs.ref());
}
else
{
totherFcs = srcPatch.faceCentres();
}
const pointField& otherFcs = totherFcs();
//forAll(uncoveredFaces, i)
//{
// const label facei = uncoveredFaces[i];
// const labelList& covered = nearestCoveredFaces[i];
// Pout<< "Uncovered face:" << facei
// << " low weihgt:" << tgtWeights()[facei]
// << " at:" << tgtPatch.faceCentres()[facei]
// << " has local donors:" << endl;
// for (const label coveredSlot : covered)
// {
// Pout<< " fc:" << stencilFcs[coveredSlot] << nl
// << " which has remotes:"
// << stencilAddressing[coveredSlot] << nl
// << " at:"
// << UIndirectList<point>
// (
// otherFcs,
// stencilAddressing[coveredSlot]
// ) << nl
// << " with weights:"
// << stencilWeights[coveredSlot] << nl
// << endl;
// }
//}
// Re-do interpolation on uncovered faces
forAll(uncoveredFaces, i)
{
const label facei = uncoveredFaces[i];
calculateStencil
(
facei,
tgtPatch.faceCentres()[facei],
nearestCoveredFaces[i],
stencilAddressing,
tgtAddress_,
tgtWeights_,
tgtWeightsSum_,
otherFcs
);
}
//OBJstream os
//(
// mesh.time().path()
// /mesh.boundaryMesh()[srcPatchi].name()+"_tgt.obj"
//);
//forAll(uncoveredFaces, i)
//{
// const label facei = uncoveredFaces[i];
// const point& fc = tgtPatch.faceCentres()[facei];
// const labelList& slots = tgtAddress_[facei];
// for (const label sloti : slots)
// {
// os.write(linePointRef(fc, otherFcs[sloti]));
// }
//}
}
////- Write face connectivity as OBJ file
//writeFaceConnectivity
//(
// srcPatch,
// tgtPatch,
// srcAddress_
//);
}
return upToDate_;
}
void Foam::lowWeightCorrection::write(Ostream& os) const
{
AMIInterpolation::write(os);
os.writeEntry("nDonors", nDonors_);
os.writeEntry("nIters", nIters_);
os.beginBlock(word(this->type() + "Coeffs"));
AMIPtr_->write(os);
os.endBlock();
}
// ************************************************************************* //

View File

@ -0,0 +1,216 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
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::lowWeightCorrection
Description
Arbitrary Mesh Interface (AMI) method wrapper that extends low-weight
faces to use donors of local neighbouring (non-low weight) faces.
It collects valid nearby faces and uses their stencil instead.
Example of the boundary specification:
\verbatim
<patchName>
{
type cyclicAMI;
AMIMethod lowWeightCorrection;
// Low weight value; which faces to override
lowWeightCorrection 0.1;
// How far to walk out from valid faces. 0 means no change.
nIters 2;
// Limit amount of valid faces to store.
nDonors 2;
// Underlying AMI method to use.
lowWeightCorrectionCoeffs
{
AMIMethod faceAreaWeightAMI;
}
}
\endverbatim
This version walks out the local faces and then collects all their
remote donors. Hence even for e.g. 2 local faces (nDonors) the stencil
might still include lots of remote faces.
The number of iterations to walk out might be set to infinite if one
wants a value on all faces.
SourceFiles
lowWeightCorrection.C
\*---------------------------------------------------------------------------*/
#ifndef lowWeightCorrection_H
#define lowWeightCorrection_H
#include "AMIInterpolation.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class lowWeightCorrection Declaration
\*---------------------------------------------------------------------------*/
class lowWeightCorrection
:
public AMIInterpolation
{
private:
// Private Data
//- Threshold weight below which interpolation switches to distance
// weighted
const scalar lowWeightCorrection_;
//- Number of donors to weigh
const label nDonors_;
//- Number of iterations to walk out to search for donors
const label nIters_;
//- Uncorrected AMI
autoPtr<AMIInterpolation> AMIPtr_;
// Private Member Functions
//- Find (local) donor faces for all uncovered faces
void findNearest
(
const polyMesh& mesh,
const primitivePatch& pp,
const globalIndex& globalFaces,
const labelList& uncoveredFaces,
labelListList& nearestCoveredFaces
) const;
//- Determine stencil for single face. Gets list of (local) faces
// whose addressing can be used for calculating the stencil
void calculateStencil
(
const label facei,
const point& sample,
const labelList& nearestCovered,
const labelListList& stencilAddressing,
labelListList& addressing,
scalarListList& weights,
scalarField& sumWeights,
const pointField& faceCentres
) const;
//- No copy assignment
void operator=(const lowWeightCorrection&) = delete;
public:
//- Runtime type information
TypeName("lowWeightCorrection");
// Constructors
//- Construct from dictionary
lowWeightCorrection
(
const dictionary& dict,
const bool reverseTarget = false
);
//- Construct from components
lowWeightCorrection
(
const bool requireMatch = true,
const bool reverseTarget = false,
const scalar lowWeightCorrection = -1
);
//- Construct as copy
lowWeightCorrection(const lowWeightCorrection& ami);
//- Construct and return a clone
virtual autoPtr<AMIInterpolation> clone() const
{
return autoPtr<AMIInterpolation>(new lowWeightCorrection(*this));
}
//- Destructor
virtual ~lowWeightCorrection() = default;
// Member Functions
//- Update addressing and weights
virtual bool calculate
(
const primitivePatch& srcPatch,
const primitivePatch& tgtPatch,
const autoPtr<searchableSurface>& surfPtr = nullptr
)
{
NotImplemented;
return false;
}
//- Update addressing and weights
virtual bool calculate
(
const polyMesh& mesh,
const label srcPatchi,
const primitivePatch& srcPatch,
const label tgtPatchi,
const primitivePatch& tgtPatch,
const autoPtr<searchableSurface>& surfPtr = nullptr
);
// I-O
//- Write
virtual void write(Ostream& os) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -435,7 +435,15 @@ void Foam::cyclicAMIPolyPatch::resetAMI(const UList<point>& points) const
// Construct/apply AMI interpolation to determine addressing and weights
AMIPtr_->upToDate() = false;
AMIPtr_->calculate(patch0, nbrPatch0, surfPtr());
AMIPtr_->calculate
(
boundaryMesh().mesh(),
index(),
patch0,
nbr.index(),
nbrPatch0,
surfPtr()
);
if (debug)
{

View File

@ -293,6 +293,7 @@ $(AMI)/GAMG/interfaceFields/cyclicACMIGAMGInterfaceField/cyclicACMIGAMGInterface
$(AMI)/triangle2D/triangle2D.C
$(AMI)/AMIInterpolation/faceAreaWeightAMI2D/faceAreaWeightAMI2D.C
$(AMI)/AMIInterpolation/lowWeightCorrection/lowWeightCorrection.C
AMICycPatches=$(AMI)/patches/cyclicAMI