IPAL52827: Very few filtering criteria available in filter dialog called from Remove Elements dialog

IPAL52831: Veeery long awaiting for Mesh Information dialog appearance
IPAL52822: Find Elements By Point does not find coincident nodes
IPAL52821: Find Elements By Point dialog: no types available for search if there are only nodes in the mesh
IPAL52823: mesh.GetSubMeshElementsId( subShape ) works wrong if subShape is retrieved indirectly
This commit is contained in:
eap 2015-07-16 20:52:26 +03:00
parent 8297100f36
commit 24fe0efaab
16 changed files with 207 additions and 86 deletions

View File

@ -11,35 +11,40 @@ import SMESH, SALOMEDS
from salome.smesh import smeshBuilder
smesh = smeshBuilder.New(salome.myStudy)
# create two faces of the box
box1 = geompy.MakeBox(0., 0., 0., 20., 20., 15.)
facesList1 = geompy.SubShapeAll(box1, geompy.ShapeType["FACE"])
face1 = facesList1[2]
# make two not sewed quadranges
OY0 = geompy.MakeVectorDXDYDZ(0, 1, 0)
OY1 = geompy.MakeTranslation( OY0, 1, 0, 0, theName="OY1" )
OY2 = geompy.MakeTranslation( OY0, 1.01, 0, 0, theName="OY2" )
OY3 = geompy.MakeTranslation( OY0, 2, 0, 0 )
q1 = geompy.MakeQuad2Edges( OY0, OY1 )
q2 = geompy.MakeQuad2Edges( OY2, OY3 )
box2 = geompy.MakeBox(0., 5., 0., 20., 20., 15.)
facesList2 = geompy.SubShapeAll(box2, geompy.ShapeType["FACE"])
face2 = facesList2[1]
edgesList = geompy.SubShapeAll(face2, geompy.ShapeType["EDGE"])
edge1 = edgesList[2]
aComp = geompy.MakeCompound([face1, face2])
geompy.addToStudy(aComp, "Two faces")
# create a mesh on two faces
mesh = smesh.Mesh(aComp, "Two faces : quadrangle mesh")
algo1D = mesh.Segment()
algo1D.NumberOfSegments(4)
algo2D = mesh.Quadrangle()
algo_local = mesh.Segment(edge1)
algo_local.Arithmetic1D(1, 4)
algo_local.Propagation()
shape = geompy.MakeCompound( [q1,q2], theName='shape' )
# make a non-uniform quadrangle mesh on two faces
mesh = smesh.Mesh(shape, "Two faces : quadrangle mesh")
mesh.Segment().Arithmetic1D( 0.1, 0.4 )
mesh.Segment(q1).NumberOfSegments( 5 )
mesh.Quadrangle()
mesh.Compute()
# sew free borders
# FirstNodeID1, SecondNodeID1, LastNodeID1,
# FirstNodeID2, SecondNodeID2, LastNodeID2, CreatePolygons, CreatePolyedrs
mesh.SewFreeBorders(6, 21, 5, 1, 12, 3, 0, 0)
segs1 = mesh.GetSubMeshElementsId( OY1 ) # mesh segments generated on borders
segs2 = mesh.GetSubMeshElementsId( OY2 )
FirstNodeID1 = mesh.GetElemNode( segs1[0], 0 )
SecondNodeID1 = mesh.GetElemNode( segs1[0], 1 )
LastNodeID1 = mesh.GetElemNode( segs1[-1], 1 )
FirstNodeID2 = mesh.GetElemNode( segs2[0], 0 )
SecondNodeID2 = mesh.GetElemNode( segs2[0], 1 )
LastNodeID2 = mesh.GetElemNode( segs2[-1], 1 )
CreatePolygons = True
CreatePolyedrs = False
res = mesh.SewFreeBorders(FirstNodeID1, SecondNodeID1, LastNodeID1,
FirstNodeID2, SecondNodeID2, LastNodeID2,
CreatePolygons, CreatePolyedrs )
print res
print "nb polygons:", mesh.NbPolygons()

View File

@ -20,7 +20,7 @@ aComp = geompy.MakeCompound([box1, box2])
geompy.addToStudy(aComp, "Two boxes")
# create a mesh on two boxes
mesh = smesh.Mesh(aComp, "Two faces : quadrangle mesh")
mesh = smesh.Mesh(aComp, "Sew Side Elements")
algo1D = mesh.Segment()
algo1D.NumberOfSegments(2)
@ -33,6 +33,31 @@ algo_local.Propagation()
mesh.Compute()
# sew side elements
# IDsOfSide1Elements, IDsOfSide2Elements,
# NodeID1OfSide1ToMerge, NodeID1OfSide2ToMerge, NodeID2OfSide1ToMerge, NodeID2OfSide2ToMerge
mesh.SewSideElements([69, 70, 71, 72], [91, 92, 89, 90], 8, 38, 23, 58)
# find elements to sew
face1 = geompy.GetFaceNearPoint( aComp, geompy.MakeVertex( 5, 10, 5 ))
IDsOfSide1Elements = mesh.GetSubMeshElementsId( face1 )
print "side faces 1:",IDsOfSide1Elements
face1Translated = geompy.MakeTranslation( face1, 0,5,0 )
faceFilter = smesh.GetFilter( SMESH.FACE, SMESH.FT_BelongToGeom,'=', face1Translated )
IDsOfSide2Elements = mesh.GetIdsFromFilter( faceFilter )
print "side faces 2:",IDsOfSide2Elements
# find corresponding nodes on sides
edge1 = geompy.GetEdgeNearPoint( aComp, geompy.MakeVertex( 0, 10, 5 ))
segs1 = mesh.GetSubMeshElementsId( edge1 ) # mesh segments generated on edge1
NodeID1OfSide1ToMerge = mesh.GetElemNode( segs1[0], 0 )
NodeID2OfSide1ToMerge = mesh.GetElemNode( segs1[0], 1 )
print "nodes of side1:", [NodeID1OfSide1ToMerge,NodeID2OfSide1ToMerge]
edge2 = geompy.GetEdgeNearPoint( aComp, geompy.MakeVertex( 0, 15, 5 ))
segs2 = mesh.GetSubMeshElementsId( edge2 ) # mesh segments generated on edge2
NodeID1OfSide2ToMerge = mesh.GetElemNode( segs2[0], 0 )
NodeID2OfSide2ToMerge = mesh.GetElemNode( segs2[0], 1 )
print "nodes of side2:", [NodeID1OfSide2ToMerge,NodeID2OfSide2ToMerge]
res = mesh.SewSideElements(IDsOfSide1Elements, IDsOfSide2Elements,
NodeID1OfSide1ToMerge, NodeID1OfSide2ToMerge,
NodeID2OfSide1ToMerge, NodeID2OfSide2ToMerge)
print res

View File

@ -50,6 +50,12 @@ corresponding end nodes of two borders will be merged. Intermediate
nodes of two borders will be either merged or inserted into faces of
the opposite border.
In practice the borders to sew often coincide and in this case it is
difficult to specify the first and the last nodes of a border since
they coincide with the first and the last nodes of the other
border. To cope with this, manually \ref merging_nodes_page to fuse
each pair of coincident nodes into one.
The sewing algorithm is as follows:
<ol>
<li>The parameter (U) of each node within a border is computed. So

View File

@ -2891,6 +2891,11 @@ void SMESHGUI_FilterDlg::Init (const int type, const bool setInViewer)
//=======================================================================
void SMESHGUI_FilterDlg::Init (const QList<int>& theTypes, const bool setInViewer)
{
if ( theTypes.empty() )
{
Init( SMESH::ALL, setInViewer );
return;
}
mySourceWg = 0;
myTypes = theTypes;
myMesh = SMESH::SMESH_Mesh::_nil();

View File

@ -665,6 +665,7 @@ void SMESHGUI_HypothesisDlg::setCustomFrame( QFrame* f )
void SMESHGUI_HypothesisDlg::accept()
{
SUIT_OverrideCursor wc; // some creators temporary set params to a hyp which can be long
QString msg;
if ( myCreator && !myCreator->checkParams( msg ) )
{

View File

@ -2258,8 +2258,11 @@ void SMESHGUI_TreeElemInfo::saveInfo( QTextStream &out )
/*!
\brief Contructor
*/
GrpComputor::GrpComputor( SMESH::SMESH_GroupBase_ptr grp, QTreeWidgetItem* item, QObject* parent )
: QObject( parent ), myItem( item )
GrpComputor::GrpComputor( SMESH::SMESH_GroupBase_ptr grp,
QTreeWidgetItem* item,
QObject* parent,
bool toComputeSize)
: QObject( parent ), myItem( item ), myToComputeSize( toComputeSize )
{
myGroup = SMESH::SMESH_GroupBase::_narrow( grp );
}
@ -2272,9 +2275,9 @@ void GrpComputor::compute()
if ( !CORBA::is_nil( myGroup ) && myItem ) {
QTreeWidgetItem* item = myItem;
myItem = 0;
int nbNodes = myGroup->GetNumberOfNodes();
int nb = myToComputeSize ? myGroup->Size() : myGroup->GetNumberOfNodes();
item->treeWidget()->removeItemWidget( item, 1 );
item->setText( 1, QString::number( nbNodes ));
item->setText( 1, QString::number( nb ));
}
}
@ -2506,10 +2509,28 @@ void SMESHGUI_AddInfo::groupInfo( SMESH::SMESH_GroupBase_ptr grp, QTreeWidgetIte
etypeItem->setText( 1, etype );
}
// size
SMESH::SMESH_Mesh_var mesh = grp->GetMesh();
bool meshLoaded = mesh->IsLoaded();
// size. Don't call grp->Size() for GroupOnFilter - issue IPAL52831
int groupSize = -1;
if ( grp->IsNodeInfoAvailable() || CORBA::is_nil( aFltGroup ))
groupSize = grp->Size();
QTreeWidgetItem* sizeItem = createItem( parent, Bold );
sizeItem->setText( 0, tr( "SIZE" ) );
sizeItem->setText( 1, QString::number( grp->Size() ) );
if ( groupSize > -1 ) {
sizeItem->setText( 1, QString::number( groupSize ) );
}
else {
QPushButton* btn = new QPushButton( tr( meshLoaded ? "COMPUTE" : "LOAD"), this );
setItemWidget( sizeItem, 1, btn );
GrpComputor* comp = new GrpComputor( grp, sizeItem, this, /*size=*/true );
connect( btn, SIGNAL( clicked() ), comp, SLOT( compute() ) );
myComputors.append( comp );
if ( !meshLoaded )
connect( btn, SIGNAL( clicked() ), this, SLOT( changeLoadToCompute() ) );
}
// color
SALOMEDS::Color color = grp->GetColor();
@ -2522,9 +2543,7 @@ void SMESHGUI_AddInfo::groupInfo( SMESH::SMESH_GroupBase_ptr grp, QTreeWidgetIte
QTreeWidgetItem* nodesItem = createItem( parent, Bold );
nodesItem->setText( 0, tr( "NB_NODES" ) );
int nbNodesLimit = SMESHGUI::resourceMgr()->integerValue( "SMESH", "info_groups_nodes_limit", 100000 );
SMESH::SMESH_Mesh_var mesh = grp->GetMesh();
bool meshLoaded = mesh->IsLoaded();
bool toShowNodes = ( grp->IsNodeInfoAvailable() || nbNodesLimit <= 0 || grp->Size() <= nbNodesLimit );
bool toShowNodes = groupSize >= 0 ? ( grp->IsNodeInfoAvailable() || nbNodesLimit <= 0 || groupSize <= nbNodesLimit ) : false;
if ( toShowNodes && meshLoaded ) {
// already calculated and up-to-date
nodesItem->setText( 1, QString::number( grp->GetNumberOfNodes() ) );

View File

@ -250,7 +250,7 @@ class GrpComputor: public QObject
Q_OBJECT;
public:
GrpComputor( SMESH::SMESH_GroupBase_ptr, QTreeWidgetItem*, QObject* );
GrpComputor( SMESH::SMESH_GroupBase_ptr, QTreeWidgetItem*, QObject*, bool = false);
QTreeWidgetItem* getItem() { return myItem; }
public slots:
@ -259,6 +259,7 @@ public slots:
private:
SMESH::SMESH_GroupBase_var myGroup;
QTreeWidgetItem* myItem;
bool myToComputeSize;
};
class SMESHGUI_EXPORT SMESHGUI_AddInfo : public QTreeWidget

View File

@ -522,6 +522,12 @@ void SMESHGUI_RemoveElementsDlg::setFilters()
if ( !myFilterDlg )
myFilterDlg = new SMESHGUI_FilterDlg( mySMESHGUI, SMESH::ALL );
QList<int> types;
if ( myMesh->NbEdges() ) types << SMESH::EDGE;
if ( myMesh->NbFaces() ) types << SMESH::FACE;
if ( myMesh->NbVolumes() ) types << SMESH::VOLUME;
myFilterDlg->Init( types );
myFilterDlg->SetSelection();
myFilterDlg->SetMesh( myMesh );
myFilterDlg->SetSourceWg( LineEditC1A1 );

View File

@ -168,6 +168,18 @@ struct SMESH_NodeSearcherImpl: public SMESH_NodeSearcher
return closestNode;
}
//---------------------------------------------------------------------
/*!
* \brief Finds nodes located within a tolerance near a point
*/
int FindNearPoint(const gp_Pnt& point,
const double tolerance,
std::vector< const SMDS_MeshNode* >& foundNodes)
{
myOctreeNode->NodesAround( point.Coord(), foundNodes, tolerance );
return foundNodes.size();
}
//---------------------------------------------------------------------
/*!
* \brief Destructor
@ -696,23 +708,23 @@ FindElementsByPoint(const gp_Pnt& point,
if ( !_nodeSearcher )
_nodeSearcher = new SMESH_NodeSearcherImpl( _mesh );
const SMDS_MeshNode* closeNode = _nodeSearcher->FindClosestTo( point );
if ( !closeNode ) return foundElements.size();
if ( point.Distance( SMESH_TNodeXYZ( closeNode )) > tolerance )
return foundElements.size(); // to far from any node
std::vector< const SMDS_MeshNode* > foundNodes;
_nodeSearcher->FindNearPoint( point, tolerance, foundNodes );
if ( type == SMDSAbs_Node )
{
foundElements.push_back( closeNode );
foundElements.assign( foundNodes.begin(), foundNodes.end() );
}
else
{
SMDS_ElemIteratorPtr elemIt = closeNode->GetInverseElementIterator( type );
for ( size_t i = 0; i < foundNodes.size(); ++i )
{
SMDS_ElemIteratorPtr elemIt = foundNodes[i]->GetInverseElementIterator( type );
while ( elemIt->more() )
foundElements.push_back( elemIt->next() );
}
}
}
// =================================================================================
else // elements more complex than 0D
{

View File

@ -55,6 +55,9 @@ struct SMESHUtils_EXPORT SMESH_NodeSearcher
{
virtual const SMDS_MeshNode* FindClosestTo( const gp_Pnt& pnt ) = 0;
virtual void MoveNode( const SMDS_MeshNode* node, const gp_Pnt& toPnt ) = 0;
virtual int FindNearPoint(const gp_Pnt& point,
const double tolerance,
std::vector< const SMDS_MeshNode* >& foundNodes) = 0;
};
//=======================================================================

View File

@ -189,6 +189,7 @@ void SMESH_OctreeNode::NodesAround (const SMDS_MeshNode * Node,
//================================================================================
/*!
* \brief Return in dist2Nodes nodes mapped to their square distance from Node
* Tries to find a closest node.
* \param node - node to find nodes closest to
* \param dist2Nodes - map of found nodes and their distances
* \param precision - radius of a sphere to check nodes inside
@ -241,6 +242,44 @@ bool SMESH_OctreeNode::NodesAround(const gp_XYZ &node,
return false;
}
//================================================================================
/*!
* \brief Return a list of nodes close to a point
* \param [in] point - point
* \param [out] nodes - found nodes
* \param [in] precision - allowed distance from \a point
*/
//================================================================================
void SMESH_OctreeNode::NodesAround(const gp_XYZ& point,
std::vector<const SMDS_MeshNode*>& nodes,
double precision)
{
if ( isInside( point, precision ))
{
if ( isLeaf() && NbNodes() )
{
double minDist2 = precision * precision;
TIDSortedNodeSet::iterator nIt = myNodes.begin();
for ( ; nIt != myNodes.end(); ++nIt )
{
SMESH_TNodeXYZ p2( *nIt );
double dist2 = ( point - p2 ).SquareModulus();
if ( dist2 <= minDist2 )
nodes.push_back( p2._node );
}
}
else if ( myChildren )
{
for (int i = 0; i < 8; i++)
{
SMESH_OctreeNode* myChild = dynamic_cast<SMESH_OctreeNode*> (myChildren[i]);
myChild->NodesAround( point, nodes, precision);
}
}
}
}
//=============================
/*!
* \brief Return in theGroupsOfNodes a list of group of nodes close to each other within theTolerance

View File

@ -30,16 +30,17 @@
#ifndef _SMESH_OCTREENODE_HXX_
#define _SMESH_OCTREENODE_HXX_
#include "SMESH_Utils.hxx"
#include "SMESH_Octree.hxx"
#include <gp_Pnt.hxx>
#include "SMDS_ElemIterator.hxx"
#include "SMDS_MeshNode.hxx"
#include "SMESH_Octree.hxx"
#include "SMESH_Utils.hxx"
#include <gp_Pnt.hxx>
#include <list>
#include <set>
#include <map>
#include "SMDS_ElemIterator.hxx"
#include <vector>
//forward declaration
class SMDS_MeshNode;
@ -49,34 +50,35 @@ typedef SMDS_Iterator<SMESH_OctreeNode*> SMESH_OctreeNodeIterator;
typedef boost::shared_ptr<SMESH_OctreeNodeIterator> SMESH_OctreeNodeIteratorPtr;
typedef std::set< const SMDS_MeshNode*, TIDCompare > TIDSortedNodeSet;
class SMESHUtils_EXPORT SMESH_OctreeNode : public SMESH_Octree {
public:
class SMESHUtils_EXPORT SMESH_OctreeNode : public SMESH_Octree
{
public:
// Constructor
SMESH_OctreeNode (const TIDSortedNodeSet& theNodes, const int maxLevel = 8,
const int maxNbNodes = 5, const double minBoxSize = 0.);
//=============================
/*!
* \brief Empty destructor
*/
//=============================
// destructor
virtual ~SMESH_OctreeNode () {};
// Tells us if Node is inside the current box with the precision "precision"
virtual const bool isInside(const gp_XYZ& p, const double precision = 0.);
// Return in Result a list of Nodes potentials to be near Node
void NodesAround(const SMDS_MeshNode * Node,
std::list<const SMDS_MeshNode*>* Result,
void NodesAround(const SMDS_MeshNode * node,
std::list<const SMDS_MeshNode*>* result,
const double precision = 0.);
// Return in dist2Nodes nodes mapped to their square distance from Node
bool NodesAround(const gp_XYZ& node,
bool NodesAround(const gp_XYZ& point,
std::map<double, const SMDS_MeshNode*>& dist2Nodes,
double precision);
// Return a list of Nodes close to a point
void NodesAround(const gp_XYZ& point,
std::vector<const SMDS_MeshNode*>& nodes,
double precision);
// Return in theGroupsOfNodes a list of group of nodes close to each other within theTolerance
// Search for all the nodes in nodes
void FindCoincidentNodes ( TIDSortedNodeSet* nodes,

View File

@ -5141,6 +5141,8 @@ SMESH::array_of_ElementType* SMESH_Mesh_i::GetTypes()
if (_impl->NbVolumes()) types[nbTypes++] = SMESH::VOLUME;
if (_impl->Nb0DElements()) types[nbTypes++] = SMESH::ELEM0D;
if (_impl->NbBalls()) types[nbTypes++] = SMESH::BALL;
if (_impl->NbNodes() &&
nbTypes == 0 ) types[nbTypes++] = SMESH::NODE;
types->length( nbTypes );
return types._retn();

View File

@ -2505,8 +2505,8 @@ class Mesh:
# @return the list of integer values
# @ingroup l1_meshinfo
def GetSubMeshElementsId(self, Shape):
if ( isinstance( Shape, geomBuilder.GEOM._objref_GEOM_Object)):
ShapeID = Shape.GetSubShapeIndices()[0]
if isinstance( Shape, geomBuilder.GEOM._objref_GEOM_Object):
ShapeID = self.geompyD.GetSubShapeID( self.geom, Shape )
else:
ShapeID = Shape
return self.mesh.GetSubMeshElementsId(ShapeID)
@ -2518,7 +2518,7 @@ class Mesh:
# @return the list of integer values
# @ingroup l1_meshinfo
def GetSubMeshNodesId(self, Shape, all):
if ( isinstance( Shape, geomBuilder.GEOM._objref_GEOM_Object)):
if isinstance( Shape, geomBuilder.GEOM._objref_GEOM_Object):
ShapeID = self.geompyD.GetSubShapeID( self.geom, Shape )
else:
ShapeID = Shape
@ -2530,8 +2530,8 @@ class Mesh:
# @return element type
# @ingroup l1_meshinfo
def GetSubMeshElementType(self, Shape):
if ( isinstance( Shape, geomBuilder.GEOM._objref_GEOM_Object)):
ShapeID = Shape.GetSubShapeIndices()[0]
if isinstance( Shape, geomBuilder.GEOM._objref_GEOM_Object):
ShapeID = self.geompyD.GetSubShapeID( self.geom, Shape )
else:
ShapeID = Shape
return self.mesh.GetSubMeshElementType(ShapeID)
@ -2920,7 +2920,7 @@ class Mesh:
# @ingroup l2_modif_add
def SetNodeOnVertex(self, NodeID, Vertex):
if ( isinstance( Vertex, geomBuilder.GEOM._objref_GEOM_Object)):
VertexID = Vertex.GetSubShapeIndices()[0]
VertexID = self.geompyD.GetSubShapeID( self.geom, Vertex )
else:
VertexID = Vertex
try:
@ -2938,7 +2938,7 @@ class Mesh:
# @ingroup l2_modif_add
def SetNodeOnEdge(self, NodeID, Edge, paramOnEdge):
if ( isinstance( Edge, geomBuilder.GEOM._objref_GEOM_Object)):
EdgeID = Edge.GetSubShapeIndices()[0]
EdgeID = self.geompyD.GetSubShapeID( self.geom, Edge )
else:
EdgeID = Edge
try:
@ -2956,7 +2956,7 @@ class Mesh:
# @ingroup l2_modif_add
def SetNodeOnFace(self, NodeID, Face, u, v):
if ( isinstance( Face, geomBuilder.GEOM._objref_GEOM_Object)):
FaceID = Face.GetSubShapeIndices()[0]
FaceID = self.geompyD.GetSubShapeID( self.geom, Face )
else:
FaceID = Face
try:
@ -2972,7 +2972,7 @@ class Mesh:
# @ingroup l2_modif_add
def SetNodeInVolume(self, NodeID, Solid):
if ( isinstance( Solid, geomBuilder.GEOM._objref_GEOM_Object)):
SolidID = Solid.GetSubShapeIndices()[0]
SolidID = self.geompyD.GetSubShapeID( self.geom, Solid )
else:
SolidID = Solid
try:
@ -2988,7 +2988,7 @@ class Mesh:
# @ingroup l2_modif_add
def SetMeshElementOnShape(self, ElementID, Shape):
if ( isinstance( Shape, geomBuilder.GEOM._objref_GEOM_Object)):
ShapeID = Shape.GetSubShapeIndices()[0]
ShapeID = self.geompyD.GetSubShapeID( self.geom, Shape )
else:
ShapeID = Shape
try:
@ -4589,7 +4589,7 @@ class Mesh:
def ClearLastCreated(self):
self.editor.ClearLastCreated()
## Creates Duplicates given elements, i.e. creates new elements based on the
## Creates duplicates of given elements, i.e. creates new elements based on the
# same nodes as the given ones.
# @param theElements - container of elements to duplicate. It can be a Mesh,
# sub-mesh, group, filter or a list of element IDs. If \a theElements is

View File

@ -288,12 +288,7 @@ class Mesh_Algorithm:
raise TypeError, "ViscousLayers are not supported by %s"%self.algo.GetName()
if faces and isinstance( faces[0], geomBuilder.GEOM._objref_GEOM_Object ):
import GEOM
faceIDs = []
for f in faces:
if self.mesh.geompyD.ShapeIdToType( f.GetType() ) == "GROUP":
faceIDs += f.GetSubShapeIndices()
else:
faceIDs += [self.mesh.geompyD.GetSubShapeID(self.mesh.geom, f)]
faceIDs = [self.mesh.geompyD.GetSubShapeID(self.mesh.geom, f) for f in faces]
faces = faceIDs
hyp = self.Hypothesis("ViscousLayers",
[thickness, numberOfLayers, stretchFactor, faces, isFacesToIgnore],