diff --git a/src/lagrangian/basic/particle/particle.C b/src/lagrangian/basic/particle/particle.C index c71559d..861b7ec 100644 --- a/src/lagrangian/basic/particle/particle.C +++ b/src/lagrangian/basic/particle/particle.C @@ -61,7 +61,8 @@ Foam::particle::particle tetFaceI_(tetFaceI), tetPtI_(tetPtI), origProc_(Pstream::myProcNo()), - origId_(getNewParticleID()) + origId_(getNewParticleID()), + markedToTeleport_(false) {} @@ -81,7 +82,8 @@ Foam::particle::particle tetFaceI_(-1), tetPtI_(-1), origProc_(Pstream::myProcNo()), - origId_(getNewParticleID()) + origId_(getNewParticleID()), + markedToTeleport_(false) { if (doCellFacePt) { @@ -100,7 +102,8 @@ Foam::particle::particle(const particle& p) tetFaceI_(p.tetFaceI_), tetPtI_(p.tetPtI_), origProc_(p.origProc_), - origId_(p.origId_) + origId_(p.origId_), + markedToTeleport_(false) {} @@ -114,7 +117,8 @@ Foam::particle::particle(const particle& p, const polyMesh& mesh) tetFaceI_(p.tetFaceI_), tetPtI_(p.tetPtI_), origProc_(p.origProc_), - origId_(p.origId_) + origId_(p.origId_), + markedToTeleport_(false) {} diff --git a/src/lagrangian/basic/particle/particle.H b/src/lagrangian/basic/particle/particle.H index 6cb1838..9d9b099 100644 --- a/src/lagrangian/basic/particle/particle.H +++ b/src/lagrangian/basic/particle/particle.H @@ -154,6 +154,8 @@ protected: //- Local particle id on originating processor label origId_; + + bool markedToTeleport_; // Private Member Functions @@ -457,6 +459,9 @@ public: //- Is the particle on the boundary/(or outside the domain)? inline bool onBoundary() const; + inline bool needsTeleportHandled() const; + inline void resetTeleportHandler(); + //- Is this global face an internal face? inline bool internalFace(const label faceI) const; @@ -511,7 +516,7 @@ public: // on entry 'stepFraction()' should be set to the fraction of the // time-step at which the tracking starts. template - scalar trackToFace(const vector& endPosition, TrackData& td); + scalar trackToFace(const vector& endPosition, TrackData& td, bool handleTeleport=false); //- Return the index of the face to be used in the interpolation // routine diff --git a/src/lagrangian/basic/particle/particleI.H b/src/lagrangian/basic/particle/particleI.H index 803a200..38a0d0b 100644 --- a/src/lagrangian/basic/particle/particleI.H +++ b/src/lagrangian/basic/particle/particleI.H @@ -850,7 +850,19 @@ inline void Foam::particle::initCellFacePt() inline bool Foam::particle::onBoundary() const { - return faceI_ != -1 && faceI_ >= mesh_.nInternalFaces(); + return faceI_ != -1 && faceI_ >= mesh_.nInternalFaces() && !markedToTeleport_; +} + + +inline bool Foam::particle::needsTeleportHandled() const +{ + return markedToTeleport_; +} + + +inline void Foam::particle::resetTeleportHandler() +{ + markedToTeleport_ = false; } diff --git a/src/lagrangian/basic/particle/particleTemplates.C b/src/lagrangian/basic/particle/particleTemplates.C index 280165d..7e44e58 100644 --- a/src/lagrangian/basic/particle/particleTemplates.C +++ b/src/lagrangian/basic/particle/particleTemplates.C @@ -203,7 +203,8 @@ template Foam::scalar Foam::particle::trackToFace ( const vector& endPosition, - TrackData& td + TrackData& td, + bool handleTeleport ) { typedef typename TrackData::cloudType cloudType; @@ -448,6 +449,9 @@ Foam::scalar Foam::particle::trackToFace // This must be a cell face crossing faceI_ = tetFaceI_; + Info << "tetFaceI_: " << tetFaceI_ << endl; + Info << "tri, lambdaMin: " << triI << " " << lambdaMin << endl; + // Set the faceHitTetIs to those for the current tet in case a // wall interaction is required with the cell face faceHitTetIs = tetIndices @@ -576,14 +580,24 @@ Foam::scalar Foam::particle::trackToFace patchI = patch(faceI_); } + Info << "origFaceI: " << origFaceI << endl; + Info << "faceI_: " << faceI_ << endl; + const polyPatch& patch = mesh_.boundaryMesh()[patchI]; if (isA(patch)) { - p.hitWedgePatch - ( - static_cast(patch), td - ); + if(!handleTeleport || markedToTeleport_) + { + p.hitWedgePatch + ( + static_cast(patch), td + ); + } + else + { + markedToTeleport_ = true; + } } else if (isA(patch)) { @@ -601,26 +615,47 @@ Foam::scalar Foam::particle::trackToFace } else if (isA(patch)) { - p.hitCyclicPatch - ( - static_cast(patch), td - ); + if(!handleTeleport || markedToTeleport_) + { + p.hitCyclicPatch + ( + static_cast(patch), td + ); + } + else + { + markedToTeleport_ = true; + } } else if (isA(patch)) { - p.hitCyclicAMIPatch - ( - static_cast(patch), - td, - endPosition - position_ - ); + if(!handleTeleport || markedToTeleport_) + { + p.hitCyclicAMIPatch + ( + static_cast(patch), + td, + endPosition - position_ + ); + } + else + { + markedToTeleport_ = true; + } } else if (isA(patch)) { - p.hitProcessorPatch - ( - static_cast(patch), td - ); + if(!handleTeleport || markedToTeleport_) + { + p.hitProcessorPatch + ( + static_cast(patch), td + ); + } + else + { + markedToTeleport_ = true; + } } else if (isA(patch)) { @@ -633,9 +668,14 @@ Foam::scalar Foam::particle::trackToFace { p.hitPatch(patch, td); } + + Info << "faceI_ After: " << faceI_ << endl; + } } + Info << "A faceI_: " << faceI_ << endl; + if (lambdaMin < SMALL) { // Apply tracking correction towards tet centre. @@ -708,6 +748,8 @@ Foam::scalar Foam::particle::trackToFace cloud.trackingRescue(); } + Info << "V faceI_: " << faceI_ << endl; + return trackFraction; } diff --git a/src/sampling/sampledSet/face/faceOnlySet.C b/src/sampling/sampledSet/face/faceOnlySet.C index 85fb7bc..8d0d19b 100644 --- a/src/sampling/sampledSet/face/faceOnlySet.C +++ b/src/sampling/sampledSet/face/faceOnlySet.C @@ -60,14 +60,19 @@ bool Foam::faceOnlySet::trackToBoundary // Alias const point& trackPt = singleParticle.position(); + bool mustResetTeleport = false; while(true) { point oldPoint = trackPt; - singleParticle.trackToFace(end_, trackData); - - if (singleParticle.face() != -1 && mag(oldPoint - trackPt) > smallDist) + singleParticle.trackToFace(end_, trackData, true); + + if ( + singleParticle.face() != -1 && + (mag(oldPoint - trackPt) > smallDist + || singleParticle.needsTeleportHandled()) + ) { // Reached face. Sample. samplingPts.append(trackPt); @@ -86,6 +91,18 @@ bool Foam::faceOnlySet::trackToBoundary // Boundary reached. return true; } + if(singleParticle.needsTeleportHandled()) + { + if(mustResetTeleport) + { + singleParticle.resetTeleportHandler(); + mustResetTeleport = false; + } + else + { + mustResetTeleport = true; + } + } } } diff --git a/src/sampling/sampledSet/midPoint/midPointSet.C b/src/sampling/sampledSet/midPoint/midPointSet.C index 168d0bb..7e76bfe 100644 --- a/src/sampling/sampledSet/midPoint/midPointSet.C +++ b/src/sampling/sampledSet/midPoint/midPointSet.C @@ -62,11 +62,18 @@ void Foam::midPointSet::genSamples() { midPoints[midI] = 0.5*(operator[](sampleI) + operator[](sampleI+1)); + + Info << "faces_[sampleI]: " << faces_[sampleI] << endl; + Info << "faces_[sampleI+1]: " << faces_[sampleI+1] << endl; label cell1 = getCell(faces_[sampleI], midPoints[midI]); label cell2 = getCell(faces_[sampleI+1], midPoints[midI]); - if (cell1 != cell2) + if (cell1 == -1 || cell2 == -1) + { + Info << "Ignoring intermediate jump" << endl; + } + else if (cell1 != cell2) { FatalErrorIn("midPointSet::genSamples()") << " sampleI:" << sampleI @@ -80,12 +87,14 @@ void Foam::midPointSet::genSamples() << " cell2:" << cell2 << abort(FatalError); } + else + { + midCells[midI] = cell1; + midSegments[midI] = segments_[sampleI]; + midCurveDist[midI] = mag(midPoints[midI] - start()); + midI++; + } - midCells[midI] = cell1; - midSegments[midI] = segments_[sampleI]; - midCurveDist[midI] = mag(midPoints[midI] - start()); - - midI++; sampleI++; } diff --git a/src/sampling/sampledSet/sampledSet/sampledSet.C b/src/sampling/sampledSet/sampledSet/sampledSet.C index 1442ca3..2283610 100644 --- a/src/sampling/sampledSet/sampledSet/sampledSet.C +++ b/src/sampling/sampledSet/sampledSet/sampledSet.C @@ -30,6 +30,11 @@ License #include "writer.H" #include "particle.H" +#include "cyclicPolyPatch.H" +#include "cyclicAMIPolyPatch.H" +#include "processorPolyPatch.H" +#include "wedgePolyPatch.H" + // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // namespace Foam @@ -64,18 +69,33 @@ Foam::label Foam::sampledSet::getCell << abort(FatalError); } - if (faceI >= mesh().nInternalFaces()) + if (!mesh().isInternalFace(faceI)) { label cellI = getBoundaryCell(faceI); if (!mesh().pointInCell(sample, cellI, searchEngine_.decompMode())) { - FatalErrorIn - ( - "sampledSet::getCell(const label, const point&)" - ) << "Found cell " << cellI << " using face " << faceI - << ". But cell does not contain point " << sample - << abort(FatalError); + label patchI = mesh_.boundaryMesh().whichPatch(faceI); + const polyPatch& patch = mesh_.boundaryMesh()[patchI]; + if( + isA(patch) || + isA(patch) || + isA(patch) || + isA(patch) + ) + { + return -1; + } + else + { + + FatalErrorIn + ( + "sampledSet::getCell(const label, const point&)" + ) << "Found cell " << cellI << " using face " << faceI + << ". But cell does not contain point " << sample + << abort(FatalError); + } } return cellI; }