-// Copyright (C) 2007-2010 CEA/DEN, EDF R&D, OPEN CASCADE
+// Copyright (C) 2007-2011 CEA/DEN, EDF R&D, OPEN CASCADE
//
-// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
-// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
+// Copyright (C) 2003-2007 OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
+// CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
//
-// This library is free software; you can redistribute it and/or
-// modify it under the terms of the GNU Lesser General Public
-// License as published by the Free Software Foundation; either
-// version 2.1 of the License.
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License.
//
-// This library is distributed in the hope that it will be useful,
-// but WITHOUT ANY WARRANTY; without even the implied warranty of
-// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
-// Lesser General Public License for more details.
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+// Lesser General Public License for more details.
//
-// You should have received a copy of the GNU Lesser General Public
-// License along with this library; if not, write to the Free Software
-// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
//
-// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
//
// File : SMDS_VolumeTool.cxx
#include "SMDS_MeshElement.hxx"
#include "SMDS_MeshNode.hxx"
-#include "SMDS_PolyhedralVolumeOfNodes.hxx"
+#include "SMDS_VtkVolume.hxx"
#include "SMDS_Mesh.hxx"
#include "utilities.h"
#include <map>
-#include <float.h>
-#include <math.h>
+#include <limits>
+#include <cmath>
using namespace std;
inline XYZ Crossed( const XYZ& other );
inline double Dot( const XYZ& other );
inline double Magnitude();
+ inline double SquareMagnitude();
};
inline XYZ XYZ::operator-( const XYZ& Right ) {
return XYZ(x - Right.x, y - Right.y, z - Right.z);
inline double XYZ::Magnitude() {
return sqrt (x * x + y * y + z * z);
}
+inline double XYZ::SquareMagnitude() {
+ return (x * x + y * y + z * z);
+}
}
//=======================================================================
myFaceNodeIndices( NULL ),
myFaceNodes( NULL )
{
+ //MESSAGE("******************************************************** SMDS_VolumeToo");
}
//=======================================================================
myFaceNodeIndices( NULL ),
myFaceNodes( NULL )
{
+ //MESSAGE("******************************************************** SMDS_VolumeToo");
Set( theVolume );
}
}
if (myVolume->IsPoly()) {
- myPolyedre = static_cast<const SMDS_PolyhedralVolumeOfNodes*>( myVolume );
+ myPolyedre = dynamic_cast<const SMDS_VtkVolume*>( myVolume );
if (!myPolyedre) {
MESSAGE("Warning: bad volumic element");
return false;
if ( !myPolyedre )
return 0.;
+ SMDS_Mesh *mesh = SMDS_Mesh::_meshList[myPolyedre->getMeshId()];
// split a polyhedron into tetrahedrons
SMDS_VolumeTool* me = const_cast< SMDS_VolumeTool* > ( this );
XYZ baryCenter;
me->GetBaryCenter(baryCenter.x, baryCenter.y, baryCenter.z);
- SMDS_MeshNode bcNode ( baryCenter.x, baryCenter.y, baryCenter.z );
+ SMDS_MeshNode *bcNode = mesh->AddNode( baryCenter.x, baryCenter.y, baryCenter.z );
for ( int f = 0; f < NbFaces(); ++f )
{
double Vn = getTetraVolume( myFaceNodes[ 0 ],
myFaceNodes[ n-1 ],
myFaceNodes[ n ],
- & bcNode );
+ bcNode );
/// cout <<"++++ " << Vn << " nodes " <<myFaceNodes[ 0 ]->GetID() << " " <<myFaceNodes[ n-1 ]->GetID() << " " <<myFaceNodes[ n ]->GetID() << " < " << V << endl;
V += externalFace ? -Vn : Vn;
}
}
+ mesh->RemoveNode(bcNode);
}
else
{
*/
//================================================================================
-bool SMDS_VolumeTool::IsOut(double X, double Y, double Z, double tol)
+bool SMDS_VolumeTool::IsOut(double X, double Y, double Z, double tol) const
{
// LIMITATION: for convex volumes only
XYZ p( X,Y,Z );
//purpose : Return number of nodes in the array of face nodes
//=======================================================================
-int SMDS_VolumeTool::NbFaceNodes( int faceIndex )
+int SMDS_VolumeTool::NbFaceNodes( int faceIndex ) const
{
if ( !setFace( faceIndex ))
return 0;
// the last node == the first one.
//=======================================================================
-const SMDS_MeshNode** SMDS_VolumeTool::GetFaceNodes( int faceIndex )
+const SMDS_MeshNode** SMDS_VolumeTool::GetFaceNodes( int faceIndex ) const
{
if ( !setFace( faceIndex ))
return 0;
// the last node index == the first one.
//=======================================================================
-const int* SMDS_VolumeTool::GetFaceNodesIndices( int faceIndex )
+const int* SMDS_VolumeTool::GetFaceNodesIndices( int faceIndex ) const
{
if ( !setFace( faceIndex ))
return 0;
if (myVolume->IsPoly())
{
- myPolyIndices.resize( myFaceNbNodes + 1 );
- myFaceNodeIndices = & myPolyIndices[0];
+ SMDS_VolumeTool* me = const_cast< SMDS_VolumeTool* > ( this );
+ me->myPolyIndices.resize( myFaceNbNodes + 1 );
+ me->myFaceNodeIndices = & me->myPolyIndices[0];
for ( int i = 0; i <= myFaceNbNodes; ++i )
- myFaceNodeIndices[i] = myVolume->GetNodeIndex( myFaceNodes[i] );
+ me->myFaceNodeIndices[i] = myVolume->GetNodeIndex( myFaceNodes[i] );
}
return myFaceNodeIndices;
}
//=======================================================================
bool SMDS_VolumeTool::GetFaceNodes (int faceIndex,
- set<const SMDS_MeshNode*>& theFaceNodes )
+ set<const SMDS_MeshNode*>& theFaceNodes ) const
{
if ( !setFace( faceIndex ))
return false;
//=======================================================================
//function : IsFaceExternal
-//purpose : Check normal orientation of a returned face
+//purpose : Check normal orientation of a given face
//=======================================================================
-bool SMDS_VolumeTool::IsFaceExternal( int faceIndex )
+bool SMDS_VolumeTool::IsFaceExternal( int faceIndex ) const
{
if ( myExternalFaces || !myVolume )
return true;
//purpose : Return a normal to a face
//=======================================================================
-bool SMDS_VolumeTool::GetFaceNormal (int faceIndex, double & X, double & Y, double & Z)
+bool SMDS_VolumeTool::GetFaceNormal (int faceIndex, double & X, double & Y, double & Z) const
{
if ( !setFace( faceIndex ))
return false;
}
double size = cross.Magnitude();
- if ( size <= DBL_MIN )
+ if ( size <= numeric_limits<double>::min() )
return false;
X = cross.x / size;
*/
//================================================================================
-bool SMDS_VolumeTool::GetFaceBaryCenter (int faceIndex, double & X, double & Y, double & Z)
+bool SMDS_VolumeTool::GetFaceBaryCenter (int faceIndex, double & X, double & Y, double & Z) const
{
if ( !setFace( faceIndex ))
return false;
//purpose : Return face area
//=======================================================================
-double SMDS_VolumeTool::GetFaceArea( int faceIndex )
+double SMDS_VolumeTool::GetFaceArea( int faceIndex ) const
{
if (myVolume->IsPoly()) {
MESSAGE("Warning: attempt to obtain area of a face of polyhedral volume");
if ( minInd >= firstMediumInd )
return false; // medium nodes are not linked
if ( maxInd < firstMediumInd ) // both nodes are corners
+ {
if ( theIgnoreMediumNodes )
type = SMDSAbs_EntityType( int(type)-1 ); // check linkage of corner nodes
else
return false; // corner nodes are not linked directly in a quadratic cell
+ }
}
switch ( type ) {
*/
//================================================================================
-int SMDS_VolumeTool::GetAllExistingFaces(vector<const SMDS_MeshElement*> & faces)
+int SMDS_VolumeTool::GetAllExistingFaces(vector<const SMDS_MeshElement*> & faces) const
{
faces.clear();
faces.reserve( NbFaces() );
return edges.size();
}
-//=======================================================================
-//function : IsFreeFace
-//purpose : check that only one volume is build on the face nodes
-//=======================================================================
+//================================================================================
+/*!
+ * \brief Return minimal square distance between connected corner nodes
+ */
+//================================================================================
-bool SMDS_VolumeTool::IsFreeFace( int faceIndex )
+double SMDS_VolumeTool::MinLinearSize2() const
{
- const int free = true;
+ double minSize = 1e+100;
+ int iQ = myVolume->IsQuadratic() ? 2 : 1;
+
+ // store current face data
+ int curFace = myCurFace, nbN = myFaceNbNodes;
+ int* ind = myFaceNodeIndices;
+ myFaceNodeIndices = NULL;
+ const SMDS_MeshNode** nodes = myFaceNodes;
+ myFaceNodes = NULL;
+
+ // it seems that compute distance twice is faster than organization of a sole computing
+ myCurFace = -1;
+ for ( int iF = 0; iF < myNbFaces; ++iF )
+ {
+ setFace( iF );
+ for ( int iN = 0; iN < myFaceNbNodes; iN += iQ )
+ {
+ XYZ n1( myFaceNodes[ iN ]);
+ XYZ n2( myFaceNodes[(iN + iQ) % myFaceNbNodes]);
+ minSize = std::min( minSize, (n1 - n2).SquareMagnitude());
+ }
+ }
+ // restore current face data
+ myCurFace = curFace;
+ myFaceNbNodes = nbN;
+ myFaceNodeIndices = ind;
+ delete [] myFaceNodes; myFaceNodes = nodes;
+
+ return minSize;
+}
+
+//================================================================================
+/*!
+ * \brief check that only one volume is build on the face nodes
+ *
+ * If a face is shared by one of <ignoreVolumes>, it is considered free
+ */
+//================================================================================
+
+bool SMDS_VolumeTool::IsFreeFace( int faceIndex, const SMDS_MeshElement** otherVol/*=0*/ ) const
+{
+ const bool isFree = true;
if (!setFace( faceIndex ))
- return !free;
+ return !isFree;
const SMDS_MeshNode** nodes = GetFaceNodes( faceIndex );
- int nbFaceNodes = myFaceNbNodes;
+ const int nbFaceNodes = myFaceNbNodes;
// evaluate nb of face nodes shared by other volume
int maxNbShared = -1;
TElemIntMap::iterator vNbIt;
for ( int iNode = 0; iNode < nbFaceNodes; iNode++ ) {
const SMDS_MeshNode* n = nodes[ iNode ];
- SMDS_ElemIteratorPtr eIt = n->GetInverseElementIterator();
+ SMDS_ElemIteratorPtr eIt = n->GetInverseElementIterator( SMDSAbs_Volume );
while ( eIt->more() ) {
const SMDS_MeshElement* elem = eIt->next();
- if ( elem != myVolume && elem->GetType() == SMDSAbs_Volume ) {
- int nbShared = 1;
- vNbIt = volNbShared.find( elem );
- if ( vNbIt == volNbShared.end() ) {
- volNbShared.insert ( TElemIntMap::value_type( elem, nbShared ));
- }
- else {
- nbShared = ++(*vNbIt).second;
- }
- if ( nbShared > maxNbShared )
- maxNbShared = nbShared;
+ if ( elem != myVolume ) {
+ vNbIt = volNbShared.insert( make_pair( elem, 0 )).first;
+ (*vNbIt).second++;
+ if ( vNbIt->second > maxNbShared )
+ maxNbShared = vNbIt->second;
}
}
}
if ( maxNbShared < 3 )
- return free; // is free
+ return isFree; // is free
// find volumes laying on the opposite side of the face
// and sharing all nodes
if ( IsFaceExternal( faceIndex ))
intNormal = XYZ( -intNormal.x, -intNormal.y, -intNormal.z );
XYZ p0 ( nodes[0] ), baryCenter;
- for ( vNbIt = volNbShared.begin(); vNbIt != volNbShared.end(); vNbIt++ ) {
- int nbShared = (*vNbIt).second;
+ for ( vNbIt = volNbShared.begin(); vNbIt != volNbShared.end(); ) {
+ const int& nbShared = (*vNbIt).second;
if ( nbShared >= 3 ) {
SMDS_VolumeTool volume( (*vNbIt).first );
volume.GetBaryCenter( baryCenter.x, baryCenter.y, baryCenter.z );
XYZ intNormal2( baryCenter - p0 );
- if ( intNormal.Dot( intNormal2 ) < 0 )
- continue; // opposite side
+ if ( intNormal.Dot( intNormal2 ) < 0 ) {
+ // opposite side
+ if ( nbShared >= nbFaceNodes )
+ {
+ // a volume shares the whole facet
+ if ( otherVol ) *otherVol = vNbIt->first;
+ return !isFree;
+ }
+ ++vNbIt;
+ continue;
+ }
}
// remove a volume from volNbShared map
- volNbShared.erase( vNbIt-- );
+ volNbShared.erase( vNbIt++ );
}
- // here volNbShared contains only volumes laying on the
- // opposite side of the face
- if ( volNbShared.empty() ) {
- return free; // is free
+ // here volNbShared contains only volumes laying on the opposite side of
+ // the face and sharing 3 or more but not all face nodes with myVolume
+ if ( volNbShared.size() < 2 ) {
+ return isFree; // is free
}
// check if the whole area of a face is shared
- bool isShared[] = { false, false, false, false }; // 4 triangle parts of a quadrangle
- for ( vNbIt = volNbShared.begin(); vNbIt != volNbShared.end(); vNbIt++ ) {
- SMDS_VolumeTool volume( (*vNbIt).first );
- bool prevLinkShared = false;
- int nbSharedLinks = 0;
- for ( int iNode = 0; iNode < nbFaceNodes; iNode++ ) {
- bool linkShared = volume.IsLinked( nodes[ iNode ], nodes[ iNode + 1] );
- if ( linkShared )
- nbSharedLinks++;
- if ( linkShared && prevLinkShared &&
- volume.IsLinked( nodes[ iNode - 1 ], nodes[ iNode + 1] ))
- isShared[ iNode ] = true;
- prevLinkShared = linkShared;
- }
- if ( nbSharedLinks == nbFaceNodes )
- return !free; // is not free
- if ( nbFaceNodes == 4 ) {
- // check traingle parts 1 & 3
- if ( isShared[1] && isShared[3] )
- return !free; // is not free
- // check triangle parts 0 & 2;
- // 0 part could not be checked in the loop; check it here
- if ( isShared[2] && prevLinkShared &&
- volume.IsLinked( nodes[ 0 ], nodes[ 1 ] ) &&
- volume.IsLinked( nodes[ 1 ], nodes[ 3 ] ) )
- return !free; // is not free
- }
+ for ( int iNode = 0; iNode < nbFaceNodes; iNode++ )
+ {
+ const SMDS_MeshNode* n = nodes[ iNode ];
+ // check if n is shared by one of volumes of volNbShared
+ bool isShared = false;
+ SMDS_ElemIteratorPtr eIt = n->GetInverseElementIterator( SMDSAbs_Volume );
+ while ( eIt->more() && !isShared )
+ isShared = volNbShared.count( eIt->next() );
+ if ( !isShared )
+ return isFree;
}
- return free;
+ if ( otherVol ) *otherVol = volNbShared.begin()->first;
+ return !isFree;
+
+// if ( !myVolume->IsPoly() )
+// {
+// bool isShared[] = { false, false, false, false }; // 4 triangle parts of a quadrangle
+// for ( vNbIt = volNbShared.begin(); vNbIt != volNbShared.end(); vNbIt++ ) {
+// SMDS_VolumeTool volume( (*vNbIt).first );
+// bool prevLinkShared = false;
+// int nbSharedLinks = 0;
+// for ( int iNode = 0; iNode < nbFaceNodes; iNode++ ) {
+// bool linkShared = volume.IsLinked( nodes[ iNode ], nodes[ iNode + 1] );
+// if ( linkShared )
+// nbSharedLinks++;
+// if ( linkShared && prevLinkShared &&
+// volume.IsLinked( nodes[ iNode - 1 ], nodes[ iNode + 1] ))
+// isShared[ iNode ] = true;
+// prevLinkShared = linkShared;
+// }
+// if ( nbSharedLinks == nbFaceNodes )
+// return !free; // is not free
+// if ( nbFaceNodes == 4 ) {
+// // check traingle parts 1 & 3
+// if ( isShared[1] && isShared[3] )
+// return !free; // is not free
+// // check triangle parts 0 & 2;
+// // 0 part could not be checked in the loop; check it here
+// if ( isShared[2] && prevLinkShared &&
+// volume.IsLinked( nodes[ 0 ], nodes[ 1 ] ) &&
+// volume.IsLinked( nodes[ 1 ], nodes[ 3 ] ) )
+// return !free; // is not free
+// }
+// }
+// }
+// return free;
}
//=======================================================================
//purpose : Return index of a face formed by theFaceNodes
//=======================================================================
-int SMDS_VolumeTool::GetFaceIndex( const set<const SMDS_MeshNode*>& theFaceNodes )
+int SMDS_VolumeTool::GetFaceIndex( const set<const SMDS_MeshNode*>& theFaceNodes ) const
{
for ( int iFace = 0; iFace < myNbFaces; iFace++ ) {
const SMDS_MeshNode** nodes = GetFaceNodes( iFace );
//purpose :
//=======================================================================
-bool SMDS_VolumeTool::setFace( int faceIndex )
+bool SMDS_VolumeTool::setFace( int faceIndex ) const
{
if ( !myVolume )
return false;
return false;
}
- // check orientation
- bool isGoodOri = true;
- if (myExternalFaces)
- isGoodOri = IsFaceExternal( faceIndex );
-
// set face nodes
int iNode;
myFaceNbNodes = myPolyedre->NbFaceNodes(faceIndex + 1);
myFaceNodes = new const SMDS_MeshNode* [myFaceNbNodes + 1];
- if (isGoodOri) {
- for ( iNode = 0; iNode < myFaceNbNodes; iNode++ )
- myFaceNodes[ iNode ] = myPolyedre->GetFaceNode(faceIndex + 1, iNode + 1);
- } else {
- for ( iNode = 0; iNode < myFaceNbNodes; iNode++ )
- myFaceNodes[ iNode ] = myPolyedre->GetFaceNode(faceIndex + 1, myFaceNbNodes - iNode);
- }
+ for ( iNode = 0; iNode < myFaceNbNodes; iNode++ )
+ myFaceNodes[ iNode ] = myPolyedre->GetFaceNode(faceIndex + 1, iNode + 1);
myFaceNodes[ myFaceNbNodes ] = myFaceNodes[ 0 ]; // last = first
+ // check orientation
+ if (myExternalFaces)
+ {
+ myCurFace = faceIndex; // avoid infinite recursion in IsFaceExternal()
+ myExternalFaces = false; // force normal computation by IsFaceExternal()
+ if ( !IsFaceExternal( faceIndex ))
+ for ( int i = 0, j = myFaceNbNodes; i < j; ++i, --j )
+ std::swap( myFaceNodes[i], myFaceNodes[j] );
+ myExternalFaces = true;
+ }
}
else {
// choose face node indices