46 bool Foam::uniformSet::nextSample
48 const point& currentPt,
50 const scalar smallDist,
55 bool pointFound =
false;
57 const vector normOffset = offset/
mag(offset);
62 for (; sampleI < nPoints_; sampleI++)
64 scalar
s = (samplePt - currentPt) & normOffset;
80 bool Foam::uniformSet::trackToBoundary
82 passiveParticleCloud& particleCloud,
83 passiveParticle& singleParticle,
86 DynamicList<point>& samplingPts,
87 DynamicList<label>& samplingCells,
88 DynamicList<label>& samplingFaces,
89 DynamicList<scalar>& samplingCurveDist
93 const vector offset = (end_ - start_)/(nPoints_ - 1);
94 const vector smallVec = tol*offset;
95 const scalar smallDist =
mag(smallVec);
98 const point& trackPt = singleParticle.position();
100 particle::TrackingData<passiveParticleCloud> trackData(particleCloud);
105 if (!nextSample(trackPt, offset, smallDist, samplePt, sampleI))
110 Pout<<
"trackToBoundary : Reached end : samplePt now:" 111 << samplePt <<
" sampleI now:" << sampleI <<
endl;
116 if (
mag(samplePt - trackPt) < smallDist)
121 Pout<<
"trackToBoundary : samplePt corresponds to trackPt : " 122 <<
" trackPt:" << trackPt <<
" samplePt:" << samplePt
126 samplingPts.append(trackPt);
127 samplingCells.append(singleParticle.cell());
128 samplingFaces.append(-1);
129 samplingCurveDist.append(
mag(trackPt - start_));
132 if (!nextSample(trackPt, offset, smallDist, samplePt, sampleI))
137 Pout<<
"trackToBoundary : Reached end : " 138 <<
" samplePt now:" << samplePt
139 <<
" sampleI now:" << sampleI
150 Pout<<
"Searching along trajectory from " 151 <<
" trackPt:" << trackPt
152 <<
" trackCelli:" << singleParticle.cell()
153 <<
" to:" << samplePt <<
endl;
156 point oldPos = trackPt;
160 singleParticle.stepFraction() = 0;
161 singleParticle.track(samplePt, trackData);
165 Pout<<
"Result of tracking " 166 <<
" trackPt:" << trackPt
167 <<
" trackCelli:" << singleParticle.cell()
168 <<
" trackFacei:" << singleParticle.face()
169 <<
" onBoundary:" << singleParticle.onBoundary()
170 <<
" samplePt:" << samplePt
171 <<
" smallDist:" << smallDist
177 !singleParticle.onBoundary()
178 && (
mag(trackPt - oldPos) < smallDist)
181 if (singleParticle.onBoundary())
184 if (
mag(trackPt - samplePt) < smallDist)
189 samplingPts.append(trackPt);
190 samplingCells.append(singleParticle.cell());
191 samplingFaces.append(facei);
192 samplingCurveDist.append(
mag(trackPt - start_));
200 samplingPts.append(trackPt);
201 samplingCells.append(singleParticle.cell());
202 samplingFaces.append(-1);
203 samplingCurveDist.append(
mag(trackPt - start_));
210 void Foam::uniformSet::calcSamples
212 DynamicList<point>& samplingPts,
213 DynamicList<label>& samplingCells,
214 DynamicList<label>& samplingFaces,
215 DynamicList<label>& samplingSegments,
216 DynamicList<scalar>& samplingCurveDist
220 if ((nPoints_ < 2) || (
mag(end_ - start_) < SMALL))
223 <<
"Incorrect sample specification. Either too few points or" 224 <<
" start equals end point." <<
endl 225 <<
"nPoints:" << nPoints_
226 <<
" start:" << start_
231 const vector offset = (end_ - start_)/(nPoints_ - 1);
232 const vector normOffset = offset/
mag(offset);
233 const vector smallVec = tol*offset;
234 const scalar smallDist =
mag(smallVec);
237 const bool oldMoving =
const_cast<polyMesh&
>(
mesh()).moving(
false);
238 passiveParticleCloud particleCloud(
mesh());
247 point bPoint(GREAT, GREAT, GREAT);
252 bPoint = bHits[0].hitPoint();
253 bFacei = bHits[0].index();
259 label trackCelli = -1;
260 label trackFacei = -1;
275 if (trackCelli == -1)
286 samplingPts.append(start_);
287 samplingCells.append(trackCelli);
288 samplingFaces.append(trackFacei);
289 samplingCurveDist.append(0.0);
300 label startSegmentI = 0;
303 point samplePt = start_;
311 passiveParticle singleParticle(
mesh(), trackPt, trackCelli);
313 bool reachedBoundary = trackToBoundary
326 for (
label i = samplingPts.size() - 1; i >= startSegmentI; --i)
328 samplingSegments.append(segmentI);
332 if (!reachedBoundary)
336 Pout<<
"calcSamples : Reached end of samples: " 337 <<
" samplePt now:" << samplePt
338 <<
" sampleI now:" << sampleI
345 bool foundValidB =
false;
347 while (bHitI < bHits.size())
350 (bHits[bHitI].hitPoint() - singleParticle.position())
355 Pout<<
"Finding next boundary : " 356 <<
"bPoint:" << bHits[bHitI].hitPoint()
357 <<
" tracking:" << singleParticle.position()
362 if (dist > smallDist)
382 trackPt = pushIn(bPoint, trackFacei);
383 trackCelli = getBoundaryCell(trackFacei);
387 startSegmentI = samplingPts.size();
390 const_cast<polyMesh&
>(
mesh()).moving(oldMoving);
394 void Foam::uniformSet::genSamples()
397 DynamicList<point> samplingPts;
398 DynamicList<label> samplingCells;
399 DynamicList<label> samplingFaces;
400 DynamicList<label> samplingSegments;
401 DynamicList<scalar> samplingCurveDist;
412 samplingPts.shrink();
413 samplingCells.shrink();
414 samplingFaces.shrink();
415 samplingSegments.shrink();
416 samplingCurveDist.shrink();
465 start_(dict.
lookup(
"start")),
Various (local, not parallel) searches on polyMesh; uses (demand driven) octree to search...
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
errorManipArg< error, int > exit(error &err, const int errNo=1)
List< pointIndexHit > intersections(const point &pStart, const point &pEnd) const
Find all intersections of boundary within segment pStart .. pEnd.
const double e
Elementary charge.
A list of keyword definitions, which are a keyword followed by any number of values (e...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Vector< scalar > vector
A scalar version of the templated Vector.
Macros for easy insertion into run-time selection tables.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Holds list of sampling points which is filled at construction time. Various implementations of this b...
A class for handling words, derived from string.
label readLabel(Istream &is)
prefixOSstream Pout(cout,"Pout")
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
defineTypeNameAndDebug(combustionModel, 0)
vector point
Point is a vector.
dimensioned< scalar > mag(const dimensioned< Type > &)
Mesh consisting of general polyhedral cells.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.