793 lines
22 KiB
C++
793 lines
22 KiB
C++
/*---------------------------------------------------------------------------*\
|
|
========= |
|
|
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
|
|
\\ / O peration | Website: https://openfoam.org
|
|
\\ / A nd | Copyright (C) 2011-2024 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 "Cloud.H"
|
|
#include "processorPolyPatch.H"
|
|
#include "globalMeshData.H"
|
|
#include "meshToMesh.H"
|
|
#include "PstreamCombineReduceOps.H"
|
|
#include "polyTopoChangeMap.H"
|
|
#include "polyMeshMap.H"
|
|
#include "polyDistributionMap.H"
|
|
#include "Time.H"
|
|
#include "OFstream.H"
|
|
#include "wallPolyPatch.H"
|
|
#include "nonConformalCyclicPolyPatch.H"
|
|
#include "cpuLoad.H"
|
|
|
|
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
|
|
|
|
namespace Foam
|
|
{
|
|
|
|
template<class Type>
|
|
struct IDLListAppendEqOp
|
|
{
|
|
void operator()(IDLList<Type>& x, const IDLList<Type>& y) const
|
|
{
|
|
if (y.size())
|
|
{
|
|
forAllConstIter(typename IDLList<Type>, y, iter)
|
|
{
|
|
x.append(new Type(iter()));
|
|
}
|
|
}
|
|
}
|
|
};
|
|
|
|
}
|
|
|
|
|
|
// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
|
|
|
|
template<class ParticleType>
|
|
Foam::labelList Foam::Cloud<ParticleType>::patchNbrProc
|
|
(
|
|
const polyMesh& pMesh
|
|
)
|
|
{
|
|
const polyBoundaryMesh& pbm = pMesh.boundaryMesh();
|
|
|
|
labelList result(pbm.size(), -1);
|
|
|
|
if (Pstream::parRun())
|
|
{
|
|
forAll(pbm, patchi)
|
|
{
|
|
if (isA<processorPolyPatch>(pbm[patchi]))
|
|
{
|
|
const processorPolyPatch& ppp =
|
|
refCast<const processorPolyPatch>(pbm[patchi]);
|
|
|
|
result[patchi] = ppp.neighbProcNo();
|
|
}
|
|
}
|
|
}
|
|
|
|
return result;
|
|
}
|
|
|
|
|
|
template<class ParticleType>
|
|
Foam::labelList Foam::Cloud<ParticleType>::patchNbrProcPatch
|
|
(
|
|
const polyMesh& pMesh
|
|
)
|
|
{
|
|
const polyBoundaryMesh& pbm = pMesh.boundaryMesh();
|
|
|
|
labelList result(pbm.size(), -1);
|
|
|
|
if (Pstream::parRun())
|
|
{
|
|
PstreamBuffers pBufs(Pstream::commsTypes::nonBlocking);
|
|
|
|
forAll(pbm, patchi)
|
|
{
|
|
if (isA<processorPolyPatch>(pbm[patchi]))
|
|
{
|
|
const processorPolyPatch& ppp =
|
|
refCast<const processorPolyPatch>(pbm[patchi]);
|
|
|
|
UOPstream(ppp.neighbProcNo(), pBufs)()
|
|
<< ppp.index();
|
|
}
|
|
}
|
|
|
|
pBufs.finishedSends();
|
|
|
|
forAll(pbm, patchi)
|
|
{
|
|
if (isA<processorPolyPatch>(pbm[patchi]))
|
|
{
|
|
const processorPolyPatch& ppp =
|
|
refCast<const processorPolyPatch>(pbm[patchi]);
|
|
|
|
UIPstream(ppp.neighbProcNo(), pBufs)()
|
|
>> result[patchi];
|
|
}
|
|
}
|
|
}
|
|
|
|
return result;
|
|
}
|
|
|
|
|
|
template<class ParticleType>
|
|
Foam::labelListList Foam::Cloud<ParticleType>::patchNonConformalCyclicPatches
|
|
(
|
|
const polyMesh& pMesh
|
|
)
|
|
{
|
|
const polyBoundaryMesh& pbm = pMesh.boundaryMesh();
|
|
|
|
labelListList result(pbm.size(), labelList());
|
|
|
|
forAll(pbm, patchi)
|
|
{
|
|
if (isA<nonConformalCyclicPolyPatch>(pbm[patchi]))
|
|
{
|
|
const nonConformalCyclicPolyPatch& nccPp =
|
|
refCast<const nonConformalCyclicPolyPatch>(pbm[patchi]);
|
|
|
|
result[nccPp.origPatchIndex()].append(patchi);
|
|
}
|
|
}
|
|
|
|
return result;
|
|
}
|
|
|
|
|
|
template<class ParticleType>
|
|
void Foam::Cloud<ParticleType>::storeRays() const
|
|
{
|
|
const polyBoundaryMesh& pbm = pMesh_.boundaryMesh();
|
|
|
|
forAll(patchNonConformalCyclicPatches_, patchi)
|
|
{
|
|
forAll(patchNonConformalCyclicPatches_[patchi], i)
|
|
{
|
|
const label nccPatchi =
|
|
patchNonConformalCyclicPatches_[patchi][i];
|
|
|
|
const nonConformalCyclicPolyPatch& nccPp =
|
|
refCast<const nonConformalCyclicPolyPatch>(pbm[nccPatchi]);
|
|
|
|
if (nccPp.owner()) nccPp.rays();
|
|
}
|
|
}
|
|
}
|
|
|
|
|
|
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
|
|
|
|
template<class ParticleType>
|
|
Foam::Cloud<ParticleType>::Cloud
|
|
(
|
|
const polyMesh& pMesh,
|
|
const word& cloudName,
|
|
const IDLList<ParticleType>& particles
|
|
)
|
|
:
|
|
cloud(pMesh, cloudName),
|
|
IDLList<ParticleType>(),
|
|
pMesh_(pMesh),
|
|
patchNbrProc_(patchNbrProc(pMesh)),
|
|
patchNbrProcPatch_(patchNbrProcPatch(pMesh)),
|
|
patchNonConformalCyclicPatches_(patchNonConformalCyclicPatches(pMesh)),
|
|
globalPositionsPtr_(),
|
|
timeIndex_(-1)
|
|
{
|
|
// 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.
|
|
pMesh_.tetBasePtIs();
|
|
pMesh_.oldCellCentres();
|
|
|
|
if (particles.size())
|
|
{
|
|
IDLList<ParticleType>::operator=(particles);
|
|
}
|
|
}
|
|
|
|
|
|
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
|
|
|
|
template<class ParticleType>
|
|
void Foam::Cloud<ParticleType>::addParticle(ParticleType* pPtr)
|
|
{
|
|
this->append(pPtr);
|
|
}
|
|
|
|
|
|
template<class ParticleType>
|
|
void Foam::Cloud<ParticleType>::deleteParticle(ParticleType& p)
|
|
{
|
|
delete(this->remove(&p));
|
|
}
|
|
|
|
|
|
template<class ParticleType>
|
|
void Foam::Cloud<ParticleType>::deleteLostParticles()
|
|
{
|
|
label lostCount = 0;
|
|
|
|
forAllIter(typename Cloud<ParticleType>, *this, pIter)
|
|
{
|
|
if (pIter().cell() == -1)
|
|
{
|
|
deleteParticle(pIter());
|
|
lostCount ++;
|
|
}
|
|
}
|
|
|
|
reduce(lostCount, sumOp<label>());
|
|
if (lostCount != 0)
|
|
{
|
|
WarningInFunction
|
|
<< "Cloud " << this->name()
|
|
<< " deleted " << lostCount << " lost particles" << endl;
|
|
}
|
|
}
|
|
|
|
|
|
template<class ParticleType>
|
|
void Foam::Cloud<ParticleType>::cloudReset(const Cloud<ParticleType>& c)
|
|
{
|
|
// Reset particle count and particles only
|
|
// - not changing the cloud object registry or reference to the polyMesh
|
|
ParticleType::particleCount_ = 0;
|
|
IDLList<ParticleType>::operator=(c);
|
|
}
|
|
|
|
|
|
template<class ParticleType>
|
|
void Foam::Cloud<ParticleType>::changeTimeStep()
|
|
{
|
|
forAllIter(typename Cloud<ParticleType>, *this, pIter)
|
|
{
|
|
pIter().reset(0);
|
|
}
|
|
|
|
timeIndex_ = pMesh_.time().timeIndex();
|
|
}
|
|
|
|
|
|
template<class ParticleType>
|
|
template<class TrackCloudType>
|
|
void Foam::Cloud<ParticleType>::move
|
|
(
|
|
TrackCloudType& cloud,
|
|
typename ParticleType::trackingData& td
|
|
)
|
|
{
|
|
// If the time has changed, modify the particles accordingly
|
|
if (timeIndex_ != pMesh_.time().timeIndex())
|
|
{
|
|
changeTimeStep();
|
|
}
|
|
|
|
// Clear the global positions as these are about to change
|
|
globalPositionsPtr_.clear();
|
|
|
|
// Ensure rays are available for non conformal transfers
|
|
storeRays();
|
|
|
|
// Create transfer buffers
|
|
PstreamBuffers pBufs(Pstream::commsTypes::nonBlocking);
|
|
|
|
// Create lists of particles and patch indices to transfer
|
|
List<IDLList<ParticleType>> sendParticles(Pstream::nProcs());
|
|
List<DynamicList<label>> sendPatchIndices(Pstream::nProcs());
|
|
|
|
optionalCpuLoad& cloudCpuTime
|
|
(
|
|
optionalCpuLoad::New(name() + ":cpuLoad", pMesh_, cloud.cpuLoad())
|
|
);
|
|
|
|
// While there are particles to transfer
|
|
while (true)
|
|
{
|
|
// Clear the transfer lists
|
|
forAll(sendParticles, proci)
|
|
{
|
|
sendParticles[proci].clear();
|
|
sendPatchIndices[proci].clear();
|
|
}
|
|
|
|
if (cloud.cpuLoad())
|
|
{
|
|
cloudCpuTime.resetCpuTime();
|
|
}
|
|
|
|
// Loop over all particles
|
|
forAllIter(typename Cloud<ParticleType>, *this, pIter)
|
|
{
|
|
ParticleType& p = pIter();
|
|
|
|
// Move the particle
|
|
const bool keepParticle = p.move(cloud, td);
|
|
|
|
if (cloud.cpuLoad())
|
|
{
|
|
cloudCpuTime.cpuTimeIncrement(p.cell());
|
|
}
|
|
|
|
// If the particle is to be kept
|
|
if (keepParticle)
|
|
{
|
|
if (td.sendToProc != -1)
|
|
{
|
|
#ifdef FULLDEBUG
|
|
if (!Pstream::parRun() || !p.onBoundaryFace(pMesh_))
|
|
{
|
|
FatalErrorInFunction
|
|
<< "Switch processor flag is true when no parallel "
|
|
<< "transfer is possible. This is a bug."
|
|
<< exit(FatalError);
|
|
}
|
|
#endif
|
|
|
|
p.prepareForParallelTransfer(cloud, td);
|
|
|
|
sendParticles[td.sendToProc].append(this->remove(&p));
|
|
|
|
sendPatchIndices[td.sendToProc].append(td.sendToPatch);
|
|
}
|
|
}
|
|
else
|
|
{
|
|
deleteParticle(p);
|
|
}
|
|
}
|
|
|
|
// If running in serial then everything has been moved, so finish
|
|
if (!Pstream::parRun())
|
|
{
|
|
break;
|
|
}
|
|
|
|
// Clear transfer buffers
|
|
pBufs.clear();
|
|
|
|
// Stream into send buffers
|
|
forAll(sendParticles, proci)
|
|
{
|
|
if (sendParticles[proci].size())
|
|
{
|
|
UOPstream particleStream(proci, pBufs);
|
|
|
|
particleStream
|
|
<< sendPatchIndices[proci]
|
|
<< sendParticles[proci];
|
|
}
|
|
}
|
|
|
|
// Start sending. Sets number of bytes transferred.
|
|
labelList receiveSizes(Pstream::nProcs());
|
|
pBufs.finishedSends(receiveSizes);
|
|
|
|
// Determine if any particles were transferred. If not, then finish.
|
|
bool transferred = false;
|
|
forAll(receiveSizes, proci)
|
|
{
|
|
if (receiveSizes[proci])
|
|
{
|
|
transferred = true;
|
|
break;
|
|
}
|
|
}
|
|
reduce(transferred, orOp<bool>());
|
|
if (!transferred)
|
|
{
|
|
break;
|
|
}
|
|
|
|
// Retrieve from receive buffers and add into the cloud
|
|
forAll(receiveSizes, proci)
|
|
{
|
|
if (receiveSizes[proci])
|
|
{
|
|
UIPstream particleStream(proci, pBufs);
|
|
|
|
const labelList receivePatchIndices(particleStream);
|
|
|
|
IDLList<ParticleType> newParticles(particleStream);
|
|
|
|
label i = 0;
|
|
|
|
forAllIter(typename Cloud<ParticleType>, newParticles, iter)
|
|
{
|
|
const label patchi = receivePatchIndices[i ++];
|
|
|
|
ParticleType& p = iter();
|
|
|
|
td.sendToPatch = patchi;
|
|
|
|
p.correctAfterParallelTransfer(cloud, td);
|
|
|
|
addParticle(newParticles.remove(&p));
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
// Warn about any approximate locates
|
|
Pstream::listCombineGather(td.patchNLocateBoundaryHits, plusEqOp<label>());
|
|
if (Pstream::master())
|
|
{
|
|
forAll(td.patchNLocateBoundaryHits, patchi)
|
|
{
|
|
if (td.patchNLocateBoundaryHits[patchi] != 0)
|
|
{
|
|
WarningInFunction
|
|
<< "Cloud " << name() << " did not accurately locate "
|
|
<< td.patchNLocateBoundaryHits[patchi]
|
|
<< " particles that transferred to patch "
|
|
<< pMesh_.boundaryMesh()[patchi].name() << nl;
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
|
|
template<class ParticleType>
|
|
void Foam::Cloud<ParticleType>::topoChange(const polyTopoChangeMap& map)
|
|
{
|
|
if (map.reverseCellMap().empty()) return;
|
|
|
|
// Ask for the tetBasePtIs to trigger all processors to build
|
|
// them, otherwise, if some processors have no particles then
|
|
// there is a comms mismatch.
|
|
pMesh_.tetBasePtIs();
|
|
pMesh_.oldCellCentres();
|
|
|
|
if (!globalPositionsPtr_.valid())
|
|
{
|
|
FatalErrorInFunction
|
|
<< "Global positions are not available. "
|
|
<< "Cloud::storeGlobalPositions has not been called."
|
|
<< exit(FatalError);
|
|
}
|
|
|
|
const vectorField& positions = globalPositionsPtr_();
|
|
|
|
label lostCount = 0;
|
|
|
|
label particlei = 0;
|
|
forAllIter(typename Cloud<ParticleType>, *this, iter)
|
|
{
|
|
const point& pos = positions[particlei ++];
|
|
|
|
const label celli = map.reverseCellMap()[iter().cell()];
|
|
|
|
if (!iter().locate(pMesh_, pos, celli))
|
|
{
|
|
this->remove(iter);
|
|
lostCount ++;
|
|
}
|
|
}
|
|
|
|
reduce(lostCount, sumOp<label>());
|
|
if (lostCount != 0)
|
|
{
|
|
WarningInFunction
|
|
<< "Topology change of cloud " << this->name()
|
|
<< " lost " << lostCount << " particles" << endl;
|
|
}
|
|
}
|
|
|
|
|
|
template<class ParticleType>
|
|
void Foam::Cloud<ParticleType>::mapMesh(const polyMeshMap& map)
|
|
{
|
|
// Ask for the tetBasePtIs to trigger all processors to build
|
|
// them, otherwise, if some processors have no particles then
|
|
// there is a comms mismatch.
|
|
pMesh_.tetBasePtIs();
|
|
pMesh_.oldCellCentres();
|
|
|
|
// Update cached mesh indexing
|
|
patchNbrProc_ = patchNbrProc(pMesh_);
|
|
patchNbrProcPatch_ = patchNbrProcPatch(pMesh_);
|
|
patchNonConformalCyclicPatches_ = patchNonConformalCyclicPatches(pMesh_);
|
|
|
|
if (!globalPositionsPtr_.valid())
|
|
{
|
|
FatalErrorInFunction
|
|
<< "Global positions are not available. "
|
|
<< "Cloud::storeGlobalPositions has not been called."
|
|
<< exit(FatalError);
|
|
}
|
|
|
|
const vectorField& positions = globalPositionsPtr_();
|
|
|
|
label lostCount = 0;
|
|
|
|
// Loop the particles. Map those that remain on this processor, and
|
|
// transfer others into send arrays.
|
|
List<DynamicList<label>> sendCellIndices(Pstream::nProcs());
|
|
List<DynamicList<point>> sendPositions(Pstream::nProcs());
|
|
List<IDLList<ParticleType>> sendParticles(Pstream::nProcs());
|
|
{
|
|
label particlei = 0;
|
|
forAllIter(typename Cloud<ParticleType>, *this, iter)
|
|
{
|
|
const point& pos = positions[particlei ++];
|
|
|
|
const remote tgtProcCell =
|
|
map.mapper().srcToTgtPoint(iter().cell(), pos);
|
|
const label proci = tgtProcCell.proci;
|
|
const label celli = tgtProcCell.elementi;
|
|
|
|
if (tgtProcCell == remote())
|
|
{
|
|
this->remove(iter);
|
|
lostCount ++;
|
|
}
|
|
else if (proci == Pstream::myProcNo())
|
|
{
|
|
if (!iter().locate(pMesh_, pos, celli))
|
|
{
|
|
this->remove(iter);
|
|
lostCount ++;
|
|
}
|
|
}
|
|
else
|
|
{
|
|
sendCellIndices[proci].append(celli);
|
|
sendPositions[proci].append(pos);
|
|
sendParticles[proci].append(this->remove(iter));
|
|
}
|
|
}
|
|
}
|
|
|
|
// If parallel then send and receive particles that move processes and map
|
|
// those sent to this process
|
|
if (Pstream::parRun())
|
|
{
|
|
// Create transfer buffers
|
|
PstreamBuffers pBufs(Pstream::commsTypes::nonBlocking);
|
|
|
|
// Stream into send buffers
|
|
forAll(sendParticles, proci)
|
|
{
|
|
if (sendParticles[proci].size())
|
|
{
|
|
UOPstream particleStream(proci, pBufs);
|
|
|
|
particleStream
|
|
<< sendCellIndices[proci]
|
|
<< sendPositions[proci]
|
|
<< sendParticles[proci];
|
|
}
|
|
}
|
|
|
|
// Finish sending
|
|
labelList receiveSizes(Pstream::nProcs());
|
|
pBufs.finishedSends(receiveSizes);
|
|
|
|
// Retrieve from receive buffers and map into the new mesh
|
|
forAll(sendParticles, proci)
|
|
{
|
|
if (receiveSizes[proci])
|
|
{
|
|
UIPstream particleStream(proci, pBufs);
|
|
|
|
const labelList receiveCellIndices(particleStream);
|
|
const List<point> receivePositions(particleStream);
|
|
IDLList<ParticleType> receiveParticles(particleStream);
|
|
|
|
label particlei = 0;
|
|
forAllIter(typename Cloud<ParticleType>, receiveParticles, iter)
|
|
{
|
|
const label celli = receiveCellIndices[particlei];
|
|
const vector& pos = receivePositions[particlei ++];
|
|
|
|
if (iter().locate(pMesh_, pos, celli))
|
|
{
|
|
this->append(receiveParticles.remove(iter));
|
|
}
|
|
else
|
|
{
|
|
receiveParticles.remove(iter);
|
|
lostCount ++;
|
|
}
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
reduce(lostCount, sumOp<label>());
|
|
if (lostCount != 0)
|
|
{
|
|
WarningInFunction
|
|
<< "Mesh-to-mesh mapping of cloud " << this->name()
|
|
<< " lost " << lostCount << " particles" << endl;
|
|
}
|
|
}
|
|
|
|
|
|
template<class ParticleType>
|
|
void Foam::Cloud<ParticleType>::distribute(const polyDistributionMap& map)
|
|
{
|
|
// Ask for the tetBasePtIs to trigger all processors to build
|
|
// them, otherwise, if some processors have no particles then
|
|
// there is a comms mismatch.
|
|
pMesh_.tetBasePtIs();
|
|
pMesh_.oldCellCentres();
|
|
|
|
// Update cached mesh indexing
|
|
patchNbrProc_ = patchNbrProc(pMesh_);
|
|
patchNbrProcPatch_ = patchNbrProcPatch(pMesh_);
|
|
patchNonConformalCyclicPatches_ = patchNonConformalCyclicPatches(pMesh_);
|
|
|
|
if (!globalPositionsPtr_.valid())
|
|
{
|
|
FatalErrorInFunction
|
|
<< "Global positions are not available. "
|
|
<< "Cloud::storeGlobalPositions has not been called."
|
|
<< exit(FatalError);
|
|
}
|
|
|
|
const vectorField& positions = globalPositionsPtr_();
|
|
|
|
// Distribute the global positions
|
|
List<List<point>> cellParticlePositions(map.nOldCells());
|
|
{
|
|
labelList cellParticleis(map.nOldCells(), 0);
|
|
forAllIter(typename Cloud<ParticleType>, *this, iter)
|
|
{
|
|
cellParticleis[iter().cell()] ++;
|
|
}
|
|
forAll(cellParticlePositions, celli)
|
|
{
|
|
cellParticlePositions[celli].resize(cellParticleis[celli]);
|
|
}
|
|
|
|
label particlei = 0;
|
|
cellParticleis = 0;
|
|
forAllIter(typename Cloud<ParticleType>, *this, iter)
|
|
{
|
|
const label celli = iter().cell();
|
|
label& cellParticlei = cellParticleis[celli];
|
|
|
|
cellParticlePositions[celli][cellParticlei ++] =
|
|
positions[particlei ++];
|
|
}
|
|
}
|
|
distributionMapBase::distribute
|
|
(
|
|
Pstream::commsTypes::nonBlocking,
|
|
List<labelPair>(),
|
|
pMesh_.nCells(),
|
|
map.cellMap().subMap(),
|
|
false,
|
|
map.cellMap().constructMap(),
|
|
false,
|
|
cellParticlePositions,
|
|
ListAppendEqOp<point>(),
|
|
flipOp(),
|
|
List<point>()
|
|
);
|
|
|
|
// Distribute the particles
|
|
List<IDLList<ParticleType>> cellParticles(map.nOldCells());
|
|
forAllIter(typename Cloud<ParticleType>, *this, iter)
|
|
{
|
|
cellParticles[iter().cell()].append(this->remove(iter));
|
|
}
|
|
distributionMapBase::distribute
|
|
(
|
|
Pstream::commsTypes::nonBlocking,
|
|
List<labelPair>(),
|
|
pMesh_.nCells(),
|
|
map.cellMap().subMap(),
|
|
false,
|
|
map.cellMap().constructMap(),
|
|
false,
|
|
cellParticles,
|
|
IDLListAppendEqOp<ParticleType>(),
|
|
flipOp(),
|
|
IDLList<ParticleType>()
|
|
);
|
|
|
|
label lostCount = 0;
|
|
|
|
// Locate the particles within the new mesh
|
|
forAll(cellParticles, celli)
|
|
{
|
|
label cellParticlei = 0;
|
|
forAllIter(typename IDLList<ParticleType>, cellParticles[celli], iter)
|
|
{
|
|
const point& pos = cellParticlePositions[celli][cellParticlei++];
|
|
|
|
if (iter().locate(pMesh_, pos, celli))
|
|
{
|
|
this->append(cellParticles[celli].remove(iter));
|
|
}
|
|
else
|
|
{
|
|
cellParticles[celli].remove(iter);
|
|
lostCount ++;
|
|
}
|
|
}
|
|
}
|
|
|
|
reduce(lostCount, sumOp<label>());
|
|
if (lostCount != 0)
|
|
{
|
|
WarningInFunction
|
|
<< "Distribution of cloud " << this->name()
|
|
<< " lost " << lostCount << " particles" << endl;
|
|
}
|
|
}
|
|
|
|
|
|
template<class ParticleType>
|
|
void Foam::Cloud<ParticleType>::writePositions() const
|
|
{
|
|
OFstream pObj
|
|
(
|
|
this->db().time().path()/this->name() + "_positions.obj"
|
|
);
|
|
|
|
forAllConstIter(typename Cloud<ParticleType>, *this, pIter)
|
|
{
|
|
const point pos = pIter().position(pMesh_);
|
|
|
|
pObj<< "v " << pos.x() << " " << pos.y() << " " << pos.z() << nl;
|
|
}
|
|
|
|
pObj.flush();
|
|
}
|
|
|
|
|
|
template<class ParticleType>
|
|
void Foam::Cloud<ParticleType>::storeGlobalPositions() const
|
|
{
|
|
// Store the global positions for later use by mapping functions. It would
|
|
// be preferable not to need this. If the objects passed to had a copy of
|
|
// the old mesh then the global positions could be recovered within the
|
|
// mapping functions, and this pre-processing would not be necessary.
|
|
|
|
globalPositionsPtr_.reset(new vectorField(this->size()));
|
|
|
|
vectorField& positions = globalPositionsPtr_();
|
|
|
|
label particlei = 0;
|
|
forAllConstIter(typename Cloud<ParticleType>, *this, iter)
|
|
{
|
|
positions[particlei++] = iter().position(pMesh_);
|
|
}
|
|
}
|
|
|
|
|
|
// * * * * * * * * * * * * * * * * IOStream operators * * * * * * * * * * * //
|
|
|
|
#include "CloudIO.C"
|
|
|
|
// ************************************************************************* //
|