v0.15.0
Loading...
Searching...
No Matches
ADOLCPlasticity::OpCalculateStress< DIM, LARGE_STRAIN > Struct Template Reference
Inheritance diagram for ADOLCPlasticity::OpCalculateStress< DIM, LARGE_STRAIN >:
[legend]
Collaboration diagram for ADOLCPlasticity::OpCalculateStress< DIM, LARGE_STRAIN >:
[legend]

Public Member Functions

MoFEMErrorCode doWork (int side, EntityType type, DataForcesAndSourcesCore::EntData &data) override
 
- Public Member Functions inherited from ADOLCPlasticity::OpCalculateStressImpl< DIM, LARGE_STRAIN >
 OpCalculateStressImpl (MoFEM::Interface &m_field, boost::shared_ptr< CommonData > common_data_ptr, boost::shared_ptr< ClosestPointProjection > cp_ptr, bool calc_lhs)
 

Additional Inherited Members

- Protected Member Functions inherited from ADOLCPlasticity::OpCalculateStressImpl< DIM, LARGE_STRAIN >
MoFEMErrorCode getTags ()
 
MoFEMErrorCode setTagsData (const EntityHandle tet, const int nb_gauss_pts, const int nb_internal_variables)
 
- Protected Attributes inherited from ADOLCPlasticity::OpCalculateStressImpl< DIM, LARGE_STRAIN >
MoFEM::InterfacemField
 
boost::shared_ptr< CommonDatacommonDataPtr
 
boost::shared_ptr< ClosestPointProjectioncpPtr
 
bool calcLhs
 
Tag thPlasticStrain
 
Tag thInternalVariables
 

Detailed Description

template<int DIM>
struct ADOLCPlasticity::OpCalculateStress< DIM, LARGE_STRAIN >

Definition at line 283 of file ADOLCPlasticity.cpp.

Member Function Documentation

◆ doWork()

template<int DIM>
MoFEMErrorCode ADOLCPlasticity::OpCalculateStress< DIM, LARGE_STRAIN >::doWork ( int side,
EntityType type,
DataForcesAndSourcesCore::EntData & data )
override

Definition at line 376 of file ADOLCPlasticity.cpp.

377 {
379
380 auto do_work_impl = [this](auto commonDataPtr, auto cpPtr) {
382
383 int nb_gauss_pts = this->getGaussPts().size2();
384
385 int nb_internal_variables = cpPtr->nbInternalVariables;
386 auto tet = this->getNumeredEntFiniteElementPtr()->getEnt();
387 CHKERR this->setTagsData(tet, nb_gauss_pts, nb_internal_variables);
388
389 commonDataPtr->activeVariablesW.resize(
390 nb_gauss_pts, 12 + cpPtr->nbInternalVariables, false);
391 commonDataPtr->gradientW.resize(nb_gauss_pts,
392 12 + cpPtr->nbInternalVariables, false);
393 commonDataPtr->materialTangentOperator.resize(36, nb_gauss_pts, false);
394 commonDataPtr->deltaGamma.resize(nb_gauss_pts);
395
396 FTensor::Index<'i', 3> i;
397 FTensor::Index<'j', 3> j;
398 FTensor::Index<'k', 3> k;
399 FTensor::Index<'l', 3> l;
400 FTensor::Index<'Z', 6> Z;
401 FTensor::Index<'Y', 6> Y;
402
403 // Code uses trial stress, until first solve of linear system of equations
404 int iter = 1;
405 if (this->getFEMethod()->snes_ctx != SnesMethod::CTX_SNESNONE) {
406 CHKERR SNESGetLinearSolveIterations(this->getFEMethod()->snes, &iter);
407 }
408
409 double ave_tr_strain = 0;
410
411 auto t_voight_stress_op = voight_to_stress_op();
412
413 auto t_grad =
415
416 auto t_Cep =
417 getFTensor4DdgFromMat<3, 3, 1>(commonDataPtr->materialTangentOperator);
418 for (auto gg = 0; gg != nb_gauss_pts; ++gg) {
419
420 auto tmp_active =
421 getVectorAdaptor(&(commonDataPtr->activeVariablesW(gg, 0)),
422 12 + cpPtr->nbInternalVariables);
423 cpPtr->activeVariablesW.swap(tmp_active);
424 auto tmp_gradient = getVectorAdaptor(&(commonDataPtr->gradientW(gg, 0)),
425 12 + cpPtr->nbInternalVariables);
426 this->cpPtr->gradientW.swap(tmp_gradient);
427
428 auto t_voigt_strain =
429 calcStrain(commonDataPtr->bBar, ave_tr_strain, t_grad,
430 getFTensor1FromPtr<6>(&(cpPtr->activeVariablesW[6])));
431 ++t_grad;
432
433 // Copy plastic strains and internal variables from tags
434 auto copy_plastic_strain_and_internal_variables = [&]() {
436
437 // Get data stored on Tags
438 VectorAdaptor plastic_strain_tags =
439 VectorAdaptor(6, ublas::shallow_array_adaptor<double>(
440 6, &commonDataPtr->plasticStrainPtr[gg * 6]));
441 VectorAdaptor internal_variables_tags = VectorAdaptor(
442 nb_internal_variables,
443 ublas::shallow_array_adaptor<double>(
444 nb_internal_variables,
446 ->internalVariablesPtr[gg * nb_internal_variables]));
447
448 // Set values to closest point projection
449 cpPtr->plasticStrain0.resize(6, false);
450 noalias(cpPtr->plasticStrain0) = plastic_strain_tags;
451 cpPtr->internalVariables0.resize(nb_internal_variables, false);
452 noalias(cpPtr->internalVariables0) = internal_variables_tags;
453
454 auto cp_plastic_strain = cpPtr->getPlasticStrain();
455 auto cp_internal_variables = cpPtr->getInternalVariables();
456 noalias(cp_plastic_strain) = plastic_strain_tags;
457 noalias(cp_internal_variables) = internal_variables_tags;
459 };
460
461 CHKERR copy_plastic_strain_and_internal_variables();
462
463 auto closest_point_projection = [&](auto &recalculate_elastic_tangent) {
467 CHKERR cpPtr->setParams(t, recalculate_elastic_tangent);
468 }
469 CHKERR cpPtr->playW_NoHessian();
470 CHKERR cpPtr->playY_NoGradient();
471 double y = cpPtr->y;
472 cpPtr->deltaGamma = 0; // zero lagrange multiplier
473 if (iter > 0 && y > std::numeric_limits<double>::epsilon()) {
474 CHKERR cpPtr->solveClosestProjection();
475 }
476 commonDataPtr->deltaGamma[gg] = cpPtr->deltaGamma;
478 };
479
480 auto evaluate_consistent_tangent_matrix =
481 [&](auto &recalculate_elastic_tangent) {
483 if (recalculate_elastic_tangent)
484 CHKERR cpPtr->playW();
485
486 if (iter > 0 &&
487 cpPtr->deltaGamma > std::numeric_limits<double>::epsilon()) {
488
489 CHKERR cpPtr->consistentTangent();
490
491 auto &m = cpPtr->Cep;
492 // for generic case of non-associated plasticity consistent
493 // tangent matrix is non-symmetric
494 auto t_voight_cep =
496 &m(0, 0), &m(0, 1), &m(0, 2), &m(0, 3), &m(0, 4),
497 &m(0, 5),
498
499 &m(1, 0), &m(1, 1), &m(1, 2), &m(1, 3), &m(1, 4),
500 &m(1, 5),
501
502 &m(2, 0), &m(2, 1), &m(2, 2), &m(2, 3), &m(2, 4),
503 &m(2, 5),
504
505 &m(3, 0), &m(3, 1), &m(3, 2), &m(3, 3), &m(3, 4),
506 &m(3, 5),
507
508 &m(4, 0), &m(4, 1), &m(4, 2), &m(4, 3), &m(4, 4),
509 &m(4, 5),
510
511 &m(5, 0), &m(5, 1), &m(5, 2), &m(5, 3), &m(5, 4), &m(5, 5)
512
513 );
514
515 t_Cep(i, j, k, l) =
516 t_voight_stress_op(i, j, Z) *
517 (t_voight_cep(Z, Y) * t_voight_stress_op(k, l, Y));
518
519 } else {
520
521 auto &m = cpPtr->C;
522 auto t_voight_cep =
524 &m(0, 0), &m(0, 1), &m(0, 2), &m(0, 3), &m(0, 4),
525 &m(0, 5),
526
527 &m(1, 1), &m(1, 2), &m(1, 3), &m(1, 4), &m(1, 5),
528
529 &m(2, 2), &m(2, 3), &m(2, 4), &m(2, 5),
530
531 &m(3, 3), &m(3, 4), &m(3, 5),
532
533 &m(4, 4), &m(4, 5),
534
535 &m(5, 5));
536
537 t_Cep(i, j, k, l) =
538 t_voight_stress_op(i, j, Z) *
539 (t_voight_cep(Z, Y) * t_voight_stress_op(k, l, Y));
540 }
541
542 ++t_Cep;
544 };
545
546 // Always recalculate elastic tangent if first element
547 bool recalculate_elastic_tangent =
548 (this->getNinTheLoop() == 0) ? true : false;
549 CHKERR closest_point_projection(recalculate_elastic_tangent);
550 // always calculate consistent tangent matrix for large strain
551 CHKERR evaluate_consistent_tangent_matrix(recalculate_elastic_tangent);
552 }
553
554 auto calcs_stress_matrix = [&]() {
556 auto t_voight_stress_op = voight_to_stress_op();
557 commonDataPtr->stressMatrix.resize(
558 6, this->commonDataPtr->gradientW.size1(), false);
559 auto t_stess =
561 auto t_voight_stress = commonDataPtr->getFTensor1StressAtGaussPts();
562 for (auto gg = 0; gg != commonDataPtr->gradientW.size1(); ++gg) {
563 t_stess(i, j) = t_voight_stress_op(i, j, Z) * t_voight_stress(Z);
564 ++t_voight_stress;
565 ++t_stess;
566 }
568 };
569
570 auto calcs_plastic_strain_matrix = [&]() {
572 auto t_voight_strain_op = voight_to_stress_op();
573 commonDataPtr->plasticStrainMatrix.resize(
574 6, commonDataPtr->activeVariablesW.size1(), false);
575 auto t_plastic_strain =
576 getFTensor2SymmetricFromMat<3>(commonDataPtr->plasticStrainMatrix);
577 auto t_voight_plastic_strain =
578 commonDataPtr->getFTensor1PlasticStrainAtGaussPts();
579 for (auto gg = 0; gg != commonDataPtr->activeVariablesW.size1(); ++gg) {
580 t_plastic_strain(i, j) =
581 t_voight_strain_op(i, j, Z) * t_voight_plastic_strain(Z);
582 ++t_voight_plastic_strain;
583 ++t_plastic_strain;
584 }
586 };
587
588 CHKERR calcs_stress_matrix();
589 CHKERR calcs_plastic_strain_matrix();
591 };
592
593 CHKERR do_work_impl(this->commonDataPtr, this->cpPtr);
594
596}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'l', 3 > l
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
FTensor::Dg< double, 3, 6 > voight_to_stress_op()
Op convert Vight stress vector to stress tensor.
FTensor::Tensor1< T2, DIM3 > calcStrain(bool b_bar, double ave_tr_strain, FTensor::Tensor2< T1, DIM1, DIM2 > &t_grad, FTensor::Tensor1< T2, DIM3 > &&t_voight_strain)
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition Types.hpp:115
auto getVectorAdaptor(T1 ptr, const size_t n)
Get Vector adaptor.
Definition Templates.hpp:31
static FTensor::Ddg< FTensor::PackPtr< T *, 1 >, Tensor_Dim01, Tensor_Dim23 > getFTensor4DdgFromMat(ublas::matrix< T, L, A > &data)
Get symmetric tensor rank 4 on first two and last indices from form data matrix.
static auto getFTensor2SymmetricFromMat(ublas::matrix< T, L, A > &data)
Get symmetric tensor rank 2 (matrix) form data matrix.
FTensor::Tensor1< FTensor::PackPtr< double *, S >, DIM > getFTensor1FromPtr(double *ptr)
Make Tensor1 from pointer.
constexpr double t
plate stiffness
Definition plate.cpp:58
FTensor::Index< 'm', 3 > m
MoFEMErrorCode setTagsData(const EntityHandle tet, const int nb_gauss_pts, const int nb_internal_variables)
boost::shared_ptr< ClosestPointProjection > cpPtr

The documentation for this struct was generated from the following file: