47 bool Foam::uniformSet::nextSample
49 const point& currentPt,
51 const scalar smallDist,
56 bool pointFound =
false;
58 const vector normOffset(offset/
mag(offset));
63 for (; sampleI < nPoints_; ++sampleI)
65 const scalar
s = (samplePt - currentPt) & normOffset;
81 bool Foam::uniformSet::trackToBoundary
83 passiveParticleCloud& particleCloud,
84 passiveParticle& singleParticle,
87 DynamicList<point>& samplingPts,
88 DynamicList<label>& samplingCells,
89 DynamicList<label>& samplingFaces,
90 DynamicList<scalar>& samplingCurveDist
94 const vector offset((end_ - start_)/(nPoints_ - 1));
95 const scalar smallDist =
mag(tol_*offset);
97 point trackPt = singleParticle.position();
102 if (!nextSample(trackPt, offset, smallDist, samplePt, sampleI))
107 Pout<<
"trackToBoundary : Reached end : samplePt now:" 108 << samplePt <<
" sampleI now:" << sampleI <<
endl;
113 if (
mag(samplePt - trackPt) < smallDist)
118 Pout<<
"trackToBoundary : samplePt corresponds to trackPt : " 119 <<
" trackPt:" << trackPt <<
" samplePt:" << samplePt
123 samplingPts.append(trackPt);
124 samplingCells.append(singleParticle.cell());
125 samplingFaces.append(-1);
126 samplingCurveDist.append(
mag(trackPt - start_));
129 if (!nextSample(trackPt, offset, smallDist, samplePt, sampleI))
134 Pout<<
"trackToBoundary : Reached end : " 135 <<
" samplePt now:" << samplePt
136 <<
" sampleI now:" << sampleI
147 Pout<<
"Searching along trajectory from " 148 <<
" trackPt:" << trackPt
149 <<
" trackCelli:" << singleParticle.cell()
150 <<
" to:" << samplePt <<
endl;
153 singleParticle.track(samplePt - trackPt, 0);
154 trackPt = singleParticle.position();
156 if (singleParticle.onBoundaryFace())
159 if (
mag(trackPt - samplePt) < smallDist)
164 samplingPts.append(trackPt);
165 samplingCells.append(singleParticle.cell());
166 samplingFaces.append(singleParticle.face());
167 samplingCurveDist.append(
mag(trackPt - start_));
175 samplingPts.append(trackPt);
176 samplingCells.append(singleParticle.cell());
177 samplingFaces.append(-1);
178 samplingCurveDist.append(
mag(trackPt - start_));
185 void Foam::uniformSet::calcSamples
187 DynamicList<point>& samplingPts,
188 DynamicList<label>& samplingCells,
189 DynamicList<label>& samplingFaces,
190 DynamicList<label>& samplingSegments,
191 DynamicList<scalar>& samplingCurveDist
195 if ((nPoints_ < 2) || (
mag(end_ - start_) < SMALL))
198 <<
"Incorrect sample specification. Either too few points or" 199 <<
" start equals end point." <<
endl 200 <<
"nPoints:" << nPoints_
201 <<
" start:" << start_
206 const vector offset((end_ - start_)/(nPoints_ - 1));
207 const vector normOffset(offset/
mag(offset));
208 const vector smallVec(tol_*offset);
209 const scalar smallDist =
mag(smallVec);
212 const bool oldMoving =
const_cast<polyMesh&
>(
mesh()).moving(
false);
213 passiveParticleCloud particleCloud(
mesh());
222 point bPoint(GREAT, GREAT, GREAT);
227 bPoint = bHits[0].hitPoint();
228 bFacei = bHits[0].index();
234 label trackCelli = -1;
235 label trackFacei = -1;
237 const bool isSample =
250 if (trackCelli == -1)
256 const_cast<polyMesh&
>(
mesh()).moving(oldMoving);
263 samplingPts.append(start_);
264 samplingCells.append(trackCelli);
265 samplingFaces.append(trackFacei);
266 samplingCurveDist.append(0.0);
277 label startSegmentI = 0;
280 point samplePt = start_;
288 passiveParticle singleParticle(
mesh(), trackPt, trackCelli);
290 bool reachedBoundary = trackToBoundary
303 for (label i = samplingPts.size() - 1; i >= startSegmentI; --i)
305 samplingSegments.append(segmentI);
309 if (!reachedBoundary)
313 Pout<<
"calcSamples : Reached end of samples: " 314 <<
" samplePt now:" << samplePt
315 <<
" sampleI now:" << sampleI
322 bool foundValidB =
false;
324 while (bHitI < bHits.size())
327 (bHits[bHitI].hitPoint() - singleParticle.position())
332 Pout<<
"Finding next boundary : " 333 <<
"bPoint:" << bHits[bHitI].point()
334 <<
" tracking:" << singleParticle.position()
339 if (dist > smallDist)
342 bPoint = bHits[bHitI].hitPoint();
343 bFacei = bHits[bHitI].index();
361 trackPt = pushIn(bPoint, trackFacei);
362 trackCelli = getBoundaryCell(trackFacei);
366 startSegmentI = samplingPts.size();
369 const_cast<polyMesh&
>(
mesh()).moving(oldMoving);
373 void Foam::uniformSet::genSamples()
376 DynamicList<point> samplingPts;
377 DynamicList<label> samplingCells;
378 DynamicList<label> samplingFaces;
379 DynamicList<label> samplingSegments;
380 DynamicList<scalar> samplingCurveDist;
391 samplingPts.shrink();
392 samplingCells.shrink();
393 samplingFaces.shrink();
394 samplingSegments.shrink();
395 samplingCurveDist.shrink();
400 std::move(samplingPts),
401 std::move(samplingCells),
402 std::move(samplingFaces),
403 std::move(samplingSegments),
404 std::move(samplingCurveDist)
419 const polyMesh&
mesh,
420 const meshSearch& searchEngine,
427 sampledSet(
name,
mesh, searchEngine, axis),
449 nPoints_(
dict.
get<label>(
"nPoints"))
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
List of values generated by applying the access operation to each list item.
Various (local, not parallel) searches on polyMesh; uses (demand driven) octree to search...
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Macros for easy insertion into run-time selection tables.
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
const dimensionedScalar e
Elementary charge.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Holds list of sampling points which is filled at construction time. Various implementations of this b...
A class for handling words, derived from Foam::string.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
int debug
Static debugging option.
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
defineTypeNameAndDebug(combustionModel, 0)
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
vector point
Point is a vector.
Mesh consisting of general polyhedral cells.
List< pointIndexHit > intersections(const point &pStart, const point &pEnd) const
Find all intersections of boundary within segment pStart .. pEnd.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.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))
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.