1 // Copyright (C) 2007-2023 CEA, EDF
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, or (at your option) any later version.
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 (EDF R&D)
21 #include "ParaMEDFileMesh.hxx"
23 #include "MEDFileMesh.hxx"
24 #include "MEDFileMeshLL.hxx"
25 #include "MEDLoader.hxx"
26 #include "MEDFileField1TS.hxx"
27 #include "MEDFileUtilities.hxx"
28 #include "MEDFileEntities.hxx"
33 // From MEDLOader.cxx TU
34 extern med_geometry_type typmai3[INTERP_KERNEL::NORM_MAXTYPE];
36 using namespace MEDCoupling;
38 void getSingleGeometricType(const std::string& fileName, const std::string& mName, INTERP_KERNEL::NormalizedCellType& geoType)
40 int meshDim, spaceDim;
41 mcIdType numberOfNodes;
42 std::vector< std::vector< std::pair<INTERP_KERNEL::NormalizedCellType,int> > > typesDistrib(GetUMeshGlobalInfo(fileName,mName,meshDim,spaceDim,numberOfNodes));
43 std::size_t numberOfTypesMaxDimension = typesDistrib[0].size();
44 if(numberOfTypesMaxDimension != 1)
45 throw INTERP_KERNEL::Exception("ParaMEDFileMesh : only mesh with single geometrical type are supported with given distribution !");
46 geoType = typesDistrib[0][0].first;
49 void checkDistribution(const MPI_Comm& com, mcIdType totalNumberOfElements, const std::vector<mcIdType>& distrib)
51 mcIdType nbEltsInDistribLoc = distrib.size();
52 mcIdType nbEltsInDistribTot = -1;
54 MPI_Allreduce(&nbEltsInDistribLoc, &nbEltsInDistribTot, 1, MPI_LONG, MPI_SUM, com);
56 throw INTERP_KERNEL::Exception("not(HAVE_MPI) incompatible with MPI_World_Size>1");
58 if(nbEltsInDistribTot != totalNumberOfElements)
60 if(nbEltsInDistribTot > totalNumberOfElements)
61 throw INTERP_KERNEL::Exception("ParaMEDFileMesh : Some of your partitions overlap each other ! Each element in your distribution vector must appear only once ! ");
63 throw INTERP_KERNEL::Exception("ParaMEDFileMesh : The distribution does not cover the whole mesh ! Each element of the mesh must appear once in your distribution vector ");
68 MEDFileMesh *ParaMEDFileMesh::New(int iPart, int nbOfParts, const std::string& fileName, const std::string& mName, int dt, int it, MEDFileMeshReadSelector *mrs)
70 MEDFileUtilities::CheckFileForRead(fileName);
71 MEDCoupling::MEDCouplingMeshType meshType;
72 MEDFileUtilities::AutoFid fid(MEDfileOpen(fileName.c_str(),MED_ACC_RDONLY));
75 MEDCoupling::MEDCouplingAxisType dummy3;
76 MEDFileMeshL2::GetMeshIdFromName(fid,mName,meshType,dummy3,dummy0,dummy1,dummy2);
81 return ParaMEDFileUMesh::New(iPart,nbOfParts,fileName,mName,dt,it,mrs);
84 throw INTERP_KERNEL::Exception("ParaMEDFileMesh::New : only unstructured mesh supported for the moment !");
88 MEDFileMesh *ParaMEDFileMesh::ParaNew(int iPart, int nbOfParts, const MPI_Comm& com, const MPI_Info& nfo, const std::string& fileName, const std::string& mName, int dt, int it, MEDFileMeshReadSelector *mrs)
90 MEDFileUtilities::CheckFileForRead(fileName);
91 MEDCoupling::MEDCouplingMeshType meshType;
92 MEDFileUtilities::AutoFid fid(MEDfileOpen(fileName.c_str(),MED_ACC_RDONLY));
95 MEDCoupling::MEDCouplingAxisType dummy3;
96 MEDFileMeshL2::GetMeshIdFromName(fid,mName,meshType,dummy3,dummy0,dummy1,dummy2);
101 return ParaMEDFileUMesh::ParaNew(iPart,nbOfParts,com,nfo,fileName,mName,dt,it,mrs);
104 throw INTERP_KERNEL::Exception("ParaMEDFileMesh::ParaNew : only unstructured mesh supported for the moment !");
108 MEDFileUMesh *ParaMEDFileUMesh::New(int iPart, int nbOfParts, const std::string& fileName, const std::string& mName, int dt, int it, MEDFileMeshReadSelector *mrs)
110 MEDFileUtilities::CheckFileForRead(fileName);
111 MEDFileUtilities::AutoFid fid(MEDfileOpen(fileName.c_str(),MED_ACC_RDONLY));
112 return ParaMEDFileUMesh::NewPrivate(fid,iPart,nbOfParts,fileName,mName,dt,it,mrs);
116 * Opens the given file in parallel so that each processor can load a specific part of the mesh \a mName.
117 * Each processor will load the cells contained in the vector \a distrib (only nodes lying on those cells will be loaded),
118 * in order to read the entire mesh in parallel (memory consumption is thus distributed among all the processes).
119 * \param [in] distrib - map defining for each geometric type, the corresponding vector of cells we want to load with c-type indexing (starting from zero).
120 * Each vector in this map has an independant numerotation, which means on one processor, vectors of different types may contain the same numbers, they will not refer to the same cells
121 * (the i-th cell of a type A does not correspond to the i-th cell of type B)
122 * However they have to differ from one processor to another, as to ensure that:
123 * 1) each processor only loads a unique part of the mesh
124 * 2) the combined distribution vectors cover the entire mesh
125 * \param [in] com - group of MPI processes that will read the file
126 * \param [in] nfo- MPI info object (used to manage MPI routines)
127 * \param [in] filename - name of the file we want to read
128 * \param [in] mName - name of the mesh we want to read
129 * \param [in] dt - order at which to read the mesh
130 * \param [in] it - iteration at which to read the mesh
131 * \param [in] mrs - object used to read additional low-level information
132 * \return MEDFileUMesh* - a new instance of MEDFileUMesh. The
133 * caller is to delete this mesh using decrRef() as it is no more needed.
134 * \throw exception if the mesh contains multiple types of cells
135 * \throw exception if the partition of the mesh cells defined by \a distrib does not cover the whole mesh
137 MEDFileUMesh *ParaMEDFileUMesh::ParaNew(const std::map<INTERP_KERNEL::NormalizedCellType,std::vector<mcIdType>> &distrib, const MPI_Comm& com, const MPI_Info& nfo, const std::string& fileName, const std::string& mName, int dt, int it, MEDFileMeshReadSelector *mrs)
139 MEDFileUtilities::CheckFileForRead(fileName);
140 #ifdef HDF5_IS_PARALLEL
141 MEDFileUtilities::AutoFid fid(MEDparFileOpen(fileName.c_str(),MED_ACC_RDONLY,com,nfo));
143 MEDFileUtilities::AutoFid fid(MEDfileOpen(fileName.c_str(),MED_ACC_RDONLY));
145 return ParaMEDFileUMesh::NewPrivate(fid,com,distrib,fileName,mName,dt,it,mrs);
150 * Opens the given file in parallel so that each processor can load its part of the mesh \a mName.
151 * The mesh will be equally and linearly distributed among all processes:
152 * the list of cells will be divided into \a nbOfParts slices and only slice \a iPart (cells and nodes lying on those cells) will be loaded by the current processor.
153 * The entire mesh is thus read in parallel and memory consumption is divided among the group of processes.
154 * \param [in] iPart - part of the mesh that will be loaded
155 * \param [in] nbOfParts - total number of parts in which to divide the mesh
156 * \param [in] com - group of MPI processes that will read the file
157 * \param [in] nfo- MPI info object (used to manage MPI routines)
158 * \param [in] filename - name of the file we want to read
159 * \param [in] mName - name of the mesh we want to read
160 * \param [in] dt - Time order at which to read the mesh
161 * \param [in] it - Time iteration at which to read the mesh
162 * \param [in] mrs - object used to read additional low-level information
163 * \return MEDFileUMesh* - a new instance of MEDFileUMesh. The
164 * caller is to delete this mesh using decrRef() as it is no more needed.
166 MEDFileUMesh *ParaMEDFileUMesh::ParaNew(int iPart, int nbOfParts, const MPI_Comm& com, const MPI_Info& nfo, const std::string& fileName, const std::string& mName, int dt, int it, MEDFileMeshReadSelector *mrs)
168 MEDFileUtilities::CheckFileForRead(fileName);
169 #ifdef HDF5_IS_PARALLEL
170 MEDFileUtilities::AutoFid fid(MEDparFileOpen(fileName.c_str(),MED_ACC_RDONLY,com,nfo)); // MPI_COMM_WORLD, MPI_INFO_NULL
172 MEDFileUtilities::AutoFid fid(MEDfileOpen(fileName.c_str(),MED_ACC_RDONLY));
174 return ParaMEDFileUMesh::NewPrivate(fid,iPart,nbOfParts,fileName,mName,dt,it,mrs);
178 * Loads mesh \a mName in parallel using a custom partition of the mesh cells among the processes.
179 * See ParaMEDFileUMesh::ParaNew for detailed description.
181 MEDFileUMesh *ParaMEDFileUMesh::NewPrivate(med_idt fid, const MPI_Comm& com, const std::map<INTERP_KERNEL::NormalizedCellType,std::vector<mcIdType>>& distrib, const std::string& fileName, const std::string& mName, int dt, int it, MEDFileMeshReadSelector *mrs)
183 MCAuto<MEDFileUMesh> ret;
184 for(std::map<INTERP_KERNEL::NormalizedCellType,std::vector<mcIdType>>::const_iterator iter=distrib.begin(); iter!= distrib.end(); iter++)
186 med_geometry_type geoMedType(typmai3[iter->first /*current geometric type*/]);
187 med_bool changement,transformation;
188 med_int totalNumberOfElements(MEDmeshnEntity(fid,mName.c_str(),dt,it,MED_CELL,geoMedType,MED_CONNECTIVITY,MED_NODAL,&changement,&transformation));
189 checkDistribution(com,totalNumberOfElements,iter->second /*distrib over this geometric type*/);
191 ret=MEDFileUMesh::LoadPartOfFromUserDistrib(fid,mName,distrib,dt,it,mrs);
196 * Loads mesh \a mName in parallel using a slice partition of the mesh cells among the processes
197 * See ParaMEDFileUMesh::ParaNew for detailed description.
199 MEDFileUMesh *ParaMEDFileUMesh::NewPrivate(med_idt fid, int iPart, int nbOfParts, const std::string& fileName, const std::string& mName, int dt, int it, MEDFileMeshReadSelector *mrs)
201 MCAuto<MEDFileUMesh> ret;
202 int meshDim, spaceDim;
203 mcIdType numberOfNodes;
204 std::vector< std::vector< std::pair<INTERP_KERNEL::NormalizedCellType,int> > > typesDistrib(GetUMeshGlobalInfo(fileName,mName,meshDim,spaceDim,numberOfNodes));
205 std::vector<INTERP_KERNEL::NormalizedCellType> types;
206 std::vector<mcIdType> distrib;
207 for(std::vector< std::vector< std::pair<INTERP_KERNEL::NormalizedCellType,int> > >::const_iterator it0=typesDistrib.begin();it0!=typesDistrib.end();it0++)
208 for(std::vector< std::pair<INTERP_KERNEL::NormalizedCellType,int> >::const_iterator it1=(*it0).begin();it1!=(*it0).end();it1++)
210 types.push_back((*it1).first);
212 DataArray::GetSlice(0,(*it1).second,1,iPart,nbOfParts,tmp[0],tmp[1]);
214 distrib.insert(distrib.end(),tmp,tmp+3);
216 ret=MEDFileUMesh::LoadPartOf(fid,mName,types,distrib,dt,it,mrs);
221 * Loads field \a fName laying on mesh \a mName from the filename \a fileName in parallel:
222 * each processor will load their portion of the field (ie the portion laying on the cells in the vector \a distrib given in the parameters).
223 * WARNING : this will only load the array of values of the field, additionnal information of the local field such as the number of its tuples might be incorrect
224 * \param [in] com - group of MPI processes that will read the file
225 * \param [in] nfo- MPI info object (used to manage MPI routines)
226 * \param [in] fileName - name of the file containing the field
227 * \param [in] fName - name of the field we want to load
228 * \param [in] mName - name of the mesh on which the field is defined
229 * \param [in] distrib - vector of cells on which we want to load the field \a fName (with c-type indexing, so starting from zero).
230 * \param [in] dt - Time order at which to read the field
231 * \param [in] it - Time iteration at which to read the field
232 * \return MEDFileField1TS* - a new instance of MEDFileField1TS. The
233 * caller is to delete it using decrRef() as it is no more needed.
234 * \throw exception if the field is not of type FLOAT64
235 * \throw exception if the mesh contains more than one geometric type
236 * \throw exception if the given distribution does not cover the entire mesh on which the field is defined
238 MEDFileField1TS *ParaMEDFileField1TS::ParaNew(const MPI_Comm& com, const MPI_Info& nfo, const std::string& fileName, const std::string& fName, const std::string& mName, const std::vector<mcIdType>& distrib, TypeOfField loc, int dt, int it)
240 MEDFileUtilities::CheckFileForRead(fileName);
241 #ifdef HDF5_IS_PARALLEL
242 MEDFileUtilities::AutoFid fid(MEDparFileOpen(fileName.c_str(),MED_ACC_RDONLY,com,nfo));
244 MEDFileUtilities::AutoFid fid(MEDfileOpen(fileName.c_str(),MED_ACC_RDONLY));
246 return ParaMEDFileField1TS::NewPrivate(fid,com,fName,mName,distrib,loc,dt,it);
250 * Loads field \a fName in parallel using a custom partition of the mesh cells on which the field is defined among the processes.
251 * See ParaMEDFileField1TS::ParaNew for detailed description.
253 MEDFileField1TS *ParaMEDFileField1TS::NewPrivate(med_idt fid, const MPI_Comm& com, const std::string& fName, const std::string& mName, const std::vector<mcIdType>& distrib, TypeOfField loc, int dt, int it)
255 INTERP_KERNEL::NormalizedCellType geoType;
256 getSingleGeometricType(MEDFileWritable::FileNameFromFID(fid),mName,geoType);
257 if(loc==ON_CELLS) //if distribution is on nodes, no fast way to check it (as a node can be shared by multiple processors)
259 med_geometry_type geoMedType(typmai3[geoType]);
260 med_bool changement,transformation;
261 med_int totalNumberOfElements(MEDmeshnEntity(fid,mName.c_str(),dt,it,MED_CELL,geoMedType,MED_FAMILY_NUMBER,MED_NODAL,&changement,&transformation));
262 checkDistribution(com,totalNumberOfElements,distrib);
264 std::vector<std::pair<TypeOfField,INTERP_KERNEL::NormalizedCellType>> tmp={ {loc, geoType} };
265 INTERP_KERNEL::AutoCppPtr<MEDFileEntities> entities(MEDFileEntities::BuildFrom(&tmp));
266 MCAuto<MEDFileAnyTypeField1TS> partFile(MEDFileAnyTypeField1TS::NewAdv(fid,fName,dt,it,entities,distrib));
268 MCAuto<MEDFileField1TS> ret(MEDCoupling::DynamicCast<MEDFileAnyTypeField1TS,MEDFileField1TS>(partFile));
272 throw INTERP_KERNEL::Exception("ParaMEDFileField1TS::ParaNew : only FLOAT64 field supported for the moment !");
276 MEDFileMeshes *ParaMEDFileMeshes::New(int iPart, int nbOfParts, const std::string& fileName)
278 std::vector<std::string> ms(GetMeshNames(fileName));
279 MCAuto<MEDFileMeshes> ret(MEDFileMeshes::New());
280 for(std::vector<std::string>::const_iterator it=ms.begin();it!=ms.end();it++)
282 MCAuto<MEDFileMesh> mesh(ParaMEDFileMesh::New(iPart,nbOfParts,fileName,(*it)));
288 MEDFileMeshes *ParaMEDFileMeshes::ParaNew(int iPart, int nbOfParts, const MPI_Comm& com, const MPI_Info& nfo, const std::string& fileName)
290 std::vector<std::string> ms(GetMeshNames(fileName));
291 MCAuto<MEDFileMeshes> ret(MEDFileMeshes::New());
292 for(std::vector<std::string>::const_iterator it=ms.begin();it!=ms.end();it++)
294 MCAuto<MEDFileMesh> mesh(ParaMEDFileMesh::ParaNew(iPart,nbOfParts,com,nfo,fileName,(*it)));