v0.13.1
elastic.cpp
Go to the documentation of this file.
1/**
2 * @file elastic.cpp
3 * @brief elastic example
4 * @version 0.13.2
5 * @date 2022-09-19
6 *
7 * @copyright Copyright (c) 2022
8 *
9 */
10
11#include <MoFEM.hpp>
12
13using namespace MoFEM;
14
15template <int DIM> struct ElementsAndOps {};
16
17template <> struct ElementsAndOps<2> {
20};
21
22template <> struct ElementsAndOps<3> {
25};
26
27//! [Define dimension]
28constexpr int SPACE_DIM =
29 EXECUTABLE_DIMENSION; //< Space dimension of problem, mesh
30constexpr AssemblyType A = AssemblyType::PETSC; //< selected assembly type
31constexpr IntegrationType I =
32 IntegrationType::GAUSS; //< selected integration type
33
40
42 I>::OpGradSymTensorGrad<1, SPACE_DIM, SPACE_DIM, 0>;
44 I>::OpGradTimesSymTensor<1, SPACE_DIM, SPACE_DIM>;
45
46struct DomainBCs {};
47struct BoundaryBCs {};
48
55
56#include <ElasticSpring.hpp>
57#include <NaturalDomainBC.hpp>
58#include <NaturalBoundaryBC.hpp>
59
60constexpr double young_modulus = 1;
61constexpr double poisson_ratio = 0.3;
62constexpr double bulk_modulus_K = young_modulus / (3 * (1 - 2 * poisson_ratio));
63constexpr double shear_modulus_G = young_modulus / (2 * (1 + poisson_ratio));
64
65struct Example {
66
67 Example(MoFEM::Interface &m_field) : mField(m_field) {}
68
70
71private:
73
81
83 boost::ptr_vector<ForcesAndSourcesCore::UserDataOperator> &pipeline,
84 std::string field_name, std::string block_name,
85 boost::shared_ptr<MatrixDouble> mat_D_Ptr, Sev sev);
86};
87
89 boost::ptr_vector<ForcesAndSourcesCore::UserDataOperator> &pipeline,
90 std::string field_name, std::string block_name,
91 boost::shared_ptr<MatrixDouble> mat_D_Ptr, Sev sev) {
93
94 struct OpMatBlock : public DomainEleOp {
95 OpMatBlock(std::string field_name, boost::shared_ptr<MatrixDouble> m,
96 Range &&ents, double bulk_modulus_K, double shear_modulus_G)
98 bulkModulusK(bulk_modulus_K), shearModulusG(shear_modulus_G) {
99 std::fill(&(doEntities[MBEDGE]), &(doEntities[MBMAXTYPE]), false);
100 }
101
102 MoFEMErrorCode doWork(int side, EntityType type,
105 if (!feEnts.empty()) {
106 if (feEnts.find(getFEEntityHandle()) != feEnts.end()) {
107 CHKERR getMatDPtr(matDPtr, bulkModulusK, shearModulusG);
108 }
109 } else {
110 CHKERR getMatDPtr(matDPtr, bulkModulusK, shearModulusG);
111 }
113 }
114
115 private:
116 boost::shared_ptr<MatrixDouble> matDPtr;
117
118 Range feEnts;
119 double bulkModulusK;
120 double shearModulusG;
121
122 MoFEMErrorCode getMatDPtr(boost::shared_ptr<MatrixDouble> mat_D_ptr,
123 double bulk_modulus_K, double shear_modulus_G) {
125 //! [Calculate elasticity tensor]
126 auto set_material_stiffness = [&]() {
132 double A = (SPACE_DIM == 2)
133 ? 2 * shear_modulus_G /
134 (bulk_modulus_K + (4. / 3.) * shear_modulus_G)
135 : 1;
136 auto t_D = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, 0>(*mat_D_ptr);
137 t_D(i, j, k, l) =
138 2 * shear_modulus_G * ((t_kd(i, k) ^ t_kd(j, l)) / 4.) +
139 A * (bulk_modulus_K - (2. / 3.) * shear_modulus_G) * t_kd(i, j) *
140 t_kd(k, l);
141 };
142 //! [Calculate elasticity tensor]
143 constexpr auto size_symm = (SPACE_DIM * (SPACE_DIM + 1)) / 2;
144 mat_D_ptr->resize(size_symm * size_symm, 1);
145 set_material_stiffness();
147 }
148 };
149
150 pipeline.push_back(new OpMatBlock(field_name, mat_D_Ptr, Range(),
152
153 auto add_op = [&](auto &&meshset_vec_ptr) {
155 for (auto m : meshset_vec_ptr) {
156 MOFEM_TAG_AND_LOG("WORLD", sev, "MatBlock") << *m;
157 std::vector<double> block_data;
158 CHKERR m->getAttributes(block_data);
159 if (block_data.size() != 2) {
160 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
161 "Expected that block has two attribute");
162 }
163 auto get_block_ents = [&]() {
164 Range ents;
165 CHKERR
166 mField.get_moab().get_entities_by_handle(m->meshset, ents, true);
167 return ents;
168 };
169
170 double young_modulus = block_data[0];
171 double poisson_ratio = block_data[1];
172 double bulk_modulus_K = young_modulus / (3 * (1 - 2 * poisson_ratio));
173 double shear_modulus_G = young_modulus / (2 * (1 + poisson_ratio));
174
175 MOFEM_TAG_AND_LOG("WORLD", sev, "MatBlock")
176 << "E = " << young_modulus << " nu = " << poisson_ratio;
177
178 pipeline.push_back(new OpMatBlock(field_name, mat_D_Ptr, get_block_ents(),
180 }
181 MOFEM_LOG_CHANNEL("WORLD");
183 };
184
185 CHKERR add_op(
186
188
189 (boost::format("%s(.*)") % block_name).str()
190
191 ))
192
193 );
194
196}
197
198//! [Run problem]
209}
210//! [Run problem]
211
212//! [Read mesh]
216 CHKERR simple->getOptions();
217 CHKERR simple->loadFile();
219}
220//! [Read mesh]
221
222//! [Set up problem]
226 // Add field
227 CHKERR simple->addDomainField("U", H1, AINSWORTH_LEGENDRE_BASE, SPACE_DIM);
228 CHKERR simple->addBoundaryField("U", H1, AINSWORTH_LEGENDRE_BASE, SPACE_DIM);
229 int order = 3;
230 CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-order", &order, PETSC_NULL);
231
232 CHKERR simple->addDataField("GEOMETRY", H1, AINSWORTH_LEGENDRE_BASE,
233 SPACE_DIM);
234
235 CHKERR simple->setFieldOrder("U", order);
236 CHKERR simple->setFieldOrder("GEOMETRY", 2);
237 CHKERR simple->setUp();
238
239 auto project_ho_geometry = [&]() {
240 Projection10NodeCoordsOnField ent_method(mField, "GEOMETRY");
241 return mField.loop_dofs("GEOMETRY", ent_method);
242 };
243 CHKERR project_ho_geometry();
244
246}
247//! [Set up problem]
248
249//! [Boundary condition]
252 auto *pipeline_mng = mField.getInterface<PipelineManager>();
254 auto bc_mng = mField.getInterface<BcManager>();
255
257 pipeline_mng->getOpDomainRhsPipeline(), {H1}, "GEOMETRY");
259 pipeline_mng->getOpBoundaryRhsPipeline(), {NOSPACE}, "GEOMETRY");
260
261 // Infernal forces
262 auto mat_grad_ptr = boost::make_shared<MatrixDouble>();
263 auto mat_strain_ptr = boost::make_shared<MatrixDouble>();
264 auto mat_stress_ptr = boost::make_shared<MatrixDouble>();
265 pipeline_mng->getOpDomainRhsPipeline().push_back(
267 mat_grad_ptr));
268 pipeline_mng->getOpDomainRhsPipeline().push_back(
269 new OpSymmetrizeTensor<SPACE_DIM>("U", mat_grad_ptr, mat_strain_ptr));
270
271 auto mat_D_ptr = boost::make_shared<MatrixDouble>();
272 CHKERR addMatBlockOps(pipeline_mng->getOpDomainRhsPipeline(), "U",
273 "MAT_ELASTIC", mat_D_ptr, Sev::inform);
274 pipeline_mng->getOpDomainRhsPipeline().push_back(
276 "U", mat_strain_ptr, mat_stress_ptr, mat_D_ptr));
277 pipeline_mng->getOpDomainRhsPipeline().push_back(new OpInternalForce(
278 "U", mat_stress_ptr,
279 [](double, double, double) constexpr { return -1; }));
280
281 // Body forces
283 pipeline_mng->getOpDomainRhsPipeline(), mField, "U", Sev::inform);
284 // Add force boundary condition
286 pipeline_mng->getOpBoundaryRhsPipeline(), mField, "U", -1, Sev::inform);
288 pipeline_mng->getOpBoundaryLhsPipeline(), mField, "U", Sev::verbose);
289
290 // Essential boundary condition
291 CHKERR bc_mng->removeBlockDOFsOnEntities<DisplacementCubitBcData>(
292 simple->getProblemName(), "U");
293
294 auto get_pre_proc_hook = [&]() {
296 mField, pipeline_mng->getDomainRhsFE(), {});
297 };
298
299 pipeline_mng->getDomainRhsFE()->preProcessHook = get_pre_proc_hook();
300 pipeline_mng->getBoundaryRhsFE()->preProcessHook = get_pre_proc_hook();
301
303}
304//! [Boundary condition]
305
306//! [Push operators to pipeline]
310
312 pipeline_mng->getOpDomainLhsPipeline(), {H1}, "GEOMETRY");
313
314 auto mat_D_ptr = boost::make_shared<MatrixDouble>();
315 CHKERR addMatBlockOps(pipeline_mng->getOpDomainLhsPipeline(), "U",
316 "MAT_ELASTIC", mat_D_ptr, Sev::verbose);
317 pipeline_mng->getOpDomainLhsPipeline().push_back(
318 new OpK("U", "U", mat_D_ptr));
319
320 auto integration_rule = [](int, int, int approx_order) {
321 return 2 * approx_order;
322 };
327
329}
330//! [Push operators to pipeline]
331
332//! [Solve]
335 auto *simple = mField.getInterface<Simple>();
336 auto *pipeline_mng = mField.getInterface<PipelineManager>();
337 auto solver = pipeline_mng->createKSP();
338 CHKERR KSPSetFromOptions(solver);
339 CHKERR KSPSetUp(solver);
340
341 auto dm = simple->getDM();
342 auto D = smartCreateDMVector(dm);
343 auto F = smartVectorDuplicate(D);
344
345 CHKERR KSPSolve(solver, F, D);
346 CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
347 CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
348 CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
350}
351//! [Solve]
352
353//! [Postprocess results]
357 auto det_ptr = boost::make_shared<VectorDouble>();
358 auto jac_ptr = boost::make_shared<MatrixDouble>();
359 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
360 pipeline_mng->getDomainRhsFE().reset();
361 pipeline_mng->getDomainLhsFE().reset();
362 pipeline_mng->getBoundaryRhsFE().reset();
363 pipeline_mng->getBoundaryLhsFE().reset();
364
365 auto post_proc_fe = boost::make_shared<PostProcEle>(mField);
366
368 post_proc_fe->getOpPtrVector(), {H1}, "GEOMETRY");
369
370 auto mat_grad_ptr = boost::make_shared<MatrixDouble>();
371 auto mat_strain_ptr = boost::make_shared<MatrixDouble>();
372 auto mat_stress_ptr = boost::make_shared<MatrixDouble>();
373
374 post_proc_fe->getOpPtrVector().push_back(
376 mat_grad_ptr));
377 post_proc_fe->getOpPtrVector().push_back(
378 new OpSymmetrizeTensor<SPACE_DIM>("U", mat_grad_ptr, mat_strain_ptr));
379
380 auto mat_D_ptr = boost::make_shared<MatrixDouble>();
381 CHKERR addMatBlockOps(post_proc_fe->getOpPtrVector(), "U", "MAT_ELASTIC",
382 mat_D_ptr, Sev::verbose);
383 post_proc_fe->getOpPtrVector().push_back(
385 "U", mat_strain_ptr, mat_stress_ptr, mat_D_ptr));
386
387 auto u_ptr = boost::make_shared<MatrixDouble>();
388 post_proc_fe->getOpPtrVector().push_back(
390 auto x_ptr = boost::make_shared<MatrixDouble>();
391 post_proc_fe->getOpPtrVector().push_back(
392 new OpCalculateVectorFieldValues<SPACE_DIM>("GEOMETRY", x_ptr));
393
395
396 post_proc_fe->getOpPtrVector().push_back(
397
398 new OpPPMap(
399
400 post_proc_fe->getPostProcMesh(), post_proc_fe->getMapGaussPts(),
401
402 {},
403
404 {{"U", u_ptr}, {"GEOMETRY", x_ptr}},
405
406 {},
407
408 {{"STRAIN", mat_strain_ptr}, {"STRESS", mat_stress_ptr}}
409
410 )
411
412 );
413
414 pipeline_mng->getBoundaryRhsFE().reset();
415 pipeline_mng->getDomainRhsFE() = post_proc_fe;
416 CHKERR pipeline_mng->loopFiniteElements();
417 CHKERR post_proc_fe->writeFile("out_elastic.h5m");
419}
420//! [Postprocessing results]
421
422//! [Check]
424 MOFEM_LOG_CHANNEL("WORLD");
428 pipeline_mng->getDomainRhsFE().reset();
429 pipeline_mng->getDomainLhsFE().reset();
430 pipeline_mng->getBoundaryRhsFE().reset();
431 pipeline_mng->getBoundaryLhsFE().reset();
432
433 auto integration_rule = [](int, int, int p_data) { return 2 * p_data; };
434 CHKERR pipeline_mng->setDomainRhsIntegrationRule(integration_rule);
435 CHKERR pipeline_mng->setBoundaryRhsIntegrationRule(integration_rule);
436
438 pipeline_mng->getOpDomainRhsPipeline(), {H1}, "GEOMETRY");
440 pipeline_mng->getOpBoundaryRhsPipeline(), {NOSPACE}, "GEOMETRY");
441
442 auto mat_grad_ptr = boost::make_shared<MatrixDouble>();
443 auto mat_strain_ptr = boost::make_shared<MatrixDouble>();
444 auto mat_stress_ptr = boost::make_shared<MatrixDouble>();
445
446 pipeline_mng->getOpDomainRhsPipeline().push_back(
448 mat_grad_ptr));
449 pipeline_mng->getOpDomainRhsPipeline().push_back(
450 new OpSymmetrizeTensor<SPACE_DIM>("U", mat_grad_ptr, mat_strain_ptr));
451
452 auto mat_D_ptr = boost::make_shared<MatrixDouble>();
453 CHKERR addMatBlockOps(pipeline_mng->getOpDomainRhsPipeline(), "U",
454 "MAT_ELASTIC", mat_D_ptr, Sev::verbose);
455 pipeline_mng->getOpDomainRhsPipeline().push_back(
457 "U", mat_strain_ptr, mat_stress_ptr, mat_D_ptr));
458 pipeline_mng->getOpDomainRhsPipeline().push_back(
459 new OpInternalForce("U", mat_stress_ptr));
460
462 pipeline_mng->getOpDomainRhsPipeline(), mField, "U", Sev::verbose);
464 pipeline_mng->getOpBoundaryRhsPipeline(), mField, "U", 1, Sev::verbose);
465
466 auto dm = simple->getDM();
467 auto res = smartCreateDMVector(dm);
468 pipeline_mng->getDomainRhsFE()->ksp_f = res;
469 pipeline_mng->getBoundaryRhsFE()->ksp_f = res;
470
471 CHKERR VecZeroEntries(res);
472
473 CHKERR mField.getInterface<FieldBlas>()->fieldScale(-1, "U");
474 CHKERR pipeline_mng->loopFiniteElements();
475 CHKERR mField.getInterface<FieldBlas>()->fieldScale(-1, "U");
476
477 CHKERR VecGhostUpdateBegin(res, ADD_VALUES, SCATTER_REVERSE);
478 CHKERR VecGhostUpdateEnd(res, ADD_VALUES, SCATTER_REVERSE);
479 CHKERR VecAssemblyBegin(res);
480 CHKERR VecAssemblyEnd(res);
481
482 double nrm2;
483 CHKERR VecNorm(res, NORM_2, &nrm2);
484 MOFEM_LOG_CHANNEL("WORLD");
485 MOFEM_LOG_C("WORLD", Sev::inform, "residual = %3.4e\n", nrm2);
486 constexpr double eps = 1e-8;
487 if (nrm2 > eps)
488 SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "Residual is not zero");
489
491}
492//! [Check]
493
494static char help[] = "...\n\n";
495
496int main(int argc, char *argv[]) {
497
498 // Initialisation of MoFEM/PETSc and MOAB data structures
499 const char param_file[] = "param_file.petsc";
501
502 try {
503
504 //! [Register MoFEM discrete manager in PETSc]
505 DMType dm_name = "DMMOFEM";
506 CHKERR DMRegister_MoFEM(dm_name);
507 //! [Register MoFEM discrete manager in PETSc
508
509 //! [Create MoAB]
510 moab::Core mb_instance; ///< mesh database
511 moab::Interface &moab = mb_instance; ///< mesh database interface
512 //! [Create MoAB]
513
514 //! [Create MoFEM]
515 MoFEM::Core core(moab); ///< finite element database
516 MoFEM::Interface &m_field = core; ///< finite element database interface
517 //! [Create MoFEM]
518
519 //! [Example]
520 Example ex(m_field);
521 CHKERR ex.runProblem();
522 //! [Example]
523 }
525
527}
std::string param_file
Implementation of elastic spring bc.
ForcesAndSourcesCore::UserDataOperator UserDataOperator
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:304
#define MOFEM_TAG_AND_LOG(channel, severity, tag)
Tag and log in channel.
Definition: LogManager.hpp:345
Implementation of natural boundary conditions.
Boundary conditions in domain, i.e. body forces.
void simple(double P1[], double P2[], double P3[], double c[], const int N)
Definition: acoustic.cpp:69
MoFEM::FaceElementForcesAndSourcesCore FaceEle
static const double eps
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
ElementsAndOps< SPACE_DIM >::BoundaryEle BoundaryEle
Kronecker Delta class symmetric.
MoFEM::EdgeElementForcesAndSourcesCore EdgeEle
#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 MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
#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
int main(int argc, char *argv[])
Definition: elastic.cpp:496
EntitiesFieldData::EntData EntData
Definition: elastic.cpp:34
static char help[]
[Check]
Definition: elastic.cpp:494
FormsIntegrators< DomainEleOp >::Assembly< A >::LinearForm< I >::OpGradTimesSymTensor< 1, SPACE_DIM, SPACE_DIM > OpInternalForce
Definition: elastic.cpp:44
constexpr int SPACE_DIM
[Define dimension]
Definition: elastic.cpp:28
constexpr double poisson_ratio
Definition: elastic.cpp:61
constexpr double shear_modulus_G
Definition: elastic.cpp:63
constexpr IntegrationType I
Definition: elastic.cpp:31
DomainEle::UserDataOperator DomainEleOp
Definition: elastic.cpp:36
constexpr double bulk_modulus_K
Definition: elastic.cpp:62
constexpr AssemblyType A
Definition: elastic.cpp:30
FormsIntegrators< DomainEleOp >::Assembly< A >::BiLinearForm< I >::OpGradSymTensorGrad< 1, SPACE_DIM, SPACE_DIM, 0 > OpK
Definition: elastic.cpp:42
constexpr double young_modulus
Definition: elastic.cpp:60
FTensor::Index< 'm', SPACE_DIM > m
auto integration_rule
constexpr auto t_kd
PetscErrorCode DMoFEMMeshToLocalVector(DM dm, Vec l, InsertMode mode, ScatterMode scatter_mode)
set local (or ghosted) vector values on mesh for partition only
Definition: DMMMoFEM.cpp:470
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMMoFEM.cpp:47
auto smartCreateDMVector(DM dm)
Get smart vector from DM.
Definition: DMMoFEM.hpp:965
MoFEMErrorCode loopFiniteElements(SmartPetscObj< DM > dm=nullptr)
Iterate finite elements.
boost::ptr_vector< UserDataOperator > & getOpDomainLhsPipeline()
Get the Op Domain Lhs Pipeline object.
IntegrationType
Form integrator integration types.
AssemblyType
[Storage and set boundary conditions]
SeverityLevel
Severity levels.
Definition: LogManager.hpp:33
#define MOFEM_LOG_CHANNEL(channel)
Set and reset channel.
Definition: LogManager.hpp:277
virtual MoFEMErrorCode loop_dofs(const Problem *problem_ptr, const std::string &field_name, RowColData rc, DofMethod &method, int lower_rank, int upper_rank, int verb=DEFAULT_VERBOSITY)=0
Make a loop over dofs.
MoFEMErrorCode getCubitMeshsetPtr(const int ms_id, const CubitBCType cubit_bc_type, const CubitMeshSets **cubit_meshset_ptr) const
get cubit meshset
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: MoFEM.hpp:24
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
SmartPetscObj< Vec > smartVectorDuplicate(SmartPetscObj< Vec > &vec)
Create duplicate vector of smart vector.
CoreTmp< 0 > Core
Definition: Core.hpp:1086
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1955
const double D
diffusivity
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
#define EXECUTABLE_DIMENSION
Definition: plastic.cpp:10
constexpr auto field_name
static constexpr int approx_order
[Example]
Definition: plastic.cpp:120
MoFEMErrorCode boundaryCondition()
MoFEMErrorCode assembleSystem()
MoFEMErrorCode readMesh()
MoFEMErrorCode addMatBlockOps(boost::ptr_vector< ForcesAndSourcesCore::UserDataOperator > &pipeline, std::string field_name, std::string block_name, boost::shared_ptr< MatrixDouble > mat_D_Ptr, Sev sev)
Definition: elastic.cpp:88
MoFEMErrorCode checkResults()
MoFEMErrorCode solveSystem()
Example(MoFEM::Interface &m_field)
Definition: elastic.cpp:67
MoFEMErrorCode runProblem()
MoFEM::Interface & mField
Definition: plastic.cpp:127
MoFEMErrorCode setupProblem()
MoFEMErrorCode outputResults()
boost::shared_ptr< MatrixDouble > matDPtr
Add operators pushing bases from local to physical configuration.
Simple interface for fast problem set-up.
Definition: BcManager.hpp:23
virtual moab::Interface & get_moab()=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.
Definition of the displacement bc data structure.
Definition: BCData.hpp:72
Data on single entity (This is passed as argument to DataOperator::doWork)
Basic algebra on fields.
Definition: FieldBlas.hpp:21
@ OPROW
operator doWork function is executed on FE rows
Interface for managing meshsets containing materials and boundary conditions.
Assembly methods.
Definition: Natural.hpp:57
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
Get values at integration pts for tensor filed rank 1, i.e. vector field.
Post post-proc data at points from hash maps.
PipelineManager interface.
boost::shared_ptr< FEMethod > & getDomainRhsFE()
boost::shared_ptr< FEMethod > & getDomainLhsFE()
boost::shared_ptr< FEMethod > & getBoundaryLhsFE()
MoFEMErrorCode setDomainRhsIntegrationRule(RuleHookFun rule)
MoFEMErrorCode setBoundaryLhsIntegrationRule(RuleHookFun rule)
MoFEMErrorCode setBoundaryRhsIntegrationRule(RuleHookFun rule)
boost::shared_ptr< FEMethod > & getBoundaryRhsFE()
MoFEMErrorCode setDomainLhsIntegrationRule(RuleHookFun rule)
Projection of edge entities with one mid-node on hierarchical basis.
Simple interface for fast problem set-up.
Definition: Simple.hpp:26
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.