#define NB_TETRA_SIDES 4
#define NB_TETRA_NODES 4
+//#define DMP_UNITTETRAINTERSECTIONBARY
+
using namespace std;
namespace INTERP_KERNEL
inline bool samePoint( const double* p1, const double* p2 )
{
- return ( p1[0]==p2[0] && p1[1]==p2[1] && p1[2]==p2[2]);
+ //return ( p1[0]==p2[0] && p1[1]==p2[1] && p1[2]==p2[2]);
+ return ( epsilonEqual( p1[0], p2[0]) && epsilonEqual( p1[1], p2[1]) && epsilonEqual( p1[2], p2[2]));
}
//================================================================================
*/
//================================================================================
- UnitTetraIntersectionBary::UnitTetraIntersectionBary():TransformedTriangle()
+ UnitTetraIntersectionBary::UnitTetraIntersectionBary(bool isTetraInversed)
+ :TransformedTriangle(),_int_volume(0),_isTetraInversed( isTetraInversed )
{
- _int_volume = 0;
//init();
}
//================================================================================
*/
//================================================================================
- void UnitTetraIntersectionBary::init()
+ void UnitTetraIntersectionBary::init(bool isTetraInversed)
{
_int_volume = 0;
+ _isTetraInversed = isTetraInversed;
_faces.clear();
+ _polyNormals.clear();
}
//================================================================================
}
crossprod<3>( p1, p2, p3, polyNormal );
bool reverse = ( dotprod<3>( triNormal, polyNormal ) < 0.0 );
+ if (_isTetraInversed) reverse = !reverse;
// store polygon
_faces.push_back( vector< double* > () );
copyVector3( *polyF, faceCorner[i] = new double[3] );
else
--i;
+ polyNormal[0] *= -1.;
+ polyNormal[1] *= -1.;
+ polyNormal[2] *= -1.;
}
else
{
{
if ( i < (int)pPolygonA->size() )
faceCorner.resize( i );
+
+ if ( _polyNormals.empty() )
+ _polyNormals.reserve(4);
+ _polyNormals.push_back( vector< double >( polyNormal, polyNormal+3 ));
}
-#ifdef DMP_ADDSIDE
+#ifdef DMP_UNITTETRAINTERSECTIONBARY
cout << "**** addSide() " << _faces.size() << endl;
for ( int i = 0; i < faceCorner.size(); ++i )
{
baryCenter[1] /= _int_volume;
baryCenter[2] /= _int_volume;
+#ifdef DMP_UNITTETRAINTERSECTIONBARY
+ cout << "**** Barycenter " << baryCenter[0] <<", "<< baryCenter[1] <<", "<< baryCenter[2]
+ << "\t **** Volume " << _int_volume << endl;
+#endif
return true;
}
int nbAddedSides = 0;
list< vector< double* > >::iterator f = _faces.begin(), fEnd = _faces.end();
for ( ; f != fEnd; ++f )
+ {
+ vector< double* >& polygon = *f;
+ double coordSum[3] = {0,0,0};
+ for ( int i = 0; i < (int)polygon.size(); ++i )
{
- vector< double* >& polygon = *f;
- double coordSum[3] = {0,0,0};
- for ( int i = 0; i < (int)polygon.size(); ++i )
- {
- double* p = polygon[i];
- coordSum[0] += p[0];
- coordSum[1] += p[1];
- coordSum[2] += p[2];
- }
- for ( int j = 0; j < 3 && !sideAdded[j]; ++j )
- {
- if ( coordSum[j] == 0 )
- sideAdded[j] = bool( ++nbAddedSides );
- }
- if ( !sideAdded[3] &&
- ( epsilonEqual( (coordSum[0]+coordSum[1]+coordSum[2]) / polygon.size(), 1. )))
- sideAdded[3] = bool( ++nbAddedSides );
+ double* p = polygon[i];
+ coordSum[0] += p[0];
+ coordSum[1] += p[1];
+ coordSum[2] += p[2];
}
+ for ( int j = 0; j < 3 && !sideAdded[j]; ++j )
+ {
+ if ( coordSum[j] == 0 )
+ sideAdded[j] = bool( ++nbAddedSides );
+ }
+ if ( !sideAdded[3] &&
+ ( epsilonEqual( (coordSum[0]+coordSum[1]+coordSum[2]) / polygon.size(), 1. )))
+ sideAdded[3] = bool( ++nbAddedSides );
+ }
if ( nbAddedSides == NB_TETRA_SIDES )
return nbAddedSides;
vector< double* > * sideFaces[ 4 ]; // future polygons on sides of tetra
for ( int i = 0; i < NB_TETRA_SIDES; ++i )
+ {
+ sideFaces[ i ]=0;
+ if ( !sideAdded[ i ] )
{
- sideFaces[ i ]=0;
- if ( !sideAdded[ i ] )
- {
- _faces.push_back( vector< double* > () );
- sideFaces[ i ] = &_faces.back();
- }
+ _faces.push_back( vector< double* > () );
+ sideFaces[ i ] = &_faces.back();
}
+ }
f = _faces.begin(), fEnd = _faces.end();
for ( int iF = 0; iF < nbIntersectPolygs; ++f, ++iF ) // loop on added intersection polygons
+ {
+ vector< double* >& polygon = *f;
+ for ( int i = 0; i < (int)polygon.size(); ++i )
{
- vector< double* >& polygon = *f;
- for ( int i = 0; i < (int)polygon.size(); ++i )
+ // segment ends
+ double* p1 = polygon[i];
+ double* p2 = polygon[(i+1)%polygon.size()];
+ bool p1OnSide, p2OnSide;//, onZeroSide = false;
+ for ( int j = 0; j < 3; ++j )
+ {
+ if ( !sideFaces[ j ] )
+ continue;
+ p1OnSide = ( p1[j] == 0. );
+ p2OnSide = ( p2[j] == 0. );
+ if ( p1OnSide && p2OnSide )
{
- // segment ends
- double* p1 = polygon[i];
- double* p2 = polygon[(i+1)%polygon.size()];
- bool p1OnSide, p2OnSide;//, onZeroSide = false;
- for ( int j = 0; j < 3; ++j )
- {
- if ( !sideFaces[ j ] )
- continue;
- p1OnSide = ( p1[j] == 0 );
- p2OnSide = ( p2[j] == 0 );
- if ( p1OnSide && p2OnSide )
- {
- // segment p1-p2 is on j-th orthogonal side of tetra
- sideFaces[j]->push_back( new double[3] );
- copyVector3( p1, sideFaces[j]->back() );
- sideFaces[j]->push_back( new double[3] );
- copyVector3( p2, sideFaces[j]->back() );
- //break;
- }
- }
- // check if the segment p1-p2 is on the inclined side
- if ( sideFaces[3] &&
- epsilonEqual( p1[_X] + p1[_Y] + p1[_Z], 1.0 ) &&
- epsilonEqual( p2[_X] + p2[_Y] + p2[_Z], 1.0 ))
- {
- sideFaces[3]->push_back( new double[3] );
- copyVector3( p1, sideFaces[3]->back() );
- sideFaces[3]->push_back( new double[3] );
- copyVector3( p2, sideFaces[3]->back() );
- }
+ // segment p1-p2 is on j-th orthogonal side of tetra
+ sideFaces[j]->push_back( new double[3] );
+ copyVector3( p1, sideFaces[j]->back() );
+ sideFaces[j]->push_back( new double[3] );
+ copyVector3( p2, sideFaces[j]->back() );
+ //break;
}
+ }
+ // check if the segment p1-p2 is on the inclined side
+ if ( sideFaces[3] &&
+ epsilonEqual( p1[_X] + p1[_Y] + p1[_Z], 1.0 ) &&
+ epsilonEqual( p2[_X] + p2[_Y] + p2[_Z], 1.0 ))
+ {
+ sideFaces[3]->push_back( new double[3] );
+ copyVector3( p1, sideFaces[3]->back() );
+ sideFaces[3]->push_back( new double[3] );
+ copyVector3( p2, sideFaces[3]->back() );
+ }
}
-#ifdef DMP_ADDSIDEFACES
+ }
+#ifdef DMP_UNITTETRAINTERSECTIONBARY
cout << "**** after Add segments to sides " << endl;
for ( int i = 0; i < NB_TETRA_SIDES; ++i )
+ {
+ cout << "\t Side " << i << endl;
+ if ( !sideFaces[i] )
{
- cout << "\t Side " << i << endl;
- if ( !sideFaces[i] )
- {
- cout << "\t cut by triagle" << endl;
- }
- else
- {
- vector< double* >& sideFace = *sideFaces[i];
- for ( int i = 0; i < sideFace.size(); ++i )
- {
- double* p = sideFace[i];
- cout << "\t" << i << ": ( " << p[0] << ", " << p[1] << ", " << p[2] << " )" << endl;
- }
- }
+ cout << "\t cut by triagle" << endl;
+ }
+ else
+ {
+ vector< double* >& sideFace = *sideFaces[i];
+ for ( int i = 0; i < sideFace.size(); ++i )
+ {
+ double* p = sideFace[i];
+ cout << "\t" << i << ": ( " << p[0] << ", " << p[1] << ", " << p[2] << " )" << endl;
+ }
}
+ }
#endif
// ---------------------------------------------------------------------------
// corners are coded like this: index = 1*X + 2*Y + 3*Z
// (0,0,0) -> index == 0; (0,0,1) -> index == 3
- vector< int > cutOffCorners(NB_TETRA_NODES, false), passedCorners(NB_TETRA_NODES, false);
+ int cutOffCorners[NB_TETRA_NODES] = { false, false, false, false };
+ int passedCorners[NB_TETRA_NODES] = { false, false, false, false };
+ // find cutOffCorners by normals of intersection polygons
int nbCutOffCorners = 0;
- for ( int i = 0; i < 3; ++i ) // loop on orthogonal faces of the unit tetra
+ for ( int ic = 0; ic < NB_TETRA_NODES; ++ic )
+ {
+ //if ( !passedCorners[ ic ] && !cutOffCorners[ ic ])
{
- if ( !sideFaces[i] ) continue;
- vector< double* >& sideFace = *sideFaces[i];
+// double corner[3] = { 0., 0., 0. };
+// if ( ic ) corner[ ic-1 ] = 1.0;
- int nbPoints = sideFace.size();
- if ( nbPoints == 0 )
- continue; // not intersected face at all - no cut off corners can be detected
+ f = _faces.begin(), fEnd = _faces.end();
+ for ( int iF = 0; iF < nbIntersectPolygs; ++f, ++iF ) // loop on added intersection polygons
+ {
+ vector< double* >& polygon = *f;
- int ind1 = (i+1)%3, ind2 = (i+2)%3; // indices of coords on i-th tetra side
+ double corner2Poly[3] = { polygon[0][0], polygon[0][1], polygon[0][2] };
+ if ( ic ) corner2Poly[ ic-1 ] -= 1.0;
- int nbCutOnSide = 0;
- bool isSegmentOnEdge;
- for ( int ip = 0; ip < nbPoints; ++ip )
+ // _polyNormals are outside of a tetrahedron
+ double dot = dotprod<3>( corner2Poly, &_polyNormals[iF][0] );
+ if ( dot < -DEFAULT_ABS_TOL )
{
- int isSegmentEnd = ( ip % 2 );
-
- double* p = sideFace[ ip ];
- double* p2 = isSegmentEnd ? 0 : sideFace[ip+1];
-
- if ( !isSegmentEnd )
- isSegmentOnEdge = false; // initialize
-
- int cutOffIndex = -1, passThIndex = -1;// no cut off neither pass through
- int pCut[] = { 0,0,0 }, pPass[] = { 0,0,0 };
-
- if ( p[ind1]==0.)
- {
- // point is in on orthogonal edge
- if ( !isSegmentEnd && p2[ind1]==0 )
- isSegmentOnEdge = true;
-
- if ( !isSegmentOnEdge )
- { // segment ends are on different edges
- pCut[ind2] = isSegmentEnd; // believe that cutting triangles are well oriented
- cutOffIndex = pCut[0] + 2*pCut[1] + 3*pCut[2];
- }
- if ( p[ind2]==0. || p[ind2]==1.)
- {
- pPass[ind2] = int(p[ind2]);
- passThIndex = pPass[0] + 2*pPass[1] + 3*pPass[2];
- }
- }
- else if ( p[ind2]==0.)
- {
- // point is on orthogonal edge
- if ( !isSegmentEnd && p2[ind2]==0 )
- isSegmentOnEdge = true;
- if ( !isSegmentEnd )
- {// segment ends are on different edges
- pCut[ind1] = 1-isSegmentEnd;
- cutOffIndex = pCut[0] + 2*pCut[1] + 3*pCut[2];
- }
- if ( p[ind1]==0. || p[ind1]==1.)
- {
- pPass[ind1] = int(p[ind1]);
- passThIndex = pPass[0] + 2*pPass[1] + 3*pPass[2];
- }
- }
- else if ( epsilonEqual(p[ind1] + p[ind2], 1.0 ))
- {
- // point is on inclined edge
- if ( !isSegmentEnd && epsilonEqual(p2[ind1] + p2[ind2], 1.0 ))
- isSegmentOnEdge = true;
- if ( !isSegmentOnEdge )
- { //segment ends are on different edges
- pCut[ind1] = isSegmentEnd;
- pCut[ind2] = 1-isSegmentEnd;
- cutOffIndex = pCut[0] + 2*pCut[1] + 3*pCut[2];
- }
- }
- else
- {
- continue;
- }
- // remember cut off and passed through points
- if ( passThIndex >= 0. )
- {
- passedCorners[ passThIndex ] = true;
- if ( cutOffCorners[ passThIndex ] )
- {
- nbCutOffCorners--;
- cutOffCorners[ passThIndex ] = false;
- }
- }
- if ( cutOffIndex >= 0. )
- {
- nbCutOnSide++;
- if ( !passedCorners[ cutOffIndex ] && !cutOffCorners[ cutOffIndex ] )
- {
- nbCutOffCorners++;
- cutOffCorners[ cutOffIndex ] = true;
- }
- }
- } // loop on points on a unit tetra side
+ cutOffCorners[ ic ] = true;
+ nbCutOffCorners++;
+ break;
+ }
+ }
+ }
+ }
+
+ for ( int i = 0; i < 3; ++i ) // loop on orthogonal faces of the unit tetra
+ {
+ if ( !sideFaces[i] ) continue;
+ vector< double* >& sideFace = *sideFaces[i];
+
+ int nbPoints = sideFace.size();
+ if ( nbPoints == 0 )
+ continue; // not intersected face at all - no cut off corners can be detected
+
+ int ind1 = (i+1)%3, ind2 = (i+2)%3; // indices of coords on i-th tetra side
+
+ int nbCutOnSide = 0;
+ bool isSegmentOnEdge;
+ for ( int ip = 0; ip < nbPoints; ++ip )
+ {
+ int isSegmentEnd = ( ip % 2 );
+
+ double* p = sideFace[ ip ];
+ double* p2 = isSegmentEnd ? 0 : sideFace[ip+1];
- if ( nbCutOnSide == 0 && nbPoints <= 2 )
- continue; // one segment laying on edge at most
+ if ( !isSegmentEnd )
+ isSegmentOnEdge = false; // initialize
- if ( nbCutOffCorners == NB_TETRA_NODES )
- break; // all tetra corners are cut off
+ int cutOffIndex = -1, passThIndex = -1;// no cut off neither pass through
+ int pCut[] = { 0,0,0 }, pPass[] = { 0,0,0 };
- if ( /*nbCutOnSide <= 2 &&*/ nbPoints >= 6 )
+ if ( p[ind1]==0.)
+ {
+ // point is on orthogonal edge
+ if ( !isSegmentEnd && p2[ind1]==0. )
+ isSegmentOnEdge = true;
+
+ if ( !isSegmentOnEdge )
+ { // segment ends are on different edges
+ pCut[ind2] = isSegmentEnd; // believe that cutting triangles are well oriented
+ cutOffIndex = pCut[0] + 2*pCut[1] + 3*pCut[2];
+ }
+ if ( p[ind2]==0. || p[ind2]==1.)
+ {
+ pPass[ind2] = int(p[ind2]);
+ passThIndex = pPass[0] + 2*pPass[1] + 3*pPass[2];
+ }
+ }
+ else if ( p[ind2]==0.)
+ {
+ // point is on orthogonal edge
+ if ( !isSegmentEnd && p2[ind2]==0. )
+ isSegmentOnEdge = true;
+ if ( !isSegmentEnd )
+ {// segment ends are on different edges
+ pCut[ind1] = 1-isSegmentEnd;
+ cutOffIndex = pCut[0] + 2*pCut[1] + 3*pCut[2];
+ }
+ if ( p[ind1]==0. || p[ind1]==1.)
+ {
+ pPass[ind1] = int(p[ind1]);
+ passThIndex = pPass[0] + 2*pPass[1] + 3*pPass[2];
+ }
+ }
+ else if ( epsilonEqual(p[ind1] + p[ind2], 1.0 ))
+ {
+ // point is on inclined edge
+ if ( !isSegmentEnd && epsilonEqual(p2[ind1] + p2[ind2], 1.0 ))
+ isSegmentOnEdge = true;
+ if ( !isSegmentOnEdge )
+ { //segment ends are on different edges
+ pCut[ind1] = isSegmentEnd;
+ pCut[ind2] = 1-isSegmentEnd;
+ cutOffIndex = pCut[0] + 2*pCut[1] + 3*pCut[2];
+ }
+ }
+ else
+ {
+ continue;
+ }
+ // remember cut off and passed through points
+ if ( passThIndex >= 0. )
+ {
+ passedCorners[ passThIndex ] = true;
+ if ( cutOffCorners[ passThIndex ] )
+ {
+ nbCutOffCorners--;
+ cutOffCorners[ passThIndex ] = false;
+ }
+ }
+ if ( cutOffIndex >= 0. )
+ {
+ nbCutOnSide++;
+ if ( !passedCorners[ cutOffIndex ] && !cutOffCorners[ cutOffIndex ] )
{
- // at least 3 segments - all corners of a side are cut off
- for (int cutIndex = 0; cutIndex < NB_TETRA_NODES; ++cutIndex )
- if ( cutIndex != i+1 && !passedCorners[ cutIndex ] && !cutOffCorners[ cutIndex ])
- cutOffCorners[ cutIndex ] = bool ( ++nbCutOffCorners );
+ nbCutOffCorners++;
+ cutOffCorners[ cutOffIndex ] = true;
}
+ }
+ } // loop on points on a unit tetra side
- } // loop on orthogonal faces of tetra
+ if ( nbCutOnSide == 0 && nbPoints <= 2 )
+ continue; // one segment laying on edge at most
- // check if all corners are cut off on the inclined tetra side
- if ( sideFaces[ XYZ ] && sideFaces[ XYZ ]->size() >= 6 )
+ if ( nbCutOffCorners == NB_TETRA_NODES )
+ break; // all tetra corners are cut off
+
+ if ( /*nbCutOnSide <= 2 &&*/ nbPoints >= 6 )
{
- for (int cutIndex = 1; cutIndex < NB_TETRA_NODES; ++cutIndex )
- if ( !passedCorners[ cutIndex ] && !cutOffCorners[ cutIndex ])
+ // at least 3 segments - all corners of a side are cut off
+ for (int cutIndex = 0; cutIndex < NB_TETRA_NODES; ++cutIndex )
+ if ( cutIndex != i+1 && !passedCorners[ cutIndex ] && !cutOffCorners[ cutIndex ])
cutOffCorners[ cutIndex ] = bool ( ++nbCutOffCorners );
}
+ }
+ // loop on orthogonal faces of tetra
+
+ // check if all corners are cut off on the inclined tetra side
+ if ( sideFaces[ XYZ ] && sideFaces[ XYZ ]->size() >= 6 )
+ {
+ for (int cutIndex = 1; cutIndex < NB_TETRA_NODES; ++cutIndex )
+ if ( !passedCorners[ cutIndex ] && !cutOffCorners[ cutIndex ])
+ cutOffCorners[ cutIndex ] = bool ( ++nbCutOffCorners );
+ }
+
// Add to faces on tetra sides the corners not cut off by segments of intersection polygons
// ----------------------------------------------------------------------------------
if ( nbCutOffCorners > 0 )
+ {
+ for ( int i = 0; i < NB_TETRA_SIDES; ++i )
{
- for ( int i = 0; i < NB_TETRA_SIDES; ++i )
+ if ( !sideFaces[ i ] ) continue;
+ vector< double* >& sideFace = *sideFaces[i];
+
+ int excludeCorner = (i + 1) % NB_TETRA_NODES;
+ for ( int ic = 0; ic < NB_TETRA_NODES; ++ic )
+ {
+ if ( !cutOffCorners[ ic ] && ic != excludeCorner )
{
- if ( !sideFaces[ i ] ) continue;
- vector< double* >& sideFace = *sideFaces[i];
-
- int excludeCorner = (i + 1) % NB_TETRA_NODES;
- for ( int ic = 0; ic < NB_TETRA_NODES; ++ic )
- {
- if ( !cutOffCorners[ ic ] && ic != excludeCorner )
- {
- sideFace.push_back( new double[3] );
- copyVector3( origin, sideFace.back() );
- if ( ic )
- sideFace.back()[ ic-1 ] = 1.0;
- }
- }
+ sideFace.push_back( new double[3] );
+ copyVector3( origin, sideFace.back() );
+ if ( ic )
+ sideFace.back()[ ic-1 ] = 1.0;
}
+ }
}
+ }
-#ifdef DMP_ADDSIDEFACES
+#ifdef DMP_UNITTETRAINTERSECTIONBARY
cout << "**** after Add corners to sides " << endl;
for ( int i = 0; i < NB_TETRA_SIDES; ++i )
+ {
+ cout << "\t Side " << i << endl;
+ if ( !sideFaces[i] ) {
+ cout << "\t cut by triagle" << endl;
+ }
+ else
{
- cout << "\t Side " << i << endl;
- if ( !sideFaces[i] ) {
- cout << "\t cut by triagle" << endl;
+ vector< double* >& sideFace = *sideFaces[i];
+ for ( int i = 0; i < sideFace.size(); ++i )
+ {
+ double* p = sideFace[i];
+ cout << "\t" << i << ": ( " << p[0] << ", " << p[1] << ", " << p[2] << " )" << endl;
}
- else
- {
- vector< double* >& sideFace = *sideFaces[i];
- for ( int i = 0; i < sideFace.size(); ++i )
- {
- double* p = sideFace[i];
- cout << "\t" << i << ": ( " << p[0] << ", " << p[1] << ", " << p[2] << " )" << endl;
- }
- }
}
+ }
cout << "Cut off corners: ";
if ( nbCutOffCorners == 0 )
cout << "NO";
int iF = 0;
for ( f = _faces.begin(); f != fEnd; ++f, ++iF )
+ {
+ vector< double* >& face = *f;
+ if ( face.size() >= 3 )
{
- vector< double* >& face = *f;
- if ( face.size() >= 3 )
+ clearPolygons(); // free memory of _polygonA
+ _polygonA = face;
+ face.clear();
+ face.reserve( _polygonA.size() );
+ if ( iF >= nbIntersectPolygs )
+ { // sort points of side faces
+ calculatePolygonBarycenter( A, _barycenterA );
+ setTriangleOnSide( iF - nbIntersectPolygs );
+ sortIntersectionPolygon( A, _barycenterA );
+ }
+ // exclude equal points
+ vector< double* >::iterator v = _polygonA.begin(), vEnd = _polygonA.end();
+ face.push_back( *v );
+ *v = 0;
+ for ( ++v; v != vEnd; ++v )
+ {
+ double* pPrev = face.back();
+ double* p = *v;
+ if ( !samePoint( p, pPrev ))
{
- clearPolygons(); // free memory of _polygonA
- _polygonA = face;
- face.clear();
- face.reserve( _polygonA.size() );
- if ( iF >= nbIntersectPolygs )
- { // sort points of side faces
- calculatePolygonBarycenter( A, _barycenterA );
- setTriangleOnSide( iF - nbIntersectPolygs );
- sortIntersectionPolygon( A, _barycenterA );
- }
- // exclude equal points
- vector< double* >::iterator v = _polygonA.begin(), vEnd = _polygonA.end();
- face.push_back( *v );
+ face.push_back( p );
*v = 0;
- for ( ++v; v != vEnd; ++v )
- {
- double* pPrev = face.back();
- double* p = *v;
- if ( !samePoint( p, pPrev ))
- {
- face.push_back( p );
- *v = 0;
- }
- }
- }
- if ( face.size() < 3 )
- { // size could decrease
- clearPolygons(); // free memory of _polygonA
- _polygonA = face;
- face.clear();
- }
- else
- {
- nbPolyhedraFaces++;
}
+ }
+ }
+ if ( face.size() < 3 )
+ { // size could decrease
+ clearPolygons(); // free memory of _polygonA
+ _polygonA = face;
+ face.clear();
+ }
+ else
+ {
+ nbPolyhedraFaces++;
}
-#ifdef DMP_ADDSIDEFACES
+ }
+#ifdef DMP_UNITTETRAINTERSECTIONBARY
cout << "**** after HEALING all faces " << endl;
for (iF=0, f = _faces.begin(); f != fEnd; ++f, ++iF )
+ {
+ cout << "\t Side " << iF << endl;
+ vector< double* >& sideFace = *f;
+ for ( int i = 0; i < sideFace.size(); ++i )
{
- cout << "\t Side " << iF << endl;
- vector< double* >& sideFace = *f;
- for ( int i = 0; i < sideFace.size(); ++i )
- {
- double* p = sideFace[i];
- cout << "\t" << i << ": ( " << p[0] << ", " << p[1] << ", " << p[2] << " )" << endl;
- }
+ double* p = sideFace[i];
+ cout << "\t" << i << ": ( " << p[0] << ", " << p[1] << ", " << p[2] << " )" << endl;
}
+ }
#endif
return nbPolyhedraFaces;
}
