29#include "MapLagrangianFields.H"
38static const scalar perturbFactor = 1
e-6;
43static label findCell(
const Cloud<passiveParticle>& cloud,
const point& pt)
49 const polyMesh&
mesh = cloud.pMesh();
62 meshSearch meshSearcher
68 label facei = meshSearcher.findNearestBoundaryFace(pt);
74 const point perturbPt = (1-perturbFactor)*pt+perturbFactor*cc;
91 const polyMesh& meshSource = interp.srcRegion();
92 const polyMesh& meshTarget = interp.tgtRegion();
93 const labelListList& sourceToTarget = interp.srcToTgtCellAddr();
104 for (
const fileName& cloudDir : cloudDirs)
110 meshSource.time().timeName(),
118 (objects.found(
"coordinates") || objects.found(
"positions")),
124 Info<<
nl <<
" processing cloud " << cloudDir <<
endl;
127 passiveParticleCloud sourceParcels
133 Info<<
" read " << sourceParcels.size()
134 <<
" parcels from source mesh." <<
endl;
137 passiveParticleCloud targetParcels
141 IDLList<passiveParticle>()
144 passiveParticle::trackingData td(targetParcels);
146 label sourceParticleI = 0;
149 DynamicList<label> addParticles(sourceParcels.size());
160 for (
const passiveParticle&
p : sourceParcels)
162 bool foundCell =
false;
168 sourceToTarget[
p.cell()];
174 for (
const label targetCell : targetCells)
177 autoPtr<passiveParticle> newPtr
184 meshTarget.cells()[targetCell][0],
188 passiveParticle& newP = newPtr();
190 newP.track(
p.position() - newP.position(), 0);
196 addParticles.append(sourceParticleI);
197 targetParcels.addParticle(newPtr.ptr());
206 unmappedSource.insert(sourceParticleI);
212 Info<<
" after meshToMesh addressing found "
213 << targetParcels.size()
214 <<
" parcels in target mesh." <<
endl;
220 if (unmappedSource.size())
224 for (passiveParticle&
p : sourceParcels)
226 if (unmappedSource.found(sourceParticleI))
228 const label targetCell =
229 findCell(targetParcels,
p.position());
233 unmappedSource.erase(sourceParticleI);
234 addParticles.append(sourceParticleI);
235 targetParcels.addParticle
244 sourceParcels.remove(&
p);
250 addParticles.shrink();
252 Info<<
" after additional mesh searching found "
253 << targetParcels.size() <<
" parcels in target mesh." <<
endl;
255 if (addParticles.size())
257 IOPosition<passiveParticleCloud>(targetParcels).
write();
265 MapLagrangianFields<label>
272 MapLagrangianFields<scalar>
279 MapLagrangianFields<vector>
286 MapLagrangianFields<sphericalTensor>
293 MapLagrangianFields<symmTensor>
300 MapLagrangianFields<tensor>
virtual bool write(const token &tok)=0
Write token to stream or otherwise handle it.
static const word prefix
The prefix to local: lagrangian.
void findCellFacePt(const point &p, label &celli, label &tetFacei, label &tetPti) const
Find the cell, tetFacei and tetPti for point p.
virtual const labelList & faceOwner() const
Return face owner.
const vectorField & cellCentres() const
List< fileName > fileNameList
A List of fileNames.
Barycentric< scalar > barycentric
A scalar version of the templated Barycentric.
List< label > labelList
A List of labels.
messageStream Info
Information stream (stdout output on master, null elsewhere)
vector point
Point is a vector.
Ostream & endl(Ostream &os)
Add newline and flush stream.
List< labelList > labelListList
A List of labelList.
fileNameList readDir(const fileName &directory, const fileName::Type type=fileName::FILE, const bool filtergz=true, const bool followLink=true)
Read a directory and return the entries as a fileName List.
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce (copy) and return value.
void mapLagrangian(const meshToMesh0 &meshToMesh0Interp)
Maps lagrangian positions and fields.
constexpr char nl
The newline '\n' character (0x0a)