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

FiniteStrainUObasedCP 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 <FiniteStrainUObasedCP.h>

Inheritance diagram for FiniteStrainUObasedCP:
[legend]

Public Member Functions

 FiniteStrainUObasedCP (const InputParameters &parameters)
 

Protected Member Functions

virtual void computeQpStress ()
 updates the stress at a quadrature point. More...
 
virtual void initQpStatefulProperties ()
 initializes the stateful properties such as stress, plastic deformation gradient, slip system resistances, etc. More...
 
virtual void calcResidJacob ()
 calls the residual and jacobian functions used in the stress update algorithm. More...
 
virtual void updateSlipSystemResistanceAndStateVariable ()
 updates the slip system resistances and state variables. More...
 
virtual void preSolveQp ()
 set variables for stress and internal variable solve. More...
 
virtual void solveQp ()
 solve stress and internal variables. More...
 
virtual void postSolveQp ()
 update stress and internal variable after solve. More...
 
virtual void preSolveStatevar ()
 set variables for internal variable solve. More...
 
virtual void solveStatevar ()
 solve internal variables. More...
 
virtual void postSolveStatevar ()
 update internal variable after solve. More...
 
virtual void preSolveStress ()
 set variables for stress solve. More...
 
virtual void solveStress ()
 solves for stress, updates plastic deformation gradient. More...
 
virtual void postSolveStress ()
 update stress and plastic deformation gradient after solve. More...
 
virtual void calcResidual ()
 calculate stress residual. More...
 
virtual void calcJacobian ()
 calculate jacobian. More...
 
virtual void getSlipRates ()
 updates the slip rates. More...
 
virtual void calcTangentModuli ()
 calculate the tangent moduli for preconditioner. More...
 
virtual void elasticTangentModuli ()
 calculate the elastic tangent moduli for preconditioner. More...
 
virtual void elastoPlasticTangentModuli ()
 calculate the exact tangent moduli for preconditioner. More...
 
bool lineSearchUpdate (const Real rnorm_prev, const RankTwoTensor)
 performs the line search update More...
 
virtual bool isStateVariablesConverged ()
 evaluates convergence of state variables. More...
 
virtual void computeQpProperties () override
 

Protected Attributes

std::vector< const CrystalPlasticitySlipRate * > _uo_slip_rates
 User objects that define the slip rate. More...
 
std::vector< const CrystalPlasticitySlipResistance * > _uo_slip_resistances
 User objects that define the slip resistance. More...
 
std::vector< const CrystalPlasticityStateVariable * > _uo_state_vars
 User objects that define the state variable. More...
 
std::vector< const CrystalPlasticityStateVarRateComponent * > _uo_state_var_evol_rate_comps
 User objects that define the state variable evolution rate component. More...
 
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_slip_rates
 Slip rates material property. More...
 
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_slip_resistances
 Slip resistance material property. More...
 
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars
 State variable material property. More...
 
std::vector< const MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars_old
 Old state variable material property. More...
 
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_var_evol_rate_comps
 State variable evolution rate component material property. More...
 
unsigned int _num_uo_slip_rates
 Number of slip rate user objects. More...
 
unsigned int _num_uo_slip_resistances
 Number of slip resistance user objects. More...
 
unsigned int _num_uo_state_vars
 Number of state variable user objects. More...
 
unsigned int _num_uo_state_var_evol_rate_comps
 Number of state variable evolution rate component user objects. More...
 
std::vector< std::vector< Real > > _state_vars_old
 Local state variable. More...
 
std::vector< std::vector< Real > > _state_vars_old_stored
 Local stored state variable (for sub-stepping) More...
 
std::vector< std::vector< Real > > _state_vars_prev
 Local old state variable. More...
 
Real _rtol
 Stress residual equation relative tolerance. More...
 
Real _abs_tol
 Stress residual equation absolute tolerance. More...
 
Real _stol
 Internal variable update equation tolerance. More...
 
Real _zero_tol
 Residual tolerance when variable value is zero. Default 1e-12. More...
 
RankTwoTensor _resid
 Residual tensor. More...
 
RankFourTensor _jac
 Jacobian tensor. More...
 
unsigned int _maxiter
 Maximum number of iterations for stress update. More...
 
unsigned int _maxiterg
 Maximum number of iterations for internal variable update. More...
 
MooseEnum _tan_mod_type
 Type of tangent moduli calculation. 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
 Line search method. More...
 
MaterialProperty< RankTwoTensor > & _fp
 
const MaterialProperty< RankTwoTensor > & _fp_old
 
MaterialProperty< RankTwoTensor > & _pk2
 
const MaterialProperty< RankTwoTensor > & _pk2_old
 
MaterialProperty< RankTwoTensor > & _lag_e
 
MaterialProperty< RankTwoTensor > & _update_rot
 
const MaterialProperty< RankTwoTensor > & _update_rot_old
 
const MaterialProperty< RankTwoTensor > & _deformation_gradient
 
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
 
const MaterialProperty< RankTwoTensor > & _crysrot
 Crystal rotation. More...
 
RankTwoTensor _dfgrd_tmp
 
RankTwoTensor _fe
 
RankTwoTensor _fp_old_inv
 
RankTwoTensor _fp_inv
 
DenseVector< Real > _tau
 
std::vector< MaterialProperty< std::vector< RankTwoTensor > > * > _flow_direction
 
bool _err_tol
 Flag to check whether convergence is achieved. More...
 
RankTwoTensor _delta_dfgrd
 Used for substepping; Uniformly divides the increment in deformation gradient. More...
 
RankTwoTensor _dfgrd_tmp_old
 
Real _dfgrd_scale_factor
 Scales the substepping increment to obtain deformation gradient at a substep iteration. More...
 
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< RankFourTensor > & _elasticity_tensor
 
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

FiniteStrainUObasedCP 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.

Involves 4 different types of user objects that calculates: State variables - update state variable (derive from CrystalPlasticityStateVariable) State variable evolution compoment - individual component of state variable incremental rate (derive from CrystalPlasticityStateVariableEvolutionRateComponent) Slip resistance - calcuate slip resistances (derive from CrystalPlasticitySlipResistances) Slip rates - calcuate flow direction and slip rates (derive from CrystalPlasticitySlipRates)

Definition at line 39 of file FiniteStrainUObasedCP.h.

Constructor & Destructor Documentation

◆ FiniteStrainUObasedCP()

FiniteStrainUObasedCP::FiniteStrainUObasedCP ( const InputParameters &  parameters)

Definition at line 65 of file FiniteStrainUObasedCP.C.

66  : ComputeStressBase(parameters),
67  _num_uo_slip_rates(parameters.get<std::vector<UserObjectName>>("uo_slip_rates").size()),
69  parameters.get<std::vector<UserObjectName>>("uo_slip_resistances").size()),
70  _num_uo_state_vars(parameters.get<std::vector<UserObjectName>>("uo_state_vars").size()),
72  parameters.get<std::vector<UserObjectName>>("uo_state_var_evol_rate_comps").size()),
73  _rtol(getParam<Real>("rtol")),
74  _abs_tol(getParam<Real>("abs_tol")),
75  _stol(getParam<Real>("stol")),
76  _zero_tol(getParam<Real>("zero_tol")),
77  _maxiter(getParam<unsigned int>("maxiter")),
78  _maxiterg(getParam<unsigned int>("maxiter_state_variable")),
79  _tan_mod_type(getParam<MooseEnum>("tan_mod_type")),
80  _max_substep_iter(getParam<unsigned int>("maximum_substep_iteration")),
81  _use_line_search(getParam<bool>("use_line_search")),
82  _min_lsrch_step(getParam<Real>("min_line_search_step_size")),
83  _lsrch_tol(getParam<Real>("line_search_tol")),
84  _lsrch_max_iter(getParam<unsigned int>("line_search_maxiter")),
85  _lsrch_method(getParam<MooseEnum>("line_search_method")),
86  _fp(declareProperty<RankTwoTensor>("fp")), // Plastic deformation gradient
87  _fp_old(getMaterialPropertyOld<RankTwoTensor>(
88  "fp")), // Plastic deformation gradient of previous increment
89  _pk2(declareProperty<RankTwoTensor>("pk2")), // 2nd Piola Kirchoff Stress
90  _pk2_old(getMaterialPropertyOld<RankTwoTensor>(
91  "pk2")), // 2nd Piola Kirchoff Stress of previous increment
92  _lag_e(declareProperty<RankTwoTensor>("lage")), // Lagrangian strain
93  _update_rot(declareProperty<RankTwoTensor>(
94  "update_rot")), // Rotation tensor considering material rotation and crystal orientation
95  _update_rot_old(getMaterialPropertyOld<RankTwoTensor>("update_rot")),
96  _deformation_gradient(getMaterialProperty<RankTwoTensor>("deformation_gradient")),
97  _deformation_gradient_old(getMaterialPropertyOld<RankTwoTensor>("deformation_gradient")),
98  _crysrot(getMaterialProperty<RankTwoTensor>("crysrot"))
99 {
100  _err_tol = false;
101 
102  _delta_dfgrd.zero();
103 
104  // resize the material properties for each userobject
110 
111  // resize the flow direction
113 
114  // resize local state variables
118 
119  // resize user objects
124 
125  // assign the user objects
126  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
127  {
128  _uo_slip_rates[i] = &getUserObjectByName<CrystalPlasticitySlipRate>(
129  parameters.get<std::vector<UserObjectName>>("uo_slip_rates")[i]);
130  _mat_prop_slip_rates[i] = &declareProperty<std::vector<Real>>(
131  parameters.get<std::vector<UserObjectName>>("uo_slip_rates")[i]);
132  _flow_direction[i] = &declareProperty<std::vector<RankTwoTensor>>(
133  parameters.get<std::vector<UserObjectName>>("uo_slip_rates")[i] + "_flow_direction");
134  }
135 
136  for (unsigned int i = 0; i < _num_uo_slip_resistances; ++i)
137  {
138  _uo_slip_resistances[i] = &getUserObjectByName<CrystalPlasticitySlipResistance>(
139  parameters.get<std::vector<UserObjectName>>("uo_slip_resistances")[i]);
140  _mat_prop_slip_resistances[i] = &declareProperty<std::vector<Real>>(
141  parameters.get<std::vector<UserObjectName>>("uo_slip_resistances")[i]);
142  }
143 
144  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
145  {
146  _uo_state_vars[i] = &getUserObjectByName<CrystalPlasticityStateVariable>(
147  parameters.get<std::vector<UserObjectName>>("uo_state_vars")[i]);
148  _mat_prop_state_vars[i] = &declareProperty<std::vector<Real>>(
149  parameters.get<std::vector<UserObjectName>>("uo_state_vars")[i]);
150  _mat_prop_state_vars_old[i] = &getMaterialPropertyOld<std::vector<Real>>(
151  parameters.get<std::vector<UserObjectName>>("uo_state_vars")[i]);
152  }
153 
154  for (unsigned int i = 0; i < _num_uo_state_var_evol_rate_comps; ++i)
155  {
156  _uo_state_var_evol_rate_comps[i] = &getUserObjectByName<CrystalPlasticityStateVarRateComponent>(
157  parameters.get<std::vector<UserObjectName>>("uo_state_var_evol_rate_comps")[i]);
158  _mat_prop_state_var_evol_rate_comps[i] = &declareProperty<std::vector<Real>>(
159  parameters.get<std::vector<UserObjectName>>("uo_state_var_evol_rate_comps")[i]);
160  }
161 }
const MaterialProperty< RankTwoTensor > & _pk2_old
std::vector< std::vector< Real > > _state_vars_old
Local state variable.
bool _err_tol
Flag to check whether convergence is achieved.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_slip_rates
Slip rates material property.
Real _stol
Internal variable update equation tolerance.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_var_evol_rate_comps
State variable evolution rate component material property.
std::vector< const CrystalPlasticitySlipResistance * > _uo_slip_resistances
User objects that define the slip resistance.
MooseEnum _tan_mod_type
Type of tangent moduli calculation.
unsigned int _lsrch_max_iter
Line search bisection method maximum iteration number.
Real _zero_tol
Residual tolerance when variable value is zero. Default 1e-12.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_slip_resistances
Slip resistance material property.
std::vector< const CrystalPlasticityStateVarRateComponent * > _uo_state_var_evol_rate_comps
User objects that define the state variable evolution rate component.
std::vector< std::vector< Real > > _state_vars_old_stored
Local stored state variable (for sub-stepping)
RankTwoTensor _delta_dfgrd
Used for substepping; Uniformly divides the increment in deformation gradient.
std::vector< std::vector< Real > > _state_vars_prev
Local old state variable.
unsigned int _maxiter
Maximum number of iterations for stress update.
std::vector< const CrystalPlasticityStateVariable * > _uo_state_vars
User objects that define the state variable.
const MaterialProperty< RankTwoTensor > & _update_rot_old
MaterialProperty< RankTwoTensor > & _fp
Real _lsrch_tol
Line search bisection method tolerance.
unsigned int _num_uo_slip_rates
Number of slip rate user objects.
unsigned int _maxiterg
Maximum number of iterations for internal variable update.
const MaterialProperty< RankTwoTensor > & _deformation_gradient
std::vector< const CrystalPlasticitySlipRate * > _uo_slip_rates
User objects that define the slip rate.
unsigned int _num_uo_state_vars
Number of state variable user objects.
unsigned int _num_uo_slip_resistances
Number of slip resistance user objects.
MaterialProperty< RankTwoTensor > & _lag_e
unsigned int _num_uo_state_var_evol_rate_comps
Number of state variable evolution rate component user objects.
Real _rtol
Stress residual equation relative tolerance.
const MaterialProperty< RankTwoTensor > & _crysrot
Crystal rotation.
ComputeStressBase(const InputParameters &parameters)
std::vector< const MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars_old
Old state variable material property.
MaterialProperty< RankTwoTensor > & _update_rot
bool _use_line_search
Flag to activate line serach.
const MaterialProperty< RankTwoTensor > & _fp_old
unsigned int _max_substep_iter
Maximum number of substep iterations.
Real _abs_tol
Stress residual equation absolute tolerance.
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars
State variable material property.
std::vector< MaterialProperty< std::vector< RankTwoTensor > > * > _flow_direction
Real _min_lsrch_step
Minimum line search step size.
MooseEnum _lsrch_method
Line search method.
MaterialProperty< RankTwoTensor > & _pk2

Member Function Documentation

◆ calcJacobian()

void FiniteStrainUObasedCP::calcJacobian ( )
protectedvirtual

calculate jacobian.

Definition at line 528 of file FiniteStrainUObasedCP.C.

Referenced by calcResidJacob().

529 {
530  RankFourTensor dfedfpinv, deedfe, dfpinvdpk2;
531 
532  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
533  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
534  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
535  dfedfpinv(i, j, k, j) = _dfgrd_tmp(i, k);
536 
537  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
538  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
539  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
540  {
541  deedfe(i, j, k, i) = deedfe(i, j, k, i) + _fe(k, j) * 0.5;
542  deedfe(i, j, k, j) = deedfe(i, j, k, j) + _fe(k, i) * 0.5;
543  }
544 
545  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
546  {
547  unsigned int nss = _uo_slip_rates[i]->variableSize();
548  std::vector<RankTwoTensor> dtaudpk2(nss), dfpinvdslip(nss);
549  std::vector<Real> dslipdtau;
550  dslipdtau.resize(nss);
551  _uo_slip_rates[i]->calcSlipRateDerivative(_qp, _dt, dslipdtau);
552  for (unsigned int j = 0; j < nss; j++)
553  {
554  dtaudpk2[j] = (*_flow_direction[i])[_qp][j];
555  dfpinvdslip[j] = -_fp_old_inv * (*_flow_direction[i])[_qp][j];
556  }
557 
558  for (unsigned int j = 0; j < nss; j++)
559  dfpinvdpk2 += (dfpinvdslip[j] * dslipdtau[j] * _dt).outerProduct(dtaudpk2[j]);
560  }
561  _jac =
562  RankFourTensor::IdentityFour() - (_elasticity_tensor[_qp] * deedfe * dfedfpinv * dfpinvdpk2);
563 }
unsigned int _num_uo_slip_rates
Number of slip rate user objects.
std::vector< const CrystalPlasticitySlipRate * > _uo_slip_rates
User objects that define the slip rate.
RankFourTensor _jac
Jacobian tensor.
const MaterialProperty< RankFourTensor > & _elasticity_tensor
std::vector< MaterialProperty< std::vector< RankTwoTensor > > * > _flow_direction

◆ calcResidJacob()

void FiniteStrainUObasedCP::calcResidJacob ( )
protectedvirtual

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

Definition at line 480 of file FiniteStrainUObasedCP.C.

Referenced by solveStress().

481 {
482  calcResidual();
483  if (_err_tol)
484  return;
485  calcJacobian();
486 }
bool _err_tol
Flag to check whether convergence is achieved.
virtual void calcResidual()
calculate stress residual.
virtual void calcJacobian()
calculate jacobian.

◆ calcResidual()

void FiniteStrainUObasedCP::calcResidual ( )
protectedvirtual

calculate stress residual.

Definition at line 502 of file FiniteStrainUObasedCP.C.

Referenced by calcResidJacob(), and lineSearchUpdate().

503 {
504  RankTwoTensor iden(RankTwoTensor::initIdentity), ce, ee, ce_pk2, eqv_slip_incr, pk2_new;
505 
506  getSlipRates();
507  if (_err_tol)
508  return;
509 
510  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
511  for (unsigned int j = 0; j < _uo_slip_rates[i]->variableSize(); ++j)
512  eqv_slip_incr += (*_flow_direction[i])[_qp][j] * (*_mat_prop_slip_rates[i])[_qp][j] * _dt;
513 
514  eqv_slip_incr = iden - eqv_slip_incr;
515  _fp_inv = _fp_old_inv * eqv_slip_incr;
516  _fe = _dfgrd_tmp * _fp_inv;
517 
518  ce = _fe.transpose() * _fe;
519  ee = ce - iden;
520  ee *= 0.5;
521 
522  pk2_new = _elasticity_tensor[_qp] * ee;
523 
524  _resid = _pk2[_qp] - pk2_new;
525 }
bool _err_tol
Flag to check whether convergence is achieved.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_slip_rates
Slip rates material property.
virtual void getSlipRates()
updates the slip rates.
unsigned int _num_uo_slip_rates
Number of slip rate user objects.
std::vector< const CrystalPlasticitySlipRate * > _uo_slip_rates
User objects that define the slip rate.
RankTwoTensor _resid
Residual tensor.
const MaterialProperty< RankFourTensor > & _elasticity_tensor
std::vector< MaterialProperty< std::vector< RankTwoTensor > > * > _flow_direction
MaterialProperty< RankTwoTensor > & _pk2

◆ calcTangentModuli()

void FiniteStrainUObasedCP::calcTangentModuli ( )
protectedvirtual

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 566 of file FiniteStrainUObasedCP.C.

Referenced by postSolveQp().

567 {
568  switch (_tan_mod_type)
569  {
570  case 0:
572  break;
573  default:
575  }
576 }
MooseEnum _tan_mod_type
Type of tangent moduli calculation.
virtual void elastoPlasticTangentModuli()
calculate the exact tangent moduli for preconditioner.
virtual void elasticTangentModuli()
calculate the elastic tangent moduli for preconditioner.

◆ 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 FiniteStrainUObasedCP::computeQpStress ( )
protectedvirtual

updates the stress at a quadrature point.

Solves stress residual equation using NR.

Updates slip system resistances iteratively.

Implements ComputeStressBase.

Definition at line 204 of file FiniteStrainUObasedCP.C.

205 {
206  // Userobject based crystal plasticity does not support face/boundary material property
207  // calculation.
208  if (isBoundaryMaterial())
209  return;
210  // Depth of substepping; Limited to maximum substep iteration
211  unsigned int substep_iter = 1;
212  // Calculated from substep_iter as 2^substep_iter
213  unsigned int num_substep = 1;
214  // Store original _dt; Reset at the end of solve
215  Real dt_original = _dt;
216 
218  if (_dfgrd_tmp_old.det() == 0)
219  _dfgrd_tmp_old.addIa(1.0);
220 
222 
223  // Saves the old stateful properties that is modified during sub stepping
224  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
226 
227  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
228  _uo_slip_rates[i]->calcFlowDirection(_qp, (*_flow_direction[i])[_qp]);
229 
230  do
231  {
232  _err_tol = false;
233 
234  preSolveQp();
235 
236  _dt = dt_original / num_substep;
237 
238  for (unsigned int istep = 0; istep < num_substep; ++istep)
239  {
240  _dfgrd_tmp = (static_cast<Real>(istep) + 1) / num_substep * _delta_dfgrd + _dfgrd_tmp_old;
241 
242  solveQp();
243 
244  if (_err_tol)
245  {
246  substep_iter++;
247  num_substep *= 2;
248  break;
249  }
250  }
251  if (substep_iter > _max_substep_iter && _err_tol)
252  throw MooseException("FiniteStrainUObasedCP: Constitutive failure.");
253  } while (_err_tol);
254 
255  _dt = dt_original;
256 
257  postSolveQp();
258 }
std::vector< std::vector< Real > > _state_vars_old
Local state variable.
bool _err_tol
Flag to check whether convergence is achieved.
RankTwoTensor _delta_dfgrd
Used for substepping; Uniformly divides the increment in deformation gradient.
virtual void postSolveQp()
update stress and internal variable after solve.
unsigned int _num_uo_slip_rates
Number of slip rate user objects.
const MaterialProperty< RankTwoTensor > & _deformation_gradient
std::vector< const CrystalPlasticitySlipRate * > _uo_slip_rates
User objects that define the slip rate.
unsigned int _num_uo_state_vars
Number of state variable user objects.
virtual void preSolveQp()
set variables for stress and internal variable solve.
std::vector< const MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars_old
Old state variable material property.
virtual void solveQp()
solve stress and internal variables.
unsigned int _max_substep_iter
Maximum number of substep iterations.
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
std::vector< MaterialProperty< std::vector< RankTwoTensor > > * > _flow_direction

◆ elasticTangentModuli()

void FiniteStrainUObasedCP::elasticTangentModuli ( )
protectedvirtual

calculate the elastic tangent moduli for preconditioner.

Definition at line 622 of file FiniteStrainUObasedCP.C.

Referenced by calcTangentModuli().

623 {
624  // update jacobian_mult
626 }
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
const MaterialProperty< RankFourTensor > & _elasticity_tensor

◆ elastoPlasticTangentModuli()

void FiniteStrainUObasedCP::elastoPlasticTangentModuli ( )
protectedvirtual

calculate the exact tangent moduli for preconditioner.

Definition at line 579 of file FiniteStrainUObasedCP.C.

Referenced by calcTangentModuli().

580 {
581  RankFourTensor tan_mod;
582  RankTwoTensor pk2fet, fepk2;
583  RankFourTensor deedfe, dsigdpk2dfe, dfedf;
584 
585  // Fill in the matrix stiffness material property
586  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
587  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
588  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
589  {
590  deedfe(i, j, k, i) = deedfe(i, j, k, i) + _fe(k, j) * 0.5;
591  deedfe(i, j, k, j) = deedfe(i, j, k, j) + _fe(k, i) * 0.5;
592  }
593 
594  dsigdpk2dfe = _fe.mixedProductIkJl(_fe) * _elasticity_tensor[_qp] * deedfe;
595 
596  pk2fet = _pk2[_qp] * _fe.transpose();
597  fepk2 = _fe * _pk2[_qp];
598 
599  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
600  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
601  for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
602  {
603  tan_mod(i, j, i, l) += pk2fet(l, j);
604  tan_mod(i, j, j, l) += fepk2(i, l);
605  }
606 
607  tan_mod += dsigdpk2dfe;
608 
609  Real je = _fe.det();
610  if (je > 0.0)
611  tan_mod /= je;
612 
613  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
614  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
615  for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
616  dfedf(i, j, i, l) = _fp_inv(l, j);
617 
618  _Jacobian_mult[_qp] = tan_mod * dfedf;
619 }
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
const MaterialProperty< RankFourTensor > & _elasticity_tensor
MaterialProperty< RankTwoTensor > & _pk2

◆ getSlipRates()

void FiniteStrainUObasedCP::getSlipRates ( )
protectedvirtual

updates the slip rates.

Definition at line 489 of file FiniteStrainUObasedCP.C.

Referenced by calcResidual().

490 {
491  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
492  {
493  if (!_uo_slip_rates[i]->calcSlipRate(_qp, _dt, (*_mat_prop_slip_rates[i])[_qp]))
494  {
495  _err_tol = true;
496  return;
497  }
498  }
499 }
bool _err_tol
Flag to check whether convergence is achieved.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_slip_rates
Slip rates material property.
unsigned int _num_uo_slip_rates
Number of slip rate user objects.
std::vector< const CrystalPlasticitySlipRate * > _uo_slip_rates
User objects that define the slip rate.

◆ initQpStatefulProperties()

void FiniteStrainUObasedCP::initQpStatefulProperties ( )
protectedvirtual

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

Reimplemented from ComputeStressBase.

Definition at line 164 of file FiniteStrainUObasedCP.C.

165 {
166  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
167  {
168  (*_mat_prop_slip_rates[i])[_qp].resize(_uo_slip_rates[i]->variableSize());
169  (*_flow_direction[i])[_qp].resize(_uo_slip_rates[i]->variableSize());
170  }
171 
172  for (unsigned int i = 0; i < _num_uo_slip_resistances; ++i)
173  (*_mat_prop_slip_resistances[i])[_qp].resize(_uo_slip_resistances[i]->variableSize());
174 
175  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
176  {
177  (*_mat_prop_state_vars[i])[_qp].resize(_uo_state_vars[i]->variableSize());
178  _state_vars_old[i].resize(_uo_state_vars[i]->variableSize());
179  _state_vars_old_stored[i].resize(_uo_state_vars[i]->variableSize());
180  _state_vars_prev[i].resize(_uo_state_vars[i]->variableSize());
181  }
182 
183  for (unsigned int i = 0; i < _num_uo_state_var_evol_rate_comps; ++i)
184  (*_mat_prop_state_var_evol_rate_comps[i])[_qp].resize(
185  _uo_state_var_evol_rate_comps[i]->variableSize());
186 
187  _stress[_qp].zero();
188  _pk2[_qp].zero();
189  _lag_e[_qp].zero();
190 
191  _fp[_qp].setToIdentity();
192  _update_rot[_qp].setToIdentity();
193 
194  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
195  // Initializes slip system related properties
196  _uo_state_vars[i]->initSlipSysProps((*_mat_prop_state_vars[i])[_qp], _q_point[_qp]);
197 }
std::vector< std::vector< Real > > _state_vars_old
Local state variable.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_slip_rates
Slip rates material property.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_var_evol_rate_comps
State variable evolution rate component material property.
std::vector< const CrystalPlasticitySlipResistance * > _uo_slip_resistances
User objects that define the slip resistance.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_slip_resistances
Slip resistance material property.
std::vector< const CrystalPlasticityStateVarRateComponent * > _uo_state_var_evol_rate_comps
User objects that define the state variable evolution rate component.
std::vector< std::vector< Real > > _state_vars_old_stored
Local stored state variable (for sub-stepping)
MaterialProperty< RankTwoTensor > & _stress
std::vector< std::vector< Real > > _state_vars_prev
Local old state variable.
std::vector< const CrystalPlasticityStateVariable * > _uo_state_vars
User objects that define the state variable.
MaterialProperty< RankTwoTensor > & _fp
unsigned int _num_uo_slip_rates
Number of slip rate user objects.
std::vector< const CrystalPlasticitySlipRate * > _uo_slip_rates
User objects that define the slip rate.
unsigned int _num_uo_state_vars
Number of state variable user objects.
unsigned int _num_uo_slip_resistances
Number of slip resistance user objects.
MaterialProperty< RankTwoTensor > & _lag_e
unsigned int _num_uo_state_var_evol_rate_comps
Number of state variable evolution rate component user objects.
MaterialProperty< RankTwoTensor > & _update_rot
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars
State variable material property.
std::vector< MaterialProperty< std::vector< RankTwoTensor > > * > _flow_direction
MaterialProperty< RankTwoTensor > & _pk2

◆ isStateVariablesConverged()

bool FiniteStrainUObasedCP::isStateVariablesConverged ( )
protectedvirtual

evaluates convergence of state variables.

Definition at line 347 of file FiniteStrainUObasedCP.C.

Referenced by solveStatevar().

348 {
349  Real diff;
350 
351  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
352  {
353  unsigned int n = (*_mat_prop_state_vars[i])[_qp].size();
354  for (unsigned j = 0; j < n; j++)
355  {
356  diff = std::abs((*_mat_prop_state_vars[i])[_qp][j] -
357  _state_vars_prev[i][j]); // Calculate increment size
358  if (std::abs(_state_vars_old_stored[i][j]) < _zero_tol && diff > _zero_tol)
359  return true;
360  if (std::abs(_state_vars_old_stored[i][j]) > _zero_tol &&
361  diff > _stol * std::abs(_state_vars_old_stored[i][j]))
362  return true;
363  }
364  }
365  return false;
366 }
Real _stol
Internal variable update equation tolerance.
Real _zero_tol
Residual tolerance when variable value is zero. Default 1e-12.
std::vector< std::vector< Real > > _state_vars_old_stored
Local stored state variable (for sub-stepping)
std::vector< std::vector< Real > > _state_vars_prev
Local old state variable.
unsigned int _num_uo_state_vars
Number of state variable user objects.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars
State variable material property.

◆ lineSearchUpdate()

bool FiniteStrainUObasedCP::lineSearchUpdate ( const Real  rnorm_prev,
const RankTwoTensor  dpk2 
)
protected

performs the line search update

Definition at line 629 of file FiniteStrainUObasedCP.C.

Referenced by solveStress().

630 {
631  switch (_lsrch_method)
632  {
633  case 0: // CUT_HALF
634  {
635  Real rnorm;
636  Real step = 1.0;
637 
638  do
639  {
640  _pk2[_qp] = _pk2[_qp] - step * dpk2;
641  step /= 2.0;
642  _pk2[_qp] = _pk2[_qp] + step * dpk2;
643 
644  calcResidual();
645  rnorm = _resid.L2norm();
646  } while (rnorm > rnorm_prev && step > _min_lsrch_step);
647 
648  // has norm improved or is the step still above minumum search step size?
649  return (rnorm <= rnorm_prev || step > _min_lsrch_step);
650  }
651 
652  case 1: // BISECTION
653  {
654  unsigned int count = 0;
655  Real step_a = 0.0;
656  Real step_b = 1.0;
657  Real step = 1.0;
658  Real s_m = 1000.0;
659  Real rnorm = 1000.0;
660 
661  calcResidual();
662  Real s_b = _resid.doubleContraction(dpk2);
663  Real rnorm1 = _resid.L2norm();
664  _pk2[_qp] = _pk2[_qp] - dpk2;
665  calcResidual();
666  Real s_a = _resid.doubleContraction(dpk2);
667  Real rnorm0 = _resid.L2norm();
668  _pk2[_qp] = _pk2[_qp] + dpk2;
669 
670  if ((rnorm1 / rnorm0) < _lsrch_tol || s_a * s_b > 0)
671  {
672  calcResidual();
673  return true;
674  }
675 
676  while ((rnorm / rnorm0) > _lsrch_tol && count < _lsrch_max_iter)
677  {
678  _pk2[_qp] = _pk2[_qp] - step * dpk2;
679  step = 0.5 * (step_b + step_a);
680  _pk2[_qp] = _pk2[_qp] + step * dpk2;
681  calcResidual();
682  s_m = _resid.doubleContraction(dpk2);
683  rnorm = _resid.L2norm();
684 
685  if (s_m * s_a < 0.0)
686  {
687  step_b = step;
688  s_b = s_m;
689  }
690  if (s_m * s_b < 0.0)
691  {
692  step_a = step;
693  s_a = s_m;
694  }
695  count++;
696  }
697 
698  // below tolerance and max iterations?
699  return ((rnorm / rnorm0) < _lsrch_tol && count < _lsrch_max_iter);
700  }
701 
702  default:
703  mooseError("Line search method is not provided.");
704  }
705 }
unsigned int _lsrch_max_iter
Line search bisection method maximum iteration number.
virtual void calcResidual()
calculate stress residual.
Real _lsrch_tol
Line search bisection method tolerance.
RankTwoTensor _resid
Residual tensor.
Real _min_lsrch_step
Minimum line search step size.
MooseEnum _lsrch_method
Line search method.
MaterialProperty< RankTwoTensor > & _pk2

◆ postSolveQp()

void FiniteStrainUObasedCP::postSolveQp ( )
protectedvirtual

update stress and internal variable after solve.

Definition at line 281 of file FiniteStrainUObasedCP.C.

Referenced by computeQpStress().

282 {
283  _stress[_qp] = _fe * _pk2[_qp] * _fe.transpose() / _fe.det();
284 
285  // Calculate jacobian for preconditioner
287 
288  RankTwoTensor iden(RankTwoTensor::initIdentity);
289 
290  _lag_e[_qp] = _deformation_gradient[_qp].transpose() * _deformation_gradient[_qp] - iden;
291  _lag_e[_qp] = _lag_e[_qp] * 0.5;
292 
293  RankTwoTensor rot;
294  // Calculate material rotation
295  _deformation_gradient[_qp].getRUDecompositionRotation(rot);
296  _update_rot[_qp] = rot * _crysrot[_qp];
297 }
MaterialProperty< RankTwoTensor > & _stress
virtual void calcTangentModuli()
calculate the tangent moduli for preconditioner.
const MaterialProperty< RankTwoTensor > & _deformation_gradient
MaterialProperty< RankTwoTensor > & _lag_e
const MaterialProperty< RankTwoTensor > & _crysrot
Crystal rotation.
MaterialProperty< RankTwoTensor > & _update_rot
MaterialProperty< RankTwoTensor > & _pk2

◆ postSolveStatevar()

void FiniteStrainUObasedCP::postSolveStatevar ( )
protectedvirtual

update internal variable after solve.

Definition at line 369 of file FiniteStrainUObasedCP.C.

Referenced by solveQp().

370 {
371  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
373 
375 }
std::vector< std::vector< Real > > _state_vars_old_stored
Local stored state variable (for sub-stepping)
unsigned int _num_uo_state_vars
Number of state variable user objects.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars
State variable material property.

◆ postSolveStress()

void FiniteStrainUObasedCP::postSolveStress ( )
protectedvirtual

update stress and plastic deformation gradient after solve.

Definition at line 452 of file FiniteStrainUObasedCP.C.

Referenced by solveStatevar().

453 {
454  _fp[_qp] = _fp_inv.inverse();
455 }
MaterialProperty< RankTwoTensor > & _fp

◆ preSolveQp()

void FiniteStrainUObasedCP::preSolveQp ( )
protectedvirtual

set variables for stress and internal variable solve.

Definition at line 261 of file FiniteStrainUObasedCP.C.

Referenced by computeQpStress().

262 {
263  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
265 
266  _pk2[_qp] = _pk2_old[_qp];
267  _fp_old_inv = _fp_old[_qp].inverse();
268 }
const MaterialProperty< RankTwoTensor > & _pk2_old
std::vector< std::vector< Real > > _state_vars_old
Local state variable.
std::vector< std::vector< Real > > _state_vars_old_stored
Local stored state variable (for sub-stepping)
unsigned int _num_uo_state_vars
Number of state variable user objects.
const MaterialProperty< RankTwoTensor > & _fp_old
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars
State variable material property.
MaterialProperty< RankTwoTensor > & _pk2

◆ preSolveStatevar()

void FiniteStrainUObasedCP::preSolveStatevar ( )
protectedvirtual

set variables for internal variable solve.

Definition at line 300 of file FiniteStrainUObasedCP.C.

Referenced by solveQp().

301 {
302  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
304 
305  for (unsigned int i = 0; i < _num_uo_slip_resistances; ++i)
306  _uo_slip_resistances[i]->calcSlipResistance(_qp, (*_mat_prop_slip_resistances[i])[_qp]);
307 
309 }
std::vector< const CrystalPlasticitySlipResistance * > _uo_slip_resistances
User objects that define the slip resistance.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_slip_resistances
Slip resistance material property.
std::vector< std::vector< Real > > _state_vars_old_stored
Local stored state variable (for sub-stepping)
unsigned int _num_uo_state_vars
Number of state variable user objects.
unsigned int _num_uo_slip_resistances
Number of slip resistance user objects.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars
State variable material property.

◆ preSolveStress()

void FiniteStrainUObasedCP::preSolveStress ( )
protectedvirtual

set variables for stress solve.

Definition at line 378 of file FiniteStrainUObasedCP.C.

Referenced by solveStatevar().

379 {
380 }

◆ solveQp()

void FiniteStrainUObasedCP::solveQp ( )
protectedvirtual

solve stress and internal variables.

Definition at line 271 of file FiniteStrainUObasedCP.C.

Referenced by computeQpStress().

272 {
274  solveStatevar();
275  if (_err_tol)
276  return;
278 }
virtual void postSolveStatevar()
update internal variable after solve.
bool _err_tol
Flag to check whether convergence is achieved.
virtual void preSolveStatevar()
set variables for internal variable solve.
virtual void solveStatevar()
solve internal variables.

◆ solveStatevar()

void FiniteStrainUObasedCP::solveStatevar ( )
protectedvirtual

solve internal variables.

Definition at line 312 of file FiniteStrainUObasedCP.C.

Referenced by solveQp().

313 {
314  unsigned int iterg;
315  bool iter_flag = true;
316 
317  iterg = 0;
318  // Check for slip system resistance update tolerance
319  while (iter_flag && iterg < _maxiterg)
320  {
321  preSolveStress();
322  solveStress();
323  if (_err_tol)
324  return;
325  postSolveStress();
326 
327  // Update slip system resistance and state variable
329 
330  if (_err_tol)
331  return;
332 
333  iter_flag = isStateVariablesConverged();
334  iterg++;
335  }
336 
337  if (iterg == _maxiterg)
338  {
339 #ifdef DEBUG
340  mooseWarning("FiniteStrainUObasedCP: Hardness Integration error\n");
341 #endif
342  _err_tol = true;
343  }
344 }
virtual bool isStateVariablesConverged()
evaluates convergence of state variables.
bool _err_tol
Flag to check whether convergence is achieved.
virtual void updateSlipSystemResistanceAndStateVariable()
updates the slip system resistances and state variables.
virtual void postSolveStress()
update stress and plastic deformation gradient after solve.
virtual void preSolveStress()
set variables for stress solve.
unsigned int _maxiterg
Maximum number of iterations for internal variable update.
virtual void solveStress()
solves for stress, updates plastic deformation gradient.

◆ solveStress()

void FiniteStrainUObasedCP::solveStress ( )
protectedvirtual

solves for stress, updates plastic deformation gradient.

Definition at line 383 of file FiniteStrainUObasedCP.C.

Referenced by solveStatevar().

384 {
385  unsigned int iter = 0;
386  RankTwoTensor dpk2;
387  Real rnorm, rnorm0, rnorm_prev;
388 
389  // Calculate stress residual
390  calcResidJacob();
391  if (_err_tol)
392  {
393 #ifdef DEBUG
394  mooseWarning("FiniteStrainUObasedCP: Slip increment exceeds tolerance - Element number ",
395  _current_elem->id(),
396  " Gauss point = ",
397  _qp);
398 #endif
399  return;
400  }
401 
402  rnorm = _resid.L2norm();
403  rnorm0 = rnorm;
404 
405  // Check for stress residual tolerance
406  while (rnorm > _rtol * rnorm0 && rnorm0 > _abs_tol && iter < _maxiter)
407  {
408  // Calculate stress increment
409  dpk2 = -_jac.invSymm() * _resid;
410  _pk2[_qp] = _pk2[_qp] + dpk2;
411  calcResidJacob();
412 
413  if (_err_tol)
414  {
415 #ifdef DEBUG
416  mooseWarning("FiniteStrainUObasedCP: Slip increment exceeds tolerance - Element number ",
417  _current_elem->id(),
418  " Gauss point = ",
419  _qp);
420 #endif
421  return;
422  }
423 
424  rnorm_prev = rnorm;
425  rnorm = _resid.L2norm();
426 
427  if (_use_line_search && rnorm > rnorm_prev && !lineSearchUpdate(rnorm_prev, dpk2))
428  {
429 #ifdef DEBUG
430  mooseWarning("FiniteStrainUObasedCP: Failed with line search");
431 #endif
432  _err_tol = true;
433  return;
434  }
435 
436  if (_use_line_search)
437  rnorm = _resid.L2norm();
438 
439  iter++;
440  }
441 
442  if (iter >= _maxiter)
443  {
444 #ifdef DEBUG
445  mooseWarning("FiniteStrainUObasedCP: Stress Integration error rmax = ", rnorm);
446 #endif
447  _err_tol = true;
448  }
449 }
bool _err_tol
Flag to check whether convergence is achieved.
unsigned int _maxiter
Maximum number of iterations for stress update.
virtual void calcResidJacob()
calls the residual and jacobian functions used in the stress update algorithm.
RankFourTensor _jac
Jacobian tensor.
Real _rtol
Stress residual equation relative tolerance.
bool lineSearchUpdate(const Real rnorm_prev, const RankTwoTensor)
performs the line search update
RankTwoTensor _resid
Residual tensor.
bool _use_line_search
Flag to activate line serach.
Real _abs_tol
Stress residual equation absolute tolerance.
MaterialProperty< RankTwoTensor > & _pk2

◆ updateSlipSystemResistanceAndStateVariable()

void FiniteStrainUObasedCP::updateSlipSystemResistanceAndStateVariable ( )
protectedvirtual

updates the slip system resistances and state variables.

override to modify slip system resistance and state variable evolution.

Definition at line 458 of file FiniteStrainUObasedCP.C.

Referenced by solveStatevar().

459 {
460  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
461  _state_vars_prev[i] = (*_mat_prop_state_vars[i])[_qp];
462 
463  for (unsigned int i = 0; i < _num_uo_state_var_evol_rate_comps; ++i)
464  _uo_state_var_evol_rate_comps[i]->calcStateVariableEvolutionRateComponent(
465  _qp, (*_mat_prop_state_var_evol_rate_comps[i])[_qp]);
466 
467  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
468  {
469  if (!_uo_state_vars[i]->updateStateVariable(
470  _qp, _dt, (*_mat_prop_state_vars[i])[_qp], _state_vars_old_stored[i]))
471  _err_tol = true;
472  }
473 
474  for (unsigned int i = 0; i < _num_uo_slip_resistances; ++i)
475  _uo_slip_resistances[i]->calcSlipResistance(_qp, (*_mat_prop_slip_resistances[i])[_qp]);
476 }
bool _err_tol
Flag to check whether convergence is achieved.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_var_evol_rate_comps
State variable evolution rate component material property.
std::vector< const CrystalPlasticitySlipResistance * > _uo_slip_resistances
User objects that define the slip resistance.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_slip_resistances
Slip resistance material property.
std::vector< const CrystalPlasticityStateVarRateComponent * > _uo_state_var_evol_rate_comps
User objects that define the state variable evolution rate component.
std::vector< std::vector< Real > > _state_vars_old_stored
Local stored state variable (for sub-stepping)
std::vector< std::vector< Real > > _state_vars_prev
Local old state variable.
std::vector< const CrystalPlasticityStateVariable * > _uo_state_vars
User objects that define the state variable.
unsigned int _num_uo_state_vars
Number of state variable user objects.
unsigned int _num_uo_slip_resistances
Number of slip resistance user objects.
unsigned int _num_uo_state_var_evol_rate_comps
Number of state variable evolution rate component user objects.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars
State variable material property.

Member Data Documentation

◆ _abs_tol

Real FiniteStrainUObasedCP::_abs_tol
protected

Stress residual equation absolute tolerance.

Definition at line 207 of file FiniteStrainUObasedCP.h.

Referenced by solveStress().

◆ _base_name

const std::string ComputeStressBase::_base_name
protectedinherited

◆ _crysrot

const MaterialProperty<RankTwoTensor>& FiniteStrainUObasedCP::_crysrot
protected

Crystal rotation.

Definition at line 257 of file FiniteStrainUObasedCP.h.

Referenced by postSolveQp().

◆ _deformation_gradient

const MaterialProperty<RankTwoTensor>& FiniteStrainUObasedCP::_deformation_gradient
protected

Definition at line 253 of file FiniteStrainUObasedCP.h.

Referenced by computeQpStress(), and postSolveQp().

◆ _deformation_gradient_old

const MaterialProperty<RankTwoTensor>& FiniteStrainUObasedCP::_deformation_gradient_old
protected

Definition at line 254 of file FiniteStrainUObasedCP.h.

Referenced by computeQpStress().

◆ _delta_dfgrd

RankTwoTensor FiniteStrainUObasedCP::_delta_dfgrd
protected

Used for substepping; Uniformly divides the increment in deformation gradient.

Definition at line 268 of file FiniteStrainUObasedCP.h.

Referenced by computeQpStress(), and FiniteStrainUObasedCP().

◆ _dfgrd_scale_factor

Real FiniteStrainUObasedCP::_dfgrd_scale_factor
protected

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

Definition at line 270 of file FiniteStrainUObasedCP.h.

◆ _dfgrd_tmp

RankTwoTensor FiniteStrainUObasedCP::_dfgrd_tmp
protected

Definition at line 259 of file FiniteStrainUObasedCP.h.

Referenced by calcJacobian(), calcResidual(), and computeQpStress().

◆ _dfgrd_tmp_old

RankTwoTensor FiniteStrainUObasedCP::_dfgrd_tmp_old
protected

Definition at line 268 of file FiniteStrainUObasedCP.h.

Referenced by computeQpStress().

◆ _elastic_strain

MaterialProperty<RankTwoTensor>& ComputeStressBase::_elastic_strain
protectedinherited

◆ _elasticity_tensor

const MaterialProperty<RankFourTensor>& ComputeStressBase::_elasticity_tensor
protectedinherited

◆ _elasticity_tensor_name

const std::string ComputeStressBase::_elasticity_tensor_name
protectedinherited

◆ _err_tol

bool FiniteStrainUObasedCP::_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 FiniteStrainUObasedCP::_fe
protected

◆ _flow_direction

std::vector<MaterialProperty<std::vector<RankTwoTensor> > *> FiniteStrainUObasedCP::_flow_direction
protected

◆ _fp

MaterialProperty<RankTwoTensor>& FiniteStrainUObasedCP::_fp
protected

Definition at line 245 of file FiniteStrainUObasedCP.h.

Referenced by initQpStatefulProperties(), and postSolveStress().

◆ _fp_inv

RankTwoTensor FiniteStrainUObasedCP::_fp_inv
protected

◆ _fp_old

const MaterialProperty<RankTwoTensor>& FiniteStrainUObasedCP::_fp_old
protected

Definition at line 246 of file FiniteStrainUObasedCP.h.

Referenced by preSolveQp().

◆ _fp_old_inv

RankTwoTensor FiniteStrainUObasedCP::_fp_old_inv
protected

◆ _initial_stress_fcn

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

initial stress components

Definition at line 50 of file ComputeStressBase.h.

◆ _jac

RankFourTensor FiniteStrainUObasedCP::_jac
protected

Jacobian tensor.

Definition at line 217 of file FiniteStrainUObasedCP.h.

Referenced by calcJacobian(), and solveStress().

◆ _Jacobian_mult

MaterialProperty<RankFourTensor>& ComputeStressBase::_Jacobian_mult
protectedinherited

◆ _lag_e

MaterialProperty<RankTwoTensor>& FiniteStrainUObasedCP::_lag_e
protected

Definition at line 249 of file FiniteStrainUObasedCP.h.

Referenced by initQpStatefulProperties(), and postSolveQp().

◆ _lsrch_max_iter

unsigned int FiniteStrainUObasedCP::_lsrch_max_iter
protected

Line search bisection method maximum iteration number.

Definition at line 240 of file FiniteStrainUObasedCP.h.

Referenced by lineSearchUpdate().

◆ _lsrch_method

MooseEnum FiniteStrainUObasedCP::_lsrch_method
protected

Line search method.

Definition at line 243 of file FiniteStrainUObasedCP.h.

Referenced by lineSearchUpdate().

◆ _lsrch_tol

Real FiniteStrainUObasedCP::_lsrch_tol
protected

Line search bisection method tolerance.

Definition at line 237 of file FiniteStrainUObasedCP.h.

Referenced by lineSearchUpdate().

◆ _mat_prop_slip_rates

std::vector<MaterialProperty<std::vector<Real> > *> FiniteStrainUObasedCP::_mat_prop_slip_rates
protected

Slip rates material property.

Definition at line 169 of file FiniteStrainUObasedCP.h.

Referenced by calcResidual(), FiniteStrainUObasedCP(), getSlipRates(), and initQpStatefulProperties().

◆ _mat_prop_slip_resistances

std::vector<MaterialProperty<std::vector<Real> > *> FiniteStrainUObasedCP::_mat_prop_slip_resistances
protected

Slip resistance material property.

Definition at line 172 of file FiniteStrainUObasedCP.h.

Referenced by FiniteStrainUObasedCP(), initQpStatefulProperties(), preSolveStatevar(), and updateSlipSystemResistanceAndStateVariable().

◆ _mat_prop_state_var_evol_rate_comps

std::vector<MaterialProperty<std::vector<Real> > *> FiniteStrainUObasedCP::_mat_prop_state_var_evol_rate_comps
protected

State variable evolution rate component material property.

Definition at line 181 of file FiniteStrainUObasedCP.h.

Referenced by FiniteStrainUObasedCP(), initQpStatefulProperties(), and updateSlipSystemResistanceAndStateVariable().

◆ _mat_prop_state_vars

std::vector<MaterialProperty<std::vector<Real> > *> FiniteStrainUObasedCP::_mat_prop_state_vars
protected

◆ _mat_prop_state_vars_old

std::vector<const MaterialProperty<std::vector<Real> > *> FiniteStrainUObasedCP::_mat_prop_state_vars_old
protected

Old state variable material property.

Definition at line 178 of file FiniteStrainUObasedCP.h.

Referenced by computeQpStress(), and FiniteStrainUObasedCP().

◆ _max_substep_iter

unsigned int FiniteStrainUObasedCP::_max_substep_iter
protected

Maximum number of substep iterations.

Definition at line 228 of file FiniteStrainUObasedCP.h.

Referenced by computeQpStress().

◆ _maxiter

unsigned int FiniteStrainUObasedCP::_maxiter
protected

Maximum number of iterations for stress update.

Definition at line 220 of file FiniteStrainUObasedCP.h.

Referenced by solveStress().

◆ _maxiterg

unsigned int FiniteStrainUObasedCP::_maxiterg
protected

Maximum number of iterations for internal variable update.

Definition at line 222 of file FiniteStrainUObasedCP.h.

Referenced by solveStatevar().

◆ _mechanical_strain

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_mechanical_strain
protectedinherited

◆ _min_lsrch_step

Real FiniteStrainUObasedCP::_min_lsrch_step
protected

Minimum line search step size.

Definition at line 234 of file FiniteStrainUObasedCP.h.

Referenced by lineSearchUpdate().

◆ _num_uo_slip_rates

unsigned int FiniteStrainUObasedCP::_num_uo_slip_rates
protected

Number of slip rate user objects.

Definition at line 184 of file FiniteStrainUObasedCP.h.

Referenced by calcJacobian(), calcResidual(), computeQpStress(), FiniteStrainUObasedCP(), getSlipRates(), and initQpStatefulProperties().

◆ _num_uo_slip_resistances

unsigned int FiniteStrainUObasedCP::_num_uo_slip_resistances
protected

Number of slip resistance user objects.

Definition at line 187 of file FiniteStrainUObasedCP.h.

Referenced by FiniteStrainUObasedCP(), initQpStatefulProperties(), preSolveStatevar(), and updateSlipSystemResistanceAndStateVariable().

◆ _num_uo_state_var_evol_rate_comps

unsigned int FiniteStrainUObasedCP::_num_uo_state_var_evol_rate_comps
protected

Number of state variable evolution rate component user objects.

Definition at line 193 of file FiniteStrainUObasedCP.h.

Referenced by FiniteStrainUObasedCP(), initQpStatefulProperties(), and updateSlipSystemResistanceAndStateVariable().

◆ _num_uo_state_vars

unsigned int FiniteStrainUObasedCP::_num_uo_state_vars
protected

◆ _pk2

MaterialProperty<RankTwoTensor>& FiniteStrainUObasedCP::_pk2
protected

◆ _pk2_old

const MaterialProperty<RankTwoTensor>& FiniteStrainUObasedCP::_pk2_old
protected

Definition at line 248 of file FiniteStrainUObasedCP.h.

Referenced by preSolveQp().

◆ _resid

RankTwoTensor FiniteStrainUObasedCP::_resid
protected

Residual tensor.

Definition at line 214 of file FiniteStrainUObasedCP.h.

Referenced by calcResidual(), lineSearchUpdate(), and solveStress().

◆ _rtol

Real FiniteStrainUObasedCP::_rtol
protected

Stress residual equation relative tolerance.

Definition at line 205 of file FiniteStrainUObasedCP.h.

Referenced by solveStress().

◆ _state_vars_old

std::vector<std::vector<Real> > FiniteStrainUObasedCP::_state_vars_old
protected

Local state variable.

Definition at line 196 of file FiniteStrainUObasedCP.h.

Referenced by computeQpStress(), FiniteStrainUObasedCP(), initQpStatefulProperties(), and preSolveQp().

◆ _state_vars_old_stored

std::vector<std::vector<Real> > FiniteStrainUObasedCP::_state_vars_old_stored
protected

◆ _state_vars_prev

std::vector<std::vector<Real> > FiniteStrainUObasedCP::_state_vars_prev
protected

◆ _stol

Real FiniteStrainUObasedCP::_stol
protected

Internal variable update equation tolerance.

Definition at line 209 of file FiniteStrainUObasedCP.h.

Referenced by isStateVariablesConverged().

◆ _stress

MaterialProperty<RankTwoTensor>& ComputeStressBase::_stress
protectedinherited

◆ _tan_mod_type

MooseEnum FiniteStrainUObasedCP::_tan_mod_type
protected

Type of tangent moduli calculation.

Definition at line 225 of file FiniteStrainUObasedCP.h.

Referenced by calcTangentModuli().

◆ _tau

DenseVector<Real> FiniteStrainUObasedCP::_tau
protected

Definition at line 261 of file FiniteStrainUObasedCP.h.

◆ _uo_slip_rates

std::vector<const CrystalPlasticitySlipRate *> FiniteStrainUObasedCP::_uo_slip_rates
protected

User objects that define the slip rate.

Definition at line 157 of file FiniteStrainUObasedCP.h.

Referenced by calcJacobian(), calcResidual(), computeQpStress(), FiniteStrainUObasedCP(), getSlipRates(), and initQpStatefulProperties().

◆ _uo_slip_resistances

std::vector<const CrystalPlasticitySlipResistance *> FiniteStrainUObasedCP::_uo_slip_resistances
protected

User objects that define the slip resistance.

Definition at line 160 of file FiniteStrainUObasedCP.h.

Referenced by FiniteStrainUObasedCP(), initQpStatefulProperties(), preSolveStatevar(), and updateSlipSystemResistanceAndStateVariable().

◆ _uo_state_var_evol_rate_comps

std::vector<const CrystalPlasticityStateVarRateComponent *> FiniteStrainUObasedCP::_uo_state_var_evol_rate_comps
protected

User objects that define the state variable evolution rate component.

Definition at line 166 of file FiniteStrainUObasedCP.h.

Referenced by FiniteStrainUObasedCP(), initQpStatefulProperties(), and updateSlipSystemResistanceAndStateVariable().

◆ _uo_state_vars

std::vector<const CrystalPlasticityStateVariable *> FiniteStrainUObasedCP::_uo_state_vars
protected

User objects that define the state variable.

Definition at line 163 of file FiniteStrainUObasedCP.h.

Referenced by FiniteStrainUObasedCP(), initQpStatefulProperties(), and updateSlipSystemResistanceAndStateVariable().

◆ _update_rot

MaterialProperty<RankTwoTensor>& FiniteStrainUObasedCP::_update_rot
protected

Definition at line 250 of file FiniteStrainUObasedCP.h.

Referenced by initQpStatefulProperties(), and postSolveQp().

◆ _update_rot_old

const MaterialProperty<RankTwoTensor>& FiniteStrainUObasedCP::_update_rot_old
protected

Definition at line 251 of file FiniteStrainUObasedCP.h.

◆ _use_line_search

bool FiniteStrainUObasedCP::_use_line_search
protected

Flag to activate line serach.

Definition at line 231 of file FiniteStrainUObasedCP.h.

Referenced by solveStress().

◆ _zero_tol

Real FiniteStrainUObasedCP::_zero_tol
protected

Residual tolerance when variable value is zero. Default 1e-12.

Definition at line 211 of file FiniteStrainUObasedCP.h.

Referenced by isStateVariablesConverged().


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