v0.13.0
Functions | Variables
mofem_part.cpp File Reference
#include <MoFEM.hpp>

Go to the source code of this file.

Functions

int main (int argc, char *argv[])
 

Variables

static char help [] = "...\n\n"
 

Function Documentation

◆ main()

int main ( int  argc,
char *  argv[] 
)

Definition at line 26 of file mofem_part.cpp.

26  {
27 
28  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
29 
30  try {
31 
32  moab::Core mb_instance;
33  moab::Interface &moab = mb_instance;
34  ParallelComm *pcomm = ParallelComm::get_pcomm(&moab, MYPCOMM_INDEX);
35  if (pcomm == NULL)
36  pcomm = new ParallelComm(&moab, PETSC_COMM_WORLD);
37 
38  // global variables
39  char mesh_file_name[255];
40  char mesh_out_file[255] = "out.h5m";
41  PetscBool flg_file = PETSC_FALSE;
42  PetscBool flg_n_part = PETSC_FALSE;
43  PetscBool flg_part = PETSC_FALSE;
44  PetscBool only_tags = PETSC_FALSE;
45  PetscInt n_partas = 1;
46  PetscBool create_lower_dim_ents = PETSC_TRUE;
47  PetscInt dim = 3;
48  PetscInt adj_dim = 2;
49 
50  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "", "none", "none");
51  CHKERRQ(ierr);
52 
53  CHKERR PetscOptionsString("-my_file", "mesh file name", "", "mesh.h5m",
54  mesh_file_name, 255, &flg_file);
55  if (flg_file != PETSC_TRUE)
56  CHKERR PetscOptionsString("-file_name", "mesh file name", "", "mesh.h5m",
57  mesh_file_name, 255, &flg_file);
58  if (flg_file != PETSC_TRUE)
59  SETERRQ(PETSC_COMM_SELF, 1,
60  "*** ERROR -my_file (-file_name) (MESH FILE NEEDED)");
61 
62  const char *option;
63  option = "";
64  CHKERR moab.load_file(mesh_file_name, 0, option);
65 
66  CHKERR PetscOptionsString("-output_file", "output mesh file name", "",
67  "out.h5m", mesh_out_file, 255, PETSC_NULL);
68  CHKERR PetscOptionsInt("-my_nparts", "number of parts", "", n_partas,
69  &n_partas, &flg_n_part);
70  CHKERR PetscOptionsInt("-nparts", "number of parts", "", n_partas,
71  &n_partas, &flg_part);
72 
73  if (!flg_n_part && !flg_part)
74  SETERRQ(PETSC_COMM_SELF, 1,
75  "*** ERROR partitioning number not given (-nparts)");
76 
77  auto get_nb_ents_by_dim = [&](const int dim) {
78  int nb;
79  CHKERR moab.get_number_entities_by_dimension(0, dim, nb);
80  return nb;
81  };
82  for (; dim >= 0; dim--) {
83  if (get_nb_ents_by_dim(dim))
84  break;
85  }
86 
87  if (!dim)
88  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
89  "Dimension of entities to partition not found");
90 
91  CHKERR PetscOptionsInt("-dim", "entities dim", "", dim, &dim, PETSC_NULL);
92  adj_dim = dim - 1;
93  CHKERR PetscOptionsInt("-adj_dim", "adjacency dim", "", adj_dim, &adj_dim,
94  PETSC_NULL);
95  CHKERR PetscOptionsBool(
96  "-my_create_lower_dim_ents", "if tru create lower dimension entireties",
97  "", create_lower_dim_ents, &create_lower_dim_ents, PETSC_NULL);
98  CHKERR PetscOptionsBool("-block_tags", "obly block and meshsests tags", "",
99  only_tags, &only_tags, PETSC_NULL);
100 
101  ierr = PetscOptionsEnd();
102  CHKERRQ(ierr);
103 
104  // Create MoFEM database
105  MoFEM::Core core(moab);
106  MoFEM::Interface &m_field = core;
107 
108  auto meshsets_interface_ptr = m_field.getInterface<MeshsetsManager>();
109  CHKERR meshsets_interface_ptr->setMeshsetFromFile();
110 
111  MOFEM_LOG_CHANNEL("WORLD");
112  MOFEM_LOG_TAG("WORLD", "mofem_part")
113  for (CubitMeshSet_multiIndex::iterator cit =
114  meshsets_interface_ptr->getBegin();
115  cit != meshsets_interface_ptr->getEnd(); cit++)
116  MOFEM_LOG("WORLD", Sev::inform) << *cit;
117  MOFEM_LOG_CHANNEL("WORLD");
118 
119  {
120  Range ents_dim;
121  CHKERR moab.get_entities_by_dimension(0, dim, ents_dim, false);
122  if (create_lower_dim_ents) {
123  if (dim == 3) {
124  Range faces;
125  CHKERR moab.get_adjacencies(ents_dim, 2, true, faces,
126  moab::Interface::UNION);
127  }
128  if (dim > 2) {
129  Range edges;
130  CHKERR moab.get_adjacencies(ents_dim, 1, true, edges,
131  moab::Interface::UNION);
132  }
133  }
134  CommInterface *comm_interafce_ptr = m_field.getInterface<CommInterface>();
135  CHKERR comm_interafce_ptr->partitionMesh(
136  ents_dim, dim, adj_dim, n_partas, nullptr, nullptr, nullptr, VERBOSE);
137  }
138 
139  auto get_tag_list = [&]() {
140  std::vector<Tag> tags_list;
141  auto meshsets_mng = m_field.getInterface<MeshsetsManager>();
142  auto &list = meshsets_mng->getMeshsetsMultindex();
143  for (auto &m : list) {
144  auto meshset = m.getMeshset();
145  std::vector<Tag> tmp_tags_list;
146  CHKERR m_field.get_moab().tag_get_tags_on_entity(meshset,
147  tmp_tags_list);
148  for (auto t : tmp_tags_list) {
149  tags_list.push_back(t);
150  }
151  }
152 
153  return tags_list;
154  };
155 
156  EntityHandle root_mesh = 0;
157  if (m_field.get_comm_rank() == 0) {
158  if (only_tags) {
159  auto tags_list = get_tag_list();
160  std::sort(tags_list.begin(), tags_list.end());
161  auto new_end = std::unique(tags_list.begin(), tags_list.end());
162  tags_list.resize(std::distance(tags_list.begin(), new_end));
163  tags_list.push_back(pcomm->part_tag());
164  // tags_list.push_back(m_field.get_moab().globalId_tag());
165  CHKERR moab.write_file(mesh_out_file, "MOAB", "", &root_mesh, 1,
166  &*tags_list.begin(), tags_list.size());
167  } else {
168  CHKERR moab.write_file(mesh_out_file, "MOAB", "");
169  }
170  }
171  }
172  CATCH_ERRORS;
173 
175  CHKERRQ(ierr);
176 
177  return 0;
178 }
@ VERBOSE
Definition: definitions.h:222
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:228
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
#define CHKERR
Inline error check.
Definition: definitions.h:548
const int dim
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:311
SeverityLevel
Severity levels.
Definition: LogManager.hpp:43
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:342
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Definition: LogManager.hpp:287
MoFEMErrorCode partitionMesh(const Range &ents, const int dim, const int adj_dim, const int n_parts, Tag *th_vertex_weights=nullptr, Tag *th_edge_weights=nullptr, Tag *th_part_weights=nullptr, int verb=VERBOSE, const bool debug=false)
Set partition tag to each finite element in the problem.
char mesh_file_name[255]
static char help[]
Definition: mofem_part.cpp:24
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
multi_index_container< CubitMeshSets, indexed_by< hashed_unique< tag< Meshset_mi_tag >, member< CubitMeshSets, EntityHandle, &CubitMeshSets::meshset > >, ordered_non_unique< tag< CubitMeshsetType_mi_tag >, const_mem_fun< CubitMeshSets, unsigned long int, &CubitMeshSets::getBcTypeULong > >, ordered_non_unique< tag< CubitMeshsetMaskedType_mi_tag >, const_mem_fun< CubitMeshSets, unsigned long int, &CubitMeshSets::getMaksedBcTypeULong > >, ordered_non_unique< tag< CubitMeshsets_name >, const_mem_fun< CubitMeshSets, std::string, &CubitMeshSets::getName > >, hashed_unique< tag< Composite_Cubit_msId_And_MeshsetType_mi_tag >, composite_key< CubitMeshSets, const_mem_fun< CubitMeshSets, int, &CubitMeshSets::getMeshsetId >, const_mem_fun< CubitMeshSets, unsigned long int, &CubitMeshSets::getMaksedBcTypeULong > > > > > CubitMeshSet_multiIndex
Stores data about meshsets (see CubitMeshSets) storing data about boundary conditions,...
CoreTmp< 0 > Core
Definition: Core.hpp:1096
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
constexpr double t
plate stiffness
Definition: plate.cpp:76
FTensor::Index< 'm', 3 > m
Managing BitRefLevels.
virtual moab::Interface & get_moab()=0
virtual int get_comm_rank() const =0
Core (interface) class.
Definition: Core.hpp:92
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:85
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:125
Deprecated interface functions.
Interface for managing meshsets containing materials and boundary conditions.
CubitMeshSet_multiIndex & getMeshsetsMultindex()
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.

Variable Documentation

◆ help

char help[] = "...\n\n"
static

Definition at line 24 of file mofem_part.cpp.