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

This class solves visco plastic model based on isotropically damaged stress The damage parameter is obtained from phase field fracture kernel Computes undamaged elastic strain energy and associated tensors used in phase field fracture kernel. More...

#include <HyperElasticPhaseFieldIsoDamage.h>

Inheritance diagram for HyperElasticPhaseFieldIsoDamage:
[legend]

Public Member Functions

 HyperElasticPhaseFieldIsoDamage (const InputParameters &parameters)
 

Protected Member Functions

virtual void computePK2StressAndDerivative ()
 This function computes PK2 stress. More...
 
virtual void computeDamageStress ()
 This function computes PK2 stress modified to account for damage Computes numerical stiffness if flag is true Computes undamaged elastic strain energy and associated tensors. More...
 
virtual void computeNumStiffness ()
 This function computes numerical stiffness. More...
 
virtual void computeQpJacobian ()
 This function computes tensors used to construct diagonal and off-diagonal Jacobian. More...
 
virtual void initUOVariables ()
 This function initializes the properties, stateful properties and user objects The properties and stateful properties associated with user objects are only initialized here The properties have the same name as the user object name. More...
 
void initNumUserObjects (const std::vector< UserObjectName > &, unsigned int &)
 This function calculates the number of each user object type. More...
 
template<typename T >
void initProp (const std::vector< UserObjectName > &, unsigned int, std::vector< MaterialProperty< T > *> &)
 This function initializes properties for each user object. More...
 
template<typename T >
void initPropOld (const std::vector< UserObjectName > &, unsigned int, std::vector< const MaterialProperty< T > *> &)
 This function initializes old for stateful properties associated with user object Only user objects that update internal variables have an associated old property. More...
 
template<typename T >
void initUserObjects (const std::vector< UserObjectName > &, unsigned int, std::vector< const T *> &)
 This function initializes user objects. More...
 
virtual void initJacobianVariables ()
 This function initialize variables required for Jacobian calculation. More...
 
virtual void initQpStatefulProperties ()
 Initializes state. More...
 
virtual void computeQpStress ()
 This function computes the Cauchy stress. More...
 
virtual void saveOldState ()
 This function saves the old stateful properties that is modified during sub stepping. More...
 
virtual void preSolveQp ()
 Sets state for solve. More...
 
virtual bool solveQp ()
 Solve state. More...
 
virtual void postSolveQp ()
 Update state for output (Outside substepping) More...
 
virtual void recoverOldState ()
 This function restores the the old stateful properties after a successful solve. More...
 
virtual void preSolveFlowrate ()
 Sets state for solve (Inside substepping) More...
 
virtual bool solveFlowrate ()
 Solve for flow rate and state. More...
 
virtual void postSolveFlowrate ()
 Update state for output (Inside substepping) More...
 
virtual bool computeFlowRateFunction ()
 Calls user objects to compute flow rates. More...
 
virtual bool computeFlowDirection ()
 Calls user objects to compute flow directions. More...
 
virtual void computeElasticRightCauchyGreenTensor ()
 Computes elastic Right Cauchy Green Tensor. More...
 
virtual void computeElasticStrain ()
 Computes elastic Lagrangian strain. More...
 
virtual void computeDeeDce ()
 Computes derivative of elastic strain w.r.t elastic Right Cauchy Green Tensor. More...
 
virtual bool computeFlowRateResidual ()
 Computes flow rate residual vector. More...
 
virtual void computeFlowRateJacobian ()
 Computes flow rate Jacobian matrix. More...
 
virtual void computeElasticPlasticDeformGrad ()
 Computes elastic and plastic deformation gradients. More...
 
virtual Real computeNorm (const std::vector< Real > &)
 Computes norm of residual vector. More...
 
virtual void updateFlowRate ()
 Update flow rate. More...
 
virtual void computeDpk2Dfpinv ()
 Computes derivative of PK2 stress wrt inverse of plastic deformation gradient. More...
 
virtual bool computeIntVarRates ()
 This function call user objects to calculate rate of internal variables. More...
 
virtual bool computeIntVar ()
 This function call user objects to integrate internal variables. More...
 
virtual bool computeStrength ()
 This function call user objects to compute strength. More...
 
virtual void computeIntVarRateDerivatives ()
 This function call user objects to compute dintvar_rate/dintvar and dintvarrate/dflowrate. More...
 
virtual void computeIntVarDerivatives ()
 This function call user objects to compute dintvar/dintvar_rate and dintvar/dflowrate. More...
 
void computeStrengthDerivatives ()
 This function call user objects to compute dstrength/dintvar. More...
 
virtual void computeQpProperties () override
 

Protected Attributes

bool _num_stiffness
 Flag to compute numerical stiffness. More...
 
Real _kdamage
 Small stiffness of completely damaged material point. More...
 
bool _use_current_hist
 Use current value of history variable. More...
 
const MaterialProperty< Real > & _l
 Material property defining crack width, declared elsewhere. More...
 
const MaterialProperty< Real > & _gc
 Material property defining gc parameter, declared elsewhere. More...
 
Real _zero_tol
 Used in numerical stiffness calculation to check near zero values. More...
 
Real _zero_pert
 Perturbation value for near zero or zero strain components. More...
 
Real _pert_val
 Perturbation value for strain components. More...
 
const VariableValue & _c
 Compupled damage variable. More...
 
bool _save_state
 Flag to save couple material properties. More...
 
MaterialProperty< RankTwoTensor > & _dstress_dc
 
RankTwoTensor _pk2_tmp
 
RankTwoTensor _dG0_dee
 
RankTwoTensor _dpk2_dc
 
RankFourTensor _dpk2_dee
 
std::vector< RankTwoTensor > _etens
 
MaterialProperty< Real > & _F
 Elastic energy and derivatives, declared in this material. More...
 
MaterialProperty< Real > & _dFdc
 
MaterialProperty< Real > & _d2Fdc2
 
MaterialProperty< RankTwoTensor > & _d2Fdcdstrain
 
MaterialProperty< Real > & _hist
 History variable that prevents crack healing, declared in this material. More...
 
const MaterialProperty< Real > & _hist_old
 Old value of history variable. More...
 
Real _resid_abs_tol
 Absolute tolerance for residual convergence check. More...
 
Real _resid_rel_tol
 Relative tolerance for residual convergence check. More...
 
unsigned int _maxiters
 Maximum number of iterations. More...
 
unsigned int _max_substep_iter
 Maximum number of substep iterations. More...
 
std::vector< UserObjectName > _flow_rate_uo_names
 Names of flow rate user objects. More...
 
std::vector< UserObjectName > _strength_uo_names
 Names of strength user objects. More...
 
std::vector< UserObjectName > _int_var_uo_names
 Names of internal variable user objects. More...
 
std::vector< UserObjectName > _int_var_rate_uo_names
 Names of internal variable rate user objects. More...
 
unsigned int _num_flow_rate_uos
 Number of flow rate user objects. More...
 
unsigned int _num_strength_uos
 Number of strength user objects. More...
 
unsigned int _num_int_var_uos
 Number of internal variable user objects. More...
 
unsigned int _num_int_var_rate_uos
 Number of internal variable rate user objects. More...
 
std::vector< const HEVPFlowRateUOBase * > _flow_rate_uo
 Flow rate user objects. More...
 
std::vector< const HEVPStrengthUOBase * > _strength_uo
 Strength user objects. More...
 
std::vector< const HEVPInternalVarUOBase * > _int_var_uo
 Internal variable user objects. More...
 
std::vector< const HEVPInternalVarRateUOBase * > _int_var_rate_uo
 Internal variable rate user objects. More...
 
std::string _pk2_prop_name
 
MaterialProperty< RankTwoTensor > & _pk2
 
MaterialProperty< RankTwoTensor > & _fp
 
const MaterialProperty< RankTwoTensor > & _fp_old
 
MaterialProperty< RankTwoTensor > & _ce
 
const MaterialProperty< RankTwoTensor > & _deformation_gradient
 
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
 
const MaterialProperty< RankTwoTensor > & _rotation_increment
 
std::vector< MaterialProperty< Real > * > _flow_rate_prop
 
std::vector< MaterialProperty< Real > * > _strength_prop
 
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop
 
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old
 
std::vector< MaterialProperty< Real > * > _int_var_rate_prop
 
std::vector< Real > _int_var_old
 
RankTwoTensor _dfgrd_tmp
 
RankTwoTensor _fp_tmp_inv
 
RankTwoTensor _fp_tmp_old_inv
 
RankTwoTensor _fe
 
RankTwoTensor _ee
 
RankTwoTensor _pk2_fet
 
RankTwoTensor _fe_pk2
 
RankFourTensor _dpk2_dce
 
RankFourTensor _dpk2_dfe
 
RankFourTensor _dfe_dfpinv
 
RankFourTensor _dpk2_dfpinv
 
RankFourTensor _dee_dce
 
RankFourTensor _dce_dfe
 
RankFourTensor _dfe_df
 
RankFourTensor _tan_mod
 
RankFourTensor _df_dstretch_inc
 
std::vector< RankTwoTensor > _flow_dirn
 
std::vector< RankTwoTensor > _dflowrate_dpk2
 
std::vector< RankTwoTensor > _dpk2_dflowrate
 
std::vector< RankTwoTensor > _dfpinv_dflowrate
 
DenseVector< Real > _dflow_rate
 
DenseVector< Real > _flow_rate
 
DenseVector< Real > _resid
 
std::vector< DenseVector< Real > > _dintvarrate_dflowrate
 Jacobian variables. More...
 
std::vector< DenseVector< Real > > _dintvar_dflowrate_tmp
 
DenseMatrix< Real > _dintvarrate_dintvar
 
DenseMatrix< Real > _dintvar_dintvarrate
 
DenseMatrix< Real > _dintvar_dintvar
 
DenseMatrix< Real > _dintvar_dflowrate
 
DenseMatrix< Real > _dstrength_dintvar
 
DenseMatrix< Real > _dflowrate_dstrength
 
DenseVector< Real > _dintvar_dintvar_x
 
DenseMatrix< Real > _jac
 
Real _dt_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< 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

This class solves visco plastic model based on isotropically damaged stress The damage parameter is obtained from phase field fracture kernel Computes undamaged elastic strain energy and associated tensors used in phase field fracture kernel.

Definition at line 26 of file HyperElasticPhaseFieldIsoDamage.h.

Constructor & Destructor Documentation

◆ HyperElasticPhaseFieldIsoDamage()

HyperElasticPhaseFieldIsoDamage::HyperElasticPhaseFieldIsoDamage ( const InputParameters &  parameters)

Definition at line 37 of file HyperElasticPhaseFieldIsoDamage.C.

39  _num_stiffness(getParam<bool>("numerical_stiffness")),
40  _kdamage(getParam<Real>("damage_stiffness")),
41  _use_current_hist(getParam<bool>("use_current_history_variable")),
42  _l(getMaterialProperty<Real>("l")),
43  _gc(getMaterialProperty<Real>("gc_prop")),
44  _zero_tol(getParam<Real>("zero_tol")),
45  _zero_pert(getParam<Real>("zero_perturb")),
46  _pert_val(getParam<Real>("perturbation_scale_factor")),
47  _c(coupledValue("c")),
48  _save_state(false),
50  declarePropertyDerivative<RankTwoTensor>(_base_name + "stress", getVar("c", 0)->name())),
51  _etens(LIBMESH_DIM),
52  _F(declareProperty<Real>(getParam<MaterialPropertyName>("F_name"))),
53  _dFdc(declarePropertyDerivative<Real>(getParam<MaterialPropertyName>("F_name"),
54  getVar("c", 0)->name())),
55  _d2Fdc2(declarePropertyDerivative<Real>(
56  getParam<MaterialPropertyName>("F_name"), getVar("c", 0)->name(), getVar("c", 0)->name())),
57  _d2Fdcdstrain(declareProperty<RankTwoTensor>("d2Fdcdstrain")),
58  _hist(declareProperty<Real>("hist")),
59  _hist_old(getMaterialPropertyOld<Real>("hist"))
60 {
61 }
const MaterialProperty< Real > & _hist_old
Old value of history variable.
bool _use_current_hist
Use current value of history variable.
Real _zero_pert
Perturbation value for near zero or zero strain components.
const std::string name
Definition: Setup.h:22
bool _save_state
Flag to save couple material properties.
const VariableValue & _c
Compupled damage variable.
Real _pert_val
Perturbation value for strain components.
Real _zero_tol
Used in numerical stiffness calculation to check near zero values.
const MaterialProperty< Real > & _l
Material property defining crack width, declared elsewhere.
const std::string _base_name
MaterialProperty< RankTwoTensor > & _dstress_dc
MaterialProperty< RankTwoTensor > & _d2Fdcdstrain
bool _num_stiffness
Flag to compute numerical stiffness.
MaterialProperty< Real > & _F
Elastic energy and derivatives, declared in this material.
FiniteStrainHyperElasticViscoPlastic(const InputParameters &parameters)
Real _kdamage
Small stiffness of completely damaged material point.
MaterialProperty< Real > & _hist
History variable that prevents crack healing, declared in this material.
const MaterialProperty< Real > & _gc
Material property defining gc parameter, declared elsewhere.

Member Function Documentation

◆ computeDamageStress()

void HyperElasticPhaseFieldIsoDamage::computeDamageStress ( )
protectedvirtual

This function computes PK2 stress modified to account for damage Computes numerical stiffness if flag is true Computes undamaged elastic strain energy and associated tensors.

Definition at line 92 of file HyperElasticPhaseFieldIsoDamage.C.

Referenced by computeNumStiffness(), and computePK2StressAndDerivative().

93 {
94  Real lambda = _elasticity_tensor[_qp](0, 0, 1, 1);
95  Real mu = _elasticity_tensor[_qp](0, 1, 0, 1);
96 
97  Real c = _c[_qp];
98  Real xfac = Utility::pow<2>(1.0 - c) + _kdamage;
99 
100  std::vector<Real> eigval;
101  RankTwoTensor evec;
102  _ee.symmetricEigenvaluesEigenvectors(eigval, evec);
103 
104  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
105  _etens[i].vectorOuterProduct(evec.column(i), evec.column(i));
106 
107  Real etr = 0.0;
108  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
109  etr += eigval[i];
110 
111  Real etrpos = (std::abs(etr) + etr) / 2.0;
112  Real etrneg = (std::abs(etr) - etr) / 2.0;
113 
114  RankTwoTensor pk2pos, pk2neg;
115 
116  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
117  {
118  pk2pos += _etens[i] * (lambda * etrpos + 2.0 * mu * (std::abs(eigval[i]) + eigval[i]) / 2.0);
119  pk2neg += _etens[i] * (lambda * etrneg + 2.0 * mu * (std::abs(eigval[i]) - eigval[i]) / 2.0);
120  }
121 
122  _pk2_tmp = pk2pos * xfac - pk2neg;
123 
124  if (_save_state)
125  {
126  std::vector<Real> epos(LIBMESH_DIM), eneg(LIBMESH_DIM);
127  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
128  {
129  epos[i] = (std::abs(eigval[i]) + eigval[i]) / 2.0;
130  eneg[i] = (std::abs(eigval[i]) - eigval[i]) / 2.0;
131  }
132 
133  // sum squares of epos and eneg
134  Real pval(0.0), nval(0.0);
135  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
136  {
137  pval += epos[i] * epos[i];
138  nval += eneg[i] * eneg[i];
139  }
140 
141  // Energy with positive principal strains
142  const Real G0_pos = lambda * etrpos * etrpos / 2.0 + mu * pval;
143  const Real G0_neg = lambda * etrneg * etrneg / 2.0 + mu * nval;
144 
145  // Assign history variable and derivative
146  if (G0_pos > _hist_old[_qp])
147  _hist[_qp] = G0_pos;
148  else
149  _hist[_qp] = _hist_old[_qp];
150 
151  Real hist_variable = _hist_old[_qp];
152  if (_use_current_hist)
153  hist_variable = _hist[_qp];
154 
155  // Elastic free energy density
156  _F[_qp] = hist_variable * xfac - G0_neg + _gc[_qp] / (2 * _l[_qp]) * c * c;
157 
158  // derivative of elastic free energy density wrt c
159  _dFdc[_qp] = -hist_variable * 2.0 * (1.0 - c) * (1 - _kdamage) + _gc[_qp] / _l[_qp] * c;
160 
161  // 2nd derivative of elastic free energy density wrt c
162  _d2Fdc2[_qp] = hist_variable * 2.0 * (1 - _kdamage) + _gc[_qp] / _l[_qp];
163 
164  _dG0_dee = pk2pos;
165 
166  _dpk2_dc = -pk2pos * 2.0 * (1.0 - c);
167  }
168 }
const MaterialProperty< Real > & _hist_old
Old value of history variable.
bool _use_current_hist
Use current value of history variable.
bool _save_state
Flag to save couple material properties.
const VariableValue & _c
Compupled damage variable.
const MaterialProperty< Real > & _l
Material property defining crack width, declared elsewhere.
const MaterialProperty< RankFourTensor > & _elasticity_tensor
MaterialProperty< Real > & _F
Elastic energy and derivatives, declared in this material.
Real _kdamage
Small stiffness of completely damaged material point.
MaterialProperty< Real > & _hist
History variable that prevents crack healing, declared in this material.
const MaterialProperty< Real > & _gc
Material property defining gc parameter, declared elsewhere.

◆ computeDeeDce()

void FiniteStrainHyperElasticViscoPlastic::computeDeeDce ( )
protectedvirtualinherited

Computes derivative of elastic strain w.r.t elastic Right Cauchy Green Tensor.

Definition at line 490 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::initJacobianVariables().

491 {
492  _dee_dce.zero();
493 
494  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
495  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
496  _dee_dce(i, j, i, j) = 0.5;
497 }

◆ computeDpk2Dfpinv()

void FiniteStrainHyperElasticViscoPlastic::computeDpk2Dfpinv ( )
protectedvirtualinherited

Computes derivative of PK2 stress wrt inverse of plastic deformation gradient.

Definition at line 521 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::computeFlowRateJacobian().

522 {
523  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
524  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
525  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
526  _dfe_dfpinv(i, j, k, j) = _dfgrd_tmp(i, k);
527 
529 }

◆ computeElasticPlasticDeformGrad()

void FiniteStrainHyperElasticViscoPlastic::computeElasticPlasticDeformGrad ( )
protectedvirtualinherited

Computes elastic and plastic deformation gradients.

Definition at line 506 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::solveFlowrate().

507 {
508  RankTwoTensor iden;
509  iden.addIa(1.0);
510 
511  RankTwoTensor val;
512  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
513  val += _flow_rate(i) * _flow_dirn[i] * _dt_substep;
514 
515  _fp_tmp_inv = _fp_tmp_old_inv * (iden - val);
516  _fp_tmp_inv = std::pow(_fp_tmp_inv.det(), -1.0 / 3.0) * _fp_tmp_inv;
518 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)

◆ computeElasticRightCauchyGreenTensor()

void FiniteStrainHyperElasticViscoPlastic::computeElasticRightCauchyGreenTensor ( )
protectedvirtualinherited

Computes elastic Right Cauchy Green Tensor.

Definition at line 500 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::computeFlowRateResidual().

501 {
502  _ce[_qp] = _fe.transpose() * _fe;
503 }

◆ computeElasticStrain()

void FiniteStrainHyperElasticViscoPlastic::computeElasticStrain ( )
protectedvirtualinherited

Computes elastic Lagrangian strain.

Definition at line 482 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by computePK2StressAndDerivative(), and FiniteStrainHyperElasticViscoPlastic::computePK2StressAndDerivative().

483 {
484  RankTwoTensor iden;
485  iden.addIa(1.0);
486  _ee = 0.5 * (_ce[_qp] - iden);
487 }

◆ computeFlowDirection()

bool FiniteStrainHyperElasticViscoPlastic::computeFlowDirection ( )
protectedvirtualinherited

Calls user objects to compute flow directions.

Definition at line 439 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::computeFlowRateResidual().

440 {
441  for (unsigned i = 0; i < _num_flow_rate_uos; ++i)
442  {
443  if (!_flow_rate_uo[i]->computeDirection(_qp, _flow_dirn[i]))
444  return false;
445  }
446  return true;
447 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< const HEVPFlowRateUOBase * > _flow_rate_uo
Flow rate user objects.

◆ computeFlowRateFunction()

bool FiniteStrainHyperElasticViscoPlastic::computeFlowRateFunction ( )
protectedvirtualinherited

Calls user objects to compute flow rates.

Definition at line 450 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::computeFlowRateResidual().

451 {
452  Real val = 0;
453  for (unsigned i = 0; i < _num_flow_rate_uos; ++i)
454  {
455  if (_flow_rate_uo[i]->computeValue(_qp, val))
456  _resid(i) = -val;
457  else
458  return false;
459  }
460  return true;
461 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< const HEVPFlowRateUOBase * > _flow_rate_uo
Flow rate user objects.

◆ computeFlowRateJacobian()

void FiniteStrainHyperElasticViscoPlastic::computeFlowRateJacobian ( )
protectedvirtualinherited

Computes flow rate Jacobian matrix.

Definition at line 401 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::solveFlowrate().

402 {
406 
407  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
408  for (unsigned int j = 0; j < _num_strength_uos; ++j)
409  _flow_rate_uo[i]->computeDerivative(_qp, _strength_uo_names[j], _dflowrate_dstrength(i, j));
410 
411  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
412  _flow_rate_uo[i]->computeTensorDerivative(_qp, _pk2_prop_name, _dflowrate_dpk2[i]);
413 
415 
416  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
417  {
420  }
421 
422  DenseMatrix<Real> dflowrate_dflowrate;
423  dflowrate_dflowrate = _dflowrate_dstrength;
424  dflowrate_dflowrate.right_multiply(_dstrength_dintvar);
425  dflowrate_dflowrate.right_multiply(_dintvar_dflowrate);
426 
427  _jac.zero();
428  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
429  for (unsigned int j = 0; j < _num_flow_rate_uos; ++j)
430  {
431  if (i == j)
432  _jac(i, j) = 1;
433  _jac(i, j) -= dflowrate_dflowrate(i, j);
434  _jac(i, j) -= _dflowrate_dpk2[i].doubleContraction(_dpk2_dflowrate[j]);
435  }
436 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< UserObjectName > _strength_uo_names
Names of strength user objects.
virtual void computeDpk2Dfpinv()
Computes derivative of PK2 stress wrt inverse of plastic deformation gradient.
void computeStrengthDerivatives()
This function call user objects to compute dstrength/dintvar.
virtual void computeIntVarRateDerivatives()
This function call user objects to compute dintvar_rate/dintvar and dintvarrate/dflowrate.
virtual void computeIntVarDerivatives()
This function call user objects to compute dintvar/dintvar_rate and dintvar/dflowrate.
unsigned int _num_strength_uos
Number of strength user objects.
std::vector< const HEVPFlowRateUOBase * > _flow_rate_uo
Flow rate user objects.

◆ computeFlowRateResidual()

bool FiniteStrainHyperElasticViscoPlastic::computeFlowRateResidual ( )
protectedvirtualinherited

Computes flow rate residual vector.

Definition at line 375 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::solveFlowrate().

376 {
377  if (!computeIntVarRates())
378  return false;
379 
380  if (!computeIntVar())
381  return false;
382 
383  if (!computeStrength())
384  return false;
385 
388 
390  return false;
391 
392  if (!computeFlowDirection())
393  return false;
394 
395  _resid += _flow_rate;
396 
397  return true;
398 }
virtual bool computeFlowDirection()
Calls user objects to compute flow directions.
virtual void computeElasticRightCauchyGreenTensor()
Computes elastic Right Cauchy Green Tensor.
virtual void computePK2StressAndDerivative()
Computes PK2 stress and derivative w.r.t elastic Right Cauchy Green Tensor.
virtual bool computeFlowRateFunction()
Calls user objects to compute flow rates.
virtual bool computeIntVar()
This function call user objects to integrate internal variables.
virtual bool computeStrength()
This function call user objects to compute strength.
virtual bool computeIntVarRates()
This function call user objects to calculate rate of internal variables.

◆ computeIntVar()

bool FiniteStrainHyperElasticViscoPlastic::computeIntVar ( )
protectedvirtualinherited

This function call user objects to integrate internal variables.

Definition at line 597 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::computeFlowRateResidual().

598 {
599  Real val = 0;
600  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
601  {
602  if (_int_var_uo[i]->computeValue(_qp, _dt_substep, val))
603  (*_int_var_stateful_prop[i])[_qp] = val;
604  else
605  return false;
606  }
607  return true;
608 }
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const HEVPInternalVarUOBase * > _int_var_uo
Internal variable user objects.
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop

◆ computeIntVarDerivatives()

void FiniteStrainHyperElasticViscoPlastic::computeIntVarDerivatives ( )
protectedvirtualinherited

This function call user objects to compute dintvar/dintvar_rate and dintvar/dflowrate.

Definition at line 638 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::computeFlowRateJacobian().

639 {
640  Real val = 0;
641 
642  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
643  for (unsigned int j = 0; j < _num_int_var_rate_uos; ++j)
644  {
645  _int_var_uo[i]->computeDerivative(_qp, _dt_substep, _int_var_rate_uo_names[j], val);
646  _dintvar_dintvarrate(i, j) = val;
647  }
648 
649  _dintvar_dintvar.zero();
650 
651  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
652  for (unsigned int j = 0; j < _num_int_var_uos; ++j)
653  {
654  if (i == j)
655  _dintvar_dintvar(i, j) = 1;
656  for (unsigned int k = 0; k < _num_int_var_rate_uos; ++k)
658  }
659 
660  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
662 
663  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
664  {
665  _dintvar_dintvar_x.zero();
667  for (unsigned int j = 0; j < _num_int_var_uos; ++j)
669  }
670 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
unsigned int _num_int_var_rate_uos
Number of internal variable rate user objects.
std::vector< DenseVector< Real > > _dintvar_dflowrate_tmp
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< DenseVector< Real > > _dintvarrate_dflowrate
Jacobian variables.
std::vector< const HEVPInternalVarUOBase * > _int_var_uo
Internal variable user objects.
std::vector< UserObjectName > _int_var_rate_uo_names
Names of internal variable rate user objects.

◆ computeIntVarRateDerivatives()

void FiniteStrainHyperElasticViscoPlastic::computeIntVarRateDerivatives ( )
protectedvirtualinherited

This function call user objects to compute dintvar_rate/dintvar and dintvarrate/dflowrate.

Definition at line 625 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::computeFlowRateJacobian().

626 {
627  Real val = 0;
628 
629  for (unsigned int i = 0; i < _num_int_var_rate_uos; ++i)
630  for (unsigned int j = 0; j < _num_flow_rate_uos; ++j)
631  {
632  _int_var_rate_uo[i]->computeDerivative(_qp, _flow_rate_uo_names[j], val);
633  _dintvarrate_dflowrate[j](i) = val;
634  }
635 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
unsigned int _num_int_var_rate_uos
Number of internal variable rate user objects.
std::vector< DenseVector< Real > > _dintvarrate_dflowrate
Jacobian variables.
std::vector< UserObjectName > _flow_rate_uo_names
Names of flow rate user objects.
std::vector< const HEVPInternalVarRateUOBase * > _int_var_rate_uo
Internal variable rate user objects.

◆ computeIntVarRates()

bool FiniteStrainHyperElasticViscoPlastic::computeIntVarRates ( )
protectedvirtualinherited

This function call user objects to calculate rate of internal variables.

Definition at line 583 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::computeFlowRateResidual().

584 {
585  Real val = 0;
586  for (unsigned int i = 0; i < _num_int_var_rate_uos; ++i)
587  {
588  if (_int_var_rate_uo[i]->computeValue(_qp, val))
589  (*_int_var_rate_prop[i])[_qp] = val;
590  else
591  return false;
592  }
593  return true;
594 }
unsigned int _num_int_var_rate_uos
Number of internal variable rate user objects.
std::vector< const HEVPInternalVarRateUOBase * > _int_var_rate_uo
Internal variable rate user objects.
std::vector< MaterialProperty< Real > * > _int_var_rate_prop

◆ computeNorm()

Real FiniteStrainHyperElasticViscoPlastic::computeNorm ( const std::vector< Real > &  var)
protectedvirtualinherited

Computes norm of residual vector.

Definition at line 532 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::solveFlowrate().

533 {
534  Real val = 0.0;
535  for (unsigned int i = 0; i < var.size(); ++i)
536  val += Utility::pow<2>(var[i]);
537  return std::sqrt(val);
538 }

◆ computeNumStiffness()

void HyperElasticPhaseFieldIsoDamage::computeNumStiffness ( )
protectedvirtual

This function computes numerical stiffness.

Definition at line 171 of file HyperElasticPhaseFieldIsoDamage.C.

Referenced by computePK2StressAndDerivative().

172 {
173  RankTwoTensor ee_tmp;
174 
175  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
176  for (unsigned int j = i; j < LIBMESH_DIM; ++j)
177  {
178  Real ee_pert = _zero_pert;
179  if (std::abs(_ee(i, j)) > _zero_tol)
180  ee_pert = _pert_val * std::abs(_ee(i, j));
181 
182  ee_tmp = _ee;
183  _ee(i, j) += ee_pert;
184 
186 
187  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
188  for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
189  {
190  _dpk2_dee(k, l, i, j) = (_pk2_tmp(k, l) - _pk2[_qp](k, l)) / ee_pert;
191  _dpk2_dee(k, l, j, i) = (_pk2_tmp(k, l) - _pk2[_qp](k, l)) / ee_pert;
192  }
193  _ee = ee_tmp;
194  }
195 }
virtual void computeDamageStress()
This function computes PK2 stress modified to account for damage Computes numerical stiffness if flag...
Real _zero_pert
Perturbation value for near zero or zero strain components.
Real _pert_val
Perturbation value for strain components.
Real _zero_tol
Used in numerical stiffness calculation to check near zero values.

◆ computePK2StressAndDerivative()

void HyperElasticPhaseFieldIsoDamage::computePK2StressAndDerivative ( )
protectedvirtual

This function computes PK2 stress.

Reimplemented from FiniteStrainHyperElasticViscoPlastic.

Definition at line 64 of file HyperElasticPhaseFieldIsoDamage.C.

65 {
67 
68  _save_state = true;
70  _pk2[_qp] = _pk2_tmp;
71 
72  _save_state = false;
73  if (_num_stiffness)
75 
76  if (_num_stiffness)
78 
79  _dce_dfe.zero();
80  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
81  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
82  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
83  {
84  _dce_dfe(i, j, k, i) = _dce_dfe(i, j, k, i) + _fe(k, j);
85  _dce_dfe(i, j, k, j) = _dce_dfe(i, j, k, j) + _fe(k, i);
86  }
87 
89 }
virtual void computeElasticStrain()
Computes elastic Lagrangian strain.
virtual void computeDamageStress()
This function computes PK2 stress modified to account for damage Computes numerical stiffness if flag...
virtual void computeNumStiffness()
This function computes numerical stiffness.
bool _save_state
Flag to save couple material properties.
bool _num_stiffness
Flag to compute numerical stiffness.

◆ computeQpJacobian()

void HyperElasticPhaseFieldIsoDamage::computeQpJacobian ( )
protectedvirtual

This function computes tensors used to construct diagonal and off-diagonal Jacobian.

Reimplemented from FiniteStrainHyperElasticViscoPlastic.

Definition at line 198 of file HyperElasticPhaseFieldIsoDamage.C.

199 {
201 
202  RankTwoTensor dG0_dce = _dee_dce.innerProductTranspose(_dG0_dee);
203  RankTwoTensor dG0_dfe = _dce_dfe.innerProductTranspose(dG0_dce);
204  RankTwoTensor dG0_df = _dfe_df.innerProductTranspose(dG0_dfe);
205 
206  // 2nd derivative wrt c and strain = 0.0 if we used the previous step's history varible
207  if (_use_current_hist)
208  _d2Fdcdstrain[_qp] =
209  -_df_dstretch_inc.innerProductTranspose(dG0_df) * 2.0 * (1.0 - _c[_qp]) * (1 - _kdamage);
210 
211  _dstress_dc[_qp] = _fe.mixedProductIkJl(_fe) * _dpk2_dc;
212 }
bool _use_current_hist
Use current value of history variable.
virtual void computeQpJacobian()
This function computes the Jacobian.
const VariableValue & _c
Compupled damage variable.
MaterialProperty< RankTwoTensor > & _dstress_dc
MaterialProperty< RankTwoTensor > & _d2Fdcdstrain
Real _kdamage
Small stiffness of completely damaged material point.

◆ 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 FiniteStrainHyperElasticViscoPlastic::computeQpStress ( )
protectedvirtualinherited

This function computes the Cauchy stress.

Implements ComputeStressBase.

Definition at line 210 of file FiniteStrainHyperElasticViscoPlastic.C.

211 {
212  bool converge;
213  RankTwoTensor delta_dfgrd = _deformation_gradient[_qp] - _deformation_gradient_old[_qp];
214  unsigned int num_substep = 1;
215  unsigned int substep_iter = 1;
216 
217  saveOldState();
218 
219  do
220  {
221  preSolveQp();
222 
223  converge = true;
224  _dt_substep = _dt / num_substep;
225 
226  for (unsigned int istep = 0; istep < num_substep; ++istep)
227  {
228  _dfgrd_tmp = (istep + 1.0) * delta_dfgrd / num_substep + _deformation_gradient_old[_qp];
229  if (!solveQp())
230  {
231  converge = false;
232  substep_iter++;
233  num_substep *= 2;
234  break;
235  }
236  }
237 
238  if (substep_iter > _max_substep_iter)
239  mooseError("Constitutive failure with substepping at quadrature point ",
240  _q_point[_qp](0),
241  " ",
242  _q_point[_qp](1),
243  " ",
244  _q_point[_qp](2));
245  } while (!converge);
246 
247  postSolveQp();
248 }
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
unsigned int _max_substep_iter
Maximum number of substep iterations.
virtual void postSolveQp()
Update state for output (Outside substepping)
virtual void saveOldState()
This function saves the old stateful properties that is modified during sub stepping.
const MaterialProperty< RankTwoTensor > & _deformation_gradient

◆ computeStrength()

bool FiniteStrainHyperElasticViscoPlastic::computeStrength ( )
protectedvirtualinherited

This function call user objects to compute strength.

Definition at line 611 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::computeFlowRateResidual().

612 {
613  Real val = 0;
614  for (unsigned int i = 0; i < _num_strength_uos; ++i)
615  {
616  if (_strength_uo[i]->computeValue(_qp, val))
617  (*_strength_prop[i])[_qp] = val;
618  else
619  return false;
620  }
621  return true;
622 }
std::vector< const HEVPStrengthUOBase * > _strength_uo
Strength user objects.
unsigned int _num_strength_uos
Number of strength user objects.
std::vector< MaterialProperty< Real > * > _strength_prop

◆ computeStrengthDerivatives()

void FiniteStrainHyperElasticViscoPlastic::computeStrengthDerivatives ( )
protectedinherited

This function call user objects to compute dstrength/dintvar.

Definition at line 673 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::computeFlowRateJacobian().

674 {
675  Real val = 0;
676 
677  for (unsigned int i = 0; i < _num_strength_uos; ++i)
678  for (unsigned int j = 0; j < _num_int_var_uos; ++j)
679  {
680  _strength_uo[i]->computeDerivative(_qp, _int_var_uo_names[j], val);
681  _dstrength_dintvar(i, j) = val;
682  }
683 }
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< UserObjectName > _int_var_uo_names
Names of internal variable user objects.
std::vector< const HEVPStrengthUOBase * > _strength_uo
Strength user objects.
unsigned int _num_strength_uos
Number of strength user objects.

◆ initJacobianVariables()

void FiniteStrainHyperElasticViscoPlastic::initJacobianVariables ( )
protectedvirtualinherited

This function initialize variables required for Jacobian calculation.

Definition at line 154 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::FiniteStrainHyperElasticViscoPlastic().

155 {
157  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
159 
161  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
163 
171 
175 
177 
178  computeDeeDce();
179 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
unsigned int _num_int_var_rate_uos
Number of internal variable rate user objects.
std::vector< DenseVector< Real > > _dintvar_dflowrate_tmp
virtual void computeDeeDce()
Computes derivative of elastic strain w.r.t elastic Right Cauchy Green Tensor.
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< DenseVector< Real > > _dintvarrate_dflowrate
Jacobian variables.
unsigned int _num_strength_uos
Number of strength user objects.

◆ initNumUserObjects()

void FiniteStrainHyperElasticViscoPlastic::initNumUserObjects ( const std::vector< UserObjectName > &  uo_names,
unsigned int &  uo_num 
)
protectedinherited

This function calculates the number of each user object type.

Definition at line 109 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::initUOVariables().

111 {
112  uo_num = uo_names.size();
113 }

◆ initProp()

template<typename T >
void FiniteStrainHyperElasticViscoPlastic::initProp ( const std::vector< UserObjectName > &  uo_names,
unsigned int  uo_num,
std::vector< MaterialProperty< T > *> &  uo_prop 
)
protectedinherited

This function initializes properties for each user object.

Definition at line 117 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::initUOVariables().

120 {
121  uo_prop.resize(uo_num);
122  for (unsigned int i = 0; i < uo_num; ++i)
123  uo_prop[i] = &declareProperty<T>(uo_names[i]);
124 }

◆ initPropOld()

template<typename T >
void FiniteStrainHyperElasticViscoPlastic::initPropOld ( const std::vector< UserObjectName > &  uo_names,
unsigned int  uo_num,
std::vector< const MaterialProperty< T > *> &  uo_prop_old 
)
protectedinherited

This function initializes old for stateful properties associated with user object Only user objects that update internal variables have an associated old property.

Definition at line 128 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::initUOVariables().

132 {
133  uo_prop_old.resize(uo_num);
134  for (unsigned int i = 0; i < uo_num; ++i)
135  uo_prop_old[i] = &getMaterialPropertyOld<T>(uo_names[i]);
136 }

◆ initQpStatefulProperties()

void FiniteStrainHyperElasticViscoPlastic::initQpStatefulProperties ( )
protectedvirtualinherited

Initializes state.

Reimplemented from ComputeStressBase.

Definition at line 182 of file FiniteStrainHyperElasticViscoPlastic.C.

183 {
184  _stress[_qp].zero();
185 
186  _fp[_qp].zero();
187  _fp[_qp].addIa(1.0);
188  _ce[_qp].zero();
189 
190  _pk2[_qp].zero();
191 
192  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
193  (*_flow_rate_prop[i])[_qp] = 0.0;
194 
195  for (unsigned int i = 0; i < _num_strength_uos; ++i)
196  (*_strength_prop[i])[_qp] = 0.0;
197 
198  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
199  {
200  (*_int_var_stateful_prop[i])[_qp] = 0.0;
201  // TODO: remove this nasty const_cast if you can figure out how
202  const_cast<MaterialProperty<Real> &>(*_int_var_stateful_prop_old[i])[_qp] = 0.0;
203  }
204 
205  for (unsigned int i = 0; i < _num_int_var_rate_uos; ++i)
206  (*_int_var_rate_prop[i])[_qp] = 0.0;
207 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< MaterialProperty< Real > * > _flow_rate_prop
unsigned int _num_int_var_rate_uos
Number of internal variable rate user objects.
MaterialProperty< RankTwoTensor > & _stress
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old
unsigned int _num_strength_uos
Number of strength user objects.
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop
std::vector< MaterialProperty< Real > * > _int_var_rate_prop
std::vector< MaterialProperty< Real > * > _strength_prop

◆ initUOVariables()

void FiniteStrainHyperElasticViscoPlastic::initUOVariables ( )
protectedvirtualinherited

This function initializes the properties, stateful properties and user objects The properties and stateful properties associated with user objects are only initialized here The properties have the same name as the user object name.

Definition at line 86 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::FiniteStrainHyperElasticViscoPlastic().

87 {
92 
97 
99 
104 
105  _int_var_old.resize(_num_int_var_uos, 0.0);
106 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< MaterialProperty< Real > * > _flow_rate_prop
void initPropOld(const std::vector< UserObjectName > &, unsigned int, std::vector< const MaterialProperty< T > *> &)
This function initializes old for stateful properties associated with user object Only user objects t...
unsigned int _num_int_var_rate_uos
Number of internal variable rate user objects.
void initNumUserObjects(const std::vector< UserObjectName > &, unsigned int &)
This function calculates the number of each user object type.
std::vector< UserObjectName > _strength_uo_names
Names of strength user objects.
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< UserObjectName > _flow_rate_uo_names
Names of flow rate user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old
std::vector< UserObjectName > _int_var_uo_names
Names of internal variable user objects.
std::vector< const HEVPInternalVarUOBase * > _int_var_uo
Internal variable user objects.
std::vector< const HEVPStrengthUOBase * > _strength_uo
Strength user objects.
unsigned int _num_strength_uos
Number of strength user objects.
std::vector< UserObjectName > _int_var_rate_uo_names
Names of internal variable rate user objects.
std::vector< const HEVPInternalVarRateUOBase * > _int_var_rate_uo
Internal variable rate user objects.
void initProp(const std::vector< UserObjectName > &, unsigned int, std::vector< MaterialProperty< T > *> &)
This function initializes properties for each user object.
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop
std::vector< const HEVPFlowRateUOBase * > _flow_rate_uo
Flow rate user objects.
void initUserObjects(const std::vector< UserObjectName > &, unsigned int, std::vector< const T *> &)
This function initializes user objects.
std::vector< MaterialProperty< Real > * > _int_var_rate_prop
std::vector< MaterialProperty< Real > * > _strength_prop

◆ initUserObjects()

template<typename T >
void FiniteStrainHyperElasticViscoPlastic::initUserObjects ( const std::vector< UserObjectName > &  uo_names,
unsigned int  uo_num,
std::vector< const T *> &  uo 
)
protectedinherited

This function initializes user objects.

Definition at line 140 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::initUOVariables().

143 {
144  uo.resize(uo_num);
145 
146  if (uo_num == 0)
147  mooseError("Specify atleast one user object of type", typeid(T).name());
148 
149  for (unsigned int i = 0; i < uo_num; ++i)
150  uo[i] = &getUserObjectByName<T>(uo_names[i]);
151 }
const std::string name
Definition: Setup.h:22

◆ postSolveFlowrate()

void FiniteStrainHyperElasticViscoPlastic::postSolveFlowrate ( )
protectedvirtualinherited

Update state for output (Inside substepping)

Definition at line 364 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::solveQp().

365 {
367 
368  // TODO: remove this nasty const_cast if you can figure out how
369  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
370  const_cast<MaterialProperty<Real> &>(*_int_var_stateful_prop_old[i])[_qp] =
371  (*_int_var_stateful_prop[i])[_qp];
372 }
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop

◆ postSolveQp()

void FiniteStrainHyperElasticViscoPlastic::postSolveQp ( )
protectedvirtualinherited

Update state for output (Outside substepping)

Definition at line 282 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::computeQpStress().

283 {
284  recoverOldState();
285 
286  _stress[_qp] = _fe * _pk2[_qp] * _fe.transpose() / _fe.det();
287  _fp[_qp] = _fp_tmp_inv.inverse();
288 
290 }
virtual void recoverOldState()
This function restores the the old stateful properties after a successful solve.
MaterialProperty< RankTwoTensor > & _stress
virtual void computeQpJacobian()
This function computes the Jacobian.

◆ preSolveFlowrate()

void FiniteStrainHyperElasticViscoPlastic::preSolveFlowrate ( )
protectedvirtualinherited

Sets state for solve (Inside substepping)

Definition at line 301 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::solveQp().

302 {
303  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
304  {
305  _flow_rate(i) = 0.0;
306  (*_flow_rate_prop[i])[_qp] = 0.0;
307  }
308 
309  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
310  (*_int_var_stateful_prop[i])[_qp] = (*_int_var_stateful_prop_old[i])[_qp];
311 
314 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< MaterialProperty< Real > * > _flow_rate_prop
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop

◆ preSolveQp()

void FiniteStrainHyperElasticViscoPlastic::preSolveQp ( )
protectedvirtualinherited

Sets state for solve.

Definition at line 258 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::computeQpStress().

259 {
260  _fp_tmp_old_inv = _fp_old[_qp].inverse();
261 
262  // TODO: remove this nasty const_cast if you can figure out how
263  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
264  (*_int_var_stateful_prop[i])[_qp] =
265  const_cast<MaterialProperty<Real> &>(*_int_var_stateful_prop_old[i])[_qp] = _int_var_old[i];
266 
268 }
const MaterialProperty< RankTwoTensor > & _fp_old
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop
const MaterialProperty< RankFourTensor > & _elasticity_tensor

◆ recoverOldState()

void FiniteStrainHyperElasticViscoPlastic::recoverOldState ( )
protectedvirtualinherited

This function restores the the old stateful properties after a successful solve.

Definition at line 293 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::postSolveQp().

294 {
295  // TODO: remove this nasty const_cast if you can figure out how
296  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
297  const_cast<MaterialProperty<Real> &>(*_int_var_stateful_prop_old[i])[_qp] = _int_var_old[i];
298 }
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old

◆ saveOldState()

void FiniteStrainHyperElasticViscoPlastic::saveOldState ( )
protectedvirtualinherited

This function saves the old stateful properties that is modified during sub stepping.

Definition at line 251 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::computeQpStress().

252 {
253  for (unsigned int i = 0; i < _num_int_var_uos; ++i)
255 }
unsigned int _num_int_var_uos
Number of internal variable user objects.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old

◆ solveFlowrate()

bool FiniteStrainHyperElasticViscoPlastic::solveFlowrate ( )
protectedvirtualinherited

Solve for flow rate and state.

Definition at line 317 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::solveQp().

318 {
319  Real resid0, rnorm;
320  unsigned int iter = 0;
321 
322 #ifdef DEBUG
323  std::vector<Real> rnormst(_maxiters + 1), flowratest(_maxiters + 1);
324 #endif
325 
327  return false;
328 
329  rnorm = computeNorm(_resid.get_values());
330  resid0 = rnorm;
331 
332 #ifdef DEBUG
333  rnormst[iter] = rnorm;
334  flowratest[iter] = computeNorm(_flow_rate.get_values());
335 #endif
336 
337  while (rnorm > _resid_abs_tol && rnorm > _resid_rel_tol * resid0 && iter < _maxiters)
338  {
340 
341  updateFlowRate();
342 
344 
346  return false;
347 
348  rnorm = computeNorm(_resid.get_values());
349  iter++;
350 
351 #ifdef DEBUG
352  rnormst[iter] = rnorm;
353  flowratest[iter] = computeNorm(_flow_rate.get_values());
354 #endif
355  }
356 
357  if (iter == _maxiters && rnorm > _resid_abs_tol && rnorm > _resid_rel_tol * resid0)
358  return false;
359 
360  return true;
361 }
unsigned int _maxiters
Maximum number of iterations.
Real _resid_abs_tol
Absolute tolerance for residual convergence check.
virtual bool computeFlowRateResidual()
Computes flow rate residual vector.
virtual void computeFlowRateJacobian()
Computes flow rate Jacobian matrix.
virtual Real computeNorm(const std::vector< Real > &)
Computes norm of residual vector.
virtual void computeElasticPlasticDeformGrad()
Computes elastic and plastic deformation gradients.
Real _resid_rel_tol
Relative tolerance for residual convergence check.

◆ solveQp()

bool FiniteStrainHyperElasticViscoPlastic::solveQp ( )
protectedvirtualinherited

Solve state.

Definition at line 271 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::computeQpStress().

272 {
274  if (!solveFlowrate())
275  return false;
277 
278  return true;
279 }
virtual bool solveFlowrate()
Solve for flow rate and state.
virtual void preSolveFlowrate()
Sets state for solve (Inside substepping)
virtual void postSolveFlowrate()
Update state for output (Inside substepping)

◆ updateFlowRate()

void FiniteStrainHyperElasticViscoPlastic::updateFlowRate ( )
protectedvirtualinherited

Update flow rate.

Definition at line 541 of file FiniteStrainHyperElasticViscoPlastic.C.

Referenced by FiniteStrainHyperElasticViscoPlastic::solveFlowrate().

542 {
543  _jac.lu_solve(_resid, _dflow_rate);
545 
546  for (unsigned int i = 0; i < _num_flow_rate_uos; ++i)
547  (*_flow_rate_prop[i])[_qp] = _flow_rate(i);
548 }
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
std::vector< MaterialProperty< Real > * > _flow_rate_prop

Member Data Documentation

◆ _base_name

const std::string ComputeStressBase::_base_name
protectedinherited

◆ _c

const VariableValue& HyperElasticPhaseFieldIsoDamage::_c
protected

Compupled damage variable.

Definition at line 62 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeDamageStress(), and computeQpJacobian().

◆ _ce

MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_ce
protectedinherited

◆ _d2Fdc2

MaterialProperty<Real>& HyperElasticPhaseFieldIsoDamage::_d2Fdc2
protected

Definition at line 75 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeDamageStress().

◆ _d2Fdcdstrain

MaterialProperty<RankTwoTensor>& HyperElasticPhaseFieldIsoDamage::_d2Fdcdstrain
protected

Definition at line 76 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeQpJacobian().

◆ _dce_dfe

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dce_dfe
protectedinherited

◆ _dee_dce

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dee_dce
protectedinherited

◆ _deformation_gradient

const MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_deformation_gradient
protectedinherited

◆ _deformation_gradient_old

const MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_deformation_gradient_old
protectedinherited

◆ _df_dstretch_inc

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_df_dstretch_inc
protectedinherited

◆ _dFdc

MaterialProperty<Real>& HyperElasticPhaseFieldIsoDamage::_dFdc
protected

Definition at line 74 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeDamageStress().

◆ _dfe_df

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dfe_df
protectedinherited

◆ _dfe_dfpinv

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dfe_dfpinv
protectedinherited

◆ _dfgrd_tmp

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_dfgrd_tmp
protectedinherited

◆ _dflow_rate

DenseVector<Real> FiniteStrainHyperElasticViscoPlastic::_dflow_rate
protectedinherited

◆ _dflowrate_dpk2

std::vector<RankTwoTensor> FiniteStrainHyperElasticViscoPlastic::_dflowrate_dpk2
protectedinherited

◆ _dflowrate_dstrength

DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dflowrate_dstrength
protectedinherited

◆ _dfpinv_dflowrate

std::vector<RankTwoTensor> FiniteStrainHyperElasticViscoPlastic::_dfpinv_dflowrate
protectedinherited

◆ _dG0_dee

RankTwoTensor HyperElasticPhaseFieldIsoDamage::_dG0_dee
protected

Definition at line 68 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeDamageStress(), and computeQpJacobian().

◆ _dintvar_dflowrate

DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dintvar_dflowrate
protectedinherited

◆ _dintvar_dflowrate_tmp

std::vector<DenseVector<Real> > FiniteStrainHyperElasticViscoPlastic::_dintvar_dflowrate_tmp
protectedinherited

◆ _dintvar_dintvar

DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dintvar_dintvar
protectedinherited

◆ _dintvar_dintvar_x

DenseVector<Real> FiniteStrainHyperElasticViscoPlastic::_dintvar_dintvar_x
protectedinherited

◆ _dintvar_dintvarrate

DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dintvar_dintvarrate
protectedinherited

◆ _dintvarrate_dflowrate

std::vector<DenseVector<Real> > FiniteStrainHyperElasticViscoPlastic::_dintvarrate_dflowrate
protectedinherited

◆ _dintvarrate_dintvar

DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dintvarrate_dintvar
protectedinherited

◆ _dpk2_dc

RankTwoTensor HyperElasticPhaseFieldIsoDamage::_dpk2_dc
protected

Definition at line 68 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeDamageStress(), and computeQpJacobian().

◆ _dpk2_dce

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dpk2_dce
protectedinherited

◆ _dpk2_dee

RankFourTensor HyperElasticPhaseFieldIsoDamage::_dpk2_dee
protected

◆ _dpk2_dfe

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dpk2_dfe
protectedinherited

◆ _dpk2_dflowrate

std::vector<RankTwoTensor> FiniteStrainHyperElasticViscoPlastic::_dpk2_dflowrate
protectedinherited

◆ _dpk2_dfpinv

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_dpk2_dfpinv
protectedinherited

◆ _dstrength_dintvar

DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_dstrength_dintvar
protectedinherited

◆ _dstress_dc

MaterialProperty<RankTwoTensor>& HyperElasticPhaseFieldIsoDamage::_dstress_dc
protected

Definition at line 66 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeQpJacobian().

◆ _dt_substep

Real FiniteStrainHyperElasticViscoPlastic::_dt_substep
protectedinherited

◆ _ee

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_ee
protectedinherited

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

◆ _etens

std::vector<RankTwoTensor> HyperElasticPhaseFieldIsoDamage::_etens
protected

Definition at line 70 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeDamageStress().

◆ _extra_stress

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_extra_stress
protectedinherited

Extra stress tensor.

Definition at line 47 of file ComputeStressBase.h.

Referenced by ComputeStressBase::computeQpProperties().

◆ _F

MaterialProperty<Real>& HyperElasticPhaseFieldIsoDamage::_F
protected

Elastic energy and derivatives, declared in this material.

Definition at line 73 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeDamageStress().

◆ _fe

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_fe
protectedinherited

◆ _fe_pk2

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_fe_pk2
protectedinherited

◆ _flow_dirn

std::vector<RankTwoTensor> FiniteStrainHyperElasticViscoPlastic::_flow_dirn
protectedinherited

◆ _flow_rate

DenseVector<Real> FiniteStrainHyperElasticViscoPlastic::_flow_rate
protectedinherited

◆ _flow_rate_prop

std::vector<MaterialProperty<Real> *> FiniteStrainHyperElasticViscoPlastic::_flow_rate_prop
protectedinherited

◆ _flow_rate_uo

std::vector<const HEVPFlowRateUOBase *> FiniteStrainHyperElasticViscoPlastic::_flow_rate_uo
protectedinherited

◆ _flow_rate_uo_names

std::vector<UserObjectName> FiniteStrainHyperElasticViscoPlastic::_flow_rate_uo_names
protectedinherited

◆ _fp

MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_fp
protectedinherited

◆ _fp_old

const MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_fp_old
protectedinherited

◆ _fp_tmp_inv

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_fp_tmp_inv
protectedinherited

◆ _fp_tmp_old_inv

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_fp_tmp_old_inv
protectedinherited

◆ _gc

const MaterialProperty<Real>& HyperElasticPhaseFieldIsoDamage::_gc
protected

Material property defining gc parameter, declared elsewhere.

Definition at line 54 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeDamageStress().

◆ _hist

MaterialProperty<Real>& HyperElasticPhaseFieldIsoDamage::_hist
protected

History variable that prevents crack healing, declared in this material.

Definition at line 79 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeDamageStress().

◆ _hist_old

const MaterialProperty<Real>& HyperElasticPhaseFieldIsoDamage::_hist_old
protected

Old value of history variable.

Definition at line 82 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeDamageStress().

◆ _initial_stress_fcn

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

initial stress components

Definition at line 50 of file ComputeStressBase.h.

◆ _int_var_old

std::vector<Real> FiniteStrainHyperElasticViscoPlastic::_int_var_old
protectedinherited

◆ _int_var_rate_prop

std::vector<MaterialProperty<Real> *> FiniteStrainHyperElasticViscoPlastic::_int_var_rate_prop
protectedinherited

◆ _int_var_rate_uo

std::vector<const HEVPInternalVarRateUOBase *> FiniteStrainHyperElasticViscoPlastic::_int_var_rate_uo
protectedinherited

◆ _int_var_rate_uo_names

std::vector<UserObjectName> FiniteStrainHyperElasticViscoPlastic::_int_var_rate_uo_names
protectedinherited

◆ _int_var_stateful_prop

std::vector<MaterialProperty<Real> *> FiniteStrainHyperElasticViscoPlastic::_int_var_stateful_prop
protectedinherited

◆ _int_var_stateful_prop_old

std::vector<const MaterialProperty<Real> *> FiniteStrainHyperElasticViscoPlastic::_int_var_stateful_prop_old
protectedinherited

◆ _int_var_uo

std::vector<const HEVPInternalVarUOBase *> FiniteStrainHyperElasticViscoPlastic::_int_var_uo
protectedinherited

◆ _int_var_uo_names

std::vector<UserObjectName> FiniteStrainHyperElasticViscoPlastic::_int_var_uo_names
protectedinherited

◆ _jac

DenseMatrix<Real> FiniteStrainHyperElasticViscoPlastic::_jac
protectedinherited

◆ _Jacobian_mult

MaterialProperty<RankFourTensor>& ComputeStressBase::_Jacobian_mult
protectedinherited

◆ _kdamage

Real HyperElasticPhaseFieldIsoDamage::_kdamage
protected

Small stiffness of completely damaged material point.

Definition at line 48 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeDamageStress(), and computeQpJacobian().

◆ _l

const MaterialProperty<Real>& HyperElasticPhaseFieldIsoDamage::_l
protected

Material property defining crack width, declared elsewhere.

Definition at line 52 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeDamageStress().

◆ _max_substep_iter

unsigned int FiniteStrainHyperElasticViscoPlastic::_max_substep_iter
protectedinherited

Maximum number of substep iterations.

Definition at line 165 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by FiniteStrainHyperElasticViscoPlastic::computeQpStress().

◆ _maxiters

unsigned int FiniteStrainHyperElasticViscoPlastic::_maxiters
protectedinherited

Maximum number of iterations.

Definition at line 163 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by FiniteStrainHyperElasticViscoPlastic::solveFlowrate().

◆ _mechanical_strain

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_mechanical_strain
protectedinherited

◆ _num_flow_rate_uos

unsigned int FiniteStrainHyperElasticViscoPlastic::_num_flow_rate_uos
protectedinherited

◆ _num_int_var_rate_uos

unsigned int FiniteStrainHyperElasticViscoPlastic::_num_int_var_rate_uos
protectedinherited

◆ _num_int_var_uos

unsigned int FiniteStrainHyperElasticViscoPlastic::_num_int_var_uos
protectedinherited

◆ _num_stiffness

bool HyperElasticPhaseFieldIsoDamage::_num_stiffness
protected

Flag to compute numerical stiffness.

Definition at line 46 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computePK2StressAndDerivative().

◆ _num_strength_uos

unsigned int FiniteStrainHyperElasticViscoPlastic::_num_strength_uos
protectedinherited

◆ _pert_val

Real HyperElasticPhaseFieldIsoDamage::_pert_val
protected

Perturbation value for strain components.

Definition at line 60 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeNumStiffness().

◆ _pk2

MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_pk2
protectedinherited

◆ _pk2_fet

RankTwoTensor FiniteStrainHyperElasticViscoPlastic::_pk2_fet
protectedinherited

◆ _pk2_prop_name

std::string FiniteStrainHyperElasticViscoPlastic::_pk2_prop_name
protectedinherited

◆ _pk2_tmp

RankTwoTensor HyperElasticPhaseFieldIsoDamage::_pk2_tmp
protected

◆ _resid

DenseVector<Real> FiniteStrainHyperElasticViscoPlastic::_resid
protectedinherited

◆ _resid_abs_tol

Real FiniteStrainHyperElasticViscoPlastic::_resid_abs_tol
protectedinherited

Absolute tolerance for residual convergence check.

Definition at line 159 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by FiniteStrainHyperElasticViscoPlastic::solveFlowrate().

◆ _resid_rel_tol

Real FiniteStrainHyperElasticViscoPlastic::_resid_rel_tol
protectedinherited

Relative tolerance for residual convergence check.

Definition at line 161 of file FiniteStrainHyperElasticViscoPlastic.h.

Referenced by FiniteStrainHyperElasticViscoPlastic::solveFlowrate().

◆ _rotation_increment

const MaterialProperty<RankTwoTensor>& FiniteStrainHyperElasticViscoPlastic::_rotation_increment
protectedinherited

◆ _save_state

bool HyperElasticPhaseFieldIsoDamage::_save_state
protected

Flag to save couple material properties.

Definition at line 64 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeDamageStress(), and computePK2StressAndDerivative().

◆ _strength_prop

std::vector<MaterialProperty<Real> *> FiniteStrainHyperElasticViscoPlastic::_strength_prop
protectedinherited

◆ _strength_uo

std::vector<const HEVPStrengthUOBase *> FiniteStrainHyperElasticViscoPlastic::_strength_uo
protectedinherited

◆ _strength_uo_names

std::vector<UserObjectName> FiniteStrainHyperElasticViscoPlastic::_strength_uo_names
protectedinherited

◆ _stress

MaterialProperty<RankTwoTensor>& ComputeStressBase::_stress
protectedinherited

Definition at line 41 of file ComputeStressBase.h.

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

◆ _tan_mod

RankFourTensor FiniteStrainHyperElasticViscoPlastic::_tan_mod
protectedinherited

◆ _use_current_hist

bool HyperElasticPhaseFieldIsoDamage::_use_current_hist
protected

Use current value of history variable.

Definition at line 50 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeDamageStress(), and computeQpJacobian().

◆ _zero_pert

Real HyperElasticPhaseFieldIsoDamage::_zero_pert
protected

Perturbation value for near zero or zero strain components.

Definition at line 58 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeNumStiffness().

◆ _zero_tol

Real HyperElasticPhaseFieldIsoDamage::_zero_tol
protected

Used in numerical stiffness calculation to check near zero values.

Definition at line 56 of file HyperElasticPhaseFieldIsoDamage.h.

Referenced by computeNumStiffness().


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