v0.14.0
Public Member Functions | Protected Attributes | List of all members
PlasticOps::OpCalculatePlasticityImpl< DIM, GAUSS, DomainEleOp > Struct Template Reference

#include <tutorials/adv-0/src/PlasticOpsGeneric.hpp>

Inheritance diagram for PlasticOps::OpCalculatePlasticityImpl< DIM, GAUSS, DomainEleOp >:
[legend]
Collaboration diagram for PlasticOps::OpCalculatePlasticityImpl< DIM, GAUSS, DomainEleOp >:
[legend]

Public Member Functions

 OpCalculatePlasticityImpl (const std::string field_name, boost::shared_ptr< CommonData > common_data_ptr, boost::shared_ptr< MatrixDouble > m_D_ptr)
 
MoFEMErrorCode doWork (int side, EntityType type, EntData &data)
 

Protected Attributes

boost::shared_ptr< CommonData > commonDataPtr
 
boost::shared_ptr< MatrixDouble > mDPtr
 

Detailed Description

template<int DIM, typename DomainEleOp>
struct PlasticOps::OpCalculatePlasticityImpl< DIM, GAUSS, DomainEleOp >

Definition at line 469 of file PlasticOpsGeneric.hpp.

Constructor & Destructor Documentation

◆ OpCalculatePlasticityImpl()

template<int DIM, typename DomainEleOp >
PlasticOps::OpCalculatePlasticityImpl< DIM, GAUSS, DomainEleOp >::OpCalculatePlasticityImpl ( const std::string  field_name,
boost::shared_ptr< CommonData >  common_data_ptr,
boost::shared_ptr< MatrixDouble >  m_D_ptr 
)
Examples
PlasticOpsGeneric.hpp.

Definition at line 481 of file PlasticOpsGeneric.hpp.

484  : DomainEleOp(field_name, DomainEleOp::OPROW),
485  commonDataPtr(common_data_ptr), mDPtr(m_D_ptr) {
486  // Opetor is only executed for vertices
487  std::fill(&DomainEleOp::doEntities[MBEDGE],
488  &DomainEleOp::doEntities[MBMAXTYPE], false);
489 }

Member Function Documentation

◆ doWork()

template<int DIM, typename DomainEleOp >
MoFEMErrorCode PlasticOps::OpCalculatePlasticityImpl< DIM, GAUSS, DomainEleOp >::doWork ( int  side,
EntityType  type,
EntData &  data 
)

< material parameters

Examples
PlasticOpsGeneric.hpp.

Definition at line 492 of file PlasticOpsGeneric.hpp.

493  {
494  MoFEMFunctionBegin;
495 
496  FTensor::Index<'i', DIM> i;
497  FTensor::Index<'j', DIM> j;
498  FTensor::Index<'k', DIM> k;
499  FTensor::Index<'l', DIM> l;
500  FTensor::Index<'m', DIM> m;
501  FTensor::Index<'n', DIM> n;
502 
503  auto &params = commonDataPtr->blockParams; ///< material parameters
504 
505  auto nb_gauss_pts = DomainEleOp::getGaussPts().size2();
506  auto t_w = DomainEleOp::getFTensor0IntegrationWeight();
507  auto t_tau = getFTensor0FromVec(commonDataPtr->plasticTau);
508  auto t_tau_dot = getFTensor0FromVec(commonDataPtr->plasticTauDot);
509  auto t_f = getFTensor0FromVec(commonDataPtr->plasticSurface);
510  auto t_flow = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->plasticFlow);
511  auto t_plastic_strain =
512  getFTensor2SymmetricFromMat<DIM>(commonDataPtr->plasticStrain);
513  auto t_plastic_strain_dot =
514  getFTensor2SymmetricFromMat<DIM>(commonDataPtr->plasticStrainDot);
515  auto t_stress =
516  getFTensor2SymmetricFromMat<DIM>(*(commonDataPtr->mStressPtr));
517 
518  auto t_D = getFTensor4DdgFromMat<DIM, DIM, 0>(*commonDataPtr->mDPtr);
519  auto t_D_Op = getFTensor4DdgFromMat<DIM, DIM, 0>(*mDPtr);
520 
521  auto t_diff_plastic_strain = diff_tensor(FTensor::Number<DIM>());
522  auto t_diff_deviator = diff_deviator(diff_tensor(FTensor::Number<DIM>()));
523 
524  FTensor::Ddg<double, DIM, DIM> t_flow_dir_dstress;
525  FTensor::Ddg<double, DIM, DIM> t_flow_dir_dstrain;
526  t_flow_dir_dstress(i, j, k, l) =
527  1.5 * (t_diff_deviator(M, N, i, j) * t_diff_deviator(M, N, k, l));
528  t_flow_dir_dstrain(i, j, k, l) =
529  t_flow_dir_dstress(i, j, m, n) * t_D_Op(m, n, k, l);
530 
531 
532  auto t_alpha_dir =
533  kinematic_hardening_dplastic_strain<DIM>(params[CommonData::C1_k]);
534 
535  commonDataPtr->resC.resize(nb_gauss_pts, false);
536  commonDataPtr->resCdTau.resize(nb_gauss_pts, false);
537  commonDataPtr->resCdStrain.resize(size_symm, nb_gauss_pts, false);
538  commonDataPtr->resCdPlasticStrain.resize(size_symm, nb_gauss_pts, false);
539  commonDataPtr->resFlow.resize(size_symm, nb_gauss_pts, false);
540  commonDataPtr->resFlowDtau.resize(size_symm, nb_gauss_pts, false);
541  commonDataPtr->resFlowDstrain.resize(size_symm * size_symm, nb_gauss_pts,
542  false);
543  commonDataPtr->resFlowDstrainDot.resize(size_symm * size_symm, nb_gauss_pts,
544  false);
545 
546  commonDataPtr->resC.clear();
547  commonDataPtr->resCdTau.clear();
548  commonDataPtr->resCdStrain.clear();
549  commonDataPtr->resCdPlasticStrain.clear();
550  commonDataPtr->resFlow.clear();
551  commonDataPtr->resFlowDtau.clear();
552  commonDataPtr->resFlowDstrain.clear();
553  commonDataPtr->resFlowDstrainDot.clear();
554 
555  auto t_res_c = getFTensor0FromVec(commonDataPtr->resC);
556  auto t_res_c_dtau = getFTensor0FromVec(commonDataPtr->resCdTau);
557  auto t_res_c_dstrain =
558  getFTensor2SymmetricFromMat<DIM>(commonDataPtr->resCdStrain);
559  auto t_res_c_plastic_strain =
560  getFTensor2SymmetricFromMat<DIM>(commonDataPtr->resCdPlasticStrain);
561  auto t_res_flow = getFTensor2SymmetricFromMat<DIM>(commonDataPtr->resFlow);
562  auto t_res_flow_dtau =
563  getFTensor2SymmetricFromMat<DIM>(commonDataPtr->resFlowDtau);
564  auto t_res_flow_dstrain =
565  getFTensor4DdgFromMat<DIM, DIM>(commonDataPtr->resFlowDstrain);
566  auto t_res_flow_dplastic_strain =
567  getFTensor4DdgFromMat<DIM, DIM>(commonDataPtr->resFlowDstrainDot);
568 
569  auto next = [&]() {
570  ++t_tau;
571  ++t_tau_dot;
572  ++t_f;
573  ++t_flow;
574  ++t_plastic_strain;
575  ++t_plastic_strain_dot;
576  ++t_stress;
577  ++t_res_c;
578  ++t_res_c_dtau;
579  ++t_res_c_dstrain;
580  ++t_res_c_plastic_strain;
581  ++t_res_flow;
582  ++t_res_flow_dtau;
583  ++t_res_flow_dstrain;
584  ++t_res_flow_dplastic_strain;
585  ++t_w;
586  };
587 
588  auto get_avtive_pts = [&]() {
589  int nb_points_avtive_on_elem = 0;
590  int nb_points_on_elem = 0;
591 
592  auto t_tau = getFTensor0FromVec(commonDataPtr->plasticTau);
593  auto t_tau_dot = getFTensor0FromVec(commonDataPtr->plasticTauDot);
594  auto t_f = getFTensor0FromVec(commonDataPtr->plasticSurface);
595  auto t_plastic_strain_dot =
596  getFTensor2SymmetricFromMat<SPACE_DIM>(commonDataPtr->plasticStrainDot);
597 
598  auto dt = this->getTStimeStep();
599 
600  for (auto gg = 0; gg != nb_gauss_pts; ++gg) {
601  auto eqiv = equivalent_strain_dot(t_plastic_strain_dot);
602  const auto ww = w(
603  eqiv, t_tau_dot, t_f,
604  iso_hardening(t_tau, params[CommonData::H], params[CommonData::QINF],
605  params[CommonData::BISO], params[CommonData::SIGMA_Y]),
606  params[CommonData::SIGMA_Y]);
607  const auto sign_ww = constrian_sign(ww, dt);
608 
609  ++nb_points_on_elem;
610  if (sign_ww > 0) {
611  ++nb_points_avtive_on_elem;
612  }
613 
614  ++t_tau;
615  ++t_tau_dot;
616  ++t_f;
617  ++t_plastic_strain_dot;
618  }
619 
620  int &active_points = PlasticOps::CommonData::activityData[0];
621  int &avtive_full_elems = PlasticOps::CommonData::activityData[1];
622  int &avtive_elems = PlasticOps::CommonData::activityData[2];
623  int &nb_points = PlasticOps::CommonData::activityData[3];
624  int &nb_elements = PlasticOps::CommonData::activityData[4];
625 
626  ++nb_elements;
627  nb_points += nb_points_on_elem;
628  if (nb_points_avtive_on_elem > 0) {
629  ++avtive_elems;
630  active_points += nb_points_avtive_on_elem;
631  if (nb_points_avtive_on_elem == nb_points_on_elem) {
632  ++avtive_full_elems;
633  }
634  }
635 
636  if (nb_points_avtive_on_elem != nb_points_on_elem)
637  return 1;
638  else
639  return 0;
640  };
641 
642  if (DomainEleOp::getTSCtx() == TSMethod::TSContext::CTX_TSSETIJACOBIAN) {
643  get_avtive_pts();
644  }
645 
646  auto dt = this->getTStimeStep();
647  for (auto gg = 0; gg != nb_gauss_pts; ++gg) {
648 
649  auto eqiv = equivalent_strain_dot(t_plastic_strain_dot);
650  auto t_diff_eqiv = diff_equivalent_strain_dot(eqiv, t_plastic_strain_dot,
651  t_diff_plastic_strain,
652  FTensor::Number<DIM>());
653 
654  const auto sigma_y =
655  iso_hardening(t_tau, params[CommonData::H], params[CommonData::QINF],
656  params[CommonData::BISO], params[CommonData::SIGMA_Y]);
657  const auto d_sigma_y =
658  iso_hardening_dtau(t_tau, params[CommonData::H],
659  params[CommonData::QINF], params[CommonData::BISO]);
660 
661  auto ww = w(eqiv, t_tau_dot, t_f, sigma_y, params[CommonData::SIGMA_Y]);
662  auto abs_ww = constrain_abs(ww, dt);
663  auto sign_ww = constrian_sign(ww, dt);
664 
665  auto c = constraint(eqiv, t_tau_dot, t_f, sigma_y, abs_ww,
666  params[CommonData::VIS_H], params[CommonData::SIGMA_Y]);
667  auto c_dot_tau = diff_constrain_ddot_tau(sign_ww, eqiv, t_tau_dot,
668  params[CommonData::VIS_H],
669  params[CommonData::SIGMA_Y]);
670  auto c_equiv = diff_constrain_deqiv(sign_ww, eqiv, t_tau_dot,
671  params[CommonData::SIGMA_Y]);
672  auto c_sigma_y = diff_constrain_dsigma_y(sign_ww);
673  auto c_f = diff_constrain_df(sign_ww);
674 
675  auto t_dev_stress = deviator(
676 
677  t_stress, trace(t_stress),
678 
679  kinematic_hardening(t_plastic_strain, params[CommonData::C1_k])
680 
681  );
682 
683  FTensor::Tensor2_symmetric<double, DIM> t_flow_dir;
684  t_flow_dir(k, l) = 1.5 * (t_dev_stress(I, J) * t_diff_deviator(I, J, k, l));
685  FTensor::Tensor2_symmetric<double, DIM> t_flow_dstrain;
686  t_flow_dstrain(i, j) = t_flow(k, l) * t_D_Op(k, l, i, j);
687 
688  auto get_res_c = [&]() { return c; };
689 
690  auto get_res_c_dstrain = [&](auto &t_diff_res) {
691  t_diff_res(i, j) = c_f * t_flow_dstrain(i, j);
692  };
693 
694  auto get_res_c_dplastic_strain = [&](auto &t_diff_res) {
695  t_diff_res(i, j) = (this->getTSa() * c_equiv) * t_diff_eqiv(i, j);
696  t_diff_res(k, l) -= c_f * t_flow(i, j) * t_alpha_dir(i, j, k, l);
697  };
698 
699  auto get_res_c_dtau = [&]() {
700  return this->getTSa() * c_dot_tau + c_sigma_y * d_sigma_y;
701  };
702 
703  auto get_res_c_plastic_strain = [&](auto &t_diff_res) {
704  t_diff_res(k, l) = -c_f * t_flow(i, j) * t_alpha_dir(i, j, k, l);
705  };
706 
707  auto get_res_flow = [&](auto &t_res_flow) {
708  const auto a = sigma_y;
709  const auto b = t_tau_dot;
710  t_res_flow(k, l) = a * t_plastic_strain_dot(k, l) - b * t_flow_dir(k, l);
711  };
712 
713  auto get_res_flow_dtau = [&](auto &t_res_flow_dtau) {
714  const auto da = d_sigma_y;
715  const auto db = this->getTSa();
716  t_res_flow_dtau(k, l) =
717  da * t_plastic_strain_dot(k, l) - db * t_flow_dir(k, l);
718  };
719 
720  auto get_res_flow_dstrain = [&](auto &t_res_flow_dstrain) {
721  const auto b = t_tau_dot;
722  t_res_flow_dstrain(m, n, k, l) = -t_flow_dir_dstrain(m, n, k, l) * b;
723  };
724 
725  auto get_res_flow_dplastic_strain = [&](auto &t_res_flow_dplastic_strain) {
726  const auto a = sigma_y;
727  t_res_flow_dplastic_strain(m, n, k, l) =
728  (a * this->getTSa()) * t_diff_plastic_strain(m, n, k, l);
729  const auto b = t_tau_dot;
730  t_res_flow_dplastic_strain(m, n, i, j) +=
731  (t_flow_dir_dstrain(m, n, k, l) * t_alpha_dir(k, l, i, j)) * b;
732  };
733 
734  t_res_c = get_res_c();
735  get_res_flow(t_res_flow);
736 
737  if (this->getTSCtx() == TSMethod::TSContext::CTX_TSSETIJACOBIAN) {
738  t_res_c_dtau = get_res_c_dtau();
739  get_res_c_dstrain(t_res_c_dstrain);
740  get_res_c_dplastic_strain(t_res_c_plastic_strain);
741  get_res_flow_dtau(t_res_flow_dtau);
742  get_res_flow_dstrain(t_res_flow_dstrain);
743  get_res_flow_dplastic_strain(t_res_flow_dplastic_strain);
744  }
745 
746  next();
747  }
748 
749  MoFEMFunctionReturn(0);
750 }

Member Data Documentation

◆ commonDataPtr

template<int DIM, typename DomainEleOp >
boost::shared_ptr<CommonData> PlasticOps::OpCalculatePlasticityImpl< DIM, GAUSS, DomainEleOp >::commonDataPtr
protected

Definition at line 476 of file PlasticOpsGeneric.hpp.

◆ mDPtr

template<int DIM, typename DomainEleOp >
boost::shared_ptr<MatrixDouble> PlasticOps::OpCalculatePlasticityImpl< DIM, GAUSS, DomainEleOp >::mDPtr
protected

Definition at line 477 of file PlasticOpsGeneric.hpp.

The documentation for this struct was generated from the following file:
PlasticOps::CommonData::BISO
@ BISO
Definition: PlasticOps.hpp:52
PlasticOps::OpCalculatePlasticityImpl< DIM, GAUSS, DomainEleOp >::commonDataPtr
boost::shared_ptr< CommonData > commonDataPtr
Definition: PlasticOpsGeneric.hpp:476
PlasticOps::w
double w(double eqiv, double dot_tau, double f, double sigma_y, double sigma_Y)
Definition: PlasticOpsGeneric.hpp:276
PlasticOps::J
FTensor::Index< 'J', 3 > J
Definition: PlasticOps.hpp:116
PlasticOps::diff_constrain_ddot_tau
double diff_constrain_ddot_tau(double sign, double eqiv, double dot_tau, double vis_H, double sigma_Y)
Definition: PlasticOpsGeneric.hpp:312
PlasticOps::diff_deviator
auto diff_deviator(FTensor::Ddg< double, DIM, DIM > &&t_diff_stress, FTensor::Number< DIM >)
Definition: PlasticOpsGeneric.hpp:125
PlasticOps::diff_constrain_df
auto diff_constrain_df(double sign)
Definition: PlasticOpsGeneric.hpp:322
PlasticOps::diff_tensor
auto diff_tensor(FTensor::Number< DIM >)
[Lambda functions]
Definition: PlasticOpsGeneric.hpp:10
PlasticOps::CommonData::C1_k
@ C1_k
Definition: PlasticOps.hpp:53
PlasticOps::equivalent_strain_dot
auto equivalent_strain_dot(FTensor::Tensor2_symmetric< T, DIM > &t_plastic_strain_dot)
Definition: PlasticOpsGeneric.hpp:350
PlasticOps::diff_equivalent_strain_dot
auto diff_equivalent_strain_dot(const T1 eqiv, T2 &t_plastic_strain_dot, T3 &t_diff_plastic_strain, FTensor::Number< DIM >)
Definition: PlasticOpsGeneric.hpp:361
PlasticOps::constrian_sign
double constrian_sign(double x, double dt)
Definition: PlasticOpsGeneric.hpp:251
iso_hardening_dtau
double iso_hardening_dtau(double tau, double H, double Qinf, double b_iso)
Definition: plastic.cpp:77
FTensor::Tensor2_symmetric
Definition: Tensor2_symmetric_value.hpp:13
PlasticOps::N
FTensor::Index< 'N', 3 > N
Definition: PlasticOps.hpp:118
c
const double c
speed of light (cm/ns)
Definition: initial_diffusion.cpp:39
PlasticOps::M
FTensor::Index< 'M', 3 > M
Definition: PlasticOps.hpp:117
PlasticOps::CommonData::H
@ H
Definition: PlasticOps.hpp:49
PlasticOps::diff_constrain_deqiv
double diff_constrain_deqiv(double sign, double eqiv, double dot_tau, double sigma_Y)
Definition: PlasticOpsGeneric.hpp:317
FTensor::Number
Definition: Number.hpp:11
DomainEleOp
DomainEle::UserDataOperator DomainEleOp
Finire element operator type.
Definition: child_and_parent.cpp:36
a
constexpr double a
Definition: approx_sphere.cpp:30
PlasticOps::CommonData::activityData
static std::array< int, 5 > activityData
Definition: PlasticOps.hpp:107
kinematic_hardening
auto kinematic_hardening(FTensor::Tensor2_symmetric< T, DIM > &t_plastic_strain, double C1_k)
Definition: plastic.cpp:91
MoFEM::getFTensor0FromVec
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
Definition: Templates.hpp:135
PlasticOps::deviator
auto deviator(FTensor::Tensor2_symmetric< T, DIM > &t_stress, double trace, FTensor::Tensor2_symmetric< double, DIM > &t_alpha, FTensor::Number< DIM >)
Definition: PlasticOpsGeneric.hpp:92
size_symm
constexpr auto size_symm
Definition: plastic.cpp:42
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
field_name
constexpr auto field_name
Definition: poisson_2d_homogeneous.cpp:13
FTensor::Index< 'i', DIM >
PlasticOps::CommonData::QINF
@ QINF
Definition: PlasticOps.hpp:51
convert.n
n
Definition: convert.py:82
PlasticOps::OpCalculatePlasticityImpl< DIM, GAUSS, DomainEleOp >::mDPtr
boost::shared_ptr< MatrixDouble > mDPtr
Definition: PlasticOpsGeneric.hpp:477
iso_hardening
double iso_hardening(double tau, double H, double Qinf, double b_iso, double sigmaY)
Definition: plastic.cpp:72
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
FTensor::Ddg
Definition: Ddg_value.hpp:7
PlasticOps::constrain_abs
double constrain_abs(double x, double dt)
Definition: PlasticOpsGeneric.hpp:265
PlasticOps::constraint
double constraint(double eqiv, double dot_tau, double f, double sigma_y, double abs_w, double vis_H, double sigma_Y)
Definition: PlasticOpsGeneric.hpp:296
PlasticOps::CommonData::VIS_H
@ VIS_H
Definition: PlasticOps.hpp:50
m
FTensor::Index< 'm', 3 > m
Definition: shallow_wave.cpp:80
PlasticOps::trace
double trace(FTensor::Tensor2_symmetric< T, 2 > &t_stress)
Definition: PlasticOpsGeneric.hpp:80
PlasticOps::I
FTensor::Index< 'I', 3 > I
[Common data]
Definition: PlasticOps.hpp:115
dt
double dt
Definition: heat_method.cpp:26
k
FTensor::Index< 'k', 3 > k
Definition: matrix_function.cpp:20
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
PlasticOps::diff_constrain_dsigma_y
auto diff_constrain_dsigma_y(double sign)
Definition: PlasticOpsGeneric.hpp:324
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
l
FTensor::Index< 'l', 3 > l
Definition: matrix_function.cpp:21
PlasticOps::CommonData::SIGMA_Y
@ SIGMA_Y
Definition: PlasticOps.hpp:48
Generated by Doxygen 1.8.17 and hosted at University of Glasgow