ENH: Don't always fatalError on a bad tet.

Some commented out debug code.  Smaller tolerance on tet quality.
This commit is contained in:
graham
2010-11-25 16:35:50 +00:00
parent 95bc057204
commit da8e668d0a
4 changed files with 59 additions and 11 deletions

View File

@ -29,8 +29,9 @@ License
// Note: the use of this tolerance is ad-hoc, there may be extreme // Note: the use of this tolerance is ad-hoc, there may be extreme
// cases where the resulting tetrahedra still have particle tracking // cases where the resulting tetrahedra still have particle tracking
// problems. // problems, or tets with lower quality may track OK.
const Foam::scalar Foam::polyMeshTetDecomposition::minTetQuality = SMALL; // const Foam::scalar Foam::polyMeshTetDecomposition::minTetQuality = SMALL;
const Foam::scalar Foam::polyMeshTetDecomposition::minTetQuality = ROOTVSMALL;
// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * // // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
@ -104,11 +105,10 @@ Foam::label Foam::polyMeshTetDecomposition::findSharedBasePoint
{ {
return faceBasePtI; return faceBasePtI;
} }
} }
// If a base point hasn't triggered a return by now, then there is // If a base point hasn't triggered a return by now, then there is
// non that can produce a good decomposition // none that can produce a good decomposition
return -1; return -1;
} }
@ -554,15 +554,17 @@ Foam::List<Foam::tetIndices> Foam::polyMeshTetDecomposition::faceTetIndices
if (tetBasePtI == -1) if (tetBasePtI == -1)
{ {
FatalErrorIn WarningIn
( (
"Foam::List<Foam::FixedList<Foam::label, 4> >" "Foam::List<Foam::FixedList<Foam::label, 4> >"
"Foam::Cloud<ParticleType>::" "Foam::Cloud<ParticleType>::"
"faceTetIndices(label fI, label cI) const" "faceTetIndices(label fI, label cI) const"
) )
<< "No base point for face " << fI << ", " << f << "No base point for face " << fI << ", " << f
<< ", produces a valid tet decomposition." << ", produces a valid tet decomposition."
<< abort(FatalError); << endl;
tetBasePtI = 0;
} }
for (label tetPtI = 1; tetPtI < f.size() - 1; tetPtI++) for (label tetPtI = 1; tetPtI < f.size() - 1; tetPtI++)

View File

@ -253,7 +253,11 @@ public:
//- Increment the nTrackingRescues counter //- Increment the nTrackingRescues counter
void trackingRescue() const void trackingRescue() const
{ {
nTrackingRescues_++; if (++nTrackingRescues_ % size() == 0)
{
Info<< " " << nTrackingRescues_
<< " tracking rescues " << endl;
}
} }
//- Whether each cell has any wall faces (demand driven data) //- Whether each cell has any wall faces (demand driven data)

View File

@ -264,6 +264,17 @@ Foam::scalar Foam::Particle<ParticleType>::trackToFace
faceI_ = -1; faceI_ = -1;
// Pout<< "Particle " << origId_ << " " << origProc_
// << " Tracking from " << position_
// << " to " << endPosition
// << endl;
// Pout<< "stepFraction " << stepFraction_ << nl
// << "cellI " << cellI_ << nl
// << "tetFaceI " << tetFaceI_ << nl
// << "tetPtI " << tetPtI_
// << endl;
scalar trackFraction = 0.0; scalar trackFraction = 0.0;
// Minimum tetrahedron decomposition of each cell of the mesh into // Minimum tetrahedron decomposition of each cell of the mesh into
@ -405,6 +416,8 @@ Foam::scalar Foam::Particle<ParticleType>::trackToFace
triI = -1; triI = -1;
lambdaMin = VGREAT; lambdaMin = VGREAT;
// Pout<< "tris " << tris << endl;
// Sets a value for lambdaMin and faceI_ if a wall face is hit // Sets a value for lambdaMin and faceI_ if a wall face is hit
// by the track. // by the track.
hitWallFaces(position_, endPosition, lambdaMin, faceHitTetIs); hitWallFaces(position_, endPosition, lambdaMin, faceHitTetIs);
@ -469,6 +482,35 @@ Foam::scalar Foam::Particle<ParticleType>::trackToFace
faceI_ = -1; faceI_ = -1;
} }
// Pout<< "track loop " << position_ << " " << endPosition << nl
// << " " << cellI_
// << " " << faceI_
// << " " << tetFaceI_
// << " " << tetPtI_
// << " " << triI
// << " " << lambdaMin
// << " " << trackFraction
// << endl;
// Pout<< "# Tracking loop tet "
// << origId_ << " " << origProc_<< nl
// << "# face: " << tetFaceI_ << nl
// << "# tetPtI: " << tetPtI_ << nl
// << "# tetBasePtI: " << mesh.tetBasePtIs()[tetFaceI_] << nl
// << "# tet.mag(): " << tet.mag() << nl
// << "# tet.quality(): " << tet.quality()
// << endl;
// meshTools::writeOBJ(Pout, tet.a());
// meshTools::writeOBJ(Pout, tet.b());
// meshTools::writeOBJ(Pout, tet.c());
// meshTools::writeOBJ(Pout, tet.d());
// Pout<< "f 1 3 2" << nl
// << "f 2 3 4" << nl
// << "f 1 4 3" << nl
// << "f 1 2 4" << endl;
// The particle can be 'outside' the tet. This will yield a // The particle can be 'outside' the tet. This will yield a
// lambda larger than 1, or smaller than 0. For values < 0, // lambda larger than 1, or smaller than 0. For values < 0,
// the particle travels away from the tet and we don't move // the particle travels away from the tet and we don't move

View File

@ -343,9 +343,6 @@ inline void Foam::Particle<ParticleType>::tetNeighbour(label triI)
label tetBasePtI = mesh.tetBasePtIs()[tetFaceI_]; label tetBasePtI = mesh.tetBasePtIs()[tetFaceI_];
label facePtI = (tetPtI_ + tetBasePtI) % f.size();
label otherFacePtI = f.fcIndex(facePtI);
if (tetBasePtI == -1) if (tetBasePtI == -1)
{ {
FatalErrorIn FatalErrorIn
@ -357,6 +354,9 @@ inline void Foam::Particle<ParticleType>::tetNeighbour(label triI)
<< abort(FatalError); << abort(FatalError);
} }
label facePtI = (tetPtI_ + tetBasePtI) % f.size();
label otherFacePtI = f.fcIndex(facePtI);
switch (triI) switch (triI)
{ {
case 0: case 0: