v0.14.0
hello_world.cpp

Prints basic information about users data operator. See more details in FUN-0: Hello world

/**
* \file hello_world.cpp
* \ingroup mofem_simple_interface
* \example hello_world.cpp
*
* Prints basic information about users data operator.
* See more details in \ref hello_world_tut1
*
*/
#include <MoFEM.hpp>
using namespace MoFEM;
static char help[] = "...\n\n";
#define HelloFunctionBegin \
MoFEMFunctionBegin; \
MOFEM_LOG_CHANNEL("SYNC"); \
MOFEM_LOG_FUNCTION(); \
MOFEM_LOG_TAG("WORLD", "HelloWorld");
OpRow(const std::string &field_name)
MoFEMErrorCode doWork(int side, EntityType type,
if (type == MBVERTEX) {
// get number of evaluated element in the loop
MOFEM_LOG("SYNC", Sev::inform) << "**** " << getNinTheLoop() << " ****";
MOFEM_LOG("SYNC", Sev::inform) << "**** Operators ****";
}
MOFEM_LOG("SYNC", Sev::inform)
<< "Hello Operator OpRow:"
<< " field name " << rowFieldName << " side " << side << " type "
<< CN::EntityTypeName(type) << " nb dofs on entity "
<< data.getIndices().size();
}
};
OpRowCol(const std::string row_field, const std::string col_field,
const bool symm)
: ForcesAndSourcesCore::UserDataOperator(row_field, col_field, OPROWCOL,
symm) {}
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
EntityType col_type,
MOFEM_LOG("SYNC", Sev::inform)
<< "Hello Operator OpRowCol:"
<< " row field name " << rowFieldName << " row side " << row_side
<< " row type " << CN::EntityTypeName(row_type)
<< " nb dofs on row entity " << row_data.getIndices().size() << " : "
<< " col field name " << colFieldName << " col side " << col_side
<< " col type " << CN::EntityTypeName(col_type)
<< " nb dofs on col entity " << col_data.getIndices().size();
}
};
OpVolume(const std::string &field_name)
field_name, OPROW) {
}
MoFEMErrorCode doWork(int side, EntityType type,
if (type == MBVERTEX) {
MOFEM_LOG("SYNC", Sev::inform) << "Hello Operator OpVolume:"
<< " volume " << getVolume();
}
}
};
OpFace(const std::string &field_name)
MoFEMErrorCode doWork(int side, EntityType type,
if (type == MBVERTEX) {
MOFEM_LOG("SYNC", Sev::inform) << "Hello Operator OpFace:"
<< " normal " << getNormal();
}
}
};
boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> &feSidePtr;
const std::string &field_name,
boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> &fe_side_ptr)
feSidePtr(fe_side_ptr) {}
MoFEMErrorCode doWork(int side, EntityType type,
if (type == MBVERTEX) {
MOFEM_LOG("SYNC", Sev::inform) << "Hello Operator OpSideFace";
CHKERR loopSideVolumes("dFE", *feSidePtr);
}
}
};
OpVolumeSide(const std::string &field_name)
field_name, field_name, OPROW) {}
MoFEMErrorCode doWork(int side, EntityType type,
if (type == MBVERTEX) {
MOFEM_LOG("SYNC", Sev::inform)
<< "Hello Operator OpVolumeSide:"
<< " volume " << getVolume() << " normal " << getNormal();
}
}
};
int main(int argc, char *argv[]) {
// initialize petsc
MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
try {
// Register DM Manager
DMType dm_name = "DMMOFEM";
// Create MoAB database
moab::Core moab_core;
moab::Interface &moab = moab_core;
// Create MoFEM database and link it to MoAB
MoFEM::Core mofem_core(moab);
MoFEM::Interface &m_field = mofem_core;
// Simple interface
// get options from command line
CHKERR simple->getOptions();
// load mesh file
CHKERR simple->loadFile();
// add fields
CHKERR simple->addDomainField("U", H1, AINSWORTH_LEGENDRE_BASE, 3);
CHKERR simple->addBoundaryField("L", H1, AINSWORTH_LEGENDRE_BASE, 3);
CHKERR simple->addSkeletonField("S", H1, AINSWORTH_LEGENDRE_BASE, 3);
// set fields order
CHKERR simple->setFieldOrder("U", 4);
CHKERR simple->setFieldOrder("L", 3);
CHKERR simple->setFieldOrder("S", 3);
// setup problem
CHKERR simple->setUp();
auto pipeline_mng = m_field.getInterface<PipelineManager>();
//! [set operator to the volume element instance]
pipeline_mng->getOpDomainRhsPipeline().push_back(new OpRow("U"));
pipeline_mng->getOpDomainRhsPipeline().push_back(new OpVolume("U"));
pipeline_mng->getOpDomainLhsPipeline().push_back(
new OpRowCol("U", "U", true));
//! [set operator to the volume element instance]
//! [set operator to the face element instance]
pipeline_mng->getOpBoundaryRhsPipeline().push_back(new OpRow("L"));
pipeline_mng->getOpBoundaryRhsPipeline().push_back(new OpFace("L"));
pipeline_mng->getOpBoundaryLhsPipeline().push_back(
new OpRowCol("U", "L", false));
//! [set operator to the face element instance]
//! [create skeleton side element and push operator to it]
boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe(
side_fe->getOpPtrVector().push_back(new OpVolumeSide("U"));
//! [create skeleton side element and push operator to it]
//! [set operator to the face element on skeleton instance]
pipeline_mng->getOpSkeletonRhsPipeline().push_back(new OpRow("S"));
pipeline_mng->getOpSkeletonRhsPipeline().push_back(
new OpFaceSide("S", side_fe));
//! [set operator to the face element on skeleton instance]
//! [executing finite elements]
CHKERR pipeline_mng->loopFiniteElements();
//! [executing finite elements]
}
// finish work cleaning memory, getting statistics, etc.
return 0;
}
MoFEM::VolumeElementForcesAndSourcesCoreOnSide::UserDataOperator
default operator for TET element
Definition: VolumeElementForcesAndSourcesCoreOnSide.hpp:97
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
MoFEM::EntitiesFieldData::EntData
Data on single entity (This is passed as argument to DataOperator::doWork)
Definition: EntitiesFieldData.hpp:128
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
H1
@ H1
continuous field
Definition: definitions.h:85
OpRowCol
Definition: hello_world.cpp:42
HelloFunctionBegin
#define HelloFunctionBegin
Definition: hello_world.cpp:16
MoFEM::CoreInterface::get_comm
virtual MPI_Comm & get_comm() const =0
help
static char help[]
Definition: hello_world.cpp:14
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::PipelineManager::getOpDomainRhsPipeline
boost::ptr_deque< UserDataOperator > & getOpDomainRhsPipeline()
Get the Op Domain Rhs Pipeline object.
Definition: PipelineManager.hpp:773
MoFEM::PipelineManager
PipelineManager interface.
Definition: PipelineManager.hpp:24
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::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
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
MoFEM::ForcesAndSourcesCore::UserDataOperator
Definition: ForcesAndSourcesCore.hpp:549
MoFEM::FaceElementForcesAndSourcesCore::UserDataOperator
default operator for TRI element
Definition: FaceElementForcesAndSourcesCore.hpp:94
simple
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition: acoustic.cpp:69
OpVolumeSide
Definition: hello_world.cpp:113
OpRow
Definition: hello_world.cpp:22
convert.type
type
Definition: convert.py:64
MOFEM_LOG_SYNCHRONISE
#define MOFEM_LOG_SYNCHRONISE(comm)
Synchronise "SYNC" channel.
Definition: LogManager.hpp:345
MoFEM::EntitiesFieldData::EntData::getIndices
const VectorInt & getIndices() const
Get global indices of dofs on entity.
Definition: EntitiesFieldData.hpp:1214
MoFEM::DMRegister_MoFEM
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMoFEM.cpp:43
OpFace
Definition: boundary_marker.cpp:19
MoFEM::FaceElementForcesAndSourcesCore
Face finite element.
Definition: FaceElementForcesAndSourcesCore.hpp:23
MoFEM::ForcesAndSourcesCore::UserDataOperator
friend class UserDataOperator
Definition: ForcesAndSourcesCore.hpp:482
MoFEM::VolumeElementForcesAndSourcesCore
Volume finite element base.
Definition: VolumeElementForcesAndSourcesCore.hpp:26
MoFEM::VolumeElementForcesAndSourcesCore::UserDataOperator
Definition: VolumeElementForcesAndSourcesCore.hpp:108
field_name
constexpr auto field_name
Definition: poisson_2d_homogeneous.cpp:13
MoFEM::ForcesAndSourcesCore
structure to get information form mofem into EntitiesFieldData
Definition: ForcesAndSourcesCore.hpp:22
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_LOG
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
CATCH_ERRORS
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:385
MoFEM::Core
CoreTmp< 0 > Core
Definition: Core.hpp:1148
main
int main(int argc, char *argv[])
Definition: hello_world.cpp:130
OpFaceSide
Definition: hello_world.cpp:94
AINSWORTH_LEGENDRE_BASE
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
MoFEM::VolumeElementForcesAndSourcesCoreOnSide
Base volume element used to integrate on skeleton.
Definition: VolumeElementForcesAndSourcesCoreOnSide.hpp:22
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
OpVolume
Definition: simple_interface.cpp:18