14#if (MED_MAJOR_NUM >= 3)
19#define med_geometrie_element med_geometry_type
20#define med_entite_maillage med_entity_type
21#define med_type_champ med_field_type
22#define MED_LECTURE MED_ACC_RDONLY
23#define MED_LECTURE_AJOUT MED_ACC_RDEXT
24#define MED_NOEUD MED_NODE
25#define MED_MAILLE MED_CELL
26#define MED_NOEUD_MAILLE MED_NODE_ELEMENT
27#define MED_NOPFL MED_NO_PROFILE
28#define MEDouvrir MEDfileOpen
29#define MEDfermer MEDfileClose
30#define MEDnChamp MEDfieldnComponent
31#define MEDnValProfil MEDprofileSizeByName
33#error "MED file is not MED_MAJOR_NUM == 3"
48 : cOre(const_cast<
Core &>(core)), flgFile(PETSC_FALSE) {
51 auto core_log = logging::core::get();
58 MOFEM_LOG(
"MEDWORLD", Sev::noisy) <<
"Mashset manager created";
65 ierr = PetscOptionsBegin(m_field.
get_comm(),
"",
"MED Interface",
"none");
67 ierr = PetscOptionsString(
"-med_file",
"med file name",
"",
"mesh.med",
70 ierr = PetscOptionsEnd();
83 med_idt fid = MEDfileOpen(file.c_str(), MED_ACC_RDONLY);
86 "Unable to open file '%s'", file.c_str());
88 med_int num_fields = MEDnField(fid);
89 for (
int index = 0; index < num_fields; index++) {
91 med_int num_comp = MEDfieldnComponent(fid, index + 1);
94 "Could not get number of components for MED field");
97 char name[MED_NAME_SIZE + 1], mesh_name[MED_NAME_SIZE + 1];
98 char dt_unit[MED_SNAME_SIZE + 1];
99 std::vector<char> comp_name(num_comp * MED_SNAME_SIZE + 1);
100 std::vector<char> comp_unit(num_comp * MED_SNAME_SIZE + 1);
101 med_int num_steps = 0;
104 if (MEDfieldInfo(fid, index + 1, name, mesh_name, &local_mesh, &type,
105 &comp_name[0], &comp_unit[0], dt_unit, &num_steps) < 0) {
107 "Could not get MED field info");
115 for (
int ff = 0; ff != num_comp; ff++) {
117 std::string(&comp_name[ff * MED_SNAME_SIZE], MED_SNAME_SIZE));
119 std::string(&comp_unit[ff * MED_SNAME_SIZE], MED_SNAME_SIZE));
127 if (MEDfileClose(fid) < 0) {
129 "Unable to close file '%s'", file.c_str());
147 med_idt fid = MEDfileOpen(file.c_str(), MED_ACC_RDONLY);
150 "Unable to open file '%s'", file.c_str());
154 MEDlibraryNumVersion(&
v[0], &
v[1], &
v[2]);
155 MEDfileNumVersionRd(fid, &vf[0], &vf[1], &vf[2]);
159 "Reading MED file V%d.%d.%d using MED library V%d.%d.%d", vf[0],
160 vf[1], vf[2],
v[0],
v[1],
v[2]);
162 if (vf[0] < 2 || (vf[0] == 2 && vf[1] < 2)) {
164 "Cannot read MED file older than V2.2");
167 for (
int i = 0;
i < MEDnMesh(fid);
i++) {
168 char mesh_name[MED_NAME_SIZE + 1], mesh_desc[MED_COMMENT_SIZE + 1];
169 bzero(mesh_name, MED_NAME_SIZE);
170 bzero(mesh_desc, MED_COMMENT_SIZE);
172 med_mesh_type mesh_type;
173 med_int mesh_dim, n_step;
174 char dt_unit[MED_SNAME_SIZE + 1];
175 char axis_name[3 * MED_SNAME_SIZE + 1], axis_unit[3 * MED_SNAME_SIZE + 1];
176 med_sorting_type sorting_type;
177 med_axis_type axis_type;
178 if (MEDmeshInfo(fid,
i + 1, mesh_name, &space_dim, &mesh_dim, &mesh_type,
179 mesh_desc, dt_unit, &sorting_type, &n_step, &axis_type,
180 axis_name, axis_unit) < 0) {
182 "Unable to read mesh information");
184 meshNames.push_back(std::string(mesh_name));
186 MOFEM_LOG_C(
"MEDWORLD", Sev::inform,
"Check mesh %s nsteps %d", mesh_name,
191 std::map<int, Range> family_elem_map;
192 std::map<string, Range> group_elem_map;
194 for (
unsigned int ii = 0; ii !=
meshNames.size(); ii++) {
200 if (MEDfileClose(fid) < 0) {
202 "Unable to close file '%s'", file.c_str());
208static std::vector<med_geometrie_element>
211 std::vector<med_geometrie_element> types;
215 types.push_back(MED_SEG2);
216 types.push_back(MED_SEG3);
219 types.push_back(MED_TRIA3);
220 types.push_back(MED_TRIA6);
223 types.push_back(MED_QUAD4);
226 types.push_back(MED_TETRA4);
227 types.push_back(MED_TETRA10);
230 types.push_back(MED_HEXA8);
233 types.push_back(MED_PENTA6);
236 types.push_back(MED_PYRA5);
239 types.push_back(MED_POINT1);
248 std::map<int, Range> &family_elem_map,
254 med_idt fid = MEDfileOpen(file.c_str(), MED_ACC_RDONLY);
257 "Unable to open file '%s'", file.c_str());
260 MEDlibraryNumVersion(&
v[0], &
v[1], &
v[2]);
261 MEDfileNumVersionRd(fid, &vf[0], &vf[1], &vf[2]);
264 "Reading MED file V%d.%d.%d using MED library V%d.%d.%d", vf[0],
265 vf[1], vf[2],
v[0],
v[1],
v[2]);
267 if (vf[0] < 2 || (vf[0] == 2 && vf[1] < 2)) {
269 "Cannot read MED file older than V2.2");
272 char mesh_name[MED_NAME_SIZE + 1], mesh_desc[MED_COMMENT_SIZE + 1];
273 bzero(mesh_name, MED_NAME_SIZE + 1);
274 bzero(mesh_desc, MED_COMMENT_SIZE + 1);
276 med_mesh_type mesh_type;
277 med_int mesh_dim, n_step;
278 char dt_unit[MED_SNAME_SIZE + 1];
279 char axis_name[3 * MED_SNAME_SIZE + 1], axis_unit[3 * MED_SNAME_SIZE + 1];
280 med_sorting_type sorting_type;
281 med_axis_type axis_type;
282 if (MEDmeshInfo(fid, index + 1, mesh_name, &space_dim, &mesh_dim, &mesh_type,
283 mesh_desc, dt_unit, &sorting_type, &n_step, &axis_type,
284 axis_name, axis_unit) < 0) {
286 "Unable to read mesh information");
289 MOFEM_LOG_C(
"MEDWORLD", Sev::inform,
"Reading mesh %s nsteps %d", mesh_name,
296 case MED_CYLINDRICAL:
299 "Curvilinear coordinate system implemented");
303 "Not implemented for space dim %d", space_dim);
312 max_id = (max_id < cit->getMeshsetId()) ? cit->getMeshsetId() : max_id;
316 std::string(mesh_name));
317 CubitMeshSet_multiIndex::index<
324 mesh_meshset = cit->getMeshset();
328 med_bool change_of_coord, geo_transform;
329 med_int num_nodes = MEDmeshnEntity(
330 fid, mesh_name, MED_NO_DT, MED_NO_IT, MED_NODE, MED_NO_GEOTYPE,
331 MED_COORDINATE, MED_NO_CMODE, &change_of_coord, &geo_transform);
334 "Could not read number of MED nodes");
336 if (num_nodes == 0) {
338 "No nodes in MED mesh");
341 MOFEM_LOG_C(
"MEDWORLD", Sev::inform,
"Read number of nodes %d", num_nodes);
343 std::vector<med_float> coord_med(space_dim * num_nodes);
344 if (MEDmeshNodeCoordinateRd(fid, mesh_name, MED_NO_DT, MED_NO_IT,
345 MED_NO_INTERLACE, &coord_med[0]) < 0) {
347 "Could not read MED node coordinates");
351 ReadUtilIface *iface;
352 vector<double *> arrays_coord;
355 CHKERR iface->get_node_coords(3, num_nodes, 0, startv, arrays_coord);
356 Range verts(startv, startv + num_nodes - 1);
357 std::copy(&coord_med[0 * num_nodes], &coord_med[1 * num_nodes],
359 std::copy(&coord_med[1 * num_nodes], &coord_med[2 * num_nodes],
361 if (space_dim == 2) {
362 std::fill(arrays_coord[2], &arrays_coord[2][num_nodes], 0.);
364 std::copy(&coord_med[2 * num_nodes], &coord_med[3 * num_nodes],
368 family_elem_map.clear();
372 std::vector<med_int> nodes_tags(num_nodes, 0);
373 if (MEDmeshEntityFamilyNumberRd(fid, mesh_name, MED_NO_DT, MED_NO_IT,
374 MED_NODE, MED_NO_GEOTYPE,
375 &nodes_tags[0]) < 0) {
383 for (
int i = 0;
i < num_nodes;
i++) {
386 family_elem_map[nodes_tags.empty() ?
i : nodes_tags[
i]].insert(verts[
i]);
391 for (
EntityType ent_type = MBVERTEX; ent_type < MBMAXTYPE; ent_type++) {
395 for (
auto type : types) {
398 med_bool change_of_coord;
399 med_bool geo_transform;
400 med_int num_ele = MEDmeshnEntity(
401 fid, mesh_name, MED_NO_DT, MED_NO_IT, MED_CELL, type,
402 MED_CONNECTIVITY, MED_NODAL, &change_of_coord, &geo_transform);
408 int num_nod_per_ele = type % 100;
411 <<
"Reading elements " << num_ele <<
" of type "
412 << moab::CN::EntityTypeName(ent_type) <<
" number of nodes "
415 std::vector<med_int> conn_med(num_ele * num_nod_per_ele);
416 if (MEDmeshElementConnectivityRd(fid, mesh_name, MED_NO_DT, MED_NO_IT,
417 MED_CELL, type, MED_NODAL,
418 MED_FULL_INTERLACE, &conn_med[0]) < 0) {
420 "Could not read MED elements");
428 if (ent_type != MBVERTEX) {
431 CHKERR iface->get_element_connect(num_ele, num_nod_per_ele, ent_type, 0,
437 for (
int ee = 0; ee != num_ele; ee++) {
439 for (
int nn = 0; nn != num_nod_per_ele; nn++) {
441 n[nn] = verts[conn_med[ii + nn] - 1];
444 std::array<int, 10> nodes_map{
453 for (
int nn = 0; nn != num_nod_per_ele; nn++, ii++) {
454 conn_moab[ii] =
n[nodes_map[nn]];
460 for (
int ee = 0; ee != num_ele; ee++) {
461 for (
int nn = 0; nn != num_nod_per_ele; nn++, ii++) {
463 conn_moab[ii] = verts[conn_med[ii] - 1];
469 CHKERR iface->update_adjacencies(starte, num_ele, num_nod_per_ele,
471 ents =
Range(starte, starte + num_ele - 1);
475 std::vector<EntityHandle> conn_moab(num_ele * num_nod_per_ele);
476 for (
int ee = 0; ee != num_ele; ++ee)
477 for (
int nn = 0; nn != num_nod_per_ele; ++nn, ++ii)
478 conn_moab[ii] = verts[conn_med[ii] - 1];
479 ents.insert_list(conn_moab.begin(), conn_moab.end());
487 std::vector<med_int> fam(num_ele, 0);
488 if (MEDmeshEntityFamilyNumberRd(fid, mesh_name, MED_NO_DT, MED_NO_IT,
489 MED_CELL, type, &fam[0]) < 0) {
491 "No family number for elements: using 0 as default family "
500 for (
int j = 0;
j < num_ele;
j++) {
502 family_elem_map[fam[
j]].insert(ents[
j]);
508 if (MEDfileClose(fid) < 0) {
510 "Unable to close file '%s'", file.c_str());
518 const std::map<int, Range> &family_elem_map,
519 std::map<string, Range> &group_elem_map,
int verb) {
524 med_idt fid = MEDfileOpen(file.c_str(), MED_ACC_RDONLY);
527 "Unable to open file '%s'", file.c_str());
530 MEDlibraryNumVersion(&
v[0], &
v[1], &
v[2]);
531 MEDfileNumVersionRd(fid, &vf[0], &vf[1], &vf[2]);
535 "Reading MED file V%d.%d.%d using MED library V%d.%d.%d", vf[0],
536 vf[1], vf[2],
v[0],
v[1],
v[2]);
538 if (vf[0] < 2 || (vf[0] == 2 && vf[1] < 2)) {
540 "Cannot read MED file older than V2.2");
544 group_elem_map.clear();
546 med_int num_families = MEDnFamily(fid,
meshNames[index].c_str());
547 if (num_families < 0) {
549 "Could not read MED families");
551 for (
int i = 0;
i < num_families;
i++) {
554 ? MEDnFamily23Attribute(fid,
meshNames[index].c_str(),
i + 1)
556 med_int num_groups = MEDnFamilyGroup(fid,
meshNames[index].c_str(),
i + 1);
557 if (num_attrib < 0 || num_groups < 0) {
559 "Could not read MED groups or attributes");
562 std::vector<med_int> attribId(num_attrib + 1);
563 std::vector<med_int> attrib_val(num_attrib + 1);
564 std::vector<char> attrib_des(MED_COMMENT_SIZE * num_attrib + 1);
565 std::vector<char> group_names(MED_LNAME_SIZE * num_groups + 1);
566 char family_name[MED_NAME_SIZE + 1];
570 if (MEDfamily23Info(fid,
meshNames[index].c_str(),
i + 1, family_name,
571 &attribId[0], &attrib_val[0], &attrib_des[0],
572 &family_num, &group_names[0]) < 0) {
574 "Could not read info for MED2 family %d",
i + 1);
577 if (MEDfamilyInfo(fid,
meshNames[index].c_str(),
i + 1, family_name,
578 &family_num, &group_names[0]) < 0) {
580 "Could not read info for MED3 family %d",
i + 1);
585 for (
int g = 0;
g != num_groups;
g++) {
587 std::string(&group_names[MED_LNAME_SIZE *
g], MED_LNAME_SIZE - 1);
588 name.resize(NAME_TAG_SIZE - 1);
589 if (family_elem_map.find(family_num) == family_elem_map.end()) {
591 "MEDWORLD", Sev::warning,
592 "Family %d not read, likely type of element is not added "
593 "to moab database. Currently only triangle, quad, tetrahedral and "
594 "hexahedral elements are read to moab database",
597 group_elem_map[name].merge(family_elem_map.at(family_num));
603 if (MEDfileClose(fid) < 0) {
605 "Unable to close file '%s'", file.c_str());
622 max_id = (max_id < cit->getMeshsetId()) ? cit->getMeshsetId() : max_id;
626 for (std::map<string, Range>::const_iterator git = group_elem_map.begin();
627 git != group_elem_map.end(); git++) {
629 CubitMeshSet_multiIndex::index<
636 if (!git->second.empty()) {
643 MOFEM_LOG(
"MEDWORLD", Sev::verbose) << *cit;
666 const bool load_series,
667 const int only_step,
int verb) {
671 med_idt fid = MEDfileOpen((
char *)file_name.c_str(), MED_LECTURE);
674 "Unable to open file '%s'", file_name.c_str());
677 med_int num_comp = MEDfieldnComponentByName(fid,
field_name.c_str());
680 "Could not get number of components for MED field");
683 char meshName[MED_NAME_SIZE + 1];
684 char dtUnit[MED_SNAME_SIZE + 1];
685 std::vector<char> compName(num_comp * MED_SNAME_SIZE + 1);
686 std::vector<char> compUnit(num_comp * MED_SNAME_SIZE + 1);
687 med_int numSteps = 0;
690 if (MEDfieldInfoByName(fid,
field_name.c_str(), meshName, &localMesh, &type,
691 &compName[0], &compUnit[0], dtUnit, &numSteps) < 0) {
693 "Could not get MED field info");
706 int num_comp_msh = (num_comp <= 1) ? 1
707 : (num_comp <= 3) ? 3
708 : (num_comp <= 9) ? 9
716 std::vector<double> def_val(num_comp_msh, 0);
718 tag_name.c_str(), num_comp_msh, MB_TYPE_DOUBLE,
th,
719 MB_TAG_CREAT | MB_TAG_SPARSE, &def_val[0]);
726 const med_entity_type entType[] = {MED_NODE, MED_CELL, MED_NODE_ELEMENT};
727 const med_geometrie_element eleType[] = {
728 MED_NONE, MED_SEG2, MED_TRIA3, MED_QUAD4, MED_TETRA4, MED_HEXA8,
729 MED_PENTA6, MED_PYRA5, MED_SEG3, MED_TRIA6, MED_QUAD9, MED_TETRA10,
730 MED_HEXA27, MED_POINT1, MED_QUAD8, MED_HEXA20, MED_PENTA15, MED_PYRA13};
734 std::vector<std::pair<int, int>> pairs;
735 for (
unsigned int i = 0;
i <
sizeof(entType) /
sizeof(entType[0]);
i++) {
736 for (
unsigned int j = 0;
j <
sizeof(eleType) /
sizeof(eleType[0]);
j++) {
737 if ((!
i && !
j) ||
j) {
738 med_int
n = numSteps;
740 pairs.push_back(std::pair<int, int>(
i,
j));
741 numSteps = std::max(numSteps,
n);
749 if (numSteps < 1 || pairs.empty()) {
751 "Nothing to import from MED file");
758 for (
int step = (only_step == -1) ? 0 : only_step; step < numSteps; step++) {
760 if (!load_series && only_step != step)
770 for (
unsigned int pair = 0; pair < pairs.size(); pair++) {
773 med_entite_maillage ent = entType[pairs[pair].first];
774 med_geometrie_element ele = eleType[pairs[pair].second];
775 med_int numdt, numit, ngauss;
777 if (MEDfieldComputingStepInfo(fid,
field_name.c_str(), step + 1, &numdt,
780 "Could not read step info");
783 char locName[MED_NAME_SIZE + 1], profileName[MED_NAME_SIZE + 1];
789 med_int numVal = MEDfieldnValueWithProfile(
790 fid,
field_name.c_str(), numdt, numit, ent, ele, 1,
791 MED_COMPACT_STMODE, profileName, &profileSize, locName, &ngauss);
806 std::vector<double> val(numVal * num_comp);
807 if (MEDfieldValueWithProfileRd(fid,
field_name.c_str(), numdt, numit, ent,
808 ele, MED_COMPACT_STMODE, profileName,
809 MED_FULL_INTERLACE, MED_ALL_CONSTITUENT,
810 (
unsigned char *)&val[0]) < 0) {
812 "Could not read field values");
836 "Not yet implemented for this cell %d", ele);
838 if (ent_type != MBMAXTYPE) {
843 double e_vals[num_comp_msh];
844 bzero(e_vals,
sizeof(
double) * num_comp_msh);
845 std::vector<double>::iterator vit = val.begin();
846 for (Range::iterator eit = ents.begin(); eit != ents.end(); eit++) {
847 for (
int ii = 0; ii != num_comp; ii++, vit++) {
856 if (ents.size() * ngauss * num_comp != val.size()) {
858 "data inconsistency");
861 double e_vals[num_comp_msh];
862 std::vector<double>::iterator vit = val.begin();
863 for (Range::iterator eit = ents.begin(); eit != ents.end(); eit++) {
864 bzero(e_vals,
sizeof(
double) * num_comp_msh);
865 for (
int gg = 0; gg != ngauss; gg++) {
866 for (
int ii = 0; ii != num_comp; ii++, vit++) {
867 e_vals[ii] += *vit / ngauss;
882 case MED_NODE_ELEMENT: {
885 CHKERR m_field.
get_moab().get_entities_by_type(meshset, ent_type, ents,
887 double e_vals[num_comp_msh];
888 bzero(e_vals,
sizeof(
double) * num_comp_msh);
889 std::vector<double>::iterator vit = val.begin();
890 for (Range::iterator eit = ents.begin(); eit != ents.end(); eit++) {
891 for (
int ii = 0; ii != num_comp; ii++, vit++) {
898 MOFEM_LOG_C(
"MEDWORLD", Sev::inform,
"Entity type %d not implemented",
909 os <<
"field name: " << field_data.
fieldName;
910 os <<
" mesh name: " << field_data.
meshName;
911 os <<
" local mesh: " << ((field_data.
localMesh) ?
"true" :
"false");
920 os <<
" componentNames:";
921 for (
unsigned int ff = 0; ff != field_data.
componentNames.size(); ff++) {
925 os <<
" componentUnits:";
926 for (
unsigned int ff = 0; ff != field_data.
componentUnits.size(); ff++) {
930 os <<
" dtUnit: " << field_data.
dtUnit << endl;
931 os <<
" number of steps: " << field_data.
ncSteps;
#define MOFEM_LOG_C(channel, severity, format,...)
Med file interface interface.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
@ MOFEM_OPERATION_UNSUCCESSFUL
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
FTensor::Index< 'n', SPACE_DIM > n
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
#define MOFEM_LOG(channel, severity)
Log.
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
MoFEMErrorCode addMeshset(const CubitBCType cubit_bc_type, const int ms_id, const std::string name="")
add cubit meshset
MoFEMErrorCode getCubitMeshsetPtr(const int ms_id, const CubitBCType cubit_bc_type, const CubitMeshSets **cubit_meshset_ptr) const
get cubit meshset
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
FTensor::Index< 'i', SPACE_DIM > i
const double v
phase velocity of light in medium (cm/ns)
FTensor::Index< 'j', 3 > j
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
std::bitset< 32 > CubitBCType
implementation of Data Operators for Forces and Sources
std::ostream & operator<<(std::ostream &os, const EntitiesFieldData::EntData &e)
static std::vector< med_geometrie_element > moab2med_element_type(const EntityType type)
constexpr auto field_name
virtual moab::Interface & get_moab()=0
virtual MPI_Comm & get_comm() const =0
this struct keeps basic methods for moab meshset about material and boundary conditions
EntityHandle getMeshset() const
get bc meshset
Deprecated interface functions.
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
static bool checkIfChannelExist(const std::string channel)
Check if channel exist.
std::vector< std::string > componentNames
std::vector< std::string > componentUnits
Interface for load MED files.
MoFEMErrorCode readMed(const string &file, int verb=1)
read MED file
std::vector< EntityHandle > meshMeshsets
meshset for each mesh
MoFEMErrorCode readFamily(const string &file, const int index, const std::map< int, Range > &family_elem_mapint, std::map< string, Range > &group_elem_map, int verb=1)
read family and groups
MoFEMErrorCode makeBlockSets(const std::map< string, Range > &group_elem_map, int verb=1)
make from groups meshsets
std::map< std::string, FieldData > fieldNames
std::string medFileName
MED file name.
MoFEMErrorCode getFileNameFromCommandLine(int verb=1)
Get MED file name from command line.
MedInterface(const MoFEM::Core &core)
std::vector< std::string > meshNames
list of meshes in MED file
MoFEM::Core & cOre
core database
MoFEMErrorCode writeMed(const string &file, int verb=1)
write MED file
MoFEMErrorCode readFields(const std::string &file_name, const std::string &field_name, const bool load_series=false, const int only_step=-1, int verb=1)
MoFEMErrorCode readMesh(const string &file, const int index, std::map< int, Range > &family_elem_map, int verb=1)
read mesh from MED file
PetscBool flgFile
true if file name given in command line
MoFEMErrorCode medGetFieldNames(const string &file, int verb=1)
Get field names in MED file.
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Interface for managing meshsets containing materials and boundary conditions.
CubitMeshSet_multiIndex & getMeshsetsMultindex()
base class for all interface classes
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.