v0.13.0
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 /* This file is part of MoFEM.
12 * MoFEM is free software: you can redistribute it and/or modify it under
13 * the terms of the GNU Lesser General Public License as published by the
14 * Free Software Foundation, either version 3 of the License, or (at your
15 * option) any later version.
16 *
17 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
18 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
19 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
20 
21 * License for more details.
22 *
23 * You should have received a copy of the GNU Lesser General Public
24 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
25 
26 #include <BasicFiniteElements.hpp>
27 using namespace MoFEM;
28 
29 static char help[] = "...\n\n";
30 
31 static map<EntityType, std::string> type_name;
32 
33 #define HelloFunctionBegin \
34  MoFEMFunctionBegin; \
35  MOFEM_LOG_CHANNEL("SYNC"); \
36  MOFEM_LOG_FUNCTION(); \
37  MOFEM_LOG_TAG("WORLD", "HelloWorld");
38 
40  OpRow(const std::string &field_name)
41  : ForcesAndSourcesCore::UserDataOperator(field_name, field_name, OPROW) {}
42  MoFEMErrorCode doWork(int side, EntityType type,
45  if (type == MBVERTEX) {
46  // get number of evaluated element in the loop
47  MOFEM_LOG("SYNC", Sev::inform) << "**** " << getNinTheLoop() << " ****";
48  MOFEM_LOG("SYNC", Sev::inform) << "**** Operators ****";
49  }
50  MOFEM_LOG("SYNC", Sev::inform)
51  << "Hello Operator OpRow:"
52  << " field name " << rowFieldName << " side " << side << " type "
53  << type_name[type] << " nb dofs on entity " << data.getIndices().size();
55  }
56 };
57 
59  OpRowCol(const std::string row_field, const std::string col_field,
60  const bool symm)
61  : ForcesAndSourcesCore::UserDataOperator(row_field, col_field, OPROWCOL,
62  symm) {}
63  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
64  EntityType col_type,
66  EntitiesFieldData::EntData &col_data) {
68  MOFEM_LOG("SYNC", Sev::inform)
69  << "Hello Operator OpRowCol:"
70  << " row field name " << rowFieldName << " row side " << row_side
71  << " row type " << type_name[row_type] << " nb dofs on row entity"
72  << row_data.getIndices().size() << " : "
73  << " col field name " << colFieldName << " col side " << col_side
74  << " col type " << type_name[col_type] << " nb dofs on col entity"
75  << col_data.getIndices().size();
77  }
78 };
79 
81  OpVolume(const std::string &field_name)
83  field_name, OPROW) {
84  }
85  MoFEMErrorCode doWork(int side, EntityType type,
88  if (type == MBVERTEX) {
89  MOFEM_LOG("SYNC", Sev::inform) << "Hello Operator OpVolume:"
90  << " volume " << getVolume();
91  }
93  }
94 };
95 
97  OpFace(const std::string &field_name)
98  : FaceElementForcesAndSourcesCore::UserDataOperator(field_name, OPROW) {}
99  MoFEMErrorCode doWork(int side, EntityType type,
102  if (type == MBVERTEX) {
103  MOFEM_LOG("SYNC", Sev::inform) << "Hello Operator OpFace:"
104  << " normal " << getNormal();
105  }
107  }
108 };
109 
111  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> &feSidePtr;
113  const std::string &field_name,
114  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> &fe_side_ptr)
116  feSidePtr(fe_side_ptr) {}
117  MoFEMErrorCode doWork(int side, EntityType type,
119 
121  if (type == MBVERTEX) {
122  MOFEM_LOG("SYNC", Sev::inform) << "Hello Operator OpSideFace";
123  CHKERR loopSideVolumes("dFE", *feSidePtr);
124  }
126  }
127 };
128 
131  OpVolumeSide(const std::string &field_name)
133  field_name, field_name, OPROW) {}
134  MoFEMErrorCode doWork(int side, EntityType type,
137  if (type == MBVERTEX) {
138  MOFEM_LOG("SYNC", Sev::inform)
139  << "Hello Operator OpVolumeSide:"
140  << " volume " << getVolume() << " normal " << getNormal();
141  }
143  }
144 };
145 
146 int main(int argc, char *argv[]) {
147 
148  // initialize petsc
149  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
150 
151  try {
152 
153  type_name[MBVERTEX] = "Vertex";
154  type_name[MBEDGE] = "Edge";
155  type_name[MBTRI] = "Triangle";
156  type_name[MBTET] = "Tetrahedra";
157 
158  // Register DM Manager
159  DMType dm_name = "DMMOFEM";
160  CHKERR DMRegister_MoFEM(dm_name);
161 
162  // Create MoAB database
163  moab::Core moab_core;
164  moab::Interface &moab = moab_core;
165 
166  // Create MoFEM database and link it to MoAB
167  MoFEM::Core mofem_core(moab);
168  MoFEM::Interface &m_field = mofem_core;
169 
170  // Simple interface
171  Simple *simple;
172  CHKERR m_field.getInterface(simple);
173 
174  // get options from command line
175  CHKERR simple->getOptions();
176  // load mesh file
177  CHKERR simple->loadFile();
178  // add fields
179  CHKERR simple->addDomainField("U", H1, AINSWORTH_LEGENDRE_BASE, 3);
180  CHKERR simple->addBoundaryField("L", H1, AINSWORTH_LEGENDRE_BASE, 3);
181  CHKERR simple->addSkeletonField("S", H1, AINSWORTH_LEGENDRE_BASE, 3);
182  // set fields order
183  CHKERR simple->setFieldOrder("U", 4);
184  CHKERR simple->setFieldOrder("L", 3);
185  CHKERR simple->setFieldOrder("S", 3);
186 
187  // setup problem
188  CHKERR simple->setUp();
189 
190  auto pipeline_mng = m_field.getInterface<PipelineManager>();
191 
192  //! [set operator to the volume element instance]
193  pipeline_mng->getOpDomainRhsPipeline().push_back(new OpRow("U"));
194  pipeline_mng->getOpDomainRhsPipeline().push_back(new OpVolume("U"));
195  pipeline_mng->getOpDomainLhsPipeline().push_back(
196  new OpRowCol("U", "U", true));
197  //! [set operator to the volume element instance]
198 
199  //! [set operator to the face element instance]
200  pipeline_mng->getOpBoundaryRhsPipeline().push_back(new OpRow("L"));
201  pipeline_mng->getOpBoundaryRhsPipeline().push_back(new OpFace("L"));
202  pipeline_mng->getOpBoundaryLhsPipeline().push_back(
203  new OpRowCol("U", "L", false));
204  //! [set operator to the face element instance]
205 
206  //! [create skeleton side element and push operator to it]
207  boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe(
209  side_fe->getOpPtrVector().push_back(new OpVolumeSide("U"));
210  //! [create skeleton side element and push operator to it]
211 
212  //! [set operator to the face element on skeleton instance]
213  pipeline_mng->getOpSkeletonRhsPipeline().push_back(new OpRow("S"));
214  pipeline_mng->getOpSkeletonRhsPipeline().push_back(
215  new OpFaceSide("S", side_fe));
216  //! [set operator to the face element on skeleton instance]
217 
218  //! [executing finite elements]
219  CHKERR pipeline_mng->loopFiniteElements();
220  //! [executing finite elements]
221 
222  MOFEM_LOG_SYNCHRONISE(m_field.get_comm());
223  }
224  CATCH_ERRORS;
225 
226  // finish work cleaning memory, getting statistics, etc.
228 
229  return 0;
230 }
ForcesAndSourcesCore::UserDataOperator UserDataOperator
#define MOFEM_LOG_SYNCHRONISE(comm)
Synchronise "SYNC" channel.
Definition: LogManager.hpp:348
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:385
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:73
@ H1
continuous field
Definition: definitions.h:98
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:59
boost::ptr_vector< UserDataOperator > & getOpDomainRhsPipeline()
Get the Op Domain Rhs Pipeline object.
VolumeElementForcesAndSourcesCoreOnSideSwitch< 0 > VolumeElementForcesAndSourcesCoreOnSide
Volume element used to integrate on skeleton.
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:311
int main(int argc, char *argv[])
#define HelloFunctionBegin
Definition: hello_world.cpp:33
static char help[]
Definition: hello_world.cpp:29
static map< EntityType, std::string > type_name
Definition: hello_world.cpp:31
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
CoreTmp< 0 > Core
Definition: Core.hpp:1096
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1965
virtual MPI_Comm & get_comm() 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.
Data on single entity (This is passed as argument to DataOperator::doWork)
const VectorInt & getIndices() const
Get global indices of dofs on entity.
structure to get information form mofem into EntitiesFieldData
PipelineManager interface.
Simple interface for fast problem set-up.
Definition: Simple.hpp:33
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
OpFace(const std::string &field_name)
Definition: hello_world.cpp:97
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:99
boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > & feSidePtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
OpFaceSide(const std::string &field_name, boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > &fe_side_ptr)
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:63
OpRowCol(const std::string row_field, const std::string col_field, const bool symm)
Definition: hello_world.cpp:59
OpRow(const std::string &field_name)
Definition: hello_world.cpp:40
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Definition: hello_world.cpp:42
OpVolume(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:85
OpVolumeSide(const std::string &field_name)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.