meshToMesh: Stabilisation

Stabilisation has been added to the mapping of fields between consistent
meshes. This means that if part of the target mesh is found not to
connect with the source mesh, then its values will be set by propagating
a value from the closest part of the target mesh did successfully
connect to the source. This propagation is achieved by means of a mesh
wave.

This stabilisation applies to both cell and patch fields, and any and
all ancillary fields that may be being stored by the patch boundary
conditions. It applies to the mapping performed by both mapFieldsPar and
the run-time mapping meshToMesh topology changer.

This fixes the previous situation in mapping between consistent meshes
in which target elements which did not connect to the source would be
given an undefined value, which would cause either a floating point
error, or (worse) incorrect operation.
This commit is contained in:
Will Bainbridge
2023-02-15 15:51:34 +00:00
parent 2c247c3e8f
commit 9175fb13fd
25 changed files with 1060 additions and 345 deletions

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2022 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2022-2023 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -98,6 +98,16 @@ public:
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
template<>
class pTraits<remote>
{
public:
typedef remote cmptType;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2022 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2022-2023 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -75,4 +75,20 @@ inline Foam::Istream& Foam::operator>>(Istream& is, remote& p)
}
// * * * * * * * * * * * * * * * Global Functions * * * * * * * * * * * * * //
namespace Foam
{
template<class Cmpt>
class Tensor;
inline remote transform(const Tensor<scalar>&, const remote& p)
{
return p;
}
}
// ************************************************************************* //

View File

@ -517,7 +517,5 @@ solver/solverNew.C
fvMeshToFvMesh/fvMeshToFvMesh.C
fvMeshToFvMesh/patchToPatchFvPatchFieldMapper.C
fvMeshToFvMesh/patchToPatchLeftOverFvPatchFieldMapper.C
fvMeshToFvMesh/patchToPatchNormalisedFvPatchFieldMapper.C
LIB = $(FOAM_LIBBIN)/libfiniteVolume

View File

@ -26,8 +26,7 @@ License
#include "fvMeshToFvMesh.H"
#include "directFvPatchFieldMapper.H"
#include "identityFvPatchFieldMapper.H"
#include "patchToPatchNormalisedFvPatchFieldMapper.H"
#include "patchToPatchLeftOverFvPatchFieldMapper.H"
#include "patchToPatchFvPatchFieldMapper.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
@ -125,10 +124,10 @@ Foam::tmp<Foam::VolField<Type>> Foam::fvMeshToFvMesh::srcToTgt
// Construct target patch fields as copies of source patch fields, but do
// not map values yet
PtrList<fvPatchField<Type>> tgtPatchFields(tgtMesh.boundary().size());
forAll(srcToTgtPatchIDs(), i)
forAll(patchIDs(), i)
{
const label srcPatchi = srcToTgtPatchIDs()[i].first();
const label tgtPatchi = srcToTgtPatchIDs()[i].second();
const label srcPatchi = patchIDs()[i].first();
const label tgtPatchi = patchIDs()[i].second();
if (!tgtPatchFields.set(tgtPatchi))
{
@ -188,18 +187,18 @@ Foam::tmp<Foam::VolField<Type>> Foam::fvMeshToFvMesh::srcToTgt
ttgtFld.ref().boundaryFieldRef();
// Mapped patches
forAll(srcToTgtPatchIDs(), i)
forAll(patchIDs(), i)
{
const label srcPatchi = srcToTgtPatchIDs()[i].first();
const label tgtPatchi = srcToTgtPatchIDs()[i].second();
const label srcPatchi = patchIDs()[i].first();
const label tgtPatchi = patchIDs()[i].second();
tgtBfld[tgtPatchi].map
(
srcFld.boundaryField()[srcPatchi],
patchToPatchNormalisedFvPatchFieldMapper
(
srcToTgtPatchToPatches()[i],
true
patchInterpolation(i),
tgtPatchStabilisation(i)
)
);
}
@ -241,10 +240,10 @@ Foam::tmp<Foam::VolField<Type>> Foam::fvMeshToFvMesh::srcToTgt
ttgtFld.ref().boundaryFieldRef();
// Mapped patches
forAll(srcToTgtPatchIDs(), i)
forAll(patchIDs(), i)
{
const label srcPatchi = srcToTgtPatchIDs()[i].first();
const label tgtPatchi = srcToTgtPatchIDs()[i].second();
const label srcPatchi = patchIDs()[i].first();
const label tgtPatchi = patchIDs()[i].second();
tgtBfld[tgtPatchi].map
(
@ -254,11 +253,7 @@ Foam::tmp<Foam::VolField<Type>> Foam::fvMeshToFvMesh::srcToTgt
tgtBfld[tgtPatchi].map
(
srcFld.boundaryField()[srcPatchi],
patchToPatchLeftOverFvPatchFieldMapper
(
srcToTgtPatchToPatches()[i],
true
)
patchToPatchLeftOverFvPatchFieldMapper(patchInterpolation(i))
);
}
@ -283,14 +278,18 @@ Foam::tmp<Foam::VolInternalField<Type>> Foam::fvMeshToFvMesh::srcToTgt
const VolInternalField<Type>& srcFld
) const
{
return
tmp<VolInternalField<Type>> ttgtFld =
VolInternalField<Type>::New
(
typedName("interpolate(" + srcFld.name() + ")"),
static_cast<const fvMesh&>(meshToMesh::tgtMesh()),
srcFld.dimensions(),
srcToTgtCellsToCells().srcToTgt(srcFld)
cellsInterpolation().srcToTgt(srcFld)
);
tgtCellsStabilisation().stabilise(ttgtFld.ref());
return ttgtFld;
}
@ -307,7 +306,7 @@ Foam::tmp<Foam::VolInternalField<Type>> Foam::fvMeshToFvMesh::srcToTgt
typedName("interpolate(" + srcFld.name() + ")"),
static_cast<const fvMesh&>(meshToMesh::tgtMesh()),
leftOverTgtFld.dimensions(),
srcToTgtCellsToCells().srcToTgt(srcFld, leftOverTgtFld)
cellsInterpolation().srcToTgt(srcFld, leftOverTgtFld)
);
}

View File

@ -48,4 +48,18 @@ Foam::patchToPatchFvPatchFieldMapper::operator()
}
FOR_ALL_FIELD_TYPES
(
IMPLEMENT_FIELD_MAPPER_OPERATOR,
patchToPatchLeftOverFvPatchFieldMapper
)
FOR_ALL_FIELD_TYPES
(
IMPLEMENT_FIELD_MAPPER_OPERATOR,
patchToPatchNormalisedFvPatchFieldMapper
)
// ************************************************************************* //

View File

@ -36,6 +36,8 @@ Description
#include "fvPatchFieldMapper.H"
#include "patchToPatch.H"
#include "patchToPatchStabilisation.H"
#include "ListOps.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -43,7 +45,7 @@ namespace Foam
{
/*---------------------------------------------------------------------------*\
Class patchToPatchFvPatchFieldMapper Declaration
Class patchToPatchFvPatchFieldMapper Declaration
\*---------------------------------------------------------------------------*/
class patchToPatchFvPatchFieldMapper
@ -54,49 +56,20 @@ protected:
// Protected Data
//- Patch-to-patch mapping engine
//- Patch-to-patch interpolation engine
const patchToPatch& pToP_;
//- Mapping direction
const bool forward_;
//- Size of the mapped field
label size_;
//- Does the mapping engine leave anything unmapped?
bool hasUnmapped_;
public:
// Constructors
//- Construct given a patch-to-patch and a mapping direction
patchToPatchFvPatchFieldMapper
(
const patchToPatch& pToP,
const bool forward
)
//- Construct given a patch-to-patch interpolation
patchToPatchFvPatchFieldMapper(const patchToPatch& pToP)
:
fvPatchFieldMapper(),
pToP_(pToP),
forward_(forward),
size_(-1),
hasUnmapped_(false)
{
const List<List<remote>> procFaces =
forward_ ? pToP_.tgtSrcProcFaces() : pToP_.srcTgtProcFaces();
size_ = procFaces.size();
forAll(procFaces, i)
{
if (procFaces[i].size() == 0)
{
hasUnmapped_ = true;
}
}
}
pToP_(pToP)
{}
//- Destructor
@ -104,18 +77,9 @@ public:
{}
// Member Functions
//- Return whether or not all faces receive a mapped value
virtual bool hasUnmapped() const
{
return hasUnmapped_;
}
// Member Operators
//- Map a label field
//- Map a label field (not implemented)
DEFINE_FIELD_MAPPER_OPERATOR(label, );
//- Map a temporary field
@ -128,6 +92,121 @@ public:
};
/*---------------------------------------------------------------------------*\
Class patchToPatchLeftOverFvPatchFieldMapper Declaration
\*---------------------------------------------------------------------------*/
class patchToPatchLeftOverFvPatchFieldMapper
:
public patchToPatchFvPatchFieldMapper
{
// Private Member Functions
//- Map from one field to another
template<class Type>
void map(Field<Type>& f, const Field<Type>& mapF) const;
//- Map a field and return the result as tmp
template<class Type>
tmp<Field<Type>> map(const Field<Type>& mapF) const;
public:
// Constructors
//- Inherit base class constructor
using patchToPatchFvPatchFieldMapper::patchToPatchFvPatchFieldMapper;
//- Destructor
virtual ~patchToPatchLeftOverFvPatchFieldMapper()
{}
// Member Functions
//- Return whether or not all faces receive a mapped value
virtual bool hasUnmapped() const
{
FatalErrorInFunction
<< "Not a valid query for this mapper, which should only be "
<< "used for modifying an existing, valid, field"
<< exit(FatalError);
return false;
}
// Member Operators
//- Map a field
FOR_ALL_FIELD_TYPES(DEFINE_FIELD_MAPPER_OPERATOR, );
};
/*---------------------------------------------------------------------------*\
Class patchToPatchNormalisedFvPatchFieldMapper Declaration
\*---------------------------------------------------------------------------*/
class patchToPatchNormalisedFvPatchFieldMapper
:
public patchToPatchFvPatchFieldMapper
{
// Private Data
//- Patch stabilisation engine
const patchToPatchStabilisation& pS_;
// Private Member Functions
//- Map from one field to another
template<class Type>
void map(Field<Type>& f, const Field<Type>& mapF) const;
//- Map a field and return the result as tmp
template<class Type>
tmp<Field<Type>> map(const Field<Type>& mapF) const;
public:
// Constructors
//- Construct given a patch-to-patch interpolation and stabilisation
patchToPatchNormalisedFvPatchFieldMapper
(
const patchToPatch& pToP,
const patchToPatchStabilisation& pS
)
:
patchToPatchFvPatchFieldMapper(pToP),
pS_(pS)
{}
//- Destructor
virtual ~patchToPatchNormalisedFvPatchFieldMapper()
{}
// Member Functions
//- Return whether or not all faces receive a mapped value
virtual bool hasUnmapped() const
{
return false;
}
// Member Operators
//- Map a field
FOR_ALL_FIELD_TYPES(DEFINE_FIELD_MAPPER_OPERATOR, );
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam

View File

@ -52,4 +52,56 @@ Foam::patchToPatchFvPatchFieldMapper::operator()
}
template<class Type>
void Foam::patchToPatchLeftOverFvPatchFieldMapper::map
(
Field<Type>& f,
const Field<Type>& mapF
) const
{
f = pToP_.srcToTgt(mapF, f);
}
template<class Type>
Foam::tmp<Foam::Field<Type>>
Foam::patchToPatchLeftOverFvPatchFieldMapper::map
(
const Field<Type>& mapF
) const
{
FatalErrorInFunction
<< "Not a valid operation for this mapper, which should only be "
<< "used for modifying an existing, valid, field"
<< exit(FatalError);
return tmp<Field<Type>>(nullptr);
}
template<class Type>
void Foam::patchToPatchNormalisedFvPatchFieldMapper::map
(
Field<Type>& f,
const Field<Type>& mapF
) const
{
f = pToP_.srcToTgt(mapF);
pS_.stabilise(f);
}
template<class Type>
Foam::tmp<Foam::Field<Type>>
Foam::patchToPatchNormalisedFvPatchFieldMapper::map
(
const Field<Type>& mapF
) const
{
tmp<Field<Type>> tf(new Field<Type>());
map(tf.ref(), mapF);
return tf;
}
// ************************************************************************* //

View File

@ -1,53 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2022-2023 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 "patchToPatchLeftOverFvPatchFieldMapper.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
template<class Type>
void Foam::patchToPatchLeftOverFvPatchFieldMapper::map
(
Field<Type>& f,
const Field<Type>& mapF
) const
{
f = forward_ ? pToP_.srcToTgt(mapF, f) : pToP_.tgtToSrc(mapF, f);
}
template<class Type>
Foam::tmp<Foam::Field<Type>>
Foam::patchToPatchLeftOverFvPatchFieldMapper::map
(
const Field<Type>& mapF
) const
{
NotImplemented;
return tmp<Field<Type>>(nullptr);
}
// ************************************************************************* //

View File

@ -1,53 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2022-2023 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 "patchToPatchNormalisedFvPatchFieldMapper.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
template<class Type>
void Foam::patchToPatchNormalisedFvPatchFieldMapper::map
(
Field<Type>& f,
const Field<Type>& mapF
) const
{
f = forward_ ? pToP_.srcToTgt(mapF) : pToP_.tgtToSrc(mapF);
}
template<class Type>
Foam::tmp<Foam::Field<Type>> Foam::patchToPatchNormalisedFvPatchFieldMapper::map
(
const Field<Type>& mapF
) const
{
tmp<Field<Type>> tf(new Field<Type>(size_));
map(tf.ref(), mapF);
return tf;
}
// ************************************************************************* //

View File

@ -215,12 +215,16 @@ patchToPatch/inverseDistance/inverseDistancePatchToPatch.C
patchToPatch/intersection/intersectionPatchToPatch.C
patchToPatch/rays/raysPatchToPatch.C
patchToPatch/patchToPatchStabilisation/patchToPatchStabilisation.C
cellsToCells/cellsToCells/cellsToCells.C
cellsToCells/cellsToCells/cellsToCellsParallelOps.C
cellsToCells/matching/matchingCellsToCells.C
cellsToCells/nearest/nearestCellsToCells.C
cellsToCells/intersection/intersectionCellsToCells.C
cellsToCells/cellsToCellsStabilisation/cellsToCellsStabilisation.C
meshToMesh/meshToMesh.C
mappedPatches/mappedPatchBase/mappedPatchBase.C

View File

@ -176,8 +176,8 @@ inline bool Foam::PatchEdgeFacePointData<Type>::updateFace
(
mesh,
patch,
edgei,
facei,
edgei,
edgeInfo,
tol,
td

View File

@ -183,6 +183,28 @@ Foam::autoPtr<Foam::cellsToCells> Foam::cellsToCells::New
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::PackedBoolList Foam::cellsToCells::srcCoupled() const
{
PackedBoolList result(srcLocalTgtCells_.size());
forAll(srcLocalTgtCells_, srcCelli)
{
result[srcCelli] = !srcLocalTgtCells_[srcCelli].empty();
}
return result;
}
Foam::PackedBoolList Foam::cellsToCells::tgtCoupled() const
{
PackedBoolList result(tgtLocalSrcCells_.size());
forAll(tgtLocalSrcCells_, tgtCelli)
{
result[tgtCelli] = !tgtLocalSrcCells_[tgtCelli].empty();
}
return result;
}
Foam::remote Foam::cellsToCells::srcToTgtPoint
(
const polyMesh& tgtMesh,

View File

@ -209,6 +209,12 @@ public:
//- Is this intersection on a single process?
inline bool isSingleProcess() const;
//- Return a list indicating which source cells are coupled
PackedBoolList srcCoupled() const;
//- Return a list indicating which target cells are coupled
PackedBoolList tgtCoupled() const;
// Interpolation

View File

@ -0,0 +1,240 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2023 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 "cellsToCellsStabilisation.H"
#include "syncTools.H"
#include "wallPoint.H"
#include "WallLocationData.H"
#include "WallInfo.H"
#include "FaceCellWave.H"
#include "globalIndex.H"
#include "OBJstream.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(cellsToCellsStabilisation, 0);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::cellsToCellsStabilisation::cellsToCellsStabilisation()
:
stabilisation_(false),
localStabilisationCells_(),
stabilisationMapPtr_(nullptr)
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::cellsToCellsStabilisation::~cellsToCellsStabilisation()
{}
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
void Foam::cellsToCellsStabilisation::update
(
const polyMesh& mesh,
const PackedBoolList& cellCoupleds
)
{
// Determine whether or not stabilisation is necessary
stabilisation_ = false;
forAll(cellCoupleds, celli)
{
if (!cellCoupleds[celli])
{
stabilisation_ = true;
break;
}
}
reduce(stabilisation_, orOp<bool>());
// Quick return if nothing is to be done
if (!stabilisation_) return;
// Get some information regarding the cells on the other side of couplings
List<remote> bFaceNbrProcCells(mesh.nFaces() - mesh.nInternalFaces());
boolList bFaceNbrIsCoupled(mesh.nFaces() - mesh.nInternalFaces());
forAll(bFaceNbrIsCoupled, bFacei)
{
const label owni = mesh.faceOwner()[bFacei + mesh.nInternalFaces()];
bFaceNbrProcCells[bFacei] = remote(Pstream::myProcNo(), owni);
bFaceNbrIsCoupled[bFacei] = cellCoupleds[owni];
}
syncTools::swapBoundaryFaceList(mesh, bFaceNbrProcCells);
syncTools::swapBoundaryFaceList(mesh, bFaceNbrIsCoupled);
// Determine the "cut" faces that separate coupled and un-coupled cellS
typedef WallInfo<WallLocationData<wallPoint, remote>> info;
DynamicList<label> cutFaces;
DynamicList<info> cutFaceInfos;
for (label facei = 0; facei < mesh.nFaces(); ++ facei)
{
const label owni = mesh.faceOwner()[facei];
const bool ownIsCoupled = cellCoupleds[owni];
if (facei < mesh.nInternalFaces())
{
const label nbri = mesh.faceNeighbour()[facei];
const bool nbrIsCoupled = cellCoupleds[nbri];
if (ownIsCoupled != nbrIsCoupled)
{
const label celli = ownIsCoupled ? owni : nbri;
cutFaces.append(facei);
cutFaceInfos.append
(
info
(
remote(Pstream::myProcNo(), celli),
mesh.faceCentres()[facei],
0
)
);
}
}
else
{
const label bFacei = facei - mesh.nInternalFaces();
const bool nbrIsCoupled = bFaceNbrIsCoupled[bFacei];
if (!ownIsCoupled && nbrIsCoupled)
{
cutFaces.append(facei);
cutFaceInfos.append
(
info
(
bFaceNbrProcCells[bFacei],
mesh.faceCentres()[facei],
0
)
);
}
}
}
// Wave the information about the cut faces' connected coupled cells into
// the un-coupled cells. Base this wave on distance to the cut face.
// Initialise coupled cells to have a distance of zero, so that we do not
// waste time waving into coupled regions of the mesh.
List<info> faceInfos(mesh.nFaces()), cellInfos(mesh.nCells());
forAll(cellCoupleds, celli)
{
if (cellCoupleds[celli])
{
cellInfos[celli] =
info
(
remote(Pstream::myProcNo(), celli),
mesh.cellCentres()[celli],
0
);
}
}
FaceCellWave<info> wave
(
mesh,
cutFaces,
cutFaceInfos,
faceInfos,
cellInfos,
mesh.globalData().nTotalCells() + 1 // max iterations
);
// Check that the wave connected to all un-coupled cells
forAll(cellCoupleds, celli)
{
if (!cellCoupleds[celli] && !cellInfos[celli].valid(wave.data()))
{
FatalErrorInFunction
<< "Un-coupled cell " << celli << " of mesh " << mesh.name()
<< " on processor " << Pstream::myProcNo() << " with centre "
<< "at " << mesh.cellCentres()[celli] << " was not connected "
<< "to a coupled cell by the stabilisation wave. This "
<< "indicates that an entire non-contiguous region of mesh "
<< "lies outside of the other mesh being mapped to. This is "
<< "not recoverable." << exit(FatalError);
}
}
// Construct the cell to local stabilisation cell map
const globalIndex cellGlobalIndex(mesh.nCells());
localStabilisationCells_.resize(mesh.nCells());
forAll(cellCoupleds, celli)
{
const remote& r = cellInfos[celli].data();
localStabilisationCells_[celli] =
cellCoupleds[celli]
? cellGlobalIndex.toGlobal(celli)
: cellGlobalIndex.toGlobal(r.proci, r.elementi);
}
// Construct the distribution map, if necessary
if (Pstream::parRun())
{
List<Map<label>> compactMap;
stabilisationMapPtr_.reset
(
new distributionMap
(
cellGlobalIndex,
localStabilisationCells_,
compactMap
)
);
}
// Write out stabilisation connections
if (debug)
{
OBJstream obj
(
typeName + "_" + mesh.name()
+ (Pstream::parRun() ? "_proc" + name(Pstream::myProcNo()) : "")
+ "_connections.obj"
);
pointField ccs(mesh.cellCentres());
stabilise(ccs);
forAll(ccs, celli)
{
const point& c = mesh.cellCentres()[celli];
if (magSqr(c - ccs[celli]) == 0) continue;
obj.write(linePointRef(ccs[celli], c));
}
}
}
// ************************************************************************* //

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2022-2023 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2023 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -22,19 +22,21 @@ License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::patchToPatchLeftOverFvPatchFieldMapper
Foam::cellsToCellsStabilisation
Description
FieldMapper which uses a patchToPatch object to map from another patch. The
source patch may be differently decomposed and/or geometrically and
topologically different from the target.
Stabilisation data and routines for cell-to-cell interpolations
SourceFiles
cellsToCellsStabilisation.C
cellsToCellsStabilisationTemplates.C
\*---------------------------------------------------------------------------*/
#ifndef patchToPatchLeftOverFvPatchFieldMapper_H
#define patchToPatchLeftOverFvPatchFieldMapper_H
#ifndef cellsToCellsStabilisation_H
#define cellsToCellsStabilisation_H
#include "patchToPatchFvPatchFieldMapper.H"
#include "cellsToCells.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -42,48 +44,51 @@ namespace Foam
{
/*---------------------------------------------------------------------------*\
Class patchToPatchLeftOverFvPatchFieldMapper Declaration
Class cellsToCellsStabilisation Declaration
\*---------------------------------------------------------------------------*/
class patchToPatchLeftOverFvPatchFieldMapper
:
public patchToPatchFvPatchFieldMapper
class cellsToCellsStabilisation
{
// Private Member Functions
// Private Data
//- Map from one field to another
template<class Type>
void map(Field<Type>& f, const Field<Type>& mapF) const;
//- Is stabilisation occurring?
bool stabilisation_;
//- Map a field and return the result as tmp
template<class Type>
tmp<Field<Type>> map(const Field<Type>& mapF) const;
//- For each cell, the local cell used to stabilise the interpolation
labelList localStabilisationCells_;
//- Map from cells to local stabilisation cells
autoPtr<distributionMap> stabilisationMapPtr_;
public:
//- Run-time type information
TypeName("cellsToCellsStabilisation");
// Constructors
//- Construct given a patch-to-patch and a mapping direction
patchToPatchLeftOverFvPatchFieldMapper
(
const patchToPatch& pToP,
const bool forward
)
:
patchToPatchFvPatchFieldMapper(pToP, forward)
{}
//- Construct null
cellsToCellsStabilisation();
//- Destructor
virtual ~patchToPatchLeftOverFvPatchFieldMapper()
{}
virtual ~cellsToCellsStabilisation();
// Member Operators
// Member Functions
//- Map a field
FOR_ALL_FIELD_TYPES(DEFINE_FIELD_MAPPER_OPERATOR, );
//- Compute the stabilisation addressing if necessary
void update
(
const polyMesh& mesh,
const PackedBoolList& cellCoupleds
);
//- Stabilise the given field if necessary
template<class Type>
void stabilise(Field<Type>& fld) const;
};
@ -94,7 +99,7 @@ public:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
#include "patchToPatchLeftOverFvPatchFieldMapperTemplates.C"
#include "cellsToCellsStabilisationTemplates.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2022-2023 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2023 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -23,15 +23,23 @@ License
\*---------------------------------------------------------------------------*/
#include "patchToPatchLeftOverFvPatchFieldMapper.H"
#include "cellsToCellsStabilisation.H"
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
FOR_ALL_FIELD_TYPES
(
IMPLEMENT_FIELD_MAPPER_OPERATOR,
patchToPatchLeftOverFvPatchFieldMapper
)
template<class Type>
void Foam::cellsToCellsStabilisation::stabilise(Field<Type>& fld) const
{
if (stabilisation_)
{
if (Pstream::parRun())
{
stabilisationMapPtr_->distribute(fld);
}
fld.map(fld, localStabilisationCells_);
}
}
// ************************************************************************* //

View File

@ -49,15 +49,19 @@ Foam::meshToMesh::meshToMesh
:
srcMesh_(srcMesh),
tgtMesh_(tgtMesh),
srcToTgtCellsToCells_(),
srcToTgtPatchIDs_(),
srcToTgtPatchToPatches_()
cellsInterpolation_(),
srcCellsStabilisation_(),
tgtCellsStabilisation_(),
patchIDs_(),
patchInterpolations_(),
srcPatchStabilisations_(),
tgtPatchStabilisations_()
{
// If no patch map was supplied, then assume a consistent pair of meshes in
// which corresponding patches have the same name
if (isNull(patchMap))
{
DynamicList<labelPair> srcToTgtPatchIDs;
DynamicList<labelPair> patchIDs;
forAll(srcMesh_.boundaryMesh(), srcPatchi)
{
@ -89,17 +93,17 @@ Foam::meshToMesh::meshToMesh
<< exit(FatalError);
}
srcToTgtPatchIDs.append(labelPair(srcPatchi, tgtPatchi));
patchIDs.append(labelPair(srcPatchi, tgtPatchi));
}
}
srcToTgtPatchIDs_.transfer(srcToTgtPatchIDs);
patchIDs_.transfer(patchIDs);
}
// If a patch mas was supplied then convert it to pairs of patch indices
else
{
srcToTgtPatchIDs_.setSize(patchMap.size());
patchIDs_.setSize(patchMap.size());
label i = 0;
forAllConstIter(HashTable<word>, patchMap, iter)
{
@ -130,7 +134,7 @@ Foam::meshToMesh::meshToMesh
<< exit(FatalError);
}
srcToTgtPatchIDs_[i ++] = labelPair(srcPatchi, tgtPatchi);
patchIDs_[i ++] = labelPair(srcPatchi, tgtPatchi);
}
}
@ -139,36 +143,41 @@ Foam::meshToMesh::meshToMesh
<< srcMesh_.name() << " and target mesh " << tgtMesh_.name()
<< " using " << engineType << endl << incrIndent;
srcToTgtCellsToCells_ = cellsToCells::New(engineType);
srcToTgtCellsToCells_->update(srcMesh_, tgtMesh_);
cellsInterpolation_ = cellsToCells::New(engineType);
cellsInterpolation_->update(srcMesh_, tgtMesh_);
srcCellsStabilisation_.clear();
tgtCellsStabilisation_.clear();
Info<< decrIndent;
// Calculate patch addressing and weights
srcToTgtPatchToPatches_.setSize(srcToTgtPatchIDs_.size());
forAll(srcToTgtPatchIDs_, i)
patchInterpolations_.setSize(patchIDs_.size());
srcPatchStabilisations_.setSize(patchIDs_.size());
tgtPatchStabilisations_.setSize(patchIDs_.size());
forAll(patchIDs_, i)
{
const label srcPatchi = srcToTgtPatchIDs_[i].first();
const label tgtPatchi = srcToTgtPatchIDs_[i].second();
const label srcPatchi = patchIDs_[i].first();
const label tgtPatchi = patchIDs_[i].second();
const polyPatch& srcPP = srcMesh_.boundaryMesh()[srcPatchi];
const polyPatch& tgtPP = tgtMesh_.boundaryMesh()[tgtPatchi];
const polyPatch& srcPp = srcMesh_.boundaryMesh()[srcPatchi];
const polyPatch& tgtPp = tgtMesh_.boundaryMesh()[tgtPatchi];
Info<< "Creating patchToPatch between source patch "
<< srcPP.name() << " and target patch " << tgtPP.name()
<< srcPp.name() << " and target patch " << tgtPp.name()
<< " using " << engineType << endl << incrIndent;
srcToTgtPatchToPatches_.set
patchInterpolations_.set
(
i,
patchToPatch::New(engineType, true)
);
srcToTgtPatchToPatches_[i].update
patchInterpolations_[i].update
(
srcPP,
PatchTools::pointNormals(srcMesh_, srcPP),
tgtPP
srcPp,
PatchTools::pointNormals(srcMesh_, srcPp),
tgtPp
);
Info<< decrIndent;
@ -190,10 +199,10 @@ bool Foam::meshToMesh::consistent() const
boolList tgtPatchIsMapped(tgtMesh_.boundaryMesh().size(), false);
// Mark anything paired as mapped
forAll(srcToTgtPatchIDs_, i)
forAll(patchIDs_, i)
{
const label srcPatchi = srcToTgtPatchIDs_[i].first();
const label tgtPatchi = srcToTgtPatchIDs_[i].second();
const label srcPatchi = patchIDs_[i].first();
const label tgtPatchi = patchIDs_[i].second();
srcPatchIsMapped[srcPatchi] = true;
tgtPatchIsMapped[tgtPatchi] = true;
@ -242,7 +251,7 @@ Foam::remote Foam::meshToMesh::srcToTgtPoint
const point& p
) const
{
return srcToTgtCellsToCells_().srcToTgtPoint(tgtMesh_, srcCelli, p);
return cellsInterpolation_().srcToTgtPoint(tgtMesh_, srcCelli, p);
}

View File

@ -37,8 +37,10 @@ SourceFiles
#define meshToMesh_H
#include "polyMesh.H"
#include "patchToPatch.H"
#include "cellsToCells.H"
#include "cellsToCellsStabilisation.H"
#include "patchToPatch.H"
#include "patchToPatchStabilisation.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -60,15 +62,25 @@ class meshToMesh
const polyMesh& tgtMesh_;
//- Interpolation engine between source and target cells
autoPtr<cellsToCells> srcToTgtCellsToCells_;
autoPtr<cellsToCells> cellsInterpolation_;
//- List of corresponding source and target patches that are to be
// mapped to each other
List<labelPair> srcToTgtPatchIDs_;
//- Stabilisation engine for the source cells
mutable autoPtr<cellsToCellsStabilisation> srcCellsStabilisation_;
//- List of patchToPatch interpolation engines between source and
// target patches
PtrList<patchToPatch> srcToTgtPatchToPatches_;
//- Stabilisation engine for the target cells
mutable autoPtr<cellsToCellsStabilisation> tgtCellsStabilisation_;
//- List of corresponding source and target patch indices
List<labelPair> patchIDs_;
//- List of interpolation engines between source and target patches
PtrList<patchToPatch> patchInterpolations_;
//- List of stabilisation engines for the source patches
mutable PtrList<patchToPatchStabilisation> srcPatchStabilisations_;
//- List of stabilisation engines for the source patches
mutable PtrList<patchToPatchStabilisation> tgtPatchStabilisations_;
protected:
@ -78,14 +90,37 @@ protected:
// Access
//- Return the interpolation engine between source and target cells
inline const cellsToCells& srcToTgtCellsToCells() const;
inline const cellsToCells& cellsInterpolation() const;
//- Return the source patch indices
inline const List<labelPair>& srcToTgtPatchIDs() const;
//- Return the stabilisation engine for the source cells
inline const cellsToCellsStabilisation&
srcCellsStabilisation() const;
//- Return the list of patchToPatch interpolation engines between
// source and target patches
inline const PtrList<patchToPatch>& srcToTgtPatchToPatches() const;
//- Return the stabilisation engine for the target cells
inline const cellsToCellsStabilisation&
tgtCellsStabilisation() const;
//- Return the list of corresponding source and target patch indices
inline const List<labelPair>& patchIDs() const;
//- Return the interpolation engine between a source and a target
// patch
inline const patchToPatch& patchInterpolation
(
const label i
) const;
//- Return the stabilisation engine for a source patch
inline const patchToPatchStabilisation& srcPatchStabilisation
(
const label i
) const;
//- Return the stabilisation engine for a target patch
inline const patchToPatchStabilisation& tgtPatchStabilisation
(
const label i
) const;
public:
@ -97,8 +132,9 @@ public:
// Constructors
//- Construct from source and target meshes. If a patchMap is supplied,
// then map between the specified patches. If not, then assume a
// consistent mesh and map all patches 1-to-1.
// then interpolate between the specified patches. If not, then assume
// a consistent mesh with consistently named patches and interpolate
// 1-to-1 between patches with the same name.
meshToMesh
(
const polyMesh& srcMesh,
@ -125,8 +161,8 @@ public:
//- Return const access to the target mesh
inline const polyMesh& tgtMesh() const;
//- Is the map between consistent meshes? I.e., are all
// (non-processor) patches mapped one-to-one?
//- Is the interpolation between consistent meshes? I.e., are all
// (non-processor) patches coupled one-to-one?
bool consistent() const;

View File

@ -27,23 +27,98 @@ License
// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
inline const Foam::cellsToCells& Foam::meshToMesh::srcToTgtCellsToCells() const
inline const Foam::cellsToCells& Foam::meshToMesh::cellsInterpolation() const
{
return srcToTgtCellsToCells_();
return cellsInterpolation_();
}
inline const Foam::List<Foam::labelPair>&
Foam::meshToMesh::srcToTgtPatchIDs() const
inline const Foam::cellsToCellsStabilisation&
Foam::meshToMesh::srcCellsStabilisation() const
{
return srcToTgtPatchIDs_;
if (!srcCellsStabilisation_.valid())
{
srcCellsStabilisation_.reset(new cellsToCellsStabilisation());
srcCellsStabilisation_->update
(
srcMesh_,
cellsInterpolation_->srcCoupled()
);
}
return srcCellsStabilisation_();
}
inline const Foam::PtrList<Foam::patchToPatch>&
Foam::meshToMesh::srcToTgtPatchToPatches() const
inline const Foam::cellsToCellsStabilisation&
Foam::meshToMesh::tgtCellsStabilisation() const
{
return srcToTgtPatchToPatches_;
if (!tgtCellsStabilisation_.valid())
{
tgtCellsStabilisation_.reset(new cellsToCellsStabilisation());
tgtCellsStabilisation_->update
(
tgtMesh_,
cellsInterpolation_->tgtCoupled()
);
}
return tgtCellsStabilisation_();
}
inline const Foam::List<Foam::labelPair>& Foam::meshToMesh::patchIDs() const
{
return patchIDs_;
}
inline const Foam::patchToPatch& Foam::meshToMesh::patchInterpolation
(
const label i
) const
{
return patchInterpolations_[i];
}
inline const Foam::patchToPatchStabilisation&
Foam::meshToMesh::srcPatchStabilisation(const label i) const
{
if (!srcPatchStabilisations_.set(i))
{
const label srcPatchi = patchIDs_[i].first();
const polyPatch& srcPp = srcMesh_.boundaryMesh()[srcPatchi];
srcPatchStabilisations_.set(i, new patchToPatchStabilisation());
srcPatchStabilisations_[i].update
(
srcPp,
patchInterpolations_[i].srcCoupled()
);
}
return srcPatchStabilisations_[i];
}
inline const Foam::patchToPatchStabilisation&
Foam::meshToMesh::tgtPatchStabilisation(const label i) const
{
if (!tgtPatchStabilisations_.set(i))
{
const label tgtPatchi = patchIDs_[i].first();
const polyPatch& tgtPp = tgtMesh_.boundaryMesh()[tgtPatchi];
tgtPatchStabilisations_.set(i, new patchToPatchStabilisation());
tgtPatchStabilisations_[i].update
(
tgtPp,
patchInterpolations_[i].tgtCoupled()
);
}
return tgtPatchStabilisations_[i];
}

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -807,6 +807,56 @@ Foam::autoPtr<Foam::patchToPatch> Foam::patchToPatch::New
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::PackedBoolList Foam::patchToPatch::srcCoupled() const
{
PackedBoolList result(srcLocalTgtFaces_.size());
forAll(srcLocalTgtFaces_, srcFacei)
{
result[srcFacei] = !srcLocalTgtFaces_[srcFacei].empty();
}
return result;
}
Foam::PackedBoolList Foam::patchToPatch::tgtCoupled() const
{
PackedBoolList result(tgtLocalSrcFaces_.size());
forAll(tgtLocalSrcFaces_, tgtFacei)
{
result[tgtFacei] = !tgtLocalSrcFaces_[tgtFacei].empty();
}
return result;
}
Foam::List<Foam::List<Foam::remote>>
Foam::patchToPatch::srcTgtProcFaces() const
{
return
isSingleProcess()
? patchToPatchTools::localToRemote(srcLocalTgtFaces_)
: patchToPatchTools::localToRemote
(
srcLocalTgtFaces_,
localTgtProcFacesPtr_()
);
}
Foam::List<Foam::List<Foam::remote>>
Foam::patchToPatch::tgtSrcProcFaces() const
{
return
isSingleProcess()
? patchToPatchTools::localToRemote(tgtLocalSrcFaces_)
: patchToPatchTools::localToRemote
(
tgtLocalSrcFaces_,
localSrcProcFacesPtr_()
);
}
void Foam::patchToPatch::update
(
const primitiveOldTimePatch& srcPatch,

View File

@ -311,11 +311,17 @@ public:
//- Is this intersection on a single process?
inline bool isSingleProcess() const;
//- Return a list indicating which source faces are coupled
PackedBoolList srcCoupled() const;
//- Return a list indicating which target faces are coupled
PackedBoolList tgtCoupled() const;
//- For each source face, the coupled target procs and faces
inline List<List<remote>> srcTgtProcFaces() const;
List<List<remote>> srcTgtProcFaces() const;
//- For each target face, the coupled source procs and faces
inline List<List<remote>> tgtSrcProcFaces() const;
List<List<remote>> tgtSrcProcFaces() const;
// Interpolation

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2022 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -46,32 +46,4 @@ inline bool Foam::patchToPatch::isSingleProcess() const
}
inline Foam::List<Foam::List<Foam::remote>>
Foam::patchToPatch::srcTgtProcFaces() const
{
return
isSingleProcess()
? patchToPatchTools::localToRemote(srcLocalTgtFaces_)
: patchToPatchTools::localToRemote
(
srcLocalTgtFaces_,
localTgtProcFacesPtr_()
);
}
inline Foam::List<Foam::List<Foam::remote>>
Foam::patchToPatch::tgtSrcProcFaces() const
{
return
isSingleProcess()
? patchToPatchTools::localToRemote(tgtLocalSrcFaces_)
: patchToPatchTools::localToRemote
(
tgtLocalSrcFaces_,
localSrcProcFacesPtr_()
);
}
// ************************************************************************* //

View File

@ -0,0 +1,207 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2023 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 "patchToPatchStabilisation.H"
#include "PatchEdgeFacePointData.H"
#include "PatchEdgeFaceWave.H"
#include "SubField.H"
#include "globalIndex.H"
#include "OBJstream.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(patchToPatchStabilisation, 0);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::patchToPatchStabilisation::patchToPatchStabilisation()
:
stabilisation_(false),
localStabilisationCells_(),
stabilisationMapPtr_(nullptr)
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::patchToPatchStabilisation::~patchToPatchStabilisation()
{}
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
void Foam::patchToPatchStabilisation::update
(
const polyPatch& patch,
const PackedBoolList& faceCoupleds
)
{
// Determine whether or not stabilisation is necessary
stabilisation_ = false;
forAll(faceCoupleds, facei)
{
if (!faceCoupleds[facei])
{
stabilisation_ = true;
break;
}
}
reduce(stabilisation_, orOp<bool>());
// Quick return if nothing is to be done
if (!stabilisation_) return;
// Construct initial edges. All edges that border a coupled face are added
// here. The wave will propagate everywhere for just the first iteration.
// Then most paths will end and subsequent iterations will propagate only
// through the uncoupled faces. This is a bit odd, but it is easier than
// doing the necessary synchronisation to determine which edges lie
// in-between coupled and non-coupled faces.
typedef PatchEdgeFacePointData<remote> info;
DynamicList<label> initialEdges(patch.nEdges());
DynamicList<info> initialEdgeInfos(patch.nEdges());
forAll(patch.edgeFaces(), edgei)
{
forAll(patch.edgeFaces()[edgei], edgeFacei)
{
const label facei = patch.edgeFaces()[edgei][edgeFacei];
if (faceCoupleds[facei])
{
initialEdges.append(edgei);
initialEdgeInfos.append
(
info
(
remote(Pstream::myProcNo(), facei),
patch.edges()[edgei].centre(patch.localPoints()),
0
)
);
break;
}
}
}
// Wave the information about the nearby coupled faces into the un-coupled
// faces. Base this wave on distance to the cut face. Initialise coupled
// faces to have a distance of zero, so that we do not waste time waving
// into coupled regions of the patch.
List<info> edgeInfos(patch.nEdges()), faceInfos(patch.size());
forAll(faceCoupleds, facei)
{
if (faceCoupleds[facei])
{
faceInfos[facei] =
info
(
remote(Pstream::myProcNo(), facei),
patch.faceCentres()[facei],
0
);
}
}
PatchEdgeFaceWave<primitivePatch, info> wave
(
patch.boundaryMesh().mesh(),
patch,
initialEdges,
initialEdgeInfos,
edgeInfos,
faceInfos,
returnReduce(patch.nEdges(), sumOp<label>())
);
// Check that the wave connected to all un-mapped faces
forAll(faceCoupleds, facei)
{
if (!faceCoupleds[facei] && !faceInfos[facei].valid(wave.data()))
{
FatalErrorInFunction
<< "Un-mapped face " << facei << " of patch " << patch.name()
<< " on processor " << Pstream::myProcNo() << " with centre "
<< "at " << patch.faceCentres()[facei] << " was not connected "
<< "to a mapped cell by the stabilisation wave. This "
<< "indicates that an entire non-contiguous region of patch "
<< "lies outside of the other patch being mapped to. This is "
<< "not recoverable." << exit(FatalError);
}
}
// Construct the cell to local stabilisation cell map
const globalIndex cellGlobalIndex(patch.size());
localStabilisationCells_.resize(patch.size());
forAll(faceCoupleds, facei)
{
const remote& r = faceInfos[facei].data();
localStabilisationCells_[facei] =
faceCoupleds[facei]
? cellGlobalIndex.toGlobal(facei)
: cellGlobalIndex.toGlobal(r.proci, r.elementi);
}
// Construct the distribution map, if necessary
if (Pstream::parRun())
{
List<Map<label>> compactMap;
stabilisationMapPtr_.reset
(
new distributionMap
(
cellGlobalIndex,
localStabilisationCells_,
compactMap
)
);
}
// Write out stabilisation connections
if (debug)
{
OBJstream obj
(
typeName + "_" + patch.name()
+ (Pstream::parRun() ? "_proc" + name(Pstream::myProcNo()) : "")
+ "_connections.obj"
);
pointField fcs(patch.faceCentres());
stabilise(fcs);
forAll(fcs, celli)
{
const point& c = patch.faceCentres()[celli];
if (magSqr(c - fcs[celli]) == 0) continue;
obj.write(linePointRef(fcs[celli], c));
}
}
}
// ************************************************************************* //

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2022-2023 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2023 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -22,19 +22,21 @@ License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::patchToPatchNormalisedFvPatchFieldMapper
Foam::patchToPatchStabilisation
Description
FieldMapper which uses a patchToPatch object to map from another patch. The
source patch may be differently decomposed and/or geometrically and
topologically different from the target.
Stabilisation data and routines for patch-to-patch interpolations
SourceFiles
patchToPatchStabilisation.C
patchToPatchStabilisationTemplates.C
\*---------------------------------------------------------------------------*/
#ifndef patchToPatchNormalisedFvPatchFieldMapper_H
#define patchToPatchNormalisedFvPatchFieldMapper_H
#ifndef patchToPatchStabilisation_H
#define patchToPatchStabilisation_H
#include "patchToPatchFvPatchFieldMapper.H"
#include "patchToPatch.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -42,48 +44,51 @@ namespace Foam
{
/*---------------------------------------------------------------------------*\
Class patchToPatchNormalisedFvPatchFieldMapper Declaration
Class patchToPatchStabilisation Declaration
\*---------------------------------------------------------------------------*/
class patchToPatchNormalisedFvPatchFieldMapper
:
public patchToPatchFvPatchFieldMapper
class patchToPatchStabilisation
{
// Private Member Functions
// Private Data
//- Map from one field to another
template<class Type>
void map(Field<Type>& f, const Field<Type>& mapF) const;
//- Is stabilisation occurring?
bool stabilisation_;
//- Map a field and return the result as tmp
template<class Type>
tmp<Field<Type>> map(const Field<Type>& mapF) const;
//- For each face, the local face used to stabilise the interpolation
labelList localStabilisationCells_;
//- Map from faces to local stabilisation faces
autoPtr<distributionMap> stabilisationMapPtr_;
public:
//- Run-time type information
TypeName("patchToPatchStabilisation");
// Constructors
//- Construct given a patch-to-patch and a mapping direction
patchToPatchNormalisedFvPatchFieldMapper
(
const patchToPatch& pToP,
const bool forward
)
:
patchToPatchFvPatchFieldMapper(pToP, forward)
{}
//- Construct null
patchToPatchStabilisation();
//- Destructor
virtual ~patchToPatchNormalisedFvPatchFieldMapper()
{}
virtual ~patchToPatchStabilisation();
// Member Operators
// Member Functions
//- Map a field
FOR_ALL_FIELD_TYPES(DEFINE_FIELD_MAPPER_OPERATOR, );
//- Compute the stabilisation addressing if necessary
void update
(
const polyPatch& patch,
const PackedBoolList& faceCoupleds
);
//- Stabilise the given field if necessary
template<class Type>
void stabilise(Field<Type>& fld) const;
};
@ -94,7 +99,7 @@ public:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
#include "patchToPatchNormalisedFvPatchFieldMapperTemplates.C"
#include "patchToPatchStabilisationTemplates.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2022-2023 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2023 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -23,15 +23,23 @@ License
\*---------------------------------------------------------------------------*/
#include "patchToPatchNormalisedFvPatchFieldMapper.H"
#include "patchToPatchStabilisation.H"
// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
FOR_ALL_FIELD_TYPES
(
IMPLEMENT_FIELD_MAPPER_OPERATOR,
patchToPatchNormalisedFvPatchFieldMapper
)
template<class Type>
void Foam::patchToPatchStabilisation::stabilise(Field<Type>& fld) const
{
if (stabilisation_)
{
if (Pstream::parRun())
{
stabilisationMapPtr_->distribute(fld);
}
fld.map(fld, localStabilisationCells_);
}
}
// ************************************************************************* //