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 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 > & _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
 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

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  _elasticity_tensor_name(_base_name + "elasticity_tensor"),
97  _elasticity_tensor(getMaterialPropertyByName<RankFourTensor>(_elasticity_tensor_name)),
98  _deformation_gradient(getMaterialProperty<RankTwoTensor>("deformation_gradient")),
99  _deformation_gradient_old(getMaterialPropertyOld<RankTwoTensor>("deformation_gradient")),
100  _crysrot(getMaterialProperty<RankTwoTensor>("crysrot"))
101 {
102  _err_tol = false;
103 
104  _delta_dfgrd.zero();
105 
106  // resize the material properties for each userobject
112 
113  // resize the flow direction
115 
116  // resize local state variables
120 
121  // resize user objects
126 
127  // assign the user objects
128  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
129  {
130  _uo_slip_rates[i] = &getUserObjectByName<CrystalPlasticitySlipRate>(
131  parameters.get<std::vector<UserObjectName>>("uo_slip_rates")[i]);
132  _mat_prop_slip_rates[i] = &declareProperty<std::vector<Real>>(
133  parameters.get<std::vector<UserObjectName>>("uo_slip_rates")[i]);
134  _flow_direction[i] = &declareProperty<std::vector<RankTwoTensor>>(
135  parameters.get<std::vector<UserObjectName>>("uo_slip_rates")[i] + "_flow_direction");
136  }
137 
138  for (unsigned int i = 0; i < _num_uo_slip_resistances; ++i)
139  {
140  _uo_slip_resistances[i] = &getUserObjectByName<CrystalPlasticitySlipResistance>(
141  parameters.get<std::vector<UserObjectName>>("uo_slip_resistances")[i]);
142  _mat_prop_slip_resistances[i] = &declareProperty<std::vector<Real>>(
143  parameters.get<std::vector<UserObjectName>>("uo_slip_resistances")[i]);
144  }
145 
146  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
147  {
148  _uo_state_vars[i] = &getUserObjectByName<CrystalPlasticityStateVariable>(
149  parameters.get<std::vector<UserObjectName>>("uo_state_vars")[i]);
150  _mat_prop_state_vars[i] = &declareProperty<std::vector<Real>>(
151  parameters.get<std::vector<UserObjectName>>("uo_state_vars")[i]);
152  _mat_prop_state_vars_old[i] = &getMaterialPropertyOld<std::vector<Real>>(
153  parameters.get<std::vector<UserObjectName>>("uo_state_vars")[i]);
154  }
155 
156  for (unsigned int i = 0; i < _num_uo_state_var_evol_rate_comps; ++i)
157  {
158  _uo_state_var_evol_rate_comps[i] = &getUserObjectByName<CrystalPlasticityStateVarRateComponent>(
159  parameters.get<std::vector<UserObjectName>>("uo_state_var_evol_rate_comps")[i]);
160  _mat_prop_state_var_evol_rate_comps[i] = &declareProperty<std::vector<Real>>(
161  parameters.get<std::vector<UserObjectName>>("uo_state_var_evol_rate_comps")[i]);
162  }
163 }
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.
const std::string _base_name
Base name prepended to all material property names to allow for multi-material systems.
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
const std::string _elasticity_tensor_name
Name of the elasticity tensor material property.
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.
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor 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 530 of file FiniteStrainUObasedCP.C.

Referenced by calcResidJacob().

531 {
532  RankFourTensor dfedfpinv, deedfe, dfpinvdpk2;
533 
534  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
535  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
536  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
537  dfedfpinv(i, j, k, j) = _dfgrd_tmp(i, k);
538 
539  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
540  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
541  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
542  {
543  deedfe(i, j, k, i) = deedfe(i, j, k, i) + _fe(k, j) * 0.5;
544  deedfe(i, j, k, j) = deedfe(i, j, k, j) + _fe(k, i) * 0.5;
545  }
546 
547  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
548  {
549  unsigned int nss = _uo_slip_rates[i]->variableSize();
550  std::vector<RankTwoTensor> dtaudpk2(nss), dfpinvdslip(nss);
551  std::vector<Real> dslipdtau;
552  dslipdtau.resize(nss);
553  _uo_slip_rates[i]->calcSlipRateDerivative(_qp, _dt, dslipdtau);
554  for (unsigned int j = 0; j < nss; j++)
555  {
556  dtaudpk2[j] = (*_flow_direction[i])[_qp][j];
557  dfpinvdslip[j] = -_fp_old_inv * (*_flow_direction[i])[_qp][j];
558  }
559 
560  for (unsigned int j = 0; j < nss; j++)
561  dfpinvdpk2 += (dfpinvdslip[j] * dslipdtau[j] * _dt).outerProduct(dtaudpk2[j]);
562  }
563  _jac =
564  RankFourTensor::IdentityFour() - (_elasticity_tensor[_qp] * deedfe * dfedfpinv * dfpinvdpk2);
565 }
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
Elasticity tensor material property.
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 482 of file FiniteStrainUObasedCP.C.

Referenced by solveStress().

483 {
484  calcResidual();
485  if (_err_tol)
486  return;
487  calcJacobian();
488 }
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 504 of file FiniteStrainUObasedCP.C.

Referenced by calcResidJacob(), and lineSearchUpdate().

505 {
506  RankTwoTensor iden(RankTwoTensor::initIdentity), ce, ee, ce_pk2, eqv_slip_incr, pk2_new;
507 
508  getSlipRates();
509  if (_err_tol)
510  return;
511 
512  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
513  for (unsigned int j = 0; j < _uo_slip_rates[i]->variableSize(); ++j)
514  eqv_slip_incr += (*_flow_direction[i])[_qp][j] * (*_mat_prop_slip_rates[i])[_qp][j] * _dt;
515 
516  eqv_slip_incr = iden - eqv_slip_incr;
517  _fp_inv = _fp_old_inv * eqv_slip_incr;
518  _fe = _dfgrd_tmp * _fp_inv;
519 
520  ce = _fe.transpose() * _fe;
521  ee = ce - iden;
522  ee *= 0.5;
523 
524  pk2_new = _elasticity_tensor[_qp] * ee;
525 
526  _resid = _pk2[_qp] - pk2_new;
527 }
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
Elasticity tensor material property.
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 568 of file FiniteStrainUObasedCP.C.

Referenced by postSolveQp().

569 {
570  switch (_tan_mod_type)
571  {
572  case 0:
574  break;
575  default:
577  }
578 }
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 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 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 206 of file FiniteStrainUObasedCP.C.

207 {
208  // Userobject based crystal plasticity does not support face/boundary material property
209  // calculation.
210  if (isBoundaryMaterial())
211  return;
212  // Depth of substepping; Limited to maximum substep iteration
213  unsigned int substep_iter = 1;
214  // Calculated from substep_iter as 2^substep_iter
215  unsigned int num_substep = 1;
216  // Store original _dt; Reset at the end of solve
217  Real dt_original = _dt;
218 
220  if (_dfgrd_tmp_old.det() == 0)
221  _dfgrd_tmp_old.addIa(1.0);
222 
224 
225  // Saves the old stateful properties that is modified during sub stepping
226  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
228 
229  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
230  _uo_slip_rates[i]->calcFlowDirection(_qp, (*_flow_direction[i])[_qp]);
231 
232  do
233  {
234  _err_tol = false;
235 
236  preSolveQp();
237 
238  _dt = dt_original / num_substep;
239 
240  for (unsigned int istep = 0; istep < num_substep; ++istep)
241  {
242  _dfgrd_tmp = (static_cast<Real>(istep) + 1) / num_substep * _delta_dfgrd + _dfgrd_tmp_old;
243 
244  solveQp();
245 
246  if (_err_tol)
247  {
248  substep_iter++;
249  num_substep *= 2;
250  break;
251  }
252  }
253  if (substep_iter > _max_substep_iter && _err_tol)
254  throw MooseException("FiniteStrainUObasedCP: Constitutive failure.");
255  } while (_err_tol);
256 
257  _dt = dt_original;
258 
259  postSolveQp();
260 }
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 624 of file FiniteStrainUObasedCP.C.

Referenced by calcTangentModuli().

625 {
626  // update jacobian_mult
628 }
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor material property.

◆ elastoPlasticTangentModuli()

void FiniteStrainUObasedCP::elastoPlasticTangentModuli ( )
protectedvirtual

calculate the exact tangent moduli for preconditioner.

Definition at line 581 of file FiniteStrainUObasedCP.C.

Referenced by calcTangentModuli().

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

◆ getSlipRates()

void FiniteStrainUObasedCP::getSlipRates ( )
protectedvirtual

updates the slip rates.

Definition at line 491 of file FiniteStrainUObasedCP.C.

Referenced by calcResidual().

492 {
493  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
494  {
495  if (!_uo_slip_rates[i]->calcSlipRate(_qp, _dt, (*_mat_prop_slip_rates[i])[_qp]))
496  {
497  _err_tol = true;
498  return;
499  }
500  }
501 }
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 166 of file FiniteStrainUObasedCP.C.

167 {
168  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
169  {
170  (*_mat_prop_slip_rates[i])[_qp].resize(_uo_slip_rates[i]->variableSize());
171  (*_flow_direction[i])[_qp].resize(_uo_slip_rates[i]->variableSize());
172  }
173 
174  for (unsigned int i = 0; i < _num_uo_slip_resistances; ++i)
175  (*_mat_prop_slip_resistances[i])[_qp].resize(_uo_slip_resistances[i]->variableSize());
176 
177  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
178  {
179  (*_mat_prop_state_vars[i])[_qp].resize(_uo_state_vars[i]->variableSize());
180  _state_vars_old[i].resize(_uo_state_vars[i]->variableSize());
181  _state_vars_old_stored[i].resize(_uo_state_vars[i]->variableSize());
182  _state_vars_prev[i].resize(_uo_state_vars[i]->variableSize());
183  }
184 
185  for (unsigned int i = 0; i < _num_uo_state_var_evol_rate_comps; ++i)
186  (*_mat_prop_state_var_evol_rate_comps[i])[_qp].resize(
187  _uo_state_var_evol_rate_comps[i]->variableSize());
188 
189  _stress[_qp].zero();
190  _pk2[_qp].zero();
191  _lag_e[_qp].zero();
192 
193  _fp[_qp].setToIdentity();
194  _update_rot[_qp].setToIdentity();
195 
196  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
197  // Initializes slip system related properties
198  _uo_state_vars[i]->initSlipSysProps((*_mat_prop_state_vars[i])[_qp], _q_point[_qp]);
199 }
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
Stress material property.
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 349 of file FiniteStrainUObasedCP.C.

Referenced by solveStatevar().

350 {
351  Real diff;
352 
353  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
354  {
355  unsigned int n = (*_mat_prop_state_vars[i])[_qp].size();
356  for (unsigned j = 0; j < n; j++)
357  {
358  diff = std::abs((*_mat_prop_state_vars[i])[_qp][j] -
359  _state_vars_prev[i][j]); // Calculate increment size
360  if (std::abs(_state_vars_old_stored[i][j]) < _zero_tol && diff > _zero_tol)
361  return true;
362  if (std::abs(_state_vars_old_stored[i][j]) > _zero_tol &&
363  diff > _stol * std::abs(_state_vars_old_stored[i][j]))
364  return true;
365  }
366  }
367  return false;
368 }
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 631 of file FiniteStrainUObasedCP.C.

Referenced by solveStress().

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

Referenced by computeQpStress().

284 {
285  _stress[_qp] = _fe * _pk2[_qp] * _fe.transpose() / _fe.det();
286 
287  // Calculate jacobian for preconditioner
289 
290  RankTwoTensor iden(RankTwoTensor::initIdentity);
291 
292  _lag_e[_qp] = _deformation_gradient[_qp].transpose() * _deformation_gradient[_qp] - iden;
293  _lag_e[_qp] = _lag_e[_qp] * 0.5;
294 
295  RankTwoTensor rot;
296  // Calculate material rotation
297  _deformation_gradient[_qp].getRUDecompositionRotation(rot);
298  _update_rot[_qp] = rot * _crysrot[_qp];
299 }
MaterialProperty< RankTwoTensor > & _stress
Stress material property.
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 371 of file FiniteStrainUObasedCP.C.

Referenced by solveQp().

372 {
373  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
375 
377 }
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 454 of file FiniteStrainUObasedCP.C.

Referenced by solveStatevar().

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

◆ preSolveQp()

void FiniteStrainUObasedCP::preSolveQp ( )
protectedvirtual

set variables for stress and internal variable solve.

Definition at line 263 of file FiniteStrainUObasedCP.C.

Referenced by computeQpStress().

264 {
265  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
267 
268  _pk2[_qp] = _pk2_old[_qp];
269  _fp_old_inv = _fp_old[_qp].inverse();
270 }
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 302 of file FiniteStrainUObasedCP.C.

Referenced by solveQp().

303 {
304  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
306 
307  for (unsigned int i = 0; i < _num_uo_slip_resistances; ++i)
308  _uo_slip_resistances[i]->calcSlipResistance(_qp, (*_mat_prop_slip_resistances[i])[_qp]);
309 
311 }
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 380 of file FiniteStrainUObasedCP.C.

Referenced by solveStatevar().

381 {
382 }

◆ solveQp()

void FiniteStrainUObasedCP::solveQp ( )
protectedvirtual

solve stress and internal variables.

Definition at line 273 of file FiniteStrainUObasedCP.C.

Referenced by computeQpStress().

274 {
276  solveStatevar();
277  if (_err_tol)
278  return;
280 }
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 314 of file FiniteStrainUObasedCP.C.

Referenced by solveQp().

315 {
316  unsigned int iterg;
317  bool iter_flag = true;
318 
319  iterg = 0;
320  // Check for slip system resistance update tolerance
321  while (iter_flag && iterg < _maxiterg)
322  {
323  preSolveStress();
324  solveStress();
325  if (_err_tol)
326  return;
327  postSolveStress();
328 
329  // Update slip system resistance and state variable
331 
332  if (_err_tol)
333  return;
334 
335  iter_flag = isStateVariablesConverged();
336  iterg++;
337  }
338 
339  if (iterg == _maxiterg)
340  {
341 #ifdef DEBUG
342  mooseWarning("FiniteStrainUObasedCP: Hardness Integration error\n");
343 #endif
344  _err_tol = true;
345  }
346 }
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 385 of file FiniteStrainUObasedCP.C.

Referenced by solveStatevar().

386 {
387  unsigned int iter = 0;
388  RankTwoTensor dpk2;
389  Real rnorm, rnorm0, rnorm_prev;
390 
391  // Calculate stress residual
392  calcResidJacob();
393  if (_err_tol)
394  {
395 #ifdef DEBUG
396  mooseWarning("FiniteStrainUObasedCP: Slip increment exceeds tolerance - Element number ",
397  _current_elem->id(),
398  " Gauss point = ",
399  _qp);
400 #endif
401  return;
402  }
403 
404  rnorm = _resid.L2norm();
405  rnorm0 = rnorm;
406 
407  // Check for stress residual tolerance
408  while (rnorm > _rtol * rnorm0 && rnorm0 > _abs_tol && iter < _maxiter)
409  {
410  // Calculate stress increment
411  dpk2 = -_jac.invSymm() * _resid;
412  _pk2[_qp] = _pk2[_qp] + dpk2;
413  calcResidJacob();
414 
415  if (_err_tol)
416  {
417 #ifdef DEBUG
418  mooseWarning("FiniteStrainUObasedCP: Slip increment exceeds tolerance - Element number ",
419  _current_elem->id(),
420  " Gauss point = ",
421  _qp);
422 #endif
423  return;
424  }
425 
426  rnorm_prev = rnorm;
427  rnorm = _resid.L2norm();
428 
429  if (_use_line_search && rnorm > rnorm_prev && !lineSearchUpdate(rnorm_prev, dpk2))
430  {
431 #ifdef DEBUG
432  mooseWarning("FiniteStrainUObasedCP: Failed with line search");
433 #endif
434  _err_tol = true;
435  return;
436  }
437 
438  if (_use_line_search)
439  rnorm = _resid.L2norm();
440 
441  iter++;
442  }
443 
444  if (iter >= _maxiter)
445  {
446 #ifdef DEBUG
447  mooseWarning("FiniteStrainUObasedCP: Stress Integration error rmax = ", rnorm);
448 #endif
449  _err_tol = true;
450  }
451 }
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 460 of file FiniteStrainUObasedCP.C.

Referenced by solveStatevar().

461 {
462  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
463  _state_vars_prev[i] = (*_mat_prop_state_vars[i])[_qp];
464 
465  for (unsigned int i = 0; i < _num_uo_state_var_evol_rate_comps; ++i)
466  _uo_state_var_evol_rate_comps[i]->calcStateVariableEvolutionRateComponent(
467  _qp, (*_mat_prop_state_var_evol_rate_comps[i])[_qp]);
468 
469  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
470  {
471  if (!_uo_state_vars[i]->updateStateVariable(
472  _qp, _dt, (*_mat_prop_state_vars[i])[_qp], _state_vars_old_stored[i]))
473  _err_tol = true;
474  }
475 
476  for (unsigned int i = 0; i < _num_uo_slip_resistances; ++i)
477  _uo_slip_resistances[i]->calcSlipResistance(_qp, (*_mat_prop_slip_resistances[i])[_qp]);
478 }
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

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>& FiniteStrainUObasedCP::_crysrot
protected

Crystal rotation.

Definition at line 261 of file FiniteStrainUObasedCP.h.

Referenced by postSolveQp().

◆ _deformation_gradient

const MaterialProperty<RankTwoTensor>& FiniteStrainUObasedCP::_deformation_gradient
protected

Definition at line 257 of file FiniteStrainUObasedCP.h.

Referenced by computeQpStress(), and postSolveQp().

◆ _deformation_gradient_old

const MaterialProperty<RankTwoTensor>& FiniteStrainUObasedCP::_deformation_gradient_old
protected

Definition at line 258 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 272 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 274 of file FiniteStrainUObasedCP.h.

◆ _dfgrd_tmp

RankTwoTensor FiniteStrainUObasedCP::_dfgrd_tmp
protected

Definition at line 263 of file FiniteStrainUObasedCP.h.

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

◆ _dfgrd_tmp_old

RankTwoTensor FiniteStrainUObasedCP::_dfgrd_tmp_old
protected

Definition at line 272 of file FiniteStrainUObasedCP.h.

Referenced by computeQpStress().

◆ _elastic_strain

MaterialProperty<RankTwoTensor>& ComputeStressBase::_elastic_strain
protectedinherited

◆ _elasticity_tensor

const MaterialProperty<RankFourTensor>& FiniteStrainUObasedCP::_elasticity_tensor
protected

Elasticity tensor material property.

Definition at line 256 of file FiniteStrainUObasedCP.h.

Referenced by calcJacobian(), calcResidual(), elasticTangentModuli(), and elastoPlasticTangentModuli().

◆ _elasticity_tensor_name

const std::string FiniteStrainUObasedCP::_elasticity_tensor_name
protected

Name of the elasticity tensor material property.

Definition at line 254 of file FiniteStrainUObasedCP.h.

◆ _err_tol

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

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(), FiniteStrainCrystalPlasticity::initQpStatefulProperties(), initQpStatefulProperties(), FiniteStrainHyperElasticViscoPlastic::initQpStatefulProperties(), ComputeMultiPlasticityStress::postReturnMap(), postSolveQp(), FiniteStrainHyperElasticViscoPlastic::postSolveQp(), FiniteStrainCrystalPlasticity::postSolveQp(), ComputeSmearedCrackingStress::updateCrackingStateAndStress(), ComputeMultipleInelasticStress::updateQpState(), and ComputeMultipleInelasticStress::updateQpStateSingleModel().

◆ _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 265 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: