zhyang-dev/OpenFOAM-12
1
0
Fork 0
Code Issues Pull Requests 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
parent 177d1d762c
commit e879da1341
2 changed files with 25 additions and 4 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)

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

@ -91,6 +91,16 @@ class particle
//- 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: