Back to index

salome-smesh  6.5.0
DriverCGNS_Read.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 // File      : DriverCGNS_Read.cxx
00023 // Created   : Thu Jun 30 10:33:31 2011
00024 // Author    : Edward AGAPOV (eap)
00025 
00026 #include "DriverCGNS_Read.hxx"
00027 
00028 #include "SMDS_MeshNode.hxx"
00029 #include "SMESHDS_Group.hxx"
00030 #include "SMESHDS_Mesh.hxx"
00031 #include "SMESH_Comment.hxx"
00032 
00033 #include <gp_XYZ.hxx>
00034 
00035 #include <cgnslib.h>
00036 
00037 #include <map>
00038 
00039 #if CGNS_VERSION < 3100
00040 # define cgsize_t int
00041 #endif
00042 
00043 #define NB_ZONE_SIZE_VAL 9
00044 #define CGNS_NAME_SIZE 33
00045 #define CGNS_STRUCT_RANGE_SZ 6
00046 
00047 using namespace std;
00048 
00049 namespace
00050 {
00051   //================================================================================
00055   struct TZoneData
00056   {
00057     int                    _id;
00058     int                    _nodeIdShift; // nb nodes in previously read zones
00059     int                    _elemIdShift; // nb faces in previously read zones
00060     int                    _nbNodes, _nbElems;
00061     int                    _meshDim;
00062     int                    _sizeX, _sizeY, _sizeZ, _nbCells; // structured
00063     cgsize_t               _sizes[NB_ZONE_SIZE_VAL];
00064     CGNS_ENUMT(ZoneType_t) _type;
00065     map< int, int >        _nodeReplacementMap;/* key:   id of node to replace (in this zone),
00066                                                   value: id of node to replace by (in another zone)
00067                                                   id values include _nodeIdShift of the zones */
00068     void SetSizeAndDim( cgsize_t* sizes, int meshDim )
00069     {
00070       _meshDim = meshDim;
00071       memcpy( _sizes, sizes, NB_ZONE_SIZE_VAL*sizeof(cgsize_t));
00072       _sizeX = _sizes[0];
00073       _sizeY = _meshDim > 1 ? _sizes[1] : 0;
00074       _sizeZ = _meshDim > 2 ? _sizes[2] : 0;
00075       _nbCells = (_sizeX - 1) * ( _meshDim > 1 ? _sizeY : 1 ) * ( _meshDim > 2 ? _sizeZ : 1 );
00076     }
00077     bool IsStructured() const { return ( _type == CGNS_ENUMV( Structured )); }
00078     int IndexSize() const { return IsStructured() ? _meshDim : 1; }
00079     string ReadZonesConnection(int file, int base, const map< string, TZoneData >& zonesByName);
00080     void ReplaceNodes( cgsize_t* ids, int nbIds, int idShift = 0 ) const;
00081 
00082     // Methods for a structured zone
00083 
00084     int NodeID( int i, int j, int k = 1 ) const
00085     {
00086       return _nodeIdShift + (k-1)*_sizeX*_sizeY + (j-1)*_sizeX + i;
00087     }
00088     int NodeID( const gp_XYZ& ijk ) const
00089     {
00090       return NodeID( int(ijk.X()), int(ijk.Y()), int(ijk.Z()));
00091     }
00092     void CellNodes( int i, int j, int k, cgsize_t* ids ) const
00093     {
00094       ids[0] = NodeID( i  , j  , k  );
00095       ids[1] = NodeID( i  , j+1, k  );
00096       ids[2] = NodeID( i+1, j+1, k  );
00097       ids[3] = NodeID( i+1, j  , k  );
00098       ids[4] = NodeID( i  , j  , k+1);
00099       ids[5] = NodeID( i  , j+1, k+1);
00100       ids[6] = NodeID( i+1, j+1, k+1);
00101       ids[7] = NodeID( i+1, j  , k+1);
00102     }
00103     void CellNodes( int i, int j, cgsize_t* ids ) const
00104     {
00105       ids[0] = NodeID( i  , j   );
00106       ids[1] = NodeID( i  , j+1 );
00107       ids[2] = NodeID( i+1, j+1 );
00108       ids[3] = NodeID( i+1, j   );
00109     }
00110     void IFaceNodes( int i, int j, int k, cgsize_t* ids ) const // face perpendiculaire to X (3D)
00111     {
00112       ids[0] = NodeID( i, j, k );
00113       ids[1] = ids[0] + _sizeX*( i==_sizeX ? 1 : _sizeY );
00114       ids[2] = ids[0] + _sizeX*( _sizeY + 1 );
00115       ids[3] = ids[0] + _sizeX*( i==_sizeX ? _sizeY : 1 );
00116     }
00117     void JFaceNodes( int i, int j, int k, cgsize_t* ids ) const
00118     {
00119       ids[0] = NodeID( i, j, k );
00120       ids[1] = ids[0] + ( j==_sizeY ? _sizeX*_sizeY : 1);
00121       ids[2] = ids[0] + _sizeX*_sizeY + 1;
00122       ids[3] = ids[0] + ( j==_sizeY ? 1 : _sizeX*_sizeY);
00123     }
00124     void KFaceNodes( int i, int j, int k, cgsize_t* ids ) const
00125     {
00126       ids[0] = NodeID( i, j, k );
00127       ids[1] = ids[0] + ( k==_sizeZ ? 1 : _sizeX);
00128       ids[2] = ids[0] + _sizeX + 1;
00129       ids[3] = ids[0] + ( k==_sizeZ ? _sizeX : 1);
00130     }
00131     void IEdgeNodes( int i, int j, int k, cgsize_t* ids ) const // edge perpendiculaire to X (2D)
00132     {
00133       ids[0] = NodeID( i, j, 0 );
00134       ids[1] = ids[0] + _sizeX;
00135     }
00136     void JEdgeNodes( int i, int j, int k, cgsize_t* ids ) const
00137     {
00138       ids[0] = NodeID( i, j, 0 );
00139       ids[1] = ids[0] + 1;
00140     }
00141 #define gpXYZ2IJK(METHOD)                                     \
00142     void METHOD( const gp_XYZ& ijk, cgsize_t* ids ) const {        \
00143       METHOD( int(ijk.X()), int(ijk.Y()), int(ijk.Z()), ids); \
00144     }
00145     gpXYZ2IJK( IFaceNodes )
00146     gpXYZ2IJK( JFaceNodes )
00147     gpXYZ2IJK( KFaceNodes )
00148     gpXYZ2IJK( IEdgeNodes )
00149     gpXYZ2IJK( JEdgeNodes )
00150   };
00151 
00152   //================================================================================
00157   class TPointRangeIterator
00158   {
00159     int _beg[3], _end[3], _cur[3], _dir[3], _dim;
00160     bool _more;
00161   public:
00162     TPointRangeIterator( const cgsize_t* range, int dim ):_dim(dim)
00163     {
00164       _more = false;
00165       for ( int i = 0; i < dim; ++i )
00166       {
00167         _beg[i] = range[i];
00168         _end[i] = range[i+dim];
00169         _dir[i] = _end[i] < _beg[i] ? -1 : 1;
00170         _end[i] += _dir[i];
00171         _cur[i] = _beg[i];
00172         if ( _end[i] - _beg[i] )
00173           _more = true;
00174       }
00175 //       for ( int i = dim; i < 3; ++i )
00176 //         _cur[i] = _beg[i] = _end[i] = _dir[i] = 0;
00177     }
00178     bool More() const
00179     {
00180       return _more;
00181     }
00182     gp_XYZ Next()
00183     {
00184       gp_XYZ res( _cur[0], _cur[1], _cur[2] );
00185       for ( int i = 0; i < _dim; ++i )
00186       {
00187         _cur[i] += _dir[i];
00188         if ( _cur[i]*_dir[i] < _end[i]*_dir[i] )
00189           break;
00190         if ( i+1 < _dim )
00191           _cur[i] = _beg[i];
00192         else
00193           _more = false;
00194       }
00195       return res;
00196     }
00197     size_t Size()  const
00198     {
00199       size_t size = 1;
00200       for ( int i = 0; i < _dim; ++i )
00201         size *= _dir[i]*(_end[i]-_beg[i]);
00202       return size;
00203     }
00204     gp_XYZ Begin() const { return gp_XYZ( _beg[0], _beg[1], _beg[2] ); }
00205     //gp_XYZ End() const { return gp_XYZ( _end[0]-1, _end[1]-1, _end[2]-1 ); }
00206   };
00207 
00208   //================================================================================
00219   //================================================================================
00220 
00221   string TZoneData::ReadZonesConnection( int                             file,
00222                                          int                             base,
00223                                          const map< string, TZoneData >& zonesByName)
00224   {
00225     string error;
00226 
00227     char connectName[ CGNS_NAME_SIZE ], donorName [ CGNS_NAME_SIZE ];
00228 
00229     // ----------------------------
00230     // read zone 1 to 1 interfaces
00231     // ----------------------------
00232     if ( IsStructured() )
00233     {
00234       int nb1to1 = 0;
00235       if ( cg_n1to1 ( file, base, _id, &nb1to1) == CG_OK )
00236       {
00237         cgsize_t range[CGNS_STRUCT_RANGE_SZ], donorRange[CGNS_STRUCT_RANGE_SZ];
00238         int transform[3] = {0,0,0};
00239 
00240         for ( int I = 1; I <= nb1to1; ++I )
00241         {
00242           if ( cg_1to1_read(file, base, _id, I, connectName,
00243                             donorName, range, donorRange, transform) == CG_OK )
00244           {
00245             map< string, TZoneData>::const_iterator n_z = zonesByName.find( donorName );
00246             if ( n_z == zonesByName.end() )
00247               continue; // donor zone not yet read
00248             const TZoneData& zone2 = n_z->second;
00249 
00250             // set up matrix to transform ijk of the zone to ijk of the zone2
00251             gp_Mat T;
00252             for ( int i = 0; i < _meshDim; ++i )
00253               if ( transform[i] )
00254               {
00255                 int row = Abs(transform[i]);
00256                 int col = i+1;
00257                 int val = transform[i] > 0 ? +1 : -1;
00258                 T( row, col ) = val;
00259               }
00260 
00261             // fill nodeReplacementMap
00262             TPointRangeIterator rangeIt1( range, _meshDim );
00263             TPointRangeIterator rangeIt2( donorRange, _meshDim );
00264             gp_XYZ begin1 = rangeIt1.Begin(), begin2 = rangeIt2.Begin(), index1, index2;
00265             if ( &zone2 == this )
00266             {
00267               // not to read twice the same interface with self
00268               TPointRangeIterator rangeIt1bis( range, _meshDim );
00269               if ( rangeIt1bis.More() )
00270               {
00271                 index1 = rangeIt1bis.Next();
00272                 index2 = T * ( index1 - begin1 ) + begin2;
00273                 int node1 = NodeID( index1 );
00274                 int node2 = zone2.NodeID( index2 );
00275                 if ( _nodeReplacementMap.count( node2 ) &&
00276                      _nodeReplacementMap[ node2 ] == node1 )
00277                   continue; // this interface already read
00278               }
00279             }
00280             while ( rangeIt1.More() )
00281             {
00282               index1 = rangeIt1.Next();
00283               index2 = T * ( index1 - begin1 ) + begin2;
00284               int node1 = NodeID( index1 );
00285               int node2 = zone2.NodeID( index2 );
00286               _nodeReplacementMap.insert( make_pair( node1, node2 ));
00287             }
00288           }
00289           else
00290           {
00291             error = cg_get_error();
00292           }
00293         }
00294       }
00295       else
00296       {
00297         error = cg_get_error();
00298       }
00299     }
00300 
00301     // ---------------------------------
00302     // read general zone connectivities
00303     // ---------------------------------
00304     int nbConn = 0;
00305     if ( cg_nconns( file, base, _id, &nbConn) == CG_OK )
00306     {
00307       cgsize_t nb, donorNb;
00308       CGNS_ENUMT(GridLocation_t) location;
00309       CGNS_ENUMT(GridConnectivityType_t) connectType;
00310       CGNS_ENUMT(PointSetType_t) ptype, donorPtype;
00311       CGNS_ENUMT(ZoneType_t) donorZonetype;
00312       CGNS_ENUMT(DataType_t) donorDatatype;
00313 
00314       for ( int I = 1; I <= nbConn; ++I )
00315       {
00316         if ( cg_conn_info(file, base, _id, I, connectName, &location, &connectType,
00317                           &ptype, &nb, donorName, &donorZonetype, &donorPtype,
00318                           &donorDatatype, &donorNb ) == CG_OK )
00319         {
00320           if ( location != CGNS_ENUMV( Vertex ))
00321             continue; // we do not support cell-to-cell connectivity
00322           if ( ptype != CGNS_ENUMV( PointList ) &&
00323                ptype != CGNS_ENUMV( PointRange ))
00324             continue;
00325           if ( donorPtype != CGNS_ENUMV( PointList ) &&
00326                donorPtype != CGNS_ENUMV( PointRange ))
00327             continue;
00328           
00329           map< string, TZoneData>::const_iterator n_z = zonesByName.find( donorName );
00330           if ( n_z == zonesByName.end() )
00331             continue; // donor zone not yet read
00332           const TZoneData& zone2 = n_z->second;
00333 
00334           vector< cgsize_t > ids( nb * IndexSize() );
00335           vector< cgsize_t > donorIds( donorNb * zone2.IndexSize() );
00336           if (cg_conn_read ( file, base, _id, I,
00337                              &ids[0], CGNS_ENUMV(Integer), &donorIds[0]) == CG_OK )
00338           {
00339             for ( int isThisZone = 0; isThisZone < 2; ++isThisZone )
00340             {
00341               const TZoneData&            zone = isThisZone ? *this : zone2;
00342               CGNS_ENUMT(PointSetType_t) type = isThisZone ? ptype : donorPtype;
00343               vector< cgsize_t >&      points = isThisZone ? ids : donorIds;
00344               if ( type == CGNS_ENUMV( PointRange ))
00345               {
00346                 TPointRangeIterator rangeIt( &points[0], zone._meshDim );
00347                 points.clear();
00348                 while ( rangeIt.More() )
00349                   points.push_back ( NodeID( rangeIt.Next() ));
00350               }
00351               else if ( zone.IsStructured() )
00352               {
00353                 vector< cgsize_t > resIDs; resIDs.reserve( points.size() / IndexSize() );
00354                 for ( size_t i = 0; i < points.size(); i += IndexSize() )
00355                   resIDs.push_back( zone.NodeID( points[i+0], points[i+1], points[i+2] ));
00356                 resIDs.swap( points );
00357               }
00358               else if ( zone._nodeIdShift > 0 )
00359               {
00360                 for ( size_t i = 0; i < points.size(); ++i )
00361                   points[i] += zone._nodeIdShift;
00362               }
00363             }
00364             for ( size_t i = 0; i < ids.size() && i < donorIds.size(); ++i )
00365               _nodeReplacementMap.insert( make_pair( ids[i], donorIds[i] ));
00366           }
00367           else
00368           {
00369             error = cg_get_error();
00370           }
00371         }
00372         else
00373         {
00374           error = cg_get_error();
00375         }
00376       }
00377     }
00378     else
00379     {
00380       error = cg_get_error();
00381     }
00382     return error;
00383   }
00384 
00385   //================================================================================
00390   //================================================================================
00391 
00392   void TZoneData::ReplaceNodes( cgsize_t* ids, int nbIds, int idShift/* = 0*/ ) const
00393   {
00394     if ( !_nodeReplacementMap.empty() )
00395     {
00396       map< int, int >::const_iterator it, end = _nodeReplacementMap.end();
00397       for ( size_t i = 0; i < nbIds; ++i )
00398         if (( it = _nodeReplacementMap.find( ids[i] + idShift)) != end )
00399           ids[i] = it->second;
00400         else
00401           ids[i] += idShift;
00402     }
00403     else if ( idShift )
00404     {
00405       for ( size_t i = 0; i < nbIds; ++i )
00406         ids[i] += idShift;
00407     }
00408   }
00409   //================================================================================
00413   SMDS_MeshElement* add_0D(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00414   {
00415     return mesh->Add0DElementWithID( ids[0], ID );
00416   }
00417   SMDS_MeshElement* add_BAR_2(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00418   {
00419     return mesh->AddEdgeWithID( ids[0], ids[1], ID );
00420   }
00421   SMDS_MeshElement* add_BAR_3(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00422   {
00423     return mesh->AddEdgeWithID( ids[0], ids[1], ids[2], ID );
00424   }
00425   SMDS_MeshElement* add_TRI_3(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00426   {
00427     return mesh->AddFaceWithID( ids[0], ids[2], ids[1], ID );
00428   }
00429   SMDS_MeshElement* add_TRI_6(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00430   {
00431     return mesh->AddFaceWithID( ids[0], ids[2], ids[1], ids[5], ids[4], ids[3], ID );
00432   }
00433   SMDS_MeshElement* add_QUAD_4(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00434   {
00435     return mesh->AddFaceWithID( ids[0], ids[3], ids[2], ids[1], ID );
00436   }
00437   SMDS_MeshElement* add_QUAD_8(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00438   {
00439     return mesh->AddFaceWithID( ids[0],ids[3],ids[2],ids[1],ids[7],ids[6],ids[5],ids[4], ID );
00440   }
00441   SMDS_MeshElement* add_QUAD_9(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00442   {
00443     return mesh->AddFaceWithID( ids[0],ids[3],ids[2],ids[1],ids[7],ids[6],ids[5],ids[4],ids[8], ID);
00444   }
00445   SMDS_MeshElement* add_TETRA_4(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00446   {
00447     return mesh->AddVolumeWithID( ids[0], ids[2], ids[1], ids[3], ID );
00448   }
00449   SMDS_MeshElement* add_TETRA_10(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00450   {
00451     return mesh->AddVolumeWithID( ids[0],ids[2],ids[1],ids[3],ids[6],
00452                                   ids[5],ids[4],ids[7],ids[9],ids[8], ID );
00453   }
00454   SMDS_MeshElement* add_PYRA_5(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00455   {
00456     return mesh->AddVolumeWithID( ids[0],ids[3],ids[2],ids[1],ids[4],ID );
00457   }
00458   SMDS_MeshElement* add_PYRA_13(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00459   {
00460     return mesh->AddVolumeWithID( ids[0],ids[3],ids[2],ids[1],ids[4],ids[8],ids[7],
00461                                   ids[6],ids[5],ids[9],ids[12],ids[11],ids[10], ID );
00462   }
00463   SMDS_MeshElement* add_PENTA_6(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00464   {
00465     return mesh->AddVolumeWithID( ids[0],ids[2],ids[1],ids[3],ids[5],ids[4], ID );
00466   }
00467   SMDS_MeshElement* add_PENTA_15(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00468   {
00469     return mesh->AddVolumeWithID( ids[0],ids[2],ids[1],ids[3],ids[5],ids[4],ids[8],ids[7],
00470                                   ids[6],ids[9],ids[11],ids[10],ids[14],ids[13],ids[12], ID );
00471   }
00472   SMDS_MeshElement* add_HEXA_8(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00473   {
00474     return mesh->AddVolumeWithID( ids[0],ids[3],ids[2],ids[1],ids[4],ids[7],ids[6],ids[5], ID );
00475   }
00476   SMDS_MeshElement* add_HEXA_20(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00477   {
00478     return mesh->AddVolumeWithID( ids[0],ids[3],ids[2],ids[1],ids[4],ids[7],ids[6],
00479                                   ids[5],ids[11],ids[10],ids[9],ids[8],ids[12],ids[15],
00480                                   ids[14],ids[13],ids[19],ids[18],ids[17],ids[16], ID );
00481   }
00482   SMDS_MeshElement* add_HEXA_27(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00483   {
00484     return mesh->AddVolumeWithID( ids[0],ids[3],ids[2],ids[1],ids[4],ids[7],ids[6],
00485                                   ids[5],ids[11],ids[10],ids[9],ids[8],ids[12],ids[15],
00486                                   ids[14],ids[13],ids[19],ids[18],ids[17],ids[16],
00487                                   ids[20],ids[24],ids[23],ids[22],ids[21],ids[25],ids[26], ID );
00488   }
00489   SMDS_MeshElement* add_NGON(cgsize_t* ids, SMESHDS_Mesh* mesh, int ID)
00490   {
00491     vector<int> idVec( ids[0] );
00492     for ( int i = 0; i < ids[0]; ++i )
00493       idVec[ i ] = (int) ids[ i + 1];
00494     return mesh->AddPolygonalFaceWithID( idVec, ID );
00495   }
00496 
00497   typedef SMDS_MeshElement* (* PAddElemFun) (cgsize_t* ids, SMESHDS_Mesh* mesh, int ID);
00498   
00499   //================================================================================
00503   //================================================================================
00504 
00505   PAddElemFun* getAddElemFunTable()
00506   {
00507     static vector< PAddElemFun > funVec;
00508     if ( funVec.empty() )
00509     {
00510       funVec.resize( NofValidElementTypes, (PAddElemFun)0 );
00511       funVec[ CGNS_ENUMV( NODE     )] = add_0D      ;
00512       funVec[ CGNS_ENUMV( BAR_2    )] = add_BAR_2   ;
00513       funVec[ CGNS_ENUMV( BAR_3    )] = add_BAR_3   ;
00514       funVec[ CGNS_ENUMV( TRI_3    )] = add_TRI_3   ;
00515       funVec[ CGNS_ENUMV( TRI_6    )] = add_TRI_6   ;
00516       funVec[ CGNS_ENUMV( QUAD_4   )] = add_QUAD_4  ;
00517       funVec[ CGNS_ENUMV( QUAD_8   )] = add_QUAD_8  ;
00518       funVec[ CGNS_ENUMV( QUAD_9   )] = add_QUAD_9  ;
00519       funVec[ CGNS_ENUMV( TETRA_4  )] = add_TETRA_4 ;
00520       funVec[ CGNS_ENUMV( TETRA_10 )] = add_TETRA_10;
00521       funVec[ CGNS_ENUMV( PYRA_5   )] = add_PYRA_5  ;
00522       funVec[ CGNS_ENUMV( PYRA_13  )] = add_PYRA_13 ;
00523       funVec[ CGNS_ENUMV( PYRA_14  )] = add_PYRA_13 ;
00524       funVec[ CGNS_ENUMV( PENTA_6  )] = add_PENTA_6 ;
00525       funVec[ CGNS_ENUMV( PENTA_15 )] = add_PENTA_15;
00526       funVec[ CGNS_ENUMV( PENTA_18 )] = add_PENTA_15;
00527       funVec[ CGNS_ENUMV( HEXA_8   )] = add_HEXA_8  ;
00528       funVec[ CGNS_ENUMV( HEXA_20  )] = add_HEXA_20 ;
00529       funVec[ CGNS_ENUMV( HEXA_27  )] = add_HEXA_27 ;
00530       funVec[ CGNS_ENUMV( NGON_n   )] = add_NGON    ;
00531     }
00532     return &funVec[0];
00533   }
00534 
00535   //================================================================================
00539   //================================================================================
00540 
00541   const SMDS_MeshElement* findElement(const cgsize_t*     nodeIDs,
00542                                       const int           nbNodes,
00543                                       const SMESHDS_Mesh* mesh)
00544   {
00545     const SMDS_MeshNode* nn[4]; // look for quad4 or seg2
00546     if (( nn[0] = mesh->FindNode( nodeIDs[0] )))
00547     {
00548       SMDSAbs_ElementType eType = nbNodes==4 ? SMDSAbs_Face : SMDSAbs_Edge;
00549       SMDS_ElemIteratorPtr eIt = nn[0]->GetInverseElementIterator( eType );
00550       if ( eIt->more() )
00551         for ( int i = 1; i < nbNodes; ++i )
00552           nn[i] = mesh->FindNode( nodeIDs[i] );
00553       while ( eIt->more() )
00554       {
00555         const SMDS_MeshElement* e = eIt->next();
00556         if ( e->NbNodes() == nbNodes )
00557         {
00558           bool elemOK = true;
00559           for ( int i = 1; i < nbNodes && elemOK; ++i )
00560             elemOK = ( e->GetNodeIndex( nn[i] ) >= 0 );
00561           if ( elemOK )
00562             return e;
00563         }
00564       } 
00565     }
00566     return 0;
00567   }
00568 
00569 } // namespace
00570 
00571 //================================================================================
00575 //================================================================================
00576 
00577 Driver_Mesh::Status DriverCGNS_Read::Perform()
00578 {
00579   myErrorMessages.clear();
00580 
00581   Status aResult;
00582   if (( aResult = open() ) != DRS_OK )
00583     return aResult;
00584 
00585   // read nb of meshes (CGNSBase_t)
00586   if ( myMeshId < 0 || myMeshId >= GetNbMeshes(aResult))
00587     return addMessage( SMESH_Comment("Invalid mesh index :") << myMeshId );
00588 
00589   // read a name and a dimension of the mesh
00590   const int cgnsBase = myMeshId + 1;
00591   char meshName[CGNS_NAME_SIZE];
00592   int meshDim, spaceDim;
00593   if ( cg_base_read( _fn, cgnsBase, meshName, &meshDim, &spaceDim) != CG_OK )
00594     return addMessage( cg_get_error() );
00595 
00596   if ( spaceDim < 1 || spaceDim > 3 )
00597     return addMessage( SMESH_Comment("Invalid space dimension: ") << spaceDim
00598                        << " in mesh '" << meshName << "'");
00599 
00600   myMeshName = meshName;
00601 
00602   // read nb of domains (Zone_t) in the mesh
00603   int nbZones = 0;
00604   if ( cg_nzones (_fn, cgnsBase, &nbZones) != CG_OK )
00605     return addMessage( cg_get_error() );
00606 
00607   if ( nbZones < 1 )
00608     return addMessage( SMESH_Comment("Empty mesh: '") << meshName << "'");
00609 
00610   // read the domains (zones)
00611   // ------------------------
00612   map< string, TZoneData > zonesByName;
00613   char name[CGNS_NAME_SIZE];
00614   cgsize_t sizes[NB_ZONE_SIZE_VAL];
00615   memset(sizes, 0, NB_ZONE_SIZE_VAL * sizeof(cgsize_t));
00616 
00617   const SMDS_MeshInfo& meshInfo = myMesh->GetMeshInfo();
00618   int groupID = myMesh->GetGroups().size();
00619 
00620   for ( int iZone = 1; iZone <= nbZones; ++iZone )
00621   {
00622     // size and name of a zone
00623     if ( cg_zone_read( _fn, cgnsBase, iZone, name, sizes) != CG_OK) {
00624       addMessage( cg_get_error() );
00625       continue;
00626     }
00627     TZoneData& zone = zonesByName[ name ];
00628     zone._id          = iZone;
00629     zone._nodeIdShift = meshInfo.NbNodes();
00630     zone._elemIdShift = meshInfo.NbElements();
00631     zone.SetSizeAndDim( sizes, meshDim );
00632 
00633     // mesh type of the zone
00634     if ( cg_zone_type ( _fn, cgnsBase, iZone, &zone._type) != CG_OK) {
00635       addMessage( cg_get_error() );
00636       continue;
00637     }
00638 
00639     switch ( zone._type )
00640     {
00641     case CGNS_ENUMV( Unstructured ):
00642     case CGNS_ENUMV( Structured ):
00643       break;
00644     case CGNS_ENUMV( ZoneTypeNull ):
00645       addMessage( "Meshes with ZoneTypeNull are not supported");
00646       continue;
00647     case CGNS_ENUMV( ZoneTypeUserDefined ):
00648       addMessage( "Meshes with ZoneTypeUserDefined are not supported");
00649       continue;
00650     default:
00651       addMessage( "Unknown ZoneType_t");
00652       continue;
00653     }
00654 
00655     // -----------
00656     // Read nodes
00657     // -----------
00658 
00659     if ( cg_ncoords( _fn, cgnsBase, iZone, &spaceDim) != CG_OK ) {
00660       addMessage( cg_get_error() );
00661       continue;
00662     }
00663     if ( spaceDim < 1 ) {
00664       addMessage( SMESH_Comment("No coordinates defined in zone ")
00665                   << iZone << " of Mesh " << myMeshId );
00666       continue;
00667     }
00668     // read coordinates
00669 
00670     cgsize_t rmin[3] = {1,1,1}; // range of nodes to read
00671     cgsize_t rmax[3] = {1,1,1};
00672     int nbNodes = rmax[0] = zone._sizes[0];
00673     if ( zone.IsStructured())
00674       for ( int i = 1; i < meshDim; ++i )
00675         nbNodes *= rmax[i] = zone._sizes[i];
00676 
00677     vector<double> coords[3];
00678     for ( int c = 1; c <= spaceDim; ++c)
00679     {
00680       coords[c-1].resize( nbNodes );
00681 
00682       CGNS_ENUMV( DataType_t ) type;
00683       if ( cg_coord_info( _fn, cgnsBase, iZone, c, &type, name) != CG_OK ||
00684            cg_coord_read( _fn, cgnsBase, iZone, name, CGNS_ENUMV(RealDouble),
00685                           rmin, rmax, (void*)&(coords[c-1][0])) != CG_OK)
00686       {
00687         addMessage( cg_get_error() );
00688         coords[c-1].clear();
00689         break;
00690       }
00691     }
00692     if ( coords[ spaceDim-1 ].empty() )
00693       continue; // there was an error while reading coordinates 
00694 
00695     // fill coords with zero if spaceDim < 3
00696     for ( int c = 2; c <= 3; ++c)
00697       if ( coords[ c-1 ].empty() )
00698         coords[ c-1 ].resize( nbNodes, 0.0 );
00699 
00700     // create nodes
00701     try {
00702       for ( int i = 0; i < nbNodes; ++i )
00703         myMesh->AddNodeWithID( coords[0][i], coords[1][i], coords[2][i], i+1+zone._nodeIdShift );
00704     }
00705     catch ( std::exception& exc ) // expect std::bad_alloc
00706     {
00707       addMessage( exc.what() );
00708       break;
00709     }
00710 
00711     // Read connectivity between zones. Nodes of the zone interface will be
00712     // replaced withing the zones read later
00713     string err = zone.ReadZonesConnection( _fn, cgnsBase, zonesByName );
00714     if ( !err.empty() )
00715       addMessage( err );
00716 
00717     // --------------
00718     // Read elements
00719     // --------------
00720     if ( zone.IsStructured())
00721     {
00722       int nbI = zone._sizeX - 1, nbJ = zone._sizeY - 1, nbK = zone._sizeZ - 1;
00723       cgsize_t nID[8];
00724       if ( meshDim > 2 && nbK > 0 )
00725       {
00726         for ( int k = 1; k <= nbK; ++k )
00727           for ( int j = 1; j <= nbJ; ++j )
00728             for ( int i = 1; i <= nbI; ++i )
00729             {
00730               zone.CellNodes( i, j, k, nID );
00731               zone.ReplaceNodes( nID, 8 );
00732               myMesh->AddVolumeWithID(nID[0],nID[1],nID[2],nID[3],nID[4],nID[5],nID[6],nID[7],
00733                                       meshInfo.NbElements()+1);
00734             }
00735       }
00736       else if ( meshDim > 1 && nbJ > 0 )
00737       {
00738         for ( int j = 1; j <= nbJ; ++j )
00739           for ( int i = 1; i <= nbI; ++i )
00740           {
00741             zone.CellNodes( i, j, nID );
00742             zone.ReplaceNodes( nID, 4 );
00743             myMesh->AddFaceWithID(nID[0],nID[1],nID[2],nID[3], meshInfo.NbElements()+1);
00744           }
00745       }
00746       else if ( meshDim > 0 && nbI > 0 )
00747       {
00748         nID[0] = zone.NodeID( 1, 0, 0 );
00749         for ( int i = 1; i <= nbI; ++i, ++nID[0] )
00750         {
00751           nID[1] = nID[0]+1;
00752           zone.ReplaceNodes( nID, 2 );
00753           myMesh->AddEdgeWithID(nID[0],nID[1], meshInfo.NbElements()+1);
00754         }
00755       }
00756     }
00757     else
00758     {
00759       // elements can be stored in different sections each dedicated to one element type
00760       int nbSections = 0;
00761       if ( cg_nsections( _fn, cgnsBase, iZone, &nbSections) != CG_OK)
00762       {
00763         addMessage( cg_get_error() );
00764         continue;
00765       }
00766       PAddElemFun* addElemFuns = getAddElemFunTable(), curAddElemFun = 0;
00767       int nbNotSuppElem = 0; // nb elements of not supported types
00768       bool polyhedError = false; // error at polyhedron creation
00769 
00770       // read element data
00771 
00772       CGNS_ENUMT( ElementType_t ) elemType;
00773       cgsize_t start, end; // range of ids of elements of a zone
00774       cgsize_t eDataSize = 0;
00775       int nbBnd, parent_flag;
00776       for ( int iSec = 1; iSec <= nbSections; ++iSec )
00777       {
00778         if ( cg_section_read( _fn, cgnsBase, iZone, iSec, name, &elemType,
00779                               &start, &end, &nbBnd, &parent_flag) != CG_OK ||
00780              cg_ElementDataSize( _fn, cgnsBase, iZone, iSec, &eDataSize ) != CG_OK )
00781         {
00782           addMessage( cg_get_error() );
00783           continue;
00784         }
00785         vector< cgsize_t > elemData( eDataSize );
00786         if ( cg_elements_read( _fn, cgnsBase, iZone, iSec, &elemData[0], NULL ) != CG_OK )
00787         {
00788           addMessage( cg_get_error() );
00789           continue;
00790         }
00791         // store elements
00792 
00793         int pos = 0, cgnsNbNodes = 0, elemID = start + zone._elemIdShift;
00794         cg_npe( elemType, &cgnsNbNodes ); // get nb nodes by element type
00795         curAddElemFun = addElemFuns[ elemType ];
00796         SMDS_MeshElement* newElem = 0;
00797         const SMDS_MeshElement* face;
00798 
00799         while ( pos < eDataSize )
00800         {
00801           CGNS_ENUMT( ElementType_t ) currentType = elemType;
00802           if ( currentType == CGNS_ENUMV( MIXED )) {
00803             //ElementConnectivity = Etype1, Node11, Node21, ... NodeN1,
00804             //                      Etype2, Node12, Node22, ... NodeN2,
00805             //                      ...
00806             //                      EtypeM, Node1M, Node2M, ... NodeNM
00807             currentType = (CGNS_ENUMT(ElementType_t)) elemData[ pos++ ];
00808             cg_npe( currentType, &cgnsNbNodes );
00809             curAddElemFun = addElemFuns[ currentType ];
00810           }
00811           if ( cgnsNbNodes < 1 ) // poly elements
00812           {
00813             if ( currentType == CGNS_ENUMV( NFACE_n )) // polyhedron
00814             {
00815               //ElementConnectivity = Nfaces1, Face11, Face21, ... FaceN1,
00816               //                      Nfaces2, Face12, Face22, ... FaceN2,
00817               //                      ...
00818               //                      NfacesM, Face1M, Face2M, ... FaceNM
00819               const int nbFaces = elemData[ pos++ ];
00820               vector<int> quantities( nbFaces );
00821               vector<const SMDS_MeshNode*> nodes, faceNodes;
00822               nodes.reserve( nbFaces * 4 );
00823               for ( int iF = 0; iF < nbFaces; ++iF )
00824               {
00825                 const int faceID = std::abs( elemData[ pos++ ]) + zone._elemIdShift; 
00826                 if (( face = myMesh->FindElement( faceID )) && face->GetType() == SMDSAbs_Face )
00827                 {
00828                   const bool reverse = ( elemData[ pos-1 ] < 0 );
00829                   const int    iQuad = face->IsQuadratic() ? 1 : 0;
00830                   SMDS_ElemIteratorPtr nIter = face->interlacedNodesElemIterator();
00831                   faceNodes.assign( SMDS_MeshElement::iterator( nIter ),
00832                                     SMDS_MeshElement::iterator());
00833                   if ( iQuad && reverse )
00834                     nodes.push_back( faceNodes[0] );
00835                   if ( reverse )
00836                     nodes.insert( nodes.end(), faceNodes.rbegin(), faceNodes.rend() - iQuad );
00837                   else
00838                     nodes.insert( nodes.end(), faceNodes.begin(), faceNodes.end() );
00839 
00840                   quantities[ iF ] = face->NbNodes();
00841                 }
00842                 else {
00843                   polyhedError = true;
00844                   break;
00845                 }
00846               }
00847               if ( quantities.back() )
00848               {
00849                 myMesh->AddPolyhedralVolumeWithID( nodes, quantities, elemID );
00850               }
00851             }
00852             else if ( currentType == CGNS_ENUMV( NGON_n )) // polygon
00853             {
00854               // ElementConnectivity = Nnodes1, Node11, Node21, ... NodeN1,
00855               //                       Nnodes2, Node12, Node22, ... NodeN2,
00856               //                       ...
00857               //                       NnodesM, Node1M, Node2M, ... NodeNM
00858               const int nbNodes = elemData[ pos ];
00859               zone.ReplaceNodes( &elemData[pos+1], nbNodes, zone._nodeIdShift );
00860               newElem = add_NGON( &elemData[pos  ], myMesh, elemID );
00861               pos += nbNodes + 1;
00862             }
00863           }
00864           else // standard elements
00865           {
00866             zone.ReplaceNodes( &elemData[pos], cgnsNbNodes, zone._nodeIdShift );
00867             newElem = curAddElemFun( &elemData[pos], myMesh, elemID );
00868             pos += cgnsNbNodes;
00869             nbNotSuppElem += int( newElem && newElem->NbNodes() != cgnsNbNodes );
00870           }
00871           elemID++;
00872 
00873         } // loop on elemData
00874       } // loop on cgns sections
00875 
00876       if ( nbNotSuppElem > 0 )
00877         addMessage( SMESH_Comment(nbNotSuppElem) << " elements of not supported types"
00878                     << " have beem converted to close types");
00879       if ( polyhedError )
00880         addMessage( "Some polyhedral elements have been skipped due to internal(?) errors" );
00881 
00882     } // reading unstructured elements
00883 
00884     zone._nbNodes = meshInfo.NbNodes() - zone._nodeIdShift;
00885     zone._nbElems = meshInfo.NbElements() - zone._elemIdShift;
00886 
00887     // -------------------------------------------
00888     // Read Boundary Conditions into SMESH groups
00889     // -------------------------------------------
00890     int nbBC = 0;
00891     if ( cg_nbocos( _fn, cgnsBase, iZone, &nbBC) == CG_OK )
00892     {
00893       CGNS_ENUMT( BCType_t ) bcType;
00894       CGNS_ENUMT( PointSetType_t ) psType;
00895       CGNS_ENUMT( DataType_t ) normDataType;
00896       cgsize_t nbPnt, normFlag;
00897       int normIndex[3], nbDS;
00898       for ( int iBC = 1; iBC <= nbBC; ++iBC )
00899       {
00900         if ( cg_boco_info( _fn, cgnsBase, iZone, iBC, name, &bcType, &psType,
00901                            &nbPnt, normIndex, &normFlag, &normDataType, &nbDS ) != CG_OK )
00902         {
00903           addMessage( cg_get_error() );
00904           continue;
00905         }
00906         vector< cgsize_t > ids( nbPnt * zone.IndexSize() );
00907         CGNS_ENUMT( GridLocation_t ) location;
00908         if ( cg_boco_read( _fn, cgnsBase, iZone, iBC, &ids[0], NULL ) != CG_OK ||
00909              cg_boco_gridlocation_read( _fn, cgnsBase, iZone, iBC, &location) != CG_OK )
00910         {
00911           addMessage( cg_get_error() );
00912           continue;
00913         }
00914         SMDSAbs_ElementType elemType = SMDSAbs_All;
00915         switch ( location ) {
00916         case CGNS_ENUMV( Vertex      ): elemType = SMDSAbs_Node; break;
00917         case CGNS_ENUMV( FaceCenter  ): elemType = SMDSAbs_Face; break;
00918         case CGNS_ENUMV( IFaceCenter ): elemType = SMDSAbs_Face; break;
00919         case CGNS_ENUMV( JFaceCenter ): elemType = SMDSAbs_Face; break;
00920         case CGNS_ENUMV( KFaceCenter ): elemType = SMDSAbs_Face; break;
00921         case CGNS_ENUMV( EdgeCenter  ): elemType = SMDSAbs_Edge; break;
00922         default:;
00923         }
00924         SMESHDS_Group* group = new SMESHDS_Group ( groupID++, myMesh, elemType );
00925         myMesh->AddGroup( group );
00926         SMESH_Comment groupName( name ); groupName << " " << cg_BCTypeName( bcType );
00927         group->SetStoreName( groupName.c_str() );
00928         SMDS_MeshGroup& groupDS = group->SMDSGroup();
00929 
00930         if ( elemType == SMDSAbs_Node )
00931         {
00932           if ( zone.IsStructured() )
00933           {
00934             vector< cgsize_t > nodeIds;
00935             if ( psType == CGNS_ENUMV( PointRange ))
00936             {
00937               // nodes are given as (ijkMin, ijkMax)
00938               TPointRangeIterator idIt( & ids[0], meshDim );
00939               nodeIds.reserve( idIt.Size() );
00940               while ( idIt.More() )
00941                 nodeIds.push_back( zone.NodeID( idIt.Next() ));
00942             }
00943             else
00944             {
00945               // nodes are given as (ijk1, ijk2, ..., ijkN)
00946               nodeIds.reserve( ids.size() / meshDim );
00947               for ( size_t i = 0; i < ids.size(); i += meshDim )
00948                 nodeIds.push_back( zone.NodeID( ids[i], ids[i+1], ids[i+2] ));
00949             }
00950             ids.swap( nodeIds );
00951           }
00952           else if ( zone._nodeIdShift )
00953           {
00954             for ( size_t i = 0; i < ids.size(); ++i )
00955               ids[i] += zone._nodeIdShift;
00956           }
00957           zone.ReplaceNodes( &ids[0], ids.size() );
00958 
00959           for ( size_t i = 0; i < ids.size(); ++i )
00960             if ( const SMDS_MeshNode* n = myMesh->FindNode( ids[i] ))
00961               groupDS.Add( n );
00962         }
00963         else // BC applied to elements
00964         {
00965           if ( zone.IsStructured() )
00966           {
00967             int axis = 0; // axis perpendiculaire to which boundary elements are oriented
00968             if ( ids.size() >= meshDim * 2 )
00969             {
00970               for ( ; axis < meshDim; ++axis )
00971                 if ( ids[axis] - ids[axis+meshDim] == 0 )
00972                   break;
00973             }
00974             else
00975             {
00976               for ( ; axis < meshDim; ++axis )
00977                 if ( normIndex[axis] != 0 )
00978                   break;
00979             }
00980             if ( axis == meshDim )
00981             {
00982               addMessage( SMESH_Comment("Invalid NormalIndex in BC ") << name );
00983               continue;
00984             }
00985             const int nbElemNodesByDim[] = { 1, 2, 4, 8 };
00986             const int nbElemNodes = nbElemNodesByDim[ meshDim ];
00987 
00988             if ( psType == CGNS_ENUMV( PointRange ) ||
00989                  psType == CGNS_ENUMV( ElementRange ))
00990             {
00991               // elements are given as (ijkMin, ijkMax)
00992               typedef void (TZoneData::*PGetNodesFun)( const gp_XYZ& ijk, cgsize_t* ids ) const;
00993               PGetNodesFun getNodesFun = 0;
00994               if ( elemType == SMDSAbs_Face  && meshDim == 3 )
00995                 switch ( axis ) {
00996                 case 0: getNodesFun = & TZoneData::IFaceNodes;
00997                 case 1: getNodesFun = & TZoneData::JFaceNodes;
00998                 case 2: getNodesFun = & TZoneData::KFaceNodes;
00999                 }
01000               else if ( elemType == SMDSAbs_Edge && meshDim == 2 )
01001                 switch ( axis ) {
01002                 case 0: getNodesFun = & TZoneData::IEdgeNodes;
01003                 case 1: getNodesFun = & TZoneData::JEdgeNodes;
01004                 }
01005               if ( !getNodesFun )
01006               {
01007                 addMessage( SMESH_Comment("Unsupported BC location in BC ") << name
01008                             << " " << cg_GridLocationName( location )
01009                             << " in " << meshDim << " mesh");
01010                 continue;
01011               }
01012               TPointRangeIterator rangeIt( & ids[0], meshDim );
01013               vector< cgsize_t > elemNodeIds( rangeIt.Size() * nbElemNodes );
01014               for ( int i = 0; rangeIt.More(); i+= nbElemNodes )
01015                 (zone.*getNodesFun)( rangeIt.Next(), &elemNodeIds[i] );
01016 
01017               ids.swap( elemNodeIds );
01018             }
01019             else
01020             {
01021               // elements are given as (ijk1, ijk2, ..., ijkN)
01022               typedef void (TZoneData::*PGetNodesFun)( int i, int j, int k, cgsize_t* ids ) const;
01023               PGetNodesFun getNodesFun = 0;
01024               if ( elemType == SMDSAbs_Face )
01025                 switch ( axis ) {
01026                 case 0: getNodesFun = & TZoneData::IFaceNodes;
01027                 case 1: getNodesFun = & TZoneData::JFaceNodes;
01028                 case 2: getNodesFun = & TZoneData::KFaceNodes;
01029                 }
01030               else if ( elemType == SMDSAbs_Edge && meshDim == 2 )
01031                 switch ( axis ) {
01032                 case 0: getNodesFun = & TZoneData::IEdgeNodes;
01033                 case 1: getNodesFun = & TZoneData::JEdgeNodes;
01034                 }
01035               if ( !getNodesFun )
01036               {
01037                 addMessage( SMESH_Comment("Unsupported BC location in BC ") << name
01038                             << " " << cg_GridLocationName( location )
01039                             << " in " << meshDim << " mesh");
01040                 continue;
01041               }
01042               vector< cgsize_t > elemNodeIds( ids.size()/meshDim * nbElemNodes );
01043               for ( size_t i = 0, j = 0; i < ids.size(); i += meshDim, j += nbElemNodes )
01044                 (zone.*getNodesFun)( ids[i], ids[i+1], ids[i+2], &elemNodeIds[j] );
01045 
01046               ids.swap( elemNodeIds );
01047             }
01048             zone.ReplaceNodes( &ids[0], ids.size() );
01049 
01050             PAddElemFun addElemFun = 0;
01051             switch ( meshDim ) {
01052             case 1: addElemFun = & add_BAR_2;
01053             case 2: addElemFun = & add_QUAD_4;
01054             case 3: addElemFun = & add_HEXA_8;
01055             }
01056             int elemID = meshInfo.NbElements();
01057             const SMDS_MeshElement* elem = 0;
01058             for ( size_t i = 0; i < ids.size(); i += nbElemNodes )
01059             {
01060               if ( iZone == 1 || !( elem = findElement( &ids[i], nbElemNodes, myMesh )))
01061                 elem = addElemFun( &ids[i], myMesh, ++elemID );
01062               groupDS.Add( elem );
01063             }
01064           }
01065           else // unstructured zone
01066           {
01067             if ( zone._elemIdShift )
01068               for ( size_t i = 0; i < ids.size(); ++i )
01069                 ids[i] += zone._elemIdShift;
01070 
01071             if ( psType == CGNS_ENUMV( PointRange ) && ids.size() == 2 )
01072             {
01073               for ( size_t i = ids[0]; i <= ids[1]; ++i )
01074                 if ( const SMDS_MeshElement* e = myMesh->FindElement( i ))
01075                   groupDS.Add( e );
01076             }
01077             else
01078             {
01079               for ( size_t i = 0; i < ids.size(); ++i )
01080                 if ( const SMDS_MeshElement* e = myMesh->FindElement( ids[i] ))
01081                   groupDS.Add( e );
01082             }
01083           }
01084         } // end "BC applied to elements"
01085 
01086         // to have group type according to a real elem type
01087         group->SetType( groupDS.GetType() );
01088 
01089       } // loop on BCs of the zone
01090     }
01091     else
01092     {
01093       addMessage( cg_get_error() );
01094     }
01095   } // loop on the zones of a mesh
01096 
01097 
01098   // ------------------------------------------------------------------------
01099   // Make groups for multiple zones and remove free nodes at zone interfaces
01100   // ------------------------------------------------------------------------
01101   map< string, TZoneData >::iterator nameZoneIt = zonesByName.begin();
01102   for ( ; nameZoneIt != zonesByName.end(); ++nameZoneIt )
01103   {
01104     TZoneData& zone = nameZoneIt->second;
01105     if ( zone._nbElems == 0 ) continue;
01106     if ( zone._nbElems == meshInfo.NbElements() ) break; // there is only one non-empty zone
01107 
01108     // make a group
01109     SMDSAbs_ElementType elemType = myMesh->GetElementType( zone._elemIdShift + 1,
01110                                                            /*iselem=*/true );
01111     SMESHDS_Group* group = new SMESHDS_Group ( groupID++, myMesh, elemType );
01112     myMesh->AddGroup( group );
01113     group->SetStoreName( nameZoneIt->first.c_str() );
01114     SMDS_MeshGroup& groupDS = group->SMDSGroup();
01115 
01116     for ( int i = 1; i <= zone._nbElems; ++i )
01117       if ( const SMDS_MeshElement* e = myMesh->FindElement( i + zone._elemIdShift ))
01118         groupDS.Add( e );
01119 
01120     // remove free nodes
01121     map< int, int >::iterator nnRmKeepIt = zone._nodeReplacementMap.begin();
01122     for ( ; nnRmKeepIt != zone._nodeReplacementMap.end(); ++nnRmKeepIt )
01123       if ( const SMDS_MeshNode* n = myMesh->FindNode( nnRmKeepIt->first ))
01124         if ( n->NbInverseElements() == 0 )
01125           myMesh->RemoveFreeNode( n, (SMESHDS_SubMesh *)0, /*fromGroups=*/false );
01126   }
01127 
01128   aResult = myErrorMessages.empty() ? DRS_OK : DRS_WARN_SKIP_ELEM;
01129 
01130   return aResult;
01131 }
01132 
01133 //================================================================================
01137 //================================================================================
01138 
01139 DriverCGNS_Read::DriverCGNS_Read()
01140 {
01141   _fn = -1;
01142 }
01143 //================================================================================
01147 //================================================================================
01148 
01149 DriverCGNS_Read::~DriverCGNS_Read()
01150 {
01151   if ( _fn > 0 )
01152     cg_close( _fn );
01153 }
01154 
01155 //================================================================================
01159 //================================================================================
01160 
01161 Driver_Mesh::Status DriverCGNS_Read::open()
01162 {
01163   if ( _fn < 0 )
01164   {
01165     
01166 #ifdef CG_MODE_READ
01167     int res = cg_open(myFile.c_str(), CG_MODE_READ, &_fn);
01168 #else
01169     int res = cg_open(myFile.c_str(), MODE_READ, &_fn);
01170 #endif
01171     if ( res != CG_OK)
01172     {
01173       addMessage( cg_get_error(), /*fatal = */true );
01174     }
01175   }
01176   return _fn >= 0 ? DRS_OK : DRS_FAIL;
01177 }
01178 
01179 //================================================================================
01183 //================================================================================
01184 
01185 int DriverCGNS_Read::GetNbMeshes(Status& theStatus)
01186 {
01187   if (( theStatus = open()) != DRS_OK )
01188     return 0;
01189 
01190   int nbases = 0;
01191   if(cg_nbases( _fn, &nbases) != CG_OK)
01192     theStatus = addMessage( cg_get_error(), /*fatal = */true );
01193 
01194   return nbases;
01195 }