1 // Copyright (C) 2007-2013 CEA/DEN, EDF R&D
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
19 // Author : Anthony Geay (CEA/DEN)
21 #include "MEDFileMeshElt.hxx"
22 #include "MEDFileMeshReadSelector.hxx"
24 #include "MEDCouplingUMesh.hxx"
26 #include "InterpKernelAutoPtr.hxx"
27 #include "CellModel.hxx"
31 extern med_geometry_type typmai3[32];
33 using namespace ParaMEDMEM;
35 MEDFileUMeshPerType *MEDFileUMeshPerType::New(med_idt fid, const char *mName, int dt, int it, int mdim, med_geometry_type geoElt, INTERP_KERNEL::NormalizedCellType geoElt2, MEDFileMeshReadSelector *mrs)
37 med_entity_type whichEntity;
38 if(!isExisting(fid,mName,dt,it,geoElt,whichEntity))
40 return new MEDFileUMeshPerType(fid,mName,dt,it,mdim,geoElt,geoElt2,whichEntity,mrs);
43 std::size_t MEDFileUMeshPerType::getHeapMemorySizeWithoutChildren() const
48 std::vector<const BigMemoryObject *> MEDFileUMeshPerType::getDirectChildren() const
50 std::vector<const BigMemoryObject *> ret;
51 if((const MEDCoupling1GTUMesh *)_m)
52 ret.push_back((const MEDCoupling1GTUMesh *)_m);
53 if((const DataArrayInt *)_num)
54 ret.push_back((const DataArrayInt *)_num);
55 if((const DataArrayInt *)_fam)
56 ret.push_back((const DataArrayInt *)_fam);
57 if((const DataArrayAsciiChar *)_names)
58 ret.push_back((const DataArrayAsciiChar *)_names);
62 bool MEDFileUMeshPerType::isExisting(med_idt fid, const char *mName, int dt, int it, med_geometry_type geoElt, med_entity_type& whichEntity)
64 static const med_entity_type entities[3]={MED_CELL,MED_DESCENDING_FACE,MED_DESCENDING_EDGE};
68 med_bool changement,transformation;
69 int tmp=MEDmeshnEntity(fid,mName,dt,it,entities[i],geoElt,MED_CONNECTIVITY,MED_NODAL,
70 &changement,&transformation);
74 whichEntity=entities[i];
76 std::cerr << "WARNING : MEDFile has been detected to be no compilant with MED 3 : Please change entity in MEDFile for geotype " << geoElt << std::endl;
82 int MEDFileUMeshPerType::getDim() const
84 return _m->getMeshDimension();
87 MEDFileUMeshPerType::MEDFileUMeshPerType(med_idt fid, const char *mName, int dt, int it, int mdim, med_geometry_type geoElt, INTERP_KERNEL::NormalizedCellType type,
88 med_entity_type entity, MEDFileMeshReadSelector *mrs):_entity(entity)
90 med_bool changement,transformation;
91 int curNbOfElem=MEDmeshnEntity(fid,mName,dt,it,entity,geoElt,MED_CONNECTIVITY,MED_NODAL,
92 &changement,&transformation);
93 const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(type);
96 loadFromStaticType(fid,mName,dt,it,mdim,curNbOfElem,geoElt,type,entity,mrs);
99 if(type==INTERP_KERNEL::NORM_POLYGON)
101 loadPolyg(fid,mName,dt,it,mdim,curNbOfElem,geoElt,entity,mrs);
104 //if(type==INTERP_KERNEL::NORM_POLYHED)
105 loadPolyh(fid,mName,dt,it,mdim,curNbOfElem,geoElt,entity,mrs);
108 void MEDFileUMeshPerType::loadFromStaticType(med_idt fid, const char *mName, int dt, int it, int mdim, int curNbOfElem, med_geometry_type geoElt, INTERP_KERNEL::NormalizedCellType type,
109 med_entity_type entity, MEDFileMeshReadSelector *mrs)
111 _m=MEDCoupling1SGTUMesh::New(mName,type);
112 MEDCoupling1SGTUMesh *mc(dynamic_cast<MEDCoupling1SGTUMesh *>((MEDCoupling1GTUMesh *)_m));
113 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> conn(DataArrayInt::New());
114 int nbOfNodesPerCell=mc->getNumberOfNodesPerCell();
115 conn->alloc(nbOfNodesPerCell*curNbOfElem,1);
116 MEDmeshElementConnectivityRd(fid,mName,dt,it,entity,geoElt,MED_NODAL,MED_FULL_INTERLACE,conn->getPointer());
117 std::transform(conn->begin(),conn->end(),conn->getPointer(),std::bind2nd(std::plus<int>(),-1));
118 mc->setNodalConnectivity(conn);
119 loadCommonPart(fid,mName,dt,it,mdim,curNbOfElem,geoElt,entity,mrs);
122 void MEDFileUMeshPerType::loadCommonPart(med_idt fid, const char *mName, int dt, int it, int mdim, int curNbOfElem, med_geometry_type geoElt,
123 med_entity_type entity, MEDFileMeshReadSelector *mrs)
125 med_bool changement,transformation;
127 if(MEDmeshnEntity(fid,mName,dt,it,entity,geoElt,MED_FAMILY_NUMBER,MED_NODAL,&changement,&transformation)>0)
129 if(!mrs || mrs->isCellFamilyFieldReading())
131 _fam=DataArrayInt::New();
132 _fam->alloc(curNbOfElem,1);
133 if(MEDmeshEntityFamilyNumberRd(fid,mName,dt,it,entity,geoElt,_fam->getPointer())!=0)
134 std::fill(_fam->getPointer(),_fam->getPointer()+curNbOfElem,0);
138 if(MEDmeshnEntity(fid,mName,dt,it,entity,geoElt,MED_NUMBER,MED_NODAL,&changement,&transformation)>0)
140 if(!mrs || mrs->isCellNumFieldReading())
142 _num=DataArrayInt::New();
143 _num->alloc(curNbOfElem,1);
144 if(MEDmeshEntityNumberRd(fid,mName,dt,it,entity,geoElt,_num->getPointer())!=0)
149 if(MEDmeshnEntity(fid,mName,dt,it,entity,geoElt,MED_NAME,MED_NODAL,&changement,&transformation)>0)
151 if(!mrs || mrs->isCellNameFieldReading())
153 _names=DataArrayAsciiChar::New();
154 _names->alloc(curNbOfElem+1,MED_SNAME_SIZE);//not a bug to avoid the memory corruption due to last \0 at the end
155 if(MEDmeshEntityNameRd(fid,mName,dt,it,entity,geoElt,_names->getPointer())!=0)
158 _names->reAlloc(curNbOfElem);//not a bug to avoid the memory corruption due to last \0 at the end
163 void MEDFileUMeshPerType::loadPolyg(med_idt fid, const char *mName, int dt, int it, int mdim, int arraySize, med_geometry_type geoElt,
164 med_entity_type entity, MEDFileMeshReadSelector *mrs)
166 med_bool changement,transformation;
167 med_int curNbOfElem=MEDmeshnEntity(fid,mName,dt,it,entity,MED_POLYGON,MED_INDEX_NODE,MED_NODAL,&changement,&transformation)-1;
168 _m=MEDCoupling1DGTUMesh::New(mName,INTERP_KERNEL::NORM_POLYGON);
169 MEDCoupling1DGTUMesh *mc(dynamic_cast<MEDCoupling1DGTUMesh *>((MEDCoupling1GTUMesh *)_m));
170 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> conn(DataArrayInt::New()),connI(DataArrayInt::New());
171 conn->alloc(arraySize,1); connI->alloc(curNbOfElem+1,1);
172 MEDmeshPolygonRd(fid,mName,dt,it,MED_CELL,MED_NODAL,connI->getPointer(),conn->getPointer());
173 std::transform(conn->begin(),conn->end(),conn->getPointer(),std::bind2nd(std::plus<int>(),-1));
174 std::transform(connI->begin(),connI->end(),connI->getPointer(),std::bind2nd(std::plus<int>(),-1));
175 mc->setNodalConnectivity(conn,connI);
176 loadCommonPart(fid,mName,dt,it,mdim,curNbOfElem,MED_POLYGON,entity,mrs);
179 void MEDFileUMeshPerType::loadPolyh(med_idt fid, const char *mName, int dt, int it, int mdim, int connFaceLgth, med_geometry_type geoElt,
180 med_entity_type entity, MEDFileMeshReadSelector *mrs)
182 med_bool changement,transformation;
183 med_int indexFaceLgth=MEDmeshnEntity(fid,mName,dt,it,MED_CELL,MED_POLYHEDRON,MED_INDEX_NODE,MED_NODAL,&changement,&transformation);
184 int curNbOfElem=MEDmeshnEntity(fid,mName,dt,it,MED_CELL,MED_POLYHEDRON,MED_INDEX_FACE,MED_NODAL,&changement,&transformation)-1;
185 _m=MEDCoupling1DGTUMesh::New(mName,INTERP_KERNEL::NORM_POLYHED);
186 MEDCoupling1DGTUMesh *mc(dynamic_cast<MEDCoupling1DGTUMesh *>((MEDCoupling1GTUMesh *)_m));
187 INTERP_KERNEL::AutoPtr<int> index=new int[curNbOfElem+1];
188 INTERP_KERNEL::AutoPtr<int> indexFace=new int[indexFaceLgth];
189 INTERP_KERNEL::AutoPtr<int> locConn=new int[connFaceLgth];
190 MEDmeshPolyhedronRd(fid,mName,dt,it,MED_CELL,MED_NODAL,index,indexFace,locConn);
191 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> conn(DataArrayInt::New()),connI(DataArrayInt::New());
192 int arraySize=connFaceLgth;
193 for(int i=0;i<curNbOfElem;i++)
194 arraySize+=index[i+1]-index[i]-1;
195 conn=DataArrayInt::New();
196 conn->alloc(arraySize,1);
197 int *wFinalConn=conn->getPointer();
198 connI->alloc(curNbOfElem+1,1);
199 int *finalIndex(connI->getPointer());
201 for(int i=0;i<curNbOfElem;i++)
203 finalIndex[i+1]=finalIndex[i]+index[i+1]-index[i]-1+indexFace[index[i+1]-1]-indexFace[index[i]-1];
204 wFinalConn=std::transform(locConn+indexFace[index[i]-1]-1,locConn+indexFace[index[i]]-1,wFinalConn,std::bind2nd(std::plus<int>(),-1));
205 for(int j=index[i];j<index[i+1]-1;j++)
208 wFinalConn=std::transform(locConn+indexFace[j]-1,locConn+indexFace[j+1]-1,wFinalConn,std::bind2nd(std::plus<int>(),-1));
211 mc->setNodalConnectivity(conn,connI);
212 loadCommonPart(fid,mName,dt,it,mdim,curNbOfElem,MED_POLYHEDRON,entity,mrs);
215 void MEDFileUMeshPerType::Write(med_idt fid, const char *mname, int mdim, const MEDCoupling1GTUMesh *m, const DataArrayInt *fam, const DataArrayInt *num, const DataArrayAsciiChar *names)
217 int nbOfCells=m->getNumberOfCells();
221 double timm=m->getTime(dt,it);
222 INTERP_KERNEL::NormalizedCellType ikt=m->getTypeOfCell(0);
223 const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(ikt);
224 med_geometry_type curMedType=typmai3[(int)ikt];
227 const MEDCoupling1SGTUMesh *m0(dynamic_cast<const MEDCoupling1SGTUMesh *>(m));
229 throw INTERP_KERNEL::Exception("MEDFileUMeshPerType::Write : internal error #1 !");
230 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> arr(m0->getNodalConnectivity()->deepCpy());
231 std::transform(arr->begin(),arr->end(),arr->getPointer(),std::bind2nd(std::plus<int>(),1));
232 MEDmeshElementConnectivityWr(fid,mname,dt,it,timm,MED_CELL,curMedType,MED_NODAL,MED_FULL_INTERLACE,nbOfCells,arr->begin());
236 const MEDCoupling1DGTUMesh *m0(dynamic_cast<const MEDCoupling1DGTUMesh *>(m));
238 throw INTERP_KERNEL::Exception("MEDFileUMeshPerType::Write : internal error #2 !");
239 if(ikt==INTERP_KERNEL::NORM_POLYGON)
241 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> arr(m0->getNodalConnectivity()->deepCpy()),arrI(m0->getNodalConnectivityIndex()->deepCpy());
242 std::transform(arr->begin(),arr->end(),arr->getPointer(),std::bind2nd(std::plus<int>(),1));
243 std::transform(arrI->begin(),arrI->end(),arrI->getPointer(),std::bind2nd(std::plus<int>(),1));
244 MEDmeshPolygonWr(fid,mname,dt,it,timm,MED_CELL,MED_NODAL,nbOfCells+1,arrI->begin(),arr->begin());
248 const int *conn(m0->getNodalConnectivity()->begin()),*connI(m0->getNodalConnectivityIndex()->begin());
249 int meshLgth=m0->getNodalConnectivityLength();
250 int nbOfFaces=std::count(conn,conn+meshLgth,-1)+nbOfCells;
251 INTERP_KERNEL::AutoPtr<int> tab1=new int[nbOfCells+1];
253 INTERP_KERNEL::AutoPtr<int> tab2=new int[nbOfFaces+1];
255 INTERP_KERNEL::AutoPtr<int> bigtab=new int[meshLgth];
257 for(int i=0;i<nbOfCells;i++,w1++)
260 for(const int *w=conn+connI[i];w!=conn+connI[i+1];w2++)
262 const int *wend=std::find(w,conn+connI[i+1],-1);
263 bt=std::transform(w,wend,bt,std::bind2nd(std::plus<int>(),1));
264 int nbOfNode=std::distance(w,wend);
265 w2[1]=w2[0]+nbOfNode;
266 if(wend!=conn+connI[i+1])
272 w1[1]=w1[0]+nbOfFaces2;
274 MEDmeshPolyhedronWr(fid,mname,dt,it,timm,MED_CELL,MED_NODAL,nbOfCells+1,tab1,nbOfFaces+1,tab2,bigtab);
278 MEDmeshEntityFamilyNumberWr(fid,mname,dt,it,MED_CELL,curMedType,nbOfCells,fam->getConstPointer());
280 MEDmeshEntityNumberWr(fid,mname,dt,it,MED_CELL,curMedType,nbOfCells,num->getConstPointer());
283 if(names->getNumberOfComponents()!=MED_SNAME_SIZE)
285 std::ostringstream oss; oss << "MEDFileUMeshPerType::write : expected a name field on cells with number of components set to " << MED_SNAME_SIZE;
286 oss << " ! The array has " << names->getNumberOfComponents() << " components !";
287 throw INTERP_KERNEL::Exception(oss.str().c_str());
289 MEDmeshEntityNameWr(fid,mname,dt,it,MED_CELL,curMedType,nbOfCells,names->getConstPointer());