Back to index

salome-smesh  6.5.0
SMDS_MeshNode.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2007-2012  CEA/DEN, EDF R&D, OPEN CASCADE
00002 //
00003 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
00004 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
00005 //
00006 // This library is free software; you can redistribute it and/or
00007 // modify it under the terms of the GNU Lesser General Public
00008 // License as published by the Free Software Foundation; either
00009 // version 2.1 of the License.
00010 //
00011 // This library is distributed in the hope that it will be useful,
00012 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00013 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00014 // Lesser General Public License for more details.
00015 //
00016 // You should have received a copy of the GNU Lesser General Public
00017 // License along with this library; if not, write to the Free Software
00018 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00019 //
00020 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00021 //
00022 
00023 //  SMESH SMDS : implementaion of Salome mesh data structure
00024 //
00025 #ifdef _MSC_VER
00026 #pragma warning(disable:4786)
00027 #endif
00028 
00029 #include "SMDS_MeshNode.hxx"
00030 #include "SMDS_SpacePosition.hxx"
00031 #include "SMDS_IteratorOfElements.hxx"
00032 #include "SMDS_Mesh.hxx"
00033 #include <vtkUnstructuredGrid.h>
00034 
00035 #include "utilities.h"
00036 #include "Utils_SALOME_Exception.hxx"
00037 #include <cassert>
00038 
00039 using namespace std;
00040 
00041 int SMDS_MeshNode::nbNodes =0;
00042 
00043 //=======================================================================
00044 //function : SMDS_MeshNode
00045 //purpose  : 
00046 //=======================================================================
00047 SMDS_MeshNode::SMDS_MeshNode() :
00048   SMDS_MeshElement(-1, -1, 0),
00049   myPosition(SMDS_SpacePosition::originSpacePosition())
00050 {
00051   nbNodes++;
00052 }
00053 
00054 SMDS_MeshNode::SMDS_MeshNode(int id, int meshId, int shapeId, double x, double y, double z):
00055   SMDS_MeshElement(id, meshId, shapeId),
00056   myPosition(SMDS_SpacePosition::originSpacePosition())
00057 {
00058   nbNodes++;
00059   init(id, meshId, shapeId, x, y ,z);
00060 }
00061 
00062 void SMDS_MeshNode::init(int id, int meshId, int shapeId, double x, double y, double z)
00063 {
00064   SMDS_MeshElement::init(id, meshId, shapeId);
00065   myVtkID = id -1;
00066   assert(myVtkID >= 0);
00067   //MESSAGE("Node " << myID << " " << myVtkID << " (" << x << ", " << y << ", " << z << ")");
00068   SMDS_Mesh* mesh = SMDS_Mesh::_meshList[myMeshId];
00069   SMDS_UnstructuredGrid * grid = mesh->getGrid();
00070   vtkPoints *points = grid->GetPoints();
00071   points->InsertPoint(myVtkID, x, y, z);
00072   SMDS_CellLinks *cellLinks = dynamic_cast<SMDS_CellLinks*>(grid->GetCellLinks());
00073   assert(cellLinks);
00074   if (myVtkID >= cellLinks->GetLinksSize())
00075           cellLinks->ResizeL(myVtkID+SMDS_Mesh::chunkSize);
00076 }
00077 
00078 SMDS_MeshNode::~SMDS_MeshNode()
00079 {
00080   nbNodes--;
00081   if ( myPosition && myPosition != SMDS_SpacePosition::originSpacePosition() )
00082     delete myPosition, myPosition = 0;
00083 }
00084 
00085 //=======================================================================
00086 //function : RemoveInverseElement
00087 //purpose  : 
00088 //=======================================================================
00089 
00090 void SMDS_MeshNode::RemoveInverseElement(const SMDS_MeshElement * parent)
00091 {
00092     //MESSAGE("RemoveInverseElement " << myID << " " << parent->GetID());
00093     const SMDS_MeshCell* cell = dynamic_cast<const SMDS_MeshCell*>(parent);
00094     MYASSERT(cell);
00095     SMDS_Mesh::_meshList[myMeshId]->getGrid()->RemoveReferenceToCell(myVtkID, cell->getVtkId());
00096 }
00097 
00098 //=======================================================================
00099 //function : Print
00100 //purpose  : 
00101 //=======================================================================
00102 
00103 void SMDS_MeshNode::Print(ostream & OS) const
00104 {
00105         OS << "Node <" << myID << "> : X = " << X() << " Y = "
00106                 << Y() << " Z = " << Z() << endl;
00107 }
00108 
00109 //=======================================================================
00110 //function : SetPosition
00111 //purpose  : 
00112 //=======================================================================
00113 
00114 void SMDS_MeshNode::SetPosition(const SMDS_PositionPtr& aPos)
00115 {
00116   if ( myPosition &&
00117        myPosition != SMDS_SpacePosition::originSpacePosition() &&
00118        myPosition != aPos )
00119     delete myPosition;
00120   myPosition = aPos;
00121 }
00122 
00123 //=======================================================================
00124 //function : GetPosition
00125 //purpose  : 
00126 //=======================================================================
00127 
00128 const SMDS_PositionPtr& SMDS_MeshNode::GetPosition() const
00129 {
00130         return myPosition;
00131 }
00132 
00133 //=======================================================================
00137 //=======================================================================
00138 
00139 class SMDS_MeshNode_MyInvIterator: public SMDS_ElemIterator
00140 {
00141 private:
00142   SMDS_Mesh* myMesh;
00143   vtkIdType* myCells;
00144   int myNcells;
00145   SMDSAbs_ElementType myType;
00146   int iter;
00147   vector<vtkIdType> cellList;
00148 
00149 public:
00150   SMDS_MeshNode_MyInvIterator(SMDS_Mesh *mesh, vtkIdType* cells, int ncells, SMDSAbs_ElementType type) :
00151     myMesh(mesh), myCells(cells), myNcells(ncells), myType(type), iter(0)
00152   {
00153     //MESSAGE("SMDS_MeshNode_MyInvIterator : ncells " << myNcells);
00154     cellList.clear();
00155     if (type == SMDSAbs_All)
00156       for (int i = 0; i < ncells; i++)
00157         cellList.push_back(cells[i]);
00158     else for (int i = 0; i < ncells; i++)
00159       {
00160         int vtkId = cells[i];
00161         int smdsId = myMesh->fromVtkToSmds(vtkId);
00162         const SMDS_MeshElement* elem = myMesh->FindElement(smdsId);
00163         if (elem->GetType() == type)
00164           {
00165             //MESSAGE("Add element vtkId " << vtkId << " " << elem->GetType())
00166             cellList.push_back(vtkId);
00167           }
00168       }
00169     myCells = &cellList[0];
00170     myNcells = cellList.size();
00171     //MESSAGE("myNcells="<<myNcells);
00172   }
00173 
00174   bool more()
00175   {
00176     //MESSAGE("iter " << iter << " ncells " << myNcells);
00177     return (iter < myNcells);
00178   }
00179 
00180   const SMDS_MeshElement* next()
00181   {
00182     int vtkId = myCells[iter];
00183     int smdsId = myMesh->fromVtkToSmds(vtkId);
00184     const SMDS_MeshElement* elem = myMesh->FindElement(smdsId);
00185     if (!elem)
00186       {
00187         MESSAGE("SMDS_MeshNode_MyInvIterator problem Null element");
00188         throw SALOME_Exception("SMDS_MeshNode_MyInvIterator problem Null element");
00189       }
00190     //MESSAGE("vtkId " << vtkId << " smdsId " << smdsId << " " << elem->GetType());
00191     iter++;
00192     return elem;
00193   }
00194 };
00195 
00196 SMDS_ElemIteratorPtr SMDS_MeshNode::
00197         GetInverseElementIterator(SMDSAbs_ElementType type) const
00198 {
00199     vtkCellLinks::Link l = SMDS_Mesh::_meshList[myMeshId]->getGrid()->GetCellLinks()->GetLink(myVtkID);
00200     //MESSAGE("myID " << myID << " ncells " << l.ncells);
00201     return SMDS_ElemIteratorPtr(new SMDS_MeshNode_MyInvIterator(SMDS_Mesh::_meshList[myMeshId], l.cells, l.ncells, type));
00202 }
00203 
00204 // Same as GetInverseElementIterator but the create iterator only return
00205 // wanted type elements.
00206 class SMDS_MeshNode_MyIterator:public SMDS_ElemIterator
00207 {
00208 private:
00209   SMDS_Mesh* myMesh;
00210   vtkIdType* myCells;
00211   int  myNcells;
00212   SMDSAbs_ElementType                                 myType;
00213   int  iter;
00214   vector<SMDS_MeshElement*> myFiltCells;
00215 
00216  public:
00217   SMDS_MeshNode_MyIterator(SMDS_Mesh *mesh,
00218                            vtkIdType* cells,
00219                            int ncells,
00220                            SMDSAbs_ElementType type):
00221     myMesh(mesh), myCells(cells), myNcells(ncells), myType(type), iter(0)
00222   {
00223         //MESSAGE("myNcells " << myNcells);
00224        for (; iter<ncells; iter++)
00225         {
00226            int vtkId = myCells[iter];
00227            int smdsId = myMesh->fromVtkToSmds(vtkId);
00228            //MESSAGE("vtkId " << vtkId << " smdsId " << smdsId);
00229            const SMDS_MeshElement* elem = myMesh->FindElement(smdsId);
00230            if (elem->GetType() == type)
00231                myFiltCells.push_back((SMDS_MeshElement*)elem);
00232         }
00233         myNcells = myFiltCells.size();
00234         //MESSAGE("myNcells " << myNcells);
00235        iter = 0;
00236         //MESSAGE("SMDS_MeshNode_MyIterator (filter) " << ncells << " " << myNcells);
00237   }
00238 
00239   bool more()
00240   {
00241       return (iter< myNcells);
00242   }
00243 
00244   const SMDS_MeshElement* next()
00245   {
00246       const SMDS_MeshElement* elem = myFiltCells[iter];
00247       iter++;
00248       return elem;
00249   }
00250 };
00251 
00252 SMDS_ElemIteratorPtr SMDS_MeshNode::
00253         elementsIterator(SMDSAbs_ElementType type) const
00254 {
00255   if(type==SMDSAbs_Node)
00256     return SMDS_MeshElement::elementsIterator(SMDSAbs_Node); 
00257   else
00258   {
00259     vtkCellLinks::Link l = SMDS_Mesh::_meshList[myMeshId]->getGrid()->GetCellLinks()->GetLink(myVtkID);
00260     return SMDS_ElemIteratorPtr(new SMDS_MeshNode_MyIterator(SMDS_Mesh::_meshList[myMeshId], l.cells, l.ncells, type));
00261   }
00262 }
00263 
00264 int SMDS_MeshNode::NbNodes() const
00265 {
00266         return 1;
00267 }
00268 
00269 double* SMDS_MeshNode::getCoord() const
00270 {
00271   return SMDS_Mesh::_meshList[myMeshId]->getGrid()->GetPoint(myVtkID);
00272 }
00273 
00274 double SMDS_MeshNode::X() const
00275 {
00276   double *coord = getCoord();
00277   return coord[0];
00278 }
00279 
00280 double SMDS_MeshNode::Y() const
00281 {
00282   double *coord = getCoord();
00283   return coord[1];
00284 }
00285 
00286 double SMDS_MeshNode::Z() const
00287 {
00288   double *coord = getCoord();
00289   return coord[2];
00290 }
00291 
00292 //================================================================================
00296 void SMDS_MeshNode::GetXYZ(double xyz[3]) const
00297 {
00298   return SMDS_Mesh::_meshList[myMeshId]->getGrid()->GetPoint(myVtkID,xyz);
00299 }
00300 
00301 //* resize the vtkPoints structure every SMDS_Mesh::chunkSize points
00302 void SMDS_MeshNode::setXYZ(double x, double y, double z)
00303 {
00304   SMDS_Mesh *mesh = SMDS_Mesh::_meshList[myMeshId];
00305   vtkPoints *points = mesh->getGrid()->GetPoints();
00306   points->InsertPoint(myVtkID, x, y, z);
00307   mesh->adjustBoundingBox(x, y, z);
00308   mesh->setMyModified();
00309 }
00310 
00311 SMDSAbs_ElementType SMDS_MeshNode::GetType() const
00312 {
00313         return SMDSAbs_Node;
00314 }
00315 
00316 vtkIdType SMDS_MeshNode::GetVtkType() const
00317 {
00318   return VTK_VERTEX;
00319 }
00320 
00321 //=======================================================================
00322 //function : AddInverseElement
00323 //purpose  :
00324 //=======================================================================
00325 void SMDS_MeshNode::AddInverseElement(const SMDS_MeshElement* ME)
00326 {
00327   const SMDS_MeshCell *cell = dynamic_cast<const SMDS_MeshCell*> (ME);
00328   assert(cell);
00329   SMDS_UnstructuredGrid* grid = SMDS_Mesh::_meshList[myMeshId]->getGrid();
00330   vtkCellLinks *Links = grid->GetCellLinks();
00331   Links->ResizeCellList(myVtkID, 1);
00332   Links->AddCellReference(cell->getVtkId(), myVtkID);
00333 }
00334 
00335 //=======================================================================
00336 //function : ClearInverseElements
00337 //purpose  :
00338 //=======================================================================
00339 void SMDS_MeshNode::ClearInverseElements()
00340 {
00341   SMDS_Mesh::_meshList[myMeshId]->getGrid()->ResizeCellList(myVtkID, 0);
00342 }
00343 
00344 bool SMDS_MeshNode::emptyInverseElements()
00345 {
00346   vtkCellLinks::Link l = SMDS_Mesh::_meshList[myMeshId]->getGrid()->GetCellLinks()->GetLink(myVtkID);
00347   return (l.ncells == 0);
00348 }
00349 
00350 //================================================================================
00354 //================================================================================
00355 
00356 int SMDS_MeshNode::NbInverseElements(SMDSAbs_ElementType type) const
00357 {
00358   vtkCellLinks::Link l = SMDS_Mesh::_meshList[myMeshId]->getGrid()->GetCellLinks()->GetLink(myVtkID);
00359 
00360   if ( type == SMDSAbs_All )
00361     return l.ncells;
00362 
00363   int nb = 0;
00364   SMDS_Mesh *mesh = SMDS_Mesh::_meshList[myMeshId];
00365   for (int i=0; i<l.ncells; i++)
00366         {
00367            const SMDS_MeshElement* elem = mesh->FindElement(mesh->fromVtkToSmds(l.cells[i]));
00368            if (elem->GetType() == type)
00369                 nb++;
00370         }
00371   return nb;
00372 }
00373 
00377 bool operator<(const SMDS_MeshNode& e1, const SMDS_MeshNode& e2)
00378 {
00379         return e1.getVtkId()<e2.getVtkId();
00380         /*if(e1.myX<e2.myX) return true;
00381         else if(e1.myX==e2.myX)
00382         {
00383                 if(e1.myY<e2.myY) return true;
00384                 else if(e1.myY==e2.myY) return (e1.myZ<e2.myZ);
00385                 else return false;
00386         }
00387         else return false;*/
00388 }
00389