Go to the documentation of this file.
39 namespace functionObjects
45 wallBoundedStreamLine,
62 const cell& cFaces = mesh_.cells()[celli];
66 scalar minDistSqr =
sqr(GREAT);
67 point nearestPt(GREAT, GREAT, GREAT);
69 for (label facei : cFaces)
71 if (isWallPatch[facei])
73 const face&
f = mesh_.faces()[facei];
74 const label fp0 = mesh_.tetBasePtIs()[facei];
75 const point& basePoint = mesh_.points()[
f[fp0]];
77 label fp =
f.fcIndex(fp0);
78 for (label i = 2; i <
f.size(); i++)
80 const point& thisPoint = mesh_.points()[
f[fp]];
81 const label nextFp =
f.fcIndex(fp);
82 const point& nextPoint = mesh_.points()[
f[nextFp]];
84 const triPointRef tri(basePoint, thisPoint, nextPoint);
88 const scalar d2 = nearInfo.
distance();
91 nearestPt = nearInfo.rawPoint();
128 return (1.0 - ROOTSMALL)*pt + ROOTSMALL*tri.
centre();
161 const label celli = seedPoints.
cells()[seedi];
165 const point& seedPt = seedPoints[seedi];
172 if (ids.
face() != -1 && isWallPatch[ids.
face()])
199 if (
trackDir_ == trackDirType::BIDIRECTIONAL)
222 Pout<<
type() <<
" : ignoring seed " << seedPt
223 <<
" since not in wall cell." <<
endl;
231 Log <<
type() <<
" : seeded " << nSeeds <<
" particles." <<
endl;
274 particles.
move(particles, td, trackTime);
317 faceSet faces(mesh_,
"lowQualityTetFaces", mesh_.nFaces()/100+1);
333 <<
"Found " << nFaces
334 <<
" faces with low quality or negative volume "
335 <<
"decomposition tets. Writing to faceSet " << faces.name()
341 for (
const cell& cFaces : mesh_.cells())
343 numFacesPerEdge.
clear();
345 for (
const label facei : cFaces)
347 const face&
f = mesh_.faces()[facei];
350 const edge e(
f[fp],
f.nextLabel(fp));
352 ++(numFacesPerEdge(
e, 0));
361 <<
"problem cell:" << cFaces
int debug
Static debugging option.
Holds list of sampling points which is filled at construction time. Various implementations of this b...
static bool checkFaceTets(const polyMesh &mesh, scalar tol=minTetQuality, const bool report=false, labelHashSet *setPtr=nullptr)
Check face-decomposition tet volume.
void addParticle(ParticleType *pPtr)
Transfer particle to cloud.
scalar trackLength_
Track length.
List< DynamicList< scalarList > > allScalars_
Per scalarField, per track, the sampled values.
static const scalar minTetQuality
Minimum tetrahedron quality.
wallBoundedStreamLine(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
label cell() const
Return the cell.
Describes the interaction of a face and a point. It carries the info of a successful hit and (if succ...
bool read(const char *buf, int32_t &val)
Same as readInt32.
const sampledSet & sampledSetPoints() const
Demand driven construction of the sampledSet.
const labelList & cells() const
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
label face() const
Return the face.
scalar distance() const noexcept
Return distance to hit.
autoPtr< indirectPrimitivePatch > wallPatch() const
Construct patch out of all wall patch faces.
Ostream & endl(Ostream &os)
Add newline and flush stream.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
virtual bool read(const dictionary &)
Read the field average data.
#define forAll(list, i)
Loop across all elements in list.
Point centre() const
Return centre (centroid)
virtual void track()
Do the actual tracking to fill the track data.
Like polyPatch but without reference to mesh. Used in boundaryMesh to hold data on patches....
A triangle primitive used to calculate face normals and swept volumes.
A Cloud of streamLine particles.
Class used to pass tracking data to the trackToEdge function.
triPointRef faceTri(const polyMesh &mesh) const
Return the geometry corresponding to the tri on the face for.
Particle class that samples fields as it passes through. Used in streamline calculation.
point pushIn(const triPointRef &tri, const point &pt) const
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
void initInterpolations(const label nSeeds, label &UIndex, PtrList< volScalarField > &vsFlds, PtrList< interpolation< scalar >> &vsInterp, PtrList< volVectorField > &vvFlds, PtrList< interpolation< vector >> &vvInterp)
Initialise fields, interpolators and track storage.
Template class for intrusive linked lists.
void move(TrackCloudType &cloud, typename ParticleType::trackingData &td, const scalar trackTime)
Move the particles.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
const wordRes fieldNames(propsDict.getOrDefault< wordRes >("fields", wordRes()))
Macros for easy insertion into run-time selection tables.
errorManip< error > abort(error &err)
List< DynamicList< vectorList > > allVectors_
Per vectorField, per track, the sampled values.
label tetPt() const
Return the characterising tetPtI.
trackDirType trackDir_
Whether to use +u or -u or both.
Storage and named access for the indices of a tet which is part of the decomposition of a cell.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
forAllConstIters(mixture.phases(), phase)
void clear()
Clear all entries from table.
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
word cloudName_
Optional specified name of particles.
virtual const word & type() const =0
Runtime type information.
dimensionedScalar sqrt(const dimensionedScalar &ds)
label lifeTime_
Maximum lifetime (= number of cells) of particle.
const fvMesh & mesh_
Reference to the fvMesh.
const dimensionedScalar e
Elementary charge.
defineTypeNameAndDebug(ObukhovLength, 0)
const T2 & second() const noexcept
Return second.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
label nFaces() const noexcept
Number of mesh faces.
A face is a list of labels corresponding to mesh vertices.
pointHit nearestPoint(const point &p) const
Return nearest point to p on triangle.
Tuple2< tetIndices, point > findNearestTet(const bitSet &isWallPatch, const point &seedPt, const label celli) const
Find wall tet on cell.
A 2-tuple for storing two objects of dissimilar types. The container is similar in purpose to std::pa...
virtual bool read(const dictionary &)
Read settings.
DynamicList< List< point > > allTracks_
All tracks. Per track the points it passed through.
const T1 & first() const noexcept
Return first.
#define WarningInFunction
Report a warning using Foam::Warning.
A cell is defined as a list of faces with extra functionality.