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 std::string _elasticity_tensor_name
 Name of the elasticity tensor material property. More...
 
const MaterialProperty< RankFourTensor > & _elasticity_tensor
 Elasticity tensor material property. More...
 
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
 Base name prepended to all material property names to allow for multi-material systems. More...
 
const MaterialProperty< RankTwoTensor > & _mechanical_strain
 Mechanical strain material property. More...
 
MaterialProperty< RankTwoTensor > & _stress
 Stress material property. More...
 
MaterialProperty< RankTwoTensor > & _elastic_strain
 Elastic strain material property. More...
 
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_name(_base_name + "elasticity_tensor"),
150  _elasticity_tensor(getMaterialPropertyByName<RankFourTensor>(_elasticity_tensor_name)),
151  _crysrot(getMaterialProperty<RankTwoTensor>("crysrot")),
152  _mo(_nss * LIBMESH_DIM),
153  _no(_nss * LIBMESH_DIM),
154  _slip_incr(_nss),
155  _tau(_nss),
156  _dslipdtau(_nss),
157  _s0(_nss),
158  _gss_tmp(_nss),
161 {
162  _err_tol = false;
163 
164  if (_num_slip_sys_props > 0)
166 
167  _pk2_tmp.zero();
168  _delta_dfgrd.zero();
169 
170  _first_step_iter = false;
171  _last_step_iter = false;
172  // Initialize variables in the first iteration of substepping
173  _first_substep = true;
174 
175  _read_from_slip_sys_file = false;
176  if (_intvar_read_type == "slip_sys_file")
178 
180  mooseError("Crystal Plasticity Error: Specify number of internal variable's initial values to "
181  "be read from slip system file");
182 
183  getSlipSystems();
184 
185  RankTwoTensor::initRandom(_rndm_seed);
186 }
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
Elasticity tensor material property.
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
const std::string _base_name
Base name prepended to all material property names to allow for multi-material systems.
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
const std::string _elasticity_tensor_name
Name of the elasticity tensor material property.
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 232 of file FiniteStrainCrystalPlasticity.C.

Referenced by initSlipSysProps().

233 {
234  _gss[_qp].resize(_nss);
235 
236  for (unsigned int i = 0; i < _nss; ++i)
237  _gss[_qp][i] = _slip_sys_props(i);
238 }
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 879 of file FiniteStrainCrystalPlasticity.C.

Referenced by solveStress().

880 {
881  calcResidual(resid);
882  if (_err_tol)
883  return;
884  calcJacobian(jac);
885 }
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 1047 of file FiniteStrainCrystalPlasticity.C.

Referenced by preSolveQp().

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

Referenced by calc_resid_jacob().

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

◆ calcResidual()

void FiniteStrainCrystalPlasticity::calcResidual ( RankTwoTensor resid)
protectedvirtual

This function calculate stress residual.

Definition at line 888 of file FiniteStrainCrystalPlasticity.C.

Referenced by calc_resid_jacob(), and line_search_update().

889 {
890  RankTwoTensor iden(RankTwoTensor::initIdentity), ce, ee, ce_pk2, eqv_slip_incr, pk2_new;
891 
892  _fe = _dfgrd_tmp * _fp_prev_inv; // _fp_inv ==> _fp_prev_inv
893 
894  ce = _fe.transpose() * _fe;
895  ce_pk2 = ce * _pk2_tmp;
896  ce_pk2 = ce_pk2 / _fe.det();
897 
898  // Calculate Schmid tensor and resolved shear stresses
899  for (unsigned int i = 0; i < _nss; ++i)
900  _tau(i) = ce_pk2.doubleContraction(_s0[i]);
901 
902  getSlipIncrements(); // Calculate dslip,dslipdtau
903 
904  if (_err_tol)
905  return;
906 
907  eqv_slip_incr.zero();
908  for (unsigned int i = 0; i < _nss; ++i)
909  eqv_slip_incr += _s0[i] * _slip_incr(i);
910 
911  eqv_slip_incr = iden - eqv_slip_incr;
912  _fp_inv = _fp_old_inv * eqv_slip_incr;
913  _fe = _dfgrd_tmp * _fp_inv;
914 
915  ce = _fe.transpose() * _fe;
916  ee = ce - iden;
917  ee *= 0.5;
918 
919  pk2_new = _elasticity_tensor[_qp] * ee;
920 
921  resid = _pk2_tmp - pk2_new;
922 }
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor material property.
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 1030 of file FiniteStrainCrystalPlasticity.C.

Referenced by postSolveQp().

1031 {
1032  RankFourTensor tan_mod;
1033 
1034  switch (_tan_mod_type)
1035  {
1036  case 0:
1037  tan_mod = elastoPlasticTangentModuli();
1038  break;
1039  default:
1040  tan_mod = elasticTangentModuli();
1041  }
1042 
1043  return tan_mod;
1044 }
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 1025 of file FiniteStrainCrystalPlasticity.C.

1026 {
1027 }

◆ computeQpProperties()

void ComputeStressBase::computeQpProperties ( )
overrideprotectedvirtualinherited

Definition at line 49 of file ComputeStressBase.C.

50 {
52 
53  // Add in extra stress
54  _stress[_qp] += _extra_stress[_qp];
55 }
virtual void computeQpStress()=0
Compute the stress and store it in the _stress material property for the current quadrature point...
MaterialProperty< RankTwoTensor > & _stress
Stress material property.
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 491 of file FiniteStrainCrystalPlasticity.C.

492 {
493  unsigned int substep_iter = 1; // Depth of substepping; Limited to maximum substep iteration
494  unsigned int num_substep = 1; // Calculated from substep_iter as 2^substep_iter
495  Real dt_original = _dt; // Stores original _dt; Reset at the end of solve
496  _first_substep = true; // Initialize variables at substep_iter = 1
497 
498  if (_max_substep_iter > 1)
499  {
501  if (_dfgrd_tmp_old.det() == 0)
502  _dfgrd_tmp_old.addIa(1.0);
503 
505  _err_tol = true; // Indicator to continue substepping
506  }
507 
508  // Substepping loop
509  while (_err_tol && _max_substep_iter > 1)
510  {
511  _dt = dt_original / num_substep;
512 
513  for (unsigned int istep = 0; istep < num_substep; ++istep)
514  {
515  _first_step_iter = false;
516  if (istep == 0)
517  _first_step_iter = true;
518 
519  _last_step_iter = false;
520  if (istep == num_substep - 1)
521  _last_step_iter = true;
522 
523  _dfgrd_scale_factor = (static_cast<Real>(istep) + 1) / num_substep;
524 
525  preSolveQp();
526  solveQp();
527 
528  if (_err_tol)
529  {
530  substep_iter++;
531  num_substep *= 2;
532  break;
533  }
534  }
535 
536  _first_substep = false; // Prevents reinitialization
537  _dt = dt_original; // Resets dt
538 
539 #ifdef DEBUG
540  if (substep_iter > _max_substep_iter)
541  mooseWarning("FiniteStrainCrystalPlasticity: Failure with substepping");
542 #endif
543 
544  if (!_err_tol || substep_iter > _max_substep_iter)
545  postSolveQp(); // Evaluate variables after successful solve or indicate failure
546  }
547 
548  // No substepping
549  if (_max_substep_iter == 1)
550  {
551  preSolveQp();
552  solveQp();
553  postSolveQp();
554  }
555 }
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 1118 of file FiniteStrainCrystalPlasticity.C.

Referenced by calcTangentModuli().

1119 {
1120  return _elasticity_tensor[_qp]; // update jacobian_mult
1121 }
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor material property.

◆ elastoPlasticTangentModuli()

RankFourTensor FiniteStrainCrystalPlasticity::elastoPlasticTangentModuli ( )
protectedvirtual

This function calculate the exact tangent moduli for preconditioner.

Definition at line 1079 of file FiniteStrainCrystalPlasticity.C.

Referenced by calcTangentModuli().

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

◆ 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 983 of file FiniteStrainCrystalPlasticity.C.

Referenced by postSolveQp().

984 {
985  return getMatRot(a);
986 }
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 342 of file FiniteStrainCrystalPlasticity.C.

Referenced by initSlipSysProps().

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

Referenced by initSlipSysProps().

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

◆ 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 260 of file FiniteStrainCrystalPlasticity.C.

Referenced by initSlipSysProps().

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

Referenced by get_current_rotation().

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

◆ 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 959 of file FiniteStrainCrystalPlasticity.C.

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

960 {
961  for (unsigned int i = 0; i < _nss; ++i)
962  {
963  _slip_incr(i) = _a0(i) * std::pow(std::abs(_tau(i) / _gss_tmp[i]), 1.0 / _xm(i)) *
964  std::copysign(1.0, _tau(i)) * _dt;
965  if (std::abs(_slip_incr(i)) > _slip_incr_tol)
966  {
967  _err_tol = true;
968 #ifdef DEBUG
969  mooseWarning("Maximum allowable slip increment exceeded ", std::abs(_slip_incr(i)));
970 #endif
971  return;
972  }
973  }
974 
975  for (unsigned int i = 0; i < _nss; ++i)
976  _dslipdtau(i) = _a0(i) / _xm(i) *
977  std::pow(std::abs(_tau(i) / _gss_tmp[i]), 1.0 / _xm(i) - 1.0) / _gss_tmp[i] *
978  _dt;
979 }
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 424 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity().

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

Referenced by initQpStatefulProperties().

483 {
484 }

◆ 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 189 of file FiniteStrainCrystalPlasticity.C.

190 {
191  _stress[_qp].zero();
192 
193  _fp[_qp].setToIdentity();
194 
195  _pk2[_qp].zero();
196  _acc_slip[_qp] = 0.0;
197  _lag_e[_qp].zero();
198 
199  _update_rot[_qp].setToIdentity();
200 
201  initSlipSysProps(); // Initializes slip system related properties
203 }
MaterialProperty< RankTwoTensor > & _lag_e
MaterialProperty< RankTwoTensor > & _pk2
MaterialProperty< RankTwoTensor > & _stress
Stress material property.
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 206 of file FiniteStrainCrystalPlasticity.C.

Referenced by initQpStatefulProperties().

207 {
208  switch (_intvar_read_type)
209  {
210  case 0:
212  break;
213  case 1:
215  break;
216  default:
218  }
219 
220  if (_slip_sys_flow_prop_file_name.length() != 0)
222  else
224 
225  if (_slip_sys_hard_prop_file_name.length() != 0)
227  else
229 }
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 1207 of file FiniteStrainCrystalPlasticity.C.

Referenced by solveStress().

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

◆ 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 1124 of file FiniteStrainCrystalPlasticity.C.

Referenced by solveStress().

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

Referenced by computeQpStress().

587 {
588  if (_err_tol)
589  {
590  _err_tol = false;
592  {
594  _rndm_scale_var = _elasticity_tensor[_qp](0, 0, 0, 0);
595 
596  _stress[_qp] = RankTwoTensor::genRandomSymmTensor(_rndm_scale_var, 1.0);
597  }
598  else
599  mooseError("FiniteStrainCrystalPlasticity: Constitutive failure");
600  }
601  else
602  {
603  _stress[_qp] = _fe * _pk2[_qp] * _fe.transpose() / _fe.det();
604 
605  _Jacobian_mult[_qp] += calcTangentModuli(); // Calculate jacobian for preconditioner
606 
607  RankTwoTensor iden(RankTwoTensor::initIdentity);
608 
609  _lag_e[_qp] = _deformation_gradient[_qp].transpose() * _deformation_gradient[_qp] - iden;
610  _lag_e[_qp] = _lag_e[_qp] * 0.5;
611 
612  RankTwoTensor rot;
613  rot = get_current_rotation(_deformation_gradient[_qp]); // Calculate material rotation
614  _update_rot[_qp] = rot * _crysrot[_qp];
615  }
616 }
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
Elasticity tensor material property.
MaterialProperty< RankTwoTensor > & _pk2
MaterialProperty< RankTwoTensor > & _stress
Stress material property.
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 681 of file FiniteStrainCrystalPlasticity.C.

Referenced by solveQp().

682 {
683  if (_max_substep_iter == 1) // No substepping
684  {
685  _gss[_qp] = _gss_tmp;
686  _acc_slip[_qp] = _accslip_tmp;
687  }
688  else
689  {
690  if (_last_step_iter)
691  {
692  _gss[_qp] = _gss_tmp;
693  _acc_slip[_qp] = _accslip_tmp;
694  }
695  else
696  {
699  }
700  }
701 }
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 800 of file FiniteStrainCrystalPlasticity.C.

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

