26 #include "MapLagrangianFields.H" 35 static const scalar perturbFactor = 1
e-6;
40 static label findCell(
const Cloud<passiveParticle>& cloud,
const point& pt)
46 const polyMesh& mesh = cloud.pMesh();
48 mesh.findCellFacePt(pt, celli, tetFacei, tetPtI);
59 meshSearch meshSearcher
65 label facei = meshSearcher.findNearestBoundaryFace(pt);
69 const point& cc = mesh.cellCentres()[mesh.faceOwner()[facei]];
71 const point perturbPt = (1-perturbFactor)*pt+perturbFactor*cc;
73 mesh.findCellFacePt(perturbPt, celli, tetFacei, tetPtI);
88 const polyMesh& meshSource = interp.srcRegion();
89 const polyMesh& meshTarget = interp.tgtRegion();
90 const labelListList& sourceToTarget = interp.srcToTgtCellAddr();
92 const pointField& targetCc = meshTarget.cellCentres();
109 meshSource.time().timeName(),
113 IOobject* positionsPtr = objects.lookup(word(
"positions"));
117 Info<<
nl <<
" processing cloud " << cloudDirs[cloudI] <<
endl;
120 passiveParticleCloud sourceParcels
126 Info<<
" read " << sourceParcels.size()
127 <<
" parcels from source mesh." <<
endl;
130 passiveParticleCloud targetParcels
134 IDLList<passiveParticle>()
137 particle::TrackingData<passiveParticleCloud> td(targetParcels);
139 label sourceParticleI = 0;
142 DynamicList<label> addParticles(sourceParcels.size());
155 bool foundCell =
false;
158 if (iter().cell() >= 0)
161 sourceToTarget[iter().cell()];
170 autoPtr<passiveParticle> newPtr
175 targetCc[targetCells[i]],
179 passiveParticle& newP = newPtr();
181 label facei = newP.track(iter().position(), td);
183 if (facei < 0 && newP.cell() >= 0)
187 addParticles.append(sourceParticleI);
188 targetParcels.addParticle(newPtr.ptr());
197 unmappedSource.insert(sourceParticleI);
203 Info<<
" after meshToMesh addressing found " 204 << targetParcels.size()
205 <<
" parcels in target mesh." <<
endl;
211 if (unmappedSource.size())
215 forAllIter(Cloud<passiveParticle>, sourceParcels, iter)
217 if (unmappedSource.found(sourceParticleI))
220 findCell(targetParcels, iter().position());
224 unmappedSource.erase(sourceParticleI);
225 addParticles.append(sourceParticleI);
226 iter().cell() = targetCell;
227 targetParcels.addParticle
229 sourceParcels.remove(&iter())
236 addParticles.shrink();
238 Info<<
" after additional mesh searching found " 239 << targetParcels.size() <<
" parcels in target mesh." <<
endl;
241 if (addParticles.size())
243 IOPosition<passiveParticleCloud>(targetParcels).
write();
251 MapLagrangianFields<label>
258 MapLagrangianFields<scalar>
265 MapLagrangianFields<vector>
272 MapLagrangianFields<sphericalTensor>
279 MapLagrangianFields<symmTensor>
286 MapLagrangianFields<tensor>
List< labelList > labelListList
A List of labelList.
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
const double e
Elementary charge.
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
Ostream & endl(Ostream &os)
Add newline and flush stream.
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
vectorField pointField
pointField is a vectorField.
List< label > labelList
A List of labels.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
static const word prefix
The prefix to local: lagrangian.
void mapLagrangian(const meshToMesh0 &meshToMesh0Interp)
Maps lagrangian positions and fields.
vector point
Point is a vector.
List< fileName > fileNameList
A List of fileNames.
fileNameList readDir(const fileName &, const fileName::Type=fileName::FILE, const bool filtergz=true)
Read a directory and return the entries as a string list.