v0.15.0
Loading...
Searching...
No Matches
elastic.cpp
Go to the documentation of this file.
1/**
2 * @file elastic.cpp
3 * @example mofem/tutorials/vec-0/elastic.cpp
4 * @brief elastic example
5 * @version 0.15.0
6 * @date 2025-10-07
7 *
8 * @copyright Anonymous authors (c) 2025 under the MIT license (see LICENSE for
9 * details)
10 */
11
12#include <MoFEM.hpp>
13
14using namespace MoFEM;
15
16constexpr int BASE_DIM = 1; //< Dimension of the base functions
17
18//! [Define dimension]
19constexpr int SPACE_DIM =
20 EXECUTABLE_DIMENSION; //< Space dimension of problem, mesh
21//! [Define dimension]
23 ? AssemblyType::BLOCK_SCHUR
24 : AssemblyType::PETSC; //< selected assembly type
25
26constexpr IntegrationType I =
27 IntegrationType::GAUSS; //< selected integration type
28
29//! [Define entities]
34using DomainEleOp = DomainEle::UserDataOperator;
35using BoundaryEleOp = BoundaryEle::UserDataOperator;
36//! [Define entities]
37
38struct DomainBCs {};
39struct BoundaryBCs {};
40
47
48template <int DIM> struct PostProcEleByDim;
49
50template <> struct PostProcEleByDim<2> {
54};
55
56template <> struct PostProcEleByDim<3> {
60};
61
65
66#include <ElasticExample.hpp>
67
68static char help[] = "...\n\n";
69
70int main(int argc, char *argv[]) {
71
72 // Initialisation of MoFEM/PETSc and MOAB data structures
73 const char param_file[] = "param_file.petsc";
74 MoFEM::Core::Initialize(&argc, &argv, param_file, help);
75
76 auto core_log = logging::core::get();
77 core_log->add_sink(
79
80 core_log->add_sink(
82 LogManager::setLog("FieldEvaluator");
83 MOFEM_LOG_TAG("FieldEvaluator", "field_eval");
84
85 try {
86
87 //! [Register MoFEM discrete manager in PETSc]
88 DMType dm_name = "DMMOFEM";
89 CHKERR DMRegister_MoFEM(dm_name);
90 DMType dm_name_mg = "DMMOFEM_MG";
92 //! [Register MoFEM discrete manager in PETSc
93
94 //! [Create MoAB]
95 moab::Core mb_instance; ///< mesh database
96 moab::Interface &moab = mb_instance; ///< mesh database interface
97 //! [Create MoAB]
98
99 //! [Create MoFEM]
100 MoFEM::Core core(moab); ///< finite element database
101 MoFEM::Interface &m_field = core; ///< finite element database interface
102 //! [Create MoFEM]
103
104 //! [Example]
105 ElasticExample ex(m_field);
106 CHKERR ex.runProblem();
107 //! [Example]
108 }
110
112}
Implementation of elastic example class.
int main()
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
ElementsAndOps< SPACE_DIM >::BoundaryEle BoundaryEle
#define SCHUR_ASSEMBLE
Definition contact.cpp:18
#define CATCH_ERRORS
Catch errors.
#define CHKERR
Inline error check.
PostProcEleByDim< SPACE_DIM >::PostProcEleDomain PostProcEleDomain
PostProcEleByDim< SPACE_DIM >::PostProcEleBdy PostProcEleBdy
static char help[]
Definition elastic.cpp:68
constexpr int SPACE_DIM
[Define dimension]
Definition elastic.cpp:19
constexpr int BASE_DIM
Definition elastic.cpp:16
constexpr IntegrationType I
Definition elastic.cpp:26
constexpr AssemblyType A
[Define dimension]
Definition elastic.cpp:22
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition DMMoFEM.cpp:43
MoFEMErrorCode DMRegister_MGViaApproxOrders(const char sname[])
Register DM for Multi-Grid via approximation orders.
IntegrationType
Form integrator integration types.
AssemblyType
[Storage and set boundary conditions]
@ PETSC
Standard PETSc assembly.
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
#define EXECUTABLE_DIMENSION
Definition plastic.cpp:13
ElementsAndOps< SPACE_DIM >::SideEle SideEle
Definition plastic.cpp:61
Boundary conditions marker.
Definition elastic.cpp:39
[Define entities]
Definition elastic.cpp:38
MoFEMErrorCode runProblem()
[Run problem]
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:118
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
static boost::shared_ptr< std::ostream > getStrmSync()
Get the strm sync object.
Assembly methods.
Definition Natural.hpp:65
Template struct for dimension-specific finite element types.
PipelineManager::ElementsAndOpsByDim< 2 >::FaceSideEle SideEle
PipelineManager::ElementsAndOpsByDim< 3 >::FaceSideEle SideEle