www.mooseframework.org
RateDepSmearCrackModel.h
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 /****************************************************************/
7 #ifndef RATEDEPSMEARCRACKMODEL_H
8 #define RATEDEPSMEARCRACKMODEL_H
9 
10 #include "ConstitutiveModel.h"
11 #include "SymmElasticityTensor.h"
12 
13 #include "petscsys.h"
14 #include "petscblaslapack.h"
15 
16 #if PETSC_VERSION_LESS_THAN(3, 5, 0)
17 extern "C" void FORTRAN_CALL(dgetri)(...); // matrix inversion routine from LAPACK
18 #endif
19 
21 
22 template <>
23 InputParameters validParams<RateDepSmearCrackModel>();
24 
31 {
32 public:
33  RateDepSmearCrackModel(const InputParameters & parameters);
34 
35  virtual ~RateDepSmearCrackModel();
36 
37 protected:
38  virtual void computeStress(const Elem & current_elem,
39  const SymmElasticityTensor & elasticity_tensor,
40  const SymmTensor & stress_old,
41  SymmTensor & strain_increment,
42  SymmTensor & stress_new);
43 
44  virtual void initQpStatefulProperties();
45 
46  virtual void initVariables();
47 
55  virtual void solve();
59  virtual void postSolveVariables();
64  virtual void postSolveStress();
65 
70  virtual void updateVariables();
71 
75  bool getConvergeVar();
76 
81  virtual void calcResidual();
82 
87  virtual void calcStateIncr();
88 
92  virtual void calcJacobian();
93 
94  int matrixInversion(std::vector<Real> & A, int n) const;
95 
97  unsigned int _nstate;
98  Real _exponent;
99  unsigned int _maxiter;
100  Real _tol;
101  Real _zero_tol;
105 
106  MaterialProperty<std::vector<Real>> & _intvar;
107  const MaterialProperty<std::vector<Real>> & _intvar_old;
108 
109  MaterialProperty<SymmTensor> & _stress_undamaged;
110  const MaterialProperty<SymmTensor> & _stress_undamaged_old;
111 
112  std::vector<Real> _intvar_incr;
113  std::vector<Real> _intvar_tmp, _intvar_old_tmp;
114  std::vector<Real> _resid;
115  std::vector<Real> _jac;
116  std::vector<Real> _dvar;
117 
121  bool _nconv;
122  bool _err_tol;
123 
124 private:
125 };
126 
127 #endif // RATEDEPSMEARCRACKMODEL
This class defines a basic set of capabilities any elasticity tensor should have. ...
void FORTRAN_CALL() dgetri(...)
virtual void initQpStatefulProperties()
Real _tol
Maximum number of Newton Raphson iteration.
Real _rndm_scale_var
Flag to specify scaling parameter to generate random stress.
InputParameters validParams< RateDepSmearCrackModel >()
virtual void computeStress(const Elem &current_elem, const SymmElasticityTensor &elasticity_tensor, const SymmTensor &stress_old, SymmTensor &strain_increment, SymmTensor &stress_new)
std::vector< Real > _intvar_tmp
std::vector< Real > _intvar_incr
RateDepSmearCrackModel(const InputParameters &parameters)
int matrixInversion(std::vector< Real > &A, int n) const
Real _exponent
Number of state variables.
virtual void solve()
This function solves the state variables.
virtual void updateVariables()
This function updates variables during solve a(i+1) = a(i) + da(i+1)
const MaterialProperty< SymmTensor > & _stress_undamaged_old
RateDepSmearCrackModel is the base class for rate dependent continuum damage model.
bool _err_tol
Convergence flag.
const MaterialProperty< std::vector< Real > > & _intvar_old
MaterialProperty< SymmTensor > & _stress_undamaged
virtual void calcResidual()
This function calculates the residual as r = v - v_old - dv.
virtual void postSolveStress()
This function updates the stress after solve.
unsigned int _nstate
reference damage rate
virtual void calcStateIncr()
This function calculated thw increment of the state variables (dv) used to form the residual...
bool _input_rndm_scale_var
Allowable relative increment size of state variables (dv)
virtual void postSolveVariables()
This function updates the internal variables after solve.
MaterialProperty< std::vector< Real > > & _intvar
Variable value.
virtual void calcJacobian()
This function calculates the Jacobian.
std::vector< Real > _intvar_old_tmp
Real _intvar_incr_tol
Tolerance for zero.
SymmElasticityTensor _elasticity
Real _zero_tol
Relative tolerance factor for convergence of the Newton Raphson solve.
bool getConvergeVar()
This function returns true if convergence is not achieved.