Skip to content
Permalink
Browse files

polyMesh/particle: Optionally store old-cell centres

The moving-mesh tracking algorithm needs the cell-centres at the
previous time-step. These were not originally stored by the polyMesh, so
they were being generated on-the-fly. In some tracking-dominated cases
this caused an unacceptable overhead.

The polyMesh now stores the old-time cell-centres on demand. They are
not stored by default (like the old-time points), so if they are needed
then the accessor should be called before any mesh motion. Typically
this will be during construction of whatever functionality requires it.
See Cloud.C for an example.

The logic for storage and update of the old-time points has also been
improved to account for the possibility of the mesh motion coming to an
end.
  • Loading branch information...
Will Bainbridge
Will Bainbridge committed Mar 4, 2019
1 parent 8f0772d commit b5e27b7a9dd91d2daa45e9922734592ecbe4508b
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -289,8 +289,10 @@ Foam::polyMesh::polyMesh(const IOobject& io)
globalMeshDataPtr_(nullptr),
moving_(false),
topoChanging_(false),
curMotionTimeIndex_(time().timeIndex()),
oldPointsPtr_(nullptr)
curMotionTimeIndex_(-1),
oldPointsPtr_(nullptr),
oldCellCentresPtr_(nullptr),
storeOldCellCentres_(false)
{
if (!owner_.headerClassName().empty())
{
@@ -471,8 +473,10 @@ Foam::polyMesh::polyMesh
globalMeshDataPtr_(nullptr),
moving_(false),
topoChanging_(false),
curMotionTimeIndex_(time().timeIndex()),
oldPointsPtr_(nullptr)
curMotionTimeIndex_(-1),
oldPointsPtr_(nullptr),
oldCellCentresPtr_(nullptr),
storeOldCellCentres_(false)
{
// Check if the faces and cells are valid
forAll(faces_, facei)
@@ -622,8 +626,10 @@ Foam::polyMesh::polyMesh
globalMeshDataPtr_(nullptr),
moving_(false),
topoChanging_(false),
curMotionTimeIndex_(time().timeIndex()),
oldPointsPtr_(nullptr)
curMotionTimeIndex_(-1),
oldPointsPtr_(nullptr),
oldCellCentresPtr_(nullptr),
storeOldCellCentres_(false)
{
// Check if faces are valid
forAll(faces_, facei)
@@ -1052,22 +1058,42 @@ const Foam::labelList& Foam::polyMesh::faceNeighbour() const

const Foam::pointField& Foam::polyMesh::oldPoints() const
{
if (oldPointsPtr_.empty())
if (!moving_)
{
if (debug)
{
WarningInFunction
<< endl;
}
return points_;
}

oldPointsPtr_.reset(new pointField(points_));
curMotionTimeIndex_ = time().timeIndex();
if (oldPointsPtr_.empty())
{
FatalErrorInFunction
<< "Old points have not been stored"
<< exit(FatalError);
}

return oldPointsPtr_();
}


const Foam::pointField& Foam::polyMesh::oldCellCentres() const
{
storeOldCellCentres_ = true;

if (!moving_)
{
return cellCentres();
}

if (oldCellCentresPtr_.empty())
{
FatalErrorInFunction
<< "Old cell centres have not been stored"
<< exit(FatalError);
}

return oldCellCentresPtr_();
}


Foam::IOobject Foam::polyMesh::points0IO
(
const polyMesh& mesh
@@ -1149,12 +1175,16 @@ Foam::tmp<Foam::scalarField> Foam::polyMesh::movePoints

moving(true);

// Pick up old points
// Pick up old points and cell centres
if (curMotionTimeIndex_ != time().timeIndex())
{
// Mesh motion in the new time step
oldPointsPtr_.clear();
oldPointsPtr_.reset(new pointField(points_));
if (storeOldCellCentres_)
{
oldCellCentresPtr_.clear();
oldCellCentresPtr_.reset(new pointField(cellCentres()));
}
curMotionTimeIndex_ = time().timeIndex();
}

@@ -1233,8 +1263,9 @@ Foam::tmp<Foam::scalarField> Foam::polyMesh::movePoints

void Foam::polyMesh::resetMotion() const
{
curMotionTimeIndex_ = 0;
curMotionTimeIndex_ = -1;
oldPointsPtr_.clear();
oldCellCentresPtr_.clear();
}


@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -186,6 +186,12 @@ private:
//- Old points (for the last mesh motion)
mutable autoPtr<pointField> oldPointsPtr_;

//- Old cell centres (for the last mesh motion)
mutable autoPtr<pointField> oldCellCentresPtr_;

//- Whether or not to store the old cell centres
mutable bool storeOldCellCentres_;


// Private Member Functions

@@ -421,6 +427,9 @@ public:
//- Return old points for mesh motion
virtual const pointField& oldPoints() const;

//- Return old points for mesh motion
virtual const pointField& oldCellCentres() const;

//- Return IO object for points0
static IOobject points0IO(const polyMesh& mesh);

@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -522,8 +522,10 @@ Foam::polyMesh::polyMesh
globalMeshDataPtr_(nullptr),
moving_(false),
topoChanging_(false),
curMotionTimeIndex_(time().timeIndex()),
oldPointsPtr_(nullptr)
curMotionTimeIndex_(-1),
oldPointsPtr_(nullptr),
oldCellCentresPtr_(nullptr),
storeOldCellCentres_(false)
{
if (debug)
{
@@ -806,8 +808,10 @@ Foam::polyMesh::polyMesh
globalMeshDataPtr_(nullptr),
moving_(false),
topoChanging_(false),
curMotionTimeIndex_(time().timeIndex()),
oldPointsPtr_(nullptr)
curMotionTimeIndex_(-1),
oldPointsPtr_(nullptr),
oldCellCentresPtr_(nullptr),
storeOldCellCentres_(false)
{
if (debug)
{
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -91,6 +91,30 @@ void Foam::polyMesh::updateMesh(const mapPolyMesh& mpm)
}
}

if (oldCellCentresPtr_.valid())
{
// Make a copy of the original cell-centres
pointField oldMotionCellCentres = oldCellCentresPtr_();

pointField& newMotionCellCentres = oldCellCentresPtr_();

// Resize the list to new size
newMotionCellCentres.setSize(cellCentres().size());

// Map the list
newMotionCellCentres.map(oldMotionCellCentres, mpm.cellMap());

// Any points created out-of-nothing get set to the current coordinate
// for lack of anything better.
forAll(mpm.cellMap(), newCelli)
{
if (mpm.cellMap()[newCelli] == -1)
{
newMotionCellCentres[newCelli] = cellCentres()[newCelli];
}
}
}

meshObject::updateMesh<polyMesh>(*this, mpm);
meshObject::updateMesh<pointMesh>(*this, mpm);

@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -78,10 +78,11 @@ Foam::Cloud<ParticleType>::Cloud
{
checkPatches();

// Ask for the tetBasePtIs to trigger all processors to build
// them, otherwise, if some processors have no particles then
// there is a comms mismatch.
// Ask for the tetBasePtIs and oldCellCentres to trigger all processors to
// build them, otherwise, if some processors have no particles then there
// is a comms mismatch.
polyMesh_.tetBasePtIs();
polyMesh_.oldCellCentres();

if (particles.size())
{
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -153,6 +153,9 @@ Foam::Cloud<ParticleType>::Cloud
{
checkPatches();

polyMesh_.tetBasePtIs();
polyMesh_.oldCellCentres();

initCloud(checkClass);
}

@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -68,13 +68,8 @@ inline void Foam::particle::movingTetGeometry
const triFace triIs(currentTetIndices().faceTriIs(mesh_));
const pointField& ptsOld = mesh_.oldPoints();
const pointField& ptsNew = mesh_.points();

// !!! <-- We would be better off using mesh_.cellCentres() here. However,
// we need to put a mesh_.oldCellCentres() method in for this to work. The
// values obtained from the mesh and those obtained from the cell do not
// necessarily match. See mantis #1993.
const vector ccOld = mesh_.cells()[celli_].centre(ptsOld, mesh_.faces());
const vector ccNew = mesh_.cells()[celli_].centre(ptsNew, mesh_.faces());
const vector ccOld = mesh_.oldCellCentres()[celli_];
const vector ccNew = mesh_.cellCentres()[celli_];

// Old and new points and cell centres are not sub-cycled. If we are sub-
// cycling, then we have to account for the timestep change here by

0 comments on commit b5e27b7

Please sign in to comment.
You can’t perform that action at this time.