Back to index

salome-med  6.5.0
Public Member Functions | Private Member Functions | Private Attributes
INTERP_KERNEL::UnitTetraIntersectionBary Class Reference

#include <UnitTetraIntersectionBary.hxx>

Inheritance diagram for INTERP_KERNEL::UnitTetraIntersectionBary:
Inheritance graph
[legend]
Collaboration diagram for INTERP_KERNEL::UnitTetraIntersectionBary:
Collaboration graph
[legend]

List of all members.

Public Member Functions

INTERPKERNEL_EXPORT UnitTetraIntersectionBary (bool isTetraInversed=false)
 Creates a ready-to-use tool.
INTERPKERNEL_EXPORT void init (bool isTetraInversed=false)
 Initializes fields.
INTERPKERNEL_EXPORT void addSide (const TransformedTriangle &triangle)
 Stores a part of triangle common with the unit tetrahedron.
INTERPKERNEL_EXPORT bool getBary (double *baryCenter)
 Computes and return coordinates of barycentre.
INTERPKERNEL_EXPORT double getVolume () const
 Returns volume of intersection.
virtual INTERPKERNEL_EXPORT ~UnitTetraIntersectionBary ()
 Destructor clears coordinates of faces.

Private Member Functions

int addSideFaces ()
 Add faces of the intersection polyhedron formed on faces of the unit tetrahedron by sides of already added faces.
void setTriangleOnSide (int i)
 set corners of inherited TransformedTriangle as corners of i-th side of the Unit tetra.
void clearPolygons (bool andFaces=false)
 Free memory of polygons.

Private Attributes

double _int_volume
 volume of intersection
std::list< std::vector< double * > > _faces
 faces of intersection polyhedron
std::vector< std::vector
< double > > 
_polyNormals
bool _isTetraInversed

Detailed Description

Definition at line 35 of file UnitTetraIntersectionBary.hxx.


Constructor & Destructor Documentation

Creates a ready-to-use tool.

Definition at line 53 of file UnitTetraIntersectionBary.cxx.

    :TransformedTriangle(),_int_volume(0),_isTetraInversed( isTetraInversed )
  {
    //init();
  }

Destructor clears coordinates of faces.

Definition at line 750 of file UnitTetraIntersectionBary.cxx.

  {
    clearPolygons(/*andFaces=*/true );
  }

Here is the call graph for this function:


Member Function Documentation

Stores a part of triangle common with the unit tetrahedron.

Parameters:
triangle- triangle side of other cell, whose calculateIntersectionVolume() must have already been called
triangle- triangle side of other cell

Definition at line 79 of file UnitTetraIntersectionBary.cxx.

  {
    _int_volume += triangle.getVolume();

    double triNormal[3], polyNormal[3];
    crossprod<3>( triangle.getCorner(P),triangle.getCorner(Q),triangle.getCorner(R), triNormal);

    const std::vector<double*> * pPolygonA = &triangle.getPolygonA();
    if ( pPolygonA->size() < 3 )
      {
        if ( !epsilonEqual( triNormal[_Z], 0 ))
          return; // not vertical triangle does not intersect the unit tetra

        // Vertical triangle. Use inherited methods of TransformedTriangle to
        // calculate intesection polygon
        *((TransformedTriangle*)this) = triangle; // copy triangle fields
        _polygonA.clear();
        _polygonB.clear();
        calculateIntersectionAndProjectionPolygons();
        if (this->_polygonA.size() < 3)
          return;
        calculatePolygonBarycenter(A, _barycenterA);
        sortIntersectionPolygon(A, _barycenterA);
        pPolygonA = & _polygonA;
      }

    // check if polygon orientation is same as the one of triangle
    std::vector<double*>::const_iterator p = pPolygonA->begin(), pEnd = pPolygonA->end();
#ifdef DMP_UNITTETRAINTERSECTIONBARY
    std::cout.precision(18);
    std::cout << "**** int polygon() " << std::endl;
    while ( p != pEnd )
    {
      double* pp = *p++;
      std::cout << pEnd - p << ": ( " << pp[0] << ", " << pp[1] << ", " << pp[2] << " )" << std::endl;
    }
    p = pPolygonA->begin();
#endif
    double* p1 = *p++;
    double* p2 = *p;
    while ( samePoint( p1, p2 ) && ++p != pEnd )
      p2 = *p;
    if ( p == pEnd )
      {
#ifdef DMP_UNITTETRAINTERSECTIONBARY
        std::cout << "All points equal" << std::endl;
#endif
        clearPolygons();
        return;
      }
    double* p3 = *p;
    while (( samePoint( p2, p3 ) || samePoint( p1, p3 )) && ++p != pEnd )
      p3 = *p;
    if ( p == pEnd )
      {
#ifdef DMP_UNITTETRAINTERSECTIONBARY
        std::cout << "Only two points differ" << std::endl;
#endif
        clearPolygons();
        return ;
      }
    crossprod<3>( p1, p2, p3, polyNormal );
    bool reverse = ( dotprod<3>( triNormal, polyNormal ) < 0.0 );
    if (_isTetraInversed) reverse = !reverse;

    // store polygon
    _faces.push_back( std::vector< double* > () );
    std::vector< double* >& faceCorner = _faces.back();
    faceCorner.resize( pPolygonA->size()/* + 1*/ );

    int i = 0;
    if ( reverse )
      {
        std::vector<double*>::const_reverse_iterator polyF = pPolygonA->rbegin(), polyEnd;
        for ( polyEnd = pPolygonA->rend(); polyF != polyEnd; ++i, ++polyF )
          if ( i==0 || !samePoint( *polyF, faceCorner[i-1] ))
            copyVector3( *polyF, faceCorner[i] = new double[3] );
          else
            --i;
        polyNormal[0] *= -1.;
        polyNormal[1] *= -1.;
        polyNormal[2] *= -1.;
      }
    else
      {
        std::vector<double*>::const_iterator polyF = pPolygonA->begin(), polyEnd;
        for ( polyEnd = pPolygonA->end(); polyF != polyEnd; ++i, ++polyF )
          if ( i==0 || !samePoint( *polyF, faceCorner[i-1] ))
            copyVector3( *polyF, faceCorner[i] = new double[3] );
          else
            --i;
      }
    if ( i < 3 )
      {
        clearPolygons(); // free memory of _polygonA
        _polygonA = faceCorner;
        _faces.pop_back();
      }
    else
      {
        if ( i < (int)pPolygonA->size() )
          faceCorner.resize( i );

        if ( _polyNormals.empty() )
          _polyNormals.reserve(4);
        _polyNormals.push_back( std::vector< double >( polyNormal, polyNormal+3 ));
      }

#ifdef DMP_UNITTETRAINTERSECTIONBARY
    std::cout << "**** addSide() " << _faces.size() << std::endl;
    for ( int i = 0; i < faceCorner.size(); ++i )
      {
        double* p = faceCorner[i];
        std::cout << i << ": ( " << p[0] << ", " << p[1] << ", " << p[2] << " )" << std::endl;
      }
    std::cout << "NORM: ( " << _polyNormals.back()[0] << ", " << _polyNormals.back()[1] << ", " << _polyNormals.back()[2] << " )" << std::endl;
#endif
    clearPolygons();
  }

Here is the call graph for this function:

Here is the caller graph for this function:

Add faces of the intersection polyhedron formed on faces of the unit tetrahedron by sides of already added faces.

Return values:
int- number of faces of intersection polyhedron

Definition at line 300 of file UnitTetraIntersectionBary.cxx.

  {
    int nbPolyhedraFaces = 0;

    if ( _faces.empty() )
      return nbPolyhedraFaces;

    // -------------------------------------------
    // Detect polygons laying on sides of a tetra
    // -------------------------------------------

    bool sideAdded[NB_TETRA_SIDES] = { false, false, false, false };
    int nbAddedSides = 0;
    std::list< std::vector< double* > >::iterator f = _faces.begin(), fEnd = _faces.end();
    for ( ; f != fEnd; ++f )
    {
      std::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 ( epsilonEqual( coordSum[j], 0.0 ))
          sideAdded[j] = ++nbAddedSides != 0 ;
      }
      if ( !sideAdded[3] &&
           ( epsilonEqual( (coordSum[0]+coordSum[1]+coordSum[2]) / (int)polygon.size(), 1. )))
        sideAdded[3] = ++nbAddedSides != 0 ;
    }
    if ( nbAddedSides == NB_TETRA_SIDES )
      return nbAddedSides;

    // ---------------------------------------------------------------------------------
    // Add segments of already added polygons to future polygonal faces on sides of tetra
    // ---------------------------------------------------------------------------------

    std::size_t nbIntersectPolygs = _faces.size();

    std::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 ] )
      {
        _faces.push_back( std::vector< double* > () );
        sideFaces[ i ] = &_faces.back();
      }
    }
    f = _faces.begin(), fEnd = _faces.end();
    for ( std::size_t iF = 0; iF < nbIntersectPolygs; ++f, ++iF ) // loop on added intersection polygons
    {
      std::vector< double* >& polygon = *f;
      for ( std::size_t i = 0; i < 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 = epsilonEqual( p1[j], 0. );
          p2OnSide = epsilonEqual( 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() );
        }
      }
    }
#ifdef DMP_UNITTETRAINTERSECTIONBARY
    std::cout << "**** after Add segments to sides " << std::endl;
    for ( int i = 0; i < NB_TETRA_SIDES; ++i )
    {
      std::cout << "\t Side " << i << std::endl;
      if ( !sideFaces[i] )
      {
        std::cout << "\t cut by triagle" << std::endl;
      }
      else
      {
        std::vector< double* >& sideFace = *sideFaces[i];
        for ( int i = 0; i < sideFace.size(); ++i )
        {
          double* p = sideFace[i];
          std::cout << "\t" << i << ": ( " << p[0] << ", " << p[1] << ", " << p[2] << " )" << std::endl;
        }
      }
    }
#endif

    // ---------------------------------------------------------------------------
    // Make closed polygons on tetra sides by adding not cut off corners of tetra
    // ---------------------------------------------------------------------------

    double origin[3] = { 0,0,0 };

    // find corners of tetra cut off by triangles of other tetra
    // ---------------------------------------------------------

    // corners are coded like this: index = 1*X + 2*Y + 3*Z
    // (0,0,0) -> index == 0; (0,0,1) -> index == 3
    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 ic = 0; ic < NB_TETRA_NODES; ++ic )
    {
      f = _faces.begin(), fEnd = _faces.end();
      for ( std::size_t iF = 0; iF < nbIntersectPolygs; ++f, ++iF ) // loop on added intersection polygons
      {
        std::vector< double* >& polygon = *f;

        double corner2Poly[3] = { polygon[0][0], polygon[0][1], polygon[0][2] };
        if ( ic ) corner2Poly[ ic-1 ] -= 1.0;

        // _polyNormals are outside of a tetrahedron
        double dot = dotprod<3>( corner2Poly, &_polyNormals[iF][0] );
        if ( dot < -DEFAULT_ABS_TOL*DEFAULT_ABS_TOL )
        {
#ifdef DMP_UNITTETRAINTERSECTIONBARY
          std::cout << "side " << iF+1 << ": cut " << ic << std::endl;
#endif
          cutOffCorners[ ic ] = true;
          nbCutOffCorners++;
          break;
        }
      }
    }

    for ( int i = 0; i < 3; ++i ) // loop on orthogonal faces of the unit tetra
    {
      if ( !sideFaces[i] ) continue;
      std::vector< double* >& sideFace = *sideFaces[i];

      std::size_t 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=false;
      for ( std::size_t ip = 0; ip < nbPoints; ++ip )
      {
        std::size_t 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 ( epsilonEqual( p[ind1], 0.))
        {
          // point is on orthogonal edge
          if ( !isSegmentEnd && epsilonEqual( p2[ind1], 0. ))
            isSegmentOnEdge = true;

          if ( !isSegmentOnEdge )
          { // segment ends are on different edges
            pCut[ind2] = (int)isSegmentEnd; // believe that cutting triangles are well oriented
            cutOffIndex = pCut[0] + 2*pCut[1] + 3*pCut[2];
          }
          if ( epsilonEqual( p[ind2], 0.) || epsilonEqual( p[ind2], 1.))
          {
            pPass[ind2] = ( p[ind2] < 0.5 ) ? 0 : 1;
            passThIndex = pPass[0] + 2*pPass[1] + 3*pPass[2];
          }
        }
        else if ( epsilonEqual( p[ind2], 0.))
        {
          // point is on orthogonal edge
          if ( !isSegmentEnd && epsilonEqual( p2[ind2], 0. ))
            isSegmentOnEdge = true;
          if ( !isSegmentEnd )
          {// segment ends are on different edges
            pCut[ind1] = 1-(int)isSegmentEnd;
            cutOffIndex = pCut[0] + 2*pCut[1] + 3*pCut[2];
          }
          if ( epsilonEqual( p[ind1], 0.) || epsilonEqual( p[ind1], 1.))
          {
            pPass[ind1] = ( p[ind1] < 0.5 ) ? 0 : 1;
            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] = (int)isSegmentEnd;
            pCut[ind2] = 1-(int)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;
#ifdef DMP_UNITTETRAINTERSECTIONBARY
            std::cout << "PASS THROUGH " << passThIndex << std::endl;
#endif
          }
        }
        if ( cutOffIndex >= 0 )
        {
          nbCutOnSide++;
          if ( !passedCorners[ cutOffIndex ] && !cutOffCorners[ cutOffIndex ] )
          {
            nbCutOffCorners++;
            cutOffCorners[ cutOffIndex ] = true;
          }
        }
      } // loop on points on a unit tetra side

      if ( nbCutOnSide == 0 && nbPoints <= 2 )
        continue; // one segment laying on edge at most

      if ( nbCutOffCorners == NB_TETRA_NODES )
        break; // all tetra corners are cut off

      if ( /*nbCutOnSide <= 2 &&*/ nbPoints >= 6 )
      {
        // 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 ] = ++nbCutOffCorners != 0 ;
      }

    }
    // 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 ] = ++nbCutOffCorners != 0 ;
    }

    // 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 )
      {
        if ( !sideFaces[ i ] ) continue;
        std::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;
          }
        }
      }
    }

#ifdef DMP_UNITTETRAINTERSECTIONBARY
    std::cout << "**** after Add corners to sides " << std::endl;
    for ( int i = 0; i < NB_TETRA_SIDES; ++i )
    {
      std::cout << "\t Side " << i << std::endl;
      if ( !sideFaces[i] ) {
        std::cout << "\t cut by triagle" << std::endl;
      }
      else 
      {
        std::vector< double* >& sideFace = *sideFaces[i];
        for ( int i = 0; i < sideFace.size(); ++i )
        {
          double* p = sideFace[i];
          std::cout << "\t" << i << ": ( " << p[0] << ", " << p[1] << ", " << p[2] << " )" << std::endl;
        }
      }
    }
    std::cout << "Cut off corners: ";
    if ( nbCutOffCorners == 0 )
      std::cout << "NO";
    else 
      for ( int ic = 0; ic < NB_TETRA_NODES; ++ic )
        std::cout << cutOffCorners[ ic ];
    std::cout << std::endl;
#endif
    // ------------------------------------------------------------------------
    // Sort corners of filled up faces on tetra sides and exclude equal points
    // ------------------------------------------------------------------------

    std::size_t iF = 0;
    for ( f = _faces.begin(); f != fEnd; ++f, ++iF )
    {
      std::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( (int)(iF-nbIntersectPolygs) );
          sortIntersectionPolygon( A, _barycenterA );
        }
        // exclude equal points
        std::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 ))
          {
            face.push_back( p );
            *v = 0;
          }
        }
      }
      if ( face.size() < 3 )
      { // size could decrease
        clearPolygons(); // free memory of _polygonA
        _polygonA = face;
        face.clear();
      }
      else
      {
        nbPolyhedraFaces++;
      }
    }
#ifdef DMP_UNITTETRAINTERSECTIONBARY
    std::cout << "**** after HEALING all faces " << std::endl;
    for (iF=0, f = _faces.begin(); f != fEnd; ++f, ++iF )
    {
      std::cout << "\t Side " << iF << std::endl;
      std::vector< double* >& sideFace = *f;
      for ( int i = 0; i < sideFace.size(); ++i )
      {
        double* p = sideFace[i];
        std::cout << "\t" << i << ": ( " << p[0] << ", " << p[1] << ", " << p[2] << " )" << std::endl;
      }
    }
#endif
    return nbPolyhedraFaces;
  }

Here is the call graph for this function:

Here is the caller graph for this function:

void INTERP_KERNEL::UnitTetraIntersectionBary::clearPolygons ( bool  andFaces = false) [private]

Free memory of polygons.

Definition at line 713 of file UnitTetraIntersectionBary.cxx.

  {
    for(std::vector<double*>::iterator it = _polygonA.begin() ; it != _polygonA.end() ; ++it)
      {  delete[] *it;
        *it = 0; 
      }
    for(std::vector<double*>::iterator it = _polygonB.begin() ; it != _polygonB.end() ; ++it)
      { 
        delete[] *it; 
        *it = 0; 
      }

    _polygonA.clear();
    _polygonB.clear();

    if ( andFaces )
      {
        std::list< std::vector< double* > >::iterator f = this->_faces.begin(), fEnd = this->_faces.end();
        for ( ; f != fEnd; ++f )
          {
            std::vector< double* >& polygon = *f;
            for(std::vector<double*>::iterator it = polygon.begin() ; it != polygon.end() ; ++it)
              { 
                delete[] *it;
                *it = 0;
              }
          }
        this->_faces.clear();
      }
  }

Here is the caller graph for this function:

Computes and return coordinates of barycentre.

Computes and returns coordinates of barycentre.

