www.mooseframework.org
FiniteStrainUObasedCP.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 FINITESTRAINUOBASEDCP_H
8 #define FINITESTRAINUOBASEDCP_H
9 
10 #include "ComputeStressBase.h"
11 
16 
18 
19 template <>
20 InputParameters validParams<FiniteStrainUObasedCP>();
21 
37 {
38 public:
39  FiniteStrainUObasedCP(const InputParameters & parameters);
40 
41 protected:
45  virtual void computeQpStress();
46 
51  virtual void initQpStatefulProperties();
52 
57  virtual void calcResidJacob();
58 
64 
68  virtual void preSolveQp();
69 
73  virtual void solveQp();
74 
78  virtual void postSolveQp();
79 
83  virtual void preSolveStatevar();
84 
88  virtual void solveStatevar();
89 
93  virtual void postSolveStatevar();
94 
98  virtual void preSolveStress();
99 
103  virtual void solveStress();
104 
108  virtual void postSolveStress();
109 
113  virtual void calcResidual();
114 
118  virtual void calcJacobian();
119 
123  virtual void getSlipRates();
124 
131  virtual void calcTangentModuli();
132 
136  virtual void elasticTangentModuli();
137 
141  virtual void elastoPlasticTangentModuli();
142 
146  bool lineSearchUpdate(const Real rnorm_prev, const RankTwoTensor);
147 
151  virtual bool isStateVariablesConverged();
152 
154  std::vector<const CrystalPlasticitySlipRate *> _uo_slip_rates;
155 
157  std::vector<const CrystalPlasticitySlipResistance *> _uo_slip_resistances;
158 
160  std::vector<const CrystalPlasticityStateVariable *> _uo_state_vars;
161 
163  std::vector<const CrystalPlasticityStateVarRateComponent *> _uo_state_var_evol_rate_comps;
164 
166  std::vector<MaterialProperty<std::vector<Real>> *> _mat_prop_slip_rates;
167 
169  std::vector<MaterialProperty<std::vector<Real>> *> _mat_prop_slip_resistances;
170 
172  std::vector<MaterialProperty<std::vector<Real>> *> _mat_prop_state_vars;
173 
175  std::vector<const MaterialProperty<std::vector<Real>> *> _mat_prop_state_vars_old;
176 
178  std::vector<MaterialProperty<std::vector<Real>> *> _mat_prop_state_var_evol_rate_comps;
179 
181  unsigned int _num_uo_slip_rates;
182 
185 
187  unsigned int _num_uo_state_vars;
188 
191 
193  std::vector<std::vector<Real>> _state_vars_old;
194 
196  std::vector<std::vector<Real>> _state_vars_prev;
197 
199  Real _rtol;
201  Real _abs_tol;
203  Real _stol;
205  Real _zero_tol;
206 
208  RankTwoTensor _resid;
209 
211  RankFourTensor _jac;
212 
214  unsigned int _maxiter;
216  unsigned int _maxiterg;
217 
219  MooseEnum _tan_mod_type;
220 
222  unsigned int _max_substep_iter;
223 
226 
229 
232 
234  unsigned int _lsrch_max_iter;
235 
237  MooseEnum _lsrch_method;
238 
239  MaterialProperty<RankTwoTensor> & _fp;
240  const MaterialProperty<RankTwoTensor> & _fp_old;
241  MaterialProperty<RankTwoTensor> & _pk2;
242  const MaterialProperty<RankTwoTensor> & _pk2_old;
243  MaterialProperty<RankTwoTensor> & _lag_e;
244  MaterialProperty<RankTwoTensor> & _update_rot;
245  const MaterialProperty<RankTwoTensor> & _update_rot_old;
246 
247  const MaterialProperty<RankTwoTensor> & _deformation_gradient;
248  const MaterialProperty<RankTwoTensor> & _deformation_gradient_old;
249 
251  const MaterialProperty<RankTwoTensor> & _crysrot;
252 
253  RankTwoTensor _dfgrd_tmp;
254  RankTwoTensor _fe, _fp_old_inv, _fp_inv;
255  DenseVector<Real> _tau;
256  std::vector<MaterialProperty<std::vector<RankTwoTensor>> *> _flow_direction;
257 
259  bool _err_tol;
260 
262  RankTwoTensor _delta_dfgrd, _dfgrd_tmp_old;
265 };
266 
267 #endif // FINITESTRAINUOBASEDCP_H
virtual bool isStateVariablesConverged()
evaluates convergence of state variables.
const MaterialProperty< RankTwoTensor > & _pk2_old
std::vector< std::vector< Real > > _state_vars_old
Local state variable.
ComputeStressBase is the base class for stress tensors.
virtual void postSolveStatevar()
update internal variable after solve.
bool _err_tol
Flag to check whether convergence is achieved.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_slip_rates
Slip rates material property.
Real _stol
Internal variable update equation tolerance.
DenseVector< Real > _tau
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_var_evol_rate_comps
State variable evolution rate component material property.
std::vector< const CrystalPlasticitySlipResistance * > _uo_slip_resistances
User objects that define the slip resistance.
MooseEnum _tan_mod_type
Type of tangent moduli calculation.
unsigned int _lsrch_max_iter
Line search bisection method maximum iteration number.
Real _zero_tol
Residual tolerance when variable value is zero. Default 1e-12.
FiniteStrainUObasedCP(const InputParameters &parameters)
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_slip_resistances
Slip resistance material property.
std::vector< const CrystalPlasticityStateVarRateComponent * > _uo_state_var_evol_rate_comps
User objects that define the state variable evolution rate component.
virtual void initQpStatefulProperties()
initializes the stateful properties such as stress, plastic deformation gradient, slip system resista...
RankTwoTensor _delta_dfgrd
Used for substepping; Uniformly divides the increment in deformation gradient.
std::vector< std::vector< Real > > _state_vars_prev
Local old state variable.
unsigned int _maxiter
Maximum number of iterations for stress update.
std::vector< const CrystalPlasticityStateVariable * > _uo_state_vars
User objects that define the state variable.
const MaterialProperty< RankTwoTensor > & _update_rot_old
virtual void calcTangentModuli()
calculate the tangent moduli for preconditioner.
Real _dfgrd_scale_factor
Scales the substepping increment to obtain deformation gradient at a substep iteration.
virtual void calcResidJacob()
calls the residual and jacobian functions used in the stress update algorithm.
virtual void updateSlipSystemResistanceAndStateVariable()
updates the slip system resistances and state variables.
virtual void calcResidual()
calculate stress residual.
virtual void computeQpStress()
updates the stress at a quadrature point.
virtual void postSolveQp()
update stress and internal variable after solve.
MaterialProperty< RankTwoTensor > & _fp
virtual void getSlipRates()
updates the slip rates.
InputParameters validParams< FiniteStrainUObasedCP >()
Real _lsrch_tol
Line search bisection method tolerance.
virtual void postSolveStress()
update stress and plastic deformation gradient after solve.
virtual void preSolveStress()
set variables for stress solve.
unsigned int _num_uo_slip_rates
Number of slip rate user objects.
unsigned int _maxiterg
Maximum number of iterations for internal variable update.
const MaterialProperty< RankTwoTensor > & _deformation_gradient
std::vector< const CrystalPlasticitySlipRate * > _uo_slip_rates
User objects that define the slip rate.
RankFourTensor _jac
Jacobian tensor.
unsigned int _num_uo_state_vars
Number of state variable user objects.
unsigned int _num_uo_slip_resistances
Number of slip resistance user objects.
MaterialProperty< RankTwoTensor > & _lag_e
unsigned int _num_uo_state_var_evol_rate_comps
Number of state variable evolution rate component user objects.
Real _rtol
Stress residual equation relative tolerance.
FiniteStrainUObasedCP uses the multiplicative decomposition of deformation gradient and solves the PK...
bool lineSearchUpdate(const Real rnorm_prev, const RankTwoTensor)
performs the line search update
virtual void preSolveQp()
set variables for stress and internal variable solve.
virtual void elastoPlasticTangentModuli()
calculate the exact tangent moduli for preconditioner.
const MaterialProperty< RankTwoTensor > & _crysrot
Crystal rotation.
virtual void elasticTangentModuli()
calculate the elastic tangent moduli for preconditioner.
std::vector< const MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars_old
Old state variable material property.
virtual void preSolveStatevar()
set variables for internal variable solve.
MaterialProperty< RankTwoTensor > & _update_rot
virtual void solveQp()
solve stress and internal variables.
RankTwoTensor _resid
Residual tensor.
bool _use_line_search
Flag to activate line serach.
const MaterialProperty< RankTwoTensor > & _fp_old
unsigned int _max_substep_iter
Maximum number of substep iterations.
Real _abs_tol
Stress residual equation absolute tolerance.
virtual void solveStatevar()
solve internal variables.
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars
State variable material property.
std::vector< MaterialProperty< std::vector< RankTwoTensor > > * > _flow_direction
Real _min_lsrch_step
Minimum line search step size.
virtual void solveStress()
solves for stress, updates plastic deformation gradient.
MooseEnum _lsrch_method
Line search method.
MaterialProperty< RankTwoTensor > & _pk2
virtual void calcJacobian()
calculate jacobian.