www.mooseframework.org
RateDepSmearIsoCrackModel.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
8 
9 template <>
10 InputParameters
12 {
13 
14  InputParameters params = validParams<RateDepSmearCrackModel>();
15 
16  params.addRequiredParam<Real>("critical_energy", "Critical Energy");
17  params.addParam<Real>("k_fail", 1e-6, "Post failure stiffness");
18  params.addParam<Real>("upper_limit_damage",
19  5.0,
20  "Upper limit of damage beyond which constitutive check is not performed");
21 
22  return params;
23 }
24 
25 RateDepSmearIsoCrackModel::RateDepSmearIsoCrackModel(const InputParameters & parameters)
26  : RateDepSmearCrackModel(parameters),
27  _crit_energy(getParam<Real>("critical_energy")),
28  _kfail(getParam<Real>("k_fail")),
29  _upper_lim_damage(getParam<Real>("upper_limit_damage")),
30  _energy(declareProperty<Real>("energy")),
31  _energy_old(getMaterialPropertyOld<Real>("energy"))
32 {
33 
34  if (_nstate != 2)
35  mooseError(" RateDpeSmearIsoCrackModel Error: Requires 2 state variables - nstate = 2");
36 }
37 
38 void
40 {
42 
43  _intvar[_qp][1] = const_cast<MaterialProperty<std::vector<Real>> &>(_intvar_old)[_qp][1] =
45 
46  _energy[_qp] = 0.0;
47 }
48 
49 void
51 {
52 
54 
56 
57  for (unsigned int i = 0; i < LIBMESH_DIM; i++)
58  {
59  _s0_diag_pos(i, i) = (std::abs(_s0_diag(i, 0)) + _s0_diag(i, 0)) / 2.0;
60  _s0_diag_neg(i, i) = (std::abs(_s0_diag(i, 0)) - _s0_diag(i, 0)) / 2.0;
61  }
62 
64 
65  for (unsigned int i = 0; i < LIBMESH_DIM; i++)
66  {
67  _dstrain_diag_pos(i, i) = (std::abs(_dstrain_diag(i, 0)) + _dstrain_diag(i, 0)) / 2.0;
68  _dstrain_diag_neg(i, i) = (std::abs(_dstrain_diag(i, 0)) - _dstrain_diag(i, 0)) / 2.0;
69  }
70 }
71 
72 Real
74 {
75 
76  Real den = 0.0;
77 
78  for (unsigned int i = 0; i < LIBMESH_DIM; i++)
79  den += _s0_diag_pos(i, i) * _dstrain_diag_pos(i, i);
80 
81  _energy[_qp] = _energy_old[_qp] + den;
82 
83  Real e = _energy[_qp];
84 
85  Real fac = e - _intvar_tmp[1];
86  Real val = 0.0;
87 
88  _ddamagerate_drs = 0.0;
89 
90  if (fac > 0.0)
91  {
92  val = std::pow(fac, _exponent);
93  val *= _ref_damage_rate * _dt;
94 
95  Real dval = -_exponent * std::pow(fac, _exponent - 1.0);
96 
97  _ddamagerate_drs = dval * _ref_damage_rate * _dt;
98  }
99 
100  return val;
101 }
102 
103 void
105 {
106 
107  Real val = damageRate();
108 
109  if (std::abs(_intvar_old_tmp[0]) < _zero_tol && std::abs(val) > _intvar_incr_tol)
110  {
111  _err_tol = true;
112  return;
113  }
114 
115  if (std::abs(_intvar_old_tmp[0]) > _zero_tol &&
116  std::abs(_intvar_old_tmp[0]) < _upper_lim_damage &&
117  std::abs(val) > _intvar_incr_tol * std::abs(_intvar_old_tmp[0]))
118  {
119  _err_tol = true;
120  return;
121  }
122 
123  for (unsigned int i = 0; i < _nstate; i++)
124  _intvar_incr[i] = val;
125 }
126 
127 void
129 {
130 
131  _jac[0] = 1.0;
132  _jac[1] = -_ddamagerate_drs;
133  _jac[2] = 0.0;
134  _jac[3] = 1 - _ddamagerate_drs;
135 }
136 
137 void
139 {
140  ColumnMajorMatrix stress_cmm = _s0_diag_pos * (std::exp(-_intvar_tmp[0]) + _kfail) - _s0_diag_neg;
141  _stress_new = _s0_evec * stress_cmm * _s0_evec.transpose();
142 }
143 
MaterialProperty< Real > & _energy
virtual void postSolveStress()
This function updates the stress after solve.
virtual void initQpStatefulProperties()
virtual void calcStateIncr()
This function calculated thw increment of the state variables (dv) used to form the residual...
virtual Real damageRate()
This function calculates rate of damage based on energy.
InputParameters validParams< RateDepSmearCrackModel >()
std::vector< Real > _intvar_tmp
std::vector< Real > _intvar_incr
ColumnMajorMatrix columnMajorMatrix() const
Definition: SymmTensor.h:423
InputParameters validParams< RateDepSmearIsoCrackModel >()
Real _upper_lim_damage
Used to avoid non-positive definiteness.
Real _exponent
Number of state variables.
Real _kfail
Required parameter.
RateDepSmearCrackModel is the base class for rate dependent continuum damage model.
bool _err_tol
Convergence flag.
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
const MaterialProperty< std::vector< Real > > & _intvar_old
virtual void calcJacobian()
This function calculates the Jacobian.
unsigned int _nstate
reference damage rate
MaterialProperty< std::vector< Real > > & _intvar
Variable value.
RateDepSmearIsoCrackModel(const InputParameters &parameters)
std::vector< Real > _intvar_old_tmp
Real _intvar_incr_tol
Tolerance for zero.
Real _zero_tol
Relative tolerance factor for convergence of the Newton Raphson solve.
const MaterialProperty< Real > & _energy_old