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

67  : ComputeStressBase(parameters),
68  _num_uo_slip_rates(parameters.get<std::vector<UserObjectName>>("uo_slip_rates").size()),
70  parameters.get<std::vector<UserObjectName>>("uo_slip_resistances").size()),
71  _num_uo_state_vars(parameters.get<std::vector<UserObjectName>>("uo_state_vars").size()),
73  parameters.get<std::vector<UserObjectName>>("uo_state_var_evol_rate_comps").size()),
74  _rtol(getParam<Real>("rtol")),
75  _abs_tol(getParam<Real>("abs_tol")),
76  _stol(getParam<Real>("stol")),
77  _zero_tol(getParam<Real>("zero_tol")),
78  _maxiter(getParam<unsigned int>("maxiter")),
79  _maxiterg(getParam<unsigned int>("maxiter_state_variable")),
80  _tan_mod_type(getParam<MooseEnum>("tan_mod_type")),
81  _max_substep_iter(getParam<unsigned int>("maximum_substep_iteration")),
82  _use_line_search(getParam<bool>("use_line_search")),
83  _min_lsrch_step(getParam<Real>("min_line_search_step_size")),
84  _lsrch_tol(getParam<Real>("line_search_tol")),
85  _lsrch_max_iter(getParam<unsigned int>("line_search_maxiter")),
86  _lsrch_method(getParam<MooseEnum>("line_search_method")),
87  _fp(declareProperty<RankTwoTensor>("fp")), // Plastic deformation gradient
88  _fp_old(getMaterialPropertyOld<RankTwoTensor>(
89  "fp")), // Plastic deformation gradient of previous increment
90  _pk2(declareProperty<RankTwoTensor>("pk2")), // 2nd Piola Kirchoff Stress
91  _pk2_old(getMaterialPropertyOld<RankTwoTensor>(
92  "pk2")), // 2nd Piola Kirchoff Stress of previous increment
93  _lag_e(declareProperty<RankTwoTensor>("lage")), // Lagrangian strain
94  _update_rot(declareProperty<RankTwoTensor>(
95  "update_rot")), // Rotation tensor considering material rotation and crystal orientation
96  _update_rot_old(getMaterialPropertyOld<RankTwoTensor>("update_rot")),
97  _deformation_gradient(getMaterialProperty<RankTwoTensor>("deformation_gradient")),
98  _deformation_gradient_old(getMaterialPropertyOld<RankTwoTensor>("deformation_gradient")),
99  _crysrot(getMaterialProperty<RankTwoTensor>("crysrot"))
100 {
101  _err_tol = false;
102 
103  _delta_dfgrd.zero();
104 
105  // resize the material properties for each userobject
111 
112  // resize the flow direction
114 
115  // 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.
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 557 of file FiniteStrainUObasedCP.C.

Referenced by calcResidJacob().

558 {
559  RankFourTensor dfedfpinv, deedfe, dfpinvdpk2;
560 
561  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
562  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
563  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
564  dfedfpinv(i, j, k, j) = _dfgrd_tmp(i, k);
565 
566  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
567  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
568  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
569  {
570  deedfe(i, j, k, i) = deedfe(i, j, k, i) + _fe(k, j) * 0.5;
571  deedfe(i, j, k, j) = deedfe(i, j, k, j) + _fe(k, i) * 0.5;
572  }
573 
574  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
575  {
576  unsigned int nss = _uo_slip_rates[i]->variableSize();
577  std::vector<RankTwoTensor> dtaudpk2(nss), dfpinvdslip(nss);
578  std::vector<Real> dslipdtau;
579  dslipdtau.resize(nss);
580  _uo_slip_rates[i]->calcSlipRateDerivative(_qp, _dt, dslipdtau);
581  for (unsigned int j = 0; j < nss; j++)
582  {
583  dtaudpk2[j] = (*_flow_direction[i])[_qp][j];
584  dfpinvdslip[j] = -_fp_old_inv * (*_flow_direction[i])[_qp][j];
585  }
586 
587  for (unsigned int j = 0; j < nss; j++)
588  dfpinvdpk2 += (dfpinvdslip[j] * dslipdtau[j] * _dt).outerProduct(dtaudpk2[j]);
589  }
590  _jac =
591  RankFourTensor::IdentityFour() - (_elasticity_tensor[_qp] * deedfe * dfedfpinv * dfpinvdpk2);
592 }
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 503 of file FiniteStrainUObasedCP.C.

Referenced by solveStress().

504 {
505  calcResidual();
506  if (_err_tol)
507  return;
508  calcJacobian();
509 }
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 525 of file FiniteStrainUObasedCP.C.

Referenced by calcResidJacob(), and lineSearchUpdate().

526 {
527  RankTwoTensor iden, ce, ee, ce_pk2, eqv_slip_incr, pk2_new;
528 
529  iden.zero();
530  iden.addIa(1.0);
531 
532  eqv_slip_incr.zero();
533 
534  getSlipRates();
535 
536  if (_err_tol)
537  return;
538 
539  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
540  for (unsigned int j = 0; j < _uo_slip_rates[i]->variableSize(); ++j)
541  eqv_slip_incr += (*_flow_direction[i])[_qp][j] * (*_mat_prop_slip_rates[i])[_qp][j] * _dt;
542 
543  eqv_slip_incr = iden - eqv_slip_incr;
544  _fp_inv = _fp_old_inv * eqv_slip_incr;
545  _fe = _dfgrd_tmp * _fp_inv;
546 
547  ce = _fe.transpose() * _fe;
548  ee = ce - iden;
549  ee *= 0.5;
550 
551  pk2_new = _elasticity_tensor[_qp] * ee;
552 
553  _resid = _pk2[_qp] - pk2_new;
554 }
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 595 of file FiniteStrainUObasedCP.C.

Referenced by postSolveQp().

596 {
597  switch (_tan_mod_type)
598  {
599  case 0:
601  break;
602  default:
604  }
605 }
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 216 of file FiniteStrainUObasedCP.C.

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

Referenced by calcTangentModuli().

652 {
653  // update jacobian_mult
655 }
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 608 of file FiniteStrainUObasedCP.C.

Referenced by calcTangentModuli().

609 {
610  RankFourTensor tan_mod;
611  RankTwoTensor pk2fet, fepk2;
612  RankFourTensor deedfe, dsigdpk2dfe, dfedf;
613 
614  // Fill in the matrix stiffness material property
615  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
616  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
617  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
618  {
619  deedfe(i, j, k, i) = deedfe(i, j, k, i) + _fe(k, j) * 0.5;
620  deedfe(i, j, k, j) = deedfe(i, j, k, j) + _fe(k, i) * 0.5;
621  }
622 
623  dsigdpk2dfe = _fe.mixedProductIkJl(_fe) * _elasticity_tensor[_qp] * deedfe;
624 
625  pk2fet = _pk2[_qp] * _fe.transpose();
626  fepk2 = _fe * _pk2[_qp];
627 
628  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
629  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
630  for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
631  {
632  tan_mod(i, j, i, l) += pk2fet(l, j);
633  tan_mod(i, j, j, l) += fepk2(i, l);
634  }
635 
636  tan_mod += dsigdpk2dfe;
637 
638  Real je = _fe.det();
639  if (je > 0.0)
640  tan_mod /= je;
641 
642  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
643  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
644  for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
645  dfedf(i, j, i, l) = _fp_inv(l, j);
646 
647  _Jacobian_mult[_qp] = tan_mod * dfedf;
648 }
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 512 of file FiniteStrainUObasedCP.C.

Referenced by calcResidual().

513 {
514  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
515  {
516  if (!_uo_slip_rates[i]->calcSlipRate(_qp, _dt, (*_mat_prop_slip_rates[i])[_qp]))
517  {
518  _err_tol = true;
519  return;
520  }
521  }
522 }
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  // TODO: remove this nasty const_cast if you can figure out how
179  const_cast<MaterialProperty<std::vector<Real>> &>(*_mat_prop_state_vars_old[i])[_qp].resize(
180  _uo_state_vars[i]->variableSize());
181  _state_vars_old[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 
191  _fp[_qp].zero();
192  _fp[_qp].addIa(1.0);
193 
194  _pk2[_qp].zero();
195 
196  _lag_e[_qp].zero();
197 
198  _update_rot[_qp].zero();
199  _update_rot[_qp].addIa(1.0);
200 
201  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
202  {
203  // Initializes slip system related properties
204  _uo_state_vars[i]->initSlipSysProps((*_mat_prop_state_vars[i])[_qp], _q_point[_qp]);
205  // TODO: remove this nasty const_cast if you can figure out how
206  const_cast<MaterialProperty<std::vector<Real>> &>(*_mat_prop_state_vars_old[i])[_qp] =
207  (*_mat_prop_state_vars[i])[_qp];
208  }
209 }
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.
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.
std::vector< const MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars_old
Old state variable material property.
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 369 of file FiniteStrainUObasedCP.C.

Referenced by solveStatevar().

370 {
371  Real diff;
372 
373  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
374  {
375  unsigned int n = (*_mat_prop_state_vars[i])[_qp].size();
376  for (unsigned j = 0; j < n; j++)
377  {
378  diff = std::abs((*_mat_prop_state_vars[i])[_qp][j] -
379  _state_vars_prev[i][j]); // Calculate increment size
380  if (std::abs((*_mat_prop_state_vars_old[i])[_qp][j]) < _zero_tol && diff > _zero_tol)
381  return true;
382  if (std::abs((*_mat_prop_state_vars_old[i])[_qp][j]) > _zero_tol &&
383  diff > _stol * std::abs((*_mat_prop_state_vars_old[i])[_qp][j]))
384  return true;
385  }
386  }
387  return false;
388 }
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_prev
Local old state variable.
unsigned int _num_uo_state_vars
Number of state variable user objects.
std::vector< const MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars_old
Old state variable material property.
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 658 of file FiniteStrainUObasedCP.C.

Referenced by solveStress().

659 {
660  switch (_lsrch_method)
661  {
662  case 0: // CUT_HALF
663  {
664  Real rnorm;
665  Real step = 1.0;
666 
667  do
668  {
669  _pk2[_qp] = _pk2[_qp] - step * dpk2;
670  step /= 2.0;
671  _pk2[_qp] = _pk2[_qp] + step * dpk2;
672 
673  calcResidual();
674  rnorm = _resid.L2norm();
675  } while (rnorm > rnorm_prev && step > _min_lsrch_step);
676 
677  // has norm improved or is the step still above minumum search step size?
678  return (rnorm <= rnorm_prev || step > _min_lsrch_step);
679  }
680 
681  case 1: // BISECTION
682  {
683  unsigned int count = 0;
684  Real step_a = 0.0;
685  Real step_b = 1.0;
686  Real step = 1.0;
687  Real s_m = 1000.0;
688  Real rnorm = 1000.0;
689 
690  calcResidual();
691  Real s_b = _resid.doubleContraction(dpk2);
692  Real rnorm1 = _resid.L2norm();
693  _pk2[_qp] = _pk2[_qp] - dpk2;
694  calcResidual();
695  Real s_a = _resid.doubleContraction(dpk2);
696  Real rnorm0 = _resid.L2norm();
697  _pk2[_qp] = _pk2[_qp] + dpk2;
698 
699  if ((rnorm1 / rnorm0) < _lsrch_tol || s_a * s_b > 0)
700  {
701  calcResidual();
702  return true;
703  }
704 
705  while ((rnorm / rnorm0) > _lsrch_tol && count < _lsrch_max_iter)
706  {
707  _pk2[_qp] = _pk2[_qp] - step * dpk2;
708  step = 0.5 * (step_b + step_a);
709  _pk2[_qp] = _pk2[_qp] + step * dpk2;
710  calcResidual();
711  s_m = _resid.doubleContraction(dpk2);
712  rnorm = _resid.L2norm();
713 
714  if (s_m * s_a < 0.0)
715  {
716  step_b = step;
717  s_b = s_m;
718  }
719  if (s_m * s_b < 0.0)
720  {
721  step_a = step;
722  s_a = s_m;
723  }
724  count++;
725  }
726 
727  // below tolerance and max iterations?
728  return ((rnorm / rnorm0) < _lsrch_tol && count < _lsrch_max_iter);
729  }
730 
731  default:
732  mooseError("Line search method is not provided.");
733  }
734 }
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 296 of file FiniteStrainUObasedCP.C.

Referenced by computeQpStress().

297 {
298  // TODO: remove this nasty const_cast if you can figure out how
299  // Restores the the old stateful properties after a successful solve
300  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
301  const_cast<MaterialProperty<std::vector<Real>> &>(*_mat_prop_state_vars_old[i])[_qp] =
302  _state_vars_old[i];
303 
304  _stress[_qp] = _fe * _pk2[_qp] * _fe.transpose() / _fe.det();
305 
306  // Calculate jacobian for preconditioner
308 
309  RankTwoTensor iden;
310  iden.addIa(1.0);
311 
312  _lag_e[_qp] = _deformation_gradient[_qp].transpose() * _deformation_gradient[_qp] - iden;
313  _lag_e[_qp] = _lag_e[_qp] * 0.5;
314 
315  RankTwoTensor rot;
316  // Calculate material rotation
317  _deformation_gradient[_qp].getRUDecompositionRotation(rot);
318  _update_rot[_qp] = rot * _crysrot[_qp];
319 }
std::vector< std::vector< Real > > _state_vars_old
Local state variable.
MaterialProperty< RankTwoTensor > & _stress
virtual void calcTangentModuli()
calculate the tangent moduli for preconditioner.
const MaterialProperty< RankTwoTensor > & _deformation_gradient
unsigned int _num_uo_state_vars
Number of state variable user objects.
MaterialProperty< RankTwoTensor > & _lag_e
const MaterialProperty< RankTwoTensor > & _crysrot
Crystal rotation.
std::vector< const MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars_old
Old state variable material property.
MaterialProperty< RankTwoTensor > & _update_rot
MaterialProperty< RankTwoTensor > & _pk2

◆ postSolveStatevar()

void FiniteStrainUObasedCP::postSolveStatevar ( )
protectedvirtual

update internal variable after solve.

Definition at line 391 of file FiniteStrainUObasedCP.C.

Referenced by solveQp().

392 {
393  // TODO: remove this nasty const_cast if you can figure out how
394  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
395  const_cast<MaterialProperty<std::vector<Real>> &>(*_mat_prop_state_vars_old[i])[_qp] =
396  (*_mat_prop_state_vars[i])[_qp];
397 
399 }
unsigned int _num_uo_state_vars
Number of state variable user objects.
std::vector< const MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars_old
Old state variable material property.
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 476 of file FiniteStrainUObasedCP.C.

Referenced by solveStatevar().

477 {
478  _fp[_qp] = _fp_inv.inverse();
479 }
MaterialProperty< RankTwoTensor > & _fp

◆ preSolveQp()

void FiniteStrainUObasedCP::preSolveQp ( )
protectedvirtual

set variables for stress and internal variable solve.

Definition at line 273 of file FiniteStrainUObasedCP.C.

Referenced by computeQpStress().

274 {
275  // TODO: remove this nasty const_cast if you can figure out how
276  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
277  (*_mat_prop_state_vars[i])[_qp] =
278  const_cast<MaterialProperty<std::vector<Real>> &>(*_mat_prop_state_vars_old[i])[_qp] =
279  _state_vars_old[i];
280 
281  _pk2[_qp] = _pk2_old[_qp];
282  _fp_old_inv = _fp_old[_qp].inverse();
283 }
const MaterialProperty< RankTwoTensor > & _pk2_old
std::vector< std::vector< Real > > _state_vars_old
Local state variable.
unsigned int _num_uo_state_vars
Number of state variable user objects.
std::vector< const MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars_old
Old state variable material property.
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 322 of file FiniteStrainUObasedCP.C.

Referenced by solveQp().

323 {
324  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
325  (*_mat_prop_state_vars[i])[_qp] = (*_mat_prop_state_vars_old[i])[_qp];
326 
327  for (unsigned int i = 0; i < _num_uo_slip_resistances; ++i)
328  _uo_slip_resistances[i]->calcSlipResistance(_qp, (*_mat_prop_slip_resistances[i])[_qp]);
329 
331 }
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.
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< const MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars_old
Old state variable material property.
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 402 of file FiniteStrainUObasedCP.C.

Referenced by solveStatevar().

403 {
404 }

◆ solveQp()

void FiniteStrainUObasedCP::solveQp ( )
protectedvirtual

solve stress and internal variables.

Definition at line 286 of file FiniteStrainUObasedCP.C.

Referenced by computeQpStress().

287 {
289  solveStatevar();
290  if (_err_tol)
291  return;
293 }
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 334 of file FiniteStrainUObasedCP.C.

Referenced by solveQp().

335 {
336  unsigned int iterg;
337  bool iter_flag = true;
338 
339  iterg = 0;
340  // Check for slip system resistance update tolerance
341  while (iter_flag && iterg < _maxiterg)
342  {
343  preSolveStress();
344  solveStress();
345  if (_err_tol)
346  return;
347  postSolveStress();
348 
349  // Update slip system resistance and state variable
351 
352  if (_err_tol)
353  return;
354 
355  iter_flag = isStateVariablesConverged();
356  iterg++;
357  }
358 
359  if (iterg == _maxiterg)
360  {
361 #ifdef DEBUG
362  mooseWarning("FiniteStrainUObasedCP: Hardness Integration error\n");
363 #endif
364  _err_tol = true;
365  }
366 }
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 407 of file FiniteStrainUObasedCP.C.

Referenced by solveStatevar().

408 {
409  unsigned int iter = 0;
410  RankTwoTensor dpk2;
411  Real rnorm, rnorm0, rnorm_prev;
412 
413  // Calculate stress residual
414  calcResidJacob();
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 = _resid.L2norm();
427  rnorm0 = rnorm;
428 
429  // Check for stress residual tolerance
430  while (rnorm > _rtol * rnorm0 && rnorm0 > _abs_tol && iter < _maxiter)
431  {
432  // Calculate stress increment
433  dpk2 = -_jac.invSymm() * _resid;
434  _pk2[_qp] = _pk2[_qp] + dpk2;
435  calcResidJacob();
436 
437  if (_err_tol)
438  {
439 #ifdef DEBUG
440  mooseWarning("FiniteStrainUObasedCP: Slip increment exceeds tolerance - Element number ",
441  _current_elem->id(),
442  " Gauss point = ",
443  _qp);
444 #endif
445  return;
446  }
447 
448  rnorm_prev = rnorm;
449  rnorm = _resid.L2norm();
450 
451  if (_use_line_search && rnorm > rnorm_prev && !lineSearchUpdate(rnorm_prev, dpk2))
452  {
453 #ifdef DEBUG
454  mooseWarning("FiniteStrainUObasedCP: Failed with line search");
455 #endif
456  _err_tol = true;
457  return;
458  }
459 
460  if (_use_line_search)
461  rnorm = _resid.L2norm();
462 
463  iter++;
464  }
465 
466  if (iter >= _maxiter)
467  {
468 #ifdef DEBUG
469  mooseWarning("FiniteStrainUObasedCP: Stress Integration error rmax = ", rnorm);
470 #endif
471  _err_tol = true;
472  }
473 }
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 482 of file FiniteStrainUObasedCP.C.

Referenced by solveStatevar().

483 {
484  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
485  _state_vars_prev[i] = (*_mat_prop_state_vars[i])[_qp];
486 
487  for (unsigned int i = 0; i < _num_uo_state_var_evol_rate_comps; ++i)
488  _uo_state_var_evol_rate_comps[i]->calcStateVariableEvolutionRateComponent(
489  _qp, (*_mat_prop_state_var_evol_rate_comps[i])[_qp]);
490 
491  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
492  {
493  if (!_uo_state_vars[i]->updateStateVariable(_qp, _dt, (*_mat_prop_state_vars[i])[_qp]))
494  _err_tol = true;
495  }
496 
497  for (unsigned int i = 0; i < _num_uo_slip_resistances; ++i)
498  _uo_slip_resistances[i]->calcSlipResistance(_qp, (*_mat_prop_slip_resistances[i])[_qp]);
499 }
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_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 204 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 254 of file FiniteStrainUObasedCP.h.

Referenced by postSolveQp().

◆ _deformation_gradient

const MaterialProperty<RankTwoTensor>& FiniteStrainUObasedCP::_deformation_gradient
protected

Definition at line 250 of file FiniteStrainUObasedCP.h.

Referenced by computeQpStress(), and postSolveQp().

◆ _deformation_gradient_old

const MaterialProperty<RankTwoTensor>& FiniteStrainUObasedCP::_deformation_gradient_old
protected

Definition at line 251 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 265 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 267 of file FiniteStrainUObasedCP.h.

◆ _dfgrd_tmp

RankTwoTensor FiniteStrainUObasedCP::_dfgrd_tmp
protected

Definition at line 256 of file FiniteStrainUObasedCP.h.

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

◆ _dfgrd_tmp_old

RankTwoTensor FiniteStrainUObasedCP::_dfgrd_tmp_old
protected

Definition at line 265 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 242 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 243 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 214 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 246 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 237 of file FiniteStrainUObasedCP.h.

Referenced by lineSearchUpdate().

◆ _lsrch_method

MooseEnum FiniteStrainUObasedCP::_lsrch_method
protected

Line search method.

Definition at line 240 of file FiniteStrainUObasedCP.h.

Referenced by lineSearchUpdate().

◆ _lsrch_tol

Real FiniteStrainUObasedCP::_lsrch_tol
protected

Line search bisection method tolerance.

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

◆ _max_substep_iter

unsigned int FiniteStrainUObasedCP::_max_substep_iter
protected

Maximum number of substep iterations.

Definition at line 225 of file FiniteStrainUObasedCP.h.

Referenced by computeQpStress().

◆ _maxiter

unsigned int FiniteStrainUObasedCP::_maxiter
protected

Maximum number of iterations for stress update.

Definition at line 217 of file FiniteStrainUObasedCP.h.

Referenced by solveStress().

◆ _maxiterg

unsigned int FiniteStrainUObasedCP::_maxiterg
protected

Maximum number of iterations for internal variable update.

Definition at line 219 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 231 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 245 of file FiniteStrainUObasedCP.h.

Referenced by preSolveQp().

◆ _resid

RankTwoTensor FiniteStrainUObasedCP::_resid
protected

Residual tensor.

Definition at line 211 of file FiniteStrainUObasedCP.h.

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

◆ _rtol

Real FiniteStrainUObasedCP::_rtol
protected

Stress residual equation relative tolerance.

Definition at line 202 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(), postSolveQp(), and preSolveQp().

◆ _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 206 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 222 of file FiniteStrainUObasedCP.h.

Referenced by calcTangentModuli().

◆ _tau

DenseVector<Real> FiniteStrainUObasedCP::_tau
protected

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

Referenced by initQpStatefulProperties(), and postSolveQp().

◆ _update_rot_old

const MaterialProperty<RankTwoTensor>& FiniteStrainUObasedCP::_update_rot_old
protected

Definition at line 248 of file FiniteStrainUObasedCP.h.

◆ _use_line_search

bool FiniteStrainUObasedCP::_use_line_search
protected

Flag to activate line serach.

Definition at line 228 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 208 of file FiniteStrainUObasedCP.h.

Referenced by isStateVariablesConverged().


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