v0.8.4
Classes | Public Member Functions | Public Attributes | Private Member Functions | Private Attributes | List of all members
MoFEM::MedInterface Struct Reference

Interface for load MED files. More...

#include <src/interfaces/MedInterface.hpp>

Inheritance diagram for MoFEM::MedInterface:
[legend]
Collaboration diagram for MoFEM::MedInterface:
[legend]

Classes

struct  FieldData
 

Public Member Functions

MoFEMErrorCode query_interface (const MOFEMuuid &uuid, UnknownInterface **iface) const
 
 MedInterface (const MoFEM::Core &core)
 
MoFEMErrorCode getFileNameFromCommandLine (int verb=1)
 Get MED file name from command line. More...
 
PetscBool getFlgFile () const
 Check if file name is given in command line. More...
 
MoFEMErrorCode medGetFieldNames (const string &file, int verb=1)
 Get field names in MED file. More...
 
MoFEMErrorCode medGetFieldNames (int verb=1)
 get field names in MED file More...
 
MoFEMErrorCode readMed (const string &file, int verb=1)
 read MED file More...
 
MoFEMErrorCode readMed (int verb=1)
 read MED file More...
 
MoFEMErrorCode writeMed (const string &file, int verb=1)
 write MED file More...
 
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)
 
- Public Member Functions inherited from MoFEM::UnknownInterface
template<class IFACE >
MoFEMErrorCode registerInterface (const MOFEMuuid &uuid, bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE , bool VERIFY = false>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, IFACE *&iface) const
 Get interface by uuid and return reference to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
virtual ~UnknownInterface ()
 
virtual MoFEMErrorCode getLibVersion (Version &version) const
 Get library version. More...
 
virtual const MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version) const
 Get database major version. More...
 
virtual MoFEMErrorCode getInterfaceVersion (Version &version) const
 Get database major version. More...
 
template<>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, UnknownInterface *&iface) const
 

Public Attributes

std::vector< std::string > meshNames
 list of meshes in MED file More...
 
std::vector< EntityHandlemeshMeshsets
 meshset for each mesh More...
 
std::map< std::string, FieldDatafieldNames
 
std::string medFileName
 MED file name. More...
 

Private Member Functions

MoFEMErrorCode readMesh (const string &file, const int index, std::map< int, Range > &family_elem_map, int verb=1)
 read mesh from MED file More...
 
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 More...
 
MoFEMErrorCode makeBlockSets (const std::map< string, Range > &group_elem_map, int verb=1)
 make from groups meshsets More...
 

Private Attributes

MoFEM::CorecOre
 core database More...
 
PetscBool flgFile
 true if file name given in command line More...
 

Additional Inherited Members

- Protected Member Functions inherited from MoFEM::UnknownInterface
boost::typeindex::type_index getClassIdx (const MOFEMuuid &uid) const
 Get type name for interface Id. More...
 
MOFEMuuid getUId (const boost::typeindex::type_index &class_idx) const
 Get interface Id for class name. More...
 

Detailed Description

Interface for load MED files.

Definition at line 35 of file MedInterface.hpp.

Constructor & Destructor Documentation

◆ MedInterface()

MoFEM::MedInterface::MedInterface ( const MoFEM::Core core)

Definition at line 61 of file MedInterface.cpp.

61  :
62  cOre(const_cast<Core&>(core)),
63  flgFile(PETSC_FALSE) {}
MoFEM::Core & cOre
core database
PetscBool flgFile
true if file name given in command line

Member Function Documentation

◆ getFileNameFromCommandLine()

MoFEMErrorCode MoFEM::MedInterface::getFileNameFromCommandLine ( int  verb = 1)

Get MED file name from command line.

Parameters
verbverbosity level
Returns
error code

Definition at line 65 of file MedInterface.cpp.

65  {
66  Interface &m_field = cOre;
67  char mesh_file_name[255];
69  ierr = PetscOptionsBegin(
70  m_field.get_comm(),"","MED Interface","none"
71  ); CHKERRG(ierr);
72  ierr = PetscOptionsString(
73  "-med_file",
74  "med file name","", "mesh.med",mesh_file_name, 255, &flgFile
75  ); CHKERRG(ierr);
76  ierr = PetscOptionsEnd(); CHKERRG(ierr);
77  if(flgFile == PETSC_TRUE) {
78  medFileName = std::string(mesh_file_name);
79  } else {
80  medFileName = std::string("mesh.med");
81  }
83  }
MoFEM::Core & cOre
core database
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:522
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:565
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:528
char mesh_file_name[255]
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
PetscBool flgFile
true if file name given in command line
std::string medFileName
MED file name.

◆ getFlgFile()

PetscBool MoFEM::MedInterface::getFlgFile ( ) const

Check if file name is given in command line.

Definition at line 50 of file MedInterface.hpp.

50 { return flgFile; }
PetscBool flgFile
true if file name given in command line

◆ makeBlockSets()

MoFEMErrorCode MoFEM::MedInterface::makeBlockSets ( const std::map< string, Range > &  group_elem_map,
int  verb = 1 
)
private

make from groups meshsets

Parameters
group_elem_mapmap of groups and elements
verbverbosity level
Returns
error code

Definition at line 567 of file MedInterface.cpp.

570  {
571 
572 
573  Interface &m_field = cOre;
575  MeshsetsManager *meshsets_manager_ptr;
576  ierr = m_field.getInterface(meshsets_manager_ptr); CHKERRG(ierr);
577 
578  int max_id = 0;
580  max_id = (max_id < cit->getMeshsetId()) ? cit->getMeshsetId() : max_id;
581  }
582  max_id++;
583 
584  // cerr << group_elem_map.size() << endl;
585  for(
586  std::map<string,Range>::const_iterator git = group_elem_map.begin();
587  git!=group_elem_map.end();git++
588  ) {
589  // cerr << "AAA\n";
590  ierr = meshsets_manager_ptr->addMeshset(BLOCKSET,max_id,git->first); CHKERRG(ierr);
591  CubitMeshSet_multiIndex::index<Composite_Cubit_msId_And_MeshSetType_mi_tag>::type::iterator cit;
592  cit = meshsets_manager_ptr->getMeshsetsMultindex().get<Composite_Cubit_msId_And_MeshSetType_mi_tag>().find(
593  boost::make_tuple(max_id,CubitBCType(BLOCKSET).to_ulong())
594  );
595  EntityHandle meshsets = cit->getMeshset();
596  if(!git->second.empty()) {
597  rval = m_field.get_moab().add_entities(meshsets,git->second); CHKERRQ_MOAB(rval);
598  }
599  max_id++;
600  // cerr << git->second << endl;
601  }
602 
603  // if(verb>0) {
605  PetscPrintf(
606  m_field.get_comm(),"%s\n",
607  static_cast<std::ostringstream&>(
608  std::ostringstream().seekp(0) << *cit
609  ).str().c_str()
610  );
611  }
612  // }
613 
615  }
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:536
MoFEM::Core & cOre
core database
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:522
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:565
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:528
std::bitset< 32 > CubitBCType
Definition: Common.hpp:200
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Common.hpp:78
#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...

◆ medGetFieldNames() [1/2]

MoFEMErrorCode MoFEM::MedInterface::medGetFieldNames ( const string &  file,
int  verb = 1 
)

Get field names in MED file.

Parameters
filefile name
verbverbosity level
Returns
error code

Definition at line 85 of file MedInterface.cpp.

85  {
86  Interface &m_field = cOre;
88  med_idt fid = MEDfileOpen(file.c_str(), MED_ACC_RDONLY);
89  if (fid < 0) {
90  SETERRQ1(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
91  "Unable to open file '%s'", file.c_str());
92  }
93  med_int num_fields = MEDnField(fid);
94  for (int index = 0; index < num_fields; index++) {
95 
96  med_int num_comp = MEDfieldnComponent(fid, index + 1);
97  if (num_comp <= 0) {
98  SETERRQ(m_field.get_comm(), MOFEM_IMPOSIBLE_CASE,
99  "Could not get number of components for MED field");
100  }
101 
102  char name[MED_NAME_SIZE + 1], mesh_name[MED_NAME_SIZE + 1];
103  char dt_unit[MED_SNAME_SIZE + 1];
104  std::vector<char> comp_name(num_comp * MED_SNAME_SIZE + 1);
105  std::vector<char> comp_unit(num_comp * MED_SNAME_SIZE + 1);
106  med_int num_steps = 0;
107  med_type_champ type;
108  med_bool local_mesh;
109  if (MEDfieldInfo(fid, index + 1, name, mesh_name, &local_mesh, &type,
110  &comp_name[0], &comp_unit[0], dt_unit, &num_steps) < 0) {
111  SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
112  "Could not get MED field info");
113  }
114 
115  std::string field_name = std::string(name);
116  fieldNames[field_name] = FieldData();
117  fieldNames[field_name].fieldName = field_name;
118  fieldNames[field_name].meshName = std::string(mesh_name);
119  fieldNames[field_name].localMesh = (local_mesh == MED_TRUE);
120  for (int ff = 0; ff != num_comp; ff++) {
121  fieldNames[field_name].componentNames.push_back(
122  std::string(&comp_name[ff * MED_SNAME_SIZE], MED_SNAME_SIZE));
123  fieldNames[field_name].componentUnits.push_back(
124  std::string(&comp_unit[ff * MED_SNAME_SIZE], MED_SNAME_SIZE));
125  }
126  fieldNames[field_name].dtUnit = std::string(&dt_unit[0]);
127  fieldNames[field_name].ncSteps = num_steps;
128 
129  if (verb > 0) {
130  std::ostringstream ss;
131  ss << fieldNames[name] << std::endl;
132  ierr = PetscPrintf(m_field.get_comm(), ss.str().c_str());
133  CHKERRG(ierr);
134  }
135  }
136  if (MEDfileClose(fid) < 0) {
137  SETERRQ1(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
138  "Unable to close file '%s'", file.c_str());
139  }
141  }
std::map< std::string, FieldData > fieldNames
MoFEM::Core & cOre
core database
double FieldData
Field data type.
Definition: Common.hpp:130
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:522
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:565
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:528
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80

◆ medGetFieldNames() [2/2]

MoFEMErrorCode MoFEM::MedInterface::medGetFieldNames ( int  verb = 1)

get field names in MED file

File name is get form command line

Parameters
verbverbosity level
Returns
error code

Definition at line 143 of file MedInterface.cpp.

143  {
145  if(medFileName.empty()) {
147  }
150  }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:565
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
MoFEMErrorCode getFileNameFromCommandLine(int verb=1)
Get MED file name from command line.
#define CHKERR
Inline error check.
Definition: definitions.h:617
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443
MoFEMErrorCode medGetFieldNames(const string &file, int verb=1)
Get field names in MED file.
std::string medFileName
MED file name.

◆ query_interface()

MoFEMErrorCode MoFEM::MedInterface::query_interface ( const MOFEMuuid uuid,
UnknownInterface **  iface 
) const
virtual

Implements MoFEM::UnknownInterface.

Definition at line 50 of file MedInterface.cpp.

50  {
52  *iface = NULL;
53  if(uuid == IDD_MOFEMMedInterface) {
54  *iface = const_cast<MedInterface*>(this);
56  }
57  SETERRQ(PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,"unknown interface");
59  }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:522
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:528
MedInterface(const MoFEM::Core &core)
static const MOFEMuuid IDD_MOFEMMedInterface

◆ readFamily()

MoFEMErrorCode MoFEM::MedInterface::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 
)
private

read family and groups

Parameters
filefile name
indexmesh index
family_elem_mapintmap of families and elements
group_elem_mapmap of groups and elements
verbverbosity level
Returns
error code

Definition at line 478 of file MedInterface.cpp.

484  {
485  //
486  Interface &m_field = cOre;
488 
489  med_idt fid = MEDfileOpen(file.c_str(),MED_ACC_RDONLY);
490  if(fid < 0) {
491  SETERRQ1(m_field.get_comm(),MOFEM_OPERATION_UNSUCCESSFUL,"Unable to open file '%s'",file.c_str());
492  }
493  med_int v[3], vf[3];
494  MEDlibraryNumVersion(&v[0], &v[1], &v[2]);
495  MEDfileNumVersionRd(fid, &vf[0], &vf[1], &vf[2]);
496  // if(verb>1) {
497  // PetscPrintf(
498  // m_field.get_comm(),
499  // "Reading MED file V%d.%d.%d using MED library V%d.%d.%d\n",
500  // vf[0], vf[1], vf[2], v[0], v[1], v[2]
501  // );
502  // }
503  if(vf[0] < 2 || (vf[0] == 2 && vf[1] < 2)){
504  SETERRQ(m_field.get_comm(),MOFEM_DATA_INCONSISTENCY,"Cannot read MED file older than V2.2");
505  }
506 
507  // clear groups
508  group_elem_map.clear();
509 
510  med_int num_families = MEDnFamily(fid,meshNames[index].c_str());
511  if(num_families < 0){
512  SETERRQ(m_field.get_comm(),MOFEM_OPERATION_UNSUCCESSFUL,"Could not read MED families");
513  }
514  for(int i = 0;i<num_families;i++) {
515  med_int num_attrib = (vf[0] == 2) ? MEDnFamily23Attribute(fid,meshNames[index].c_str(),i+1) : 0;
516  med_int num_groups = MEDnFamilyGroup(fid,meshNames[index].c_str(),i+1);
517  if(num_attrib < 0 || num_groups < 0){
518  SETERRQ(
519  m_field.get_comm(),
521  "Could not read MED groups or attributes"
522  );
523  }
524 
525  std::vector<med_int> attribId(num_attrib + 1);
526  std::vector<med_int> attrib_val(num_attrib + 1);
527  std::vector<char> attrib_des(MED_COMMENT_SIZE * num_attrib + 1);
528  std::vector<char> group_names(MED_LNAME_SIZE * num_groups + 1);
529  char family_name[MED_NAME_SIZE + 1];
530  med_int family_num;
531 
532  if(vf[0] == 2) { // MED2 file
533  if(MEDfamily23Info(
534  fid,meshNames[index].c_str(), i + 1, family_name, &attribId[0],
535  &attrib_val[0], &attrib_des[0], &family_num,
536  &group_names[0]) < 0
537  ) {
538  SETERRQ1(m_field.get_comm(),MOFEM_OPERATION_UNSUCCESSFUL,"Could not read info for MED2 family %d",i+1);
539  }
540  }
541  else{
542  if(MEDfamilyInfo(
543  fid,meshNames[index].c_str(),i+1,family_name,&family_num,
544  &group_names[0]
545  ) < 0) {
546  SETERRQ1(m_field.get_comm(),MOFEM_OPERATION_UNSUCCESSFUL,"Could not read info for MED3 family %d",i+1);
547  }
548  }
549 
550  // cerr << family_name << " " << family_num << " " << num_groups << endl;
551  for(int g = 0;g!=num_groups;g++) {
552  std::string name = std::string(&group_names[MED_LNAME_SIZE*g],MED_LNAME_SIZE-1);
553  name.resize(NAME_TAG_SIZE-1);
554  group_elem_map[name].merge(family_elem_map.at(family_num));
555  // cerr << string(&group_names[MED_LNAME_SIZE*g]) << endl;
556  }
557 
558  }
559 
560  if(MEDfileClose(fid) < 0) {
561  SETERRQ1(m_field.get_comm(),MOFEM_OPERATION_UNSUCCESSFUL,"Unable to close file '%s'",file.c_str());
562  }
563 
565  }
MoFEM::Core & cOre
core database
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:522
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:528
std::vector< std::string > meshNames
list of meshes in MED file

◆ readFields()

MoFEMErrorCode MoFEM::MedInterface::readFields ( const std::string &  file_name,
const std::string &  field_name,
const bool  load_series = false,
const int  only_step = -1,
int  verb = 1 
)

Read fields

Parameters
file_namefile name
file_indexmesh index
loadSeriesload field into series
onlyStepread only one step
Returns
error code

Definition at line 633 of file MedInterface.cpp.

639  {
640 
641 
642  Interface &m_field = cOre;
644  med_idt fid = MEDfileOpen((char*)file_name.c_str(), MED_LECTURE);
645  if(fid < 0){
646  SETERRQ1(m_field.get_comm(),MOFEM_DATA_INCONSISTENCY,"Unable to open file '%s'", file_name.c_str());
647  }
648 
649  med_int num_comp = MEDfieldnComponentByName(fid,field_name.c_str());
650  if(num_comp <= 0){
651  SETERRQ(m_field.get_comm(),MOFEM_DATA_INCONSISTENCY,"Could not get number of components for MED field");
652  }
653 
654  char meshName[MED_NAME_SIZE + 1];
655  char dtUnit[MED_SNAME_SIZE + 1];
656  std::vector<char> compName(num_comp * MED_SNAME_SIZE + 1);
657  std::vector<char> compUnit(num_comp * MED_SNAME_SIZE + 1);
658  med_int numSteps = 0;
659  med_type_champ type;
660  med_bool localMesh;
661  if(MEDfieldInfoByName(
662  fid, field_name.c_str(), meshName, &localMesh, &type,
663  &compName[0], &compUnit[0], dtUnit, &numSteps
664  ) < 0) {
665  SETERRQ(m_field.get_comm(),MOFEM_DATA_INCONSISTENCY,"Could not get MED field info");
666  }
667 
668  // Get meshset
669  MeshsetsManager *meshsets_manager_ptr;
670  ierr = m_field.getInterface(meshsets_manager_ptr); CHKERRG(ierr);
671  const CubitMeshSets *cubit_meshset_ptr;
672  ierr = meshsets_manager_ptr->getCubitMeshsetPtr(meshName,&cubit_meshset_ptr); CHKERRG(ierr);
673  EntityHandle meshset = cubit_meshset_ptr->getMeshset();
674 
675  int num_comp_msh =
676  (num_comp <= 1) ? 1 : (num_comp <= 3) ? 3 : (num_comp <= 9) ? 9 : num_comp;
677 
678  // Create tag to store nodal or cell values read form med file. Note that tag
679  // has prefix MED to avoid collison with other tags.
680  Tag th;
681  std::string tag_name = "MED_"+field_name;
682  {
683  std::vector<double> def_val(num_comp_msh,0);
684  rval = m_field.get_moab().tag_get_handle(
685  tag_name.c_str(),num_comp_msh,MB_TYPE_DOUBLE,th,MB_TAG_CREAT|MB_TAG_SPARSE,&def_val[0]
686  ); CHKERRQ_MOAB(rval);
687  }
688 
689  // Warning! The ordering of the elements in the last two lists is
690  // important: it should match the ordering of the MSH element types
691  // (when elements are saved without tags, this governs the order
692  // with which we implicitly index them in GModel::readMED)
693  const med_entity_type entType[] = { MED_NODE, MED_CELL, MED_NODE_ELEMENT };
694  const med_geometrie_element eleType[] = {
695  MED_NONE, MED_SEG2, MED_TRIA3, MED_QUAD4, MED_TETRA4, MED_HEXA8,
696  MED_PENTA6, MED_PYRA5, MED_SEG3, MED_TRIA6, MED_QUAD9, MED_TETRA10,
697  MED_HEXA27, MED_POINT1, MED_QUAD8, MED_HEXA20, MED_PENTA15, MED_PYRA13
698  };
699  //const int nodesPerEle[] = {0, 2, 3, 4, 4, 8, 6, 5, 3, 6, 9, 10, 27, 1, 8, 20, 15, 13};
700 
701  std::vector<std::pair<int, int> > pairs;
702  for(unsigned int i = 0; i < sizeof(entType) / sizeof(entType[0]); i++){
703  for(unsigned int j = 0; j < sizeof(eleType) / sizeof(eleType[0]); j++){
704  if((!i && !j) || j){
705  med_int n = numSteps;
706  if(n > 0){
707  pairs.push_back(std::pair<int, int>(i, j));
708  numSteps = std::max(numSteps, n);
709  }
710  if(!i && !j) break; // MED_NOEUD does not care about eleType
711  }
712  }
713  }
714 
715  if(numSteps < 1 || pairs.empty()){
716  SETERRQ(m_field.get_comm(),MOFEM_DATA_INCONSISTENCY,"Nothing to import from MED file");
717  }
718 
719  if(load_series) {
720  SETERRQ(m_field.get_comm(), MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
721  }
722 
723  for(
724  int step = (only_step == -1) ? 0 : only_step;
725  step < numSteps;
726  step++
727  ) {
728 
729 
730  if(!load_series && only_step != step) continue;
731 
732  // cerr << only_step << " " << step << endl;
733 
734 
735  // FIXME: in MED3 we might want to loop over all profiles instead
736  // of relying of the default one
737 
738  // FIXME: MED3 allows to store multi-step meshes; we should
739 
740  for(unsigned int pair = 0; pair < pairs.size(); pair++) {
741 
742  // get step info
743  med_entite_maillage ent = entType[pairs[pair].first];
744  med_geometrie_element ele = eleType[pairs[pair].second];
745  med_int numdt, numit, ngauss;
746  med_float dt;
747  if(MEDfieldComputingStepInfo(fid, field_name.c_str(), step + 1, &numdt, &numit, &dt) < 0){
748  SETERRQ(m_field.get_comm(),MOFEM_OPERATION_UNSUCCESSFUL,"Could not read step info");
749  }
750 
751  char locName[MED_NAME_SIZE + 1], profileName[MED_NAME_SIZE + 1];
752 
753  // get number of values in the field (numVal takes the number of
754  // Gauss points or the number of nodes per element into account,
755  // but not the number of components)
756  med_int profileSize;
757  med_int numVal = MEDfieldnValueWithProfile(
758  fid, field_name.c_str(), numdt, numit, ent, ele,
759  1, MED_COMPACT_STMODE, profileName,
760  &profileSize, locName, &ngauss
761  );
762  numVal *= ngauss;
763 
764  if(numVal <= 0) continue;
765 
766  //int mult = 1;
767  //if(ent == MED_NODE_ELEMENT) {
768  // mult = nodesPerEle[pairs[pair].second];
769  //}
770  //else if(ngauss != 1){
771  // mult = ngauss;
772  //}
773 
774  // read field data
775  std::vector<double> val(numVal * num_comp);
776  if(
777  MEDfieldValueWithProfileRd(
778  fid, field_name.c_str(), numdt, numit, ent, ele, MED_COMPACT_STMODE,
779  profileName, MED_FULL_INTERLACE, MED_ALL_CONSTITUENT,
780  (unsigned char*)&val[0]
781  ) < 0
782  ) {
783  SETERRQ(m_field.get_comm(),MOFEM_OPERATION_UNSUCCESSFUL,"Could not read field values");
784  }
785 
786  if(verb>2) {
787  // FIXME: This not looks ok for me
788  cerr << ent << " " << ele << endl;
789  cerr << string(meshName)
790  << " : " << string(profileName)
791  << " : " << string(locName) << " : "
792  << profileSize << " : " << ngauss << endl;
793  }
794 
795  switch (ent) {
796  case MED_CELL: {
797  EntityType ent_type;
798  switch (ele) {
799  case MED_TETRA4:
800  case MED_TETRA10:
801  ent_type = MBTET;
802  break;
803  case MED_HEXA8:
804  ent_type = MBHEX;
805  break;
806  default:
807  SETERRQ1(
808  m_field.get_comm(),
810  "Not yet implemented for this cell %d",ele
811  );
812  }
813  if(ngauss==1) {
814  Range ents;
815  rval = m_field.get_moab().get_entities_by_type(
816  meshset,ent_type,ents,true
817  ); CHKERRQ_MOAB(rval);
818  double e_vals[num_comp_msh];
819  bzero(e_vals,sizeof(double)*num_comp_msh);
820  std::vector<double>::iterator vit = val.begin();
821  for(Range::iterator eit = ents.begin();eit!=ents.end();eit++) {
822  for(int ii = 0;ii!=num_comp;ii++,vit++) {
823  e_vals[ii] = *vit;
824  }
825  rval = m_field.get_moab().tag_set_data(th,&*eit,1,e_vals); CHKERRQ_MOAB(rval);
826  }
827  } else {
828  Range ents;
829  rval = m_field.get_moab().get_entities_by_type(
830  meshset,ent_type,ents,true
831  ); CHKERRQ_MOAB(rval);
832  if(ents.size()*ngauss*num_comp!=val.size()) {
833  SETERRQ(m_field.get_comm(),MOFEM_DATA_INCONSISTENCY,"data inconsistency");
834  }
835  // FIXME simply average gauss values, far from perfect need fix
836  double e_vals[num_comp_msh];
837  std::vector<double>::iterator vit = val.begin();
838  for(Range::iterator eit = ents.begin();eit!=ents.end();eit++) {
839  bzero(e_vals,sizeof(double)*num_comp_msh);
840  for(int gg = 0;gg!=ngauss;gg++) {
841  for(int ii = 0;ii!=num_comp;ii++,vit++) {
842  e_vals[ii] += *vit/ngauss;
843  }
844  }
845  rval = m_field.get_moab().tag_set_data(th,&*eit,1,e_vals); CHKERRQ_MOAB(rval);
846  }
847  }
848  // SETERRQ1(
849  // m_field.get_comm(),
850  // MOFEM_NOT_IMPLEMENTED,
851  // "Not implemented for more gauss pts ngauss = %d",
852  // ngauss
853  // );
854  }
855  break;
856  case MED_NODE:
857  case MED_NODE_ELEMENT:
858  default:
859  SETERRQ(m_field.get_comm(),MOFEM_NOT_IMPLEMENTED,"Not yet implemented");
860  }
861  }
862 
863  }
864 
866  }
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:536
MoFEM::Core & cOre
core database
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:522
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:565
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:528
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Common.hpp:78

◆ readMed() [1/2]

MoFEMErrorCode MoFEM::MedInterface::readMed ( const string &  file,
int  verb = 1 
)

read MED file

Parameters
filefiled name
verbverbosity level
Returns
error code

Definition at line 152 of file MedInterface.cpp.

152  {
153  Interface &m_field = cOre;
155 
156  med_idt fid = MEDfileOpen(file.c_str(), MED_ACC_RDONLY);
157  if (fid < 0) {
158  SETERRQ1(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
159  "Unable to open file '%s'", file.c_str());
160  }
161 
162  med_int v[3], vf[3];
163  MEDlibraryNumVersion(&v[0], &v[1], &v[2]);
164  MEDfileNumVersionRd(fid, &vf[0], &vf[1], &vf[2]);
165 
166  if (verb > 0) {
167  PetscPrintf(m_field.get_comm(),
168  "Reading MED file V%d.%d.%d using MED library V%d.%d.%d\n",
169  vf[0], vf[1], vf[2], v[0], v[1], v[2]);
170  }
171  if (vf[0] < 2 || (vf[0] == 2 && vf[1] < 2)) {
172  SETERRQ(m_field.get_comm(), MOFEM_DATA_INCONSISTENCY,
173  "Cannot read MED file older than V2.2");
174  }
175 
176  for (int i = 0; i < MEDnMesh(fid); i++) {
177  char mesh_name[MED_NAME_SIZE + 1], mesh_desc[MED_COMMENT_SIZE + 1];
178  bzero(mesh_name, MED_NAME_SIZE);
179  bzero(mesh_desc, MED_COMMENT_SIZE);
180  med_int space_dim;
181  med_mesh_type mesh_type;
182  med_int mesh_dim, n_step;
183  char dt_unit[MED_SNAME_SIZE + 1];
184  char axis_name[3 * MED_SNAME_SIZE + 1], axis_unit[3 * MED_SNAME_SIZE + 1];
185  med_sorting_type sorting_type;
186  med_axis_type axis_type;
187  if (MEDmeshInfo(fid, i + 1, mesh_name, &space_dim, &mesh_dim, &mesh_type,
188  mesh_desc, dt_unit, &sorting_type, &n_step, &axis_type,
189  axis_name, axis_unit) < 0) {
190  SETERRQ(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
191  "Unable to read mesh information");
192  }
193  meshNames.push_back(std::string(mesh_name));
194  if (verb > 0) {
195  PetscPrintf(m_field.get_comm(), "Check mesh %s nsteps %d\n", mesh_name,
196  n_step);
197  }
198  }
199 
200  std::map<int, Range> family_elem_map;
201  std::map<string, Range> group_elem_map;
202 
203  for (unsigned int ii = 0; ii != meshNames.size(); ii++) {
204  CHKERR readMesh(file, ii, family_elem_map, verb);
205  CHKERR readFamily(file, ii, family_elem_map, group_elem_map, verb);
206  CHKERR makeBlockSets(group_elem_map, verb);
207  }
208 
209  if (MEDfileClose(fid) < 0) {
210  SETERRQ1(m_field.get_comm(), MOFEM_OPERATION_UNSUCCESSFUL,
211  "Unable to close file '%s'", file.c_str());
212  }
213 
215  }
MoFEMErrorCode makeBlockSets(const std::map< string, Range > &group_elem_map, int verb=1)
make from groups meshsets
MoFEMErrorCode readMesh(const string &file, const int index, std::map< int, Range > &family_elem_map, int verb=1)
read mesh from MED file
MoFEM::Core & cOre
core database
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
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
std::vector< std::string > meshNames
list of meshes in MED file
#define CHKERR
Inline error check.
Definition: definitions.h:617
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443

◆ readMed() [2/2]

MoFEMErrorCode MoFEM::MedInterface::readMed ( int  verb = 1)

read MED file

File name is form command line

Parameters
verbverbosity level
Returns
error code

Definition at line 617 of file MedInterface.cpp.

617  {
619  if (medFileName.empty()) {
621  }
622  CHKERR readMed(medFileName, verb);
624  }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:498
MoFEMErrorCode getFileNameFromCommandLine(int verb=1)
Get MED file name from command line.
MoFEMErrorCode readMed(const string &file, int verb=1)
read MED file
#define CHKERR
Inline error check.
Definition: definitions.h:617
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:443
std::string medFileName
MED file name.

◆ readMesh()

MoFEMErrorCode MoFEM::MedInterface::readMesh ( const string &  file,
const int  index,
std::map< int, Range > &  family_elem_map,
int  verb = 1 
)
private

read mesh from MED file

Parameters
filefile name
indexindex of mesh
family_elem_mapmap of families and elements
verbverbosity level
Returns
error code

Definition at line 231 of file MedInterface.cpp.

236  {
237 
238  Interface &m_field = cOre;
240 
241  med_idt fid = MEDfileOpen(file.c_str(), MED_ACC_RDONLY);
242  if(fid < 0) {
243  SETERRQ1(m_field.get_comm(),MOFEM_OPERATION_UNSUCCESSFUL,"Unable to open file '%s'",file.c_str());
244  }
245  med_int v[3], vf[3];
246  MEDlibraryNumVersion(&v[0], &v[1], &v[2]);
247  MEDfileNumVersionRd(fid, &vf[0], &vf[1], &vf[2]);
248  // if(verb>1) {
249  // PetscPrintf(
250  // m_field.get_comm(),
251  // "Reading MED file V%d.%d.%d using MED library V%d.%d.%d\n",
252  // vf[0], vf[1], vf[2], v[0], v[1], v[2]
253  // );
254  // }
255  if(vf[0] < 2 || (vf[0] == 2 && vf[1] < 2)){
256  SETERRQ(m_field.get_comm(),MOFEM_DATA_INCONSISTENCY,"Cannot read MED file older than V2.2");
257  }
258 
259  char mesh_name[MED_NAME_SIZE + 1], mesh_desc[MED_COMMENT_SIZE + 1];
260  bzero(mesh_name,MED_NAME_SIZE+1);
261  bzero(mesh_desc,MED_COMMENT_SIZE+1);
262  med_int space_dim;
263  med_mesh_type mesh_type;
264  med_int mesh_dim, n_step;
265  char dt_unit[ MED_SNAME_SIZE + 1];
266  char axis_name[3 * MED_SNAME_SIZE + 1], axis_unit[3 * MED_SNAME_SIZE + 1];
267  med_sorting_type sorting_type;
268  med_axis_type axis_type;
269  if(MEDmeshInfo(
270  fid,index+1,mesh_name,&space_dim,&mesh_dim,&mesh_type,mesh_desc,
271  dt_unit,&sorting_type,&n_step,&axis_type,axis_name,axis_unit) < 0
272  ) {
273  SETERRQ(m_field.get_comm(),MOFEM_OPERATION_UNSUCCESSFUL,"Unable to read mesh information");
274  }
275  if(verb>0) {
276  PetscPrintf(m_field.get_comm(),"Reading mesh %s nsteps %d\n",mesh_name,n_step);
277  }
278 
279 
280  switch (axis_type) {
281  case MED_CARTESIAN:
282  break;
283  case MED_SPHERICAL:
284  case MED_CYLINDRICAL:
285  default:
286  SETERRQ(m_field.get_comm(),MOFEM_NOT_IMPLEMENTED,"Curvilinear coordinate system implemented");
287  }
288  if(space_dim!=3) {
289  SETERRQ1(m_field.get_comm(),MOFEM_NOT_IMPLEMENTED,"Not implemented for space dim %d",space_dim);
290  }
291 
292  EntityHandle mesh_meshset;
293  {
294  MeshsetsManager *meshsets_manager_ptr;
295  ierr = m_field.getInterface(meshsets_manager_ptr); CHKERRG(ierr);
296  int max_id = 0;
298  max_id = (max_id < cit->getMeshsetId()) ? cit->getMeshsetId() : max_id;
299  }
300  max_id++;
301  ierr = meshsets_manager_ptr->addMeshset(BLOCKSET,max_id,std::string(mesh_name)); CHKERRG(ierr);
302  CubitMeshSet_multiIndex::index<Composite_Cubit_msId_And_MeshSetType_mi_tag>::type::iterator cit;
303  cit = meshsets_manager_ptr->getMeshsetsMultindex().get<Composite_Cubit_msId_And_MeshSetType_mi_tag>().find(
304  boost::make_tuple(max_id,CubitBCType(BLOCKSET).to_ulong())
305  );
306  max_id++;
307  mesh_meshset = cit->getMeshset();
308  meshMeshsets.push_back(mesh_meshset);
309  }
310 
311  med_bool change_of_coord, geo_transform;
312  med_int num_nodes = MEDmeshnEntity(
313  fid, mesh_name, MED_NO_DT, MED_NO_IT, MED_NODE,
314  MED_NO_GEOTYPE, MED_COORDINATE, MED_NO_CMODE,
315  &change_of_coord, &geo_transform
316  );
317  if(num_nodes < 0){
318  SETERRQ(
319  m_field.get_comm(),
321  "Could not read number of MED nodes"
322  );
323  }
324  if(num_nodes == 0) {
325  SETERRQ(m_field.get_comm(),MOFEM_OPERATION_UNSUCCESSFUL,"No nodes in MED mesh");
326  }
327  if(verb>0) {
328  PetscPrintf(m_field.get_comm(),"Read number of nodes %d\n",num_nodes);
329  }
330 
331  std::vector<med_float> coord_med(space_dim * num_nodes);
332  if(MEDmeshNodeCoordinateRd(
333  fid, mesh_name, MED_NO_DT, MED_NO_IT, MED_NO_INTERLACE,&coord_med[0]
334  ) < 0 ) {
335  SETERRQ(m_field.get_comm(),MOFEM_OPERATION_UNSUCCESSFUL,"Could not read MED node coordinates");
336  }
337 
338  // Add vertices to moab
339  ReadUtilIface* iface;
340  vector<double*> arrays_coord;
341  EntityHandle startv;
342  rval = m_field.get_moab().query_interface(iface); CHKERRQ_MOAB(rval);
343  rval = iface->get_node_coords(3,num_nodes,0,startv,arrays_coord); CHKERRQ_MOAB(rval);
344  Range verts(startv,startv+num_nodes-1);
345  std::copy(&coord_med[0*num_nodes],&coord_med[1*num_nodes],arrays_coord[0]);
346  std::copy(&coord_med[1*num_nodes],&coord_med[2*num_nodes],arrays_coord[1]);
347  std::copy(&coord_med[2*num_nodes],&coord_med[3*num_nodes],arrays_coord[2]);
348  ierr = m_field.get_moab().add_entities(mesh_meshset,verts); CHKERRG(ierr);
349  family_elem_map.clear();
350 
351  // get family for vertices
352  {
353  std::vector<med_int> nodes_tags(num_nodes,0);
354  if(MEDmeshEntityFamilyNumberRd(
355  fid,mesh_name,MED_NO_DT,MED_NO_IT,MED_NODE,MED_NO_GEOTYPE,&nodes_tags[0]
356  ) < 0) {
357  nodes_tags.clear();
358  // SETERRQ(
359  // m_field.get_comm(),
360  // MOFEM_OPERATION_UNSUCCESSFUL,
361  // "No family number for elements: using 0 as default family number"
362  // );
363  }
364  for(int i = 0; i < num_nodes; i++) {
365  // cerr << verts[i] << " " /*<< ele_tags[j] << " "*/ << nodes_tags[i] << endl;
366  family_elem_map[nodes_tags.empty() ? i : nodes_tags[i]].insert(verts[i]);
367  }
368  }
369 
370  // read elements (loop over all possible MSH element types)
371  for(EntityType ent_type = MBVERTEX; ent_type < MBMAXTYPE; ent_type++){
372 
373  med_geometrie_element type = moab2med_element_type(ent_type);
374  if(type == MED_NONE) continue;
375 
376  // get number of cells
377  med_bool change_of_coord;
378  med_bool geo_transform;
379  med_int num_ele = MEDmeshnEntity(
380  fid, mesh_name, MED_NO_DT, MED_NO_IT, MED_CELL,
381  type, MED_CONNECTIVITY, MED_NODAL, &change_of_coord,
382  &geo_transform
383  );
384 
385  // get connectivity
386  if(num_ele <= 0) continue;
387  int num_nod_per_ele = type % 100;
388  std::vector<med_int> conn_med(num_ele*num_nod_per_ele);
389  if(MEDmeshElementConnectivityRd(
390  fid, mesh_name, MED_NO_DT, MED_NO_IT, MED_CELL,
391  type, MED_NODAL, MED_FULL_INTERLACE, &conn_med[0]
392  ) < 0){
393  SETERRQ(
394  m_field.get_comm(),MOFEM_OPERATION_UNSUCCESSFUL,"Could not read MED elements"
395  );
396  }
397 
398  // cerr << "type " << ent_type << " ";
399  // cerr << "num_ele " << num_ele << " " << num_nod_per_ele << endl;;
400 
401  EntityHandle* conn_moab;
402  EntityHandle starte;
403  rval = iface->get_element_connect(
404  num_ele,num_nod_per_ele,ent_type,0,starte,conn_moab
405  ); CHKERRQ_MOAB(rval);
406  switch (ent_type) {
407  // FIXME: Some connectivity could not work, need to run and test
408  case MBTET:
409  {
410  int ii = 0;
411  for(int ee = 0;ee!=num_ele;ee++) {
412  EntityHandle n[4];
413  for(int nn = 0;nn!=num_nod_per_ele;nn++) {
414  // conn_moab[ii] = verts[conn_med[ii]-1];
415  n[nn] = verts[conn_med[ii+nn]-1];
416  }
417  EntityHandle n0 = n[0];
418  n[0] = n[1];
419  n[1] = n0;
420  for(int nn = 0;nn!=num_nod_per_ele;nn++,ii++) {
421  conn_moab[ii] = n[nn];
422  }
423  }
424  }
425  break;
426  default:
427  {
428  int ii = 0;
429  for(int ee = 0;ee!=num_ele;ee++) {
430  for(int nn = 0;nn!=num_nod_per_ele;nn++,ii++) {
431  // cerr << conn_med[ii] << " ";
432  conn_moab[ii] = verts[conn_med[ii]-1];
433  }
434  // cerr << endl;
435  }
436  }
437  }
438  rval = iface->update_adjacencies(
439  starte,num_ele,num_nod_per_ele,conn_moab
440  ); CHKERRQ_MOAB(rval);
441 
442  Range ents(starte,starte+num_ele-1);
443  ierr = m_field.get_moab().add_entities(mesh_meshset,ents); CHKERRG(ierr);
444 
445  // get family for cells
446  {
447  std::vector<med_int> fam(num_ele,0);
448  if(MEDmeshEntityFamilyNumberRd(
449  fid,mesh_name,MED_NO_DT,MED_NO_IT,MED_CELL,type,&fam[0]
450  ) < 0) {
451  SETERRQ(
452  m_field.get_comm(),
454  "No family number for elements: using 0 as default family number"
455  );
456  }
457  // std::vector<med_int> ele_tags(num_ele);
458  // if(MEDmeshEntityNumberRd(
459  // fid, mesh_name,MED_NO_DT,MED_NO_IT,MED_CELL,type,&ele_tags[0]) < 0
460  // ) {
461  // ele_tags.clear();
462  // }
463  for(int j = 0; j < num_ele; j++) {
464  // cerr << ents[j] << " " /*<< ele_tags[j] << " "*/ << fam[j] << endl;
465  family_elem_map[fam[j]].insert(ents[j]);
466  }
467  }
468 
469  }
470 
471  if(MEDfileClose(fid) < 0) {
472  SETERRQ1(m_field.get_comm(),MOFEM_OPERATION_UNSUCCESSFUL,"Unable to close file '%s'",file.c_str());
473  }
474 
476  }
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:536
MoFEM::Core & cOre
core database
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:522
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:565
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:528
std::bitset< 32 > CubitBCType
Definition: Common.hpp:200
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Common.hpp:78
static med_geometrie_element moab2med_element_type(const EntityType type)
std::vector< EntityHandle > meshMeshsets
meshset for each mesh
#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...

◆ writeMed()

MoFEMErrorCode MoFEM::MedInterface::writeMed ( const string &  file,
int  verb = 1 
)

write MED file

Parameters
filefile name
verbverbosity level
Returns
error code

Definition at line 626 of file MedInterface.cpp.

626  {
627  Interface &m_field = cOre;
629  SETERRQ(m_field.get_comm(),MOFEM_NOT_IMPLEMENTED,"Not yet implemented");
631  }
MoFEM::Core & cOre
core database
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:522
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:528

Member Data Documentation

◆ cOre

MoFEM::Core& MoFEM::MedInterface::cOre
private

core database

Definition at line 132 of file MedInterface.hpp.

◆ fieldNames

std::map<std::string,FieldData> MoFEM::MedInterface::fieldNames

Definition at line 126 of file MedInterface.hpp.

◆ flgFile

PetscBool MoFEM::MedInterface::flgFile
private

true if file name given in command line

Definition at line 134 of file MedInterface.hpp.

◆ medFileName

std::string MoFEM::MedInterface::medFileName

MED file name.

Definition at line 127 of file MedInterface.hpp.

◆ meshMeshsets

std::vector<EntityHandle> MoFEM::MedInterface::meshMeshsets

meshset for each mesh

Definition at line 125 of file MedInterface.hpp.

◆ meshNames

std::vector<std::string> MoFEM::MedInterface::meshNames

list of meshes in MED file

Definition at line 124 of file MedInterface.hpp.


The documentation for this struct was generated from the following files: