v0.14.0
Functions | Variables
partition_mesh.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[] 
)
Examples
partition_mesh.cpp.

Definition at line 15 of file partition_mesh.cpp.

15  {
16 
17  // initialize petsc
18  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
19 
20  int nb_vertices;
21 
22  try {
23 
24  PetscBool flg = PETSC_TRUE;
25  char mesh_file_name[255];
26 #if PETSC_VERSION_GE(3, 6, 4)
27  CHKERR PetscOptionsGetString(PETSC_NULL, "", "-my_file", mesh_file_name,
28  255, &flg);
29 #else
30  CHKERR PetscOptionsGetString(PETSC_NULL, PETSC_NULL, "-my_file",
31  mesh_file_name, 255, &flg);
32 #endif
33  if (flg != PETSC_TRUE)
34  SETERRQ(PETSC_COMM_SELF, 1, "*** ERROR -my_file (MESH FILE NEEDED)");
35 
36  // read mesh and create moab and mofem datastrutures
37 
38  moab::Core mb_instance;
39  moab::Interface &moab = mb_instance;
40  {
41  const char *option;
42  option = "";
43  CHKERR moab.load_file(mesh_file_name, 0, option);
44  }
45  MoFEM::Core core(moab);
46  MoFEM::Interface &m_field = core;
47 
48  CHKERR moab.get_number_entities_by_dimension(0, 0, nb_vertices, true);
49 
50  EntityHandle root_set = moab.get_root_set();
51  Range tets;
52  moab.get_entities_by_type(root_set, MBTET, tets, false);
53 
54  Tag th_vertex_weight;
55  int def_val = 1;
56  CHKERR moab.tag_get_handle("VERTEX_WEIGHT", 1, MB_TYPE_INTEGER,
57  th_vertex_weight, MB_TAG_CREAT | MB_TAG_DENSE,
58  &def_val);
59 
60  CommInterface *comm_interafce_ptr = m_field.getInterface<CommInterface>();
61  CHKERR comm_interafce_ptr->partitionMesh(
62  tets, 3, 2, m_field.get_comm_size(), &th_vertex_weight, NULL, NULL,
63  VERBOSE, false);
64 
65  if (!m_field.get_comm_rank()) {
66  CHKERR moab.write_file("partitioned_mesh.h5m");
67  }
68  }
70 
71  PetscBarrier(PETSC_NULL);
72 
73  try {
74 
75  moab::Core mb_instance2;
76  moab::Interface &moab2 = mb_instance2;
77 
78 
79  MoFEM::CoreTmp<1> core(moab2);
80  MoFEM::Interface &m_field = core;
81 
82  // Register DM Manager
83  DMType dm_name = "DMMOFEM";
84  CHKERR DMRegister_MoFEM(dm_name);
85 
86  // Test build simple problem
87  const char *option = "DEBUG_IO;"
88  "PARALLEL=READ_PART;"
89  "PARALLEL_RESOLVE_SHARED_ENTS;"
90  "PARTITION=PARALLEL_PARTITION;";
91 
92  CHKERR m_field.getInterface<Simple>()->getOptions();
93  CHKERR m_field.getInterface<Simple>()->loadFile(option,
94  "partitioned_mesh.h5m");
95  CHKERR m_field.getInterface<Simple>()->addDomainField(
96  "U", H1, AINSWORTH_LEGENDRE_BASE, 1);
97  CHKERR m_field.getInterface<Simple>()->setFieldOrder("U", 1);
98  CHKERR m_field.getInterface<Simple>()->setUp();
99 
100  auto dm = m_field.getInterface<Simple>()->getDM();
101 
102  const MoFEM::Problem *problem_ptr;
103  CHKERR DMMoFEMGetProblemPtr(dm, &problem_ptr);
104 
105  if(problem_ptr->nbDofsRow != nb_vertices)
106  SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
107  "Number of vertices and DOFs is inconstent");
108 
109  MOFEM_LOG_CHANNEL("WORLD");
110  MOFEM_TAG_AND_LOG("WORLD", Sev::inform, "Atom test")
111  << "All is good in this test";
112  }
113  CATCH_ERRORS;
114 
116 
117  return 0;
118 }

Variable Documentation

◆ help

char help[] = "...\n\n"
static
Examples
partition_mesh.cpp.

Definition at line 13 of file partition_mesh.cpp.

MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
H1
@ H1
continuous field
Definition: definitions.h:85
MOFEM_LOG_CHANNEL
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Definition: LogManager.hpp:284
EntityHandle
MoFEM::CoreInterface::get_comm_rank
virtual int get_comm_rank() const =0
MoFEM::CommInterface::partitionMesh
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.
Definition: CommInterface.cpp:868
MoFEM::CoreTmp< 0 >::Finalize
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:112
MoFEM::Simple
Simple interface for fast problem set-up.
Definition: Simple.hpp:27
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:2010
VERBOSE
@ VERBOSE
Definition: definitions.h:222
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::CoreTmp
Definition: Core.hpp:36
help
static char help[]
Definition: partition_mesh.cpp:13
MoFEM::DMRegister_MoFEM
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMoFEM.cpp:43
MoFEM::CoreInterface::get_comm_size
virtual int get_comm_size() const =0
mesh_file_name
char mesh_file_name[255]
Definition: mesh_smoothing.cpp:22
Range
MoFEM::CoreTmp< 0 >::Initialize
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:72
MOFEM_TAG_AND_LOG
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
Definition: LogManager.hpp:362
CATCH_ERRORS
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
MoFEM::CommInterface
Managing BitRefLevels.
Definition: CommInterface.hpp:21
MoFEM::Core
CoreTmp< 0 > Core
Definition: Core.hpp:1148
AINSWORTH_LEGENDRE_BASE
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
MoFEM::Problem::nbDofsRow
DofIdx nbDofsRow
Global number of DOFs in row.
Definition: ProblemsMultiIndices.hpp:65
MoFEM::PetscOptionsGetString
PetscErrorCode PetscOptionsGetString(PetscOptions *, const char pre[], const char name[], char str[], size_t size, PetscBool *set)
Definition: DeprecatedPetsc.hpp:172
MoFEM::DMMoFEMGetProblemPtr
PetscErrorCode DMMoFEMGetProblemPtr(DM dm, const MoFEM::Problem **problem_ptr)
Get pointer to problem data structure.
Definition: DMMoFEM.cpp:426
MOFEM_ATOM_TEST_INVALID
@ MOFEM_ATOM_TEST_INVALID
Definition: definitions.h:40
MoFEM::Problem
keeps basic data about problem
Definition: ProblemsMultiIndices.hpp:54