1 // Copyright (C) 2007-2012 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"
23 #include "MEDCouplingUMesh.hxx"
25 #include "InterpKernelAutoPtr.hxx"
26 #include "CellModel.hxx"
30 extern med_geometry_type typmai3[32];
32 using namespace ParaMEDMEM;
34 MEDFileUMeshPerType *MEDFileUMeshPerType::New(med_idt fid, const char *mName, int dt, int it, int mdim, med_geometry_type geoElt, INTERP_KERNEL::NormalizedCellType geoElt2)
36 med_entity_type whichEntity;
37 if(!isExisting(fid,mName,dt,it,geoElt,whichEntity))
39 return new MEDFileUMeshPerType(fid,mName,dt,it,mdim,geoElt,geoElt2,whichEntity);
42 bool MEDFileUMeshPerType::isExisting(med_idt fid, const char *mName, int dt, int it, med_geometry_type geoElt, med_entity_type& whichEntity)
44 static const med_entity_type entities[3]={MED_CELL,MED_DESCENDING_FACE,MED_DESCENDING_EDGE};
48 med_bool changement,transformation;
49 int tmp=MEDmeshnEntity(fid,mName,dt,it,entities[i],geoElt,MED_CONNECTIVITY,MED_NODAL,
50 &changement,&transformation);
54 whichEntity=entities[i];
56 std::cerr << "WARNING : MEDFile has been detected to be no compilant with MED 3 : Please change entity in MEDFile for geotype " << geoElt << std::endl;
62 int MEDFileUMeshPerType::getDim() const
64 const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(_type);
65 return cm.getDimension();
68 MEDFileUMeshPerType::MEDFileUMeshPerType(med_idt fid, const char *mName, int dt, int it, int mdim, med_geometry_type geoElt, INTERP_KERNEL::NormalizedCellType type,
69 med_entity_type entity):_type(type),_entity(entity)
71 med_bool changement,transformation;
72 int curNbOfElem=MEDmeshnEntity(fid,mName,dt,it,entity,geoElt,MED_CONNECTIVITY,MED_NODAL,
73 &changement,&transformation);
74 if(type!=INTERP_KERNEL::NORM_POLYGON && type!=INTERP_KERNEL::NORM_POLYHED)
76 loadFromStaticType(fid,mName,dt,it,mdim,curNbOfElem,geoElt,type,entity);
79 if(type==INTERP_KERNEL::NORM_POLYGON)
81 loadPolyg(fid,mName,dt,it,mdim,curNbOfElem,geoElt,entity);
84 //if(type==INTERP_KERNEL::NORM_POLYHED)
85 loadPolyh(fid,mName,dt,it,mdim,curNbOfElem,geoElt,entity);
88 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,
89 med_entity_type entity)
91 _conn=DataArrayInt::New();
92 int nbOfNodesPerCell=(geoElt%100);
93 _conn->alloc((nbOfNodesPerCell+1)*curNbOfElem,1);
94 _conn_index=DataArrayInt::New();
95 _conn_index->alloc(curNbOfElem+1,1);
96 INTERP_KERNEL::AutoPtr<int> connTab=new int[(nbOfNodesPerCell)*curNbOfElem];
97 _num=DataArrayInt::New();
98 _num->alloc(curNbOfElem,1);
99 _fam=DataArrayInt::New();
100 _fam->alloc(curNbOfElem,1);
101 med_bool changement,transformation;
102 INTERP_KERNEL::AutoPtr<char> noms=new char[MED_SNAME_SIZE*(std::size_t)curNbOfElem+1];
103 MEDmeshElementConnectivityRd(fid,mName,dt,it,entity,geoElt,MED_NODAL,MED_FULL_INTERLACE,connTab);
104 if(MEDmeshnEntity(fid,mName,dt,it,entity,geoElt,MED_FAMILY_NUMBER,MED_NODAL,&changement,&transformation)>0)
106 if(MEDmeshEntityFamilyNumberRd(fid,mName,dt,it,entity,geoElt,_fam->getPointer())!=0)
107 std::fill(_fam->getPointer(),_fam->getPointer()+curNbOfElem,0);
110 std::fill(_fam->getPointer(),_fam->getPointer()+curNbOfElem,0);
111 if(MEDmeshnEntity(fid,mName,dt,it,entity,geoElt,MED_NUMBER,MED_NODAL,&changement,&transformation)>0)
113 _num=DataArrayInt::New();
114 _num->alloc(curNbOfElem,1);
115 if(MEDmeshEntityNumberRd(fid,mName,dt,it,entity,geoElt,_num->getPointer())!=0)
120 int *w1=_conn->getPointer();
121 int *w2=_conn_index->getPointer();
123 const int *wi=connTab;
124 for(int i=0;i<curNbOfElem;i++,wi+=nbOfNodesPerCell,w2++)
127 w1=std::transform(wi,wi+nbOfNodesPerCell,w1,std::bind2nd(std::plus<int>(),-1));
128 *w2=w2[-1]+nbOfNodesPerCell+1;
132 void MEDFileUMeshPerType::loadPolyg(med_idt fid, const char *mName, int dt, int it, int mdim, int arraySize, med_geometry_type geoElt,
133 med_entity_type entity)
135 med_bool changement,transformation;
136 med_int curNbOfElem=MEDmeshnEntity(fid,mName,dt,it,entity,MED_POLYGON,MED_INDEX_NODE,MED_NODAL,&changement,&transformation)-1;
137 _conn_index=DataArrayInt::New();
138 _conn_index->alloc(curNbOfElem+1,1);
139 _conn=DataArrayInt::New();
140 _conn->alloc(arraySize+curNbOfElem,1);
141 _num=DataArrayInt::New();
142 _num->alloc(curNbOfElem,1);
143 _fam=DataArrayInt::New();
144 _fam->alloc(curNbOfElem,1);
145 INTERP_KERNEL::AutoPtr<int> locConn=new int[arraySize];
146 MEDmeshPolygonRd(fid,mName,dt,it,MED_CELL,MED_NODAL,_conn_index->getPointer(),locConn);
147 int *w1=_conn->getPointer();
148 int *w2=_conn_index->getPointer();
149 const int *wi=locConn;
150 for(int i=0;i<curNbOfElem;i++,w2++)
152 *w1++=(int)INTERP_KERNEL::NORM_POLYGON;
153 const int *wi2=wi+(w2[1]-w2[0]);
154 w1=std::transform(wi,wi2,w1,std::bind2nd(std::plus<int>(),-1));
158 *w2=*w2-1+curNbOfElem;
159 if(MEDmeshnEntity(fid,mName,dt,it,entity,MED_POLYGON,MED_FAMILY_NUMBER,MED_NODAL,&changement,&transformation)>0)
161 if(MEDmeshEntityFamilyNumberRd(fid,mName,dt,it,entity,MED_POLYGON,_fam->getPointer())!=0)
162 std::fill(_fam->getPointer(),_fam->getPointer()+curNbOfElem,0);
165 std::fill(_fam->getPointer(),_fam->getPointer()+curNbOfElem,0);
166 if(MEDmeshnEntity(fid,mName,dt,it,entity,MED_POLYGON,MED_NUMBER,MED_NODAL,&changement,&transformation)>0)
168 if(MEDmeshEntityNumberRd(fid,mName,dt,it,entity,MED_POLYGON,_num->getPointer())!=0)
175 void MEDFileUMeshPerType::loadPolyh(med_idt fid, const char *mName, int dt, int it, int mdim, int connFaceLgth, med_geometry_type geoElt,
176 med_entity_type entity)
178 med_bool changement,transformation;
179 med_int indexFaceLgth=MEDmeshnEntity(fid,mName,dt,it,MED_CELL,MED_POLYHEDRON,MED_INDEX_NODE,MED_NODAL,&changement,&transformation);
180 int curNbOfElem=MEDmeshnEntity(fid,mName,dt,it,MED_CELL,MED_POLYHEDRON,MED_INDEX_FACE,MED_NODAL,&changement,&transformation)-1;
181 INTERP_KERNEL::AutoPtr<int> index=new int[curNbOfElem+1];
182 INTERP_KERNEL::AutoPtr<int> indexFace=new int[indexFaceLgth];
183 INTERP_KERNEL::AutoPtr<int> locConn=new int[connFaceLgth];
184 _fam=DataArrayInt::New();
185 _fam->alloc(curNbOfElem,1);
186 MEDmeshPolyhedronRd(fid,mName,dt,it,MED_CELL,MED_NODAL,index,indexFace,locConn);
187 if(MEDmeshnEntity(fid,mName,dt,it,entity,MED_POLYHEDRON,MED_FAMILY_NUMBER,MED_NODAL,&changement,&transformation)>0)
189 if(MEDmeshEntityFamilyNumberRd(fid,mName,dt,it,entity,MED_POLYHEDRON,_fam->getPointer())!=0)
190 std::fill(_fam->getPointer(),_fam->getPointer()+curNbOfElem,0);
193 std::fill(_fam->getPointer(),_fam->getPointer()+curNbOfElem,0);
194 int arraySize=connFaceLgth;
195 for(int i=0;i<curNbOfElem;i++)
196 arraySize+=index[i+1]-index[i]-1;
197 _conn=DataArrayInt::New();
198 _conn->alloc(arraySize+curNbOfElem,1);
199 int *wFinalConn=_conn->getPointer();
200 _conn_index=DataArrayInt::New();
201 _conn_index->alloc(curNbOfElem+1,1);
202 int *finalIndex=_conn_index->getPointer();
204 for(int i=0;i<curNbOfElem;i++)
206 *wFinalConn++=(int)INTERP_KERNEL::NORM_POLYHED;
207 finalIndex[i+1]=finalIndex[i]+index[i+1]-index[i]-1+indexFace[index[i+1]-1]-indexFace[index[i]-1]+1;
208 wFinalConn=std::transform(locConn+indexFace[index[i]-1]-1,locConn+indexFace[index[i]]-1,wFinalConn,std::bind2nd(std::plus<int>(),-1));
209 for(int j=index[i];j<index[i+1]-1;j++)
212 wFinalConn=std::transform(locConn+indexFace[j]-1,locConn+indexFace[j+1]-1,wFinalConn,std::bind2nd(std::plus<int>(),-1));
215 _num=DataArrayInt::New();
216 _num->alloc(curNbOfElem,1);
217 if(MEDmeshnEntity(fid,mName,dt,it,entity,MED_POLYHEDRON,MED_NUMBER,MED_NODAL,&changement,&transformation)>0)
219 if(MEDmeshEntityNumberRd(fid,mName,dt,it,entity,MED_POLYHEDRON,_num->getPointer())!=0)
226 void MEDFileUMeshPerType::write(med_idt fid, const char *mname, int mdim, const MEDCouplingUMesh *m, const DataArrayInt *fam, const DataArrayInt *num)
228 int nbOfCells=m->getNumberOfCells();
232 double timm=m->getTime(dt,it);
233 INTERP_KERNEL::NormalizedCellType ikt=m->getTypeOfCell(0);
234 const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(ikt);
235 med_geometry_type curMedType=typmai3[(int)ikt];
236 const int *conn=m->getNodalConnectivity()->getConstPointer();
237 const int *connI=m->getNodalConnectivityIndex()->getConstPointer();
238 if(ikt!=INTERP_KERNEL::NORM_POLYGON && ikt!=INTERP_KERNEL::NORM_POLYHED)
240 int nbNodesPerCell=cm.getNumberOfNodes();
241 INTERP_KERNEL::AutoPtr<int> tab=new int[nbNodesPerCell*nbOfCells];
243 for(int i=0;i<nbOfCells;i++)
244 w=std::transform(conn+connI[i]+1,conn+connI[i+1],w,std::bind2nd(std::plus<int>(),1));
245 MEDmeshElementConnectivityWr(fid,mname,dt,it,timm,MED_CELL,curMedType,MED_NODAL,MED_FULL_INTERLACE,nbOfCells,tab);
249 if(ikt==INTERP_KERNEL::NORM_POLYGON)
251 INTERP_KERNEL::AutoPtr<int> tab1=new int[nbOfCells+1];
252 INTERP_KERNEL::AutoPtr<int> tab2=new int[m->getMeshLength()];
255 for(int i=0;i<nbOfCells;i++,wI++)
257 wI[1]=wI[0]+connI[i+1]-connI[i]-1;
258 w=std::transform(conn+connI[i]+1,conn+connI[i+1],w,std::bind2nd(std::plus<int>(),1));
260 MEDmeshPolygonWr(fid,mname,dt,it,timm,MED_CELL,MED_NODAL,nbOfCells+1,tab1,tab2);
264 int meshLgth=m->getMeshLength();
265 int nbOfFaces=std::count(conn,conn+meshLgth,-1)+nbOfCells;
266 INTERP_KERNEL::AutoPtr<int> tab1=new int[nbOfCells+1];
268 INTERP_KERNEL::AutoPtr<int> tab2=new int[nbOfFaces+1];
270 INTERP_KERNEL::AutoPtr<int> bigtab=new int[meshLgth-nbOfCells];
272 for(int i=0;i<nbOfCells;i++,w1++)
275 for(const int *w=conn+connI[i]+1;w!=conn+connI[i+1];w2++)
277 const int *wend=std::find(w,conn+connI[i+1],-1);
278 bt=std::transform(w,wend,bt,std::bind2nd(std::plus<int>(),1));
279 int nbOfNode=std::distance(w,wend);
280 w2[1]=w2[0]+nbOfNode;
281 if(wend!=conn+connI[i+1])
287 w1[1]=w1[0]+nbOfFaces2;
289 MEDmeshPolyhedronWr(fid,mname,dt,it,timm,MED_CELL,MED_NODAL,nbOfCells+1,tab1,nbOfFaces+1,tab2,bigtab);
293 MEDmeshEntityFamilyNumberWr(fid,mname,dt,it,MED_CELL,curMedType,nbOfCells,fam->getConstPointer());
295 MEDmeshEntityNumberWr(fid,mname,dt,it,MED_CELL,curMedType,nbOfCells,num->getConstPointer());