zhyang-dev/openfoam
1
0
Fork 0
mirror of https://develop.openfoam.com/Development/openfoam.git synced 2025-11-28 03:28:01 +00:00
Code Issues Actions Packages Projects Releases Wiki Activity

lagrangian: Fixed infinite loops

Browse Source 
Tracking through an inverted region of the mesh happens in a reversed
direction relative to a non-inverted region. Usually, this allows the
tracking to propagate normally, regardless of the sign of the space.
However, in rare cases, it is possible for a straight trajectory to form
a closed loop through both positive and negative regions. This causes
the tracking to loop indefinitely.

To fix this, the displacement through inverted regions has been
artifically increased by a small amount (1% at the moment). This has the
effect that the change in track fraction over the negative part of the
loop no longer exactly cancels the change over the positive part, and
the track therefore terminates.
This commit is contained in:
Will Bainbridge
2017-07-18 09:19:32 +01:00
committed by Andrew Heather
parent e274e08b1a
commit a4762ea6ea
2 changed files with 29 additions and 8 deletions
Download Patch File Download Diff File Expand all files Collapse all files

19
src/lagrangian/basic/particle/particle.C
View File

@ -30,6 +30,8 @@ License
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
const Foam::scalar Foam::particle::negativeSpaceDisplacementFactor = 1.01;
Foam::label Foam::particle::particleCount_ = 0;
namespace Foam
@ -685,8 +687,11 @@ Foam::scalar Foam::particle::trackToStationaryTri
<< "Start local coordinates = " << y0 << endl;
}
// Get the factor by which the displacement is increased
const scalar f = detA >= 0 ? 1 : negativeSpaceDisplacementFactor;
// Calculate the local tracking displacement
barycentric Tx1(x1 & T);
barycentric Tx1(f*x1 & T);
if (debug)
{
@ -783,16 +788,22 @@ Foam::scalar Foam::particle::trackToMovingTri
FixedList<barycentricTensor, 3> T;
movingTetReverseTransform(fraction, centre, detA, T);
// Get the factor by which the displacement is increased
const scalar f = detA[0] >= 0 ? 1 : negativeSpaceDisplacementFactor;
// Get the relative global position
const vector x0Rel = x0 - centre[0];
const vector x1Rel = x1 - centre[1];
const vector x1Rel = f*x1 - centre[1];
// Form the determinant and hit equations
cubicEqn detAEqn(sqr(detA[0])*detA[3], detA[0]*detA[2], detA[1], 1);
const barycentric yC(1, 0, 0, 0);
const barycentric hitEqnA = ((x1Rel & T[2]) + detA[3]*yC)*sqr(detA[0]);
const barycentric hitEqnA =
((x1Rel & T[2]) + detA[3]*yC)*sqr(detA[0]);
const barycentric hitEqnB =
((x1Rel & T[1]) + (x0Rel & T[2]) + detA[2]*yC)*detA[0];
const barycentric hitEqnC = (x1Rel & T[0]) + (x0Rel & T[1]) + detA[1]*yC;
const barycentric hitEqnC =
((x1Rel & T[0]) + (x0Rel & T[1]) + detA[1]*yC);
const barycentric hitEqnD = y0;
FixedList<cubicEqn, 4> hitEqn;
forAll(hitEqn, i)

18
src/lagrangian/basic/particle/particle.H
View File

@ -88,12 +88,22 @@ class particle
//- Size in bytes of the position data
static const std::size_t sizeofPosition_;
//- Size in bytes of the fields
static const std::size_t sizeofFields;
//- The factor by which the displacement is increased when passing
// through negative space. This should be slightly bigger than one.
// This is needed as a straight trajectory can form a closed loop
// through regions of overlapping positive and negative space, leading
// to a track which never ends. By increasing the rate of displacement
// through negative regions, the change in track fraction over this
// part of the loop no longer exactly cancels the change over the
// positive part, and the track therefore terminates.
static const scalar negativeSpaceDisplacementFactor;
public:
//- Size in bytes of the fields
static const std::size_t sizeofFields;
template<class CloudType>
class TrackingData
{