Back to index

salome-smesh  6.5.0
SMDS_UnstructuredGrid.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2010-2012  CEA/DEN, EDF R&D, OPEN CASCADE
00002 //
00003 // This library is free software; you can redistribute it and/or
00004 // modify it under the terms of the GNU Lesser General Public
00005 // License as published by the Free Software Foundation; either
00006 // version 2.1 of the License.
00007 //
00008 // This library is distributed in the hope that it will be useful,
00009 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00010 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00011 // Lesser General Public License for more details.
00012 //
00013 // You should have received a copy of the GNU Lesser General Public
00014 // License along with this library; if not, write to the Free Software
00015 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00016 //
00017 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00018 //
00019 
00020 #define CHRONODEF
00021 #include "SMDS_UnstructuredGrid.hxx"
00022 #include "SMDS_Mesh.hxx"
00023 #include "SMDS_MeshInfo.hxx"
00024 #include "SMDS_Downward.hxx"
00025 #include "SMDS_MeshVolume.hxx"
00026 
00027 #include "utilities.h"
00028 
00029 #include <vtkCellArray.h>
00030 #include <vtkCellLinks.h>
00031 #include <vtkIdTypeArray.h>
00032 #include <vtkUnsignedCharArray.h>
00033 
00034 #include <list>
00035 #include <climits>
00036 
00037 using namespace std;
00038 
00039 SMDS_CellLinks* SMDS_CellLinks::New()
00040 {
00041   MESSAGE("SMDS_CellLinks::New");
00042   return new SMDS_CellLinks();
00043 }
00044 
00045 vtkCellLinks::Link* SMDS_CellLinks::ResizeL(vtkIdType sz)
00046 {
00047   return vtkCellLinks::Resize(sz);
00048 }
00049 
00050 vtkIdType SMDS_CellLinks::GetLinksSize()
00051 {
00052   return this->Size;
00053 }
00054 
00055 SMDS_CellLinks::SMDS_CellLinks() :
00056   vtkCellLinks()
00057 {
00058 }
00059 
00060 SMDS_CellLinks::~SMDS_CellLinks()
00061 {
00062 }
00063 
00064 SMDS_UnstructuredGrid* SMDS_UnstructuredGrid::New()
00065 {
00066   MESSAGE("SMDS_UnstructuredGrid::New");
00067   return new SMDS_UnstructuredGrid();
00068 }
00069 
00070 SMDS_UnstructuredGrid::SMDS_UnstructuredGrid() :
00071   vtkUnstructuredGrid()
00072 {
00073   _cellIdToDownId.clear();
00074   _downTypes.clear();
00075   _downArray.clear();
00076   _mesh = 0;
00077 }
00078 
00079 SMDS_UnstructuredGrid::~SMDS_UnstructuredGrid()
00080 {
00081 }
00082 
00083 unsigned long SMDS_UnstructuredGrid::GetMTime()
00084 {
00085   unsigned long mtime = vtkUnstructuredGrid::GetMTime();
00086   MESSAGE("vtkUnstructuredGrid::GetMTime: " << mtime);
00087   return mtime;
00088 }
00089 
00090 void SMDS_UnstructuredGrid::Update()
00091 {
00092   MESSAGE("SMDS_UnstructuredGrid::Update");
00093   return vtkUnstructuredGrid::Update();
00094 }
00095 
00096 void SMDS_UnstructuredGrid::UpdateInformation()
00097 {
00098   MESSAGE("SMDS_UnstructuredGrid::UpdateInformation");
00099   return vtkUnstructuredGrid::UpdateInformation();
00100 }
00101 
00102 vtkPoints* SMDS_UnstructuredGrid::GetPoints()
00103 {
00104   // TODO erreur incomprehensible de la macro vtk GetPoints apparue avec la version paraview de fin aout 2010
00105   //MESSAGE("*********************** SMDS_UnstructuredGrid::GetPoints " << this->Points << " " << vtkUnstructuredGrid::GetPoints());
00106   return this->Points;
00107 }
00108 
00109 //#ifdef VTK_HAVE_POLYHEDRON
00110 int SMDS_UnstructuredGrid::InsertNextLinkedCell(int type, int npts, vtkIdType *pts)
00111 {
00112   if (type != VTK_POLYHEDRON)
00113     return vtkUnstructuredGrid::InsertNextLinkedCell(type, npts, pts);
00114 
00115   // --- type = VTK_POLYHEDRON
00116   //MESSAGE("InsertNextLinkedCell VTK_POLYHEDRON");
00117   int cellid = this->InsertNextCell(type, npts, pts);
00118 
00119   set<vtkIdType> setOfNodes;
00120   setOfNodes.clear();
00121   int nbfaces = npts;
00122   int i = 0;
00123   for (int nf = 0; nf < nbfaces; nf++)
00124     {
00125       int nbnodes = pts[i];
00126       i++;
00127       for (int k = 0; k < nbnodes; k++)
00128         {
00129           //MESSAGE(" cell " << cellid << " face " << nf << " node " << pts[i]);
00130           setOfNodes.insert(pts[i]);
00131           i++;
00132         }
00133     }
00134 
00135   set<vtkIdType>::iterator it = setOfNodes.begin();
00136   for (; it != setOfNodes.end(); ++it)
00137     {
00138       //MESSAGE("reverse link for node " << *it << " cell " << cellid);
00139       this->Links->ResizeCellList(*it, 1);
00140       this->Links->AddCellReference(cellid, *it);
00141     }
00142 
00143   return cellid;
00144 }
00145 //#endif
00146 
00147 void SMDS_UnstructuredGrid::setSMDS_mesh(SMDS_Mesh *mesh)
00148 {
00149   _mesh = mesh;
00150 }
00151 
00152 void SMDS_UnstructuredGrid::compactGrid(std::vector<int>& idNodesOldToNew, int newNodeSize,
00153                                         std::vector<int>& idCellsOldToNew, int newCellSize)
00154 {
00155   MESSAGE("------------------------- SMDS_UnstructuredGrid::compactGrid " << newNodeSize << " " << newCellSize);CHRONO(1);
00156   int alreadyCopied = 0;
00157 
00158   // --- if newNodeSize, create a new compacted vtkPoints
00159 
00160   vtkPoints *newPoints = vtkPoints::New();
00161   newPoints->SetDataType(VTK_DOUBLE);
00162   newPoints->SetNumberOfPoints(newNodeSize);
00163   if (newNodeSize)
00164     {
00165       MESSAGE("-------------- compactGrid, newNodeSize " << newNodeSize);
00166       // rnv: to fix bug "21125: EDF 1233 SMESH: Degradation of precision in a test case for quadratic conversion"
00167       // using double type for storing coordinates of nodes instead float.
00168       int oldNodeSize = idNodesOldToNew.size();
00169 
00170       int i = 0;
00171       while ( i < oldNodeSize )
00172       {
00173         // skip a hole if any
00174         while ( i < oldNodeSize && idNodesOldToNew[i] < 0 )
00175           ++i;
00176         int startBloc = i;
00177         // look for a block end
00178         while ( i < oldNodeSize && idNodesOldToNew[i] >= 0 )
00179           ++i;
00180         int endBloc = i;
00181         copyNodes(newPoints, idNodesOldToNew, alreadyCopied, startBloc, endBloc);
00182       }
00183       newPoints->Squeeze();
00184     }
00185 
00186   // --- create new compacted Connectivity, Locations and Types
00187 
00188   int oldCellSize = this->Types->GetNumberOfTuples();
00189 
00190   vtkCellArray *newConnectivity = vtkCellArray::New();
00191   newConnectivity->Initialize();
00192   int oldCellDataSize = this->Connectivity->GetData()->GetSize();
00193   newConnectivity->Allocate(oldCellDataSize);
00194   MESSAGE("oldCellSize="<< oldCellSize << " oldCellDataSize=" << oldCellDataSize);
00195 
00196   vtkUnsignedCharArray *newTypes = vtkUnsignedCharArray::New();
00197   newTypes->Initialize();
00198   newTypes->SetNumberOfValues(newCellSize);
00199 
00200   vtkIdTypeArray *newLocations = vtkIdTypeArray::New();
00201   newLocations->Initialize();
00202   newLocations->SetNumberOfValues(newCellSize);
00203 
00204   // TODO some polyhedron may be huge (only in some tests)
00205   vtkIdType tmpid[NBMAXNODESINCELL];
00206   vtkIdType *pointsCell = &tmpid[0]; // --- points id to fill a new cell
00207 
00208   alreadyCopied = 0;
00209   int i = 0;
00210   while ( i < oldCellSize )
00211   {
00212     // skip a hole if any
00213     while ( i < oldCellSize && this->Types->GetValue(i) == VTK_EMPTY_CELL )
00214       ++i;
00215     int startBloc = i;
00216     // look for a block end
00217     while ( i < oldCellSize && this->Types->GetValue(i) != VTK_EMPTY_CELL )
00218       ++i;
00219     int endBloc = i;
00220     if ( endBloc > startBloc )
00221       copyBloc(newTypes, idCellsOldToNew, idNodesOldToNew, newConnectivity, newLocations, pointsCell,
00222                alreadyCopied, startBloc, endBloc);
00223   }
00224 
00225   newConnectivity->Squeeze();
00226 
00227   if (1/*newNodeSize*/)
00228     {
00229       MESSAGE("------- newNodeSize, setPoints");
00230       this->SetPoints(newPoints);
00231       MESSAGE("NumberOfPoints: " << this->GetNumberOfPoints());
00232     }
00233 
00234   if (this->FaceLocations)
00235     {
00236       vtkIdTypeArray *newFaceLocations = vtkIdTypeArray::New();
00237       newFaceLocations->Initialize();
00238       newFaceLocations->Allocate(newTypes->GetSize());
00239       vtkIdTypeArray *newFaces = vtkIdTypeArray::New();
00240       newFaces->Initialize();
00241       newFaces->Allocate(this->Faces->GetSize());
00242       for (int i = 0; i < oldCellSize; i++)
00243         {
00244           if (this->Types->GetValue(i) == VTK_EMPTY_CELL)
00245             continue;
00246           int newCellId = idCellsOldToNew[i];
00247           if (newTypes->GetValue(newCellId) == VTK_POLYHEDRON)
00248             {
00249                newFaceLocations->InsertNextValue(newFaces->GetMaxId()+1);
00250                int oldFaceLoc = this->FaceLocations->GetValue(i);
00251                int nCellFaces = this->Faces->GetValue(oldFaceLoc++);
00252                newFaces->InsertNextValue(nCellFaces);
00253                for (int n=0; n<nCellFaces; n++)
00254                  {
00255                    int nptsInFace = this->Faces->GetValue(oldFaceLoc++);
00256                    newFaces->InsertNextValue(nptsInFace);
00257                    for (int k=0; k<nptsInFace; k++)
00258                      {
00259                        int oldpt = this->Faces->GetValue(oldFaceLoc++);
00260                        newFaces->InsertNextValue(idNodesOldToNew[oldpt]);
00261                      }
00262                  }
00263             }
00264           else
00265             {
00266                newFaceLocations->InsertNextValue(-1);
00267             }
00268         }
00269       newFaceLocations->Squeeze();
00270       newFaces->Squeeze();
00271       newFaceLocations->Register(this);
00272       newFaces->Register(this);
00273       this->SetCells(newTypes, newLocations, newConnectivity, newFaceLocations, newFaces);
00274       newFaceLocations->Delete();
00275       newFaces->Delete();
00276     }
00277   else
00278     this->SetCells(newTypes, newLocations, newConnectivity, FaceLocations, Faces);
00279 
00280   newPoints->Delete();
00281   newTypes->Delete();
00282   newLocations->Delete();
00283   newConnectivity->Delete();
00284   this->BuildLinks();
00285 }
00286 
00287 void SMDS_UnstructuredGrid::copyNodes(vtkPoints *newPoints, std::vector<int>& idNodesOldToNew, int& alreadyCopied,
00288                                       int start, int end)
00289 {
00290   MESSAGE("copyNodes " << alreadyCopied << " " << start << " " << end << " size: " << end - start << " total: " << alreadyCopied + end - start);
00291   void *target = newPoints->GetVoidPointer(3 * alreadyCopied);
00292   void *source = this->Points->GetVoidPointer(3 * start);
00293   int nbPoints = end - start;
00294   if (nbPoints > 0)
00295     {
00296       memcpy(target, source, 3 * sizeof(double) * nbPoints);
00297       for (int j = start; j < end; j++)
00298         idNodesOldToNew[j] = alreadyCopied++; // old vtkId --> new vtkId
00299     }
00300 }
00301 
00302 void SMDS_UnstructuredGrid::copyBloc(vtkUnsignedCharArray *newTypes, std::vector<int>& idCellsOldToNew,
00303                                      std::vector<int>& idNodesOldToNew, vtkCellArray* newConnectivity,
00304                                      vtkIdTypeArray* newLocations, vtkIdType* pointsCell, int& alreadyCopied,
00305                                      int start, int end)
00306 {
00307   MESSAGE("copyBloc " << alreadyCopied << " " << start << " " << end << " size: " << end - start << " total: " << alreadyCopied + end - start);
00308   for (int j = start; j < end; j++)
00309     {
00310       newTypes->SetValue(alreadyCopied, this->Types->GetValue(j));
00311       idCellsOldToNew[j] = alreadyCopied; // old vtkId --> new vtkId
00312       vtkIdType oldLoc = this->Locations->GetValue(j);
00313       vtkIdType nbpts;
00314       vtkIdType *oldPtsCell = 0;
00315       this->Connectivity->GetCell(oldLoc, nbpts, oldPtsCell);
00316       assert(nbpts < NBMAXNODESINCELL);
00317       //MESSAGE(j << " " << alreadyCopied << " " << (int)this->Types->GetValue(j) << " " << oldLoc << " " << nbpts );
00318       for (int l = 0; l < nbpts; l++)
00319         {
00320           int oldval = oldPtsCell[l];
00321           pointsCell[l] = idNodesOldToNew[oldval];
00322           //MESSAGE("   " << oldval << " " << pointsCell[l]);
00323         }
00324       /*int newcnt = */newConnectivity->InsertNextCell(nbpts, pointsCell);
00325       int newLoc = newConnectivity->GetInsertLocation(nbpts);
00326       //MESSAGE(newcnt << " " << newLoc);
00327       newLocations->SetValue(alreadyCopied, newLoc);
00328       alreadyCopied++;
00329     }
00330 }
00331 
00332 int SMDS_UnstructuredGrid::CellIdToDownId(int vtkCellId)
00333 {
00334   if((vtkCellId < 0) || (vtkCellId >= _cellIdToDownId.size()))
00335     {
00336       //MESSAGE("SMDS_UnstructuredGrid::CellIdToDownId structure not up to date: vtkCellId="
00337       //    << vtkCellId << " max="<< _cellIdToDownId.size());
00338       return -1;
00339     }
00340   return _cellIdToDownId[vtkCellId];
00341 }
00342 
00343 void SMDS_UnstructuredGrid::setCellIdToDownId(int vtkCellId, int downId)
00344 {
00345   // ASSERT((vtkCellId >= 0) && (vtkCellId < _cellIdToDownId.size()));
00346   _cellIdToDownId[vtkCellId] = downId;
00347 }
00348 
00349 void SMDS_UnstructuredGrid::CleanDownwardConnectivity()
00350 {
00351   for (int i = 0; i < _downArray.size(); i++)
00352     {
00353       if (_downArray[i])
00354         delete _downArray[i];
00355       _downArray[i] = 0;
00356     }
00357   _cellIdToDownId.clear();
00358 }
00359 
00364 void SMDS_UnstructuredGrid::BuildDownwardConnectivity(bool withEdges)
00365 {
00366   MESSAGE("SMDS_UnstructuredGrid::BuildDownwardConnectivity");CHRONO(2);
00367   // TODO calcul partiel sans edges
00368 
00369   // --- erase previous data if any
00370 
00371   this->CleanDownwardConnectivity();
00372 
00373   // --- create SMDS_Downward structures (in _downArray vector[vtkCellType])
00374 
00375   _downArray.resize(VTK_MAXTYPE + 1, 0);
00376 
00377   _downArray[VTK_LINE]                    = new SMDS_DownEdge(this);
00378   _downArray[VTK_QUADRATIC_EDGE]          = new SMDS_DownQuadEdge(this);
00379   _downArray[VTK_TRIANGLE]                = new SMDS_DownTriangle(this);
00380   _downArray[VTK_QUADRATIC_TRIANGLE]      = new SMDS_DownQuadTriangle(this);
00381   _downArray[VTK_QUAD]                    = new SMDS_DownQuadrangle(this);
00382   _downArray[VTK_QUADRATIC_QUAD]          = new SMDS_DownQuadQuadrangle(this);
00383   _downArray[VTK_BIQUADRATIC_QUAD]        = new SMDS_DownQuadQuadrangle(this);
00384   _downArray[VTK_TETRA]                   = new SMDS_DownTetra(this);
00385   _downArray[VTK_QUADRATIC_TETRA]         = new SMDS_DownQuadTetra(this);
00386   _downArray[VTK_PYRAMID]                 = new SMDS_DownPyramid(this);
00387   _downArray[VTK_QUADRATIC_PYRAMID]       = new SMDS_DownQuadPyramid(this);
00388   _downArray[VTK_WEDGE]                   = new SMDS_DownPenta(this);
00389   _downArray[VTK_QUADRATIC_WEDGE]         = new SMDS_DownQuadPenta(this);
00390   _downArray[VTK_HEXAHEDRON]              = new SMDS_DownHexa(this);
00391   _downArray[VTK_QUADRATIC_HEXAHEDRON]    = new SMDS_DownQuadHexa(this);
00392   _downArray[VTK_TRIQUADRATIC_HEXAHEDRON] = new SMDS_DownQuadHexa(this);
00393   _downArray[VTK_HEXAGONAL_PRISM]         = new SMDS_DownPenta(this);
00394 
00395   // --- get detailed info of number of cells of each type, allocate SMDS_downward structures
00396 
00397   const SMDS_MeshInfo &meshInfo = _mesh->GetMeshInfo();
00398 
00399   int nbLinTetra  = meshInfo.NbTetras  (ORDER_LINEAR);
00400   int nbQuadTetra = meshInfo.NbTetras  (ORDER_QUADRATIC);
00401   int nbLinPyra   = meshInfo.NbPyramids(ORDER_LINEAR);
00402   int nbQuadPyra  = meshInfo.NbPyramids(ORDER_QUADRATIC);
00403   int nbLinPrism  = meshInfo.NbPrisms  (ORDER_LINEAR);
00404   int nbQuadPrism = meshInfo.NbPrisms  (ORDER_QUADRATIC);
00405   int nbLinHexa   = meshInfo.NbHexas   (ORDER_LINEAR);
00406   int nbQuadHexa  = meshInfo.NbHexas   (ORDER_QUADRATIC);
00407   int nbHexPrism  = meshInfo.NbHexPrisms();
00408 
00409   int nbLineGuess     = int((4.0 / 3.0) * nbLinTetra + 2 * nbLinPrism + 2.5 * nbLinPyra + 3 * nbLinHexa);
00410   int nbQuadEdgeGuess = int((4.0 / 3.0) * nbQuadTetra + 2 * nbQuadPrism + 2.5 * nbQuadPyra + 3 * nbQuadHexa);
00411   int nbLinTriaGuess  = 2 * nbLinTetra + nbLinPrism + 2 * nbLinPyra;
00412   int nbQuadTriaGuess = 2 * nbQuadTetra + nbQuadPrism + 2 * nbQuadPyra;
00413   int nbLinQuadGuess  = int((2.0 / 3.0) * nbLinPrism + (1.0 / 2.0) * nbLinPyra + 3 * nbLinHexa);
00414   int nbQuadQuadGuess = int((2.0 / 3.0) * nbQuadPrism + (1.0 / 2.0) * nbQuadPyra + 3 * nbQuadHexa);
00415 
00416   int GuessSize[VTK_MAXTYPE];
00417   GuessSize[VTK_LINE]                    = nbLineGuess;
00418   GuessSize[VTK_QUADRATIC_EDGE]          = nbQuadEdgeGuess;
00419   GuessSize[VTK_TRIANGLE]                = nbLinTriaGuess;
00420   GuessSize[VTK_QUADRATIC_TRIANGLE]      = nbQuadTriaGuess;
00421   GuessSize[VTK_QUAD]                    = nbLinQuadGuess;
00422   GuessSize[VTK_QUADRATIC_QUAD]          = nbQuadQuadGuess;
00423   GuessSize[VTK_BIQUADRATIC_QUAD]        = nbQuadQuadGuess;
00424   GuessSize[VTK_TETRA]                   = nbLinTetra;
00425   GuessSize[VTK_QUADRATIC_TETRA]         = nbQuadTetra;
00426   GuessSize[VTK_PYRAMID]                 = nbLinPyra;
00427   GuessSize[VTK_QUADRATIC_PYRAMID]       = nbQuadPyra;
00428   GuessSize[VTK_WEDGE]                   = nbLinPrism;
00429   GuessSize[VTK_QUADRATIC_WEDGE]         = nbQuadPrism;
00430   GuessSize[VTK_HEXAHEDRON]              = nbLinHexa;
00431   GuessSize[VTK_QUADRATIC_HEXAHEDRON]    = nbQuadHexa;
00432   GuessSize[VTK_TRIQUADRATIC_HEXAHEDRON] = nbQuadHexa;
00433   GuessSize[VTK_HEXAGONAL_PRISM]         = nbHexPrism;
00434 
00435   _downArray[VTK_LINE]                   ->allocate(nbLineGuess);
00436   _downArray[VTK_QUADRATIC_EDGE]         ->allocate(nbQuadEdgeGuess);
00437   _downArray[VTK_TRIANGLE]               ->allocate(nbLinTriaGuess);
00438   _downArray[VTK_QUADRATIC_TRIANGLE]     ->allocate(nbQuadTriaGuess);
00439   _downArray[VTK_QUAD]                   ->allocate(nbLinQuadGuess);
00440   _downArray[VTK_QUADRATIC_QUAD]         ->allocate(nbQuadQuadGuess);
00441   _downArray[VTK_BIQUADRATIC_QUAD]       ->allocate(nbQuadQuadGuess);
00442   _downArray[VTK_TETRA]                  ->allocate(nbLinTetra);
00443   _downArray[VTK_QUADRATIC_TETRA]        ->allocate(nbQuadTetra);
00444   _downArray[VTK_PYRAMID]                ->allocate(nbLinPyra);
00445   _downArray[VTK_QUADRATIC_PYRAMID]      ->allocate(nbQuadPyra);
00446   _downArray[VTK_WEDGE]                  ->allocate(nbLinPrism);
00447   _downArray[VTK_QUADRATIC_WEDGE]        ->allocate(nbQuadPrism);
00448   _downArray[VTK_HEXAHEDRON]             ->allocate(nbLinHexa);
00449   _downArray[VTK_QUADRATIC_HEXAHEDRON]   ->allocate(nbQuadHexa);
00450   _downArray[VTK_TRIQUADRATIC_HEXAHEDRON]->allocate(nbQuadHexa);
00451   _downArray[VTK_HEXAGONAL_PRISM]        ->allocate(nbHexPrism);
00452 
00453   // --- iteration on vtkUnstructuredGrid cells, only faces
00454   //     for each vtk face:
00455   //       create a downward face entry with its downward id.
00456   //       compute vtk volumes, create downward volumes entry.
00457   //       mark face in downward volumes
00458   //       mark volumes in downward face
00459 
00460   MESSAGE("--- iteration on vtkUnstructuredGrid cells, only faces");CHRONO(20);
00461   int cellSize = this->Types->GetNumberOfTuples();
00462   _cellIdToDownId.resize(cellSize, -1);
00463 
00464   for (int i = 0; i < cellSize; i++)
00465     {
00466       int vtkFaceType = this->GetCellType(i);
00467       if (SMDS_Downward::getCellDimension(vtkFaceType) == 2)
00468         {
00469           int vtkFaceId = i;
00470           //ASSERT(_downArray[vtkFaceType]);
00471           int connFaceId = _downArray[vtkFaceType]->addCell(vtkFaceId);
00472           SMDS_Down2D* downFace = static_cast<SMDS_Down2D*> (_downArray[vtkFaceType]);
00473           downFace->setTempNodes(connFaceId, vtkFaceId);
00474           int vols[2] = { -1, -1 };
00475           int nbVolumes = downFace->computeVolumeIds(vtkFaceId, vols);
00476           //MESSAGE("nbVolumes="<< nbVolumes);
00477           for (int ivol = 0; ivol < nbVolumes; ivol++)
00478             {
00479               int vtkVolId = vols[ivol];
00480               int vtkVolType = this->GetCellType(vtkVolId);
00481               //ASSERT(_downArray[vtkVolType]);
00482               int connVolId = _downArray[vtkVolType]->addCell(vtkVolId);
00483               _downArray[vtkVolType]->addDownCell(connVolId, connFaceId, vtkFaceType);
00484               _downArray[vtkFaceType]->addUpCell(connFaceId, connVolId, vtkVolType);
00485               // MESSAGE("Face " << vtkFaceId << " belongs to volume " << vtkVolId);
00486             }
00487         }
00488     }
00489 
00490   // --- iteration on vtkUnstructuredGrid cells, only volumes
00491   //     for each vtk volume:
00492   //       create downward volumes entry if not already done
00493   //       build a temporary list of faces described with their nodes
00494   //       for each face
00495   //         compute the vtk volumes containing this face
00496   //         check if the face is already listed in the volumes (comparison of ordered list of nodes)
00497   //         if not, create a downward face entry (resizing of structure required)
00498   //         (the downward faces store a temporary list of nodes to ease the comparison)
00499   //         create downward volumes entry if not already done
00500   //         mark volumes in downward face
00501   //         mark face in downward volumes
00502 
00503   CHRONOSTOP(20);
00504   MESSAGE("--- iteration on vtkUnstructuredGrid cells, only volumes");CHRONO(21);
00505 
00506   for (int i = 0; i < cellSize; i++)
00507     {
00508       int vtkType = this->GetCellType(i);
00509       if (SMDS_Downward::getCellDimension(vtkType) == 3)
00510         {
00511           //CHRONO(31);
00512           int vtkVolId = i;
00513           // MESSAGE("vtk volume " << vtkVolId);
00514           //ASSERT(_downArray[vtkType]);
00515           /*int connVolId = */_downArray[vtkType]->addCell(vtkVolId);
00516 
00517           // --- find all the faces of the volume, describe the faces by their nodes
00518 
00519           SMDS_Down3D* downVol = static_cast<SMDS_Down3D*> (_downArray[vtkType]);
00520           ListElemByNodesType facesWithNodes;
00521           downVol->computeFacesWithNodes(vtkVolId, facesWithNodes);
00522           // MESSAGE("vtk volume " << vtkVolId << " contains " << facesWithNodes.nbElems << " faces");
00523           //CHRONOSTOP(31);
00524           for (int iface = 0; iface < facesWithNodes.nbElems; iface++)
00525             {
00526               // --- find the volumes containing the face
00527 
00528               //CHRONO(32);
00529               int vtkFaceType = facesWithNodes.elems[iface].vtkType;
00530               SMDS_Down2D* downFace = static_cast<SMDS_Down2D*> (_downArray[vtkFaceType]);
00531               int vols[2] = { -1, -1 };
00532               int *nodes = &facesWithNodes.elems[iface].nodeIds[0];
00533               int lg = facesWithNodes.elems[iface].nbNodes;
00534               int nbVolumes = downFace->computeVolumeIdsFromNodesFace(nodes, lg, vols);
00535               // MESSAGE("vtk volume " << vtkVolId << " face " << iface << " belongs to " << nbVolumes << " volumes");
00536 
00537               // --- check if face is registered in the volumes
00538               //CHRONOSTOP(32);
00539 
00540               //CHRONO(33);
00541               int connFaceId = -1;
00542               for (int ivol = 0; ivol < nbVolumes; ivol++)
00543                 {
00544                   int vtkVolId2 = vols[ivol];
00545                   int vtkVolType = this->GetCellType(vtkVolId2);
00546                   //ASSERT(_downArray[vtkVolType]);
00547                   int connVolId2 = _downArray[vtkVolType]->addCell(vtkVolId2);
00548                   SMDS_Down3D* downVol2 = static_cast<SMDS_Down3D*> (_downArray[vtkVolType]);
00549                   connFaceId = downVol2->FindFaceByNodes(connVolId2, facesWithNodes.elems[iface]);
00550                   if (connFaceId >= 0)
00551                     break; // --- face already created
00552                 }//CHRONOSTOP(33);
00553 
00554               // --- if face is not registered in the volumes, create face
00555 
00556               //CHRONO(34);
00557               if (connFaceId < 0)
00558                 {
00559                   connFaceId = _downArray[vtkFaceType]->addCell();
00560                   downFace->setTempNodes(connFaceId, facesWithNodes.elems[iface]);
00561                 }//CHRONOSTOP(34);
00562 
00563               // --- mark volumes in downward face and mark face in downward volumes
00564 
00565               //CHRONO(35);
00566               for (int ivol = 0; ivol < nbVolumes; ivol++)
00567                 {
00568                   int vtkVolId2 = vols[ivol];
00569                   int vtkVolType = this->GetCellType(vtkVolId2);
00570                   //ASSERT(_downArray[vtkVolType]);
00571                   int connVolId2 = _downArray[vtkVolType]->addCell(vtkVolId2);
00572                   _downArray[vtkVolType]->addDownCell(connVolId2, connFaceId, vtkFaceType);
00573                   _downArray[vtkFaceType]->addUpCell(connFaceId, connVolId2, vtkVolType);
00574                   // MESSAGE(" From volume " << vtkVolId << " face " << connFaceId << " belongs to volume " << vtkVolId2);
00575                 }//CHRONOSTOP(35);
00576             }
00577         }
00578     }
00579 
00580   // --- iteration on vtkUnstructuredGrid cells, only edges
00581   //     for each vtk edge:
00582   //       create downward edge entry
00583   //       store the nodes id's in downward edge (redundant with vtkUnstructuredGrid)
00584   //       find downward faces
00585   //       (from vtk faces or volumes, get downward faces, they have a temporary list of nodes)
00586   //       mark edge in downward faces
00587   //       mark faces in downward edge
00588 
00589   CHRONOSTOP(21);
00590   MESSAGE("--- iteration on vtkUnstructuredGrid cells, only edges");CHRONO(22);
00591 
00592   for (int i = 0; i < cellSize; i++)
00593     {
00594       int vtkEdgeType = this->GetCellType(i);
00595       if (SMDS_Downward::getCellDimension(vtkEdgeType) == 1)
00596         {
00597           int vtkEdgeId = i;
00598           //ASSERT(_downArray[vtkEdgeType]);
00599           int connEdgeId = _downArray[vtkEdgeType]->addCell(vtkEdgeId);
00600           SMDS_Down1D* downEdge = static_cast<SMDS_Down1D*> (_downArray[vtkEdgeType]);
00601           downEdge->setNodes(connEdgeId, vtkEdgeId);
00602           vector<int> vtkIds;
00603           int nbVtkCells = downEdge->computeVtkCells(connEdgeId, vtkIds);
00604           int downFaces[1000];
00605           unsigned char downTypes[1000];
00606           int nbDownFaces = downEdge->computeFaces(connEdgeId, &vtkIds[0], nbVtkCells, downFaces, downTypes);
00607           for (int n = 0; n < nbDownFaces; n++)
00608             {
00609               _downArray[downTypes[n]]->addDownCell(downFaces[n], connEdgeId, vtkEdgeType);
00610               _downArray[vtkEdgeType]->addUpCell(connEdgeId, downFaces[n], downTypes[n]);
00611             }
00612         }
00613     }
00614 
00615   // --- iteration on downward faces (they are all listed now)
00616   //     for each downward face:
00617   //       build a temporary list of edges with their ordered list of nodes
00618   //       for each edge:
00619   //         find all the vtk cells containing this edge
00620   //         then identify all the downward faces containing the edge, from the vtk cells
00621   //         check if the edge is already listed in the faces (comparison of ordered list of nodes)
00622   //         if not, create a downward edge entry with the node id's
00623   //         mark edge in downward faces
00624   //         mark downward faces in edge (size of list unknown, to be allocated)
00625 
00626   CHRONOSTOP(22);CHRONO(23);
00627 
00628   for (int vtkFaceType = 0; vtkFaceType < VTK_QUADRATIC_PYRAMID; vtkFaceType++)
00629     {
00630       if (SMDS_Downward::getCellDimension(vtkFaceType) != 2)
00631         continue;
00632 
00633       // --- find all the edges of the face, describe the edges by their nodes
00634 
00635       SMDS_Down2D* downFace = static_cast<SMDS_Down2D*> (_downArray[vtkFaceType]);
00636       int maxId = downFace->getMaxId();
00637       for (int faceId = 0; faceId < maxId; faceId++)
00638         {
00639           //CHRONO(40);
00640           ListElemByNodesType edgesWithNodes;
00641           downFace->computeEdgesWithNodes(faceId, edgesWithNodes);
00642           // MESSAGE("downward face type " << vtkFaceType << " num " << faceId << " contains " << edgesWithNodes.nbElems << " edges");
00643 
00644           //CHRONOSTOP(40);
00645           for (int iedge = 0; iedge < edgesWithNodes.nbElems; iedge++)
00646             {
00647 
00648               // --- check if the edge is already registered by exploration of the faces
00649 
00650               //CHRONO(41);
00651               vector<int> vtkIds;
00652               unsigned char vtkEdgeType = edgesWithNodes.elems[iedge].vtkType;
00653               int *pts = &edgesWithNodes.elems[iedge].nodeIds[0];
00654               SMDS_Down1D* downEdge = static_cast<SMDS_Down1D*> (_downArray[vtkEdgeType]);
00655               int nbVtkCells = downEdge->computeVtkCells(pts, vtkIds);
00656               //CHRONOSTOP(41);CHRONO(42);
00657               int downFaces[1000];
00658               unsigned char downTypes[1000];
00659               int nbDownFaces = downEdge->computeFaces(pts, &vtkIds[0], nbVtkCells, downFaces, downTypes);
00660               //CHRONOSTOP(42);
00661 
00662               //CHRONO(43);
00663               int connEdgeId = -1;
00664               for (int idf = 0; idf < nbDownFaces; idf++)
00665                 {
00666                   int faceId2 = downFaces[idf];
00667                   int faceType = downTypes[idf];
00668                   //ASSERT(_downArray[faceType]);
00669                   SMDS_Down2D* downFace2 = static_cast<SMDS_Down2D*> (_downArray[faceType]);
00670                   connEdgeId = downFace2->FindEdgeByNodes(faceId2, edgesWithNodes.elems[iedge]);
00671                   if (connEdgeId >= 0)
00672                     break; // --- edge already created
00673                 }//CHRONOSTOP(43);
00674 
00675               // --- if edge is not registered in the faces, create edge
00676 
00677               if (connEdgeId < 0)
00678                 {
00679                   //CHRONO(44);
00680                   connEdgeId = _downArray[vtkEdgeType]->addCell();
00681                   downEdge->setNodes(connEdgeId, edgesWithNodes.elems[iedge].nodeIds);
00682                   //CHRONOSTOP(44);
00683                 }
00684 
00685               // --- mark faces in downward edge and mark edge in downward faces
00686 
00687               //CHRONO(45);
00688               for (int idf = 0; idf < nbDownFaces; idf++)
00689                 {
00690                   int faceId2 = downFaces[idf];
00691                   int faceType = downTypes[idf];
00692                   //ASSERT(_downArray[faceType]);
00693                   //SMDS_Down2D* downFace2 = static_cast<SMDS_Down2D*> (_downArray[faceType]);
00694                   _downArray[vtkEdgeType]->addUpCell(connEdgeId, faceId2, faceType);
00695                   _downArray[faceType]->addDownCell(faceId2, connEdgeId, vtkEdgeType);
00696                   // MESSAGE(" From face t:" << vtkFaceType << " " << faceId <<
00697                   //  " edge " << connEdgeId << " belongs to face t:" << faceType << " " << faceId2);
00698                 }//CHRONOSTOP(45);
00699             }
00700         }
00701     }
00702 
00703   CHRONOSTOP(23);CHRONO(24);
00704 
00705   // compact downward connectivity structure: adjust downward arrays size, replace vector<vector int>> by a single vector<int>
00706   // 3D first then 2D and last 1D to release memory before edge upCells reorganization, (temporary memory use)
00707 
00708   for (int vtkType = VTK_QUADRATIC_PYRAMID; vtkType >= 0; vtkType--)
00709     {
00710       if (SMDS_Downward *down = _downArray[vtkType])
00711         {
00712           down->compactStorage();
00713         }
00714     }
00715 
00716   // --- Statistics
00717 
00718   for (int vtkType = 0; vtkType <= VTK_QUADRATIC_PYRAMID; vtkType++)
00719     {
00720       if (SMDS_Downward *down = _downArray[vtkType])
00721         {
00722           if (down->getMaxId())
00723             {
00724               MESSAGE("Cells of Type " << vtkType << " : number of entities, est: "
00725                   << GuessSize[vtkType] << " real: " << down->getMaxId());
00726             }
00727         }
00728     }CHRONOSTOP(24);CHRONOSTOP(2);
00729   counters::stats();
00730 }
00731 
00742 int SMDS_UnstructuredGrid::GetNeighbors(int* neighborsVtkIds, int* downIds, unsigned char* downTypes, int vtkId)
00743 {
00744   int vtkType = this->GetCellType(vtkId);
00745   int cellDim = SMDS_Downward::getCellDimension(vtkType);
00746   if (cellDim <2)
00747     return 0; // TODO voisins des edges = edges connectees
00748   int cellId = this->_cellIdToDownId[vtkId];
00749 
00750   int nbCells = _downArray[vtkType]->getNumberOfDownCells(cellId);
00751   const int *downCells = _downArray[vtkType]->getDownCells(cellId);
00752   const unsigned char* downTyp = _downArray[vtkType]->getDownTypes(cellId);
00753 
00754   // --- iteration on faces of the 3D cell (or edges on the 2D cell).
00755 
00756   int nb = 0;
00757   for (int i = 0; i < nbCells; i++)
00758     {
00759       int downId = downCells[i];
00760       int cellType = downTyp[i];
00761       int nbUp = _downArray[cellType]->getNumberOfUpCells(downId);
00762       const int *upCells = _downArray[cellType]->getUpCells(downId);
00763       const unsigned char* upTypes = _downArray[cellType]->getUpTypes(downId);
00764 
00765       // ---for a volume, max 2 upCells, one is this cell, the other is a neighbor
00766       //    for a face, number of neighbors (connected faces) not known
00767 
00768       for (int j = 0; j < nbUp; j++)
00769         {
00770           if ((upCells[j] == cellId) && (upTypes[j] == vtkType))
00771             continue;
00772           int vtkNeighbor = _downArray[upTypes[j]]->getVtkCellId(upCells[j]);
00773           neighborsVtkIds[nb] = vtkNeighbor;
00774           downIds[nb] = downId;
00775           downTypes[nb] = cellType;
00776           nb++;
00777         }
00778       if (nb >= NBMAXNEIGHBORS)
00779         assert(0);
00780     }
00781   return nb;
00782 }
00783 
00789 int SMDS_UnstructuredGrid::GetParentVolumes(int* volVtkIds, int vtkId)
00790 {
00791   int vtkType = this->GetCellType(vtkId);
00792   int dim = SMDS_Downward::getCellDimension(vtkType);
00793   int nbFaces = 0;
00794   unsigned char cellTypes[1000];
00795   int downCellId[1000];
00796   if (dim == 1)
00797     {
00798       int downId = this->CellIdToDownId(vtkId);
00799       if (downId < 0)
00800         {
00801           MESSAGE("Downward structure not up to date: new edge not taken into account");
00802           return 0;
00803         }
00804       nbFaces = _downArray[vtkType]->getNumberOfUpCells(downId);
00805       const int *upCells = _downArray[vtkType]->getUpCells(downId);
00806       const unsigned char* upTypes = _downArray[vtkType]->getUpTypes(downId);
00807       for (int i=0; i< nbFaces; i++)
00808         {
00809           cellTypes[i] = upTypes[i];
00810           downCellId[i] = upCells[i];
00811         }
00812     }
00813   else if (dim == 2)
00814     {
00815       nbFaces = 1;
00816       cellTypes[0] = this->GetCellType(vtkId);
00817       int downId = this->CellIdToDownId(vtkId);
00818       if (downId < 0)
00819         {
00820           MESSAGE("Downward structure not up to date: new face not taken into account");
00821           return 0;
00822         }
00823       downCellId[0] = downId;
00824     }
00825 
00826   int nbvol =0;
00827   for (int i=0; i<nbFaces; i++)
00828     {
00829       int vtkTypeFace = cellTypes[i];
00830       int downId = downCellId[i];
00831       int nv = _downArray[vtkTypeFace]->getNumberOfUpCells(downId);
00832       const int *upCells = _downArray[vtkTypeFace]->getUpCells(downId);
00833       const unsigned char* upTypes = _downArray[vtkTypeFace]->getUpTypes(downId);
00834        for (int j=0; j<nv; j++)
00835         {
00836           int vtkVolId = _downArray[upTypes[j]]->getVtkCellId(upCells[j]);
00837           if (vtkVolId >= 0)
00838             volVtkIds[nbvol++] = vtkVolId;
00839         }
00840     }
00841   return nbvol;
00842 }
00843 
00850 int SMDS_UnstructuredGrid::GetParentVolumes(int* volVtkIds, int downId, unsigned char downType)
00851 {
00852   int vtkType = downType;
00853   int dim = SMDS_Downward::getCellDimension(vtkType);
00854   int nbFaces = 0;
00855   unsigned char cellTypes[1000];
00856   int downCellId[1000];
00857   if (dim == 1)
00858     {
00859       nbFaces = _downArray[vtkType]->getNumberOfUpCells(downId);
00860       const int *upCells = _downArray[vtkType]->getUpCells(downId);
00861       const unsigned char* upTypes = _downArray[vtkType]->getUpTypes(downId);
00862       for (int i=0; i< nbFaces; i++)
00863         {
00864           cellTypes[i] = upTypes[i];
00865           downCellId[i] = upCells[i];
00866         }
00867     }
00868   else if (dim == 2)
00869     {
00870       nbFaces = 1;
00871       cellTypes[0] = vtkType;
00872       downCellId[0] = downId;
00873     }
00874 
00875   int nbvol =0;
00876   for (int i=0; i<nbFaces; i++)
00877     {
00878       int vtkTypeFace = cellTypes[i];
00879       int downId = downCellId[i];
00880       int nv = _downArray[vtkTypeFace]->getNumberOfUpCells(downId);
00881       const int *upCells = _downArray[vtkTypeFace]->getUpCells(downId);
00882       const unsigned char* upTypes = _downArray[vtkTypeFace]->getUpTypes(downId);
00883        for (int j=0; j<nv; j++)
00884         {
00885           int vtkVolId = _downArray[upTypes[j]]->getVtkCellId(upCells[j]);
00886           if (vtkVolId >= 0)
00887             volVtkIds[nbvol++] = vtkVolId;
00888         }
00889     }
00890   return nbvol;
00891 }
00892 
00899 void SMDS_UnstructuredGrid::GetNodeIds(std::set<int>& nodeSet, int downId, unsigned char downType)
00900 {
00901   _downArray[downType]->getNodeIds(downId, nodeSet);
00902 }
00903 
00909 void SMDS_UnstructuredGrid::ModifyCellNodes(int vtkVolId, std::map<int, int> localClonedNodeIds)
00910 {
00911   vtkIdType npts = 0;
00912   vtkIdType *pts; // will refer to the point id's of the face
00913   this->GetCellPoints(vtkVolId, npts, pts);
00914   for (int i = 0; i < npts; i++)
00915     {
00916       if (localClonedNodeIds.count(pts[i]))
00917         {
00918           vtkIdType oldpt = pts[i];
00919           pts[i] = localClonedNodeIds[oldpt];
00920           //MESSAGE(oldpt << " --> " << pts[i]);
00921           //this->RemoveReferenceToCell(oldpt, vtkVolId);
00922           //this->AddReferenceToCell(pts[i], vtkVolId);
00923         }
00924     }
00925 }
00926 
00932 int SMDS_UnstructuredGrid::getOrderedNodesOfFace(int vtkVolId, int& dim, std::vector<vtkIdType>& orderedNodes)
00933 {
00934   int vtkType = this->GetCellType(vtkVolId);
00935   dim = SMDS_Downward::getCellDimension(vtkType);
00936   if (dim == 3)
00937     {
00938       SMDS_Down3D *downvol = static_cast<SMDS_Down3D*> (_downArray[vtkType]);
00939       int downVolId = this->_cellIdToDownId[vtkVolId];
00940       downvol->getOrderedNodesOfFace(downVolId, orderedNodes);
00941     }
00942   // else nothing to do;
00943   return orderedNodes.size();
00944 }
00945 
00946 void SMDS_UnstructuredGrid::BuildLinks()
00947 {
00948   // Remove the old links if they are already built
00949   if (this->Links)
00950     {
00951     this->Links->UnRegister(this);
00952     }
00953 
00954   this->Links = SMDS_CellLinks::New();
00955   this->Links->Allocate(this->GetNumberOfPoints());
00956   this->Links->Register(this);
00957   this->Links->BuildLinks(this, this->Connectivity);
00958   this->Links->Delete();
00959 }
00960 
00975 SMDS_MeshCell* SMDS_UnstructuredGrid::extrudeVolumeFromFace(int vtkVolId,
00976                                                   int domain1,
00977                                                   int domain2,
00978                                                   std::set<int>& originalNodes,
00979                                                   std::map<int, std::map<int, int> >& nodeDomains,
00980                                                   std::map<int, std::map<long, int> >& nodeQuadDomains)
00981 {
00982   //MESSAGE("extrudeVolumeFromFace " << vtkVolId);
00983   vector<vtkIdType> orderedOriginals;
00984   orderedOriginals.clear();
00985   set<int>::const_iterator it = originalNodes.begin();
00986   for (; it != originalNodes.end(); ++it)
00987     orderedOriginals.push_back(*it);
00988 
00989   int dim = 0;
00990   int nbNodes = this->getOrderedNodesOfFace(vtkVolId, dim, orderedOriginals);
00991   vector<vtkIdType> orderedNodes;
00992 
00993   bool isQuadratic = false;
00994   switch (orderedOriginals.size())
00995   {
00996     case 3:
00997       if (dim == 2)
00998         isQuadratic = true;
00999       break;
01000     case 6:
01001     case 8:
01002       isQuadratic = true;
01003       break;
01004     default:
01005       isQuadratic = false;
01006       break;
01007   }
01008 
01009   if (isQuadratic)
01010     {
01011       long dom1 = domain1;
01012       long dom2 = domain2;
01013       long dom1_2; // for nodeQuadDomains
01014       if (domain1 < domain2)
01015         dom1_2 = dom1 + INT_MAX * dom2;
01016       else
01017         dom1_2 = dom2 + INT_MAX * dom1;
01018       //cerr << "dom1=" << dom1 << " dom2=" << dom2 << " dom1_2=" << dom1_2 << endl;
01019       int ima = orderedOriginals.size();
01020       int mid = orderedOriginals.size() / 2;
01021       //cerr << "ima=" << ima << " mid=" << mid << endl;
01022       for (int i = 0; i < mid; i++)
01023         orderedNodes.push_back(nodeDomains[orderedOriginals[i]][domain1]);
01024       for (int i = 0; i < mid; i++)
01025         orderedNodes.push_back(nodeDomains[orderedOriginals[i]][domain2]);
01026       for (int i = mid; i < ima; i++)
01027         orderedNodes.push_back(nodeDomains[orderedOriginals[i]][domain1]);
01028       for (int i = mid; i < ima; i++)
01029         orderedNodes.push_back(nodeDomains[orderedOriginals[i]][domain2]);
01030       for (int i = 0; i < mid; i++)
01031         {
01032           int oldId = orderedOriginals[i];
01033           int newId;
01034           if (nodeQuadDomains.count(oldId) && nodeQuadDomains[oldId].count(dom1_2))
01035             newId = nodeQuadDomains[oldId][dom1_2];
01036           else
01037             {
01038               double *coords = this->GetPoint(oldId);
01039               SMDS_MeshNode *newNode = _mesh->AddNode(coords[0], coords[1], coords[2]);
01040               newId = newNode->getVtkId();
01041               std::map<long, int> emptyMap;
01042               nodeQuadDomains[oldId] = emptyMap;
01043               nodeQuadDomains[oldId][dom1_2] = newId;
01044             }
01045           orderedNodes.push_back(newId);
01046         }
01047     }
01048   else
01049     {
01050       for (int i = 0; i < nbNodes; i++)
01051         orderedNodes.push_back(nodeDomains[orderedOriginals[i]][domain1]);
01052       if (dim == 3)
01053         for (int i = 0; i < nbNodes; i++)
01054           orderedNodes.push_back(nodeDomains[orderedOriginals[i]][domain2]);
01055       else
01056         for (int i = nbNodes-1; i >= 0; i--)
01057           orderedNodes.push_back(nodeDomains[orderedOriginals[i]][domain2]);
01058 
01059     }
01060 
01061   if (dim == 3)
01062     {
01063       SMDS_MeshVolume *vol = _mesh->AddVolumeFromVtkIds(orderedNodes);
01064       return vol;
01065     }
01066   else if (dim == 2)
01067     {
01068       SMDS_MeshFace *face = _mesh->AddFaceFromVtkIds(orderedNodes);
01069       return face;
01070     }
01071 
01072   // TODO update sub-shape list of elements and nodes
01073   return 0;
01074 }