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

#include <FiniteStrainCPSlipRateRes.h>

Inheritance diagram for FiniteStrainCPSlipRateRes:
[legend]

Public Member Functions

 FiniteStrainCPSlipRateRes (const InputParameters &parameters)
 

Protected Member Functions

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

Protected Attributes

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

Detailed Description

Definition at line 18 of file FiniteStrainCPSlipRateRes.h.

Constructor & Destructor Documentation

FiniteStrainCPSlipRateRes::FiniteStrainCPSlipRateRes ( const InputParameters &  parameters)

Definition at line 19 of file FiniteStrainCPSlipRateRes.C.

20  : FiniteStrainCrystalPlasticity(parameters),
21  _resid(_nss),
24  _jacob(_nss, _nss),
26 {
27 }
FiniteStrainCrystalPlasticity(const InputParameters &parameters)
const unsigned int _nss
Number of slip system resistance.
DenseMatrix< Real > _dsliprate_dsliprate

Member Function Documentation

void ComputeStressBase::addQpInitialStress ( )
protectedinherited

InitialStress Deprecation: remove this method.

Adds initial stress, if it is provided, to _stress[_qp]. This function should NOT be used if you calculate stress using

stress = stress_old + elasticity * strain_increment

because stress_old will already include initial stress. However this function SHOULD be used if your code uses

stress = elasticity * (elastic_strain_old + strain_increment) or stress = elasticity * elastic_strain

since in these cases the elastic_strain and elastic_strain_old will not include any contribution from initial stress.

Definition at line 112 of file ComputeStressBase.C.

Referenced by ComputeLinearElasticStress::computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeSmearedCrackingStress::computeQpStress(), ComputeMultipleInelasticStress::computeQpStressIntermediateConfiguration(), ComputeStressBase::initQpStatefulProperties(), ComputeMultipleInelasticStress::updateQpState(), and ComputeMultipleInelasticStress::updateQpStateSingleModel().

113 {
115  _stress[_qp] += (*_initial_stress)[_qp];
116 }
MaterialProperty< RankTwoTensor > & _stress
const bool _initial_stress_provided
Whether initial stress was provided. InitialStress Deprecation: remove this.
void FiniteStrainCrystalPlasticity::assignSlipSysRes ( )
protectedvirtualinherited

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

Definition at line 227 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::initSlipSysProps().

228 {
229  _gss[_qp].resize(_nss);
230 
231  for (unsigned int i = 0; i < _nss; ++i)
232  _gss[_qp][i] = _slip_sys_props(i);
233 }
const unsigned int _nss
Number of slip system resistance.
MaterialProperty< std::vector< Real > > & _gss
void FiniteStrainCrystalPlasticity::calc_resid_jacob ( RankTwoTensor &  resid,
RankFourTensor &  jac 
)
protectedvirtualinherited

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

Definition at line 875 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::solveStress().

876 {
877  calcResidual(resid);
878  if (_err_tol)
879  return;
880  calcJacobian(jac);
881 }
virtual void calcJacobian(RankFourTensor &)
This function calculate jacobian.
virtual void calcResidual(RankTwoTensor &)
This function calculate stress residual.
void FiniteStrainCrystalPlasticity::calc_schmid_tensor ( )
protectedinherited

This function calculate the Schmid tensor.

Definition at line 1046 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::preSolveQp().

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

This function calculates partial derivative of slip system resistances with respect to split rate.

Definition at line 218 of file FiniteStrainCPSlipRateRes.C.

Referenced by calcJacobianSlipRate().

219 {
220  for (unsigned int i = 0; i < _nss; ++i)
221  for (unsigned int j = 0; j < _nss; ++j)
223 }
const unsigned int _nss
Number of slip system resistance.
DenseMatrix< Real > _dsliprate_dsliprate
void FiniteStrainCPSlipRateRes::calcDtauDsliprate ( )
protectedvirtual

This function calculates partial derivative of resolved shear stress with respect to split rate.

Definition at line 183 of file FiniteStrainCPSlipRateRes.C.

Referenced by calcJacobianSlipRate().

184 {
185  RankFourTensor dfedfpinv, deedfe, dfpinvdpk2;
186  std::vector<RankTwoTensor> dtaudpk2(_nss), dfpinvdsliprate(_nss);
187 
188  for (unsigned int i = 0; i < _nss; ++i)
189  {
190  dtaudpk2[i] = _s0[i];
191  dfpinvdsliprate[i] = -_fp_old_inv * _s0[i] * _dt;
192  }
193 
194  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
195  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
196  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
197  dfedfpinv(i, j, k, j) = _dfgrd_tmp(i, k);
198 
199  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
200  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
201  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
202  {
203  deedfe(i, j, k, i) = deedfe(i, j, k, i) + _fe(k, j) * 0.5;
204  deedfe(i, j, k, j) = deedfe(i, j, k, j) + _fe(k, i) * 0.5;
205  }
206 
207  RankFourTensor dpk2dfpinv;
208 
209  dpk2dfpinv = _elasticity_tensor[_qp] * deedfe * dfedfpinv;
210 
211  for (unsigned int i = 0; i < _nss; ++i)
212  for (unsigned int j = 0; j < _nss; ++j)
213  _dsliprate_dsliprate(i, j) =
214  _dslipdtau(i) * dtaudpk2[i].doubleContraction(dpk2dfpinv * dfpinvdsliprate[j]);
215 }
const MaterialProperty< RankFourTensor > & _elasticity_tensor
const unsigned int _nss
Number of slip system resistance.
DenseMatrix< Real > _dsliprate_dsliprate
void FiniteStrainCrystalPlasticity::calcJacobian ( RankFourTensor &  jac)
protectedvirtualinherited

This function calculate jacobian.

Definition at line 924 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::calc_resid_jacob().

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

This function calculates jacobian of slip rate.

Definition at line 166 of file FiniteStrainCPSlipRateRes.C.

Referenced by calcResidJacobSlipRate(), and solveStress().

167 {
168  //_dsliprate_dsliprate not reinitialized to zero, hence order is important
171 
172  for (unsigned int i = 0; i < _nss; ++i)
173  for (unsigned int j = 0; j < _nss; ++j)
174  {
175  _jacob(i, j) = 0.0;
176  if (i == j)
177  _jacob(i, j) += 1.0;
178  _jacob(i, j) -= _dsliprate_dsliprate(i, j);
179  }
180 }
const unsigned int _nss
Number of slip system resistance.
DenseMatrix< Real > _dsliprate_dsliprate
virtual void calcDtauDsliprate()
This function calculates partial derivative of resolved shear stress with respect to split rate...
virtual void calcDgssDsliprate()
This function calculates partial derivative of slip system resistances with respect to split rate...
Real FiniteStrainCPSlipRateRes::calcResidDotProdUpdate ( const DenseVector< Real > &  update)
protected

This function calculates the dot product of residual and update.

Definition at line 364 of file FiniteStrainCPSlipRateRes.C.

Referenced by lineSearchUpdateSlipRate().

365 {
366  Real dotprod = 0.0;
367  for (unsigned int i = 0; i < _nss; ++i)
368  dotprod += _resid(i) * update(i);
369  return dotprod;
370 }
const unsigned int _nss
Number of slip system resistance.
void FiniteStrainCPSlipRateRes::calcResidJacobSlipRate ( )
protectedvirtual

This function calculates residual and jacobian of slip rate.

Definition at line 121 of file FiniteStrainCPSlipRateRes.C.

Referenced by solveStress().

122 {
124  if (_err_tol)
125  return;
127 }
virtual void calcJacobianSlipRate()
This function calculates jacobian of slip rate.
virtual void calcResidualSlipRate()
This function calculates residual of slip rate.
Real FiniteStrainCPSlipRateRes::calcResidNorm ( )
protectedvirtual

This function calculates the residual norm.

Definition at line 256 of file FiniteStrainCPSlipRateRes.C.

Referenced by lineSearchUpdateSlipRate(), and solveStress().

257 {
258  Real rnorm = 0.0;
259  for (unsigned int i = 0; i < _nss; ++i)
260  rnorm += Utility::pow<2>(_resid(i));
261  rnorm = std::sqrt(rnorm) / _nss;
262 
263  return rnorm;
264 }
const unsigned int _nss
Number of slip system resistance.
void FiniteStrainCrystalPlasticity::calcResidual ( RankTwoTensor &  resid)
protectedvirtualinherited

This function calculate stress residual.

Definition at line 884 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::calc_resid_jacob(), and FiniteStrainCrystalPlasticity::line_search_update().

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

This function calculates residual of slip rate.

Definition at line 130 of file FiniteStrainCPSlipRateRes.C.

Referenced by calcResidJacobSlipRate(), lineSearchUpdateSlipRate(), and solveStress().

131 {
132  RankTwoTensor eqv_slip_incr, ce, ee;
133  const RankTwoTensor iden(RankTwoTensor::initIdentity);
134 
136  _slip_incr *= _dt;
137 
138  for (unsigned int i = 0; i < _nss; ++i)
139  eqv_slip_incr += _s0[i] * _slip_incr(i);
140 
141  eqv_slip_incr = iden - eqv_slip_incr;
142 
143  _fp_inv = _fp_old_inv * eqv_slip_incr;
144  _fe = _dfgrd_tmp * _fp_inv;
145 
146  ce = _fe.transpose() * _fe;
147  ee = ce - iden;
148  ee *= 0.5;
149 
150  _pk2_tmp = _elasticity_tensor[_qp] * ee;
151 
152  for (unsigned int i = 0; i < _nss; ++i)
153  _tau(i) = _pk2_tmp.doubleContraction(_s0[i]);
154 
157 
158  if (_err_tol)
159  return;
160 
161  for (unsigned int i = 0; i < _nss; ++i)
162  _resid(i) = _slip_rate(i) - _slip_incr(i) / _dt;
163 }
const MaterialProperty< RankFourTensor > & _elasticity_tensor
virtual void getSlipIncrements()
This function updates the slip system resistances.
virtual void update_slip_system_resistance()
This function updates the slip system resistances.
const unsigned int _nss
Number of slip system resistance.
RankFourTensor FiniteStrainCrystalPlasticity::calcTangentModuli ( )
protectedvirtualinherited

