v0.14.0
forces_and_sources_getting_orders_indices_atom_test.cpp
Go to the documentation of this file.
1 /** \file forces_and_sources_getting_orders_indices_atom_test.cpp
2  \brief Atome test for getting orders on entities
3 */
4 
5 
6 
7 #include <MoFEM.hpp>
8 
9 namespace bio = boost::iostreams;
10 using bio::stream;
11 using bio::tee_device;
12 
13 using namespace MoFEM;
14 
15 static char help[] = "...\n\n";
16 
17 int main(int argc, char *argv[]) {
18 
19  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
20 
21  try {
22 
23  moab::Core mb_instance;
24  moab::Interface &moab = mb_instance;
25  int rank;
26  MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
27 
28  PetscBool flg = PETSC_TRUE;
29  char mesh_file_name[255];
30 #if PETSC_VERSION_GE(3, 6, 4)
31  CHKERR PetscOptionsGetString(PETSC_NULL, "", "-my_file", mesh_file_name,
32  255, &flg);
33 #else
34  CHKERR PetscOptionsGetString(PETSC_NULL, PETSC_NULL, "-my_file",
35  mesh_file_name, 255, &flg);
36 #endif
37  if (flg != PETSC_TRUE) {
38  SETERRQ(PETSC_COMM_SELF, 1, "*** ERROR -my_file (MESH FILE NEEDED)");
39  }
40 
41  // Create MoFEM (Joseph) database
42  MoFEM::Core core(moab);
43  MoFEM::Interface &m_field = core;
44 
45  const char *option;
46  option = "";
47  CHKERR moab.load_file(mesh_file_name, 0, option);
48 
49  // set entitities bit level
50  BitRefLevel bit_level0;
51  bit_level0.set(0);
52  EntityHandle meshset_level0;
53  CHKERR moab.create_meshset(MESHSET_SET, meshset_level0);
54  CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevelByDim(
55  0, 3, bit_level0);
56 
57  // Fields
58  CHKERR m_field.add_field("FIELD1", H1, AINSWORTH_LEGENDRE_BASE, 1);
59  CHKERR m_field.add_field("FIELD2", H1, AINSWORTH_LEGENDRE_BASE, 3);
60 
61  // FE
62  CHKERR m_field.add_finite_element("TEST_FE");
63 
64  // Define rows/cols and element data
65  CHKERR m_field.modify_finite_element_add_field_row("TEST_FE", "FIELD1");
66  CHKERR m_field.modify_finite_element_add_field_col("TEST_FE", "FIELD2");
67  CHKERR m_field.modify_finite_element_add_field_data("TEST_FE", "FIELD1");
68  CHKERR m_field.modify_finite_element_add_field_data("TEST_FE", "FIELD2");
69 
70  // Problem
71  CHKERR m_field.add_problem("TEST_PROBLEM");
72 
73  // set finite elements for problem
74  CHKERR m_field.modify_problem_add_finite_element("TEST_PROBLEM", "TEST_FE");
75  // set refinement level for problem
76  CHKERR m_field.modify_problem_ref_level_add_bit("TEST_PROBLEM", bit_level0);
77 
78  // meshset consisting all entities in mesh
79  EntityHandle root_set = moab.get_root_set();
80  // add entities to field
81  CHKERR m_field.add_ents_to_field_by_type(root_set, MBTET, "FIELD1");
82  CHKERR m_field.add_ents_to_field_by_type(root_set, MBTET, "FIELD2");
83  // add entities to finite element
84  CHKERR m_field.add_ents_to_finite_element_by_type(root_set, MBTET,
85  "TEST_FE");
86 
87  // set app. order
88  // see Hierarchic Finite Element Bases on Unstructured Tetrahedral Meshes
89  // (Mark Ainsworth & Joe Coyle)
90  int order = 5;
91  CHKERR m_field.set_field_order(root_set, MBTET, "FIELD1", order);
92  CHKERR m_field.set_field_order(root_set, MBTRI, "FIELD1", order);
93  CHKERR m_field.set_field_order(root_set, MBEDGE, "FIELD1", order);
94  CHKERR m_field.set_field_order(root_set, MBVERTEX, "FIELD1", 1);
95  CHKERR m_field.set_field_order(root_set, MBTET, "FIELD2", order);
96  CHKERR m_field.set_field_order(root_set, MBTRI, "FIELD2", order);
97  CHKERR m_field.set_field_order(root_set, MBEDGE, "FIELD2", order);
98  CHKERR m_field.set_field_order(root_set, MBVERTEX, "FIELD2", 1);
99 
100  /****/
101  // build database
102  // build field
103  CHKERR m_field.build_fields();
104  // build finite elemnts
105  CHKERR m_field.build_finite_elements();
106  // build adjacencies
107  CHKERR m_field.build_adjacencies(bit_level0);
108 
109  /****/
110  // build problem
111  ProblemsManager *prb_mng_ptr;
112  CHKERR m_field.getInterface(prb_mng_ptr);
113  CHKERR prb_mng_ptr->buildProblem("TEST_PROBLEM", false);
114  // mesh partitioning
115  // partition
116  CHKERR prb_mng_ptr->partitionSimpleProblem("TEST_PROBLEM");
117  CHKERR prb_mng_ptr->partitionFiniteElements("TEST_PROBLEM");
118  // what are ghost nodes, see Petsc Manual
119  CHKERR prb_mng_ptr->partitionGhostDofs("TEST_PROBLEM");
120 
121  struct ForcesAndSourcesCore_TestFE : public ForcesAndSourcesCore {
122 
123  typedef tee_device<std::ostream, std::ofstream> TeeDevice;
124  typedef stream<TeeDevice> TeeStream;
125 
126  std::ofstream ofs;
127  TeeDevice my_tee;
128  TeeStream my_split;
129 
130  boost::shared_ptr<EntitiesFieldData> dataPtr;
131  boost::shared_ptr<DerivedEntitiesFieldData> derivedDataPtr;
132 
133  ForcesAndSourcesCore_TestFE(MoFEM::Interface &_m_field)
134  : ForcesAndSourcesCore(_m_field),
135  ofs("forces_and_sources_getting_orders_indices_atom_test.txt"),
136  my_tee(std::cout, ofs), my_split(my_tee),
137  dataPtr(new EntitiesFieldData(MBTET)),
138  derivedDataPtr(new DerivedEntitiesFieldData(dataPtr)){};
139 
140  MoFEMErrorCode preProcess() {
143  }
144 
145  MoFEMErrorCode operator()() {
147 
148  EntitiesFieldData &data = *dataPtr;
149  DerivedEntitiesFieldData &derived_data = *derivedDataPtr;
150 
151  my_split << "\n\nNEXT ELEM\n\n";
152 
153  CHKERR getSpacesAndBaseOnEntities(data);
154 
155  CHKERR getEntitySense<MBEDGE>(data);
156  CHKERR getEntitySense<MBTRI>(data);
157  CHKERR getEntityDataOrder<MBEDGE>(data, H1);
158  CHKERR getEntityDataOrder<MBTRI>(data, H1);
159  CHKERR getEntityDataOrder<MBTET>(data, H1);
160  CHKERR getFaceNodes(data);
161 
162  MatrixDouble gauss_pts(4, 4);
163  for (int gg = 0; gg < 4; gg++) {
164  gauss_pts(0, gg) = G_TET_X4[gg];
165  gauss_pts(1, gg) = G_TET_Y4[gg];
166  gauss_pts(2, gg) = G_TET_Z4[gg];
167  gauss_pts(3, gg) = G_TET_W4[gg];
168  }
170  gauss_pts,
171  boost::shared_ptr<BaseFunctionCtx>(new EntPolynomialBaseCtx(
173 
174  const auto bn1 = mField.get_field_bit_number("FIELD1");
175  CHKERR getRowNodesIndices(data, bn1);
176  CHKERR getEntityFieldData(data, bn1, MBEDGE);
177  CHKERR getEntityRowIndices(data, bn1, MBEDGE);
178  CHKERR getNodesFieldData(data, bn1);
179 
180  data.dataOnEntities[MBVERTEX][0].getFieldData().resize(0);
181  my_split << "FIELD1:\n";
182  my_split << data << std::endl;
183 
184  derived_data.dataOnEntities[MBVERTEX][0].getBase() =
186  const auto bn2 = mField.get_field_bit_number("FIELD2");
187  CHKERR getColNodesIndices(derived_data, bn2);
188  CHKERR getEntityFieldData(derived_data, bn2, MBEDGE);
189  CHKERR getEntityColIndices(derived_data, bn2, MBEDGE);
190  CHKERR getNodesFieldData(derived_data, bn2);
191 
192  derived_data.dataOnEntities[MBVERTEX][0].getFieldData().resize(0);
193 
194  my_split << "FIELD2:\n";
195  my_split << derived_data << std::endl;
196 
198  }
199 
200  MoFEMErrorCode postProcess() {
202 
203  my_split.close();
204 
206  }
207  };
208 
209  ForcesAndSourcesCore_TestFE fe1(m_field);
210  CHKERR m_field.loop_finite_elements("TEST_PROBLEM", "TEST_FE", fe1);
211 
212  }
213  CATCH_ERRORS;
214 
216 
217  return 0;
218 }
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
main
int main(int argc, char *argv[])
Definition: forces_and_sources_getting_orders_indices_atom_test.cpp:17
TeeStream
stream< TeeDevice > TeeStream
Definition: forces_and_sources_testing_contact_prism_element.cpp:10
MoFEM::CoreInterface::loop_finite_elements
virtual MoFEMErrorCode loop_finite_elements(const std::string problem_name, const std::string &fe_name, FEMethod &method, boost::shared_ptr< NumeredEntFiniteElement_multiIndex > fe_ptr=nullptr, MoFEMTypes bh=MF_EXIST, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr(), int verb=DEFAULT_VERBOSITY)=0
Make a loop over finite elements.
MoFEM::ProblemsManager::buildProblem
MoFEMErrorCode buildProblem(const std::string name, const bool square_matrix, int verb=VERBOSE)
build problem data structures
Definition: ProblemsManager.cpp:87
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
H1
@ H1
continuous field
Definition: definitions.h:85
EntityHandle
MoFEM::ProblemsManager
Problem manager is used to build and partition problems.
Definition: ProblemsManager.hpp:21
MoFEM::EntPolynomialBaseCtx
Class used to pass element data to calculate base functions on tet,triangle,edge.
Definition: EntPolynomialBaseCtx.hpp:22
G_TET_Z4
static const double G_TET_Z4[]
Definition: fem_tools.h:1120
MoFEM::CoreInterface::modify_finite_element_add_field_row
virtual MoFEMErrorCode modify_finite_element_add_field_row(const std::string &fe_name, const std::string name_row)=0
set field row which finite element use
help
static char help[]
Definition: forces_and_sources_getting_orders_indices_atom_test.cpp:15
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::Types::MatrixDouble
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
G_TET_Y4
static const double G_TET_Y4[]
Definition: fem_tools.h:1118
MoFEM::DerivedEntitiesFieldData
this class derive data form other data structure
Definition: EntitiesFieldData.hpp:110
MoFEM.hpp
MoFEM::CoreTmp< 0 >::Finalize
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:112
MoFEM::CoreInterface::add_ents_to_field_by_type
virtual MoFEMErrorCode add_ents_to_field_by_type(const Range &ents, const EntityType type, const std::string &name, int verb=DEFAULT_VERBOSITY)=0
Add entities to field meshset.
order
constexpr int order
Definition: dg_projection.cpp:18
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
TeeDevice
tee_device< std::ostream, std::ofstream > TeeDevice
Definition: forces_and_sources_testing_contact_prism_element.cpp:9
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:2010
MoFEM::CoreInterface::add_ents_to_finite_element_by_type
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::CoreInterface::add_finite_element
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
MoFEM::CoreInterface::modify_finite_element_add_field_col
virtual MoFEMErrorCode modify_finite_element_add_field_col(const std::string &fe_name, const std::string name_row)=0
set field col which finite element use
MoFEM::CoreInterface::build_finite_elements
virtual MoFEMErrorCode build_finite_elements(int verb=DEFAULT_VERBOSITY)=0
Build finite elements.
MoFEM::CoreInterface::modify_finite_element_add_field_data
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string name_field)=0
set finite element field data
MoFEM::ProblemsManager::partitionFiniteElements
MoFEMErrorCode partitionFiniteElements(const std::string name, bool part_from_moab=false, int low_proc=-1, int hi_proc=-1, int verb=VERBOSE)
partition finite elements
Definition: ProblemsManager.cpp:2167
MoFEM::CoreInterface::modify_problem_ref_level_add_bit
virtual MoFEMErrorCode modify_problem_ref_level_add_bit(const std::string &name_problem, const BitRefLevel &bit)=0
add ref level to problem
mesh_file_name
char mesh_file_name[255]
Definition: mesh_smoothing.cpp:22
MoFEM::ForcesAndSourcesCore
structure to get information form mofem into EntitiesFieldData
Definition: ForcesAndSourcesCore.hpp:22
G_TET_X4
static const double G_TET_X4[]
Definition: fem_tools.h:1116
G_TET_W4
static const double G_TET_W4[]
Definition: fem_tools.h:1122
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
CATCH_ERRORS
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
MoFEM::Core
CoreTmp< 0 > Core
Definition: Core.hpp:1148
MoFEM::ProblemsManager::partitionSimpleProblem
MoFEMErrorCode partitionSimpleProblem(const std::string name, int verb=VERBOSE)
partition problem dofs
Definition: ProblemsManager.cpp:1537
AINSWORTH_LEGENDRE_BASE
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
MoFEM::PetscOptionsGetString
PetscErrorCode PetscOptionsGetString(PetscOptions *, const char pre[], const char name[], char str[], size_t size, PetscBool *set)
Definition: DeprecatedPetsc.hpp:172
MoFEM::CoreInterface::build_fields
virtual MoFEMErrorCode build_fields(int verb=DEFAULT_VERBOSITY)=0
MoFEM::TetPolynomialBase::getValue
MoFEMErrorCode getValue(MatrixDouble &pts, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
Definition: TetPolynomialBase.cpp:2021
MoFEM::EntitiesFieldData::dataOnEntities
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
Definition: EntitiesFieldData.hpp:57
MoFEM::BitRefManager
Managing BitRefLevels.
Definition: BitRefManager.hpp:21
MoFEM::CoreInterface::modify_problem_add_finite_element
virtual MoFEMErrorCode modify_problem_add_finite_element(const std::string name_problem, const std::string &fe_name)=0
add finite element to problem, this add entities assigned to finite element to a particular problem
MoFEM::Types::BitRefLevel
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
MoFEM::CoreInterface::build_adjacencies
virtual MoFEMErrorCode build_adjacencies(const Range &ents, int verb=DEFAULT_VERBOSITY)=0
build adjacencies
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
MoFEM::CoreInterface::set_field_order
virtual MoFEMErrorCode set_field_order(const EntityHandle meshset, const EntityType type, const std::string &name, const ApproximationOrder order, int verb=DEFAULT_VERBOSITY)=0
Set order approximation of the entities in the field.
MoFEM::EntitiesFieldData
data structure for finite element entity
Definition: EntitiesFieldData.hpp:40
MoFEM::ProblemsManager::partitionGhostDofs
MoFEMErrorCode partitionGhostDofs(const std::string name, int verb=VERBOSE)
determine ghost nodes
Definition: ProblemsManager.cpp:2339
CONTINUOUS
@ CONTINUOUS
Regular field.
Definition: definitions.h:100
MoFEM::CoreInterface::add_problem
virtual MoFEMErrorCode add_problem(const std::string &name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
Add problem.
MoFEM::CoreInterface::add_field
virtual MoFEMErrorCode add_field(const std::string name, const FieldSpace space, const FieldApproximationBase base, const FieldCoefficientsNumber nb_of_coefficients, const TagType tag_type=MB_TAG_SPARSE, const enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
Add field.
MoFEM::TetPolynomialBase
Calculate base functions on tetrahedral.
Definition: TetPolynomialBase.hpp:17