v0.13.2
Loading...
Searching...
No Matches
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 OpMatBlocks : public DomainEleOp {
95 OpMatBlocks(std::string field_name, boost::shared_ptr<MatrixDouble> m,
96 double bulk_modulus_K, double shear_modulus_G,
97 MoFEM::Interface &m_field, Sev sev,
98 std::vector<const CubitMeshSets *> meshset_vec_ptr)
100 bulkModulusKDefault(bulk_modulus_K),
101 shearModulusGDefault(shear_modulus_G) {
102 std::fill(&(doEntities[MBEDGE]), &(doEntities[MBMAXTYPE]), false);
103 CHK_THROW_MESSAGE(extractBlockData(m_field, meshset_vec_ptr, sev),
104 "Can not get data from block");
105 }
106
107 MoFEMErrorCode doWork(int side, EntityType type,
110
111 for (auto &b : blockData) {
112
113 if (b.blockEnts.find(getFEEntityHandle()) != b.blockEnts.end()) {
114 CHKERR getMatDPtr(matDPtr, b.bulkModulusK, b.shearModulusG);
116 }
117 }
118
119 CHKERR getMatDPtr(matDPtr, bulkModulusKDefault, shearModulusGDefault);
121 }
122
123 private:
124 boost::shared_ptr<MatrixDouble> matDPtr;
125
126 struct BlockData {
127 double bulkModulusK;
128 double shearModulusG;
129 Range blockEnts;
130 };
131
132 double bulkModulusKDefault;
133 double shearModulusGDefault;
134 std::vector<BlockData> blockData;
135
137 extractBlockData(MoFEM::Interface &m_field,
138 std::vector<const CubitMeshSets *> meshset_vec_ptr,
139 Sev sev) {
141
142 for (auto m : meshset_vec_ptr) {
143 MOFEM_TAG_AND_LOG("WORLD", sev, "MatBlock") << *m;
144 std::vector<double> block_data;
145 CHKERR m->getAttributes(block_data);
146 if (block_data.size() != 2) {
147 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
148 "Expected that block has two attribute");
149 }
150 auto get_block_ents = [&]() {
151 Range ents;
152 CHKERR
153 m_field.get_moab().get_entities_by_handle(m->meshset, ents, true);
154 return ents;
155 };
156
157 double young_modulus = block_data[0];
158 double poisson_ratio = block_data[1];
159 double bulk_modulus_K = young_modulus / (3 * (1 - 2 * poisson_ratio));
160 double shear_modulus_G = young_modulus / (2 * (1 + poisson_ratio));
161
162 MOFEM_TAG_AND_LOG("WORLD", sev, "MatBlock")
163 << "E = " << young_modulus << " nu = " << poisson_ratio;
164
165 blockData.push_back(
166 {bulk_modulus_K, shear_modulus_G, get_block_ents()});
167 }
168 MOFEM_LOG_CHANNEL("WORLD");
170 }
171
172 MoFEMErrorCode getMatDPtr(boost::shared_ptr<MatrixDouble> mat_D_ptr,
173 double bulk_modulus_K, double shear_modulus_G) {
175 //! [Calculate elasticity tensor]
176 auto set_material_stiffness = [&]() {
182 double A = (SPACE_DIM == 2)
183 ? 2 * shear_modulus_G /
184 (bulk_modulus_K + (4. / 3.) * shear_modulus_G)
185 : 1;
186 auto t_D = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, 0>(*mat_D_ptr);
187 t_D(i, j, k, l) =
188 2 * shear_modulus_G * ((t_kd(i, k) ^ t_kd(j, l)) / 4.) +
189 A * (bulk_modulus_K - (2. / 3.) * shear_modulus_G) * t_kd(i, j) *
190 t_kd(k, l);
191 };
192 //! [Calculate elasticity tensor]
193 constexpr auto size_symm = (SPACE_DIM * (SPACE_DIM + 1)) / 2;
194 mat_D_ptr->resize(size_symm * size_symm, 1);
195 set_material_stiffness();
197 }
198 };
199
200 pipeline.push_back(new OpMatBlocks(
202
203 // Get blockset using regular expression
205
206 (boost::format("%s(.*)") % block_name).str()
207
208 ))
209
210 ));
211
213}
214
215//! [Run problem]
226}
227//! [Run problem]
228
229//! [Read mesh]
233 CHKERR simple->getOptions();
234 CHKERR simple->loadFile();
236}
237//! [Read mesh]
238
239//! [Set up problem]
243 // Add field
244 CHKERR simple->addDomainField("U", H1, AINSWORTH_LEGENDRE_BASE, SPACE_DIM);
245 CHKERR simple->addBoundaryField("U", H1, AINSWORTH_LEGENDRE_BASE, SPACE_DIM);
246 int order = 3;
247 CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-order", &order, PETSC_NULL);
248
249 CHKERR simple->addDataField("GEOMETRY", H1, AINSWORTH_LEGENDRE_BASE,
250 SPACE_DIM);
251
252 CHKERR simple->setFieldOrder("U", order);
253 CHKERR simple->setFieldOrder("GEOMETRY", 2);
254 CHKERR simple->setUp();
255
256 auto project_ho_geometry = [&]() {
257 Projection10NodeCoordsOnField ent_method(mField, "GEOMETRY");
258 return mField.loop_dofs("GEOMETRY", ent_method);
259 };
260 CHKERR project_ho_geometry();
261
263}
264//! [Set up problem]
265
266//! [Boundary condition]
269 auto *pipeline_mng = mField.getInterface<PipelineManager>();
271 auto bc_mng = mField.getInterface<BcManager>();
272
273 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "REMOVE_X",
274 "U", 0, 0);
275 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "REMOVE_Y",
276 "U", 1, 1);
277 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "REMOVE_Z",
278 "U", 2, 2);
279 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(),
280 "REMOVE_ALL", "U", 0, 3);
281
283 pipeline_mng->getOpDomainRhsPipeline(), {H1}, "GEOMETRY");
285 pipeline_mng->getOpBoundaryRhsPipeline(), {NOSPACE}, "GEOMETRY");
286
287 // Infernal forces
288 auto mat_grad_ptr = boost::make_shared<MatrixDouble>();
289 auto mat_strain_ptr = boost::make_shared<MatrixDouble>();
290 auto mat_stress_ptr = boost::make_shared<MatrixDouble>();
291 pipeline_mng->getOpDomainRhsPipeline().push_back(
293 mat_grad_ptr));
294 pipeline_mng->getOpDomainRhsPipeline().push_back(
295 new OpSymmetrizeTensor<SPACE_DIM>("U", mat_grad_ptr, mat_strain_ptr));
296
297 auto mat_D_ptr = boost::make_shared<MatrixDouble>();
298 CHKERR addMatBlockOps(pipeline_mng->getOpDomainRhsPipeline(), "U",
299 "MAT_ELASTIC", mat_D_ptr, Sev::inform);
300 pipeline_mng->getOpDomainRhsPipeline().push_back(
302 "U", mat_strain_ptr, mat_stress_ptr, mat_D_ptr));
303 pipeline_mng->getOpDomainRhsPipeline().push_back(
304 new OpInternalForce("U", mat_stress_ptr,
305 [](double, double, double) constexpr { return -1; }));
306
307 // Body forces
309 pipeline_mng->getOpDomainRhsPipeline(), mField, "U", Sev::inform);
310 // Add force boundary condition
312 pipeline_mng->getOpBoundaryRhsPipeline(), mField, "U", -1, Sev::inform);
314 pipeline_mng->getOpBoundaryLhsPipeline(), mField, "U", Sev::verbose);
315
316 // Essential boundary condition
317 CHKERR bc_mng->removeBlockDOFsOnEntities<DisplacementCubitBcData>(
318 simple->getProblemName(), "U");
319
320 auto get_pre_proc_hook = [&]() {
322 mField, pipeline_mng->getDomainRhsFE(), {});
323 };
324
325 pipeline_mng->getDomainRhsFE()->preProcessHook = get_pre_proc_hook();
326 pipeline_mng->getBoundaryRhsFE()->preProcessHook = get_pre_proc_hook();
327
329}
330//! [Boundary condition]
331
332//! [Push operators to pipeline]
336
338 pipeline_mng->getOpDomainLhsPipeline(), {H1}, "GEOMETRY");
339
340 auto mat_D_ptr = boost::make_shared<MatrixDouble>();
341 CHKERR addMatBlockOps(pipeline_mng->getOpDomainLhsPipeline(), "U",
342 "MAT_ELASTIC", mat_D_ptr, Sev::verbose);
343 pipeline_mng->getOpDomainLhsPipeline().push_back(
344 new OpK("U", "U", mat_D_ptr));
345
346 auto integration_rule = [](int, int, int approx_order) {
347 return 2 * approx_order;
348 };
353
355}
356//! [Push operators to pipeline]
357
358//! [Solve]
361 auto *simple = mField.getInterface<Simple>();
362 auto *pipeline_mng = mField.getInterface<PipelineManager>();
363 auto solver = pipeline_mng->createKSP();
364 CHKERR KSPSetFromOptions(solver);
365 CHKERR KSPSetUp(solver);
366
367 auto dm = simple->getDM();
368 auto D = smartCreateDMVector(dm);
369 auto F = smartVectorDuplicate(D);
370
371 CHKERR KSPSolve(solver, F, D);
372 CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
373 CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
374 CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
376}
377//! [Solve]
378
379//! [Postprocess results]
383 auto det_ptr = boost::make_shared<VectorDouble>();
384 auto jac_ptr = boost::make_shared<MatrixDouble>();
385 auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
386 pipeline_mng->getDomainRhsFE().reset();
387 pipeline_mng->getDomainLhsFE().reset();
388 pipeline_mng->getBoundaryRhsFE().reset();
389 pipeline_mng->getBoundaryLhsFE().reset();
390
391 auto post_proc_fe = boost::make_shared<PostProcEle>(mField);
392
394 post_proc_fe->getOpPtrVector(), {H1}, "GEOMETRY");
395
396 auto mat_grad_ptr = boost::make_shared<MatrixDouble>();
397 auto mat_strain_ptr = boost::make_shared<MatrixDouble>();
398 auto mat_stress_ptr = boost::make_shared<MatrixDouble>();
399
400 post_proc_fe->getOpPtrVector().push_back(
402 mat_grad_ptr));
403 post_proc_fe->getOpPtrVector().push_back(
404 new OpSymmetrizeTensor<SPACE_DIM>("U", mat_grad_ptr, mat_strain_ptr));
405
406 auto mat_D_ptr = boost::make_shared<MatrixDouble>();
407 CHKERR addMatBlockOps(post_proc_fe->getOpPtrVector(), "U", "MAT_ELASTIC",
408 mat_D_ptr, Sev::verbose);
409 post_proc_fe->getOpPtrVector().push_back(
411 "U", mat_strain_ptr, mat_stress_ptr, mat_D_ptr));
412
413 auto u_ptr = boost::make_shared<MatrixDouble>();
414 post_proc_fe->getOpPtrVector().push_back(
416 auto x_ptr = boost::make_shared<MatrixDouble>();
417 post_proc_fe->getOpPtrVector().push_back(
418 new OpCalculateVectorFieldValues<SPACE_DIM>("GEOMETRY", x_ptr));
419
421
422 post_proc_fe->getOpPtrVector().push_back(
423
424 new OpPPMap(
425
426 post_proc_fe->getPostProcMesh(), post_proc_fe->getMapGaussPts(),
427
428 {},
429
430 {{"U", u_ptr}, {"GEOMETRY", x_ptr}},
431
432 {},
433
434 {{"STRAIN", mat_strain_ptr}, {"STRESS", mat_stress_ptr}}
435
436 )
437
438 );
439
440 pipeline_mng->getBoundaryRhsFE().reset();
441 pipeline_mng->getDomainRhsFE() = post_proc_fe;
442 CHKERR pipeline_mng->loopFiniteElements();
443 CHKERR post_proc_fe->writeFile("out_elastic.h5m");
445}
446//! [Postprocessing results]
447
448//! [Check]
450 MOFEM_LOG_CHANNEL("WORLD");
454 pipeline_mng->getDomainRhsFE().reset();
455 pipeline_mng->getDomainLhsFE().reset();
456 pipeline_mng->getBoundaryRhsFE().reset();
457 pipeline_mng->getBoundaryLhsFE().reset();
458
459 auto integration_rule = [](int, int, int p_data) { return 2 * p_data; };
460 CHKERR pipeline_mng->setDomainRhsIntegrationRule(integration_rule);
461 CHKERR pipeline_mng->setBoundaryRhsIntegrationRule(integration_rule);
462
464 pipeline_mng->getOpDomainRhsPipeline(), {H1}, "GEOMETRY");
466 pipeline_mng->getOpBoundaryRhsPipeline(), {NOSPACE}, "GEOMETRY");
467
468 auto mat_grad_ptr = boost::make_shared<MatrixDouble>();
469 auto mat_strain_ptr = boost::make_shared<MatrixDouble>();
470 auto mat_stress_ptr = boost::make_shared<MatrixDouble>();
471
472 pipeline_mng->getOpDomainRhsPipeline().push_back(
474 mat_grad_ptr));
475 pipeline_mng->getOpDomainRhsPipeline().push_back(
476 new OpSymmetrizeTensor<SPACE_DIM>("U", mat_grad_ptr, mat_strain_ptr));
477
478 auto mat_D_ptr = boost::make_shared<MatrixDouble>();
479 CHKERR addMatBlockOps(pipeline_mng->getOpDomainRhsPipeline(), "U",
480 "MAT_ELASTIC", mat_D_ptr, Sev::verbose);
481 pipeline_mng->getOpDomainRhsPipeline().push_back(
483 "U", mat_strain_ptr, mat_stress_ptr, mat_D_ptr));
484 pipeline_mng->getOpDomainRhsPipeline().push_back(
485 new OpInternalForce("U", mat_stress_ptr));
486
488 pipeline_mng->getOpDomainRhsPipeline(), mField, "U", Sev::verbose);
490 pipeline_mng->getOpBoundaryRhsPipeline(), mField, "U", 1, Sev::verbose);
491
492 auto dm = simple->getDM();
493 auto res = smartCreateDMVector(dm);
494 pipeline_mng->getDomainRhsFE()->ksp_f = res;
495 pipeline_mng->getBoundaryRhsFE()->ksp_f = res;
496
497 CHKERR VecZeroEntries(res);
498
499 CHKERR mField.getInterface<FieldBlas>()->fieldScale(-1, "U");
500 CHKERR pipeline_mng->loopFiniteElements();
501 CHKERR mField.getInterface<FieldBlas>()->fieldScale(-1, "U");
502
503 CHKERR VecGhostUpdateBegin(res, ADD_VALUES, SCATTER_REVERSE);
504 CHKERR VecGhostUpdateEnd(res, ADD_VALUES, SCATTER_REVERSE);
505 CHKERR VecAssemblyBegin(res);
506 CHKERR VecAssemblyEnd(res);
507
508 double nrm2;
509 CHKERR VecNorm(res, NORM_2, &nrm2);
510 MOFEM_LOG_CHANNEL("WORLD");
511 MOFEM_LOG_C("WORLD", Sev::inform, "residual = %3.4e\n", nrm2);
512 constexpr double eps = 1e-8;
513 if (nrm2 > eps)
514 SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "Residual is not zero");
515
517}
518//! [Check]
519
520static char help[] = "...\n\n";
521
522int main(int argc, char *argv[]) {
523
524 // Initialisation of MoFEM/PETSc and MOAB data structures
525 const char param_file[] = "param_file.petsc";
527
528 try {
529
530 //! [Register MoFEM discrete manager in PETSc]
531 DMType dm_name = "DMMOFEM";
532 CHKERR DMRegister_MoFEM(dm_name);
533 //! [Register MoFEM discrete manager in PETSc
534
535 //! [Create MoAB]
536 moab::Core mb_instance; ///< mesh database
537 moab::Interface &moab = mb_instance; ///< mesh database interface
538 //! [Create MoAB]
539
540 //! [Create MoFEM]
541 MoFEM::Core core(moab); ///< finite element database
542 MoFEM::Interface &m_field = core; ///< finite element database interface
543 //! [Create MoFEM]
544
545 //! [Example]
546 Example ex(m_field);
547 CHKERR ex.runProblem();
548 //! [Example]
549 }
551
553}
std::string param_file
Implementation of elastic spring bc.
#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
int main()
Definition: adol-c_atom.cpp:46
static const double eps
constexpr int SPACE_DIM
DomainEle::UserDataOperator DomainEleOp
Finire element operator type.
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
ElementsAndOps< SPACE_DIM >::BoundaryEle BoundaryEle
Kronecker Delta class symmetric.
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:372
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:595
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
@ 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
static char help[]
[Check]
Definition: elastic.cpp:520
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
constexpr double bulk_modulus_K
Definition: elastic.cpp:62
constexpr AssemblyType A
Definition: elastic.cpp:30
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: DMMoFEM.cpp:491
PetscErrorCode DMRegister_MoFEM(const char sname[])
Register MoFEM problem.
Definition: DMMoFEM.cpp:47
auto smartCreateDMVector(DM dm)
Get smart vector from DM.
Definition: DMMoFEM.hpp:978
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
double D
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: Common.hpp:10
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.
constexpr AssemblyType A
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
double young_modulus
Definition: plastic.cpp:96
#define EXECUTABLE_DIMENSION
Definition: plastic.cpp:10
constexpr auto size_symm
Definition: plastic.cpp:33
constexpr auto field_name
static constexpr int approx_order
[Example]
Definition: plastic.cpp:138
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:145
MoFEMErrorCode setupProblem()
MoFEMErrorCode outputResults()
Add operators pushing bases from local to physical configuration.
Simple interface for fast problem set-up.
Definition: BcManager.hpp:24
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)
Class (Function) to enforce essential constrains.
Definition: Essential.hpp:39
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()
MoFEM::FaceElementForcesAndSourcesCore FaceEle
boost::shared_ptr< FEMethod > & getBoundaryLhsFE()
MoFEMErrorCode setDomainRhsIntegrationRule(RuleHookFun rule)
MoFEMErrorCode setBoundaryLhsIntegrationRule(RuleHookFun rule)
MoFEMErrorCode setBoundaryRhsIntegrationRule(RuleHookFun rule)
MoFEM::EdgeElementForcesAndSourcesCore EdgeEle
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:27
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.