v0.14.0
Loading...
Searching...
No Matches
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
12using namespace MoFEM;
13
14static char help[] = "...\n\n";
15
16#define HelloFunctionBegin \
17 MoFEMFunctionBegin; \
18 MOFEM_LOG_CHANNEL("SYNC"); \
19 MOFEM_LOG_FUNCTION(); \
20 MOFEM_LOG_TAG("WORLD", "HelloWorld");
21
23 OpRow(const std::string &field_name)
28 if (type == MBVERTEX) {
29 // get number of evaluated element in the loop
30 MOFEM_LOG("SYNC", Sev::inform) << "**** " << getNinTheLoop() << " ****";
31 MOFEM_LOG("SYNC", Sev::inform) << "**** Operators ****";
32 }
33 MOFEM_LOG("SYNC", Sev::inform)
34 << "Hello Operator OpRow:"
35 << " field name " << rowFieldName << " side " << side << " type "
36 << CN::EntityTypeName(type) << " nb dofs on entity "
37 << data.getIndices().size();
39 }
40};
41
43 OpRowCol(const std::string row_field, const std::string col_field,
44 const bool symm)
45 : ForcesAndSourcesCore::UserDataOperator(row_field, col_field, OPROWCOL,
46 symm) {}
47 MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
48 EntityType col_type,
52 MOFEM_LOG("SYNC", Sev::inform)
53 << "Hello Operator OpRowCol:"
54 << " row field name " << rowFieldName << " row side " << row_side
55 << " row type " << CN::EntityTypeName(row_type)
56 << " nb dofs on row entity " << row_data.getIndices().size() << " : "
57 << " col field name " << colFieldName << " col side " << col_side
58 << " col type " << CN::EntityTypeName(col_type)
59 << " nb dofs on col entity " << col_data.getIndices().size();
61 }
62};
63
65 OpVolume(const std::string &field_name)
67 field_name, OPROW) {
68 }
72 if (type == MBVERTEX) {
73 MOFEM_LOG("SYNC", Sev::inform) << "Hello Operator OpVolume:"
74 << " volume " << getVolume();
75 }
77 }
78};
79
81 OpFace(const std::string &field_name)
86 if (type == MBVERTEX) {
87 MOFEM_LOG("SYNC", Sev::inform) << "Hello Operator OpFace:"
88 << " normal " << getNormal();
89 }
91 }
92};
93
95 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> &feSidePtr;
97 const std::string &field_name,
98 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> &fe_side_ptr)
100 feSidePtr(fe_side_ptr) {}
103
105 if (type == MBVERTEX) {
106 MOFEM_LOG("SYNC", Sev::inform) << "Hello Operator OpSideFace";
107 CHKERR loopSideVolumes("dFE", *feSidePtr);
108 }
110 }
111};
112
115 OpVolumeSide(const std::string &field_name)
117 field_name, field_name, OPROW) {}
121 if (type == MBVERTEX) {
122 MOFEM_LOG("SYNC", Sev::inform)
123 << "Hello Operator OpVolumeSide:"
124 << " volume " << getVolume() << " normal " << getNormal();
125 }
127 }
128};
129
130int main(int argc, char *argv[]) {
131
132 // initialize petsc
133 MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
134
135 try {
136
137 // Register DM Manager
138 DMType dm_name = "DMMOFEM";
139 CHKERR DMRegister_MoFEM(dm_name);
140
141 // Create MoAB database
142 moab::Core moab_core;
143 moab::Interface &moab = moab_core;
144
145 // Create MoFEM database and link it to MoAB
146 MoFEM::Core mofem_core(moab);
147 MoFEM::Interface &m_field = mofem_core;
148
149 // Simple interface
150 Simple *simple;
151 CHKERR m_field.getInterface(simple);
152
153 // get options from command line
154 CHKERR simple->getOptions();
155 // load mesh file
156 CHKERR simple->loadFile();
157 // add fields
158 CHKERR simple->addDomainField("U", H1, AINSWORTH_LEGENDRE_BASE, 3);
159 CHKERR simple->addBoundaryField("L", H1, AINSWORTH_LEGENDRE_BASE, 3);
160 CHKERR simple->addSkeletonField("S", H1, AINSWORTH_LEGENDRE_BASE, 3);
161 // set fields order
162 CHKERR simple->setFieldOrder("U", 4);
163 CHKERR simple->setFieldOrder("L", 3);
164 CHKERR simple->setFieldOrder("S", 3);
165
166 // setup problem
167 CHKERR simple->setUp();
168
169 auto pipeline_mng = m_field.getInterface<PipelineManager>();
170
171 //! [set operator to the volume element instance]
172 pipeline_mng->getOpDomainRhsPipeline().push_back(new OpRow("U"));
173 pipeline_mng->getOpDomainRhsPipeline().push_back(new OpVolume("U"));
174 pipeline_mng->getOpDomainLhsPipeline().push_back(
175 new OpRowCol("U", "U", true));
176 //! [set operator to the volume element instance]
177
178 //! [set operator to the face element instance]
179 pipeline_mng->getOpBoundaryRhsPipeline().push_back(new OpRow("L"));
180 pipeline_mng->getOpBoundaryRhsPipeline().push_back(new OpFace("L"));
181 pipeline_mng->getOpBoundaryLhsPipeline().push_back(
182 new OpRowCol("U", "L", false));
183 //! [set operator to the face element instance]
184
185 //! [create skeleton side element and push operator to it]
186 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe(
188 side_fe->getOpPtrVector().push_back(new OpVolumeSide("U"));
189 //! [create skeleton side element and push operator to it]
190
191 //! [set operator to the face element on skeleton instance]
192 pipeline_mng->getOpSkeletonRhsPipeline().push_back(new OpRow("S"));
193 pipeline_mng->getOpSkeletonRhsPipeline().push_back(
194 new OpFaceSide("S", side_fe));
195 //! [set operator to the face element on skeleton instance]
196
197 //! [executing finite elements]
198 CHKERR pipeline_mng->loopFiniteElements();
199 //! [executing finite elements]
200
202 }
204
205 // finish work cleaning memory, getting statistics, etc.
207
208 return 0;
209}
ForcesAndSourcesCore::UserDataOperator UserDataOperator
#define MOFEM_LOG_SYNCHRONISE(comm)
Synchronise "SYNC" channel.
Definition: LogManager.hpp:345
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition: acoustic.cpp:69
int main()
Definition: adol-c_atom.cpp:46
#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: DMMoFEM.cpp:47
boost::ptr_deque< UserDataOperator > & getOpDomainRhsPipeline()
Get the Op Domain Rhs Pipeline object.
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
#define HelloFunctionBegin
Definition: hello_world.cpp:16
static char help[]
Definition: hello_world.cpp:14
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
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:27
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:81
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:83
boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > & feSidePtr
Definition: hello_world.cpp:95
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:96
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:47
OpRowCol(const std::string row_field, const std::string col_field, const bool symm)
Definition: hello_world.cpp:43
OpRow(const std::string &field_name)
Definition: hello_world.cpp:23
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Definition: hello_world.cpp:25
OpVolume(const std::string &field_name)
Definition: hello_world.cpp:65
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Definition: hello_world.cpp:69
OpVolumeSide(const std::string &field_name)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)