Definition at line 205 of file UnitTetraIntersectionBary.cxx.

  {
    baryCenter[0] = baryCenter[1] = baryCenter[2] = -1.0;
    if ( addSideFaces() < NB_TETRA_SIDES )
      {
        // tetra is not intersected
        if ( fabs(_int_volume) > 1e-10 )
          {
            // tetra is fully inside the other cell
            baryCenter[0] = baryCenter[1] = baryCenter[2] = 0.25;
            _int_volume = 0.16666666666666666;
            return true;
          }
        return false;
      }
    // Algo:
    // - pick up one point P among the summits of the polyhedron
    // - for each face of the polyhedron which does not contain the point :
    //   - compute the barycenter of the volume obtained by forming the "pyramid" with
    //     the face as a base and point P as a summit
    //   - compute the volume of the "pyramid"
    // - Add up all barycenter positions weighting them with the volumes.

    baryCenter[0] = baryCenter[1] = baryCenter[2] = 0.;

    std::list< std::vector< double* > >::iterator f = _faces.begin(), fEnd = _faces.end();
    double * PP = f->at(0);

    for ( ++f; f != fEnd; ++f )
      {
        std::vector< double* >& polygon = *f;
        if ( polygon.empty() )
          continue;

        bool pBelongsToPoly = false;
        std::vector<double*>::iterator v = polygon.begin(), vEnd = polygon.end();
        for ( ; !pBelongsToPoly && v != vEnd; ++v )
          pBelongsToPoly = samePoint( PP, *v );
        if ( pBelongsToPoly )
          continue;

        // Compute the barycenter of the volume. Barycenter of pyramid is on line
        // ( barycenter of polygon -> PP ) with 1/4 of pyramid height from polygon.

        double bary[] = { 0, 0, 0 };

        // base polygon bary
        for ( v = polygon.begin(); v != vEnd ; ++v )
          {
            double* p = *v;
            bary[0] += p[0];
            bary[1] += p[1];
            bary[2] += p[2];
          }
        bary[0] /= (int)polygon.size();
        bary[1] /= (int)polygon.size();
        bary[2] /= (int)polygon.size();

        // pyramid volume
        double vol = 0;
        for ( int i = 0; i < (int)polygon.size(); ++i )
          {
            double* p1 = polygon[i];
            double* p2 = polygon[(i+1)%polygon.size()];
            vol += std::fabs( calculateVolumeForTetra( p1, p2, bary, PP ));
          }

        // put bary on the line ( barycenter of polygon -> PP ) and multiply by volume
        baryCenter[0] += ( bary[0] * 0.75 + PP[0] * 0.25 ) * vol;
        baryCenter[1] += ( bary[1] * 0.75 + PP[1] * 0.25 ) * vol;
        baryCenter[2] += ( bary[2] * 0.75 + PP[2] * 0.25 ) * vol;
      }
    if ( _int_volume < 0. )
      _int_volume = -_int_volume;
    baryCenter[0] /= _int_volume;
    baryCenter[1] /= _int_volume;
    baryCenter[2] /= _int_volume;

#ifdef DMP_UNITTETRAINTERSECTIONBARY
    std::cout.precision(5);
    std::cout << "**** Barycenter " << baryCenter[0] <<", "<< baryCenter[1] <<", "<< baryCenter[2]
         << "\t **** Volume " << _int_volume << std::endl;
    std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << std::endl;
#endif
    return true;
  }

Here is the call graph for this function:

Returns volume of intersection.

Return values:
double-

Definition at line 57 of file UnitTetraIntersectionBary.hxx.

{ return _int_volume; }
void INTERP_KERNEL::UnitTetraIntersectionBary::init ( bool  isTetraInversed = false)

Initializes fields.

Definition at line 64 of file UnitTetraIntersectionBary.cxx.

  {
    _int_volume = 0;
    _isTetraInversed = isTetraInversed;
    _faces.clear();
    _polyNormals.clear();
  }

Here is the caller graph for this function:

set corners of inherited TransformedTriangle as corners of i-th side of the Unit tetra.

It is necessary to sort points of faces on sides of the unit tetrahedron using sortIntersectionPolygon(A)

Definition at line 695 of file UnitTetraIntersectionBary.cxx.

  {
    if ( iSide >= 3 )
      iSide = 0;
    for(int i = 0 ; i < 3 ; ++i) 
      {
        _coords[5*i] = _coords[5*i + 1] = _coords[5*i + 2] = 0.;
        if ( i != iSide )
          _coords[5*i + i] = 1.;
      }
  }

Here is the caller graph for this function:


Member Data Documentation

std::list< std::vector< double* > > INTERP_KERNEL::UnitTetraIntersectionBary::_faces [private]

faces of intersection polyhedron

Definition at line 73 of file UnitTetraIntersectionBary.hxx.

volume of intersection

Definition at line 70 of file UnitTetraIntersectionBary.hxx.

Definition at line 76 of file UnitTetraIntersectionBary.hxx.

std::vector< std::vector< double > > INTERP_KERNEL::UnitTetraIntersectionBary::_polyNormals [private]

Definition at line 74 of file UnitTetraIntersectionBary.hxx.


The documentation for this class was generated from the following files: