Files
openfoam/src/finiteVolume/fvMesh/wallDist/wallDistAddressing/wallDistAddressing.C
2024-05-30 15:15:56 +00:00

504 lines
14 KiB
C

/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2023-2024 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 "FaceCellWave.H"
#include "wallDistAddressing.H"
#include "wallPointAddressing.H"
#include "addToRunTimeSelectionTable.H"
#include "surfaceFields.H"
#include "wallPolyPatch.H"
#include "patchDistMethod.H"
#include "OBJstream.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(wallDistAddressing, 0);
}
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
Foam::label Foam::wallDistAddressing::getValues
(
const FaceCellWave<wallPointAddressing>& wave,
const List<wallPointAddressing>& allCellInfo,
const List<wallPointAddressing>& allFaceInfo,
volScalarField& y
) const
{
// Extract volField. See patchWave::getValues
label nIllegal = 0;
forAll(allCellInfo, celli)
{
const scalar dist = allCellInfo[celli].distSqr();
if (allCellInfo[celli].valid(wave.data()))
{
y[celli] = Foam::sqrt(dist);
}
else
{
y[celli] = dist;
++nIllegal;
}
}
// Copy boundary values
auto& bfld = y.boundaryFieldRef();
for (auto& pfld : bfld)
{
scalarField patchField(pfld.size(), Zero);
forAll(pfld, patchFacei)
{
const label meshFacei = pfld.patch().start() + patchFacei;
const scalar dist = allFaceInfo[meshFacei].distSqr();
if (allFaceInfo[meshFacei].valid(wave.data()))
{
// Adding SMALL to avoid problems with /0 in the turbulence
// models
patchField[patchFacei] = Foam::sqrt(dist) + SMALL;
}
else
{
patchField[patchFacei] = dist;
++nIllegal;
}
}
pfld = patchField;
}
return nIllegal;
}
void Foam::wallDistAddressing::addItem
(
const label item,
const labelPair& data,
label& untransformi,
label& transformi,
labelPairList& transformedWallInfo
)
{
const auto& gt = mesh_.globalData().globalTransforms();
if (gt.transformIndex(data) == gt.nullTransformIndex())
{
// Store item and global index of untransformed data
const label proci = gt.processor(data);
const label index = gt.index(data);
untransformedItems_[untransformi] = item;
untransformedSlots_[untransformi++] =
globalWallsPtr_().toGlobal(proci, index);
}
else
{
// Store item and transformed data
transformedItems_[transformi] = item;
transformedWallInfo[transformi++] = data;
}
}
void Foam::wallDistAddressing::correct(volScalarField& y)
{
y = dimensionedScalar("yWall", dimLength, GREAT);
const auto& C = mesh_.C();
const auto& gt = mesh_.globalData().globalTransforms();
label nWalls = 0;
for (const label patchi : patchIDs_)
{
nWalls += C.boundaryField()[patchi].size();
}
globalWallsPtr_.reset(new globalIndex(nWalls));
const globalIndex& globalWalls = globalWallsPtr_();
DynamicList<label> seedFaces(nWalls);
DynamicList<wallPointAddressing> seedInfo(nWalls);
nWalls = 0;
for (const label patchi : patchIDs_)
{
const auto& fc = C.boundaryField()[patchi];
const auto& areaFraction = fc.patch().patch().areaFraction();
forAll(fc, patchFacei)
{
if
(
!areaFraction
|| (areaFraction()[patchFacei] > 0.5) // mostly wall
)
{
seedFaces.append(fc.patch().start()+patchFacei);
seedInfo.append
(
wallPointAddressing
(
fc[patchFacei],
gt.encode
(
Pstream::myProcNo(),
nWalls,
gt.nullTransformIndex()
),
scalar(0.0)
)
);
}
nWalls++;
}
}
List<wallPointAddressing> allCellInfo(mesh_.nCells());
// Initialise the passive data (since no mesh reference in the constructor
// of wallPointAddressing)
forAll(allCellInfo, celli)
{
allCellInfo[celli].data() = gt.encode
(
Pstream::myProcNo(),
celli,
gt.nullTransformIndex()
);
}
List<wallPointAddressing> allFaceInfo(mesh_.nFaces());
forAll(allFaceInfo, facei)
{
allFaceInfo[facei].data() = gt.encode
(
Pstream::myProcNo(),
facei,
gt.nullTransformIndex()
);
}
// Propagate information inwards
FaceCellWave<wallPointAddressing> wave
(
mesh_,
seedFaces,
seedInfo,
allFaceInfo,
allCellInfo,
mesh_.globalData().nTotalCells()+1
);
// Extract wall distance
// ~~~~~~~~~~~~~~~~~~~~~
getValues(wave, allCellInfo, allFaceInfo, y);
const labelHashSet patchSet(patchIDs_);
// Correct wall cells for true distance
// - store for corrected cells the nearest wall face
// - update distance field (y) for these cells
Map<label> cellToWallFace;
if (correctWalls_)
{
cellToWallFace.reserve(nWalls);
correctBoundaryFaceCells
(
patchSet,
y,
cellToWallFace
);
correctBoundaryPointCells
(
patchSet,
y,
cellToWallFace
);
}
// Make sure boundary values are up to date
y.correctBoundaryConditions();
// Extract all addressing
// ~~~~~~~~~~~~~~~~~~~~~~
// - untransformed : -local cell/boundary, -nearest (global) wall
// - transformed : -local cell/boundary, -nearest (global) wall
label untransformi = 0;
untransformedItems_.resize(mesh_.nCells()+mesh_.nBoundaryFaces());
untransformedSlots_.resize(untransformedItems_.size());
label transformi = 0;
transformedItems_.resize(mesh_.nCells()+mesh_.nBoundaryFaces());
labelPairList transformedWallInfo(transformedItems_.size());
for (label celli = 0; celli < mesh_.nCells(); celli++)
{
const auto cellFnd = cellToWallFace.find(celli);
if (cellFnd)
{
const label wallFacei = cellFnd();
const auto& data = allFaceInfo[wallFacei].data();
addItem(celli, data, untransformi, transformi, transformedWallInfo);
}
else if (allCellInfo[celli].valid(wave.data()))
{
const auto& data = allCellInfo[celli].data();
addItem(celli, data, untransformi, transformi, transformedWallInfo);
}
}
untransformedPatchStarts_.resize(mesh_.boundary().size()+1);
transformedPatchStarts_.resize(mesh_.boundary().size()+1);
for (const auto& pp : mesh_.boundary())
{
untransformedPatchStarts_[pp.index()] = untransformi;
transformedPatchStarts_[pp.index()] = transformi;
forAll(pp, patchFacei)
{
const label facei = pp.start()+patchFacei;
if (allFaceInfo[facei].valid(wave.data()))
{
const auto& data = allFaceInfo[facei].data();
addItem
(
facei,
data,
untransformi,
transformi,
transformedWallInfo
);
}
}
}
untransformedItems_.resize(untransformi);
untransformedSlots_.resize(untransformi);
untransformedPatchStarts_.back() = untransformi;
transformedItems_.resize(transformi);
transformedWallInfo.resize(transformi);
transformedPatchStarts_.back() = transformi;
if (debug)
{
Pout<< typeName
<< " : untransformed:" << untransformi
<< " transformed:" << transformi
<< " out of nWalls:" << nWalls
<< " untransformStart:" << flatOutput(untransformedPatchStarts_)
<< " transformStart:" << flatOutput(transformedPatchStarts_)
<< endl;
}
List<Map<label>> compactMap;
mapPtr_.reset
(
new mapDistribute
(
globalWalls, // globalIndex
untransformedSlots_, // untransformedElements
gt, // globalIndexAndTransform
transformedWallInfo, // transformedElements
transformedSlots_, // transformedIndices
compactMap
)
);
if (debug & 2)
{
// Use the obtained map to write nearest
{
OBJstream os(mesh_.polyMesh::path()/"nearest.obj");
Pout<< typeName
<< " : writing line from cell/face to nearest wall face to "
<< os.name() << endl;
// Seed all wall faces with centroid and override all other
// cells/faces with that of nearest wall
volVectorField wallCentre(mesh_.C());
this->map(wallCentre, mapDistribute::transformPosition());
// Seed all wall faces with normal and override all other
// cells/faces with that of nearest wall
volVectorField wallNormal
(
IOobject
(
"n" & patchTypeName_,
mesh_.time().timeName(),
mesh_,
IOobjectOption::NO_REGISTER
),
mesh_,
dimensionedVector(dimless, Zero),
patchDistMethod::patchTypes<vector>(mesh_, patchSet)
);
// Calculate normal
for (const label patchi : patchIDs_)
{
auto& pnf = wallNormal.boundaryFieldRef()[patchi];
pnf == pnf.patch().nf();
}
this->map(wallNormal, mapDistribute::transform());
forAll(wallCentre, celli)
{
const point& cc = mesh_.C()[celli];
const point& wallC = wallCentre[celli];
const vector& wallN = wallNormal[celli];
os.write(linePointRef(cc, wallC), wallN, wallN);
}
}
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::wallDistAddressing::wallDistAddressing
(
const fvMesh& mesh,
const bool correctWalls,
const label updateInterval
)
:
// Register as "wallDistAddressing"
MeshObject_type(mesh),
cellDistFuncs(mesh),
patchIDs_(mesh.boundaryMesh().findPatchIDs<wallPolyPatch>().sortedToc()),
patchTypeName_("wall"),
updateInterval_(updateInterval),
correctWalls_(correctWalls),
requireUpdate_(true),
y_
(
IOobject
(
"y" & patchTypeName_,
mesh.time().timeName(),
mesh.thisDb(),
IOobjectOption::NO_REGISTER
),
mesh,
dimensionedScalar("y", dimLength, GREAT) // Same as wallDistData
)
{
movePoints();
}
Foam::wallDistAddressing::wallDistAddressing
(
const word& patchTypeName,
const fvMesh& mesh,
const labelList& patchIDs,
const bool correctWalls,
const label updateInterval
)
:
MeshObject_type(patchTypeName, mesh),
cellDistFuncs(mesh),
patchIDs_(patchIDs),
patchTypeName_(patchTypeName),
updateInterval_(updateInterval),
correctWalls_(correctWalls),
requireUpdate_(true),
y_
(
IOobject
(
"y" & patchTypeName_,
mesh.time().timeName(),
mesh.thisDb(),
IOobjectOption::NO_REGISTER
),
mesh,
dimensionedScalar("y", dimLength, GREAT) // Same as wallDistData
)
{
movePoints();
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::wallDistAddressing::~wallDistAddressing()
{} // mapDistribute was forward declared
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::wallDistAddressing::movePoints()
{
if
(
(updateInterval_ > 0)
&& ((mesh_.time().timeIndex() % updateInterval_) == 0)
)
{
requireUpdate_ = true;
}
if (requireUpdate_)
{
DebugInfo<< "Updating wall distance" << endl;
requireUpdate_ = false;
correct(y_);
return true;
}
return false;
}
void Foam::wallDistAddressing::updateMesh(const mapPolyMesh& mpm)
{
// Force update if performing topology change
// Note: needed?
// - field would have been mapped, so if using updateInterval option (!= 1)
// live with error associated of not updating and use mapped values?
requireUpdate_ = true;
movePoints();
}
// ************************************************************************* //