v0.14.0
Public Member Functions | List of all members
ADOLCPlasticity::J2Plasticity< 2 > Struct Reference

[J2 2D] More...

#include <users_modules/adolc-plasticity/src/ADOLCPlasticityMaterialModels.hpp>

Inheritance diagram for ADOLCPlasticity::J2Plasticity< 2 >:
[legend]
Collaboration diagram for ADOLCPlasticity::J2Plasticity< 2 >:
[legend]

Public Member Functions

MoFEMErrorCode freeHemholtzEnergy ()
 [Free energy J2] More...
 
- Public Member Functions inherited from ADOLCPlasticity::J2Plasticity< 3 >
 J2Plasticity ()
 
MoFEMErrorCode getDefaultMaterialParameters ()
 
MoFEMErrorCode addMatBlockOps (MoFEM::Interface &m_field, boost::ptr_deque< ForcesAndSourcesCore::UserDataOperator > &pip, std::string block_name, Sev sev)
 Get material parameters from block. More...
 
MoFEMErrorCode setParams (short tag, bool &recalculate_elastic_tangent)
 Set material parameters at integration point. More...
 
MoFEMErrorCode evalF ()
 [Yield function J2] More...
 
MoFEMErrorCode yieldFunction ()
 Set yield function. More...
 
MoFEMErrorCode flowPotential ()
 [Yield function J2] More...
 
- Public Member Functions inherited from ADOLCPlasticity::ClosestPointProjection
VectorAdaptor getPlasticStrain ()
 
VectorAdaptor getTotalStrain ()
 
VectorAdaptor getInternalVariables ()
 
VectorAdaptor getActiveVariablesYH ()
 
VectorAdaptor getStress ()
 
VectorAdaptor getInternalFluxes ()
 
 ClosestPointProjection ()
 
MoFEMErrorCode recordW ()
 Record strain energy. More...
 
MoFEMErrorCode recordY ()
 Record yield function. More...
 
MoFEMErrorCode recordH ()
 Record flow potential. More...
 
MoFEMErrorCode playW ()
 
MoFEMErrorCode playW_NoHessian ()
 
MoFEMErrorCode playY ()
 
MoFEMErrorCode playY_NoGradient ()
 
MoFEMErrorCode playH ()
 
MoFEMErrorCode playH_NoHessian ()
 
MoFEMErrorCode createMatAVecR ()
 
MoFEMErrorCode evaluatePotentials ()
 
MoFEMErrorCode playPotentials ()
 
MoFEMErrorCode playPotentials_NoHessian ()
 
MoFEMErrorCode calculateR (Vec R)
 
MoFEMErrorCode calculateA ()
 
MoFEMErrorCode snesCreate ()
 Create nested snes. More...
 
MoFEMErrorCode solveClosetProjection ()
 Solve nonlinear system of equations to find stress, internal fluxes, and Lagrange plastic multiplier. More...
 
MoFEMErrorCode consistentTangent ()
 Calculate consistent tangent matrix. More...
 
MoFEMErrorCode recordTapes ()
 Record tapes. More...
 

Additional Inherited Members

- Public Types inherited from ADOLCPlasticity::J2Plasticity< 3 >
enum  ModelParams {
  LAMBDA, MU, ISOH, KINK,
  PHI, SIGMA_Y, LAST_PARAM
}
 
using ParamsArray = std::array< double, LAST_PARAM >
 
- Public Types inherited from ADOLCPlasticity::ClosestPointProjection
enum  TypesTags { W = 0, Y, H, LAST_TAPE }
 
- Public Attributes inherited from ADOLCPlasticity::J2Plasticity< 3 >
FTensor::Index< 'i', 3 > i
 
FTensor::Index< 'j', 3 > j
 
FTensor::Index< 'Z', 6 > Z
 
FTensor::Tensor2_symmetric< adouble, 3 > tTotalStrain
 
FTensor::Tensor2_symmetric< adouble, 3 > tPlasticStrain
 
FTensor::Tensor2_symmetric< adouble, 3 > tElasticStrain
 
FTensor::Tensor2_symmetric< adouble, 3 > tInternalTensor
 
FTensor::Tensor2_symmetric< adouble, 3 > tStress
 
FTensor::Tensor2_symmetric< adouble, 3 > tBackStress
 
FTensor::Tensor2_symmetric< adouble, 3 > tDeviator
 
adouble tR
 
double mu
 
double lambda
 
double H
 
double K
 
double phi
 
double sigmaY
 
boost::shared_ptr< ParamsArraydefaultParamsArrayPtr = nullptr
 
boost::shared_ptr< ParamsArrayparamsArrayPtr = nullptr
 
boost::shared_ptr< ParamsArrayoldParamsArrayPtr = nullptr
 
adouble f
 [Free energy J2] More...
 
- Public Attributes inherited from ADOLCPlasticity::ClosestPointProjection
int nbInternalVariables
 
boost::function< int(int, int, int)> integrationRule
 
VectorDouble internalVariables0
 
VectorDouble plasticStrain0
 
std::array< int, LAST_TAPEtapesTags
 
VectorAdaptor activeVariablesW
 
VectorAdaptor gradientW
 
double w
 
double y
 
double h
 
double deltaGamma
 
MatrixDouble Ep
 
MatrixDouble Cp
 
MatrixDouble Cep
 
ublas::symmetric_matrix< double, ublas::lower > C
 
ublas::symmetric_matrix< double, ublas::lower > D
 
MatrixDouble hessianW
 
VectorDouble gradientY
 
VectorDouble gradientH
 
MatrixDouble hessianH
 
MatrixDouble partialWStrainPlasticStrain
 
VectorAdaptor partialYSigma
 
VectorAdaptor partialYFlux
 
VectorAdaptor partialHSigma
 
VectorAdaptor partialHFlux
 
ublas::symmetric_matrix< double, ublas::lower > partial2HSigma
 
ublas::symmetric_matrix< double, ublas::lower > partial2HFlux
 
MatrixDouble partial2HSigmaFlux
 
SmartPetscObj< Mat > A
 
SmartPetscObj< Vec > R
 
SmartPetscObj< Vec > Chi
 
SmartPetscObj< Vec > dChi
 
ublas::matrix< double, ublas::column_major > dataA
 
SmartPetscObj< SNES > sNes
 
ublas::vector< adoublea_sTrain
 
ublas::vector< adoublea_plasticStrain
 
ublas::vector< adoublea_internalVariables
 
ublas::vector< adoublea_sTress
 
ublas::vector< adoublea_internalFluxes
 
adouble a_w
 
adouble a_y
 
adouble a_h
 

Detailed Description

[J2 2D]

Definition at line 410 of file ADOLCPlasticityMaterialModels.hpp.

Member Function Documentation

◆ freeHemholtzEnergy()

MoFEMErrorCode ADOLCPlasticity::J2Plasticity< 2 >::freeHemholtzEnergy ( )
inlinevirtual

[Free energy J2]

[Plane stress]

[Plane stress]

Reimplemented from ADOLCPlasticity::J2Plasticity< 3 >.

Definition at line 413 of file ADOLCPlasticityMaterialModels.hpp.

413  {
415 
416  auto p_lambda = mkparam(lambda); // 0
417  auto p_mu = mkparam(mu); // 1
418  auto p_H = mkparam(H); // 2
419  auto p_K = mkparam(K); // 3
420  auto p_phi = mkparam(phi); // 4
421  auto p_sigma_y = mkparam(sigmaY); // 5
422 
423  auto t_voight_strain_op = voight_to_strain_op();
424  auto t_voight_total_strain = getFTensor1FromPtr<6>(&a_sTrain[0]);
425  auto t_voight_plastic_strain = getFTensor1FromPtr<6>(&a_plasticStrain[0]);
426  auto t_voight_internal_tensor =
427  getFTensor1FromPtr<6>(&a_internalVariables[1]);
428 
429  tTotalStrain(i, j) = t_voight_strain_op(i, j, Z) * t_voight_total_strain(Z);
430  tPlasticStrain(i, j) =
431  t_voight_strain_op(i, j, Z) * t_voight_plastic_strain(Z);
432  tInternalTensor(i, j) =
433  t_voight_strain_op(i, j, Z) * t_voight_internal_tensor(Z);
434 
436 
437  //! [Plane stress]
438  // Calculate st rain at ezz enforcing plane strain case
439  tElasticStrain(2, 2) =
440  (-tElasticStrain(0, 0) * p_lambda - tElasticStrain(1, 1) * p_lambda) /
441  (p_lambda + 2 * p_mu);
442  //! [Plane stress]
443 
444  tR = tElasticStrain(i, i);
445 
446  a_w = 0.5 * p_lambda * tR * tR +
447  p_mu * (tElasticStrain(i, j) * tElasticStrain(i, j));
448 
449  // plastic part
450  // isotropic
451  a_w +=
452  a_internalVariables[0] * p_sigma_y +
453  0.5 * p_phi * p_H * (a_internalVariables[0] * a_internalVariables[0]);
454  // kinematic
455  a_w += (0.5 * (1 - p_phi) * p_K) *
457 
459  }

The documentation for this struct was generated from the following file:
ADOLCPlasticity::J2Plasticity< 3 >::lambda
double lambda
Definition: ADOLCPlasticityMaterialModels.hpp:95
ADOLCPlasticity::ClosestPointProjection::a_plasticStrain
ublas::vector< adouble > a_plasticStrain
Definition: ADOLCPlasticity.hpp:325
ADOLCPlasticity::ClosestPointProjection::a_sTrain
ublas::vector< adouble > a_sTrain
Definition: ADOLCPlasticity.hpp:324
ADOLCPlasticity::J2Plasticity< 3 >::i
FTensor::Index< 'i', 3 > i
Definition: ADOLCPlasticityMaterialModels.hpp:80
ADOLCPlasticity::J2Plasticity< 3 >::tR
adouble tR
Definition: ADOLCPlasticityMaterialModels.hpp:92
ADOLCPlasticity::J2Plasticity< 3 >::Z
FTensor::Index< 'Z', 6 > Z
Definition: ADOLCPlasticityMaterialModels.hpp:82
ADOLCPlasticity::J2Plasticity< 3 >::j
FTensor::Index< 'j', 3 > j
Definition: ADOLCPlasticityMaterialModels.hpp:81
ADOLCPlasticity::voight_to_strain_op
FTensor::Dg< double, 3, 6 > voight_to_strain_op()
Op convert Vight strain vector to strain tensor.
Definition: ADOLCPlasticity.hpp:29
ADOLCPlasticity::J2Plasticity< 3 >::K
double K
Definition: ADOLCPlasticityMaterialModels.hpp:97
ADOLCPlasticity::ClosestPointProjection::a_w
adouble a_w
Definition: ADOLCPlasticity.hpp:328
ADOLCPlasticity::J2Plasticity< 3 >::tElasticStrain
FTensor::Tensor2_symmetric< adouble, 3 > tElasticStrain
Definition: ADOLCPlasticityMaterialModels.hpp:86
ADOLCPlasticity::J2Plasticity< 3 >::mu
double mu
Definition: ADOLCPlasticityMaterialModels.hpp:94
ADOLCPlasticity::J2Plasticity< 3 >::tInternalTensor
FTensor::Tensor2_symmetric< adouble, 3 > tInternalTensor
Definition: ADOLCPlasticityMaterialModels.hpp:87
ADOLCPlasticity::ClosestPointProjection::a_internalVariables
ublas::vector< adouble > a_internalVariables
Definition: ADOLCPlasticity.hpp:326
ADOLCPlasticity::J2Plasticity< 3 >::phi
double phi
Definition: ADOLCPlasticityMaterialModels.hpp:98
ADOLCPlasticity::J2Plasticity< 3 >::tPlasticStrain
FTensor::Tensor2_symmetric< adouble, 3 > tPlasticStrain
Definition: ADOLCPlasticityMaterialModels.hpp:85
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
ADOLCPlasticity::J2Plasticity< 3 >::tTotalStrain
FTensor::Tensor2_symmetric< adouble, 3 > tTotalStrain
Definition: ADOLCPlasticityMaterialModels.hpp:84
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
ADOLCPlasticity::J2Plasticity< 3 >::sigmaY
double sigmaY
Definition: ADOLCPlasticityMaterialModels.hpp:99
ADOLCPlasticity::J2Plasticity< 3 >::H
double H
Definition: ADOLCPlasticityMaterialModels.hpp:96