v0.15.0
Loading...
Searching...
No Matches
nonlinear_elastic.cpp
Go to the documentation of this file.
1/**
2 * @file nonlinear_elastic.cpp
3 * @example nonlinear_elastic.cpp
4 * @date 2025-10-5
5 *
6 * @copyright Anonymous authors (c) 2025 under the MIT license (see LICENSE for
7 * details)
8 */
9
10#include <MoFEM.hpp>
11#include <MatrixFunction.hpp>
12
13using namespace MoFEM;
14
15constexpr int SPACE_DIM =
16 EXECUTABLE_DIMENSION; //< Space dimension of problem, mesh
17
22using DomainEleOp = DomainEle::UserDataOperator;
23using BoundaryEleOp = BoundaryEle::UserDataOperator;
24
29
33
34template <int DIM> struct PostProcEleByDim;
35
36template <> struct PostProcEleByDim<2> {
40};
41
42template <> struct PostProcEleByDim<3> {
46};
47
53
54#include <HenckyOps.hpp>
55using namespace HenckyOps;
56#include <Monitor.hpp>
57
60
61static char help[] = "...\n\n";
62
63int main(int argc, char *argv[]) {
64
65 // Initialisation of MoFEM/PETSc and MOAB data structures
66 const char param_file[] = "param_file.petsc";
67 MoFEM::Core::Initialize(&argc, &argv, param_file, help);
68
69 // Add logging channel for example
70 auto core_log = logging::core::get();
71 core_log->add_sink(
73 LogManager::setLog("EXAMPLE");
74 MOFEM_LOG_TAG("EXAMPLE", "example");
75
76 try {
77
78 //! [Register MoFEM discrete manager in PETSc]
79 DMType dm_name = "DMMOFEM";
80 CHKERR DMRegister_MoFEM(dm_name);
81 //! [Register MoFEM discrete manager in PETSc
82
83 //! [Create MoAB]
84 moab::Core mb_instance; ///< mesh database
85 moab::Interface &moab = mb_instance; ///< mesh database interface
86 //! [Create MoAB]
87
88 //! [Create MoFEM]
89 MoFEM::Core core(moab); ///< finite element database
90 MoFEM::Interface &m_field = core; ///< finite element database interface
91 //! [Create MoFEM]
92
93 //! [Example]
94 ExampleNonlinearElastic ex(m_field);
95 CHKERR ex.runProblem();
96 //! [Example]
97 }
99
101}
int main()
ElementsAndOps< SPACE_DIM >::BoundaryEle BoundaryEle
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
[Define dimension]
#define CATCH_ERRORS
Catch errors.
#define CHKERR
Inline error check.
PostProcEleByDim< SPACE_DIM >::PostProcEleDomain PostProcEleDomain
PostProcEleByDim< SPACE_DIM >::PostProcEleBdy PostProcEleBdy
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition DMMoFEM.cpp:43
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
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.
Assembly methods.
Definition Natural.hpp:65
Force scale operator for reading two columns.
PipelineManager::ElementsAndOpsByDim< 2 >::FaceSideEle SideEle
PipelineManager::ElementsAndOpsByDim< 3 >::FaceSideEle SideEle
static char help[]
constexpr int SPACE_DIM