v0.13.1
hello_world.cpp
Go to the documentation of this file.
1/**
2 * \file hello_world.cpp
3 * \ingroup mofem_simple_interface
4 * \example hello_world.cpp
5 *
6 * Prints basic information about users data operator.
7 * See more details in \ref hello_world_tut1
8 *
9 */
10
11#include <BasicFiniteElements.hpp>
12using namespace MoFEM;
13
14static char help[] = "...\n\n";
15
16static map<EntityType, std::string> type_name;
17
18#define HelloFunctionBegin \
19 MoFEMFunctionBegin; \
20 MOFEM_LOG_CHANNEL("SYNC"); \
21 MOFEM_LOG_FUNCTION(); \
22 MOFEM_LOG_TAG("WORLD", "HelloWorld");
23
25 OpRow(const std::string &field_name)
30 if (type == MBVERTEX) {
31 // get number of evaluated element in the loop
32 MOFEM_LOG("SYNC", Sev::inform) << "**** " << getNinTheLoop() << " ****";
33 MOFEM_LOG("SYNC", Sev::inform) << "**** Operators ****";
34 }
35 MOFEM_LOG("SYNC", Sev::inform)
36 << "Hello Operator OpRow:"
37 << " field name " << rowFieldName << " side " << side << " type "
38 << type_name[type] << " nb dofs on entity " << data.getIndices().size();
40 }
41};
42
44 OpRowCol(const std::string row_field, const std::string col_field,
45 const bool symm)
46 : ForcesAndSourcesCore::UserDataOperator(row_field, col_field, OPROWCOL,
47 symm) {}
48 MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
49 EntityType col_type,
53 MOFEM_LOG("SYNC", Sev::inform)
54 << "Hello Operator OpRowCol:"
55 << " row field name " << rowFieldName << " row side " << row_side
56 << " row type " << type_name[row_type] << " nb dofs on row entity"
57 << row_data.getIndices().size() << " : "
58 << " col field name " << colFieldName << " col side " << col_side
59 << " col type " << type_name[col_type] << " nb dofs on col entity"
60 << col_data.getIndices().size();
62 }
63};
64
66 OpVolume(const std::string &field_name)
68 field_name, OPROW) {
69 }
73 if (type == MBVERTEX) {
74 MOFEM_LOG("SYNC", Sev::inform) << "Hello Operator OpVolume:"
75 << " volume " << getVolume();
76 }
78 }
79};
80
82 OpFace(const std::string &field_name)
87 if (type == MBVERTEX) {
88 MOFEM_LOG("SYNC", Sev::inform) << "Hello Operator OpFace:"
89 << " normal " << getNormal();
90 }
92 }
93};
94
96 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> &feSidePtr;
98 const std::string &field_name,
99 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> &fe_side_ptr)
101 feSidePtr(fe_side_ptr) {}
104
106 if (type == MBVERTEX) {
107 MOFEM_LOG("SYNC", Sev::inform) << "Hello Operator OpSideFace";
108 CHKERR loopSideVolumes("dFE", *feSidePtr);
109 }
111 }
112};
113
116 OpVolumeSide(const std::string &field_name)
118 field_name, field_name, OPROW) {}
122 if (type == MBVERTEX) {
123 MOFEM_LOG("SYNC", Sev::inform)
124 << "Hello Operator OpVolumeSide:"
125 << " volume " << getVolume() << " normal " << getNormal();
126 }
128 }
129};
130
131int main(int argc, char *argv[]) {
132
133 // initialize petsc
134 MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
135
136 try {
137
138 type_name[MBVERTEX] = "Vertex";
139 type_name[MBEDGE] = "Edge";
140 type_name[MBTRI] = "Triangle";
141 type_name[MBTET] = "Tetrahedra";
142
143 // Register DM Manager
144 DMType dm_name = "DMMOFEM";
145 CHKERR DMRegister_MoFEM(dm_name);
146
147 // Create MoAB database
148 moab::Core moab_core;
149 moab::Interface &moab = moab_core;
150
151 // Create MoFEM database and link it to MoAB
152 MoFEM::Core mofem_core(moab);
153 MoFEM::Interface &m_field = mofem_core;
154
155 // Simple interface
156 Simple *simple;
157 CHKERR m_field.getInterface(simple);
158
159 // get options from command line
160 CHKERR simple->getOptions();
161 // load mesh file
162 CHKERR simple->loadFile();
163 // add fields
164 CHKERR simple->addDomainField("U", H1, AINSWORTH_LEGENDRE_BASE, 3);
165 CHKERR simple->addBoundaryField("L", H1, AINSWORTH_LEGENDRE_BASE, 3);
166 CHKERR simple->addSkeletonField("S", H1, AINSWORTH_LEGENDRE_BASE, 3);
167 // set fields order
168 CHKERR simple->setFieldOrder("U", 4);
169 CHKERR simple->setFieldOrder("L", 3);
170 CHKERR simple->setFieldOrder("S", 3);
171
172 // setup problem
173 CHKERR simple->setUp();
174
175 auto pipeline_mng = m_field.getInterface<PipelineManager>();
176
177 //! [set operator to the volume element instance]
178 pipeline_mng->getOpDomainRhsPipeline().push_back(new OpRow("U"));
179 pipeline_mng->getOpDomainRhsPipeline().push_back(new OpVolume("U"));
180 pipeline_mng->getOpDomainLhsPipeline().push_back(
181 new OpRowCol("U", "U", true));
182 //! [set operator to the volume element instance]
183
184 //! [set operator to the face element instance]
185 pipeline_mng->getOpBoundaryRhsPipeline().push_back(new OpRow("L"));
186 pipeline_mng->getOpBoundaryRhsPipeline().push_back(new OpFace("L"));
187 pipeline_mng->getOpBoundaryLhsPipeline().push_back(
188 new OpRowCol("U", "L", false));
189 //! [set operator to the face element instance]
190
191 //! [create skeleton side element and push operator to it]
192 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe(
194 side_fe->getOpPtrVector().push_back(new OpVolumeSide("U"));
195 //! [create skeleton side element and push operator to it]
196
197 //! [set operator to the face element on skeleton instance]
198 pipeline_mng->getOpSkeletonRhsPipeline().push_back(new OpRow("S"));
199 pipeline_mng->getOpSkeletonRhsPipeline().push_back(
200 new OpFaceSide("S", side_fe));
201 //! [set operator to the face element on skeleton instance]
202
203 //! [executing finite elements]
204 CHKERR pipeline_mng->loopFiniteElements();
205 //! [executing finite elements]
206
208 }
210
211 // finish work cleaning memory, getting statistics, etc.
213
214 return 0;
215}
ForcesAndSourcesCore::UserDataOperator UserDataOperator
#define MOFEM_LOG_SYNCHRONISE(comm)
Synchronise "SYNC" channel.
Definition: LogManager.hpp:338
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition: acoustic.cpp:69
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:372
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
@ H1
continuous field
Definition: definitions.h:85
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:47
boost::ptr_vector< UserDataOperator > & getOpDomainRhsPipeline()
Get the Op Domain Rhs Pipeline object.
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:301
int main(int argc, char *argv[])
#define HelloFunctionBegin
Definition: hello_world.cpp:18
static char help[]
Definition: hello_world.cpp:14
static map< EntityType, std::string > type_name
Definition: hello_world.cpp:16
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: MoFEM.hpp:24
CoreTmp< 0 > Core
Definition: Core.hpp:1086
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1955
FaceSide::UserDataOperator OpFaceSide
constexpr auto field_name
virtual MPI_Comm & get_comm() const =0
Core (interface) class.
Definition: Core.hpp:82
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
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:112
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
const VectorInt & getIndices() const
Get global indices of dofs on entity.
@ OPROW
operator doWork function is executed on FE rows
structure to get information form mofem into EntitiesFieldData
PipelineManager interface.
Simple interface for fast problem set-up.
Definition: Simple.hpp:26
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Base volume element used to integrate on skeleton.
OpFace(const std::string &field_name)
Definition: hello_world.cpp:82
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: hello_world.cpp:84
boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > & feSidePtr
Definition: hello_world.cpp:96
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
OpFaceSide(const std::string &field_name, boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > &fe_side_ptr)
Definition: hello_world.cpp:97
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
Definition: hello_world.cpp:48
OpRowCol(const std::string row_field, const std::string col_field, const bool symm)
Definition: hello_world.cpp:44
OpRow(const std::string &field_name)
Definition: hello_world.cpp:25
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Definition: hello_world.cpp:27
OpVolume(const std::string &field_name)
Definition: hello_world.cpp:66
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Definition: hello_world.cpp:70
OpVolumeSide(const std::string &field_name)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)