v0.15.0
Loading...
Searching...
No Matches
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 ()
 Set Hemholtz energy.
 
MoFEMErrorCode codedHessianW (vector< double * > hessian)
 
- 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.
 
MoFEMErrorCode setParams (short tag, bool &recalculate_elastic_tangent)
 Set material parameters at integration point.
 
MoFEMErrorCode freeHemholtzEnergy ()
 [Free energy J2]
 
MoFEMErrorCode evalF ()
 [Yield function J2]
 
MoFEMErrorCode yieldFunction ()
 Set yield function.
 
MoFEMErrorCode flowPotential ()
 [Yield function J2]
 
MoFEMErrorCode codedHessianW (vector< double * > hessian)
 
- 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.
 
MoFEMErrorCode recordY ()
 Record yield function.
 
MoFEMErrorCode recordH ()
 Record flow potential.
 
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.
 
MoFEMErrorCode solveClosestProjection ()
 Solve nonlinear system of equations to find stress, internal fluxes, and Lagrange plastic multiplier.
 
MoFEMErrorCode consistentTangent ()
 Calculate consistent tangent matrix.
 
MoFEMErrorCode recordTapes ()
 Record tapes.
 

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]
 
- 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
 
PetscBool implHessianW
 
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 467 of file ADOLCPlasticityMaterialModels.hpp.

Member Function Documentation

◆ codedHessianW()

MoFEMErrorCode ADOLCPlasticity::J2Plasticity< 2 >::codedHessianW ( vector< double * > hessian)
inlinevirtual

Reimplemented from ADOLCPlasticity::ClosestPointProjection.

Definition at line 518 of file ADOLCPlasticityMaterialModels.hpp.

518 {
520 // Not implemented yet
521 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
522 "2D Energy Hessian not yet implemented");
524 }
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_NOT_IMPLEMENTED
Definition definitions.h:32
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()

◆ freeHemholtzEnergy()

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

Set Hemholtz energy.

[Plane stress]

[Plane stress]

Implements ADOLCPlasticity::ClosestPointProjection.

Definition at line 470 of file ADOLCPlasticityMaterialModels.hpp.

470 {
472
473 auto p_lambda = mkparam(lambda); // 0
474 auto p_mu = mkparam(mu); // 1
475 auto p_H = mkparam(H); // 2
476 auto p_K = mkparam(K); // 3
477 auto p_phi = mkparam(phi); // 4
478 auto p_sigma_y = mkparam(sigmaY); // 5
479
480 auto t_voight_strain_op = voight_to_strain_op();
481 auto t_voight_total_strain = getFTensor1FromPtr<6>(&a_sTrain[0]);
482 auto t_voight_plastic_strain = getFTensor1FromPtr<6>(&a_plasticStrain[0]);
483 auto t_voight_internal_tensor =
484 getFTensor1FromPtr<6>(&a_internalVariables[1]);
485
486 tTotalStrain(i, j) = t_voight_strain_op(i, j, Z) * t_voight_total_strain(Z);
487 tPlasticStrain(i, j) =
488 t_voight_strain_op(i, j, Z) * t_voight_plastic_strain(Z);
490 t_voight_strain_op(i, j, Z) * t_voight_internal_tensor(Z);
491
493
494 //! [Plane stress]
495 // Calculate st rain at ezz enforcing plane strain case
496 tElasticStrain(2, 2) =
497 (-tElasticStrain(0, 0) * p_lambda - tElasticStrain(1, 1) * p_lambda) /
498 (p_lambda + 2 * p_mu);
499 //! [Plane stress]
500
501 tR = tElasticStrain(i, i);
502
503 a_w = 0.5 * p_lambda * tR * tR +
504 p_mu * (tElasticStrain(i, j) * tElasticStrain(i, j));
505
506 // plastic part
507 // isotropic
508 a_w +=
509 a_internalVariables[0] * p_sigma_y +
510 0.5 * p_phi * p_H * (a_internalVariables[0] * a_internalVariables[0]);
511 // kinematic
512 a_w += (0.5 * (1 - p_phi) * p_K) *
514
516 }
FTensor::Dg< double, 3, 6 > voight_to_strain_op()
Op convert Vight strain vector to strain tensor.
FTensor::Tensor2_symmetric< adouble, 3 > tInternalTensor
FTensor::Tensor2_symmetric< adouble, 3 > tPlasticStrain
FTensor::Tensor2_symmetric< adouble, 3 > tTotalStrain
FTensor::Tensor2_symmetric< adouble, 3 > tElasticStrain

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