ForkMaster/smesh
1
0
Fork 0
mirror of https://git.salome-platform.org/gitpub/modules/smesh.git synced 2025-02-09 23:53:08 +05:00
Code Issues Packages Projects Releases Wiki Activity

IPAL52781: Orientation of a wall face created by revolution is wrong

Browse Source 
IPAL52850: ConvertToQuadratic creates invalid triangles
Regression of imps_09/K0
This commit is contained in:
eap 2015-08-28 20:03:25 +03:00
parent 0db2de9312
commit 1dd2f82c6d
12 changed files with 531 additions and 279 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

BIN
doc/salome/gui/SMESH/images/revolutionsn1.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 2.7 KiB

After

Width:  |  Height:  |  Size: 1.1 KiB

BIN
doc/salome/gui/SMESH/images/revolutionsn2.png
View File

Binary file not shown.
Side by Side Swipe Overlay

Before

Width:  |  Height:  |  Size: 2.6 KiB

After

Width:  |  Height:  |  Size: 1.9 KiB

5
doc/salome/gui/SMESH/input/extrusion.doc
View File

@ -48,6 +48,11 @@ The following dialog will appear:
<li>In this dialog: <li>In this dialog:
<ul> <ul>
<li>Use \a Selection button to specify what you are going to
select at a given moment, \b Nodes, \b Edges or \b Faces.
\image html image120.png
<center><em>"Selection" button</em></center>
</li>
<li>Specify \b Nodes, \b Edges and \b Faces, which will be extruded, by one <li>Specify \b Nodes, \b Edges and \b Faces, which will be extruded, by one
of following means: of following means:
<ul> <ul>

12
doc/salome/gui/SMESH/input/revolution.doc
View File

@ -68,16 +68,16 @@ The following dialog will appear:
<ul> <ul>
<li> <b>Angle by Step</b> - the elements are revolved by the <li> <b>Angle by Step</b> - the elements are revolved by the
specified angle at each step (i.e. for Angle=30 and Number of specified angle at each step (i.e. for Angle=30 and Number of
Steps=2, the elements will be extruded by 30 degrees twice for a Steps=3, the elements will be extruded by 30 degrees twice for a
total of 30*2=60) total of 30*3=90)
\image html revolutionsn2.png "Example of Revolution with Angle by Step" \image html revolutionsn2.png "Example of Revolution with Angle by Step. Angle=30 and Number of Steps=3"
</li> </li>
<li> <b>Total Angle</b> - the elements are revolved by the <li> <b>Total Angle</b> - the elements are revolved by the
specified angle only once and the number of steps defines the specified angle only once and the number of steps defines the
number of iterations (i.e. for Angle=30 and Number of Steps=2, number of iterations (i.e. for Angle=30 and Number of Steps=3,
the elements will be revolved by 30/2=15 degrees twice for a the elements will be revolved by 30/3=10 degrees twice for a
total of 30). total of 30).
\image html revolutionsn1.png "Example of Revolution with Total Angle" \image html revolutionsn1.png "Example of Revolution with Total Angle. Angle=30 and Number of Steps=3"
</li> </li>
</ul> </ul>
</li> </li>

611
src/SMDS/SMDS_VolumeTool.cxx
View File

@ -40,6 +40,8 @@
#include <map> #include <map>
#include <limits> #include <limits>
#include <cmath> #include <cmath>
#include <numeric>
#include <algorithm>
using namespace std; using namespace std;
@ -410,14 +412,33 @@ inline double XYZ::SquareMagnitude() {
} // namespace } // namespace
//================================================================================
/*!
* \brief Saver/restorer of a SMDS_VolumeTool::myCurFace
*/
//================================================================================
struct SMDS_VolumeTool::SaveFacet
{
SMDS_VolumeTool::Facet mySaved;
SMDS_VolumeTool::Facet& myToRestore;
SaveFacet( SMDS_VolumeTool::Facet& facet ): myToRestore( facet )
{
mySaved = facet;
}
~SaveFacet()
{
if ( myToRestore.myIndex != mySaved.myIndex )
myToRestore = mySaved;
}
};
//======================================================================= //=======================================================================
//function : SMDS_VolumeTool //function : SMDS_VolumeTool
//purpose : //purpose :
//======================================================================= //=======================================================================
SMDS_VolumeTool::SMDS_VolumeTool () SMDS_VolumeTool::SMDS_VolumeTool ()
: myVolumeNodes( NULL ),
myFaceNodes( NULL )
{ {
Set( 0 ); Set( 0 );
} }
@ -429,8 +450,6 @@ SMDS_VolumeTool::SMDS_VolumeTool ()
SMDS_VolumeTool::SMDS_VolumeTool (const SMDS_MeshElement* theVolume, SMDS_VolumeTool::SMDS_VolumeTool (const SMDS_MeshElement* theVolume,
const bool ignoreCentralNodes) const bool ignoreCentralNodes)
: myVolumeNodes( NULL ),
myFaceNodes( NULL )
{ {
Set( theVolume, ignoreCentralNodes ); Set( theVolume, ignoreCentralNodes );
} }
@ -442,11 +461,7 @@ SMDS_VolumeTool::SMDS_VolumeTool (const SMDS_MeshElement* theVolume,
SMDS_VolumeTool::~SMDS_VolumeTool() SMDS_VolumeTool::~SMDS_VolumeTool()
{ {
if ( myVolumeNodes != NULL ) delete [] myVolumeNodes; myCurFace.myNodeIndices = NULL;
if ( myFaceNodes != NULL ) delete [] myFaceNodes;
myFaceNodeIndices = NULL;
myVolumeNodes = myFaceNodes = NULL;
} }
//======================================================================= //=======================================================================
@ -464,42 +479,37 @@ bool SMDS_VolumeTool::Set (const SMDS_MeshElement* theVolume,
myVolForward = true; myVolForward = true;
myNbFaces = 0; myNbFaces = 0;
myVolumeNbNodes = 0; myVolumeNodes.clear();
if (myVolumeNodes != NULL) {
delete [] myVolumeNodes;
myVolumeNodes = NULL;
}
myPolyIndices.clear(); myPolyIndices.clear();
myPolyQuantities.clear();
myPolyFacetOri.clear();
myFwdLinks.clear();
myExternalFaces = false; myExternalFaces = false;
myAllFacesNodeIndices_F = 0; myAllFacesNodeIndices_F = 0;
//myAllFacesNodeIndices_FE = 0;
myAllFacesNodeIndices_RE = 0; myAllFacesNodeIndices_RE = 0;
myAllFacesNbNodes = 0; myAllFacesNbNodes = 0;
myCurFace = -1; myCurFace.myIndex = -1;
myFaceNbNodes = 0; myCurFace.myNodeIndices = NULL;
myFaceNodeIndices = NULL; myCurFace.myNodes.clear();
if (myFaceNodes != NULL) {
delete [] myFaceNodes;
myFaceNodes = NULL;
}
// set volume data // set volume data
if ( !theVolume || theVolume->GetType() != SMDSAbs_Volume ) if ( !theVolume || theVolume->GetType() != SMDSAbs_Volume )
return false; return false;
myVolume = theVolume; myVolume = theVolume;
if (myVolume->IsPoly())
myPolyedre = dynamic_cast<const SMDS_VtkVolume*>( myVolume );
myNbFaces = theVolume->NbFaces(); myNbFaces = theVolume->NbFaces();
myVolumeNbNodes = theVolume->NbNodes(); if ( myVolume->IsPoly() )
{
myPolyedre = dynamic_cast<const SMDS_VtkVolume*>( myVolume );
myPolyFacetOri.resize( myNbFaces, 0 );
}
// set nodes // set nodes
int iNode = 0; int iNode = 0;
myVolumeNodes = new const SMDS_MeshNode* [myVolumeNbNodes]; myVolumeNodes.resize( myVolume->NbNodes() );
SMDS_ElemIteratorPtr nodeIt = myVolume->nodesIterator(); SMDS_ElemIteratorPtr nodeIt = myVolume->nodesIterator();
while ( nodeIt->more() ) while ( nodeIt->more() )
myVolumeNodes[ iNode++ ] = static_cast<const SMDS_MeshNode*>( nodeIt->next() ); myVolumeNodes[ iNode++ ] = static_cast<const SMDS_MeshNode*>( nodeIt->next() );
@ -512,7 +522,8 @@ bool SMDS_VolumeTool::Set (const SMDS_MeshElement* theVolume,
{ {
// define volume orientation // define volume orientation
XYZ botNormal; XYZ botNormal;
GetFaceNormal( 0, botNormal.x, botNormal.y, botNormal.z ); if ( GetFaceNormal( 0, botNormal.x, botNormal.y, botNormal.z ))
{
const SMDS_MeshNode* botNode = myVolumeNodes[ 0 ]; const SMDS_MeshNode* botNode = myVolumeNodes[ 0 ];
int topNodeIndex = myVolume->NbCornerNodes() - 1; int topNodeIndex = myVolume->NbCornerNodes() - 1;
while ( !IsLinked( 0, topNodeIndex, /*ignoreMediumNodes=*/true )) --topNodeIndex; while ( !IsLinked( 0, topNodeIndex, /*ignoreMediumNodes=*/true )) --topNodeIndex;
@ -521,9 +532,9 @@ bool SMDS_VolumeTool::Set (const SMDS_MeshElement* theVolume,
topNode->Y() - botNode->Y(), topNode->Y() - botNode->Y(),
topNode->Z() - botNode->Z() ); topNode->Z() - botNode->Z() );
myVolForward = ( botNormal.Dot( upDir ) < 0 ); myVolForward = ( botNormal.Dot( upDir ) < 0 );
}
if ( !myVolForward ) if ( !myVolForward )
myCurFace = -1; // previous setFace(0) didn't take myVolForward into account myCurFace.myIndex = -1; // previous setFace(0) didn't take myVolForward into account
} }
return true; return true;
} }
@ -549,10 +560,10 @@ void SMDS_VolumeTool::Inverse ()
} }
myVolForward = !myVolForward; myVolForward = !myVolForward;
myCurFace = -1; myCurFace.myIndex = -1;
// inverse top and bottom faces // inverse top and bottom faces
switch ( myVolumeNbNodes ) { switch ( myVolumeNodes.size() ) {
case 4: case 4:
SWAP_NODES( myVolumeNodes, 1, 2 ); SWAP_NODES( myVolumeNodes, 1, 2 );
break; break;
@ -626,7 +637,7 @@ SMDS_VolumeTool::VolumeType SMDS_VolumeTool::GetVolumeType() const
if ( myPolyedre ) if ( myPolyedre )
return POLYHEDA; return POLYHEDA;
switch( myVolumeNbNodes ) { switch( myVolumeNodes.size() ) {
case 4: return TETRA; case 4: return TETRA;
case 5: return PYRAM; case 5: return PYRAM;
case 6: return PENTA; case 6: return PENTA;
@ -653,28 +664,18 @@ static double getTetraVolume(const SMDS_MeshNode* n1,
const SMDS_MeshNode* n3, const SMDS_MeshNode* n3,
const SMDS_MeshNode* n4) const SMDS_MeshNode* n4)
{ {
double X1 = n1->X(); double p1[3], p2[3], p3[3], p4[3];
double Y1 = n1->Y(); n1->GetXYZ( p1 );
double Z1 = n1->Z(); n2->GetXYZ( p2 );
n3->GetXYZ( p3 );
n4->GetXYZ( p4 );
double X2 = n2->X(); double Q1 = -(p1[ 0 ]-p2[ 0 ])*(p3[ 1 ]*p4[ 2 ]-p4[ 1 ]*p3[ 2 ]);
double Y2 = n2->Y(); double Q2 = (p1[ 0 ]-p3[ 0 ])*(p2[ 1 ]*p4[ 2 ]-p4[ 1 ]*p2[ 2 ]);
double Z2 = n2->Z(); double R1 = -(p1[ 0 ]-p4[ 0 ])*(p2[ 1 ]*p3[ 2 ]-p3[ 1 ]*p2[ 2 ]);
double R2 = -(p2[ 0 ]-p3[ 0 ])*(p1[ 1 ]*p4[ 2 ]-p4[ 1 ]*p1[ 2 ]);
double X3 = n3->X(); double S1 = (p2[ 0 ]-p4[ 0 ])*(p1[ 1 ]*p3[ 2 ]-p3[ 1 ]*p1[ 2 ]);
double Y3 = n3->Y(); double S2 = -(p3[ 0 ]-p4[ 0 ])*(p1[ 1 ]*p2[ 2 ]-p2[ 1 ]*p1[ 2 ]);
double Z3 = n3->Z();
double X4 = n4->X();
double Y4 = n4->Y();
double Z4 = n4->Z();
double Q1 = -(X1-X2)*(Y3*Z4-Y4*Z3);
double Q2 = (X1-X3)*(Y2*Z4-Y4*Z2);
double R1 = -(X1-X4)*(Y2*Z3-Y3*Z2);
double R2 = -(X2-X3)*(Y1*Z4-Y4*Z1);
double S1 = (X2-X4)*(Y1*Z3-Y3*Z1);
double S2 = -(X3-X4)*(Y1*Z2-Y2*Z1);
return (Q1+Q2+R1+R2+S1+S2)/6.0; return (Q1+Q2+R1+R2+S1+S2)/6.0;
} }
@ -697,23 +698,21 @@ double SMDS_VolumeTool::GetSize() const
// split a polyhedron into tetrahedrons // split a polyhedron into tetrahedrons
int saveCurFace = myCurFace; SaveFacet savedFacet( myCurFace );
SMDS_VolumeTool* me = const_cast< SMDS_VolumeTool* > ( this ); SMDS_VolumeTool* me = const_cast< SMDS_VolumeTool* > ( this );
for ( int f = 0; f < NbFaces(); ++f ) for ( int f = 0; f < NbFaces(); ++f )
{ {
me->setFace( f ); me->setFace( f );
XYZ area (0,0,0), p1( myFaceNodes[0] ); XYZ area (0,0,0), p1( myCurFace.myNodes[0] );
for ( int n = 0; n < myFaceNbNodes; ++n ) for ( int n = 0; n < myCurFace.myNbNodes; ++n )
{ {
XYZ p2( myFaceNodes[ n+1 ]); XYZ p2( myCurFace.myNodes[ n+1 ]);
area = area + p1.Crossed( p2 ); area = area + p1.Crossed( p2 );
p1 = p2; p1 = p2;
} }
V += p1.Dot( area ); V += p1.Dot( area );
} }
V /= 6; V /= 6;
if ( saveCurFace > -1 && saveCurFace != myCurFace )
me->setFace( myCurFace );
} }
else else
{ {
@ -844,14 +843,14 @@ bool SMDS_VolumeTool::GetBaryCenter(double & X, double & Y, double & Z) const
if ( !myVolume ) if ( !myVolume )
return false; return false;
for ( int i = 0; i < myVolumeNbNodes; i++ ) { for ( int i = 0; i < myVolumeNodes.size(); i++ ) {
X += myVolumeNodes[ i ]->X(); X += myVolumeNodes[ i ]->X();
Y += myVolumeNodes[ i ]->Y(); Y += myVolumeNodes[ i ]->Y();
Z += myVolumeNodes[ i ]->Z(); Z += myVolumeNodes[ i ]->Z();
} }
X /= myVolumeNbNodes; X /= myVolumeNodes.size();
Y /= myVolumeNbNodes; Y /= myVolumeNodes.size();
Z /= myVolumeNbNodes; Z /= myVolumeNodes.size();
return true; return true;
} }
@ -875,7 +874,7 @@ bool SMDS_VolumeTool::IsOut(double X, double Y, double Z, double tol) const
if ( !IsFaceExternal( iF )) if ( !IsFaceExternal( iF ))
faceNormal = XYZ() - faceNormal; // reverse faceNormal = XYZ() - faceNormal; // reverse
XYZ face2p( p - XYZ( myFaceNodes[0] )); XYZ face2p( p - XYZ( myCurFace.myNodes[0] ));
if ( face2p.Dot( faceNormal ) > tol ) if ( face2p.Dot( faceNormal ) > tol )
return true; return true;
} }
@ -890,7 +889,7 @@ bool SMDS_VolumeTool::IsOut(double X, double Y, double Z, double tol) const
void SMDS_VolumeTool::SetExternalNormal () void SMDS_VolumeTool::SetExternalNormal ()
{ {
myExternalFaces = true; myExternalFaces = true;
myCurFace = -1; myCurFace.myIndex = -1;
} }
//======================================================================= //=======================================================================
@ -902,7 +901,7 @@ int SMDS_VolumeTool::NbFaceNodes( int faceIndex ) const
{ {
if ( !setFace( faceIndex )) if ( !setFace( faceIndex ))
return 0; return 0;
return myFaceNbNodes; return myCurFace.myNbNodes;
} }
//======================================================================= //=======================================================================
@ -917,7 +916,7 @@ const SMDS_MeshNode** SMDS_VolumeTool::GetFaceNodes( int faceIndex ) const
{ {
if ( !setFace( faceIndex )) if ( !setFace( faceIndex ))
return 0; return 0;
return myFaceNodes; return &myCurFace.myNodes[0];
} }
//======================================================================= //=======================================================================
@ -933,15 +932,7 @@ const int* SMDS_VolumeTool::GetFaceNodesIndices( int faceIndex ) const
if ( !setFace( faceIndex )) if ( !setFace( faceIndex ))
return 0; return 0;
if (myPolyedre) return myCurFace.myNodeIndices;
{
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 )
me->myFaceNodeIndices[i] = myVolume->GetNodeIndex( myFaceNodes[i] );
}
return myFaceNodeIndices;
} }
//======================================================================= //=======================================================================
@ -956,11 +947,44 @@ bool SMDS_VolumeTool::GetFaceNodes (int faceIndex,
return false; return false;
theFaceNodes.clear(); theFaceNodes.clear();
theFaceNodes.insert( myFaceNodes, myFaceNodes + myFaceNbNodes ); theFaceNodes.insert( myCurFace.myNodes.begin(), myCurFace.myNodes.end() );
return true; return true;
} }
namespace
{
struct NLink : public std::pair<int,int>
{
int myOri;
NLink(const SMDS_MeshNode* n1=0, const SMDS_MeshNode* n2=0, int ori=1 )
{
if ( n1 )
{
if (( myOri = ( n1->GetID() < n2->GetID() )))
{
first = n1->GetID();
second = n2->GetID();
}
else
{
myOri = -1;
first = n2->GetID();
second = n1->GetID();
}
myOri *= ori;
}
else
{
myOri = first = second = 0;
}
}
//int Node1() const { return myOri == -1 ? second : first; }
//bool IsSameOri( const std::pair<int,int>& link ) const { return link.first == Node1(); }
};
}
//======================================================================= //=======================================================================
//function : IsFaceExternal //function : IsFaceExternal
//purpose : Check normal orientation of a given face //purpose : Check normal orientation of a given face
@ -971,39 +995,179 @@ bool SMDS_VolumeTool::IsFaceExternal( int faceIndex ) const
if ( myExternalFaces || !myVolume ) if ( myExternalFaces || !myVolume )
return true; return true;
if (myVolume->IsPoly()) { if ( !myPolyedre ) // all classical volumes have external facet normals
XYZ aNormal, baryCenter, p0 (myPolyedre->GetFaceNode(faceIndex + 1, 1));
GetFaceNormal(faceIndex, aNormal.x, aNormal.y, aNormal.z);
GetBaryCenter(baryCenter.x, baryCenter.y, baryCenter.z);
XYZ insideVec (baryCenter - p0);
if (insideVec.Dot(aNormal) > 0)
return false;
return true; return true;
SMDS_VolumeTool* me = const_cast< SMDS_VolumeTool* >( this );
if ( myPolyFacetOri[ faceIndex ])
return myPolyFacetOri[ faceIndex ] > 0;
int ori = 0; // -1-in, +1-out, 0-undef
double minProj, maxProj;
if ( projectNodesToNormal( faceIndex, minProj, maxProj ))
{
// all nodes are on the same side of the facet
ori = ( minProj < 0 ? +1 : -1 );
me->myPolyFacetOri[ faceIndex ] = ori;
if ( !me->myFwdLinks.empty() ) // concave polyhedron; collect oriented links
for ( int i = 0; i < myCurFace.myNbNodes; ++i )
{
NLink link( myCurFace.myNodes[i], myCurFace.myNodes[i+1], ori );
me->myFwdLinks.insert( make_pair( link, link.myOri ));
}
return ori > 0;
} }
// switch ( myVolumeNbNodes ) { SaveFacet savedFacet( myCurFace );
// case 4:
// case 5: // concave polyhedron
// case 10:
// case 13: if ( me->myFwdLinks.empty() ) // get links of the least ambiguously oriented facet
// // only the bottom of a reversed tetrahedron can be internal {
// return ( myVolForward || faceIndex != 0 ); for ( size_t i = 0; i < myPolyFacetOri.size() && !ori; ++i )
// case 6: ori = me->myPolyFacetOri[ i ];
// case 15:
// case 12: if ( !ori ) // none facet is oriented yet
// // in a forward prism, the top is internal, in a reversed one - bottom {
// return ( myVolForward ? faceIndex != 1 : faceIndex != 0 ); // find the least ambiguously oriented facet
// case 8: int faceMostConvex = -1;
// case 20: std::map< double, int > convexity2face;
// case 27: { for ( size_t iF = 0; iF < myPolyFacetOri.size() && faceMostConvex < 0; ++iF )
// // in a forward hexahedron, even face normal is external, odd - internal {
// bool odd = faceIndex % 2; if ( projectNodesToNormal( iF, minProj, maxProj ))
// return ( myVolForward ? !odd : odd ); {
// all nodes are on the same side of the facet
me->myPolyFacetOri[ iF ] = ( minProj < 0 ? +1 : -1 );
faceMostConvex = iF;
}
else
{
ori = ( -minProj < maxProj ? -1 : +1 );
double convexity = std::min( -minProj, maxProj ) / std::max( -minProj, maxProj );
convexity2face.insert( make_pair( convexity, iF * ori ));
}
}
if ( faceMostConvex < 0 ) // none facet has nodes on the same side
{
// use the least ambiguous facet
faceMostConvex = convexity2face.begin()->second;
ori = ( faceMostConvex < 0 ? -1 : +1 );
faceMostConvex = std::abs( faceMostConvex );
me->myPolyFacetOri[ faceMostConvex ] = ori;
}
}
// collect links of the oriented facets in me->myFwdLinks
for ( size_t iF = 0; iF < myPolyFacetOri.size(); ++iF )
{
ori = me->myPolyFacetOri[ iF ];
if ( !ori ) continue;
setFace( iF );
for ( int i = 0; i < myCurFace.myNbNodes; ++i )
{
NLink link( myCurFace.myNodes[i], myCurFace.myNodes[i+1], ori );
me->myFwdLinks.insert( make_pair( link, link.myOri ));
}
}
}
// compare orientation of links of the facet with myFwdLinks
ori = 0;
setFace( faceIndex );
vector< NLink > links( myCurFace.myNbNodes ), links2;
for ( int i = 0; i < myCurFace.myNbNodes && !ori; ++i )
{
NLink link( myCurFace.myNodes[i], myCurFace.myNodes[i+1] );
std::map<Link, int>::iterator l2o = me->myFwdLinks.find( link );
if ( l2o != me->myFwdLinks.end() )
ori = link.myOri * l2o->second * -1;
links[ i ] = link;
}
while ( !ori ) // the facet has no common links with already oriented facets
{
// orient and collect links of other non-oriented facets
for ( size_t iF = 0; iF < myPolyFacetOri.size(); ++iF )
{
if ( me->myPolyFacetOri[ iF ] ) continue; // already oriented
setFace( iF );
links2.clear();
ori = 0;
for ( int i = 0; i < myCurFace.myNbNodes && !ori; ++i )
{
NLink link( myCurFace.myNodes[i], myCurFace.myNodes[i+1] );
std::map<Link, int>::iterator l2o = me->myFwdLinks.find( link );
if ( l2o != me->myFwdLinks.end() )
ori = link.myOri * l2o->second * -1;
links2.push_back( link );
}
if ( ori ) // one more facet oriented
{
me->myPolyFacetOri[ iF ] = ori;
for ( size_t i = 0; i < links2.size(); ++i )
me->myFwdLinks.insert( make_pair( links2[i], links2[i].myOri * ori ));
break;
}
}
if ( !ori )
return false; // error in algorithm: infinite loop
// try to orient the facet again
ori = 0;
for ( size_t i = 0; i < links.size() && !ori; ++i )
{
std::map<Link, int>::iterator l2o = me->myFwdLinks.find( links[i] );
if ( l2o != me->myFwdLinks.end() )
ori = links[i].myOri * l2o->second * -1;
}
me->myPolyFacetOri[ faceIndex ] = ori;
}
return ori > 0;
}
//=======================================================================
//function : projectNodesToNormal
//purpose : compute min and max projections of all nodes to normal of a facet.
//=======================================================================
bool SMDS_VolumeTool::projectNodesToNormal( int faceIndex,
double& minProj,
double& maxProj ) const
{
minProj = std::numeric_limits<double>::max();
maxProj = std::numeric_limits<double>::min();
XYZ normal;
if ( !GetFaceNormal( faceIndex, normal.x, normal.y, normal.z ))
return false;
XYZ p0 ( myCurFace.myNodes[0] );
for ( size_t i = 0; i < myVolumeNodes.size(); ++i )
{
if ( std::find( myCurFace.myNodes.begin() + 1,
myCurFace.myNodes.end(),
myVolumeNodes[ i ] ) != myCurFace.myNodes.end() )
continue; // node of the faceIndex-th facet
double proj = normal.Dot( XYZ( myVolumeNodes[ i ]) - p0 );
if ( proj < minProj ) minProj = proj;
if ( proj > maxProj ) maxProj = proj;
}
const double tol = 1e-7;
minProj += tol;
maxProj -= tol;
bool diffSize = ( minProj * maxProj < 0 );
// if ( diffSize )
// {
// minProj = -minProj;
// } // }
// default:; // else if ( minProj < 0 )
// {
// minProj = -minProj;
// maxProj = -maxProj;
// } // }
// return false;
return true; return !diffSize; // ? 0 : (minProj >= 0);
} }
//======================================================================= //=======================================================================
@ -1016,16 +1180,16 @@ bool SMDS_VolumeTool::GetFaceNormal (int faceIndex, double & X, double & Y, doub
if ( !setFace( faceIndex )) if ( !setFace( faceIndex ))
return false; return false;
const int iQuad = ( myFaceNbNodes > 6 && !myPolyedre ) ? 2 : 1; const int iQuad = ( !myPolyedre && myCurFace.myNbNodes > 6 ) ? 2 : 1;
XYZ p1 ( myFaceNodes[0*iQuad] ); XYZ p1 ( myCurFace.myNodes[0*iQuad] );
XYZ p2 ( myFaceNodes[1*iQuad] ); XYZ p2 ( myCurFace.myNodes[1*iQuad] );
XYZ p3 ( myFaceNodes[2*iQuad] ); XYZ p3 ( myCurFace.myNodes[2*iQuad] );
XYZ aVec12( p2 - p1 ); XYZ aVec12( p2 - p1 );
XYZ aVec13( p3 - p1 ); XYZ aVec13( p3 - p1 );
XYZ cross = aVec12.Crossed( aVec13 ); XYZ cross = aVec12.Crossed( aVec13 );
if ( myFaceNbNodes >3*iQuad ) { if ( myCurFace.myNbNodes >3*iQuad ) {
XYZ p4 ( myFaceNodes[3*iQuad] ); XYZ p4 ( myCurFace.myNodes[3*iQuad] );
XYZ aVec14( p4 - p1 ); XYZ aVec14( p4 - p1 );
XYZ cross2 = aVec13.Crossed( aVec14 ); XYZ cross2 = aVec13.Crossed( aVec14 );
cross = cross + cross2; cross = cross + cross2;
@ -1054,11 +1218,11 @@ bool SMDS_VolumeTool::GetFaceBaryCenter (int faceIndex, double & X, double & Y,
return false; return false;
X = Y = Z = 0.0; X = Y = Z = 0.0;
for ( int i = 0; i < myFaceNbNodes; ++i ) for ( int i = 0; i < myCurFace.myNbNodes; ++i )
{ {
X += myFaceNodes[i]->X() / myFaceNbNodes; X += myCurFace.myNodes[i]->X() / myCurFace.myNbNodes;
Y += myFaceNodes[i]->Y() / myFaceNbNodes; Y += myCurFace.myNodes[i]->Y() / myCurFace.myNbNodes;
Z += myFaceNodes[i]->Z() / myFaceNbNodes; Z += myCurFace.myNodes[i]->Z() / myCurFace.myNbNodes;
} }
return true; return true;
} }
@ -1070,27 +1234,36 @@ bool SMDS_VolumeTool::GetFaceBaryCenter (int faceIndex, double & X, double & Y,
double SMDS_VolumeTool::GetFaceArea( int faceIndex ) const double SMDS_VolumeTool::GetFaceArea( int faceIndex ) const
{ {
if (myVolume->IsPoly()) { double area = 0;
MESSAGE("Warning: attempt to obtain area of a face of polyhedral volume");
return 0;
}
if ( !setFace( faceIndex )) if ( !setFace( faceIndex ))
return 0; return area;
XYZ p1 ( myFaceNodes[0] ); XYZ p1 ( myCurFace.myNodes[0] );
XYZ p2 ( myFaceNodes[1] ); XYZ p2 ( myCurFace.myNodes[1] );
XYZ p3 ( myFaceNodes[2] ); XYZ p3 ( myCurFace.myNodes[2] );
XYZ aVec12( p2 - p1 ); XYZ aVec12( p2 - p1 );
XYZ aVec13( p3 - p1 ); XYZ aVec13( p3 - p1 );
double area = aVec12.Crossed( aVec13 ).Magnitude() * 0.5; area += aVec12.Crossed( aVec13 ).Magnitude();
if ( myFaceNbNodes == 4 ) { if (myVolume->IsPoly())
XYZ p4 ( myFaceNodes[3] ); {
XYZ aVec14( p4 - p1 ); for ( int i = 3; i < myCurFace.myNbNodes; ++i )
area += aVec14.Crossed( aVec13 ).Magnitude() * 0.5; {
XYZ pI ( myCurFace.myNodes[i] );
XYZ aVecI( pI - p1 );
area += aVec13.Crossed( aVecI ).Magnitude();
aVec13 = aVecI;
} }
return area; }
else
{
if ( myCurFace.myNbNodes == 4 ) {
XYZ p4 ( myCurFace.myNodes[3] );
XYZ aVec14( p4 - p1 );
area += aVec14.Crossed( aVec13 ).Magnitude();
}
}
return area / 2;
} }
//================================================================================ //================================================================================
@ -1101,7 +1274,7 @@ double SMDS_VolumeTool::GetFaceArea( int faceIndex ) const
int SMDS_VolumeTool::GetCenterNodeIndex( int faceIndex ) const int SMDS_VolumeTool::GetCenterNodeIndex( int faceIndex ) const
{ {
if ( myAllFacesNbNodes && myVolumeNbNodes == 27 ) // classic element with 27 nodes if ( myAllFacesNbNodes && myVolumeNodes.size() == 27 ) // classic element with 27 nodes
{ {
switch ( faceIndex ) { switch ( faceIndex ) {
case 0: return 20; case 0: return 20;
@ -1129,7 +1302,7 @@ int SMDS_VolumeTool::GetOppFaceIndex( int faceIndex ) const
const int nbHoriFaces = 2; const int nbHoriFaces = 2;
if ( faceIndex >= 0 && faceIndex < NbFaces() ) { if ( faceIndex >= 0 && faceIndex < NbFaces() ) {
switch ( myVolumeNbNodes ) { switch ( myVolumeNodes.size() ) {
case 6: case 6:
case 15: case 15:
if ( faceIndex == 0 || faceIndex == 1 ) if ( faceIndex == 0 || faceIndex == 1 )
@ -1182,31 +1355,42 @@ bool SMDS_VolumeTool::IsLinked (const SMDS_MeshNode* theNode1,
MESSAGE("Warning: bad volumic element"); MESSAGE("Warning: bad volumic element");
return false; return false;
} }
bool isLinked = false; if ( !myAllFacesNbNodes ) {
int iface; SMDS_VolumeTool* me = const_cast< SMDS_VolumeTool* >( this );
for (iface = 1; iface <= myNbFaces && !isLinked; iface++) { me->myPolyQuantities = myPolyedre->GetQuantities();
int inode, nbFaceNodes = myPolyedre->NbFaceNodes(iface); myAllFacesNbNodes = &myPolyQuantities[0];
}
for (inode = 1; inode <= nbFaceNodes && !isLinked; inode++) { int from, to = 0, d1 = 1, d2 = 2;
const SMDS_MeshNode* curNode = myPolyedre->GetFaceNode(iface, inode); if ( myPolyedre->IsQuadratic() ) {
if ( theIgnoreMediumNodes ) {
if (curNode == theNode1 || curNode == theNode2) { d1 = 2; d2 = 0;
int inextnode = (inode == nbFaceNodes) ? 1 : inode + 1; }
const SMDS_MeshNode* nextNode = myPolyedre->GetFaceNode(iface, inextnode); } else {
d2 = 0;
if ((curNode == theNode1 && nextNode == theNode2) || }
(curNode == theNode2 && nextNode == theNode1)) { vector<const SMDS_MeshNode*>::const_iterator i;
isLinked = true; for (int iface = 0; iface < myNbFaces; iface++)
{
from = to;
to += myPolyQuantities[iface];
i = std::find( myVolumeNodes.begin() + from, myVolumeNodes.begin() + to, theNode1 );
if ( i != myVolumeNodes.end() )
{
if (( theNode2 == *( i-d1 ) ||
theNode2 == *( i+d1 )))
return true;
if (( d2 ) &&
(( theNode2 == *( i-d2 ) ||
theNode2 == *( i+d2 ))))
return true;
} }
} }
} return false;
}
return isLinked;
} }
// find nodes indices // find nodes indices
int i1 = -1, i2 = -1, nbFound = 0; int i1 = -1, i2 = -1, nbFound = 0;
for ( int i = 0; i < myVolumeNbNodes && nbFound < 2; i++ ) for ( int i = 0; i < myVolumeNodes.size() && nbFound < 2; i++ )
{ {
if ( myVolumeNodes[ i ] == theNode1 ) if ( myVolumeNodes[ i ] == theNode1 )
i1 = i, ++nbFound; i1 = i, ++nbFound;
@ -1234,7 +1418,7 @@ bool SMDS_VolumeTool::IsLinked (const int theNode1Index,
int minInd = min( theNode1Index, theNode2Index ); int minInd = min( theNode1Index, theNode2Index );
int maxInd = max( theNode1Index, theNode2Index ); int maxInd = max( theNode1Index, theNode2Index );
if ( minInd < 0 || maxInd > myVolumeNbNodes - 1 || maxInd == minInd ) if ( minInd < 0 || maxInd > myVolumeNodes.size() - 1 || maxInd == minInd )
return false; return false;
VolumeType type = GetVolumeType(); VolumeType type = GetVolumeType();
@ -1351,7 +1535,7 @@ bool SMDS_VolumeTool::IsLinked (const int theNode1Index,
int SMDS_VolumeTool::GetNodeIndex(const SMDS_MeshNode* theNode) const int SMDS_VolumeTool::GetNodeIndex(const SMDS_MeshNode* theNode) const
{ {
if ( myVolume ) { if ( myVolume ) {
for ( int i = 0; i < myVolumeNbNodes; i++ ) { for ( int i = 0; i < myVolumeNodes.size(); i++ ) {
if ( myVolumeNodes[ i ] == theNode ) if ( myVolumeNodes[ i ] == theNode )
return i; return i;
} }
@ -1370,6 +1554,14 @@ int SMDS_VolumeTool::GetNodeIndex(const SMDS_MeshNode* theNode) const
int SMDS_VolumeTool::GetAllExistingFaces(vector<const SMDS_MeshElement*> & faces) const int SMDS_VolumeTool::GetAllExistingFaces(vector<const SMDS_MeshElement*> & faces) const
{ {
faces.clear(); faces.clear();
SaveFacet savedFacet( myCurFace );
if ( IsPoly() )
for ( int iF = 0; iF < NbFaces(); ++iF ) {
if ( setFace( iF ))
if ( const SMDS_MeshElement* face = SMDS_Mesh::FindFace( myCurFace.myNodes ))
faces.push_back( face );
}
else
for ( int iF = 0; iF < NbFaces(); ++iF ) { for ( int iF = 0; iF < NbFaces(); ++iF ) {
const SMDS_MeshFace* face = 0; const SMDS_MeshFace* face = 0;
const SMDS_MeshNode** nodes = GetFaceNodes( iF ); const SMDS_MeshNode** nodes = GetFaceNodes( iF );
@ -1403,9 +1595,9 @@ int SMDS_VolumeTool::GetAllExistingFaces(vector<const SMDS_MeshElement*> & faces
int SMDS_VolumeTool::GetAllExistingEdges(vector<const SMDS_MeshElement*> & edges) const int SMDS_VolumeTool::GetAllExistingEdges(vector<const SMDS_MeshElement*> & edges) const
{ {
edges.clear(); edges.clear();
edges.reserve( myVolumeNbNodes * 2 ); edges.reserve( myVolumeNodes.size() * 2 );
for ( int i = 0; i < myVolumeNbNodes-1; ++i ) { for ( int i = 0; i < myVolumeNodes.size()-1; ++i ) {
for ( int j = i + 1; j < myVolumeNbNodes; ++j ) { for ( int j = i + 1; j < myVolumeNodes.size(); ++j ) {
if ( IsLinked( i, j )) { if ( IsLinked( i, j )) {
const SMDS_MeshElement* edge = const SMDS_MeshElement* edge =
SMDS_Mesh::FindEdge( myVolumeNodes[i], myVolumeNodes[j] ); SMDS_Mesh::FindEdge( myVolumeNodes[i], myVolumeNodes[j] );
@ -1428,30 +1620,20 @@ double SMDS_VolumeTool::MinLinearSize2() const
double minSize = 1e+100; double minSize = 1e+100;
int iQ = myVolume->IsQuadratic() ? 2 : 1; int iQ = myVolume->IsQuadratic() ? 2 : 1;
// store current face data SaveFacet savedFacet( myCurFace );
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 // it seems that compute distance twice is faster than organization of a sole computing
myCurFace = -1; myCurFace.myIndex = -1;
for ( int iF = 0; iF < myNbFaces; ++iF ) for ( int iF = 0; iF < myNbFaces; ++iF )
{ {
setFace( iF ); setFace( iF );
for ( int iN = 0; iN < myFaceNbNodes; iN += iQ ) for ( int iN = 0; iN < myCurFace.myNbNodes; iN += iQ )
{ {
XYZ n1( myFaceNodes[ iN ]); XYZ n1( myCurFace.myNodes[ iN ]);
XYZ n2( myFaceNodes[(iN + iQ) % myFaceNbNodes]); XYZ n2( myCurFace.myNodes[(iN + iQ) % myCurFace.myNbNodes]);
minSize = std::min( minSize, (n1 - n2).SquareMagnitude()); minSize = std::min( minSize, (n1 - n2).SquareMagnitude());
} }
} }
// restore current face data
myCurFace = curFace;
myFaceNbNodes = nbN;
myFaceNodeIndices = ind;
delete [] myFaceNodes; myFaceNodes = nodes;
return minSize; return minSize;
} }
@ -1467,30 +1649,20 @@ double SMDS_VolumeTool::MaxLinearSize2() const
double maxSize = -1e+100; double maxSize = -1e+100;
int iQ = myVolume->IsQuadratic() ? 2 : 1; int iQ = myVolume->IsQuadratic() ? 2 : 1;
// store current face data SaveFacet savedFacet( myCurFace );
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 // it seems that compute distance twice is faster than organization of a sole computing
myCurFace = -1; myCurFace.myIndex = -1;
for ( int iF = 0; iF < myNbFaces; ++iF ) for ( int iF = 0; iF < myNbFaces; ++iF )
{ {
setFace( iF ); setFace( iF );
for ( int iN = 0; iN < myFaceNbNodes; iN += iQ ) for ( int iN = 0; iN < myCurFace.myNbNodes; iN += iQ )
{ {
XYZ n1( myFaceNodes[ iN ]); XYZ n1( myCurFace.myNodes[ iN ]);
XYZ n2( myFaceNodes[(iN + iQ) % myFaceNbNodes]); XYZ n2( myCurFace.myNodes[(iN + iQ) % myCurFace.myNbNodes]);
maxSize = std::max( maxSize, (n1 - n2).SquareMagnitude()); maxSize = std::max( maxSize, (n1 - n2).SquareMagnitude());
} }
} }
// restore current face data
myCurFace = curFace;
myFaceNbNodes = nbN;
myFaceNodeIndices = ind;
delete [] myFaceNodes; myFaceNodes = nodes;
return maxSize; return maxSize;
} }
@ -1512,7 +1684,7 @@ bool SMDS_VolumeTool::IsFreeFace( int faceIndex, const SMDS_MeshElement** otherV
const SMDS_MeshNode** nodes = GetFaceNodes( faceIndex ); const SMDS_MeshNode** nodes = GetFaceNodes( faceIndex );
const int di = myVolume->IsQuadratic() ? 2 : 1; const int di = myVolume->IsQuadratic() ? 2 : 1;
const int nbN = ( myFaceNbNodes/di <= 4 && !IsPoly()) ? 3 : myFaceNbNodes/di; // nb nodes to check const int nbN = ( myCurFace.myNbNodes/di <= 4 && !IsPoly()) ? 3 : myCurFace.myNbNodes/di; // nb nodes to check
SMDS_ElemIteratorPtr eIt = nodes[0]->GetInverseElementIterator( SMDSAbs_Volume ); SMDS_ElemIteratorPtr eIt = nodes[0]->GetInverseElementIterator( SMDSAbs_Volume );
while ( eIt->more() ) while ( eIt->more() )
@ -1528,7 +1700,7 @@ bool SMDS_VolumeTool::IsFreeFace( int faceIndex, const SMDS_MeshElement** otherV
{ {
// if ( vol->IsPoly() || vol->NbFaces() > 6 ) // vol is polyhed or hex prism // if ( vol->IsPoly() || vol->NbFaces() > 6 ) // vol is polyhed or hex prism
// { // {
// int nb = myFaceNbNodes; // int nb = myCurFace.myNbNodes;
// if ( myVolume->GetEntityType() != vol->GetEntityType() ) // if ( myVolume->GetEntityType() != vol->GetEntityType() )
// nb -= ( GetCenterNodeIndex(0) > 0 ); // nb -= ( GetCenterNodeIndex(0) > 0 );
// set<const SMDS_MeshNode*> faceNodes( nodes, nodes + nb ); // set<const SMDS_MeshNode*> faceNodes( nodes, nodes + nb );
@ -1557,7 +1729,7 @@ bool SMDS_VolumeTool::IsFreeFaceAdv( int faceIndex, const SMDS_MeshElement** oth
return !isFree; return !isFree;
const SMDS_MeshNode** nodes = GetFaceNodes( faceIndex ); const SMDS_MeshNode** nodes = GetFaceNodes( faceIndex );
const int nbFaceNodes = myFaceNbNodes; const int nbFaceNodes = myCurFace.myNbNodes;
// evaluate nb of face nodes shared by other volumes // evaluate nb of face nodes shared by other volumes
int maxNbShared = -1; int maxNbShared = -1;
@ -1736,19 +1908,14 @@ bool SMDS_VolumeTool::setFace( int faceIndex ) const
if ( !myVolume ) if ( !myVolume )
return false; return false;
if ( myCurFace == faceIndex ) if ( myCurFace.myIndex == faceIndex )
return true; return true;
myCurFace = -1; myCurFace.myIndex = -1;
if ( faceIndex < 0 || faceIndex >= NbFaces() ) if ( faceIndex < 0 || faceIndex >= NbFaces() )
return false; return false;
if (myFaceNodes != NULL) {
delete [] myFaceNodes;
myFaceNodes = NULL;
}
if (myVolume->IsPoly()) if (myVolume->IsPoly())
{ {
if (!myPolyedre) { if (!myPolyedre) {
@ -1757,21 +1924,31 @@ bool SMDS_VolumeTool::setFace( int faceIndex ) const
} }
// set face nodes // set face nodes
int iNode; SMDS_VolumeTool* me = const_cast< SMDS_VolumeTool* >( this );
myFaceNbNodes = myPolyedre->NbFaceNodes(faceIndex + 1); if ( !myAllFacesNbNodes ) {
myFaceNodes = new const SMDS_MeshNode* [myFaceNbNodes + 1]; me->myPolyQuantities = myPolyedre->GetQuantities();
for ( iNode = 0; iNode < myFaceNbNodes; iNode++ ) myAllFacesNbNodes = &myPolyQuantities[0];
myFaceNodes[ iNode ] = myPolyedre->GetFaceNode(faceIndex + 1, iNode + 1); }
myFaceNodes[ myFaceNbNodes ] = myFaceNodes[ 0 ]; // last = first myCurFace.myNbNodes = myAllFacesNbNodes[ faceIndex ];
myCurFace.myNodes.resize( myCurFace.myNbNodes + 1 );
me->myPolyIndices.resize( myCurFace.myNbNodes + 1 );
myCurFace.myNodeIndices = & me->myPolyIndices[0];
int shift = std::accumulate( myAllFacesNbNodes, myAllFacesNbNodes+faceIndex, 0 );
for ( int iNode = 0; iNode < myCurFace.myNbNodes; iNode++ )
{
myCurFace.myNodes [ iNode ] = myVolumeNodes[ shift + iNode ];
myCurFace.myNodeIndices[ iNode ] = shift + iNode;
}
myCurFace.myNodes [ myCurFace.myNbNodes ] = myCurFace.myNodes[ 0 ]; // last = first
myCurFace.myNodeIndices[ myCurFace.myNbNodes ] = myCurFace.myNodeIndices[ 0 ];
// check orientation // check orientation
if (myExternalFaces) if (myExternalFaces)
{ {
myCurFace = faceIndex; // avoid infinite recursion in IsFaceExternal() myCurFace.myIndex = faceIndex; // avoid infinite recursion in IsFaceExternal()
myExternalFaces = false; // force normal computation by IsFaceExternal() myExternalFaces = false; // force normal computation by IsFaceExternal()
if ( !IsFaceExternal( faceIndex )) if ( !IsFaceExternal( faceIndex ))
for ( int i = 0, j = myFaceNbNodes; i < j; ++i, --j ) std::reverse( myCurFace.myNodes.begin(), myCurFace.myNodes.end() );
std::swap( myFaceNodes[i], myFaceNodes[j] );
myExternalFaces = true; myExternalFaces = true;
} }
} }
@ -1780,7 +1957,7 @@ bool SMDS_VolumeTool::setFace( int faceIndex ) const
if ( !myAllFacesNodeIndices_F ) if ( !myAllFacesNodeIndices_F )
{ {
// choose data for an element type // choose data for an element type
switch ( myVolumeNbNodes ) { switch ( myVolumeNodes.size() ) {
case 4: case 4:
myAllFacesNodeIndices_F = &Tetra_F [0][0]; myAllFacesNodeIndices_F = &Tetra_F [0][0];
//myAllFacesNodeIndices_FE = &Tetra_F [0][0]; //myAllFacesNodeIndices_FE = &Tetra_F [0][0];
@ -1837,7 +2014,7 @@ bool SMDS_VolumeTool::setFace( int faceIndex ) const
myAllFacesNodeIndices_RE = &QuadHexa_RE[0][0]; myAllFacesNodeIndices_RE = &QuadHexa_RE[0][0];
myAllFacesNbNodes = QuadHexa_nbN; myAllFacesNbNodes = QuadHexa_nbN;
myMaxFaceNbNodes = sizeof(QuadHexa_F[0])/sizeof(QuadHexa_F[0][0]); myMaxFaceNbNodes = sizeof(QuadHexa_F[0])/sizeof(QuadHexa_F[0][0]);
if ( !myIgnoreCentralNodes && myVolumeNbNodes == 27 ) if ( !myIgnoreCentralNodes && myVolumeNodes.size() == 27 )
{ {
myAllFacesNodeIndices_F = &TriQuadHexa_F [0][0]; myAllFacesNodeIndices_F = &TriQuadHexa_F [0][0];
//myAllFacesNodeIndices_FE = &TriQuadHexa_FE[0][0]; //myAllFacesNodeIndices_FE = &TriQuadHexa_FE[0][0];
@ -1857,21 +2034,21 @@ bool SMDS_VolumeTool::setFace( int faceIndex ) const
return false; return false;
} }
} }
myFaceNbNodes = myAllFacesNbNodes[ faceIndex ]; myCurFace.myNbNodes = myAllFacesNbNodes[ faceIndex ];
// if ( myExternalFaces ) // if ( myExternalFaces )
// myFaceNodeIndices = (int*)( myVolForward ? myAllFacesNodeIndices_FE + faceIndex*myMaxFaceNbNodes : myAllFacesNodeIndices_RE + faceIndex*myMaxFaceNbNodes ); // myCurFace.myNodeIndices = (int*)( myVolForward ? myAllFacesNodeIndices_FE + faceIndex*myMaxFaceNbNodes : myAllFacesNodeIndices_RE + faceIndex*myMaxFaceNbNodes );
// else // else
// myFaceNodeIndices = (int*)( myAllFacesNodeIndices_F + faceIndex*myMaxFaceNbNodes ); // myCurFace.myNodeIndices = (int*)( myAllFacesNodeIndices_F + faceIndex*myMaxFaceNbNodes );
myFaceNodeIndices = (int*)( myVolForward ? myAllFacesNodeIndices_F + faceIndex*myMaxFaceNbNodes : myAllFacesNodeIndices_RE + faceIndex*myMaxFaceNbNodes ); myCurFace.myNodeIndices = (int*)( myVolForward ? myAllFacesNodeIndices_F + faceIndex*myMaxFaceNbNodes : myAllFacesNodeIndices_RE + faceIndex*myMaxFaceNbNodes );
// set face nodes // set face nodes
myFaceNodes = new const SMDS_MeshNode* [myFaceNbNodes + 1]; myCurFace.myNodes.resize( myCurFace.myNbNodes + 1 );
for ( int iNode = 0; iNode < myFaceNbNodes; iNode++ ) for ( int iNode = 0; iNode < myCurFace.myNbNodes; iNode++ )
myFaceNodes[ iNode ] = myVolumeNodes[ myFaceNodeIndices[ iNode ]]; myCurFace.myNodes[ iNode ] = myVolumeNodes[ myCurFace.myNodeIndices[ iNode ]];
myFaceNodes[ myFaceNbNodes ] = myFaceNodes[ 0 ]; myCurFace.myNodes[ myCurFace.myNbNodes ] = myCurFace.myNodes[ 0 ];
} }
myCurFace = faceIndex; myCurFace.myIndex = faceIndex;
return true; return true;
} }

30
src/SMDS/SMDS_VolumeTool.hxx
View File

@ -38,6 +38,7 @@ class SMDS_MeshVolume;
#include <vector> #include <vector>
#include <set> #include <set>
#include <map>
// ========================================================================= // =========================================================================
// //
@ -90,10 +91,10 @@ class SMDS_EXPORT SMDS_VolumeTool
// top and bottom faces are reversed. // top and bottom faces are reversed.
// Result of IsForward() and methods returning nodes change // Result of IsForward() and methods returning nodes change
const SMDS_MeshNode** GetNodes() { return myVolumeNodes; } const SMDS_MeshNode** GetNodes() { return &myVolumeNodes[0]; }
// Return array of volume nodes // Return array of volume nodes
int NbNodes() { return myVolumeNbNodes; } int NbNodes() { return myVolumeNodes.size(); }
// Return array of volume nodes // Return array of volume nodes
double GetSize() const; double GetSize() const;
@ -246,15 +247,21 @@ private:
bool setFace( int faceIndex ) const; bool setFace( int faceIndex ) const;
bool projectNodesToNormal( int faceIndex, double& minProj, double& maxProj ) const;
const SMDS_MeshElement* myVolume; const SMDS_MeshElement* myVolume;
const SMDS_VtkVolume* myPolyedre; const SMDS_VtkVolume* myPolyedre;
bool myIgnoreCentralNodes; bool myIgnoreCentralNodes;
bool myVolForward; bool myVolForward;
int myNbFaces; int myNbFaces;
int myVolumeNbNodes; std::vector<const SMDS_MeshNode*> myVolumeNodes;
const SMDS_MeshNode** myVolumeNodes; std::vector< int > myPolyIndices; // of a myCurFace
std::vector< int > myPolyIndices; std::vector< int > myPolyQuantities;
std::vector< int > myPolyFacetOri; // -1-in, +1-out, 0-undef
typedef std::pair<int,int> Link;
std::map<Link, int> myFwdLinks; // used in IsFaceExternal() to find out myPolyFacetOri
mutable bool myExternalFaces; mutable bool myExternalFaces;
@ -263,10 +270,15 @@ private:
mutable const int* myAllFacesNbNodes; mutable const int* myAllFacesNbNodes;
mutable int myMaxFaceNbNodes; mutable int myMaxFaceNbNodes;
mutable int myCurFace; struct SaveFacet;
mutable int myFaceNbNodes; struct Facet
mutable int* myFaceNodeIndices; {
mutable const SMDS_MeshNode** myFaceNodes; int myIndex;
int myNbNodes;
int* myNodeIndices;
std::vector<const SMDS_MeshNode*> myNodes;
};
mutable Facet myCurFace;
}; };
#endif #endif

19
src/SMESH/SMESH_MeshEditor.cxx
View File

@ -12577,6 +12577,11 @@ int SMESH_MeshEditor::MakeBoundaryMesh(const TIDSortedElemSet& elements,
TConnectivity tgtNodes; TConnectivity tgtNodes;
ElemFeatures elemKind( missType ), elemToCopy; ElemFeatures elemKind( missType ), elemToCopy;
vector<const SMDS_MeshElement*> presentBndElems;
vector<TConnectivity> missingBndElems;
vector<int> freeFacets;
TConnectivity nodes, elemNodes;
SMDS_ElemIteratorPtr eIt; SMDS_ElemIteratorPtr eIt;
if (elements.empty()) eIt = aMesh->elementsIterator(elemType); if (elements.empty()) eIt = aMesh->elementsIterator(elemType);
else eIt = elemSetIterator( elements ); else eIt = elemSetIterator( elements );
@ -12590,18 +12595,24 @@ int SMESH_MeshEditor::MakeBoundaryMesh(const TIDSortedElemSet& elements,
// ------------------------------------------------------------------------------------ // ------------------------------------------------------------------------------------
// 1. For an elem, get present bnd elements and connectivities of missing bnd elements // 1. For an elem, get present bnd elements and connectivities of missing bnd elements
// ------------------------------------------------------------------------------------ // ------------------------------------------------------------------------------------
vector<const SMDS_MeshElement*> presentBndElems; presentBndElems.clear();
vector<TConnectivity> missingBndElems; missingBndElems.clear();
TConnectivity nodes, elemNodes; freeFacets.clear(); nodes.clear(); elemNodes.clear();
if ( vTool.Set(elem, /*ignoreCentralNodes=*/true) ) // elem is a volume -------------- if ( vTool.Set(elem, /*ignoreCentralNodes=*/true) ) // elem is a volume --------------
{ {
vTool.SetExternalNormal();
const SMDS_MeshElement* otherVol = 0; const SMDS_MeshElement* otherVol = 0;
for ( int iface = 0, n = vTool.NbFaces(); iface < n; iface++ ) for ( int iface = 0, n = vTool.NbFaces(); iface < n; iface++ )
{ {
if ( !vTool.IsFreeFace(iface, &otherVol) && if ( !vTool.IsFreeFace(iface, &otherVol) &&
( !aroundElements || elements.count( otherVol ))) ( !aroundElements || elements.count( otherVol )))
continue; continue;
freeFacets.push_back( iface );
}
if ( missType == SMDSAbs_Face )
vTool.SetExternalNormal();
for ( size_t i = 0; i < freeFacets.size(); ++i )
{
int iface = freeFacets[i];
const SMDS_MeshNode** nn = vTool.GetFaceNodes(iface); const SMDS_MeshNode** nn = vTool.GetFaceNodes(iface);
const size_t nbFaceNodes = vTool.NbFaceNodes (iface); const size_t nbFaceNodes = vTool.NbFaceNodes (iface);
if ( missType == SMDSAbs_Edge ) // boundary edges if ( missType == SMDSAbs_Edge ) // boundary edges

38
src/SMESH/SMESH_MesherHelper.cxx
View File

@ -807,6 +807,24 @@ GeomAPI_ProjectPointOnSurf& SMESH_MesherHelper::GetProjector(const TopoDS_Face&
return *( i_proj->second ); return *( i_proj->second );
} }
//=======================================================================
//function : GetSurface
//purpose : Return a cached ShapeAnalysis_Surface of a FACE
//=======================================================================
Handle(ShapeAnalysis_Surface) SMESH_MesherHelper::GetSurface(const TopoDS_Face& F )
{
Handle(Geom_Surface) surface = BRep_Tool::Surface( F );
int faceID = GetMeshDS()->ShapeToIndex( F );
TID2Surface::iterator i_surf = myFace2Surface.find( faceID );
if ( i_surf == myFace2Surface.end() && faceID )
{
Handle(ShapeAnalysis_Surface) surf( new ShapeAnalysis_Surface( surface ));
i_surf = myFace2Surface.insert( make_pair( faceID, surf )).first;
}
return i_surf->second;
}
namespace namespace
{ {
gp_XY AverageUV(const gp_XY& uv1, const gp_XY& uv2) { return ( uv1 + uv2 ) / 2.; } gp_XY AverageUV(const gp_XY& uv1, const gp_XY& uv2) { return ( uv1 + uv2 ) / 2.; }
@ -1581,6 +1599,26 @@ const SMDS_MeshNode* SMESH_MesherHelper::GetMediumNode(const SMDS_MeshNode* n1,
{ {
F = TopoDS::Face(meshDS->IndexToShape( faceID = pos.first )); F = TopoDS::Face(meshDS->IndexToShape( faceID = pos.first ));
uv[0] = GetNodeUV(F,n1,n2, force3d ? 0 : &uvOK[0]); uv[0] = GetNodeUV(F,n1,n2, force3d ? 0 : &uvOK[0]);
if ( HasDegeneratedEdges() && !force3d ) // IPAL52850 (degen VERTEX not at singularity)
{
// project middle point to a surface
SMESH_TNodeXYZ p1( n1 ), p2( n2 );
gp_Pnt pMid = 0.5 * ( p1 + p2 );
Handle(ShapeAnalysis_Surface) projector = GetSurface( F );
gp_Pnt2d uvMid;
if ( uvOK[0] )
uvMid = projector->NextValueOfUV( uv[0], pMid, BRep_Tool::Tolerance( F ));
else
uvMid = projector->ValueOfUV( pMid, getFaceMaxTol( F ));
if ( projector->Gap() * projector->Gap() < ( p1 - p2 ).SquareModulus() / 4 )
{
gp_Pnt pProj = projector->Value( uvMid );
n12 = meshDS->AddNode( pProj.X(), pProj.Y(), pProj.Z() );
meshDS->SetNodeOnFace( n12, faceID, uvMid.X(), uvMid.Y() );
myTLinkNodeMap.insert( make_pair ( link, n12 ));
return n12;
}
}
uv[1] = GetNodeUV(F,n2,n1, force3d ? 0 : &uvOK[1]); uv[1] = GetNodeUV(F,n2,n1, force3d ? 0 : &uvOK[1]);
} }
else if ( pos.second == TopAbs_EDGE ) else if ( pos.second == TopAbs_EDGE )

7
src/SMESH/SMESH_MesherHelper.hxx
View File

@ -34,6 +34,7 @@
#include <SMDS_QuadraticEdge.hxx> #include <SMDS_QuadraticEdge.hxx>
#include <Geom_Surface.hxx> #include <Geom_Surface.hxx>
#include <ShapeAnalysis_Surface.hxx>
#include <TopoDS_Face.hxx> #include <TopoDS_Face.hxx>
#include <TopoDS_Shape.hxx> #include <TopoDS_Shape.hxx>
#include <gp_Pnt2d.hxx> #include <gp_Pnt2d.hxx>
@ -528,6 +529,10 @@ public:
GeomAPI_ProjectPointOnSurf& GetProjector(const TopoDS_Face& F, GeomAPI_ProjectPointOnSurf& GetProjector(const TopoDS_Face& F,
TopLoc_Location& loc, TopLoc_Location& loc,
double tol=0 ) const; double tol=0 ) const;
/*!
* \brief Return a cached ShapeAnalysis_Surface of a FACE
*/
Handle(ShapeAnalysis_Surface) GetSurface(const TopoDS_Face& F );
/*! /*!
* \brief Check if shape is a degenerated edge or it's vertex * \brief Check if shape is a degenerated edge or it's vertex
@ -731,8 +736,10 @@ public:
std::map< int, double > myFaceMaxTol; std::map< int, double > myFaceMaxTol;
typedef std::map< int, Handle(ShapeAnalysis_Surface)> TID2Surface;
typedef std::map< int, GeomAPI_ProjectPointOnSurf* > TID2ProjectorOnSurf; typedef std::map< int, GeomAPI_ProjectPointOnSurf* > TID2ProjectorOnSurf;
typedef std::map< int, GeomAPI_ProjectPointOnCurve* > TID2ProjectorOnCurve; typedef std::map< int, GeomAPI_ProjectPointOnCurve* > TID2ProjectorOnCurve;
TID2Surface myFace2Surface;
TID2ProjectorOnSurf myFace2Projector; TID2ProjectorOnSurf myFace2Projector;
TID2ProjectorOnCurve myEdge2Projector; TID2ProjectorOnCurve myEdge2Projector;

2
src/SMESH_SWIG/smeshBuilder.py
View File

@ -4910,6 +4910,8 @@ class Mesh:
fun = self._getFunctor( funType ) fun = self._getFunctor( funType )
if fun: if fun:
if meshPart: if meshPart:
if hasattr( meshPart, "SetMesh" ):
meshPart.SetMesh( self.mesh ) # set mesh to filter
hist = fun.GetLocalHistogram( 1, False, meshPart ) hist = fun.GetLocalHistogram( 1, False, meshPart )
else: else:
hist = fun.GetHistogram( 1, False ) hist = fun.GetHistogram( 1, False )

2
src/StdMeshers/StdMeshers_FaceSide.cxx
View File

@ -410,7 +410,7 @@ const vector<UVPtStruct>& StdMeshers_FaceSide::GetUVPtStruct(bool isXConst,
if ( u2node.empty() ) return myPoints; if ( u2node.empty() ) return myPoints;
const SMDS_MeshNode* node; const SMDS_MeshNode* node;
if ( IsClosed() ) if ( IsClosed() && !proxySubMesh[0] )
node = u2node.begin()->second; node = u2node.begin()->second;
else else
{ {