This function calculate the tangent moduli for preconditioner.

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

Definition at line 1029 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::postSolveQp().

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

This function calculates and updates the residual of slip rate.

Definition at line 242 of file FiniteStrainCPSlipRateRes.C.

Referenced by solveStress().

243 {
244  DenseMatrix<Real> A = _jacob;
245  DenseVector<Real> r(_nss);
246  DenseVector<Real> x(_nss);
247 
248  r = _resid;
249 
250  A.lu_solve(r, x);
251 
252  _resid = x;
253 }
const unsigned int _nss
Number of slip system resistance.
void FiniteStrainCrystalPlasticity::computeQpElasticityTensor ( )
protectedvirtualinherited

This function updates the elasticity tensor at a quadrature point.

Presently void.

Definition at line 1024 of file FiniteStrainCrystalPlasticity.C.

1025 {
1026 }
void ComputeStressBase::computeQpProperties ( )
overrideprotectedvirtualinherited

Reimplemented in ComputeBirchMurnaghanEquationOfStress, and ComputeStressEosBase.

Definition at line 99 of file ComputeStressBase.C.

100 {
101  // InitialStress Deprecation: remove the following 2 lines
103  (*_initial_stress)[_qp] = (*_initial_stress_old)[_qp];
104 
105  computeQpStress();
106 
107  // Add in extra stress
108  _stress[_qp] += _extra_stress[_qp];
109 }
virtual void computeQpStress()=0
MaterialProperty< RankTwoTensor > & _stress
const bool _initial_stress_provided
Whether initial stress was provided. InitialStress Deprecation: remove this.
const MaterialProperty< RankTwoTensor > & _extra_stress
Extra stress tensor.
void FiniteStrainCrystalPlasticity::computeQpStress ( )
protectedvirtualinherited

This function updates the stress at a quadrature point.

Solves stress residual equation using NR.

Updates slip system resistances iteratively.

Implements ComputeStressBase.

Definition at line 486 of file FiniteStrainCrystalPlasticity.C.

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

This function calculate the elastic tangent moduli for preconditioner.

Definition at line 1117 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::calcTangentModuli().

1118 {
1119  return _elasticity_tensor[_qp]; // update jacobian_mult
1120 }
const MaterialProperty< RankFourTensor > & _elasticity_tensor
RankFourTensor FiniteStrainCrystalPlasticity::elastoPlasticTangentModuli ( )
protectedvirtualinherited

This function calculate the exact tangent moduli for preconditioner.

Definition at line 1078 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::calcTangentModuli().

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

This function perform RU decomposition to obtain the rotation tensor.

Definition at line 982 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::postSolveQp().

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

This function assign flow rate parameters - see test.

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

Definition at line 337 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::initSlipSysProps().

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

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

Definition at line 405 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::initSlipSysProps().

406 {
407  if (_hprops.size() <= 0)
408  mooseError("FiniteStrainCrystalPLasticity: Error in reading hardness properties: Specify input "
409  "in .i file or a slip_sys_hard_prop_file_name");
410 
411  _r = _hprops[0];
412  _h0 = _hprops[1];
413  _tau_init = _hprops[2];
414  _tau_sat = _hprops[3];
415 }
void FiniteStrainCrystalPlasticity::getInitSlipSysRes ( )
protectedvirtualinherited

This function assign slip system resistances - see test.

.i input file format start_slip_sys_num, end_slip_sys_num, value.

Definition at line 255 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::initSlipSysProps().

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

This function perform RU decomposition to obtain the rotation tensor.

Definition at line 989 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::get_current_rotation().

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

This function updates the slip system resistances.

Reimplemented from FiniteStrainCrystalPlasticity.

Definition at line 226 of file FiniteStrainCPSlipRateRes.C.

Referenced by calcResidualSlipRate().

227 {
229 
230  if (_err_tol)
231  return;
232 
233  _dslipdtau *= 1.0 / _dt;
234 
235  for (unsigned int i = 0; i < _nss; ++i)
236  _dsliprate_dgss(i) = -_a0(i) / _xm(i) *
237  std::pow(std::abs(_tau(i) / _gss_tmp[i]), 1.0 / _xm(i) - 1.0) * _tau(i) /
238  std::pow(_gss_tmp[i], 2.0);
239 }
const unsigned int _nss
Number of slip system resistance.
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
virtual void getSlipIncrements()
This function updates the slip increments.
void FiniteStrainCrystalPlasticity::getSlipSystems ( )
protectedvirtualinherited

This function reads slip system from file - see test.

Definition at line 419 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::FiniteStrainCrystalPlasticity().

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

This function initializes additional parameters.

Definition at line 477 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::initQpStatefulProperties().

478 {
479 }
void FiniteStrainCrystalPlasticity::initQpStatefulProperties ( )
protectedvirtualinherited

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

Reimplemented from ComputeStressBase.

Definition at line 182 of file FiniteStrainCrystalPlasticity.C.

183 {
184  _stress[_qp].zero();
185 
186  _fp[_qp].zero();
187  _fp[_qp].addIa(1.0);
188 
189  _pk2[_qp].zero();
190  _acc_slip[_qp] = 0.0;
191  _lag_e[_qp].zero();
192 
193  _update_rot[_qp].zero();
194  _update_rot[_qp].addIa(1.0);
195 
196  initSlipSysProps(); // Initializes slip system related properties
198 }
MaterialProperty< RankTwoTensor > & _lag_e
MaterialProperty< RankTwoTensor > & _pk2
MaterialProperty< RankTwoTensor > & _stress
MaterialProperty< RankTwoTensor > & _fp
virtual void initSlipSysProps()
This function initializes slip system resistances.
virtual void initAdditionalProps()
This function initializes additional parameters.
MaterialProperty< RankTwoTensor > & _update_rot
void FiniteStrainCrystalPlasticity::initSlipSysProps ( )
protectedvirtualinherited

This function initializes slip system resistances.

Definition at line 201 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::initQpStatefulProperties().

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

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

Definition at line 1206 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::solveStress().

1207 {
1208  _fp_prev_inv = _fp_inv; // update _fp_prev_inv
1209 }
bool FiniteStrainCrystalPlasticity::line_search_update ( const Real  rnorm_prev,
const RankTwoTensor  dpk2 
)
protectedinherited

This function performs the line search update.

Definition at line 1123 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::solveStress().

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

This function performs the line search update.

Definition at line 267 of file FiniteStrainCPSlipRateRes.C.

Referenced by solveStress().

269 {
270  if (_lsrch_method == "CUT_HALF")
271  {
272  Real rnorm;
273  Real step = 1.0;
274  do
275  {
276  for (unsigned int i = 0; i < update.size(); ++i)
277  _slip_rate(i) += step * update(i);
278 
279  step /= 2.0;
280 
281  for (unsigned int i = 0; i < update.size(); ++i)
282  _slip_rate(i) -= step * update(i);
283 
285  if (_err_tol)
286  return false;
287  rnorm = calcResidNorm();
288  } while (rnorm > rnorm_prev && step > _min_lsrch_step);
289 
290  if (rnorm > rnorm_prev && step <= _min_lsrch_step)
291  return false;
292 
293  return true;
294  }
295  else if (_lsrch_method == "BISECTION")
296  {
297  unsigned int count = 0;
298  Real step_a = 0.0;
299  Real step_b = 1.0;
300  Real step = 1.0;
301  Real s_m = 1000.0;
302  Real rnorm = 1000.0;
303 
304  Real s_b = calcResidDotProdUpdate(update);
305  Real rnorm1 = calcResidNorm();
306 
307  for (unsigned int i = 0; i < update.size(); ++i)
308  _slip_rate(i) += update(i);
309 
311  Real s_a = calcResidDotProdUpdate(update);
312  Real rnorm0 = calcResidNorm();
313 
314  for (unsigned int i = 0; i < update.size(); ++i)
315  _slip_rate(i) -= update(i);
316 
317  if ((rnorm1 / rnorm0) < _lsrch_tol || s_a * s_b > 0)
318  {
320  return true;
321  }
322 
323  while ((rnorm / rnorm0) > _lsrch_tol && count < _lsrch_max_iter)
324  {
325 
326  for (unsigned int i = 0; i < update.size(); ++i)
327  _slip_rate(i) += step * update(i);
328 
329  step = 0.5 * (step_a + step_b);
330 
331  for (unsigned int i = 0; i < update.size(); ++i)
332  _slip_rate(i) -= step * update(i);
333 
335  s_m = calcResidDotProdUpdate(update);
336  rnorm = calcResidNorm();
337 
338  if (s_m * s_a < 0.0)
339  {
340  step_b = step;
341  s_b = s_m;
342  }
343  if (s_m * s_b < 0.0)
344  {
345  step_a = step;
346  s_a = s_m;
347  }
348  count++;
349  }
350 
351  if ((rnorm / rnorm0) < _lsrch_tol && count < _lsrch_max_iter)
352  return true;
353 
354  return false;
355  }
356  else
357  {
358  mooseError("Line search meothod is not provided.");
359  return false;
360  }
361 }
Real calcResidDotProdUpdate(const DenseVector< Real > &)
This function calculates the dot product of residual and update.
unsigned int _lsrch_max_iter
Line search bisection method maximum iteration number.
virtual Real calcResidNorm()
This function calculates the residual norm.
Real _lsrch_tol
Line search bisection method tolerance.
virtual void calcResidualSlipRate()
This function calculates residual of slip rate.
Real _min_lsrch_step
Minimum line search step size.
void FiniteStrainCrystalPlasticity::postSolveQp ( )
protectedvirtualinherited

This function update stress and internal variable after solve.

Definition at line 581 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::computeQpStress().

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

This function update internal variable after solve.

Definition at line 677 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::solveQp().

678 {
679  if (_max_substep_iter == 1) // No substepping
680  {
681  _gss[_qp] = _gss_tmp;
682  _acc_slip[_qp] = _accslip_tmp;
683  }
684  else
685  {
686  if (_last_step_iter)
687  {
688  _gss[_qp] = _gss_tmp;
689  _acc_slip[_qp] = _accslip_tmp;
690  }
691  else
692  {
695  }
696  }
697 }
MaterialProperty< std::vector< Real > > & _gss
unsigned int _max_substep_iter
Maximum number of substep iterations.
void FiniteStrainCrystalPlasticity::postSolveStress ( )
protectedvirtualinherited

This function update stress and plastic deformation gradient after solve.

Definition at line 796 of file FiniteStrainCrystalPlasticity.C.

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

797 {
798  if (_max_substep_iter == 1) // No substepping
799  {
800  _fp[_qp] = _fp_inv.inverse();
801  _pk2[_qp] = _pk2_tmp;
802  }
803  else
804  {
805  if (_last_step_iter)
806  {
807  _fp[_qp] = _fp_inv.inverse();
808  _pk2[_qp] = _pk2_tmp;
809  }
810  else
811  {
814  }
815  }
816 }
MaterialProperty< RankTwoTensor > & _pk2
MaterialProperty< RankTwoTensor > & _fp
unsigned int _max_substep_iter
Maximum number of substep iterations.
void FiniteStrainCrystalPlasticity::preSolveQp ( )
protectedvirtualinherited

This function set variables for stress and internal variable solve.

Definition at line 553 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::computeQpStress().

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

This function set variables for internal variable solve.

Definition at line 615 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::solveQp().

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

This function sets variable for internal variable solve.

Reimplemented from FiniteStrainCrystalPlasticity.

Definition at line 40 of file FiniteStrainCPSlipRateRes.C.

Referenced by solveStatevar().

41 {
43  _slip_rate.zero();
44 }
virtual void preSolveStress()
This function set variables for stress solve.
void FiniteStrainCrystalPlasticity::readFileFlowRateParams ( )
protectedvirtualinherited

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

Definition at line 308 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::initSlipSysProps().

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

This function read hardness parameters from file.

Definition at line 399 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::initSlipSysProps().

400 {
401 }
void FiniteStrainCrystalPlasticity::readFileInitSlipSysRes ( )
protectedvirtualinherited

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

Definition at line 237 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::initSlipSysProps().

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

This function solves stress and internal variables.

Definition at line 571 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::computeQpStress().

572 {
574  solveStatevar();
575  if (_err_tol)
576  return;
578 }
virtual void solveStatevar()
This function solves internal variables.
virtual void postSolveStatevar()
This function update internal variable after solve.
virtual void preSolveStatevar()
This function set variables for internal variable solve.
void FiniteStrainCPSlipRateRes::solveStatevar ( )
protectedvirtual

This function solves internal variables.

Reimplemented from FiniteStrainCrystalPlasticity.

Definition at line 30 of file FiniteStrainCPSlipRateRes.C.

31 {
33  solveStress();
34  if (_err_tol)
35  return;
37 }
virtual void solveStress()
This function solves for stress, updates plastic deformation gradient.
virtual void preSolveStress()
This function sets variable for internal variable solve.
virtual void postSolveStress()
This function update stress and plastic deformation gradient after solve.
void FiniteStrainCPSlipRateRes::solveStress ( )
protectedvirtual

This function solves for stress, updates plastic deformation gradient.

Reimplemented from FiniteStrainCrystalPlasticity.

Definition at line 47 of file FiniteStrainCPSlipRateRes.C.

Referenced by solveStatevar().

48 {
49  Real rnorm, rnorm0, rnorm_prev;
50  unsigned int iter = 0;
51 
52 #ifdef DEBUG
53  std::vector<Real> rnormst(_maxiter + 1), slipratest(_maxiter + 1); // Use for Debugging
54 #endif
55 
57  if (_err_tol)
58  return;
59  rnorm = calcResidNorm();
60  rnorm0 = rnorm;
61 
62 #ifdef DEBUG
63  rnormst[iter] = rnorm;
64  Real slipratemax = 0.0;
65  for (unsigned int i = 0; i < _nss; ++i)
66  if (std::abs(_slip_rate(i)) > slipratemax)
67  slipratemax = std::abs(_slip_rate(i));
68  slipratest[iter] = slipratemax;
69 #endif
70 
71  while (rnorm > _rtol * rnorm0 && rnorm0 > _abs_tol && iter < _maxiter)
72  {
73  calcUpdate();
74 
75  DenseVector<Real> update = _resid;
76 
77  _slip_rate -= update;
78 
80  if (_err_tol)
81  return;
82  rnorm_prev = rnorm;
83  rnorm = calcResidNorm();
84 
85  if (_use_line_search && rnorm > rnorm_prev && !lineSearchUpdateSlipRate(rnorm_prev, update))
86  {
87 #ifdef DEBUG
88  mooseWarning("FiniteStrainCrystalPLasticity: Failed with line search");
89 #endif
90  _err_tol = true;
91  return;
92  }
93 
95 
96  if (_use_line_search)
97  rnorm = calcResidNorm();
98  iter++;
99 
100 #ifdef DEBUG
101  slipratemax = 0.0;
102  for (unsigned int i = 0; i < _nss; ++i)
103  if (std::abs(_slip_rate(i)) > slipratemax)
104  slipratemax = std::abs(_slip_rate(i));
105  rnormst[iter] = rnorm;
106  slipratest[iter] = slipratemax;
107 #endif
108  }
109 
110  if (iter == _maxiter)
111  {
112 #ifdef DEBUG
113  mooseWarning("FiniteStrainCPSlipRateRes: NR exceeds maximum iteration ", iter, " ", rnorm);
114 #endif
115  _err_tol = true;
116  return;
117  }
118 }
bool _use_line_search
Flag to activate line serach.
Real _rtol
Stress residual equation relative tolerance.
unsigned int _maxiter
Maximum number of iterations for stress update.
const unsigned int _nss
Number of slip system resistance.
void calcUpdate()
This function calculates and updates the residual of slip rate.
virtual Real calcResidNorm()
This function calculates the residual norm.
virtual void calcResidJacobSlipRate()
This function calculates residual and jacobian of slip rate.
virtual void calcJacobianSlipRate()
This function calculates jacobian of slip rate.
virtual void calcResidualSlipRate()
This function calculates residual of slip rate.
bool lineSearchUpdateSlipRate(const Real, const DenseVector< Real > &)
This function performs the line search update.
Real _abs_tol
Stress residual equation absolute tolerance.
void FiniteStrainCrystalPlasticity::update_slip_system_resistance ( )
protectedvirtualinherited

This function updates the slip system resistances.

Definition at line 820 of file FiniteStrainCrystalPlasticity.C.

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

821 {
822  updateGss();
823 }
virtual void updateGss()
This function updates the slip system resistances.
void FiniteStrainCrystalPlasticity::updateGss ( )
protectedvirtualinherited

This function updates the slip system resistances.

Old function to update slip system resistances.

Kept to avoid code break at computeQpstress

Definition at line 830 of file FiniteStrainCrystalPlasticity.C.

Referenced by FiniteStrainCrystalPlasticity::update_slip_system_resistance().

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

Member Data Documentation

DenseVector<Real> FiniteStrainCrystalPlasticity::_a0
protectedinherited
Real FiniteStrainCrystalPlasticity::_abs_tol
protectedinherited

Stress residual equation absolute tolerance.

Definition at line 248 of file FiniteStrainCrystalPlasticity.h.

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

MaterialProperty<Real>& FiniteStrainCrystalPlasticity::_acc_slip
protectedinherited
const MaterialProperty<Real>& FiniteStrainCrystalPlasticity::_acc_slip_old
protectedinherited
Real FiniteStrainCrystalPlasticity::_accslip_tmp
protectedinherited
Real FiniteStrainCrystalPlasticity::_accslip_tmp_old
protectedinherited
const std::string ComputeStressBase::_base_name
protectedinherited
const MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_crysrot
protectedinherited
const MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_deformation_gradient
protectedinherited
const MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_deformation_gradient_old
protectedinherited
RankTwoTensor FiniteStrainCrystalPlasticity::_delta_dfgrd
protectedinherited

Flag to check whether convergence is achieved.

Used for substepping; Uniformly divides the increment in deformation gradient

Definition at line 347 of file FiniteStrainCrystalPlasticity.h.

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

Real FiniteStrainCrystalPlasticity::_dfgrd_scale_factor
protectedinherited

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

Definition at line 349 of file FiniteStrainCrystalPlasticity.h.

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

RankTwoTensor FiniteStrainCrystalPlasticity::_dfgrd_tmp
protectedinherited
RankTwoTensor FiniteStrainCrystalPlasticity::_dfgrd_tmp_old
protectedinherited
DenseMatrix<Real> FiniteStrainCrystalPlasticity::_dgss_dsliprate
protectedinherited
DenseVector<Real> FiniteStrainCrystalPlasticity::_dslipdtau
protectedinherited
DenseVector<Real> FiniteStrainCPSlipRateRes::_dsliprate_dgss
protected

Definition at line 93 of file FiniteStrainCPSlipRateRes.h.

Referenced by calcDgssDsliprate(), and getSlipIncrements().

DenseMatrix<Real> FiniteStrainCPSlipRateRes::_dsliprate_dsliprate
protected
MaterialProperty<RankTwoTensor>& ComputeStressBase::_elastic_strain
protectedinherited
const MaterialProperty<RankFourTensor>& FiniteStrainCrystalPlasticity::_elasticity_tensor
protectedinherited
const std::string ComputeStressBase::_elasticity_tensor_name
protectedinherited
bool FiniteStrainCrystalPlasticity::_err_tol
protectedinherited
const MaterialProperty<RankTwoTensor>& ComputeStressBase::_extra_stress
protectedinherited
RankTwoTensor FiniteStrainCrystalPlasticity::_fe
protectedinherited
bool FiniteStrainCrystalPlasticity::_first_step_iter
protectedinherited
bool FiniteStrainCrystalPlasticity::_first_substep
protectedinherited
std::vector<Real> FiniteStrainCrystalPlasticity::_flowprops
protectedinherited
MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_fp
protectedinherited
RankTwoTensor FiniteStrainCrystalPlasticity::_fp_inv
protectedinherited
const MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_fp_old
protectedinherited
RankTwoTensor FiniteStrainCrystalPlasticity::_fp_old_inv
protectedinherited
RankTwoTensor FiniteStrainCrystalPlasticity::_fp_prev_inv
protectedinherited
bool FiniteStrainCrystalPlasticity::_gen_rndm_stress_flag
protectedinherited
std::vector<Real> FiniteStrainCrystalPlasticity::_gprops
protectedinherited
MaterialProperty<std::vector<Real> >& FiniteStrainCrystalPlasticity::_gss
protectedinherited
const MaterialProperty<std::vector<Real> >& FiniteStrainCrystalPlasticity::_gss_old
protectedinherited
std::vector<Real> FiniteStrainCrystalPlasticity::_gss_tmp
protectedinherited
std::vector<Real> FiniteStrainCrystalPlasticity::_gss_tmp_old
protectedinherited
Real FiniteStrainCrystalPlasticity::_gtol
protectedinherited

Internal variable update equation tolerance.

Definition at line 250 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity::solveStatevar().

Real FiniteStrainCrystalPlasticity::_h0
protectedinherited
std::vector<Real> FiniteStrainCrystalPlasticity::_hprops
protectedinherited
MaterialProperty<RankTwoTensor>* ComputeStressBase::_initial_stress
protectedinherited

Initial stress, if provided. InitialStress Deprecation: remove this.

Definition at line 74 of file ComputeStressBase.h.

Referenced by ComputeStressBase::initQpStatefulProperties().

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

initial stress components

Definition at line 62 of file ComputeStressBase.h.

Referenced by ComputeStressBase::ComputeStressBase(), and ComputeStressBase::initQpStatefulProperties().

const MaterialProperty<RankTwoTensor>* ComputeStressBase::_initial_stress_old
protectedinherited

Old value of initial stress, which is needed to correctly implement finite-strain rotations. InitialStress Deprecation: remove this.

Definition at line 77 of file ComputeStressBase.h.

const bool ComputeStressBase::_initial_stress_provided
protectedinherited
bool FiniteStrainCrystalPlasticity::_input_rndm_scale_var
protectedinherited

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

Definition at line 274 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity::postSolveQp().

MooseEnum FiniteStrainCrystalPlasticity::_intvar_read_type
protectedinherited

Read from options for initial values of internal variables.

Definition at line 266 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity::FiniteStrainCrystalPlasticity(), and FiniteStrainCrystalPlasticity::initSlipSysProps().

DenseMatrix<Real> FiniteStrainCPSlipRateRes::_jacob
protected

Definition at line 94 of file FiniteStrainCPSlipRateRes.h.

Referenced by calcJacobianSlipRate(), and calcUpdate().

MaterialProperty<RankFourTensor>& ComputeStressBase::_Jacobian_mult
protectedinherited
MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_lag_e
protectedinherited
const MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_lag_e_old
protectedinherited

Definition at line 305 of file FiniteStrainCrystalPlasticity.h.

bool FiniteStrainCrystalPlasticity::_last_step_iter
protectedinherited
unsigned int FiniteStrainCrystalPlasticity::_lsrch_max_iter
protectedinherited

Line search bisection method maximum iteration number.

Definition at line 295 of file FiniteStrainCrystalPlasticity.h.

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

MooseEnum FiniteStrainCrystalPlasticity::_lsrch_method
protectedinherited
Real FiniteStrainCrystalPlasticity::_lsrch_tol
protectedinherited

Line search bisection method tolerance.

Definition at line 292 of file FiniteStrainCrystalPlasticity.h.

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

unsigned int FiniteStrainCrystalPlasticity::_max_substep_iter
protectedinherited
unsigned int FiniteStrainCrystalPlasticity::_maxiter
protectedinherited

Maximum number of iterations for stress update.

Definition at line 255 of file FiniteStrainCrystalPlasticity.h.

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

unsigned int FiniteStrainCrystalPlasticity::_maxiterg
protectedinherited

Maximum number of iterations for internal variable update.

Definition at line 257 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity::solveStatevar().

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_mechanical_strain
protectedinherited
Real FiniteStrainCrystalPlasticity::_min_lsrch_step
protectedinherited

Minimum line search step size.

Definition at line 289 of file FiniteStrainCrystalPlasticity.h.

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

DenseVector<Real> FiniteStrainCrystalPlasticity::_mo
protectedinherited
DenseVector<Real> FiniteStrainCrystalPlasticity::_no
protectedinherited
const unsigned int FiniteStrainCrystalPlasticity::_nss
protectedinherited
unsigned int FiniteStrainCrystalPlasticity::_num_slip_sys_flowrate_props
protectedinherited

Number of slip system flow rate parameters.

Definition at line 260 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity::getFlowRateParams(), and FiniteStrainCrystalPlasticity::readFileFlowRateParams().

unsigned int FiniteStrainCrystalPlasticity::_num_slip_sys_props
protectedinherited

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

Definition at line 269 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity::FiniteStrainCrystalPlasticity(), and FiniteStrainCrystalPlasticity::getSlipSystems().

MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_pk2
protectedinherited
const MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_pk2_old
protectedinherited
RankTwoTensor FiniteStrainCrystalPlasticity::_pk2_tmp
protectedinherited
RankTwoTensor FiniteStrainCrystalPlasticity::_pk2_tmp_old
protectedinherited
Real FiniteStrainCrystalPlasticity::_r
protectedinherited
bool FiniteStrainCrystalPlasticity::_read_from_slip_sys_file
protectedinherited
DenseVector<Real> FiniteStrainCPSlipRateRes::_resid
protected
Real FiniteStrainCrystalPlasticity::_rndm_scale_var
protectedinherited

Scaling value.

Definition at line 277 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity::postSolveQp().

unsigned int FiniteStrainCrystalPlasticity::_rndm_seed
protectedinherited
Real FiniteStrainCrystalPlasticity::_rtol
protectedinherited

Stress residual equation relative tolerance.

Definition at line 246 of file FiniteStrainCrystalPlasticity.h.

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

std::vector<RankTwoTensor> FiniteStrainCrystalPlasticity::_s0
protectedinherited
DenseVector<Real> FiniteStrainCrystalPlasticity::_slip_incr
protectedinherited
Real FiniteStrainCrystalPlasticity::_slip_incr_tol
protectedinherited

Slip increment tolerance.

Definition at line 252 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity::getSlipIncrements().

DenseVector<Real> FiniteStrainCPSlipRateRes::_slip_rate
protected
std::string FiniteStrainCrystalPlasticity::_slip_sys_file_name
protectedinherited

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

Definition at line 231 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity::getSlipSystems().

std::string FiniteStrainCrystalPlasticity::_slip_sys_flow_prop_file_name
protectedinherited

File should contain values of the flow rate equation parameters.

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

Definition at line 240 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity::initSlipSysProps(), and FiniteStrainCrystalPlasticity::readFileFlowRateParams().

std::string FiniteStrainCrystalPlasticity::_slip_sys_hard_prop_file_name
protectedinherited

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

Definition at line 243 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity::initSlipSysProps().

DenseVector<Real> FiniteStrainCrystalPlasticity::_slip_sys_props
protectedinherited
std::string FiniteStrainCrystalPlasticity::_slip_sys_res_prop_file_name
protectedinherited

File should contain initial values of the slip system resistances.

Definition at line 234 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity::readFileInitSlipSysRes().

const bool ComputeStressBase::_store_stress_old
protectedinherited

Parameter which decides whether to store old stress. This is required for HHT time integration and Rayleigh damping.

Definition at line 68 of file ComputeStressBase.h.

MaterialProperty<RankTwoTensor>& ComputeStressBase::_stress
protectedinherited

Definition at line 53 of file ComputeStressBase.h.

Referenced by ComputeStressBase::addQpInitialStress(), ComputeMultipleInelasticCosseratStress::computeAdmissibleState(), ComputeMultipleInelasticStress::computeAdmissibleState(), ComputeStressBase::computeQpProperties(), ComputeStressEosBase::computeQpProperties(), ComputeBirchMurnaghanEquationOfStress::computeQpProperties(), ComputeLinearElasticStress::computeQpStress(), ComputeStrainIncrementBasedStress::computeQpStress(), ComputeLinearElasticPFFractureStress::computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeIsotropicLinearElasticPFFractureStress::computeQpStress(), ComputeFiniteStrainElasticStressBirchMurnaghan::computeQpStress(), FiniteStrainPlasticMaterial::computeQpStress(), ComputeMultiPlasticityStress::computeQpStress(), ComputeSmearedCrackingStress::computeQpStress(), ComputeLinearViscoelasticStress::computeQpStress(), ComputeMultipleInelasticStress::computeQpStressIntermediateConfiguration(), ComputeSmearedCrackingStress::crackingStressRotation(), ComputeMultipleInelasticStress::finiteStrainRotation(), ComputeStressBase::initQpStatefulProperties(), FiniteStrainCrystalPlasticity::initQpStatefulProperties(), FiniteStrainUObasedCP::initQpStatefulProperties(), FiniteStrainHyperElasticViscoPlastic::initQpStatefulProperties(), ComputeMultiPlasticityStress::postReturnMap(), FiniteStrainUObasedCP::postSolveQp(), FiniteStrainHyperElasticViscoPlastic::postSolveQp(), FiniteStrainCrystalPlasticity::postSolveQp(), ComputeMultipleInelasticStress::updateQpState(), ComputeMultipleInelasticStress::updateQpStateSingleModel(), and LinearIsoElasticPFDamage::updateVar().

MooseEnum FiniteStrainCrystalPlasticity::_tan_mod_type
protectedinherited

Type of tangent moduli calculation.

Definition at line 263 of file FiniteStrainCrystalPlasticity.h.

Referenced by FiniteStrainCrystalPlasticity::calcTangentModuli().

DenseVector<Real> FiniteStrainCrystalPlasticity::_tau
protectedinherited
Real FiniteStrainCrystalPlasticity::_tau_init
protectedinherited
Real FiniteStrainCrystalPlasticity::_tau_sat
protectedinherited
MaterialProperty<RankTwoTensor>& FiniteStrainCrystalPlasticity::_update_rot
protectedinherited
bool FiniteStrainCrystalPlasticity::_use_line_search
protectedinherited

Flag to activate line serach.

Definition at line 286 of file FiniteStrainCrystalPlasticity.h.

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

DenseVector<Real> FiniteStrainCrystalPlasticity::_xm
protectedinherited

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