58 <<
"Sizes not equal : "
59 <<
" points:" <<
size()
77 if (facei <
mesh().nInternalFaces())
94 const label
cells[4] =
97 getNeighbourCell(faces_[samplei]),
99 getNeighbourCell(faces_[samplei+1])
120 <<
"Could not find mid-point " <<
p
121 <<
" cell " << cellm <<
endl;
129 for (label i=0; i<4; ++i)
140 <<
"Could not find cell for mid-point" <<
nl
141 <<
" samplei: " << samplei
142 <<
" pts[samplei]: " << operator[](samplei)
143 <<
" face[samplei]: " << faces_[samplei]
144 <<
" pts[samplei+1]: " << operator[](samplei+1)
145 <<
" face[samplei+1]: " << faces_[samplei+1]
146 <<
" cellio: " <<
cells[0]
147 <<
" cellin: " <<
cells[1]
148 <<
" celljo: " <<
cells[2]
149 <<
" celljn: " <<
cells[3]
166 scalar magVec =
mag(vec);
186 const scalar smallDist
208 if (dist < smallDist)
210 return myFaces[myFacei];
226 point newPosition = facePt;
239 const scalar trackingCorrectionTol = 1
e-5;
241 if (tetFacei == -1 || tetPtI == -1)
243 newPosition = facePt;
245 label trap(1.0/trackingCorrectionTol + 1);
251 newPosition += trackingCorrectionTol*(cC - facePt);
263 }
while (tetFacei < 0 && iterNo <= trap);
269 <<
"After pushing " << facePt <<
" to " << newPosition
270 <<
" it is still outside face " << facei
272 <<
" of cell " << celli
273 <<
" at " << cC <<
endl
274 <<
"Please change your starting point"
284 const point& samplePt,
287 const scalar smallDist,
294 bool isGoodSample =
false;
299 trackCelli =
mesh().
findCell(samplePt, searchEngine_.decompMode());
304 || !
mesh().pointInCell
308 searchEngine_.decompMode()
318 isGoodSample =
false;
330 else if (
mag(samplePt - bPoint) < smallDist)
333 trackPt = pushIn(bPoint, bFacei);
335 trackCelli = getBoundaryCell(trackFacei);
341 scalar
sign = calcSign(bFacei, samplePt);
348 trackCelli =
mesh().
findCell(trackPt, searchEngine_.decompMode());
355 trackPt = pushIn(bPoint, bFacei);
357 trackCelli = getBoundaryCell(trackFacei);
359 isGoodSample =
false;
364 <<
" samplePt:" << samplePt
365 <<
" bPoint:" << bPoint
366 <<
" bFacei:" << bFacei <<
nl
367 <<
" Calculated first tracking point :"
368 <<
" trackPt:" << trackPt
369 <<
" trackCelli:" << trackCelli
370 <<
" trackFacei:" << trackFacei
371 <<
" isGoodSample:" << isGoodSample
387 setPoints(samplingPts);
388 setDistance(samplingDistance,
false);
390 segments_ = samplingSegments;
391 cells_ = samplingCells;
392 faces_ = samplingFaces;
407 setPoints(std::move(samplingPts));
408 setDistance(std::move(samplingDistance),
false);
410 segments_ = std::move(samplingSegments);
411 cells_ = std::move(samplingCells);
412 faces_ = std::move(samplingFaces);
430 searchEngine_(searchEngine),
447 searchEngine_(searchEngine),
464 searchEngine_(searchEngine),
502 auto* ctorPtr = wordConstructorTable(sampleType);
511 *wordConstructorTablePtr_
532 os <<
nl <<
"\t(celli)\t(facei)" <<
nl;
536 os <<
'\t' << cells_[samplei]
537 <<
'\t' << faces_[samplei]
Minimal example by using system/controlDict.functions:
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Describes the interaction of a face and a point. It carries the info of a successful hit and (if succ...
const point_type & missPoint() const
Return the miss point. Fatal if hit.
bool hit() const noexcept
Is there a hit.
const point_type & hitPoint() const
Return the hit point. Fatal if not hit.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
label size() const noexcept
The number of elements in the UList.
void size(const label n)
Older name for setAddressableSize.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
A cell is defined as a list of faces with extra functionality.
Holds list of sampling positions.
scalarList distance_
Cumulative distance for "distance" write specifier.
coordFormat
Enumeration defining the output format for coordinates.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
const dictionary & optionalSubDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary, otherwise return this dictionary.
A face is a list of labels corresponding to mesh vertices.
virtual bool write()
Write the output fields.
Various (local, not parallel) searches on polyMesh; uses (demand driven) octree to search.
Mesh consisting of general polyhedral cells.
virtual const faceList & faces() const
Return raw faces.
label findCell(const point &p, const cellDecomposition=CELL_TETS) const
Find cell enclosing this location and return index.
virtual const labelList & faceOwner() const
Return face owner.
void findTetFacePt(const label celli, const point &p, label &tetFacei, label &tetPti) const
Find the tetFacei and tetPti for point p in celli.
bool pointInCell(const point &p, label celli, const cellDecomposition=CELL_TETS) const
Test if point p is in the celli.
virtual const labelList & faceNeighbour() const
Return face neighbour.
const vectorField & faceCentres() const
const vectorField & cellCentres() const
const cellList & cells() const
Holds list of sampling points which is filled at construction time. Various implementations of this b...
labelList faces_
Face numbers (-1 if not known)
void setSamples(const List< point > &samplingPts, const labelList &samplingCells, const labelList &samplingFaces, const labelList &samplingSegments, const scalarList &samplingDistance)
Set sample data. Copy list contents.
label pointInCell(const point &p, const label samplei) const
Return the cell in which the point on the sample line.
label getNeighbourCell(const label) const
Returns the neighbour cell or the owner if face in on the boundary.
scalar calcSign(const label facei, const point &sample) const
Calculates inproduct of face normal and vector sample-face centre.
bool getTrackingPoint(const point &samplePt, const point &bPoint, const label bFacei, const scalar smallDist, point &trackPt, label &trackCelli, label &trackFacei) const
Calculates start of tracking given samplePt and first boundary.
labelList segments_
Segment numbers.
labelList cells_
Cell numbers.
label getBoundaryCell(const label) const
Returns cell next to boundary face.
void checkDimensions() const
Check for consistent sizing.
point pushIn(const point &sample, const label facei) const
Moves sample in direction of -n to it is 'inside' of facei.
label findNearFace(const label celli, const point &sample, const scalar smallDist) const
Returns face label (or -1) of face which is close to sample.
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
#define FatalIOErrorInLookup(ios, lookupTag, lookupName, lookupTable)
Report an error message using Foam::FatalIOError.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
OBJstream os(runTime.globalPath()/outputName)
#define WarningInFunction
Report a warning using Foam::Warning.
#define DebugInFunction
Report an information message using Foam::Info.
dimensionedScalar sign(const dimensionedScalar &ds)
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
errorManip< error > abort(error &err)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
errorManipArg< error, int > exit(error &err, const int errNo=1)
constexpr char nl
The newline '\n' character (0x0a)
#define defineRunTimeSelectionTable(baseType, argNames)
Define run-time selection table.
#define forAll(list, i)
Loop across all elements in list.