801 {
802  if (_max_substep_iter == 1) // No substepping
803  {
804  _fp[_qp] = _fp_inv.inverse();
805  _pk2[_qp] = _pk2_tmp;
806  }
807  else
808  {
809  if (_last_step_iter)
810  {
811  _fp[_qp] = _fp_inv.inverse();
812  _pk2[_qp] = _pk2_tmp;
813  }
814  else
815  {
818  }
819  }
820 }
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 558 of file FiniteStrainCrystalPlasticity.C.

Referenced by computeQpStress().

559 {
560  // Initialize variable
561  if (_first_substep)
562  {
563  _Jacobian_mult[_qp].zero(); // Initializes jacobian for preconditioner
565  }
566 
567  if (_max_substep_iter == 1)
568  _dfgrd_tmp = _deformation_gradient[_qp]; // Without substepping
569  else
571 
572  _err_tol = false;
573 }
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 619 of file FiniteStrainCrystalPlasticity.C.

Referenced by solveQp().

620 {
621  if (_max_substep_iter == 1) // No substepping
622  {
623  _gss_tmp = _gss_old[_qp];
625  }
626  else
627  {
628  if (_first_step_iter)
629  {
630  _gss_tmp = _gss_tmp_old = _gss_old[_qp];
632  }
633  else
635  }
636 }
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 704 of file FiniteStrainCrystalPlasticity.C.

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

705 {
706  if (_max_substep_iter == 1) // No substepping
707  {
708  _pk2_tmp = _pk2_old[_qp];
709  _fp_old_inv = _fp_old[_qp].inverse();
712  }
713  else
714  {
715  if (_first_step_iter)
716  {
717  _pk2_tmp = _pk2_tmp_old = _pk2_old[_qp];
718  _fp_old_inv = _fp_old[_qp].inverse();
719  }
720  else
722 
725  }
726 }
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 313 of file FiniteStrainCrystalPlasticity.C.

Referenced by initSlipSysProps().

314 {
315  _a0.resize(_nss);
316  _xm.resize(_nss);
317 
318  MooseUtils::checkFileReadable(_slip_sys_flow_prop_file_name);
319 
320  std::ifstream file;
321  file.open(_slip_sys_flow_prop_file_name.c_str());
322 
323  std::vector<Real> vec;
324  vec.resize(_num_slip_sys_flowrate_props);
325 
326  for (unsigned int i = 0; i < _nss; ++i)
327  {
328  for (unsigned int j = 0; j < _num_slip_sys_flowrate_props; ++j)
329  if (!(file >> vec[j]))
330  mooseError(
331  "Error FiniteStrainCrystalPlasticity: Premature end of slip_sys_flow_rate_param file");
332 
333  _a0(i) = vec[0];
334  _xm(i) = vec[1];
335  }
336 
337  file.close();
338 }
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 404 of file FiniteStrainCrystalPlasticity.C.

Referenced by initSlipSysProps().

405 {
406 }

◆ readFileInitSlipSysRes()

void FiniteStrainCrystalPlasticity::readFileInitSlipSysRes ( )
protectedvirtual

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

Definition at line 242 of file FiniteStrainCrystalPlasticity.C.

Referenced by initSlipSysProps().

243 {
244  _gss[_qp].resize(_nss);
245 
246  MooseUtils::checkFileReadable(_slip_sys_res_prop_file_name);
247 
248  std::ifstream file;
249  file.open(_slip_sys_res_prop_file_name.c_str());
250 
251  for (unsigned int i = 0; i < _nss; ++i)
252  if (!(file >> _gss[_qp][i]))
253  mooseError("Error FiniteStrainCrystalPlasticity: Premature end of slip_sys_res_prop file");
254 
255  file.close();
256 }
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 576 of file FiniteStrainCrystalPlasticity.C.

Referenced by computeQpStress().

577 {
579  solveStatevar();
580  if (_err_tol)
581  return;
583 }
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 639 of file FiniteStrainCrystalPlasticity.C.

Referenced by solveQp().

640 {
641  Real gmax, gdiff;
642  unsigned int iterg;
643  std::vector<Real> gss_prev(_nss);
644 
645  gmax = 1.1 * _gtol;
646  iterg = 0;
647 
648  while (gmax > _gtol && iterg < _maxiterg) // Check for slip system resistance update tolerance
649  {
650  preSolveStress();
651  solveStress();
652  if (_err_tol)
653  return;
654  postSolveStress();
655 
656  gss_prev = _gss_tmp;
657 
658  update_slip_system_resistance(); // Update slip system resistance
659 
660  gmax = 0.0;
661  for (unsigned i = 0; i < _nss; ++i)
662  {
663  gdiff = std::abs(gss_prev[i] - _gss_tmp[i]); // Calculate increment size
664 
665  if (gdiff > gmax)
666  gmax = gdiff;
667  }
668  iterg++;
669  }
670 
671  if (iterg == _maxiterg)
672  {
673 #ifdef DEBUG
674  mooseWarning("FiniteStrainCrystalPLasticity: Hardness Integration error gmax", gmax, "\n");
675 #endif
676  _err_tol = true;
677  }
678 }
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 729 of file FiniteStrainCrystalPlasticity.C.

Referenced by solveStatevar().

730 {
731  unsigned int iter = 0;
732  RankTwoTensor resid, dpk2;
733  RankFourTensor jac;
734  Real rnorm, rnorm0, rnorm_prev;
735 
736  calc_resid_jacob(resid, jac); // Calculate stress residual
737  if (_err_tol)
738  {
739 #ifdef DEBUG
740  mooseWarning(
741  "FiniteStrainCrystalPLasticity: Slip increment exceeds tolerance - Element number ",
742  _current_elem->id(),
743  " Gauss point = ",
744  _qp);
745 #endif
746  return;
747  }
748 
749  rnorm = resid.L2norm();
750  rnorm0 = rnorm;
751 
752  while (rnorm > _rtol * rnorm0 && rnorm0 > _abs_tol &&
753  iter < _maxiter) // Check for stress residual tolerance
754  {
755  dpk2 = -jac.invSymm() * resid; // Calculate stress increment
756  _pk2_tmp = _pk2_tmp + dpk2; // Update stress
757  calc_resid_jacob(resid, jac);
758  internalVariableUpdateNRiteration(); // update _fp_prev_inv
759 
760  if (_err_tol)
761  {
762 #ifdef DEBUG
763  mooseWarning(
764  "FiniteStrainCrystalPLasticity: Slip increment exceeds tolerance - Element number ",
765  _current_elem->id(),
766  " Gauss point = ",
767  _qp);
768 #endif
769  return;
770  }
771 
772  rnorm_prev = rnorm;
773  rnorm = resid.L2norm();
774 
775  if (_use_line_search && rnorm > rnorm_prev && !line_search_update(rnorm_prev, dpk2))
776  {
777 #ifdef DEBUG
778  mooseWarning("FiniteStrainCrystalPLasticity: Failed with line search");
779 #endif
780  _err_tol = true;
781  return;
782  }
783 
784  if (_use_line_search)
785  rnorm = resid.L2norm();
786 
787  iter++;
788  }
789 
790  if (iter >= _maxiter)
791  {
792 #ifdef DEBUG
793  mooseWarning("FiniteStrainCrystalPLasticity: Stress Integration error rmax = ", rnorm);
794 #endif
795  _err_tol = true;
796  }
797 }
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 824 of file FiniteStrainCrystalPlasticity.C.

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

825 {
826  updateGss();
827 }
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 834 of file FiniteStrainCrystalPlasticity.C.

Referenced by update_slip_system_resistance().

835 {
836  DenseVector<Real> hb(_nss);
837  Real qab;
838 
839  Real a = _hprops[4]; // Kalidindi
840 
842  for (unsigned int i = 0; i < _nss; ++i)
843  _accslip_tmp += std::abs(_slip_incr(i));
844 
845  // Real val = std::cosh(_h0 * _accslip_tmp / (_tau_sat - _tau_init)); // Karthik
846  // val = _h0 * std::pow(1.0/val,2.0); // Kalidindi
847 
848  for (unsigned int i = 0; i < _nss; ++i)
849  // hb(i)=val;
850  hb(i) = _h0 * std::pow(std::abs(1.0 - _gss_tmp[i] / _tau_sat), a) *
851  std::copysign(1.0, 1.0 - _gss_tmp[i] / _tau_sat);
852 
853  for (unsigned int i = 0; i < _nss; ++i)
854  {
855  if (_max_substep_iter == 1) // No substepping
856  _gss_tmp[i] = _gss_old[_qp][i];
857  else
858  _gss_tmp[i] = _gss_tmp_old[i];
859 
860  for (unsigned int j = 0; j < _nss; ++j)
861  {
862  unsigned int iplane, jplane;
863  iplane = i / 3;
864  jplane = j / 3;
865 
866  if (iplane == jplane) // Kalidindi
867  qab = 1.0;
868  else
869  qab = _r;
870 
871  _gss_tmp[i] += qab * hb(j) * std::abs(_slip_incr(j));
872  _dgss_dsliprate(i, j) = qab * hb(j) * std::copysign(1.0, _slip_incr(j)) * _dt;
873  }
874  }
875 }
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 340 of file FiniteStrainCrystalPlasticity.h.

Referenced by postSolveStatevar(), and updateGss().

◆ _accslip_tmp_old

Real FiniteStrainCrystalPlasticity::_accslip_tmp_old
protected

Definition at line 340 of file FiniteStrainCrystalPlasticity.h.

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

◆ _base_name

const std::string ComputeStressBase::_base_name
protectedinherited

Base name prepended to all material property names to allow for multi-material systems.

Definition at line 44 of file ComputeStressBase.h.

Referenced by ComputeLinearElasticStress::initialSetup(), and ComputeCosseratLinearElasticStress::initialSetup().

◆ _crysrot

const MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_crysrot
protected

Definition at line 321 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 353 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 355 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 353 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 FiniteStrainCrystalPlasticity::_elasticity_tensor_name
protected

Name of the elasticity tensor material property.

Definition at line 318 of file FiniteStrainCrystalPlasticity.h.

◆ _err_tol

bool FiniteStrainCrystalPlasticity::_err_tol
protected

◆ _extra_stress

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_extra_stress
protectedinherited

Extra stress tensor.

Definition at line 54 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 357 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 342 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 329 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 57 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 323 of file FiniteStrainCrystalPlasticity.h.

Referenced by calc_schmid_tensor(), and getSlipSystems().

◆ _no

DenseVector<Real> FiniteStrainCrystalPlasticity::_no
protected

Definition at line 324 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 339 of file FiniteStrainCrystalPlasticity.h.

Referenced by postSolveStress(), and preSolveStress().

◆ _r

Real FiniteStrainCrystalPlasticity::_r
protected

Definition at line 332 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

Stress material property.

Definition at line 49 of file ComputeStressBase.h.

Referenced by ComputeMultipleInelasticCosseratStress::computeAdmissibleState(), ComputeMultipleInelasticStress::computeAdmissibleState(), ComputeStressBase::computeQpProperties(), ComputeStrainIncrementBasedStress::computeQpStress(), ComputeLinearElasticStress::computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeDamageStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeSmearedCrackingStress::computeQpStress(), FiniteStrainPlasticMaterial::computeQpStress(), ComputeLinearElasticPFFractureStress::computeQpStress(), ComputeMultiPlasticityStress::computeQpStress(), ComputeLinearViscoelasticStress::computeQpStress(), ComputeMultipleInelasticStress::computeQpStressIntermediateConfiguration(), ComputeLinearElasticPFFractureStress::computeStrainSpectral(), ComputeLinearElasticPFFractureStress::computeStrainVolDev(), ComputeLinearElasticPFFractureStress::computeStressSpectral(), ComputeMultipleInelasticStress::finiteStrainRotation(), ComputeStressBase::initQpStatefulProperties(), initQpStatefulProperties(), FiniteStrainUObasedCP::initQpStatefulProperties(), FiniteStrainHyperElasticViscoPlastic::initQpStatefulProperties(), ComputeMultiPlasticityStress::postReturnMap(), FiniteStrainUObasedCP::postSolveQp(), FiniteStrainHyperElasticViscoPlastic::postSolveQp(), postSolveQp(), ComputeSmearedCrackingStress::updateCrackingStateAndStress(), ComputeMultipleInelasticStress::updateQpState(), and ComputeMultipleInelasticStress::updateQpStateSingleModel().

◆ _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 331 of file FiniteStrainCrystalPlasticity.h.

Referenced by getHardnessParams().

◆ _tau_sat

Real FiniteStrainCrystalPlasticity::_tau_sat
protected

Definition at line 330 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: