diff --git a/src/functionObjects/field/streamLine/streamLine.C b/src/functionObjects/field/streamLine/streamLine.C index 881330886de1d7495af89856740223a1743fe3ff..f2bbacfdf09a9822e57aba68eee7ab60e94fb850 100644 --- a/src/functionObjects/field/streamLine/streamLine.C +++ b/src/functionObjects/field/streamLine/streamLine.C @@ -56,18 +56,35 @@ void Foam::functionObjects::streamLine::track() const sampledSet& seedPoints = sampledSetPoints(); - forAll(seedPoints, i) + forAll(seedPoints, seedi) { particles.addParticle ( new streamLineParticle ( mesh_, - seedPoints[i], - seedPoints.cells()[i], + seedPoints[seedi], + seedPoints.cells()[seedi], + (trackDir_ == trackDirType::FORWARD), lifeTime_ ) ); + + if (trackDir_ == trackDirType::BIDIRECTIONAL) + { + // Add additional particle for the forward bit of the track + particles.addParticle + ( + new streamLineParticle + ( + mesh_, + seedPoints[seedi], + seedPoints.cells()[seedi], + true, + lifeTime_ + ) + ); + } } label nSeeds = returnReduce(particles.size(), sumOp<label>()); @@ -99,7 +116,6 @@ void Foam::functionObjects::streamLine::track() vsInterp, vvInterp, UIndex, // index of U in vvInterp - trackForward_, // track in +u direction? nSubCycle_, // automatic track control:step through cells in steps? trackLength_, // fixed track length diff --git a/src/functionObjects/field/streamLine/streamLine.H b/src/functionObjects/field/streamLine/streamLine.H index 0f4885a7c747ebc237e6851c77d77ab5542e53dd..9457abcf66391820a8531d1ea362e920404cde19 100644 --- a/src/functionObjects/field/streamLine/streamLine.H +++ b/src/functionObjects/field/streamLine/streamLine.H @@ -45,7 +45,8 @@ Usage setFormat vtk; U U; - trackForward yes; + //Deprecated: trackForward yes; + direction bidirectional; // or forward/backward fields ( @@ -77,6 +78,7 @@ Usage setFormat | Output data type | yes | U | Tracking velocity field name | no | U fields | Fields to sample | yes | + direction | Direction to track | yes | lifetime | Maximum number of particle tracking steps | yes | trackLength | Tracking segment length | no | nSubCycle | Number of tracking steps per cell | no| diff --git a/src/functionObjects/field/streamLine/streamLineBase.C b/src/functionObjects/field/streamLine/streamLineBase.C index c1f31c200b8b44b553e87a00eb75f608e87cf09b..055ec83bfc4c5d55ca80adfc0b18b9f03607118d 100644 --- a/src/functionObjects/field/streamLine/streamLineBase.C +++ b/src/functionObjects/field/streamLine/streamLineBase.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2015-2018 OpenCFD Ltd. + \\ / A nd | Copyright (C) 2015-2019 OpenCFD Ltd. \\/ M anipulation | ------------------------------------------------------------------------------- | Copyright (C) 2015 OpenFOAM Foundation @@ -47,6 +47,18 @@ namespace functionObjects } +const Foam::Enum +< + Foam::functionObjects::streamLineBase::trackDirType +> +Foam::functionObjects::streamLineBase::trackDirTypeNames +({ + { trackDirType::FORWARD, "forward" }, + { trackDirType::BACKWARD, "backward" }, + { trackDirType::BIDIRECTIONAL, "bidirectional" } +}); + + // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // const Foam::word& @@ -886,7 +898,27 @@ bool Foam::functionObjects::streamLineBase::read(const dictionary& dict) Info<< " Employing velocity field " << UName_ << endl; - dict.readEntry("trackForward", trackForward_); + bool trackForward; + if (dict.readIfPresent("trackForward", trackForward)) + { + trackDir_ = + ( + trackForward + ? trackDirType::FORWARD + : trackDirType::BACKWARD + ); + + if (dict.found("direction")) + { + FatalIOErrorInFunction(dict) << "Cannot specify both " + << "\"trackForward\" and \"direction\"" + << exit(FatalIOError); + } + } + else + { + trackDir_ = trackDirTypeNames.get("direction", dict); + } dict.readEntry("lifeTime", lifeTime_); if (lifeTime_ < 1) { diff --git a/src/functionObjects/field/streamLine/streamLineBase.H b/src/functionObjects/field/streamLine/streamLineBase.H index 776696f7682327bd1295ee2dc749d087f8fa90d4..9ab8ec4998c64e47b42d42efbcc17229e3376c91 100644 --- a/src/functionObjects/field/streamLine/streamLineBase.H +++ b/src/functionObjects/field/streamLine/streamLineBase.H @@ -45,6 +45,7 @@ SourceFiles #include "writer.H" #include "indirectPrimitivePatch.H" #include "interpolation.H" +#include "Enum.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -65,6 +66,22 @@ class streamLineBase : public fvMeshFunctionObject { +public: + + // Public data types + + //- Enumeration defining the track direction + enum trackDirType : char + { + FORWARD, + BACKWARD, + BIDIRECTIONAL + }; + + //- Names for the trackDir + static const Enum<trackDirType> trackDirTypeNames; + + protected: //- Seed set engine @@ -85,8 +102,8 @@ protected: //- Interpolation scheme to use word interpolationScheme_; - //- Whether to use +u or -u - bool trackForward_; + //- Whether to use +u or -u or both + trackDirType trackDir_; //- Maximum lifetime (= number of cells) of particle label lifeTime_; diff --git a/src/functionObjects/field/streamLine/streamLineParticle.C b/src/functionObjects/field/streamLine/streamLineParticle.C index 385747ca1ae87ecf3ca97c38c5a5d28de7191b5c..d7c1a7db2689859ba408d7e619b716e227085e5b 100644 --- a/src/functionObjects/field/streamLine/streamLineParticle.C +++ b/src/functionObjects/field/streamLine/streamLineParticle.C @@ -86,10 +86,12 @@ Foam::streamLineParticle::streamLineParticle const polyMesh& mesh, const vector& position, const label celli, + const bool trackForward, const label lifeTime ) : particle(mesh, position, celli), + trackForward_(trackForward), lifeTime_(lifeTime) {} @@ -98,17 +100,19 @@ Foam::streamLineParticle::streamLineParticle ( const polyMesh& mesh, Istream& is, - bool readFields + bool readFields, + bool newFormat ) : - particle(mesh, is, readFields) + particle(mesh, is, readFields, newFormat) { if (readFields) { List<scalarList> sampledScalars; List<vectorList> sampledVectors; - is >> lifeTime_ >> sampledPositions_ >> sampledScalars + is >> trackForward_ >> lifeTime_ + >> sampledPositions_ >> sampledScalars >> sampledVectors; sampledScalars_.setSize(sampledScalars.size()); @@ -133,9 +137,11 @@ Foam::streamLineParticle::streamLineParticle ) : particle(p), + trackForward_(p.trackForward_), lifeTime_(p.lifeTime_), sampledPositions_(p.sampledPositions_), - sampledScalars_(p.sampledScalars_) + sampledScalars_(p.sampledScalars_), + sampledVectors_(p.sampledVectors_) {} @@ -168,7 +174,7 @@ bool Foam::streamLineParticle::move sampledPositions_.append(position()); vector U = interpolateFields(td, position(), cell(), face()); - if (!td.trackForward_) + if (!trackForward_) { U = -U; } @@ -366,11 +372,7 @@ void Foam::streamLineParticle::hitWallPatch void Foam::streamLineParticle::readFields(Cloud<streamLineParticle>& c) { -// if (!c.size()) -// { -// return; -// } - bool valid = c.size(); + const bool valid = c.size(); particle::readFields(c); @@ -433,6 +435,7 @@ void Foam::streamLineParticle::writeFields(const Cloud<streamLineParticle>& c) Foam::Ostream& Foam::operator<<(Ostream& os, const streamLineParticle& p) { os << static_cast<const particle&>(p) + << token::SPACE << p.trackForward_ << token::SPACE << p.lifeTime_ << token::SPACE << p.sampledPositions_ << token::SPACE << p.sampledScalars_ diff --git a/src/functionObjects/field/streamLine/streamLineParticle.H b/src/functionObjects/field/streamLine/streamLineParticle.H index 3763472bd1c33fe49d589815de6cd99552996aa3..bc6acb6bb72b620bb3e878cf3ba722f0344e8fd3 100644 --- a/src/functionObjects/field/streamLine/streamLineParticle.H +++ b/src/functionObjects/field/streamLine/streamLineParticle.H @@ -79,8 +79,6 @@ public: const label UIndex_; - const bool trackForward_; - const label nSubCycle_; const scalar trackLength_; @@ -101,7 +99,6 @@ public: const PtrList<interpolation<scalar>>& vsInterp, const PtrList<interpolation<vector>>& vvInterp, const label UIndex, - const bool trackForward, const label nSubCycle, const scalar trackLength, DynamicList<List<point>>& allPositions, @@ -113,7 +110,6 @@ public: vsInterp_(vsInterp), vvInterp_(vvInterp), UIndex_(UIndex), - trackForward_(trackForward), nSubCycle_(nSubCycle), trackLength_(trackLength), allPositions_(allPositions), @@ -127,6 +123,9 @@ private: // Private data + //- Whether particle transports with +U or -U + bool trackForward_; + //- Lifetime of particle. Particle dies when reaches 0. label lifeTime_; @@ -162,6 +161,7 @@ public: const polyMesh& c, const vector& position, const label celli, + const bool trackForward, const label lifeTime ); @@ -170,7 +170,8 @@ public: ( const polyMesh& c, Istream& is, - bool readFields = true + bool readFields = true, + bool newFormat = true ); //- Construct copy diff --git a/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLine.C b/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLine.C index f4e61cbdb4f2c8a6b6e45dd42021bd4eb6154df3..3d49eb97a90d88d01ee179df584b52b20a39e121 100644 --- a/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLine.C +++ b/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLine.C @@ -51,7 +51,8 @@ namespace functionObjects // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // -Foam::tetIndices Foam::functionObjects::wallBoundedStreamLine::findNearestTet +Foam::Tuple2<Foam::tetIndices, Foam::point> +Foam::functionObjects::wallBoundedStreamLine::findNearestTet ( const bitSet& isWallPatch, const point& seedPt, @@ -63,6 +64,7 @@ Foam::tetIndices Foam::functionObjects::wallBoundedStreamLine::findNearestTet label minFacei = -1; label minTetPti = -1; scalar minDistSqr = sqr(GREAT); + point nearestPt(GREAT, GREAT, GREAT); for (label facei : cFaces) { @@ -76,14 +78,17 @@ Foam::tetIndices Foam::functionObjects::wallBoundedStreamLine::findNearestTet for (label i = 2; i < f.size(); i++) { const point& thisPoint = mesh_.points()[f[fp]]; - label nextFp = f.fcIndex(fp); + const label nextFp = f.fcIndex(fp); const point& nextPoint = mesh_.points()[f[nextFp]]; const triPointRef tri(basePoint, thisPoint, nextPoint); - scalar d2 = magSqr(tri.centre() - seedPt); + //const scalar d2 = magSqr(tri.centre() - seedPt); + const pointHit nearInfo(tri.nearestPoint(seedPt)); + const scalar d2 = nearInfo.distance(); if (d2 < minDistSqr) { + nearestPt = nearInfo.rawPoint(); minDistSqr = d2; minFacei = facei; minTetPti = i-1; @@ -93,16 +98,37 @@ Foam::tetIndices Foam::functionObjects::wallBoundedStreamLine::findNearestTet } } - // Put particle in tet - return tetIndices + // Return tet and nearest point on wall triangle + return Tuple2<Foam::tetIndices, Foam::point> ( - celli, - minFacei, - minTetPti + tetIndices + ( + celli, + minFacei, + minTetPti + ), + nearestPt ); } +Foam::point Foam::functionObjects::wallBoundedStreamLine::pushIn +( + const triPointRef& tri, + const point& pt +) const +{ + //pointHit nearPt(tri.nearestPoint(pt)); + //if (nearPt.distance() > ROOTSMALL) + //{ + // FatalErrorInFunction << "tri:" << tri + // << " seed:" << pt << exit(FatalError); + //} + + return (1.0-ROOTSMALL)*pt+ROOTSMALL*tri.centre(); +} + + void Foam::functionObjects::wallBoundedStreamLine::track() { // Determine the 'wall' patches @@ -130,14 +156,18 @@ void Foam::functionObjects::wallBoundedStreamLine::track() const sampledSet& seedPoints = sampledSetPoints(); - forAll(seedPoints, i) + forAll(seedPoints, seedi) { - const label celli = seedPoints.cells()[i]; + const label celli = seedPoints.cells()[seedi]; if (celli != -1) { - const point& seedPt = seedPoints[i]; - tetIndices ids(findNearestTet(isWallPatch, seedPt, celli)); + const point& seedPt = seedPoints[seedi]; + Tuple2<tetIndices, point> nearestId + ( + findNearestTet(isWallPatch, seedPt, celli) + ); + const tetIndices& ids = nearestId.first(); if (ids.face() != -1 && isWallPatch[ids.face()]) { @@ -154,15 +184,38 @@ void Foam::functionObjects::wallBoundedStreamLine::track() new wallBoundedStreamLineParticle ( mesh_, - ids.faceTri(mesh_).centre(), + // Perturb seed point to be inside triangle + pushIn(ids.faceTri(mesh_), nearestId.second()), ids.cell(), ids.face(), // tetFace ids.tetPt(), -1, // not on a mesh edge -1, // not on a diagonal edge + (trackDir_ == trackDirType::FORWARD), // forward? lifeTime_ // lifetime ) ); + + if (trackDir_ == trackDirType::BIDIRECTIONAL) + { + // Additional particle for other half of track + particles.addParticle + ( + new wallBoundedStreamLineParticle + ( + mesh_, + // Perturb seed point to be inside triangle + pushIn(ids.faceTri(mesh_), nearestId.second()), + ids.cell(), + ids.face(), // tetFace + ids.tetPt(), + -1, // not on a mesh edge + -1, // not on a diagonal edge + true, // forward + lifeTime_ // lifetime + ) + ); + } } else { @@ -173,7 +226,7 @@ void Foam::functionObjects::wallBoundedStreamLine::track() } } - label nSeeds = returnReduce(particles.size(), sumOp<label>()); + const label nSeeds = returnReduce(particles.size(), sumOp<label>()); Log << type() << " : seeded " << nSeeds << " particles." << endl; @@ -204,7 +257,6 @@ void Foam::functionObjects::wallBoundedStreamLine::track() vsInterp, vvInterp, UIndex, // index of U in vvInterp - trackForward_, // track in +u direction? trackLength_, // fixed track length isWallPatch, // which faces are to follow diff --git a/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLine.H b/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLine.H index 1b64c517224238f519a9274c255a7043bbb18d96..078ed93b485a0954b43d1ec6f5f922e0e11f7dc6 100644 --- a/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLine.H +++ b/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLine.H @@ -46,7 +46,8 @@ Usage setFormat vtk; U UNear; - trackForward yes; + //Deprecated: trackForward yes; + direction bidirectional; // or forward/backward fields ( @@ -76,6 +77,7 @@ Usage setFormat | Output data type | yes | U | Tracking velocity field name | yes | fields | Fields to sample | yes | + direction | Direction to track | yes | lifetime | Maximum number of particle tracking steps | yes | trackLength | Tracking segment length | no | nSubCycle | Number of tracking steps per cell | no| @@ -139,13 +141,21 @@ protected: // Protected Member Functions //- Find wall tet on cell - tetIndices findNearestTet + Tuple2<tetIndices, point> findNearestTet ( const bitSet& isWallPatch, const point& seedPt, const label celli ) const; + //- Push a point a tiny bit towards the centre of the triangle it is in + // to avoid tracking problems + point pushIn + ( + const triPointRef& tri, + const point& pt + ) const; + public: diff --git a/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLineParticle.C b/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLineParticle.C index 7e5a9253803a99903d11cc8bb70c538fddbbeaad..1e17dff0177bd232040e9fcf5e9f8227bd62df8c 100644 --- a/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLineParticle.C +++ b/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLineParticle.C @@ -96,7 +96,7 @@ Foam::vector Foam::wallBoundedStreamLineParticle::sample { vector U = interpolateFields(td, localPosition_, cell(), face()); - if (!td.trackForward_) + if (!trackForward_) { U = -U; } @@ -127,6 +127,7 @@ Foam::wallBoundedStreamLineParticle::wallBoundedStreamLineParticle const label tetPti, const label meshEdgeStart, const label diagEdge, + const bool trackForward, const label lifeTime ) : @@ -140,6 +141,7 @@ Foam::wallBoundedStreamLineParticle::wallBoundedStreamLineParticle meshEdgeStart, diagEdge ), + trackForward_(trackForward), lifeTime_(lifeTime) {} @@ -159,7 +161,7 @@ Foam::wallBoundedStreamLineParticle::wallBoundedStreamLineParticle List<scalarList> sampledScalars; List<vectorList> sampledVectors; - is >> lifeTime_ + is >> trackForward_ >> lifeTime_ >> sampledPositions_ >> sampledScalars >> sampledVectors; sampledScalars_.setSize(sampledScalars.size()); @@ -184,6 +186,7 @@ Foam::wallBoundedStreamLineParticle::wallBoundedStreamLineParticle ) : wallBoundedParticle(p), + trackForward_(p.trackForward_), lifeTime_(p.lifeTime_), sampledPositions_(p.sampledPositions_), sampledScalars_(p.sampledScalars_), @@ -269,6 +272,7 @@ Foam::Ostream& Foam::operator<< ) { os << static_cast<const wallBoundedParticle&>(p) + << token::SPACE << p.trackForward_ << token::SPACE << p.lifeTime_ << token::SPACE << p.sampledPositions_ << token::SPACE << p.sampledScalars_ diff --git a/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLineParticle.H b/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLineParticle.H index 4a0f74b30a45ae768cc50806dedfef22dd8d5e4d..fc2cef0c95bcbbd1ae492cd0df68b8fffaf80ebf 100644 --- a/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLineParticle.H +++ b/src/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLineParticle.H @@ -82,7 +82,6 @@ public: const PtrList<interpolation<scalar>>& vsInterp_; const PtrList<interpolation<vector>>& vvInterp_; const label UIndex_; - const bool trackForward_; const scalar trackLength_; DynamicList<vectorList>& allPositions_; @@ -99,7 +98,6 @@ public: const PtrList<interpolation<scalar>>& vsInterp, const PtrList<interpolation<vector>>& vvInterp, const label UIndex, - const bool trackForward, const scalar trackLength, const bitSet& isWallPatch, @@ -116,7 +114,6 @@ public: vsInterp_(vsInterp), vvInterp_(vvInterp), UIndex_(UIndex), - trackForward_(trackForward), trackLength_(trackLength), allPositions_(allPositions), @@ -132,6 +129,9 @@ protected: // Protected data + //- Track with +U or -U + bool trackForward_; + //- Lifetime of particle. Particle dies when reaches 0. label lifeTime_; @@ -172,6 +172,7 @@ public: const label tetPti, const label meshEdgeStart, const label diagEdge, + const bool trackForward, const label lifeTime );