1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | foam-extend: Open Source CFD
4 \\ / O peration | Version: 3.2
5 \\ / A nd | Web: http://www.foam-extend.org
6 \\/ M anipulation | For copyright notice see file Copyright
7 -------------------------------------------------------------------------------
9 This file is part of foam-extend.
11 foam-extend is free software: you can redistribute it and/or modify it
12 under the terms of the GNU General Public License as published by the
13 Free Software Foundation, either version 3 of the License, or (at your
14 option) any later version.
16 foam-extend is distributed in the hope that it will be useful, but
17 WITHOUT ANY WARRANTY; without even the implied warranty of
18 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19 General Public License for more details.
21 You should have received a copy of the GNU General Public License
22 along with foam-extend. If not, see <http://www.gnu.org/licenses/>.
24 \*---------------------------------------------------------------------------*/
26 #include "CloudTemplate.H"
27 #include "processorPolyPatch.H"
28 #include "globalMeshData.H"
29 #include "PstreamCombineReduceOps.H"
30 #include "mapPolyMesh.H"
34 #include "profiling.H"
36 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
38 template<class ParticleType>
39 Foam::Cloud<ParticleType>::Cloud
41 const polyMesh& pMesh,
42 const IDLList<ParticleType>& particles
46 IDLList<ParticleType>(),
50 IDLList<ParticleType>::operator=(particles);
54 template<class ParticleType>
55 Foam::Cloud<ParticleType>::Cloud
57 const polyMesh& pMesh,
58 const word& cloudName,
59 const IDLList<ParticleType>& particles
62 cloud(pMesh, cloudName),
63 IDLList<ParticleType>(),
67 IDLList<ParticleType>::operator=(particles);
71 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
73 template<class ParticleType>
74 Foam::label Foam::Cloud<ParticleType>::getNewParticleID() const
76 label id = particleCount_++;
80 WarningIn("Cloud<ParticleType>::getNewParticleID() const")
81 << "Particle counter has overflowed. This might cause problems"
82 << " when reconstructing particle tracks." << endl;
88 template<class ParticleType>
89 void Foam::Cloud<ParticleType>::addParticle(ParticleType* pPtr)
95 template<class ParticleType>
96 void Foam::Cloud<ParticleType>::deleteParticle(ParticleType& p)
98 delete(this->remove(&p));
105 class combineNsTransPs
110 void operator()(labelListList& x, const labelListList& y) const
122 } // End namespace Foam
125 template<class ParticleType>
126 template<class TrackingData>
127 void Foam::Cloud<ParticleType>::move(TrackingData& td)
129 addProfile2(moveProfile,"Cloud<ParticleType>::move_"+this->name());
131 const globalMeshData& pData = polyMesh_.globalData();
132 const labelList& processorPatches = pData.processorPatches();
133 const labelList& processorPatchIndices = pData.processorPatchIndices();
134 const labelList& processorPatchNeighbours =
135 pData.processorPatchNeighbours();
137 // Initialise the setpFraction moved for the particles
138 forAllIter(typename Cloud<ParticleType>, *this, pIter)
140 pIter().stepFraction() = 0;
143 // Assume there will be particles to transfer
144 bool transfered = true;
146 // While there are particles to transfer
149 // List of lists of particles to be transfered for all the processor
151 List<IDLList<ParticleType> > transferList(processorPatches.size());
153 // Loop over all particles
154 forAllIter(typename Cloud<ParticleType>, *this, pIter)
156 ParticleType& p = pIter();
159 bool keepParticle = p.move(td);
161 // If the particle is to be kept
162 // (i.e. it hasn't passed through an inlet or outlet)
165 // If we are running in parallel and the particle is on a
167 if (Pstream::parRun() && p.facei_ >= pMesh().nInternalFaces())
169 label patchi = pMesh().boundaryMesh().whichPatch(p.facei_);
170 label n = processorPatchIndices[patchi];
172 // ... and the face is on a processor patch
173 // prepare it for transfer
176 p.prepareForParallelTransfer(patchi, td);
177 transferList[n].append(this->remove(&p));
187 if (Pstream::parRun())
189 // List of the numbers of particles to be transfered across the
191 labelList nsTransPs(transferList.size());
193 forAll(transferList, i)
195 nsTransPs[i] = transferList[i].size();
198 // List of the numbers of particles to be transfered across the
199 // processor patches for all the processors
200 labelListList allNTrans(Pstream::nProcs());
201 allNTrans[Pstream::myProcNo()] = nsTransPs;
202 combineReduce(allNTrans, combineNsTransPs());
208 forAll(allNTrans[i], j)
223 forAll(transferList, i)
225 if (transferList[i].size())
227 OPstream particleStream
230 refCast<const processorPolyPatch>
232 pMesh().boundaryMesh()[processorPatches[i]]
236 particleStream << transferList[i];
240 forAll(processorPatches, i)
242 label patchi = processorPatches[i];
244 const processorPolyPatch& procPatch =
245 refCast<const processorPolyPatch>
246 (pMesh().boundaryMesh()[patchi]);
249 procPatch.neighbProcNo() - Pstream::masterNo();
251 label neighbProcPatchi = processorPatchNeighbours[patchi];
253 label nRecPs = allNTrans[neighbProci][neighbProcPatchi];
257 IPstream particleStream
260 procPatch.neighbProcNo()
262 IDLList<ParticleType> newParticles
265 typename ParticleType::iNew(*this)
270 typename Cloud<ParticleType>,
275 ParticleType& newp = newpIter();
276 newp.correctAfterParallelTransfer(patchi, td);
277 addParticle(newParticles.remove(&newp));
290 template<class ParticleType>
291 void Foam::Cloud<ParticleType>::autoMap(const mapPolyMesh& mapper)
295 Info<< "Cloud<ParticleType>::autoMap(const morphFieldMapper& map) "
296 "for lagrangian cloud " << cloud::name() << endl;
299 const labelList& reverseCellMap = mapper.reverseCellMap();
300 const labelList& reverseFaceMap = mapper.reverseFaceMap();
302 forAllIter(typename Cloud<ParticleType>, *this, pIter)
304 if (reverseCellMap[pIter().celli_] >= 0)
306 pIter().celli_ = reverseCellMap[pIter().celli_];
308 if (pIter().facei_ >= 0 && reverseFaceMap[pIter().facei_] >= 0)
310 pIter().facei_ = reverseFaceMap[pIter().facei_];
319 label trackStartCell = mapper.mergedCell(pIter().celli_);
321 // Tommaso, bug fix 27/02/2008
322 // when the topology of the mesh is changed
323 // it is necessary to update the cell in which the parcel is,
324 // in particular when cells are removed
326 // HJ, merge 1.5.x 20/Oct/2008
327 if (trackStartCell < 0)
329 pIter().celli_ = polyMesh_.findCell(pIter().position());
332 vector p = pIter().position();
333 const_cast<vector&>(pIter().position()) =
334 polyMesh_.cellCentres()[trackStartCell];
335 pIter().stepFraction() = 0;
342 template<class ParticleType>
343 void Foam::Cloud<ParticleType>::writePositions() const
347 this->db().time().path()/this->name() + "_positions.obj"
350 forAllConstIter(typename Cloud<ParticleType>, *this, pIter)
352 const ParticleType& p = pIter();
353 pObj<< "v " << p.position().x() << " " << p.position().y() << " "
354 << p.position().z() << nl;
361 // * * * * * * * * * * * * * * * * IOStream operators * * * * * * * * * * * //
363 #include "CloudTemplateIO.C"
365 // ************************************************************************* //