www.mooseframework.org
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
FiniteStrainCrystalPlasticity Class Reference

FiniteStrainCrystalPlasticity uses the multiplicative decomposition of deformation gradient and solves the PK2 stress residual equation at the intermediate configuration to evolve the material state. More...

#include <FiniteStrainCrystalPlasticity.h>

Inheritance diagram for FiniteStrainCrystalPlasticity:
[legend]

Public Member Functions

 FiniteStrainCrystalPlasticity (const InputParameters &parameters)
 

Protected Member Functions

virtual void computeQpStress ()
 This function updates the stress at a quadrature point. More...
 
virtual void computeQpElasticityTensor ()
 This function updates the elasticity tensor at a quadrature point. More...
 
virtual void initQpStatefulProperties ()
 This function initializes the stateful properties such as stress, plastic deformation gradient, slip system resistances, etc. More...
 
virtual void calc_resid_jacob (RankTwoTensor &, RankFourTensor &)
 This function calls the residual and jacobian functions used in the stress update algorithm. More...
 
virtual void getSlipIncrements ()
 This function updates the slip increments. More...
 
virtual void update_slip_system_resistance ()
 This function updates the slip system resistances. More...
 
virtual void updateGss ()
 This function updates the slip system resistances. More...
 
virtual void getSlipSystems ()
 This function reads slip system from file - see test. More...
 
virtual void assignSlipSysRes ()
 This function assign initial values of slip system resistances/internal variables read from getSlipSystems(). More...
 
virtual void readFileInitSlipSysRes ()
 This function read slip system resistances from file - see test. More...
 
virtual void getInitSlipSysRes ()
 This function assign slip system resistances - see test. More...
 
virtual void readFileFlowRateParams ()
 This function read flow rate parameters from file - see test. More...
 
virtual void getFlowRateParams ()
 This function assign flow rate parameters - see test. More...
 
virtual void readFileHardnessParams ()
 This function read hardness parameters from file. More...
 
virtual void getHardnessParams ()
 This function assign flow rate parameters from .i file - see test. More...
 
virtual void initSlipSysProps ()
 This function initializes slip system resistances. More...
 
virtual void initAdditionalProps ()
 This function initializes additional parameters. More...
 
virtual void preSolveQp ()
 This function set variables for stress and internal variable solve. More...
 
virtual void solveQp ()
 This function solves stress and internal variables. More...
 
virtual void postSolveQp ()
 This function update stress and internal variable after solve. More...
 
virtual void preSolveStatevar ()
 This function set variables for internal variable solve. More...
 
virtual void solveStatevar ()
 This function solves internal variables. More...
 
virtual void postSolveStatevar ()
 This function update internal variable after solve. More...
 
virtual void preSolveStress ()
 This function set variables for stress solve. More...
 
virtual void solveStress ()
 This function solves for stress, updates plastic deformation gradient. More...
 
virtual void postSolveStress ()
 This function update stress and plastic deformation gradient after solve. More...
 
virtual void calcResidual (RankTwoTensor &)
 This function calculate stress residual. More...
 
virtual void calcJacobian (RankFourTensor &)
 This function calculate jacobian. More...
 
virtual RankFourTensor calcTangentModuli ()
 This function calculate the tangent moduli for preconditioner. More...
 
virtual RankFourTensor elasticTangentModuli ()
 This function calculate the elastic tangent moduli for preconditioner. More...
 
virtual RankFourTensor elastoPlasticTangentModuli ()
 This function calculate the exact tangent moduli for preconditioner. More...
 
RankTwoTensor get_current_rotation (const RankTwoTensor &a)
 This function perform RU decomposition to obtain the rotation tensor. More...
 
RankTwoTensor getMatRot (const RankTwoTensor &a)
 This function perform RU decomposition to obtain the rotation tensor. More...
 
void calc_schmid_tensor ()
 This function calculate the Schmid tensor. More...
 
bool line_search_update (const Real rnorm_prev, const RankTwoTensor)
 This function performs the line search update. More...
 
void internalVariableUpdateNRiteration ()
 This function updates internal variables after each NewTon Raphson iteration (_fp_inv) More...
 
virtual void computeQpProperties () override
 

Protected Attributes

const unsigned int _nss
 Number of slip system resistance. More...
 
std::vector< Real > _gprops
 
std::vector< Real > _hprops
 
std::vector< Real > _flowprops
 
std::string _slip_sys_file_name
 File should contain slip plane normal and direction. See test. More...
 
std::string _slip_sys_res_prop_file_name
 File should contain initial values of the slip system resistances. More...
 
std::string _slip_sys_flow_prop_file_name
 File should contain values of the flow rate equation parameters. More...
 
std::string _slip_sys_hard_prop_file_name
 The hardening parameters in this class are read from .i file. The user can override to read from file. More...
 
Real _rtol
 Stress residual equation relative tolerance. More...
 
Real _abs_tol
 Stress residual equation absolute tolerance. More...
 
Real _gtol
 Internal variable update equation tolerance. More...
 
Real _slip_incr_tol
 Slip increment tolerance. More...
 
unsigned int _maxiter
 Maximum number of iterations for stress update. More...
 
unsigned int _maxiterg
 Maximum number of iterations for internal variable update. More...
 
unsigned int _num_slip_sys_flowrate_props
 Number of slip system flow rate parameters. More...
 
MooseEnum _tan_mod_type
 Type of tangent moduli calculation. More...
 
MooseEnum _intvar_read_type
 Read from options for initial values of internal variables. More...
 
unsigned int _num_slip_sys_props
 Number of slip system specific properties provided in the file containing slip system normals and directions. More...
 
bool _gen_rndm_stress_flag
 
bool _input_rndm_scale_var
 Input option for scaling variable to generate random stress when convergence fails. More...
 
Real _rndm_scale_var
 Scaling value. More...
 
unsigned int _rndm_seed
 Seed value. More...
 
unsigned int _max_substep_iter
 Maximum number of substep iterations. More...
 
bool _use_line_search
 Flag to activate line serach. More...
 
Real _min_lsrch_step
 Minimum line search step size. More...
 
Real _lsrch_tol
 Line search bisection method tolerance. More...
 
unsigned int _lsrch_max_iter
 Line search bisection method maximum iteration number. More...
 
MooseEnum _lsrch_method
 
MaterialProperty< RankTwoTensor > & _fp
 
const MaterialProperty< RankTwoTensor > & _fp_old
 
MaterialProperty< RankTwoTensor > & _pk2
 
const MaterialProperty< RankTwoTensor > & _pk2_old
 
MaterialProperty< RankTwoTensor > & _lag_e
 
const MaterialProperty< RankTwoTensor > & _lag_e_old
 
MaterialProperty< std::vector< Real > > & _gss
 
const MaterialProperty< std::vector< Real > > & _gss_old
 
MaterialProperty< Real > & _acc_slip
 
const MaterialProperty< Real > & _acc_slip_old
 
MaterialProperty< RankTwoTensor > & _update_rot
 
const MaterialProperty< RankTwoTensor > & _deformation_gradient
 
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
 
const MaterialProperty< RankFourTensor > & _elasticity_tensor
 
const MaterialProperty< RankTwoTensor > & _crysrot
 
DenseVector< Real > _mo
 
DenseVector< Real > _no
 
DenseVector< Real > _a0
 
DenseVector< Real > _xm
 
Real _h0
 
Real _tau_sat
 
Real _tau_init
 
Real _r
 
RankTwoTensor _dfgrd_tmp
 
RankTwoTensor _fe
 
RankTwoTensor _fp_old_inv
 
RankTwoTensor _fp_inv
 
RankTwoTensor _fp_prev_inv
 
DenseVector< Real > _slip_incr
 
DenseVector< Real > _tau
 
DenseVector< Real > _dslipdtau
 
std::vector< RankTwoTensor_s0
 
RankTwoTensor _pk2_tmp
 
RankTwoTensor _pk2_tmp_old
 
Real _accslip_tmp
 
Real _accslip_tmp_old
 
std::vector< Real > _gss_tmp
 
std::vector< Real > _gss_tmp_old
 
DenseVector< Real > _slip_sys_props
 
DenseMatrix< Real > _dgss_dsliprate
 
bool _read_from_slip_sys_file
 
bool _err_tol
 
RankTwoTensor _delta_dfgrd
 Flag to check whether convergence is achieved. More...
 
RankTwoTensor _dfgrd_tmp_old
 
Real _dfgrd_scale_factor
 Scales the substepping increment to obtain deformation gradient at a substep iteration. More...
 
bool _first_step_iter
 Flags to reset variables and reinitialize variables. More...
 
bool _last_step_iter
 
bool _first_substep
 
const std::string _base_name
 
const std::string _elasticity_tensor_name
 
const MaterialProperty< RankTwoTensor > & _mechanical_strain
 
MaterialProperty< RankTwoTensor > & _stress
 
MaterialProperty< RankTwoTensor > & _elastic_strain
 
const MaterialProperty< RankTwoTensor > & _extra_stress
 Extra stress tensor. More...
 
std::vector< Function * > _initial_stress_fcn
 initial stress components More...
 
MaterialProperty< RankFourTensor > & _Jacobian_mult
 derivative of stress w.r.t. strain (_dstress_dstrain) More...
 

Detailed Description

FiniteStrainCrystalPlasticity uses the multiplicative decomposition of deformation gradient and solves the PK2 stress residual equation at the intermediate configuration to evolve the material state.

The internal variables are updated using an interative predictor-corrector algorithm. Backward Euler integration rule is used for the rate equations.

Definition at line 27 of file FiniteStrainCrystalPlasticity.h.

Constructor & Destructor Documentation

◆ FiniteStrainCrystalPlasticity()

FiniteStrainCrystalPlasticity::FiniteStrainCrystalPlasticity ( const InputParameters &  parameters)

Definition at line 100 of file FiniteStrainCrystalPlasticity.C.

101  : ComputeStressBase(parameters),
102  _nss(getParam<int>("nss")),
103  _gprops(getParam<std::vector<Real>>("gprops")),
104  _hprops(getParam<std::vector<Real>>("hprops")),
105  _flowprops(getParam<std::vector<Real>>("flowprops")),
106  _slip_sys_file_name(getParam<FileName>("slip_sys_file_name")),
107  _slip_sys_res_prop_file_name(getParam<FileName>("slip_sys_res_prop_file_name")),
108  _slip_sys_flow_prop_file_name(getParam<FileName>("slip_sys_flow_prop_file_name")),
109  _slip_sys_hard_prop_file_name(getParam<FileName>("slip_sys_hard_prop_file_name")),
110  _rtol(getParam<Real>("rtol")),
111  _abs_tol(getParam<Real>("abs_tol")),
112  _gtol(getParam<Real>("gtol")),
113  _slip_incr_tol(getParam<Real>("slip_incr_tol")),
114  _maxiter(getParam<unsigned int>("maxiter")),
115  _maxiterg(getParam<unsigned int>("maxitergss")),
116  _num_slip_sys_flowrate_props(getParam<unsigned int>("num_slip_sys_flowrate_props")),
117  _tan_mod_type(getParam<MooseEnum>("tan_mod_type")),
118  _intvar_read_type(getParam<MooseEnum>("intvar_read_type")),
119  _num_slip_sys_props(getParam<unsigned int>("num_slip_sys_props")),
120  _gen_rndm_stress_flag(getParam<bool>("gen_random_stress_flag")),
121  _input_rndm_scale_var(getParam<bool>("input_random_scaling_var")),
122  _rndm_scale_var(getParam<Real>("random_scaling_var")),
123  _rndm_seed(getParam<unsigned int>("random_seed")),
124  _max_substep_iter(getParam<unsigned int>("maximum_substep_iteration")),
125  _use_line_search(getParam<bool>("use_line_search")),
126  _min_lsrch_step(getParam<Real>("min_line_search_step_size")),
127  _lsrch_tol(getParam<Real>("line_search_tol")),
128  _lsrch_max_iter(getParam<unsigned int>("line_search_maxiter")),
129  _lsrch_method(getParam<MooseEnum>("line_search_method")),
130  _fp(declareProperty<RankTwoTensor>("fp")), // Plastic deformation gradient
131  _fp_old(getMaterialPropertyOld<RankTwoTensor>(
132  "fp")), // Plastic deformation gradient of previous increment
133  _pk2(declareProperty<RankTwoTensor>("pk2")), // 2nd Piola Kirchoff Stress
134  _pk2_old(getMaterialPropertyOld<RankTwoTensor>(
135  "pk2")), // 2nd Piola Kirchoff Stress of previous increment
136  _lag_e(declareProperty<RankTwoTensor>("lage")), // Lagrangian strain
137  _lag_e_old(
138  getMaterialPropertyOld<RankTwoTensor>("lage")), // Lagrangian strain of previous increment
139  _gss(declareProperty<std::vector<Real>>("gss")), // Slip system resistances
140  _gss_old(getMaterialPropertyOld<std::vector<Real>>(
141  "gss")), // Slip system resistances of previous increment
142  _acc_slip(declareProperty<Real>("acc_slip")), // Accumulated slip
144  getMaterialPropertyOld<Real>("acc_slip")), // Accumulated alip of previous increment
145  _update_rot(declareProperty<RankTwoTensor>(
146  "update_rot")), // Rotation tensor considering material rotation and crystal orientation
147  _deformation_gradient(getMaterialProperty<RankTwoTensor>("deformation_gradient")),
148  _deformation_gradient_old(getMaterialPropertyOld<RankTwoTensor>("deformation_gradient")),
149  _elasticity_tensor(getMaterialProperty<RankFourTensor>("elasticity_tensor")),
150  _crysrot(getMaterialProperty<RankTwoTensor>("crysrot")),
151  _mo(_nss * LIBMESH_DIM),
152  _no(_nss * LIBMESH_DIM),
153  _slip_incr(_nss),
154  _tau(_nss),
155  _dslipdtau(_nss),
156  _s0(_nss),
157  _gss_tmp(_nss),
160 {
161  _err_tol = false;
162 
163  if (_num_slip_sys_props > 0)
165 
166  _pk2_tmp.zero();
167  _delta_dfgrd.zero();
168 
169  _first_step_iter = false;
170  _last_step_iter = false;
171  // Initialize variables in the first iteration of substepping
172  _first_substep = true;
173 
174  _read_from_slip_sys_file = false;
175  if (_intvar_read_type == "slip_sys_file")
177 
179  mooseError("Crystal Plasticity Error: Specify number of internal variable's initial values to "
180  "be read from slip system file");
181 
182  getSlipSystems();
183 
184  RankTwoTensor::initRandom(_rndm_seed);
185 }
const MaterialProperty< RankTwoTensor > & _pk2_old
MaterialProperty< RankTwoTensor > & _lag_e
const MaterialProperty< RankTwoTensor > & _deformation_gradient
bool _input_rndm_scale_var
Input option for scaling variable to generate random stress when convergence fails.
bool _use_line_search
Flag to activate line serach.
const MaterialProperty< std::vector< Real > > & _gss_old
const MaterialProperty< RankFourTensor > & _elasticity_tensor
MaterialProperty< RankTwoTensor > & _pk2
std::string _slip_sys_flow_prop_file_name
File should contain values of the flow rate equation parameters.
unsigned int _num_slip_sys_flowrate_props
Number of slip system flow rate parameters.
const MaterialProperty< RankTwoTensor > & _fp_old
Real _rtol
Stress residual equation relative tolerance.
std::string _slip_sys_hard_prop_file_name
The hardening parameters in this class are read from .i file. The user can override to read from file...
unsigned int _maxiter
Maximum number of iterations for stress update.
const unsigned int _nss
Number of slip system resistance.
RankTwoTensor _delta_dfgrd
Flag to check whether convergence is achieved.
MaterialProperty< std::vector< Real > > & _gss
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
const MaterialProperty< RankTwoTensor > & _crysrot
MaterialProperty< RankTwoTensor > & _fp
unsigned int _lsrch_max_iter
Line search bisection method maximum iteration number.
Real _gtol
Internal variable update equation tolerance.
std::string _slip_sys_res_prop_file_name
File should contain initial values of the slip system resistances.
const MaterialProperty< RankTwoTensor > & _lag_e_old
ComputeStressBase(const InputParameters &parameters)
MooseEnum _tan_mod_type
Type of tangent moduli calculation.
Real _lsrch_tol
Line search bisection method tolerance.
virtual void getSlipSystems()
This function reads slip system from file - see test.
unsigned int _num_slip_sys_props
Number of slip system specific properties provided in the file containing slip system normals and dir...
unsigned int _max_substep_iter
Maximum number of substep iterations.
bool _first_step_iter
Flags to reset variables and reinitialize variables.
MaterialProperty< RankTwoTensor > & _update_rot
Real _abs_tol
Stress residual equation absolute tolerance.
unsigned int _maxiterg
Maximum number of iterations for internal variable update.
const MaterialProperty< Real > & _acc_slip_old
MooseEnum _intvar_read_type
Read from options for initial values of internal variables.
Real _min_lsrch_step
Minimum line search step size.
std::string _slip_sys_file_name
File should contain slip plane normal and direction. See test.
Real _slip_incr_tol
Slip increment tolerance.

Member Function Documentation

◆ assignSlipSysRes()

void FiniteStrainCrystalPlasticity::assignSlipSysRes ( )
protectedvirtual

This function assign initial values of slip system resistances/internal variables read from getSlipSystems().

Definition at line 231 of file FiniteStrainCrystalPlasticity.C.

Referenced by initSlipSysProps().

232 {
233  _gss[_qp].resize(_nss);
234 
235  for (unsigned int i = 0; i < _nss; ++i)
236  _gss[_qp][i] = _slip_sys_props(i);
237 }
const unsigned int _nss
Number of slip system resistance.
MaterialProperty< std::vector< Real > > & _gss

◆ calc_resid_jacob()

void FiniteStrainCrystalPlasticity::calc_resid_jacob ( RankTwoTensor resid,
RankFourTensor jac 
)
protectedvirtual

This function calls the residual and jacobian functions used in the stress update algorithm.

Definition at line 878 of file FiniteStrainCrystalPlasticity.C.

Referenced by solveStress().

879 {
880  calcResidual(resid);
881  if (_err_tol)
882  return;
883  calcJacobian(jac);
884 }
virtual void calcJacobian(RankFourTensor &)
This function calculate jacobian.
virtual void calcResidual(RankTwoTensor &)
This function calculate stress residual.

◆ calc_schmid_tensor()

void FiniteStrainCrystalPlasticity::calc_schmid_tensor ( )
protected

This function calculate the Schmid tensor.

Definition at line 1046 of file FiniteStrainCrystalPlasticity.C.

Referenced by preSolveQp().

1047 {
1048  DenseVector<Real> mo(LIBMESH_DIM * _nss), no(LIBMESH_DIM * _nss);
1049 
1050  // Update slip direction and normal with crystal orientation
1051  for (unsigned int i = 0; i < _nss; ++i)
1052  {
1053  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
1054  {
1055  mo(i * LIBMESH_DIM + j) = 0.0;
1056  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
1057  mo(i * LIBMESH_DIM + j) =
1058  mo(i * LIBMESH_DIM + j) + _crysrot[_qp](j, k) * _mo(i * LIBMESH_DIM + k);
1059  }
1060 
1061  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
1062  {
1063  no(i * LIBMESH_DIM + j) = 0.0;
1064  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
1065  no(i * LIBMESH_DIM + j) =
1066  no(i * LIBMESH_DIM + j) + _crysrot[_qp](j, k) * _no(i * LIBMESH_DIM + k);
1067  }
1068  }
1069 
1070  // Calculate Schmid tensor and resolved shear stresses
1071  for (unsigned int i = 0; i < _nss; ++i)
1072  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
1073  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
1074  _s0[i](j, k) = mo(i * LIBMESH_DIM + j) * no(i * LIBMESH_DIM + k);
1075 }
const unsigned int _nss
Number of slip system resistance.
const MaterialProperty< RankTwoTensor > & _crysrot

◆ calcJacobian()

void FiniteStrainCrystalPlasticity::calcJacobian ( RankFourTensor jac)
protectedvirtual

This function calculate jacobian.

Definition at line 924 of file FiniteStrainCrystalPlasticity.C.

Referenced by calc_resid_jacob().

925 {
926  RankFourTensor dfedfpinv, deedfe, dfpinvdpk2;
927 
928  std::vector<RankTwoTensor> dtaudpk2(_nss), dfpinvdslip(_nss);
929 
930  for (unsigned int i = 0; i < _nss; ++i)
931  {
932  dtaudpk2[i] = _s0[i];
933  dfpinvdslip[i] = -_fp_old_inv * _s0[i];
934  }
935 
936  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
937  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
938  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
939  dfedfpinv(i, j, k, j) = _dfgrd_tmp(i, k);
940 
941  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
942  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
943  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
944  {
945  deedfe(i, j, k, i) = deedfe(i, j, k, i) + _fe(k, j) * 0.5;
946  deedfe(i, j, k, j) = deedfe(i, j, k, j) + _fe(k, i) * 0.5;
947  }
948 
949  for (unsigned int i = 0; i < _nss; ++i)
950  dfpinvdpk2 += (dfpinvdslip[i] * _dslipdtau(i)).outerProduct(dtaudpk2[i]);
951 
952  jac =
953  RankFourTensor::IdentityFour() - (_elasticity_tensor[_qp] * deedfe * dfedfpinv * dfpinvdpk2);
954 }
const MaterialProperty< RankFourTensor > & _elasticity_tensor
const unsigned int _nss
Number of slip system resistance.

◆ calcResidual()

void FiniteStrainCrystalPlasticity::calcResidual ( RankTwoTensor resid)
protectedvirtual

This function calculate stress residual.

Definition at line 887 of file FiniteStrainCrystalPlasticity.C.

Referenced by calc_resid_jacob(), and line_search_update().

888 {
889  RankTwoTensor iden(RankTwoTensor::initIdentity), ce, ee, ce_pk2, eqv_slip_incr, pk2_new;
890 
891  _fe = _dfgrd_tmp * _fp_prev_inv; // _fp_inv ==> _fp_prev_inv
892 
893  ce = _fe.transpose() * _fe;
894  ce_pk2 = ce * _pk2_tmp;
895  ce_pk2 = ce_pk2 / _fe.det();
896 
897  // Calculate Schmid tensor and resolved shear stresses
898  for (unsigned int i = 0; i < _nss; ++i)
899  _tau(i) = ce_pk2.doubleContraction(_s0[i]);
900 
901  getSlipIncrements(); // Calculate dslip,dslipdtau
902 
903  if (_err_tol)
904  return;
905 
906  eqv_slip_incr.zero();
907  for (unsigned int i = 0; i < _nss; ++i)
908  eqv_slip_incr += _s0[i] * _slip_incr(i);
909 
910  eqv_slip_incr = iden - eqv_slip_incr;
911  _fp_inv = _fp_old_inv * eqv_slip_incr;
912  _fe = _dfgrd_tmp * _fp_inv;
913 
914  ce = _fe.transpose() * _fe;
915  ee = ce - iden;
916  ee *= 0.5;
917 
918  pk2_new = _elasticity_tensor[_qp] * ee;
919 
920  resid = _pk2_tmp - pk2_new;
921 }
const MaterialProperty< RankFourTensor > & _elasticity_tensor
const unsigned int _nss
Number of slip system resistance.
virtual void getSlipIncrements()
This function updates the slip increments.

◆ calcTangentModuli()

RankFourTensor FiniteStrainCrystalPlasticity::calcTangentModuli ( )
protectedvirtual

This function calculate the tangent moduli for preconditioner.

Default is the elastic stiffness matrix. Exact jacobian is currently implemented. tan_mod_type can be modified to exact in .i file to turn it on.

Definition at line 1029 of file FiniteStrainCrystalPlasticity.C.

Referenced by postSolveQp().

1030 {
1031  RankFourTensor tan_mod;
1032 
1033  switch (_tan_mod_type)
1034  {
1035  case 0:
1036  tan_mod = elastoPlasticTangentModuli();
1037  break;
1038  default:
1039  tan_mod = elasticTangentModuli();
1040  }
1041 
1042  return tan_mod;
1043 }
virtual RankFourTensor elastoPlasticTangentModuli()
This function calculate the exact tangent moduli for preconditioner.
virtual RankFourTensor elasticTangentModuli()
This function calculate the elastic tangent moduli for preconditioner.
MooseEnum _tan_mod_type
Type of tangent moduli calculation.

◆ computeQpElasticityTensor()

void FiniteStrainCrystalPlasticity::computeQpElasticityTensor ( )
protectedvirtual

This function updates the elasticity tensor at a quadrature point.

Presently void.

Definition at line 1024 of file FiniteStrainCrystalPlasticity.C.

1025 {
1026 }

◆ computeQpProperties()

void ComputeStressBase::computeQpProperties ( )
overrideprotectedvirtualinherited

Definition at line 51 of file ComputeStressBase.C.

52 {
54 
55  // Add in extra stress
56  _stress[_qp] += _extra_stress[_qp];
57 }
virtual void computeQpStress()=0
MaterialProperty< RankTwoTensor > & _stress
const MaterialProperty< RankTwoTensor > & _extra_stress
Extra stress tensor.

◆ computeQpStress()

void FiniteStrainCrystalPlasticity::computeQpStress ( )
protectedvirtual

This function updates the stress at a quadrature point.

Solves stress residual equation using NR.

Updates slip system resistances iteratively.

Implements ComputeStressBase.

Definition at line 490 of file FiniteStrainCrystalPlasticity.C.

491 {
492  unsigned int substep_iter = 1; // Depth of substepping; Limited to maximum substep iteration
493  unsigned int num_substep = 1; // Calculated from substep_iter as 2^substep_iter
494  Real dt_original = _dt; // Stores original _dt; Reset at the end of solve
495  _first_substep = true; // Initialize variables at substep_iter = 1
496 
497  if (_max_substep_iter > 1)
498  {
500  if (_dfgrd_tmp_old.det() == 0)
501  _dfgrd_tmp_old.addIa(1.0);
502 
504  _err_tol = true; // Indicator to continue substepping
505  }
506 
507  // Substepping loop
508  while (_err_tol && _max_substep_iter > 1)
509  {
510  _dt = dt_original / num_substep;
511 
512  for (unsigned int istep = 0; istep < num_substep; ++istep)
513  {
514  _first_step_iter = false;
515  if (istep == 0)
516  _first_step_iter = true;
517 
518  _last_step_iter = false;
519  if (istep == num_substep - 1)
520  _last_step_iter = true;
521 
522  _dfgrd_scale_factor = (static_cast<Real>(istep) + 1) / num_substep;
523 
524  preSolveQp();
525  solveQp();
526 
527  if (_err_tol)
528  {
529  substep_iter++;
530  num_substep *= 2;
531  break;
532  }
533  }
534 
535  _first_substep = false; // Prevents reinitialization
536  _dt = dt_original; // Resets dt
537 
538 #ifdef DEBUG
539  if (substep_iter > _max_substep_iter)
540  mooseWarning("FiniteStrainCrystalPlasticity: Failure with substepping");
541 #endif
542 
543  if (!_err_tol || substep_iter > _max_substep_iter)
544  postSolveQp(); // Evaluate variables after successful solve or indicate failure
545  }
546 
547  // No substepping
548  if (_max_substep_iter == 1)
549  {
550  preSolveQp();
551  solveQp();
552  postSolveQp();
553  }
554 }
const MaterialProperty< RankTwoTensor > & _deformation_gradient
Real _dfgrd_scale_factor
Scales the substepping increment to obtain deformation gradient at a substep iteration.
RankTwoTensor _delta_dfgrd
Flag to check whether convergence is achieved.
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
virtual void postSolveQp()
This function update stress and internal variable after solve.
unsigned int _max_substep_iter
Maximum number of substep iterations.
bool _first_step_iter
Flags to reset variables and reinitialize variables.
virtual void solveQp()
This function solves stress and internal variables.
virtual void preSolveQp()
This function set variables for stress and internal variable solve.

◆ elasticTangentModuli()

RankFourTensor FiniteStrainCrystalPlasticity::elasticTangentModuli ( )
protectedvirtual

This function calculate the elastic tangent moduli for preconditioner.

Definition at line 1117 of file FiniteStrainCrystalPlasticity.C.

Referenced by calcTangentModuli().

1118 {
1119  return _elasticity_tensor[_qp]; // update jacobian_mult
1120 }
const MaterialProperty< RankFourTensor > & _elasticity_tensor

◆ elastoPlasticTangentModuli()

RankFourTensor FiniteStrainCrystalPlasticity::elastoPlasticTangentModuli ( )
protectedvirtual

This function calculate the exact tangent moduli for preconditioner.

Definition at line 1078 of file FiniteStrainCrystalPlasticity.C.

Referenced by calcTangentModuli().

1079 {
1080  RankFourTensor tan_mod;
1081  RankTwoTensor pk2fet, fepk2;
1082  RankFourTensor deedfe, dsigdpk2dfe;
1083 
1084  // Fill in the matrix stiffness material property
1085 
1086  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
1087  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
1088  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
1089  {
1090  deedfe(i, j, k, i) = deedfe(i, j, k, i) + _fe(k, j) * 0.5;
1091  deedfe(i, j, k, j) = deedfe(i, j, k, j) + _fe(k, i) * 0.5;
1092  }
1093 
1094  dsigdpk2dfe = _fe.mixedProductIkJl(_fe) * _elasticity_tensor[_qp] * deedfe;
1095 
1096  pk2fet = _pk2_tmp * _fe.transpose();
1097  fepk2 = _fe * _pk2_tmp;
1098 
1099  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
1100  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
1101  for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
1102  {
1103  tan_mod(i, j, i, l) = tan_mod(i, j, i, l) + pk2fet(l, j);
1104  tan_mod(i, j, j, l) = tan_mod(i, j, j, l) + fepk2(i, l);
1105  }
1106 
1107  tan_mod += dsigdpk2dfe;
1108 
1109  Real je = _fe.det();
1110  if (je > 0.0)
1111  tan_mod /= je;
1112 
1113  return tan_mod;
1114 }
const MaterialProperty< RankFourTensor > & _elasticity_tensor

◆ get_current_rotation()

RankTwoTensor FiniteStrainCrystalPlasticity::get_current_rotation ( const RankTwoTensor a)
protected

This function perform RU decomposition to obtain the rotation tensor.

Definition at line 982 of file FiniteStrainCrystalPlasticity.C.

Referenced by postSolveQp().

983 {
984  return getMatRot(a);
985 }
RankTwoTensor getMatRot(const RankTwoTensor &a)
This function perform RU decomposition to obtain the rotation tensor.

◆ getFlowRateParams()

void FiniteStrainCrystalPlasticity::getFlowRateParams ( )
protectedvirtual

This function assign flow rate parameters - see test.

.i input file format start_slip_sys_num, end_slip_sys_num, value1, value2

Definition at line 341 of file FiniteStrainCrystalPlasticity.C.

Referenced by initSlipSysProps().

342 {
343  if (_flowprops.size() <= 0)
344  mooseError("FiniteStrainCrystalPLasticity: Error in reading flow rate properties: Specify "
345  "input in .i file or a slip_sys_flow_prop_file_name");
346 
347  _a0.resize(_nss);
348  _xm.resize(_nss);
349 
350  unsigned int num_data_grp = 2 + _num_slip_sys_flowrate_props; // Number of data per group e.g.
351  // start_slip_sys, end_slip_sys,
352  // value1, value2, ..
353 
354  for (unsigned int i = 0; i < _flowprops.size() / num_data_grp; ++i)
355  {
356  Real vs, ve;
357  unsigned int is, ie;
358 
359  vs = _flowprops[i * num_data_grp];
360  ve = _flowprops[i * num_data_grp + 1];
361 
362  if (vs <= 0 || ve <= 0)
363  mooseError("FiniteStrainCrystalPLasticity: Indices in flow rate parameter read must be "
364  "positive integers: is = ",
365  vs,
366  " ie = ",
367  ve);
368 
369  if (vs != floor(vs) || ve != floor(ve))
370  mooseError("FiniteStrainCrystalPLasticity: Error in reading flow props: Values specifying "
371  "start and end number of slip system groups should be integer");
372 
373  is = static_cast<unsigned int>(vs);
374  ie = static_cast<unsigned int>(ve);
375 
376  if (is > ie)
377  mooseError("FiniteStrainCrystalPLasticity: Start index is = ",
378  is,
379  " should be greater than end index ie = ",
380  ie,
381  " in flow rate parameter read");
382 
383  for (unsigned int j = is; j <= ie; ++j)
384  {
385  _a0(j - 1) = _flowprops[i * num_data_grp + 2];
386  _xm(j - 1) = _flowprops[i * num_data_grp + 3];
387  }
388  }
389 
390  for (unsigned int i = 0; i < _nss; ++i)
391  {
392  if (!(_a0(i) > 0.0 && _xm(i) > 0.0))
393  {
394  mooseWarning(
395  "FiniteStrainCrystalPlasticity: Non-positive flow rate parameters ", _a0(i), ",", _xm(i));
396  break;
397  }
398  }
399 }
unsigned int _num_slip_sys_flowrate_props
Number of slip system flow rate parameters.
const unsigned int _nss
Number of slip system resistance.

◆ getHardnessParams()

void FiniteStrainCrystalPlasticity::getHardnessParams ( )
protectedvirtual

This function assign flow rate parameters from .i file - see test.

Definition at line 409 of file FiniteStrainCrystalPlasticity.C.

Referenced by initSlipSysProps().

410 {
411  if (_hprops.size() <= 0)
412  mooseError("FiniteStrainCrystalPLasticity: Error in reading hardness properties: Specify input "
413  "in .i file or a slip_sys_hard_prop_file_name");
414 
415  _r = _hprops[0];
416  _h0 = _hprops[1];
417  _tau_init = _hprops[2];
418  _tau_sat = _hprops[3];
419 }

◆ getInitSlipSysRes()

void FiniteStrainCrystalPlasticity::getInitSlipSysRes ( )
protectedvirtual

This function assign slip system resistances - see test.

.i input file format start_slip_sys_num, end_slip_sys_num, value.

Definition at line 259 of file FiniteStrainCrystalPlasticity.C.

Referenced by initSlipSysProps().

260 {
261  if (_gprops.size() <= 0)
262  mooseError("FiniteStrainCrystalPLasticity: Error in reading slip system resistance properties: "
263  "Specify input in .i file or in slip_sys_res_prop_file or in slip_sys_file");
264 
265  _gss[_qp].resize(_nss, 0.0);
266 
267  unsigned int num_data_grp = 3; // Number of data per group e.g. start_slip_sys, end_slip_sys,
268  // value
269 
270  for (unsigned int i = 0; i < _gprops.size() / num_data_grp; ++i)
271  {
272  Real vs, ve;
273  unsigned int is, ie;
274 
275  vs = _gprops[i * num_data_grp];
276  ve = _gprops[i * num_data_grp + 1];
277 
278  if (vs <= 0 || ve <= 0)
279  mooseError("FiniteStrainCrystalPLasticity: Indices in gss property read must be positive "
280  "integers: is = ",
281  vs,
282  " ie = ",
283  ve);
284 
285  if (vs != floor(vs) || ve != floor(ve))
286  mooseError("FiniteStrainCrystalPLasticity: Error in reading slip system resistances: Values "
287  "specifying start and end number of slip system groups should be integer");
288 
289  is = static_cast<unsigned int>(vs);
290  ie = static_cast<unsigned int>(ve);
291 
292  if (is > ie)
293  mooseError("FiniteStrainCrystalPLasticity: Start index is = ",
294  is,
295  " should be greater than end index ie = ",
296  ie,
297  " in slip system resistance property read");
298 
299  for (unsigned int j = is; j <= ie; ++j)
300  _gss[_qp][j - 1] = _gprops[i * num_data_grp + 2];
301  }
302 
303  for (unsigned int i = 0; i < _nss; ++i)
304  if (_gss[_qp][i] <= 0.0)
305  mooseError("FiniteStrainCrystalPLasticity: Value of resistance for slip system ",
306  i + 1,
307  " non positive");
308 }
const unsigned int _nss
Number of slip system resistance.
MaterialProperty< std::vector< Real > > & _gss

◆ getMatRot()

RankTwoTensor FiniteStrainCrystalPlasticity::getMatRot ( const RankTwoTensor a)
protected

This function perform RU decomposition to obtain the rotation tensor.

Definition at line 989 of file FiniteStrainCrystalPlasticity.C.

Referenced by get_current_rotation().

990 {
991  RankTwoTensor rot;
992  RankTwoTensor c, diag, evec;
993  PetscScalar cmat[LIBMESH_DIM][LIBMESH_DIM], work[10];
994  PetscReal w[LIBMESH_DIM];
995  PetscBLASInt nd = LIBMESH_DIM, lwork = 10, info;
996 
997  c = a.transpose() * a;
998 
999  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
1000  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
1001  cmat[i][j] = c(i, j);
1002 
1003  LAPACKsyev_("V", "U", &nd, &cmat[0][0], &nd, w, work, &lwork, &info);
1004 
1005  if (info != 0)
1006  mooseError("FiniteStrainCrystalPLasticity: DSYEV function call in getMatRot function failed");
1007 
1008  diag.zero();
1009 
1010  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
1011  diag(i, i) = std::sqrt(w[i]);
1012 
1013  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
1014  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
1015  evec(i, j) = cmat[i][j];
1016 
1017  rot = a * ((evec.transpose() * diag * evec).inverse());
1018 
1019  return rot;
1020 }

◆ getSlipIncrements()

void FiniteStrainCrystalPlasticity::getSlipIncrements ( )
protectedvirtual

This function updates the slip increments.

And derivative of slip w.r.t. resolved shear stress.

Reimplemented in FiniteStrainCPSlipRateRes.

Definition at line 958 of file FiniteStrainCrystalPlasticity.C.

Referenced by calcResidual(), and FiniteStrainCPSlipRateRes::getSlipIncrements().

959 {
960  for (unsigned int i = 0; i < _nss; ++i)
961  {
962  _slip_incr(i) = _a0(i) * std::pow(std::abs(_tau(i) / _gss_tmp[i]), 1.0 / _xm(i)) *
963  std::copysign(1.0, _tau(i)) * _dt;
964  if (std::abs(_slip_incr(i)) > _slip_incr_tol)
965  {
966  _err_tol = true;
967 #ifdef DEBUG
968  mooseWarning("Maximum allowable slip increment exceeded ", std::abs(_slip_incr(i)));
969 #endif
970  return;
971  }
972  }
973 
974  for (unsigned int i = 0; i < _nss; ++i)
975  _dslipdtau(i) = _a0(i) / _xm(i) *
976  std::pow(std::abs(_tau(i) / _gss_tmp[i]), 1.0 / _xm(i) - 1.0) / _gss_tmp[i] *
977  _dt;
978 }
const unsigned int _nss
Number of slip system resistance.
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
Real _slip_incr_tol
Slip increment tolerance.

◆ getSlipSystems()

void FiniteStrainCrystalPlasticity::getSlipSystems ( )
protectedvirtual

This function reads slip system from file - see test.

Definition at line 423 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity().

424 {
425  Real vec[LIBMESH_DIM];
426  std::ifstream fileslipsys;
427 
428  MooseUtils::checkFileReadable(_slip_sys_file_name);
429 
430  fileslipsys.open(_slip_sys_file_name.c_str());
431 
432  for (unsigned int i = 0; i < _nss; ++i)
433  {
434  // Read the slip normal
435  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
436  if (!(fileslipsys >> vec[j]))
437  mooseError("Crystal Plasticity Error: Premature end of file reading slip system file \n");
438 
439  // Normalize the vectors
440  Real mag;
441  mag = Utility::pow<2>(vec[0]) + Utility::pow<2>(vec[1]) + Utility::pow<2>(vec[2]);
442  mag = std::sqrt(mag);
443 
444  for (unsigned j = 0; j < LIBMESH_DIM; ++j)
445  _no(i * LIBMESH_DIM + j) = vec[j] / mag;
446 
447  // Read the slip direction
448  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
449  if (!(fileslipsys >> vec[j]))
450  mooseError("Crystal Plasticity Error: Premature end of file reading slip system file \n");
451 
452  // Normalize the vectors
453  mag = Utility::pow<2>(vec[0]) + Utility::pow<2>(vec[1]) + Utility::pow<2>(vec[2]);
454  mag = std::sqrt(mag);
455 
456  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
457  _mo(i * LIBMESH_DIM + j) = vec[j] / mag;
458 
459  mag = 0.0;
460  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
461  mag += _mo(i * LIBMESH_DIM + j) * _no(i * LIBMESH_DIM + j);
462 
463  if (std::abs(mag) > 1e-8)
464  mooseError(
465  "Crystal Plasicity Error: Slip direction and normal not orthonormal, System number = ",
466  i,
467  "\n");
468 
470  for (unsigned int j = 0; j < _num_slip_sys_props; ++j)
471  if (!(fileslipsys >> _slip_sys_props(i * _num_slip_sys_props + j)))
472  mooseError("Crystal Plasticity Error: Premature end of file reading slip system file - "
473  "check in slip system file read input options/values\n");
474  }
475 
476  fileslipsys.close();
477 }
const unsigned int _nss
Number of slip system resistance.
unsigned int _num_slip_sys_props
Number of slip system specific properties provided in the file containing slip system normals and dir...
std::string _slip_sys_file_name
File should contain slip plane normal and direction. See test.

◆ initAdditionalProps()

void FiniteStrainCrystalPlasticity::initAdditionalProps ( )
protectedvirtual

This function initializes additional parameters.

Definition at line 481 of file FiniteStrainCrystalPlasticity.C.

Referenced by initQpStatefulProperties().

482 {
483 }

◆ initQpStatefulProperties()

void FiniteStrainCrystalPlasticity::initQpStatefulProperties ( )
protectedvirtual

This function initializes the stateful properties such as stress, plastic deformation gradient, slip system resistances, etc.

Reimplemented from ComputeStressBase.

Definition at line 188 of file FiniteStrainCrystalPlasticity.C.

189 {
190  _stress[_qp].zero();
191 
192  _fp[_qp].setToIdentity();
193 
194  _pk2[_qp].zero();
195  _acc_slip[_qp] = 0.0;
196  _lag_e[_qp].zero();
197 
198  _update_rot[_qp].setToIdentity();
199 
200  initSlipSysProps(); // Initializes slip system related properties
202 }
MaterialProperty< RankTwoTensor > & _lag_e
MaterialProperty< RankTwoTensor > & _pk2
MaterialProperty< RankTwoTensor > & _stress
MaterialProperty< RankTwoTensor > & _fp
virtual void initSlipSysProps()
This function initializes slip system resistances.
virtual void initAdditionalProps()
This function initializes additional parameters.
MaterialProperty< RankTwoTensor > & _update_rot

◆ initSlipSysProps()

void FiniteStrainCrystalPlasticity::initSlipSysProps ( )
protectedvirtual

This function initializes slip system resistances.

Definition at line 205 of file FiniteStrainCrystalPlasticity.C.

Referenced by initQpStatefulProperties().

206 {
207  switch (_intvar_read_type)
208  {
209  case 0:
211  break;
212  case 1:
214  break;
215  default:
217  }
218 
219  if (_slip_sys_flow_prop_file_name.length() != 0)
221  else
223 
224  if (_slip_sys_hard_prop_file_name.length() != 0)
226  else
228 }
virtual void getHardnessParams()
This function assign flow rate parameters from .i file - see test.
std::string _slip_sys_flow_prop_file_name
File should contain values of the flow rate equation parameters.
virtual void readFileHardnessParams()
This function read hardness parameters from file.
std::string _slip_sys_hard_prop_file_name
The hardening parameters in this class are read from .i file. The user can override to read from file...
virtual void assignSlipSysRes()
This function assign initial values of slip system resistances/internal variables read from getSlipSy...
virtual void getFlowRateParams()
This function assign flow rate parameters - see test.
virtual void readFileFlowRateParams()
This function read flow rate parameters from file - see test.
virtual void getInitSlipSysRes()
This function assign slip system resistances - see test.
virtual void readFileInitSlipSysRes()
This function read slip system resistances from file - see test.
MooseEnum _intvar_read_type
Read from options for initial values of internal variables.

◆ internalVariableUpdateNRiteration()

void FiniteStrainCrystalPlasticity::internalVariableUpdateNRiteration ( )
protected

This function updates internal variables after each NewTon Raphson iteration (_fp_inv)

Definition at line 1206 of file FiniteStrainCrystalPlasticity.C.

Referenced by solveStress().

1207 {
1208  _fp_prev_inv = _fp_inv; // update _fp_prev_inv
1209 }

◆ line_search_update()

bool FiniteStrainCrystalPlasticity::line_search_update ( const Real  rnorm_prev,
const RankTwoTensor  dpk2 
)
protected

This function performs the line search update.

Definition at line 1123 of file FiniteStrainCrystalPlasticity.C.

Referenced by solveStress().

1124 {
1125  if (_lsrch_method == "CUT_HALF")
1126  {
1127  Real rnorm;
1128  RankTwoTensor resid;
1129  Real step = 1.0;
1130 
1131  do
1132  {
1133  _pk2_tmp = _pk2_tmp - step * dpk2;
1134  step /= 2.0;
1135  _pk2_tmp = _pk2_tmp + step * dpk2;
1136 
1137  calcResidual(resid);
1138  rnorm = resid.L2norm();
1139  } while (rnorm > rnorm_prev && step > _min_lsrch_step);
1140 
1141  if (rnorm > rnorm_prev && step <= _min_lsrch_step)
1142  return false;
1143 
1144  return true;
1145  }
1146  else if (_lsrch_method == "BISECTION")
1147  {
1148  unsigned int count = 0;
1149  Real step_a = 0.0;
1150  Real step_b = 1.0;
1151  Real step = 1.0;
1152  Real s_m = 1000.0;
1153  Real rnorm = 1000.0;
1154 
1155  RankTwoTensor resid;
1156  calcResidual(resid);
1157  Real s_b = resid.doubleContraction(dpk2);
1158  Real rnorm1 = resid.L2norm();
1159  _pk2_tmp = _pk2_tmp - dpk2;
1160  calcResidual(resid);
1161  Real s_a = resid.doubleContraction(dpk2);
1162  Real rnorm0 = resid.L2norm();
1163  _pk2_tmp = _pk2_tmp + dpk2;
1164 
1165  if ((rnorm1 / rnorm0) < _lsrch_tol || s_a * s_b > 0)
1166  {
1167  calcResidual(resid);
1168  return true;
1169  }
1170 
1171  while ((rnorm / rnorm0) > _lsrch_tol && count < _lsrch_max_iter)
1172  {
1173  _pk2_tmp = _pk2_tmp - step * dpk2;
1174  step = 0.5 * (step_b + step_a);
1175  _pk2_tmp = _pk2_tmp + step * dpk2;
1176  calcResidual(resid);
1177  s_m = resid.doubleContraction(dpk2);
1178  rnorm = resid.L2norm();
1179 
1180  if (s_m * s_a < 0.0)
1181  {
1182  step_b = step;
1183  s_b = s_m;
1184  }
1185  if (s_m * s_b < 0.0)
1186  {
1187  step_a = step;
1188  s_a = s_m;
1189  }
1190  count++;
1191  }
1192 
1193  if ((rnorm / rnorm0) < _lsrch_tol && count < _lsrch_max_iter)
1194  return true;
1195 
1196  return false;
1197  }
1198  else
1199  {
1200  mooseError("Line search meothod is not provided.");
1201  return false;
1202  }
1203 }
virtual void calcResidual(RankTwoTensor &)
This function calculate stress residual.
unsigned int _lsrch_max_iter
Line search bisection method maximum iteration number.
Real _lsrch_tol
Line search bisection method tolerance.
Real _min_lsrch_step
Minimum line search step size.

◆ postSolveQp()

void FiniteStrainCrystalPlasticity::postSolveQp ( )
protectedvirtual

This function update stress and internal variable after solve.

Definition at line 585 of file FiniteStrainCrystalPlasticity.C.

Referenced by computeQpStress().

586 {
587  if (_err_tol)
588  {
589  _err_tol = false;
591  {
593  _rndm_scale_var = _elasticity_tensor[_qp](0, 0, 0, 0);
594 
595  _stress[_qp] = RankTwoTensor::genRandomSymmTensor(_rndm_scale_var, 1.0);
596  }
597  else
598  mooseError("FiniteStrainCrystalPlasticity: Constitutive failure");
599  }
600  else
601  {
602  _stress[_qp] = _fe * _pk2[_qp] * _fe.transpose() / _fe.det();
603 
604  _Jacobian_mult[_qp] += calcTangentModuli(); // Calculate jacobian for preconditioner
605 
606  RankTwoTensor iden(RankTwoTensor::initIdentity);
607 
608  _lag_e[_qp] = _deformation_gradient[_qp].transpose() * _deformation_gradient[_qp] - iden;
609  _lag_e[_qp] = _lag_e[_qp] * 0.5;
610 
611  RankTwoTensor rot;
612  rot = get_current_rotation(_deformation_gradient[_qp]); // Calculate material rotation
613  _update_rot[_qp] = rot * _crysrot[_qp];
614  }
615 }
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
MaterialProperty< RankTwoTensor > & _lag_e
const MaterialProperty< RankTwoTensor > & _deformation_gradient
bool _input_rndm_scale_var
Input option for scaling variable to generate random stress when convergence fails.
const MaterialProperty< RankFourTensor > & _elasticity_tensor
MaterialProperty< RankTwoTensor > & _pk2
MaterialProperty< RankTwoTensor > & _stress
virtual RankFourTensor calcTangentModuli()
This function calculate the tangent moduli for preconditioner.
const MaterialProperty< RankTwoTensor > & _crysrot
MaterialProperty< RankTwoTensor > & _update_rot
RankTwoTensor get_current_rotation(const RankTwoTensor &a)
This function perform RU decomposition to obtain the rotation tensor.

◆ postSolveStatevar()

void FiniteStrainCrystalPlasticity::postSolveStatevar ( )
protectedvirtual

This function update internal variable after solve.

Definition at line 680 of file FiniteStrainCrystalPlasticity.C.

Referenced by solveQp().

681 {
682  if (_max_substep_iter == 1) // No substepping
683  {
684  _gss[_qp] = _gss_tmp;
685  _acc_slip[_qp] = _accslip_tmp;
686  }
687  else
688  {
689  if (_last_step_iter)
690  {
691  _gss[_qp] = _gss_tmp;
692  _acc_slip[_qp] = _accslip_tmp;
693  }
694  else
695  {
698  }
699  }
700 }
MaterialProperty< std::vector< Real > > & _gss
unsigned int _max_substep_iter
Maximum number of substep iterations.

◆ postSolveStress()

void FiniteStrainCrystalPlasticity::postSolveStress ( )
protectedvirtual

This function update stress and plastic deformation gradient after solve.

Definition at line 799 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCPSlipRateRes::solveStatevar(), and solveStatevar().

800 {
801  if (_max_substep_iter == 1) // No substepping
802  {
803  _fp[_qp] = _fp_inv.inverse();
804  _pk2[_qp] = _pk2_tmp;
805  }
806  else
807  {
808  if (_last_step_iter)
809  {
810  _fp[_qp] = _fp_inv.inverse();
811  _pk2[_qp] = _pk2_tmp;
812  }
813  else
814  {
817  }
818  }
819 }
MaterialProperty< RankTwoTensor > & _pk2
MaterialProperty< RankTwoTensor > & _fp
unsigned int _max_substep_iter
Maximum number of substep iterations.

◆ preSolveQp()

void FiniteStrainCrystalPlasticity::preSolveQp ( )
protectedvirtual

This function set variables for stress and internal variable solve.

Definition at line 557 of file FiniteStrainCrystalPlasticity.C.

Referenced by computeQpStress().

558 {
559  // Initialize variable
560  if (_first_substep)
561  {
562  _Jacobian_mult[_qp].zero(); // Initializes jacobian for preconditioner
564  }
565 
566  if (_max_substep_iter == 1)
567  _dfgrd_tmp = _deformation_gradient[_qp]; // Without substepping
568  else
570 
571  _err_tol = false;
572 }
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
const MaterialProperty< RankTwoTensor > & _deformation_gradient
Real _dfgrd_scale_factor
Scales the substepping increment to obtain deformation gradient at a substep iteration.
RankTwoTensor _delta_dfgrd
Flag to check whether convergence is achieved.
void calc_schmid_tensor()
This function calculate the Schmid tensor.
unsigned int _max_substep_iter
Maximum number of substep iterations.

◆ preSolveStatevar()

void FiniteStrainCrystalPlasticity::preSolveStatevar ( )
protectedvirtual

This function set variables for internal variable solve.

Definition at line 618 of file FiniteStrainCrystalPlasticity.C.

Referenced by solveQp().

619 {
620  if (_max_substep_iter == 1) // No substepping
621  {
622  _gss_tmp = _gss_old[_qp];
624  }
625  else
626  {
627  if (_first_step_iter)
628  {
629  _gss_tmp = _gss_tmp_old = _gss_old[_qp];
631  }
632  else
634  }
635 }
const MaterialProperty< std::vector< Real > > & _gss_old
unsigned int _max_substep_iter
Maximum number of substep iterations.
bool _first_step_iter
Flags to reset variables and reinitialize variables.
const MaterialProperty< Real > & _acc_slip_old

◆ preSolveStress()

void FiniteStrainCrystalPlasticity::preSolveStress ( )
protectedvirtual

This function set variables for stress solve.

Reimplemented in FiniteStrainCPSlipRateRes.

Definition at line 703 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCPSlipRateRes::preSolveStress(), and solveStatevar().

704 {
705  if (_max_substep_iter == 1) // No substepping
706  {
707  _pk2_tmp = _pk2_old[_qp];
708  _fp_old_inv = _fp_old[_qp].inverse();
711  }
712  else
713  {
714  if (_first_step_iter)
715  {
716  _pk2_tmp = _pk2_tmp_old = _pk2_old[_qp];
717  _fp_old_inv = _fp_old[_qp].inverse();
718  }
719  else
721 
724  }
725 }
const MaterialProperty< RankTwoTensor > & _pk2_old
const MaterialProperty< RankTwoTensor > & _fp_old
unsigned int _max_substep_iter
Maximum number of substep iterations.
bool _first_step_iter
Flags to reset variables and reinitialize variables.

◆ readFileFlowRateParams()

void FiniteStrainCrystalPlasticity::readFileFlowRateParams ( )
protectedvirtual

This function read flow rate parameters from file - see test.

Definition at line 312 of file FiniteStrainCrystalPlasticity.C.

Referenced by initSlipSysProps().

313 {
314  _a0.resize(_nss);
315  _xm.resize(_nss);
316 
317  MooseUtils::checkFileReadable(_slip_sys_flow_prop_file_name);
318 
319  std::ifstream file;
320  file.open(_slip_sys_flow_prop_file_name.c_str());
321 
322  std::vector<Real> vec;
323  vec.resize(_num_slip_sys_flowrate_props);
324 
325  for (unsigned int i = 0; i < _nss; ++i)
326  {
327  for (unsigned int j = 0; j < _num_slip_sys_flowrate_props; ++j)
328  if (!(file >> vec[j]))
329  mooseError(
330  "Error FiniteStrainCrystalPlasticity: Premature end of slip_sys_flow_rate_param file");
331 
332  _a0(i) = vec[0];
333  _xm(i) = vec[1];
334  }
335 
336  file.close();
337 }
std::string _slip_sys_flow_prop_file_name
File should contain values of the flow rate equation parameters.
unsigned int _num_slip_sys_flowrate_props
Number of slip system flow rate parameters.
const unsigned int _nss
Number of slip system resistance.

◆ readFileHardnessParams()

void FiniteStrainCrystalPlasticity::readFileHardnessParams ( )
protectedvirtual

This function read hardness parameters from file.

Definition at line 403 of file FiniteStrainCrystalPlasticity.C.

Referenced by initSlipSysProps().

404 {
405 }

◆ readFileInitSlipSysRes()

void FiniteStrainCrystalPlasticity::readFileInitSlipSysRes ( )
protectedvirtual

This function read slip system resistances from file - see test.

Definition at line 241 of file FiniteStrainCrystalPlasticity.C.

Referenced by initSlipSysProps().

242 {
243  _gss[_qp].resize(_nss);
244 
245  MooseUtils::checkFileReadable(_slip_sys_res_prop_file_name);
246 
247  std::ifstream file;
248  file.open(_slip_sys_res_prop_file_name.c_str());
249 
250  for (unsigned int i = 0; i < _nss; ++i)
251  if (!(file >> _gss[_qp][i]))
252  mooseError("Error FiniteStrainCrystalPlasticity: Premature end of slip_sys_res_prop file");
253 
254  file.close();
255 }
const unsigned int _nss
Number of slip system resistance.
MaterialProperty< std::vector< Real > > & _gss
std::string _slip_sys_res_prop_file_name
File should contain initial values of the slip system resistances.

◆ solveQp()

void FiniteStrainCrystalPlasticity::solveQp ( )
protectedvirtual

This function solves stress and internal variables.

Definition at line 575 of file FiniteStrainCrystalPlasticity.C.

Referenced by computeQpStress().

576 {
578  solveStatevar();
579  if (_err_tol)
580  return;
582 }
virtual void solveStatevar()
This function solves internal variables.
virtual void postSolveStatevar()
This function update internal variable after solve.
virtual void preSolveStatevar()
This function set variables for internal variable solve.

◆ solveStatevar()

void FiniteStrainCrystalPlasticity::solveStatevar ( )
protectedvirtual

This function solves internal variables.

Reimplemented in FiniteStrainCPSlipRateRes.

Definition at line 638 of file FiniteStrainCrystalPlasticity.C.

Referenced by solveQp().

639 {
640  Real gmax, gdiff;
641  unsigned int iterg;
642  std::vector<Real> gss_prev(_nss);
643 
644  gmax = 1.1 * _gtol;
645  iterg = 0;
646 
647  while (gmax > _gtol && iterg < _maxiterg) // Check for slip system resistance update tolerance
648  {
649  preSolveStress();
650  solveStress();
651  if (_err_tol)
652  return;
653  postSolveStress();
654 
655  gss_prev = _gss_tmp;
656 
657  update_slip_system_resistance(); // Update slip system resistance
658 
659  gmax = 0.0;
660  for (unsigned i = 0; i < _nss; ++i)
661  {
662  gdiff = std::abs(gss_prev[i] - _gss_tmp[i]); // Calculate increment size
663 
664  if (gdiff > gmax)
665  gmax = gdiff;
666  }
667  iterg++;
668  }
669 
670  if (iterg == _maxiterg)
671  {
672 #ifdef DEBUG
673  mooseWarning("FiniteStrainCrystalPLasticity: Hardness Integration error gmax", gmax, "\n");
674 #endif
675  _err_tol = true;
676  }
677 }
virtual void solveStress()
This function solves for stress, updates plastic deformation gradient.
virtual void update_slip_system_resistance()
This function updates the slip system resistances.
const unsigned int _nss
Number of slip system resistance.
Real _gtol
Internal variable update equation tolerance.
virtual void preSolveStress()
This function set variables for stress solve.
virtual void postSolveStress()
This function update stress and plastic deformation gradient after solve.
unsigned int _maxiterg
Maximum number of iterations for internal variable update.

◆ solveStress()

void FiniteStrainCrystalPlasticity::solveStress ( )
protectedvirtual

This function solves for stress, updates plastic deformation gradient.

Reimplemented in FiniteStrainCPSlipRateRes.

Definition at line 728 of file FiniteStrainCrystalPlasticity.C.

Referenced by solveStatevar().

729 {
730  unsigned int iter = 0;
731  RankTwoTensor resid, dpk2;
732  RankFourTensor jac;
733  Real rnorm, rnorm0, rnorm_prev;
734 
735  calc_resid_jacob(resid, jac); // Calculate stress residual
736  if (_err_tol)
737  {
738 #ifdef DEBUG
739  mooseWarning(
740  "FiniteStrainCrystalPLasticity: Slip increment exceeds tolerance - Element number ",
741  _current_elem->id(),
742  " Gauss point = ",
743  _qp);
744 #endif
745  return;
746  }
747 
748  rnorm = resid.L2norm();
749  rnorm0 = rnorm;
750 
751  while (rnorm > _rtol * rnorm0 && rnorm0 > _abs_tol &&
752  iter < _maxiter) // Check for stress residual tolerance
753  {
754  dpk2 = -jac.invSymm() * resid; // Calculate stress increment
755  _pk2_tmp = _pk2_tmp + dpk2; // Update stress
756  calc_resid_jacob(resid, jac);
757  internalVariableUpdateNRiteration(); // update _fp_prev_inv
758 
759  if (_err_tol)
760  {
761 #ifdef DEBUG
762  mooseWarning(
763  "FiniteStrainCrystalPLasticity: Slip increment exceeds tolerance - Element number ",
764  _current_elem->id(),
765  " Gauss point = ",
766  _qp);
767 #endif
768  return;
769  }
770 
771  rnorm_prev = rnorm;
772  rnorm = resid.L2norm();
773 
774  if (_use_line_search && rnorm > rnorm_prev && !line_search_update(rnorm_prev, dpk2))
775  {
776 #ifdef DEBUG
777  mooseWarning("FiniteStrainCrystalPLasticity: Failed with line search");
778 #endif
779  _err_tol = true;
780  return;
781  }
782 
783  if (_use_line_search)
784  rnorm = resid.L2norm();
785 
786  iter++;
787  }
788 
789  if (iter >= _maxiter)
790  {
791 #ifdef DEBUG
792  mooseWarning("FiniteStrainCrystalPLasticity: Stress Integration error rmax = ", rnorm);
793 #endif
794  _err_tol = true;
795  }
796 }
void internalVariableUpdateNRiteration()
This function updates internal variables after each NewTon Raphson iteration (_fp_inv) ...
bool _use_line_search
Flag to activate line serach.
Real _rtol
Stress residual equation relative tolerance.
unsigned int _maxiter
Maximum number of iterations for stress update.
bool line_search_update(const Real rnorm_prev, const RankTwoTensor)
This function performs the line search update.
virtual void calc_resid_jacob(RankTwoTensor &, RankFourTensor &)
This function calls the residual and jacobian functions used in the stress update algorithm...
Real _abs_tol
Stress residual equation absolute tolerance.

◆ update_slip_system_resistance()

void FiniteStrainCrystalPlasticity::update_slip_system_resistance ( )
protectedvirtual

This function updates the slip system resistances.

Definition at line 823 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCPSlipRateRes::calcResidualSlipRate(), and solveStatevar().

824 {
825  updateGss();
826 }
virtual void updateGss()
This function updates the slip system resistances.

◆ updateGss()

void FiniteStrainCrystalPlasticity::updateGss ( )
protectedvirtual

This function updates the slip system resistances.

Old function to update slip system resistances.

Kept to avoid code break at computeQpstress

Definition at line 833 of file FiniteStrainCrystalPlasticity.C.

Referenced by update_slip_system_resistance().

834 {
835  DenseVector<Real> hb(_nss);
836  Real qab;
837 
838  Real a = _hprops[4]; // Kalidindi
839 
841  for (unsigned int i = 0; i < _nss; ++i)
842  _accslip_tmp += std::abs(_slip_incr(i));
843 
844  // Real val = std::cosh(_h0 * _accslip_tmp / (_tau_sat - _tau_init)); // Karthik
845  // val = _h0 * std::pow(1.0/val,2.0); // Kalidindi
846 
847  for (unsigned int i = 0; i < _nss; ++i)
848  // hb(i)=val;
849  hb(i) = _h0 * std::pow(std::abs(1.0 - _gss_tmp[i] / _tau_sat), a) *
850  std::copysign(1.0, 1.0 - _gss_tmp[i] / _tau_sat);
851 
852  for (unsigned int i = 0; i < _nss; ++i)
853  {
854  if (_max_substep_iter == 1) // No substepping
855  _gss_tmp[i] = _gss_old[_qp][i];
856  else
857  _gss_tmp[i] = _gss_tmp_old[i];
858 
859  for (unsigned int j = 0; j < _nss; ++j)
860  {
861  unsigned int iplane, jplane;
862  iplane = i / 3;
863  jplane = j / 3;
864 
865  if (iplane == jplane) // Kalidindi
866  qab = 1.0;
867  else
868  qab = _r;
869 
870  _gss_tmp[i] += qab * hb(j) * std::abs(_slip_incr(j));
871  _dgss_dsliprate(i, j) = qab * hb(j) * std::copysign(1.0, _slip_incr(j)) * _dt;
872  }
873  }
874 }
const MaterialProperty< std::vector< Real > > & _gss_old
const unsigned int _nss
Number of slip system resistance.
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
unsigned int _max_substep_iter
Maximum number of substep iterations.

Member Data Documentation

◆ _a0

DenseVector<Real> FiniteStrainCrystalPlasticity::_a0
protected

◆ _abs_tol

Real FiniteStrainCrystalPlasticity::_abs_tol
protected

Stress residual equation absolute tolerance.

Definition at line 251 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCPSlipRateRes::solveStress(), and solveStress().

◆ _acc_slip

MaterialProperty<Real>& FiniteStrainCrystalPlasticity::_acc_slip
protected

Definition at line 311 of file FiniteStrainCrystalPlasticity.h.

Referenced by initQpStatefulProperties(), and postSolveStatevar().

◆ _acc_slip_old

const MaterialProperty<Real>& FiniteStrainCrystalPlasticity::_acc_slip_old
protected

Definition at line 312 of file FiniteStrainCrystalPlasticity.h.

Referenced by preSolveStatevar().

◆ _accslip_tmp

Real FiniteStrainCrystalPlasticity::_accslip_tmp
protected

Definition at line 337 of file FiniteStrainCrystalPlasticity.h.

Referenced by postSolveStatevar(), and updateGss().

◆ _accslip_tmp_old

Real FiniteStrainCrystalPlasticity::_accslip_tmp_old
protected

Definition at line 337 of file FiniteStrainCrystalPlasticity.h.

Referenced by postSolveStatevar(), preSolveStatevar(), and updateGss().

◆ _base_name

const std::string ComputeStressBase::_base_name
protectedinherited

◆ _crysrot

const MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_crysrot
protected

Definition at line 318 of file FiniteStrainCrystalPlasticity.h.

Referenced by calc_schmid_tensor(), and postSolveQp().

◆ _deformation_gradient

const MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_deformation_gradient
protected

Definition at line 315 of file FiniteStrainCrystalPlasticity.h.

Referenced by computeQpStress(), postSolveQp(), and preSolveQp().

◆ _deformation_gradient_old

const MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_deformation_gradient_old
protected

Definition at line 316 of file FiniteStrainCrystalPlasticity.h.

Referenced by computeQpStress().

◆ _delta_dfgrd

RankTwoTensor FiniteStrainCrystalPlasticity::_delta_dfgrd
protected

Flag to check whether convergence is achieved.

Used for substepping; Uniformly divides the increment in deformation gradient

Definition at line 350 of file FiniteStrainCrystalPlasticity.h.

Referenced by computeQpStress(), FiniteStrainCrystalPlasticity(), and preSolveQp().

◆ _dfgrd_scale_factor

Real FiniteStrainCrystalPlasticity::_dfgrd_scale_factor
protected

Scales the substepping increment to obtain deformation gradient at a substep iteration.

Definition at line 352 of file FiniteStrainCrystalPlasticity.h.

Referenced by computeQpStress(), and preSolveQp().

◆ _dfgrd_tmp

RankTwoTensor FiniteStrainCrystalPlasticity::_dfgrd_tmp
protected

◆ _dfgrd_tmp_old

RankTwoTensor FiniteStrainCrystalPlasticity::_dfgrd_tmp_old
protected

Definition at line 350 of file FiniteStrainCrystalPlasticity.h.

Referenced by computeQpStress(), and preSolveQp().

◆ _dgss_dsliprate

DenseMatrix<Real> FiniteStrainCrystalPlasticity::_dgss_dsliprate
protected

◆ _dslipdtau

DenseVector<Real> FiniteStrainCrystalPlasticity::_dslipdtau
protected

◆ _elastic_strain

MaterialProperty<RankTwoTensor>& ComputeStressBase::_elastic_strain
protectedinherited

◆ _elasticity_tensor

const MaterialProperty<RankFourTensor>& FiniteStrainCrystalPlasticity::_elasticity_tensor
protected

◆ _elasticity_tensor_name

const std::string ComputeStressBase::_elasticity_tensor_name
protectedinherited

◆ _err_tol

bool FiniteStrainCrystalPlasticity::_err_tol
protected

◆ _extra_stress

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_extra_stress
protectedinherited

Extra stress tensor.

Definition at line 47 of file ComputeStressBase.h.

Referenced by ComputeStressBase::computeQpProperties().

◆ _fe

RankTwoTensor FiniteStrainCrystalPlasticity::_fe
protected

◆ _first_step_iter

bool FiniteStrainCrystalPlasticity::_first_step_iter
protected

Flags to reset variables and reinitialize variables.

Definition at line 354 of file FiniteStrainCrystalPlasticity.h.

Referenced by computeQpStress(), FiniteStrainCrystalPlasticity(), preSolveStatevar(), and preSolveStress().

◆ _first_substep

bool FiniteStrainCrystalPlasticity::_first_substep
protected

◆ _flowprops

std::vector<Real> FiniteStrainCrystalPlasticity::_flowprops
protected

Definition at line 231 of file FiniteStrainCrystalPlasticity.h.

Referenced by getFlowRateParams().

◆ _fp

MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_fp
protected

Definition at line 303 of file FiniteStrainCrystalPlasticity.h.

Referenced by initQpStatefulProperties(), and postSolveStress().

◆ _fp_inv

RankTwoTensor FiniteStrainCrystalPlasticity::_fp_inv
protected

◆ _fp_old

const MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_fp_old
protected

Definition at line 304 of file FiniteStrainCrystalPlasticity.h.

Referenced by preSolveStress().

◆ _fp_old_inv

RankTwoTensor FiniteStrainCrystalPlasticity::_fp_old_inv
protected

◆ _fp_prev_inv

RankTwoTensor FiniteStrainCrystalPlasticity::_fp_prev_inv
protected

◆ _gen_rndm_stress_flag

bool FiniteStrainCrystalPlasticity::_gen_rndm_stress_flag
protected

Definition at line 274 of file FiniteStrainCrystalPlasticity.h.

Referenced by postSolveQp().

◆ _gprops

std::vector<Real> FiniteStrainCrystalPlasticity::_gprops
protected

Definition at line 229 of file FiniteStrainCrystalPlasticity.h.

Referenced by getInitSlipSysRes().

◆ _gss

MaterialProperty<std::vector<Real> >& FiniteStrainCrystalPlasticity::_gss
protected

◆ _gss_old

const MaterialProperty<std::vector<Real> >& FiniteStrainCrystalPlasticity::_gss_old
protected

Definition at line 310 of file FiniteStrainCrystalPlasticity.h.

Referenced by preSolveStatevar(), and updateGss().

◆ _gss_tmp

std::vector<Real> FiniteStrainCrystalPlasticity::_gss_tmp
protected

◆ _gss_tmp_old

std::vector<Real> FiniteStrainCrystalPlasticity::_gss_tmp_old
protected

Definition at line 339 of file FiniteStrainCrystalPlasticity.h.

Referenced by postSolveStatevar(), preSolveStatevar(), and updateGss().

◆ _gtol

Real FiniteStrainCrystalPlasticity::_gtol
protected

Internal variable update equation tolerance.

Definition at line 253 of file FiniteStrainCrystalPlasticity.h.

Referenced by solveStatevar().

◆ _h0

Real FiniteStrainCrystalPlasticity::_h0
protected

Definition at line 326 of file FiniteStrainCrystalPlasticity.h.

Referenced by getHardnessParams(), and updateGss().

◆ _hprops

std::vector<Real> FiniteStrainCrystalPlasticity::_hprops
protected

Definition at line 230 of file FiniteStrainCrystalPlasticity.h.

Referenced by getHardnessParams(), and updateGss().

◆ _initial_stress_fcn

std::vector<Function *> ComputeStressBase::_initial_stress_fcn
protectedinherited

initial stress components

Definition at line 50 of file ComputeStressBase.h.

◆ _input_rndm_scale_var

bool FiniteStrainCrystalPlasticity::_input_rndm_scale_var
protected

Input option for scaling variable to generate random stress when convergence fails.

Definition at line 277 of file FiniteStrainCrystalPlasticity.h.

Referenced by postSolveQp().

◆ _intvar_read_type

MooseEnum FiniteStrainCrystalPlasticity::_intvar_read_type
protected

Read from options for initial values of internal variables.

Definition at line 269 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity(), and initSlipSysProps().

◆ _Jacobian_mult

MaterialProperty<RankFourTensor>& ComputeStressBase::_Jacobian_mult
protectedinherited

◆ _lag_e

MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_lag_e
protected

Definition at line 307 of file FiniteStrainCrystalPlasticity.h.

Referenced by initQpStatefulProperties(), and postSolveQp().

◆ _lag_e_old

const MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_lag_e_old
protected

Definition at line 308 of file FiniteStrainCrystalPlasticity.h.

◆ _last_step_iter

bool FiniteStrainCrystalPlasticity::_last_step_iter
protected

◆ _lsrch_max_iter

unsigned int FiniteStrainCrystalPlasticity::_lsrch_max_iter
protected

Line search bisection method maximum iteration number.

Definition at line 298 of file FiniteStrainCrystalPlasticity.h.

Referenced by line_search_update(), and FiniteStrainCPSlipRateRes::lineSearchUpdateSlipRate().

◆ _lsrch_method

MooseEnum FiniteStrainCrystalPlasticity::_lsrch_method
protected

◆ _lsrch_tol

Real FiniteStrainCrystalPlasticity::_lsrch_tol
protected

Line search bisection method tolerance.

Definition at line 295 of file FiniteStrainCrystalPlasticity.h.

Referenced by line_search_update(), and FiniteStrainCPSlipRateRes::lineSearchUpdateSlipRate().

◆ _max_substep_iter

unsigned int FiniteStrainCrystalPlasticity::_max_substep_iter
protected

Maximum number of substep iterations.

Definition at line 286 of file FiniteStrainCrystalPlasticity.h.

Referenced by computeQpStress(), postSolveStatevar(), postSolveStress(), preSolveQp(), preSolveStatevar(), preSolveStress(), and updateGss().

◆ _maxiter

unsigned int FiniteStrainCrystalPlasticity::_maxiter
protected

Maximum number of iterations for stress update.

Definition at line 258 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCPSlipRateRes::solveStress(), and solveStress().

◆ _maxiterg

unsigned int FiniteStrainCrystalPlasticity::_maxiterg
protected

Maximum number of iterations for internal variable update.

Definition at line 260 of file FiniteStrainCrystalPlasticity.h.

Referenced by solveStatevar().

◆ _mechanical_strain

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_mechanical_strain
protectedinherited

◆ _min_lsrch_step

Real FiniteStrainCrystalPlasticity::_min_lsrch_step
protected

Minimum line search step size.

Definition at line 292 of file FiniteStrainCrystalPlasticity.h.

Referenced by line_search_update(), and FiniteStrainCPSlipRateRes::lineSearchUpdateSlipRate().

◆ _mo

DenseVector<Real> FiniteStrainCrystalPlasticity::_mo
protected

Definition at line 320 of file FiniteStrainCrystalPlasticity.h.

Referenced by calc_schmid_tensor(), and getSlipSystems().

◆ _no

DenseVector<Real> FiniteStrainCrystalPlasticity::_no
protected

Definition at line 321 of file FiniteStrainCrystalPlasticity.h.

Referenced by calc_schmid_tensor(), and getSlipSystems().

◆ _nss

const unsigned int FiniteStrainCrystalPlasticity::_nss
protected

◆ _num_slip_sys_flowrate_props

unsigned int FiniteStrainCrystalPlasticity::_num_slip_sys_flowrate_props
protected

Number of slip system flow rate parameters.

Definition at line 263 of file FiniteStrainCrystalPlasticity.h.

Referenced by getFlowRateParams(), and readFileFlowRateParams().

◆ _num_slip_sys_props

unsigned int FiniteStrainCrystalPlasticity::_num_slip_sys_props
protected

Number of slip system specific properties provided in the file containing slip system normals and directions.

Definition at line 272 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity(), and getSlipSystems().

◆ _pk2

MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_pk2
protected

◆ _pk2_old

const MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_pk2_old
protected

Definition at line 306 of file FiniteStrainCrystalPlasticity.h.

Referenced by preSolveStress().

◆ _pk2_tmp

RankTwoTensor FiniteStrainCrystalPlasticity::_pk2_tmp
protected

◆ _pk2_tmp_old

RankTwoTensor FiniteStrainCrystalPlasticity::_pk2_tmp_old
protected

Definition at line 336 of file FiniteStrainCrystalPlasticity.h.

Referenced by postSolveStress(), and preSolveStress().

◆ _r

Real FiniteStrainCrystalPlasticity::_r
protected

Definition at line 329 of file FiniteStrainCrystalPlasticity.h.

Referenced by getHardnessParams(), and updateGss().

◆ _read_from_slip_sys_file

bool FiniteStrainCrystalPlasticity::_read_from_slip_sys_file
protected

◆ _rndm_scale_var

Real FiniteStrainCrystalPlasticity::_rndm_scale_var
protected

Scaling value.

Definition at line 280 of file FiniteStrainCrystalPlasticity.h.

Referenced by postSolveQp().

◆ _rndm_seed

unsigned int FiniteStrainCrystalPlasticity::_rndm_seed
protected

Seed value.

Definition at line 283 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity().

◆ _rtol

Real FiniteStrainCrystalPlasticity::_rtol
protected

Stress residual equation relative tolerance.

Definition at line 249 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCPSlipRateRes::solveStress(), and solveStress().

◆ _s0

std::vector<RankTwoTensor> FiniteStrainCrystalPlasticity::_s0
protected

◆ _slip_incr

DenseVector<Real> FiniteStrainCrystalPlasticity::_slip_incr
protected

◆ _slip_incr_tol

Real FiniteStrainCrystalPlasticity::_slip_incr_tol
protected

Slip increment tolerance.

Definition at line 255 of file FiniteStrainCrystalPlasticity.h.

Referenced by getSlipIncrements().

◆ _slip_sys_file_name

std::string FiniteStrainCrystalPlasticity::_slip_sys_file_name
protected

File should contain slip plane normal and direction. See test.

Definition at line 234 of file FiniteStrainCrystalPlasticity.h.

Referenced by getSlipSystems().

◆ _slip_sys_flow_prop_file_name

std::string FiniteStrainCrystalPlasticity::_slip_sys_flow_prop_file_name
protected

File should contain values of the flow rate equation parameters.

Values for every slip system must be provided. Should have the same order of slip systens as in slip_sys_file. See test. The option of reading all the properties from .i is still present.

Definition at line 243 of file FiniteStrainCrystalPlasticity.h.

Referenced by initSlipSysProps(), and readFileFlowRateParams().

◆ _slip_sys_hard_prop_file_name

std::string FiniteStrainCrystalPlasticity::_slip_sys_hard_prop_file_name
protected

The hardening parameters in this class are read from .i file. The user can override to read from file.

Definition at line 246 of file FiniteStrainCrystalPlasticity.h.

Referenced by initSlipSysProps().

◆ _slip_sys_props

DenseVector<Real> FiniteStrainCrystalPlasticity::_slip_sys_props
protected

◆ _slip_sys_res_prop_file_name

std::string FiniteStrainCrystalPlasticity::_slip_sys_res_prop_file_name
protected

File should contain initial values of the slip system resistances.

Definition at line 237 of file FiniteStrainCrystalPlasticity.h.

Referenced by readFileInitSlipSysRes().

◆ _stress

MaterialProperty<RankTwoTensor>& ComputeStressBase::_stress
protectedinherited

◆ _tan_mod_type

MooseEnum FiniteStrainCrystalPlasticity::_tan_mod_type
protected

Type of tangent moduli calculation.

Definition at line 266 of file FiniteStrainCrystalPlasticity.h.

Referenced by calcTangentModuli().

◆ _tau

DenseVector<Real> FiniteStrainCrystalPlasticity::_tau
protected

◆ _tau_init

Real FiniteStrainCrystalPlasticity::_tau_init
protected

Definition at line 328 of file FiniteStrainCrystalPlasticity.h.

Referenced by getHardnessParams().

◆ _tau_sat

Real FiniteStrainCrystalPlasticity::_tau_sat
protected

Definition at line 327 of file FiniteStrainCrystalPlasticity.h.

Referenced by getHardnessParams(), and updateGss().

◆ _update_rot

MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_update_rot
protected

Definition at line 313 of file FiniteStrainCrystalPlasticity.h.

Referenced by initQpStatefulProperties(), and postSolveQp().

◆ _use_line_search

bool FiniteStrainCrystalPlasticity::_use_line_search
protected

Flag to activate line serach.

Definition at line 289 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCPSlipRateRes::solveStress(), and solveStress().

◆ _xm

DenseVector<Real> FiniteStrainCrystalPlasticity::_xm
protected

The documentation for this class was generated from the following files: