www.mooseframework.org
FiniteStrainHyperElasticViscoPlastic.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 FINITESTRAINHYPERELASTICVISCOPLASTIC_H
8 #define FINITESTRAINHYPERELASTICVISCOPLASTIC_H
9 
10 #include "ComputeStressBase.h"
11 #include "HEVPFlowRateUOBase.h"
12 #include "HEVPStrengthUOBase.h"
13 #include "HEVPInternalVarUOBase.h"
15 
17 
18 template <>
20 
32 {
33 public:
34  FiniteStrainHyperElasticViscoPlastic(const InputParameters & parameters);
35 
36 protected:
42  virtual void initUOVariables();
43 
45  void initNumUserObjects(const std::vector<UserObjectName> &, unsigned int &);
46 
48  template <typename T>
49  void
50  initProp(const std::vector<UserObjectName> &, unsigned int, std::vector<MaterialProperty<T> *> &);
51 
56  template <typename T>
57  void initPropOld(const std::vector<UserObjectName> &,
58  unsigned int,
59  std::vector<const MaterialProperty<T> *> &);
60 
62  template <typename T>
63  void initUserObjects(const std::vector<UserObjectName> &, unsigned int, std::vector<const T *> &);
64 
66  virtual void initJacobianVariables();
67 
69  virtual void initQpStatefulProperties();
70 
72  virtual void computeQpStress();
73 
75  virtual void computeQpJacobian();
76 
78  virtual void saveOldState();
79 
81  virtual void preSolveQp();
82 
84  virtual bool solveQp();
85 
87  virtual void postSolveQp();
88 
90  virtual void recoverOldState();
91 
93  virtual void preSolveFlowrate();
94 
96  virtual bool solveFlowrate();
97 
99  virtual void postSolveFlowrate();
100 
102  virtual bool computeFlowRateFunction();
103 
105  virtual bool computeFlowDirection();
106 
109 
111  virtual void computePK2StressAndDerivative();
112 
114  virtual void computeElasticStrain();
115 
117  virtual void computeDeeDce();
118 
120  virtual bool computeFlowRateResidual();
121 
123  virtual void computeFlowRateJacobian();
124 
126  virtual void computeElasticPlasticDeformGrad();
127 
129  virtual Real computeNorm(const std::vector<Real> &);
130 
132  virtual void updateFlowRate();
133 
135  virtual void computeDpk2Dfpinv();
136 
138  virtual bool computeIntVarRates();
139 
141  virtual bool computeIntVar();
142 
144  virtual bool computeStrength();
145 
147  virtual void computeIntVarRateDerivatives();
148 
150  virtual void computeIntVarDerivatives();
151 
154 
160  unsigned int _maxiters;
162  unsigned int _max_substep_iter;
163 
165  std::vector<UserObjectName> _flow_rate_uo_names;
167  std::vector<UserObjectName> _strength_uo_names;
169  std::vector<UserObjectName> _int_var_uo_names;
171  std::vector<UserObjectName> _int_var_rate_uo_names;
172 
174  unsigned int _num_flow_rate_uos;
176  unsigned int _num_strength_uos;
178  unsigned int _num_int_var_uos;
180  unsigned int _num_int_var_rate_uos;
181 
183  std::vector<const HEVPFlowRateUOBase *> _flow_rate_uo;
185  std::vector<const HEVPStrengthUOBase *> _strength_uo;
187  std::vector<const HEVPInternalVarUOBase *> _int_var_uo;
189  std::vector<const HEVPInternalVarRateUOBase *> _int_var_rate_uo;
190 
191  std::string _pk2_prop_name;
192  MaterialProperty<RankTwoTensor> & _pk2;
193  MaterialProperty<RankTwoTensor> & _fp;
194  const MaterialProperty<RankTwoTensor> & _fp_old;
195  MaterialProperty<RankTwoTensor> & _ce;
196 
197  const MaterialProperty<RankTwoTensor> & _deformation_gradient;
198  const MaterialProperty<RankTwoTensor> & _deformation_gradient_old;
199  const MaterialProperty<RankTwoTensor> & _rotation_increment;
200 
201  std::vector<MaterialProperty<Real> *> _flow_rate_prop;
202  std::vector<MaterialProperty<Real> *> _strength_prop;
203  std::vector<MaterialProperty<Real> *> _int_var_stateful_prop;
204  std::vector<const MaterialProperty<Real> *> _int_var_stateful_prop_old;
205  std::vector<MaterialProperty<Real> *> _int_var_rate_prop;
206  std::vector<Real> _int_var_old;
207 
208  RankTwoTensor _dfgrd_tmp;
209  RankTwoTensor _fp_tmp_inv, _fp_tmp_old_inv;
210  RankTwoTensor _fe, _ee;
211  RankTwoTensor _pk2_fet, _fe_pk2;
213  RankFourTensor _dee_dce, _dce_dfe, _dfe_df;
214  RankFourTensor _tan_mod, _df_dstretch_inc;
215 
216  std::vector<RankTwoTensor> _flow_dirn;
217  std::vector<RankTwoTensor> _dflowrate_dpk2;
218  std::vector<RankTwoTensor> _dpk2_dflowrate;
219  std::vector<RankTwoTensor> _dfpinv_dflowrate;
220 
221  DenseVector<Real> _dflow_rate;
222  DenseVector<Real> _flow_rate;
223  DenseVector<Real> _resid;
224 
226  std::vector<DenseVector<Real>> _dintvarrate_dflowrate;
227  std::vector<DenseVector<Real>> _dintvar_dflowrate_tmp;
228 
229  DenseMatrix<Real> _dintvarrate_dintvar;
230  DenseMatrix<Real> _dintvar_dintvarrate;
231  DenseMatrix<Real> _dintvar_dintvar;
232  DenseMatrix<Real> _dintvar_dflowrate;
233  DenseMatrix<Real> _dstrength_dintvar;
234  DenseMatrix<Real> _dflowrate_dstrength;
235  DenseVector<Real> _dintvar_dintvar_x;
236  DenseMatrix<Real> _jac;
237 
239 };
240 
241 #endif // FINITESTRAINHYPERELASTICVISCOPLASTIC_H
unsigned int _num_flow_rate_uos
Number of flow rate user objects.
unsigned int _maxiters
Maximum number of iterations.
virtual void recoverOldState()
This function restores the the old stateful properties after a successful solve.
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
This class solves the viscoplastic flow rate equations in the total form Involves 4 different types o...
ComputeStressBase is the base class for stress tensors.
std::vector< MaterialProperty< Real > * > _flow_rate_prop
void initUserObjects(const std::vector< UserObjectName > &, unsigned int, std::vector< const T * > &)
This function initializes 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 computeElasticStrain()
Computes elastic Lagrangian strain.
void initNumUserObjects(const std::vector< UserObjectName > &, unsigned int &)
This function calculates the number of each user object type.
unsigned int _max_substep_iter
Maximum number of substep iterations.
std::vector< UserObjectName > _strength_uo_names
Names of strength user objects.
virtual void initUOVariables()
This function initializes the properties, stateful properties and user objects The properties and sta...
virtual void computeDpk2Dfpinv()
Computes derivative of PK2 stress wrt inverse of plastic deformation gradient.
const MaterialProperty< RankTwoTensor > & _fp_old
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.
Real _resid_abs_tol
Absolute tolerance for residual convergence check.
virtual bool solveFlowrate()
Solve for flow rate and state.
virtual void computeQpJacobian()
This function computes the Jacobian.
virtual void postSolveQp()
Update state for output (Outside substepping)
std::vector< UserObjectName > _flow_rate_uo_names
Names of flow rate user objects.
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 initJacobianVariables()
This function initialize variables required for Jacobian calculation.
std::vector< const MaterialProperty< Real > * > _int_var_stateful_prop_old
std::vector< UserObjectName > _int_var_uo_names
Names of internal variable user objects.
virtual bool computeFlowRateResidual()
Computes flow rate residual vector.
virtual void computeFlowRateJacobian()
Computes flow rate Jacobian matrix.
std::vector< const HEVPInternalVarUOBase * > _int_var_uo
Internal variable user objects.
virtual void preSolveFlowrate()
Sets state for solve (Inside substepping)
virtual void computeIntVarDerivatives()
This function call user objects to compute dintvar/dintvar_rate and dintvar/dflowrate.
std::vector< const HEVPStrengthUOBase * > _strength_uo
Strength user objects.
InputParameters validParams< FiniteStrainHyperElasticViscoPlastic >()
virtual Real computeNorm(const std::vector< Real > &)
Computes norm of residual vector.
virtual void computeElasticPlasticDeformGrad()
Computes elastic and plastic deformation gradients.
virtual void computeQpStress()
This function computes the Cauchy stress.
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.
virtual bool computeFlowDirection()
Calls user objects to compute flow directions.
virtual void postSolveFlowrate()
Update state for output (Inside substepping)
const MaterialProperty< RankTwoTensor > & _rotation_increment
std::vector< MaterialProperty< Real > * > _int_var_stateful_prop
void initProp(const std::vector< UserObjectName > &, unsigned int, std::vector< MaterialProperty< T > * > &)
This function initializes properties for each user object.
std::vector< const HEVPFlowRateUOBase * > _flow_rate_uo
Flow rate user objects.
virtual void saveOldState()
This function saves the old stateful properties that is modified during sub stepping.
const MaterialProperty< RankTwoTensor > & _deformation_gradient
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.
std::vector< MaterialProperty< Real > * > _int_var_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...
Real _resid_rel_tol
Relative tolerance for residual convergence check.
FiniteStrainHyperElasticViscoPlastic(const InputParameters &parameters)
virtual bool computeStrength()
This function call user objects to compute strength.
virtual void initQpStatefulProperties()
Initializes state.
virtual bool computeIntVarRates()
This function call user objects to calculate rate of internal variables.
std::vector< MaterialProperty< Real > * > _strength_prop