v0.13.1
Loading...
Searching...
No Matches
plastic.cpp
Go to the documentation of this file.
1/**
2 * \file plastic.cpp
3 * \example plastic.cpp
4 *
5 * Plasticity in 2d and 3d
6 *
7 */
8
9#ifndef EXECUTABLE_DIMENSION
10#define EXECUTABLE_DIMENSION 3
11#endif
12
13#include <MoFEM.hpp>
14#include <MatrixFunction.hpp>
15#include <IntegrationRules.hpp>
16
17#ifdef __cplusplus
18extern "C" {
19#endif
20#include <cblas.h>
21#include <lapack_wrap.h>
22// #include <gm_rule.h>
23#include <quad.h>
24#ifdef __cplusplus
25}
26#endif
27
28using namespace MoFEM;
29
30template <int DIM> struct ElementsAndOps {};
31
32template <> struct ElementsAndOps<2> {
35};
36
37template <> struct ElementsAndOps<3> {
40};
41
42constexpr int SPACE_DIM =
43 EXECUTABLE_DIMENSION; //< Space dimension of problem, mesh
44constexpr auto size_symm = (SPACE_DIM * (SPACE_DIM + 1)) / 2;
45
46constexpr AssemblyType A = AssemblyType::PETSC; //< selected assembly type
47constexpr IntegrationType G =
48 IntegrationType::GAUSS; //< selected integration type
49
56
58
59//! [Only used with Hooke equation (linear material model)]
61 GAUSS>::OpGradSymTensorGrad<1, SPACE_DIM, SPACE_DIM, 0>;
64//! [Only used with Hooke equation (linear material model)]
65
66//! [Only used for dynamics]
71//! [Only used for dynamics]
72
73//! [Only used with Hencky/nonlinear material]
75 GAUSS>::OpGradTensorGrad<1, SPACE_DIM, SPACE_DIM, 1>;
78//! [Only used with Hencky/nonlinear material]
79
80//! [Essential boundary conditions]
87//! [Essential boundary conditions]
89
93
95using OpForce =
97
102
103PetscBool is_large_strains = PETSC_TRUE;
104
105double scale = 1.;
106
107double young_modulus = 206913;
108double poisson_ratio = 0.29;
109double rho = 0;
110double sigmaY = 450;
111double H = 129;
112double visH = 0;
113double cn = 1;
114double zeta = 5e-2;
115double eta = 1e2;
116double Qinf = 265;
117double b_iso = 16.93;
118
119int order = 2; ///< Order if fixed.
120int geom_order = 2; ///< Order if fixed.
121
122inline long double hardening(long double tau) {
123 return H * tau + Qinf * (1. - std::exp(-b_iso * tau)) + sigmaY;
124}
125
126inline long double hardening_dtau(long double tau) {
127 return H + Qinf * b_iso * std::exp(-b_iso * tau);
128}
129
130inline long double hardening_dtau2(long double tau) {
131 return -(Qinf * (b_iso * b_iso)) * std::exp(-b_iso * tau);
132}
133
134#include <HenckyOps.hpp>
135#include <PlasticOps.hpp>
136
137using namespace PlasticOps;
138using namespace HenckyOps;
139struct Example {
140
141 Example(MoFEM::Interface &m_field) : mField(m_field) {}
142
144
145private:
147
153
154 boost::shared_ptr<PlasticOps::CommonData> commonPlasticDataPtr;
155 boost::shared_ptr<HenckyOps::CommonData> commonHenckyDataPtr;
156 boost::shared_ptr<DomainEle> reactionFe;
157
158 std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter>> uXScatter;
159 std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter>> uYScatter;
160 std::tuple<SmartPetscObj<Vec>, SmartPetscObj<VecScatter>> uZScatter;
161
162 boost::shared_ptr<std::vector<unsigned char>> boundaryMarker;
163 boost::shared_ptr<std::vector<unsigned char>> reactionMarker;
164
167 double getScale(const double time) {
168 return scale * MoFEM::TimeScale::getScale(time);
169 };
170 };
171};
172
173//! [Run problem]
178 CHKERR bC();
179 CHKERR OPs();
180 CHKERR tsSolve();
182}
183//! [Run problem]
184
185//! [Set up problem]
189
190 Range domain_ents;
191 CHKERR mField.get_moab().get_entities_by_dimension(0, SPACE_DIM, domain_ents,
192 true);
193 auto get_ents_by_dim = [&](const auto dim) {
194 if (dim == SPACE_DIM) {
195 return domain_ents;
196 } else {
197 Range ents;
198 if (dim == 0)
199 CHKERR mField.get_moab().get_connectivity(domain_ents, ents, true);
200 else
201 CHKERR mField.get_moab().get_entities_by_dimension(0, dim, ents, true);
202 return ents;
203 }
204 };
205
206 auto get_base = [&]() {
207 auto domain_ents = get_ents_by_dim(SPACE_DIM);
208 if (domain_ents.empty())
209 CHK_THROW_MESSAGE(MOFEM_NOT_FOUND, "Empty mesh");
210 const auto type = type_from_handle(domain_ents[0]);
211 switch (type) {
212 case MBQUAD:
214 case MBHEX:
216 case MBTRI:
218 case MBTET:
220 default:
221 CHK_THROW_MESSAGE(MOFEM_NOT_FOUND, "Element type not handled");
222 }
223 return NOBASE;
224 };
225
226 const auto base = get_base();
227 MOFEM_LOG("WORLD", Sev::inform) << "Base " << ApproximationBaseNames[base];
228
229 CHKERR simple->addDomainField("U", H1, base, SPACE_DIM);
230 CHKERR simple->addDomainField("TAU", L2, base, 1);
231 CHKERR simple->addDomainField("EP", L2, base, size_symm);
232 CHKERR simple->addBoundaryField("U", H1, base, SPACE_DIM);
233
234 CHKERR simple->addDataField("GEOMETRY", H1, base, SPACE_DIM);
235
236 auto ents = get_ents_by_dim(0);
237 ents.merge(get_ents_by_dim(1));
238 // ents.merge(get_ents_by_dim(2));
239 CHKERR simple->setFieldOrder("U", order, &ents);
240 CHKERR simple->setFieldOrder("EP", order - 1);
241 CHKERR simple->setFieldOrder("TAU", order - 2);
242
243 CHKERR simple->setFieldOrder("GEOMETRY", geom_order);
244
245 CHKERR simple->setUp();
246 // CHKERR simple->addFieldToEmptyFieldBlocks("U", "TAU");
247
248 auto project_ho_geometry = [&]() {
249 Projection10NodeCoordsOnField ent_method(mField, "GEOMETRY");
250 return mField.loop_dofs("GEOMETRY", ent_method);
251 };
252 CHKERR project_ho_geometry();
253
255}
256//! [Set up problem]
257
258//! [Create common data]
261
262 auto get_command_line_parameters = [&]() {
264 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-scale", &scale, PETSC_NULL);
265 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-rho", &rho, PETSC_NULL);
266 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-young_modulus",
267 &young_modulus, PETSC_NULL);
268 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-poisson_ratio",
269 &poisson_ratio, PETSC_NULL);
270 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-hardening", &H, PETSC_NULL);
271 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-hardening_viscous", &visH,
272 PETSC_NULL);
273 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-yield_stress", &sigmaY,
274 PETSC_NULL);
275 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-cn", &cn, PETSC_NULL);
276 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-zeta", &zeta, PETSC_NULL);
277 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-eta", &zeta, PETSC_NULL);
278 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-Qinf", &Qinf, PETSC_NULL);
279 CHKERR PetscOptionsGetScalar(PETSC_NULL, "", "-b_iso", &b_iso, PETSC_NULL);
280 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-large_strains",
281 &is_large_strains, PETSC_NULL);
282
283 CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-order", &order, PETSC_NULL);
284 CHKERR PetscOptionsGetInt(PETSC_NULL, "", "-geom_order", &geom_order,
285 PETSC_NULL);
286
287 MOFEM_LOG("EXAMPLE", Sev::inform) << "Young modulus " << young_modulus;
288 MOFEM_LOG("EXAMPLE", Sev::inform) << "Poisson ratio " << poisson_ratio;
289 MOFEM_LOG("EXAMPLE", Sev::inform) << "Yield stress " << sigmaY;
290 MOFEM_LOG("EXAMPLE", Sev::inform) << "Hardening " << H;
291 MOFEM_LOG("EXAMPLE", Sev::inform) << "Viscous hardening " << visH;
292 MOFEM_LOG("EXAMPLE", Sev::inform) << "Saturation yield stress " << Qinf;
293 MOFEM_LOG("EXAMPLE", Sev::inform) << "Saturation exponent " << b_iso;
294 MOFEM_LOG("EXAMPLE", Sev::inform) << "cn " << cn;
295 MOFEM_LOG("EXAMPLE", Sev::inform) << "zeta " << zeta;
296 MOFEM_LOG("EXAMPLE", Sev::inform) << "eta " << zeta;
297
298 MOFEM_LOG("EXAMPLE", Sev::inform) << "order " << order;
299 MOFEM_LOG("EXAMPLE", Sev::inform) << "geom order " << geom_order;
300
301 PetscBool is_scale = PETSC_TRUE;
302 CHKERR PetscOptionsGetBool(PETSC_NULL, "", "-is_scale", &is_scale,
303 PETSC_NULL);
304 if (is_scale) {
307 rho *= scale;
308 sigmaY *= scale;
309 H *= scale;
310 Qinf *= scale;
311 visH *= scale;
312
313 MOFEM_LOG("EXAMPLE", Sev::inform)
314 << "Scaled Young modulus " << young_modulus;
315 MOFEM_LOG("EXAMPLE", Sev::inform)
316 << "Scaled Poisson ratio " << poisson_ratio;
317 MOFEM_LOG("EXAMPLE", Sev::inform) << "Scaled Yield stress " << sigmaY;
318 MOFEM_LOG("EXAMPLE", Sev::inform) << "Scaled Hardening " << H;
319 MOFEM_LOG("EXAMPLE", Sev::inform) << "Scaled Viscous hardening " << visH;
320 MOFEM_LOG("EXAMPLE", Sev::inform)
321 << "Scaled Saturation yield stress " << Qinf;
322 }
323
325 };
326
327 auto set_matrial_stiffness = [&]() {
334 const double bulk_modulus_K = young_modulus / (3 * (1 - 2 * poisson_ratio));
335 const double shear_modulus_G = young_modulus / (2 * (1 + poisson_ratio));
336
337 // Plane stress or when 1, plane strain or 3d
338 const double A = (SPACE_DIM == 2)
339 ? 2 * shear_modulus_G /
340 (bulk_modulus_K + (4. / 3.) * shear_modulus_G)
341 : 1;
342
343 auto t_D = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, 0>(
344 *commonPlasticDataPtr->mDPtr);
345 auto t_D_axiator = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, 0>(
346 *commonPlasticDataPtr->mDPtr_Axiator);
347 auto t_D_deviator = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, 0>(
348 *commonPlasticDataPtr->mDPtr_Deviator);
349
350 constexpr double third = boost::math::constants::third<double>();
351 t_D_axiator(i, j, k, l) = A *
352 (bulk_modulus_K - (2. / 3.) * shear_modulus_G) *
353 t_kd(i, j) * t_kd(k, l);
354 t_D_deviator(i, j, k, l) =
355 2 * shear_modulus_G * ((t_kd(i, k) ^ t_kd(j, l)) / 4.);
356 t_D(i, j, k, l) = t_D_axiator(i, j, k, l) + t_D_deviator(i, j, k, l);
357
359 };
360
361 auto make_d_mat = []() {
362 return boost::make_shared<MatrixDouble>(size_symm * size_symm, 1);
363 };
364
365 commonPlasticDataPtr = boost::make_shared<PlasticOps::CommonData>();
366 commonPlasticDataPtr->mDPtr = make_d_mat();
367 commonPlasticDataPtr->mDPtr_Axiator = make_d_mat();
368 commonPlasticDataPtr->mDPtr_Deviator = make_d_mat();
369
370 commonPlasticDataPtr->mGradPtr = boost::make_shared<MatrixDouble>();
371 commonPlasticDataPtr->mStrainPtr = boost::make_shared<MatrixDouble>();
372 commonPlasticDataPtr->mStressPtr = boost::make_shared<MatrixDouble>();
373
374 CHKERR get_command_line_parameters();
375 CHKERR set_matrial_stiffness();
376
377 if (is_large_strains) {
378 commonHenckyDataPtr = boost::make_shared<HenckyOps::CommonData>();
379 commonHenckyDataPtr->matGradPtr = commonPlasticDataPtr->mGradPtr;
381 commonHenckyDataPtr->matLogCPlastic =
382 commonPlasticDataPtr->getPlasticStrainPtr();
383 commonPlasticDataPtr->mStrainPtr = commonHenckyDataPtr->getMatLogC();
384 commonPlasticDataPtr->mStressPtr =
385 commonHenckyDataPtr->getMatHenckyStress();
386 }
387
389}
390//! [Create common data]
391
392//! [Boundary condition]
395
397 auto bc_mng = mField.getInterface<BcManager>();
398 auto prb_mng = mField.getInterface<ProblemsManager>();
399
400 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "REMOVE_X",
401 "U", 0, 0);
402 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "REMOVE_Y",
403 "U", 1, 1);
404 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(), "REMOVE_Z",
405 "U", 2, 2);
406 CHKERR bc_mng->removeBlockDOFsOnEntities(simple->getProblemName(),
407 "REMOVE_ALL", "U", 0, 3);
408
409 CHKERR bc_mng->pushMarkDOFsOnEntities<DisplacementCubitBcData>(
410 simple->getProblemName(), "U");
411
412 auto &bc_map = bc_mng->getBcMapByBlockName();
413 boundaryMarker = bc_mng->getMergedBlocksMarker(vector<string>{"FIX_"});
414
415 CHKERR bc_mng->pushMarkDOFsOnEntities(simple->getProblemName(), "REACTION",
416 "U", 0, 3);
417
418 for (auto bc : bc_map)
419 MOFEM_LOG("EXAMPLE", Sev::verbose) << "Marker " << bc.first;
420
421 // OK. We have problem with GMesh, it adding empty characters at the end of
422 // block. So first block is search by regexp. popMarkDOFsOnEntities should
423 // work with regexp.
424 std::string reaction_block_set;
425 for (auto bc : bc_map) {
426 if (bc_mng->checkBlock(bc, "REACTION")) {
427 reaction_block_set = bc.first;
428 break;
429 }
430 }
431
432 if (auto bc = bc_mng->popMarkDOFsOnEntities(reaction_block_set)) {
433 reactionMarker = bc->getBcMarkersPtr();
434
435 // Only take reaction from nodes
436 Range nodes;
437 CHKERR mField.get_moab().get_entities_by_type(0, MBVERTEX, nodes, true);
438 CHKERR prb_mng->markDofs(simple->getProblemName(), ROW,
439 ProblemsManager::MarkOP::AND, nodes,
441
442 } else {
443 MOFEM_LOG("EXAMPLE", Sev::warning) << "REACTION blockset does not exist";
444 }
445
446 if (!reactionMarker) {
447 MOFEM_LOG("EXAMPLE", Sev::warning) << "REACTION blockset does not exist";
448 }
449
451}
452//! [Boundary condition]
453
454//! [Push operators to pipeline]
457 auto pipeline_mng = mField.getInterface<PipelineManager>();
459 auto bc_mng = mField.getInterface<BcManager>();
460
461 auto add_domain_base_ops = [&](auto &pipeline) {
463
465 "GEOMETRY");
466
467 pipeline.push_back(new OpCalculateScalarFieldValuesDot(
468 "TAU", commonPlasticDataPtr->getPlasticTauDotPtr()));
469 pipeline.push_back(new OpCalculateScalarFieldValues(
470 "TAU", commonPlasticDataPtr->getPlasticTauPtr()));
472 "EP", commonPlasticDataPtr->getPlasticStrainPtr()));
474 "EP", commonPlasticDataPtr->getPlasticStrainDotPtr()));
476 "U", commonPlasticDataPtr->mGradPtr));
477
479 };
480
481 auto add_domain_stress_ops = [&](auto &pipeline, auto m_D_ptr) {
483
484 if (is_large_strains) {
485
486 if (commonPlasticDataPtr->mGradPtr != commonHenckyDataPtr->matGradPtr)
487 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
488 "Wrong pointer for grad");
489
490 pipeline.push_back(
492 pipeline.push_back(
494 pipeline.push_back(
496 pipeline.push_back(new OpCalculateHenckyPlasticStress<SPACE_DIM>(
497 "U", commonHenckyDataPtr, m_D_ptr));
498 pipeline.push_back(
500
501 } else {
502 pipeline.push_back(
504 commonPlasticDataPtr->mStrainPtr));
505 pipeline.push_back(
506 new OpPlasticStress("U", commonPlasticDataPtr, m_D_ptr, 1));
507 }
508
509 if (m_D_ptr != commonPlasticDataPtr->mDPtr_Axiator) {
510 pipeline.push_back(
512 pipeline.push_back(
514 }
515
517 };
518
519 auto add_domain_ops_lhs_mechanical = [&](auto &pipeline, auto m_D_ptr) {
521 pipeline.push_back(new OpSetBc("U", true, boundaryMarker));
522
523 if (is_large_strains) {
524 pipeline.push_back(
526 pipeline.push_back(
527 new OpKPiola("U", "U", commonHenckyDataPtr->getMatTangent()));
529 "U", "EP", commonPlasticDataPtr, commonHenckyDataPtr, m_D_ptr));
530 // pipeline.push_back(new OpCalculateArgLhs_LogStrain_dUdTau(
531 // "U", "TAU", commonPlasticDataPtr, commonHenckyDataPtr));
532 // pipeline.push_back(new OpCalculateArgLhs_LogStrain_dUdU(
533 // "U", "U", commonPlasticDataPtr, commonHenckyDataPtr));
534 } else {
535 pipeline.push_back(new OpKCauchy("U", "U", m_D_ptr));
536 pipeline.push_back(new OpCalculatePlasticInternalForceLhs_dEP(
537 "U", "EP", commonPlasticDataPtr, m_D_ptr));
538 // pipeline.push_back(
539 // new OpCalculateArgLhs_dUdTau("U", "TAU", commonPlasticDataPtr));
540 // pipeline.push_back(
541 // new OpCalculateArgLhs_dUdU("U", "U", commonPlasticDataPtr));
542 }
543
544 pipeline.push_back(new OpUnSetBc("U"));
546 };
547
548 auto add_domain_ops_rhs_mechanical = [&](auto &pipeline) {
550 pipeline.push_back(new OpSetBc("U", true, boundaryMarker));
551
553 pipeline, mField, "U", {boost::make_shared<PlasticityTimeScale>()},
554 "BODY_FORCE", Sev::inform);
555
556 // Calculate internal forces
557 if (is_large_strains) {
558 pipeline.push_back(new OpInternalForcePiola(
559 "U", commonHenckyDataPtr->getMatFirstPiolaStress()));
560 } else {
561 pipeline.push_back(
562 new OpInternalForceCauchy("U", commonPlasticDataPtr->mStressPtr));
563 }
564
565 pipeline.push_back(new OpUnSetBc("U"));
567 };
568
569 auto add_domain_ops_lhs_constrain = [&](auto &pipeline, auto m_D_ptr) {
571 pipeline.push_back(new OpSetBc("U", true, boundaryMarker));
572
573 if (is_large_strains) {
574 pipeline.push_back(new OpCalculateConstraintsLhs_LogStrain_dU(
575 "TAU", "U", commonPlasticDataPtr, commonHenckyDataPtr, m_D_ptr));
576 pipeline.push_back(new OpCalculatePlasticFlowLhs_LogStrain_dU(
577 "EP", "U", commonPlasticDataPtr, commonHenckyDataPtr, m_D_ptr));
578 } else {
579 pipeline.push_back(new OpCalculateConstraintsLhs_dU(
580 "TAU", "U", commonPlasticDataPtr, m_D_ptr));
581 pipeline.push_back(new OpCalculatePlasticFlowLhs_dU(
582 "EP", "U", commonPlasticDataPtr, m_D_ptr));
583 }
584
585 pipeline.push_back(new OpCalculatePlasticFlowLhs_dEP(
586 "EP", "EP", commonPlasticDataPtr, m_D_ptr));
587 pipeline.push_back(new OpCalculatePlasticFlowLhs_dTAU(
588 "EP", "TAU", commonPlasticDataPtr, m_D_ptr));
589 pipeline.push_back(new OpCalculateConstraintsLhs_dEP(
590 "TAU", "EP", commonPlasticDataPtr, m_D_ptr));
591 pipeline.push_back(
593
594 pipeline.push_back(new OpUnSetBc("U"));
596 };
597
598 auto add_domain_ops_rhs_constrain = [&](auto &pipeline, auto m_D_ptr) {
600 pipeline.push_back(new OpSetBc("U", true, boundaryMarker));
601
602 pipeline.push_back(
604 pipeline.push_back(
606
607 // if (is_large_strains) {
608 // pipeline.push_back(new OpCalculateArgRhs_LogStrain_dU(
609 // "U", commonPlasticDataPtr, commonHenckyDataPtr));
610 // } else {
611 // pipeline.push_back(new OpCalculateArgRhs_dU("U", commonPlasticDataPtr));
612 // }
613
614 pipeline.push_back(new OpUnSetBc("U"));
616 };
617
618 auto add_boundary_ops_lhs_mechanical = [&](auto &pipeline) {
622 mField, pipeline, simple->getProblemName(), "U");
624 };
625
626 auto add_boundary_ops_rhs_mechanical = [&](auto &pipeline) {
628
630 pipeline, {NOSPACE}, "GEOMETRY");
631
632 pipeline.push_back(new OpSetBc("U", true, boundaryMarker));
634 pipeline, mField, "U", {boost::make_shared<PlasticityTimeScale>()},
635 "FORCE", Sev::inform);
636
637 pipeline.push_back(new OpUnSetBc("U"));
638
639 auto u_mat_ptr = boost::make_shared<MatrixDouble>();
640 pipeline.push_back(
641 new OpCalculateVectorFieldValues<SPACE_DIM>("U", u_mat_ptr));
642
645 mField, pipeline, simple->getProblemName(), "U", u_mat_ptr,
646 {boost::make_shared<TimeScale>()}); // note displacements have no
647 // scaling
648
650 };
651
652 CHKERR add_domain_base_ops(pipeline_mng->getOpDomainLhsPipeline());
653 CHKERR add_domain_stress_ops(pipeline_mng->getOpDomainLhsPipeline(),
654 commonPlasticDataPtr->mDPtr);
655 CHKERR add_domain_ops_lhs_mechanical(pipeline_mng->getOpDomainLhsPipeline(),
656 commonPlasticDataPtr->mDPtr);
657 CHKERR add_domain_ops_lhs_constrain(pipeline_mng->getOpDomainLhsPipeline(),
658 commonPlasticDataPtr->mDPtr);
659 CHKERR
660 add_boundary_ops_lhs_mechanical(pipeline_mng->getOpBoundaryLhsPipeline());
661
662 CHKERR add_domain_base_ops(pipeline_mng->getOpDomainRhsPipeline());
663 CHKERR add_domain_stress_ops(pipeline_mng->getOpDomainRhsPipeline(),
664 commonPlasticDataPtr->mDPtr);
665 CHKERR add_domain_ops_rhs_mechanical(pipeline_mng->getOpDomainRhsPipeline());
666 CHKERR add_domain_ops_rhs_constrain(pipeline_mng->getOpDomainRhsPipeline(),
667 commonPlasticDataPtr->mDPtr);
668
669 // Boundary
670 CHKERR
671 add_boundary_ops_rhs_mechanical(pipeline_mng->getOpBoundaryRhsPipeline());
672
673 auto integration_rule_bc = [](int, int, int ao) { return 2 * ao; };
674
675 auto vol_rule = [](int, int, int ao) { return 2 * ao + geom_order - 1; };
676
677 CHKERR pipeline_mng->setDomainRhsIntegrationRule(vol_rule);
678 CHKERR pipeline_mng->setDomainLhsIntegrationRule(vol_rule);
679
680 CHKERR pipeline_mng->setBoundaryLhsIntegrationRule(integration_rule_bc);
681 CHKERR pipeline_mng->setBoundaryRhsIntegrationRule(integration_rule_bc);
682
683 auto create_reaction_pipeline = [&](auto &pipeline) {
685
686 if (reactionMarker) {
687
689 "GEOMETRY");
690
691 pipeline.push_back(
693 "U", commonPlasticDataPtr->mGradPtr));
695 "EP", commonPlasticDataPtr->getPlasticStrainPtr()));
696
697 if (is_large_strains) {
698
699 if (commonPlasticDataPtr->mGradPtr != commonHenckyDataPtr->matGradPtr)
700 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
701 "Wrong pointer for grad");
702
703 pipeline.push_back(
705 pipeline.push_back(
707 pipeline.push_back(
709 pipeline.push_back(new OpCalculateHenckyPlasticStress<SPACE_DIM>(
711 pipeline.push_back(
713
714 } else {
715 pipeline.push_back(new OpSymmetrizeTensor<SPACE_DIM>(
716 "U", commonPlasticDataPtr->mGradPtr,
717 commonPlasticDataPtr->mStrainPtr));
718 pipeline.push_back(new OpPlasticStress(
720 }
721
722 pipeline.push_back(new OpSetBc("U", false, reactionMarker));
723 // Calculate internal force
724 if (is_large_strains) {
725 pipeline.push_back(new OpInternalForcePiola(
726 "U", commonHenckyDataPtr->getMatFirstPiolaStress()));
727 } else {
728 pipeline.push_back(
729 new OpInternalForceCauchy("U", commonPlasticDataPtr->mStressPtr));
730 }
731 pipeline.push_back(new OpUnSetBc("U"));
732 }
733
735 };
736
737 reactionFe = boost::make_shared<DomainEle>(mField);
738 reactionFe->getRuleHook = vol_rule;
739
740 CHKERR create_reaction_pipeline(reactionFe->getOpPtrVector());
741
743}
744//! [Push operators to pipeline]
745
746//! [Solve]
749
752 ISManager *is_manager = mField.getInterface<ISManager>();
753
754 auto snes_ctx_ptr = smartGetDMSnesCtx(simple->getDM());
755
756 auto set_section_monitor = [&](auto solver) {
758 SNES snes;
759 CHKERR TSGetSNES(solver, &snes);
760 CHKERR SNESMonitorSet(snes,
761 (MoFEMErrorCode(*)(SNES, PetscInt, PetscReal,
763 (void *)(snes_ctx_ptr.get()), nullptr);
765 };
766
767 auto create_post_process_element = [&]() {
768 auto pp_fe = boost::make_shared<PostProcEle>(mField);
770
772 pp_fe->getOpPtrVector(), {H1, L2}, "GEOMETRY");
773
774 auto x_ptr = boost::make_shared<MatrixDouble>();
775 pp_fe->getOpPtrVector().push_back(
776 new OpCalculateVectorFieldValues<SPACE_DIM>("GEOMETRY", x_ptr));
777 auto u_ptr = boost::make_shared<MatrixDouble>();
778 pp_fe->getOpPtrVector().push_back(
780
781 pp_fe->getOpPtrVector().push_back(
783 "U", commonPlasticDataPtr->mGradPtr));
784 pp_fe->getOpPtrVector().push_back(new OpCalculateScalarFieldValues(
785 "TAU", commonPlasticDataPtr->getPlasticTauPtr()));
786 pp_fe->getOpPtrVector().push_back(
788 "EP", commonPlasticDataPtr->getPlasticStrainPtr()));
789
790 if (is_large_strains) {
791
792 if (commonPlasticDataPtr->mGradPtr != commonHenckyDataPtr->matGradPtr)
793 CHK_THROW_MESSAGE(MOFEM_DATA_INCONSISTENCY, "Wrong pointer for grad");
794
795 pp_fe->getOpPtrVector().push_back(
797 pp_fe->getOpPtrVector().push_back(
799 pp_fe->getOpPtrVector().push_back(
801 pp_fe->getOpPtrVector().push_back(
804 1 /*scale*/));
805 pp_fe->getOpPtrVector().push_back(
807
808 pp_fe->getOpPtrVector().push_back(
809
810 new OpPPMap(
811
812 pp_fe->getPostProcMesh(), pp_fe->getMapGaussPts(),
813
814 {},
815
816 {{"U", u_ptr}, {"GEOMETRY", x_ptr}},
817
818 {{"GRAD", commonPlasticDataPtr->mGradPtr},
819 {"FIRST_PIOLA", commonHenckyDataPtr->getMatFirstPiolaStress()}},
820
821 {}
822
823 )
824
825 );
826
827 } else {
828 pp_fe->getOpPtrVector().push_back(
830 commonPlasticDataPtr->mStrainPtr));
831 pp_fe->getOpPtrVector().push_back(new OpPlasticStress(
832 "U", commonPlasticDataPtr, commonPlasticDataPtr->mDPtr, 1 /*scale*/));
833
834 pp_fe->getOpPtrVector().push_back(
835
836 new OpPPMap(
837
838 pp_fe->getPostProcMesh(), pp_fe->getMapGaussPts(),
839
840 {},
841
842 {{"U", u_ptr}, {"GEOMETRY", x_ptr}},
843
844 {},
845
846 {{"STRAIN", commonPlasticDataPtr->mStrainPtr},
847 {"STRESS", commonPlasticDataPtr->mStressPtr}}
848
849 )
850
851 );
852 }
853
854 pp_fe->getOpPtrVector().push_back(
855 new OpCalculatePlasticSurface("U", commonPlasticDataPtr));
856
857 pp_fe->getOpPtrVector().push_back(
858
859 new OpPPMap(
860
861 pp_fe->getPostProcMesh(), pp_fe->getMapGaussPts(),
862
863 {{"PLASTIC_SURFACE", commonPlasticDataPtr->getPlasticSurfacePtr()},
864 {"PLASTIC_MULTIPLIER", commonPlasticDataPtr->getPlasticTauPtr()}},
865
866 {},
867
868 {},
869
870 {{"PLASTIC_STRAIN", commonPlasticDataPtr->getPlasticStrainPtr()},
871 {"PLASTIC_FLOW", commonPlasticDataPtr->getPlasticFlowPtr()}}
872
873 )
874
875 );
876
877 return pp_fe;
878 };
879
880 auto scatter_create = [&](auto D, auto coeff) {
882 CHKERR is_manager->isCreateProblemFieldAndRank(simple->getProblemName(),
883 ROW, "U", coeff, coeff, is);
884 int loc_size;
885 CHKERR ISGetLocalSize(is, &loc_size);
886 Vec v;
887 CHKERR VecCreateMPI(mField.get_comm(), loc_size, PETSC_DETERMINE, &v);
888 VecScatter scatter;
889 CHKERR VecScatterCreate(D, is, v, PETSC_NULL, &scatter);
890 return std::make_tuple(SmartPetscObj<Vec>(v),
892 };
893
894 auto set_time_monitor = [&](auto dm, auto solver) {
896 boost::shared_ptr<Monitor> monitor_ptr(
897 new Monitor(dm, create_post_process_element(), reactionFe, uXScatter,
898 uYScatter, uZScatter));
899 boost::shared_ptr<ForcesAndSourcesCore> null;
900 CHKERR DMMoFEMTSSetMonitor(dm, solver, simple->getDomainFEName(),
901 monitor_ptr, null, null);
903 };
904
905 auto set_fieldsplit_preconditioner = [&](auto solver) {
907
908 SNES snes;
909 CHKERR TSGetSNES(solver, &snes);
910 KSP ksp;
911 CHKERR SNESGetKSP(snes, &ksp);
912 PC pc;
913 CHKERR KSPGetPC(ksp, &pc);
914 PetscBool is_pcfs = PETSC_FALSE;
915 PetscObjectTypeCompare((PetscObject)pc, PCFIELDSPLIT, &is_pcfs);
916
917 // Setup fieldsplit (block) solver - optional: yes/no
918 if (is_pcfs == PETSC_TRUE) {
919
920 auto bc_mng = mField.getInterface<BcManager>();
921 auto name_prb = simple->getProblemName();
922 auto is_all_bc = bc_mng->getBlockIS(name_prb, "FIX_X", "U", 0, 0);
923 is_all_bc = bc_mng->getBlockIS(name_prb, "FIX_Y", "U", 1, 1, is_all_bc);
924 is_all_bc = bc_mng->getBlockIS(name_prb, "FIX_Z", "U", 2, 2, is_all_bc);
925 is_all_bc = bc_mng->getBlockIS(name_prb, "FIX_ALL", "U", 0, 2, is_all_bc);
926
927 int is_all_bc_size;
928 CHKERR ISGetSize(is_all_bc, &is_all_bc_size);
929 MOFEM_LOG("EXAMPLE", Sev::inform)
930 << "Field split block size " << is_all_bc_size;
931
932 CHKERR PCFieldSplitSetIS(pc, PETSC_NULL,
933 is_all_bc); // boundary block
934 }
935
937 };
938
939 auto dm = simple->getDM();
940 auto D = smartCreateDMVector(dm);
941 uXScatter = scatter_create(D, 0);
942 uYScatter = scatter_create(D, 1);
943 if constexpr (SPACE_DIM == 3)
944 uZScatter = scatter_create(D, 2);
945
946 auto solver = pipeline_mng->createTSIM();
947
948 auto post_rhs = [&]() {
950 auto get_iter = [&]() {
951 SNES snes;
952 CHKERR TSGetSNES(solver, &snes);
953 int iter;
954 CHKERR SNESGetIterationNumber(snes, &iter);
955 return iter;
956 };
957
959 };
960
961 mField.getInterface<PipelineManager>()->getDomainRhsFE()->postProcessHook =
962 post_rhs;
963
964 CHKERR TSSetSolution(solver, D);
965 CHKERR set_section_monitor(solver);
966 CHKERR set_time_monitor(dm, solver);
967 CHKERR TSSetSolution(solver, D);
968 CHKERR TSSetFromOptions(solver);
969 CHKERR set_fieldsplit_preconditioner(solver);
970 CHKERR TSSetUp(solver);
971 CHKERR TSSolve(solver, NULL);
972
973 CHKERR VecGhostUpdateBegin(D, INSERT_VALUES, SCATTER_FORWARD);
974 CHKERR VecGhostUpdateEnd(D, INSERT_VALUES, SCATTER_FORWARD);
975 CHKERR DMoFEMMeshToLocalVector(dm, D, INSERT_VALUES, SCATTER_REVERSE);
976
978}
979//! [Solve]
980
981static char help[] = "...\n\n";
982
983int main(int argc, char *argv[]) {
984
985 // Initialisation of MoFEM/PETSc and MOAB data structures
986 const char param_file[] = "param_file.petsc";
988
989 // Add logging channel for example
990 auto core_log = logging::core::get();
991 core_log->add_sink(
993 LogManager::setLog("EXAMPLE");
994 MOFEM_LOG_TAG("EXAMPLE", "example");
995
996 try {
997
998 //! [Register MoFEM discrete manager in PETSc]
999 DMType dm_name = "DMMOFEM";
1000 CHKERR DMRegister_MoFEM(dm_name);
1001 //! [Register MoFEM discrete manager in PETSc
1002
1003 //! [Create MoAB]
1004 moab::Core mb_instance; ///< mesh database
1005 moab::Interface &moab = mb_instance; ///< mesh database interface
1006 //! [Create MoAB]
1007
1008 //! [Create MoFEM]
1009 MoFEM::Core core(moab); ///< finite element database
1010 MoFEM::Interface &m_field = core; ///< finite element database insterface
1011 //! [Create MoFEM]
1012
1013 //! [Load mesh]
1014 Simple *simple = m_field.getInterface<Simple>();
1015 CHKERR simple->getOptions();
1016 CHKERR simple->loadFile();
1017 //! [Load mesh]
1018
1019 //! [Example]
1020 Example ex(m_field);
1021 CHKERR ex.runProblem();
1022 //! [Example]
1023 }
1025
1027}
std::string param_file
constexpr double third
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
constexpr int SPACE_DIM
ElementsAndOps< SPACE_DIM >::DomainEle DomainEle
ElementsAndOps< SPACE_DIM >::BoundaryEle BoundaryEle
Kronecker Delta class symmetric.
@ ROW
Definition: definitions.h:123
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:372
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
@ NOBASE
Definition: definitions.h:59
@ DEMKOWICZ_JACOBI_BASE
Definition: definitions.h:66
#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
@ L2
field with C-1 continuity
Definition: definitions.h:88
@ H1
continuous field
Definition: definitions.h:85
@ NOSPACE
Definition: definitions.h:83
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_NOT_FOUND
Definition: definitions.h:33
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
static const char *const ApproximationBaseNames[]
Definition: definitions.h:72
#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
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
constexpr double shear_modulus_G
Definition: elastic.cpp:63
constexpr double bulk_modulus_K
Definition: elastic.cpp:62
const int dim
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::BiLinearForm< GAUSS >::OpMass< 1, SPACE_DIM > OpMass
ElementsAndOps< SPACE_DIM >::PostProcEle PostProcEle
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
IntegrationType
Form integrator integration types.
AssemblyType
[Storage and set boundary conditions]
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
Definition: LogManager.cpp:364
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:301
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:332
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.
double D
const double v
phase velocity of light in medium (cm/ns)
FTensor::Index< 'k', 3 > k
FormsIntegrators< DomainEleOp >::Assembly< USER_ASSEMBLE >::LinearForm< GAUSS >::OpBaseTimesVector< 1, 3, 1 > OpInertiaForce
FormsIntegrators< DomainEleOp >::Assembly< USER_ASSEMBLE >::LinearForm< GAUSS >::OpSource< 1, 3 > OpBodyForce
const FTensor::Tensor2< T, Dim, Dim > Vec
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
auto type_from_handle(const EntityHandle h)
get type from entity handle
Definition: Templates.hpp:1592
PetscErrorCode DMMoFEMTSSetMonitor(DM dm, TS ts, const std::string fe_name, boost::shared_ptr< MoFEM::FEMethod > method, boost::shared_ptr< MoFEM::BasicMethod > pre_only, boost::shared_ptr< MoFEM::BasicMethod > post_only)
Set Monitor To TS solver.
Definition: DMMMoFEM.cpp:1003
PetscErrorCode PetscOptionsGetInt(PetscOptions *, const char pre[], const char name[], PetscInt *ivalue, PetscBool *set)
MoFEMErrorCode MoFEMSNESMonitorFields(SNES snes, PetscInt its, PetscReal fgnorm, SnesCtx *snes_ctx)
Sens monitor printing residual field by field.
Definition: SnesCtx.cpp:224
PetscErrorCode PetscOptionsGetBool(PetscOptions *, const char pre[], const char name[], PetscBool *bval, PetscBool *set)
PetscErrorCode PetscOptionsGetScalar(PetscOptions *, const char pre[], const char name[], PetscScalar *dval, PetscBool *set)
OpCalculateScalarFieldValuesFromPetscVecImpl< PetscData::CTX_SET_X_T > OpCalculateScalarFieldValuesDot
auto smartGetDMSnesCtx(DM dm)
Get SNES context data structure used by DM.
Definition: DMMoFEM.hpp:987
OpPostProcMapInMoab< SPACE_DIM, SPACE_DIM > OpPPMap
double young_modulus
Definition: plastic.cpp:107
double Qinf
Definition: plastic.cpp:116
static char help[]
[Solve]
Definition: plastic.cpp:981
double rho
Definition: plastic.cpp:109
#define EXECUTABLE_DIMENSION
Definition: plastic.cpp:10
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< G >::OpGradTimesTensor< 1, SPACE_DIM, SPACE_DIM > OpInternalForcePiola
Definition: plastic.cpp:77
FormsIntegrators< DomainEleOp >::Assembly< A >::BiLinearForm< GAUSS >::OpGradSymTensorGrad< 1, SPACE_DIM, SPACE_DIM, 0 > OpKCauchy
[Only used with Hooke equation (linear material model)]
Definition: plastic.cpp:61
constexpr int SPACE_DIM
Definition: plastic.cpp:42
double visH
Definition: plastic.cpp:112
double poisson_ratio
Definition: plastic.cpp:108
FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::BiLinearForm< G >::OpMass< 1, SPACE_DIM > OpBoundaryMass
[Only used with Hencky/nonlinear material]
Definition: plastic.cpp:82
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< G >::OpGradTimesSymTensor< 1, SPACE_DIM, SPACE_DIM > OpInternalForceCauchy
Definition: plastic.cpp:63
FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::LinearForm< G >::OpBaseTimesVector< 1, SPACE_DIM, 0 > OpBoundaryVec
Definition: plastic.cpp:84
double scale
Definition: plastic.cpp:105
constexpr auto size_symm
Definition: plastic.cpp:44
long double hardening_dtau2(long double tau)
Definition: plastic.cpp:130
double zeta
Definition: plastic.cpp:114
double H
Definition: plastic.cpp:111
long double hardening(long double tau)
Definition: plastic.cpp:122
constexpr AssemblyType A
Definition: plastic.cpp:46
int order
Order if fixed.
Definition: plastic.cpp:119
double b_iso
Definition: plastic.cpp:117
PetscBool is_large_strains
Definition: plastic.cpp:103
int geom_order
Order if fixed.
Definition: plastic.cpp:120
FormsIntegrators< DomainEleOp >::Assembly< A >::BiLinearForm< GAUSS >::OpGradTensorGrad< 1, SPACE_DIM, SPACE_DIM, 1 > OpKPiola
[Only used for dynamics]
Definition: plastic.cpp:75
double sigmaY
Definition: plastic.cpp:110
constexpr IntegrationType G
Definition: plastic.cpp:47
long double hardening_dtau(long double tau)
Definition: plastic.cpp:126
double eta
Definition: plastic.cpp:115
double cn
Definition: plastic.cpp:113
FormsIntegrators< BoundaryEleOp >::Assembly< PETSC >::LinearForm< G >::OpBaseTimesVector< 1, SPACE_DIM, 1 > OpBoundaryInternal
Definition: plastic.cpp:86
OpBaseImpl< PETSC, EdgeEleOp > OpBase
Definition: radiation.cpp:29
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'i', 3 > i
double getScale(const double time)
Get scaling at a given time.
Definition: plastic.cpp:167
[Example]
Definition: plastic.cpp:139
MoFEMErrorCode tsSolve()
[Push operators to pipeline]
Definition: plastic.cpp:747
FieldApproximationBase base
Definition: plot_base.cpp:68
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uYScatter
Definition: plastic.cpp:159
MoFEMErrorCode createCommonData()
[Set up problem]
Definition: plastic.cpp:259
Example(MoFEM::Interface &m_field)
Definition: plastic.cpp:141
MoFEMErrorCode OPs()
[Boundary condition]
Definition: plastic.cpp:455
MoFEMErrorCode runProblem()
[Run problem]
Definition: plastic.cpp:174
MoFEM::Interface & mField
Definition: plastic.cpp:146
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uZScatter
Definition: plastic.cpp:160
boost::shared_ptr< std::vector< unsigned char > > boundaryMarker
Definition: plastic.cpp:162
MoFEMErrorCode setupProblem()
[Run problem]
Definition: plastic.cpp:186
MoFEMErrorCode bC()
[Create common data]
Definition: plastic.cpp:393
boost::shared_ptr< std::vector< unsigned char > > reactionMarker
Definition: plastic.cpp:163
boost::shared_ptr< DomainEle > reactionFe
Definition: plastic.cpp:156
boost::shared_ptr< PlasticOps::CommonData > commonPlasticDataPtr
Definition: plastic.cpp:154
std::tuple< SmartPetscObj< Vec >, SmartPetscObj< VecScatter > > uXScatter
Definition: plastic.cpp:158
boost::shared_ptr< HenckyOps::CommonData > commonHenckyDataPtr
Definition: plastic.cpp:155
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)
Natural boundary conditions.
Definition: Essential.hpp:101
Section manager is used to create indexes and sections.
Definition: ISManager.hpp:23
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
Definition: LogManager.cpp:279
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:323
Assembly methods.
Definition: Natural.hpp:57
Get value at integration points for scalar field.
Calculate symmetric tensor field rates ant integratio pts.
Calculate symmetric tensor field values at integration pts.
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.
Enforce essential constrains on lhs.
Definition: Essential.hpp:71
Enforce essential constrains on rhs.
Definition: Essential.hpp:55
Post post-proc data at points from hash maps.
Scale base functions by inverses of measure of element.
Set indices on entities on finite element.
PipelineManager interface.
MoFEM::VolumeElementForcesAndSourcesCore VolEle
MoFEM::FaceElementForcesAndSourcesCore FaceEle
MoFEM::EdgeElementForcesAndSourcesCore EdgeEle
Problem manager is used to build and partition problems.
Projection of edge entities with one mid-node on hierarchical basis.
Simple interface for fast problem set-up.
Definition: Simple.hpp:26
intrusive_ptr for managing petsc objects
Force scale operator for reading two columns.
double getScale(const double time)
Get scaling at a given time.
TimeScale(std::string file_name="", bool error_if_file_not_given=false)
TimeScale constructor.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.