www.mooseframework.org
FiniteStrainCrystalPlasticity.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 FINITESTRAINCRYSTALPLASTICITY_H
8 #define FINITESTRAINCRYSTALPLASTICITY_H
9 
10 #include "ComputeStressBase.h"
11 
20 
21 template <>
23 
25 {
26 public:
27  FiniteStrainCrystalPlasticity(const InputParameters & parameters);
28 
29 protected:
33  virtual void computeQpStress();
34 
39  virtual void computeQpElasticityTensor();
40 
45  virtual void initQpStatefulProperties();
46 
51  virtual void calc_resid_jacob(RankTwoTensor &, RankFourTensor &);
52 
57  virtual void getSlipIncrements();
58 
59  // Override to modify slip system resistance evolution
63  virtual void update_slip_system_resistance();
64 
65  // Old function: Kept to avoid code break in computeQpStress
69  virtual void updateGss();
70 
74  virtual void getSlipSystems();
75 
80  virtual void assignSlipSysRes();
81 
85  virtual void readFileInitSlipSysRes();
86 
91  virtual void getInitSlipSysRes();
92 
96  virtual void readFileFlowRateParams();
97 
102  virtual void getFlowRateParams();
103 
107  virtual void readFileHardnessParams();
108 
112  virtual void getHardnessParams();
113 
117  virtual void initSlipSysProps();
118 
122  virtual void initAdditionalProps();
123 
127  virtual void preSolveQp();
128 
132  virtual void solveQp();
133 
137  virtual void postSolveQp();
138 
142  virtual void preSolveStatevar();
143 
147  virtual void solveStatevar();
148 
152  virtual void postSolveStatevar();
153 
157  virtual void preSolveStress();
158 
162  virtual void solveStress();
163 
167  virtual void postSolveStress();
168 
172  virtual void calcResidual(RankTwoTensor &);
173 
177  virtual void calcJacobian(RankFourTensor &);
178 
185  virtual RankFourTensor calcTangentModuli();
186 
190  virtual RankFourTensor elasticTangentModuli();
191 
195  virtual RankFourTensor elastoPlasticTangentModuli();
196 
200  RankTwoTensor get_current_rotation(const RankTwoTensor & a);
201 
203 
206  RankTwoTensor getMatRot(const RankTwoTensor & a);
207 
211  void calc_schmid_tensor();
212 
216  bool line_search_update(const Real rnorm_prev, const RankTwoTensor);
217 
222 
224  const unsigned int _nss;
225 
226  std::vector<Real> _gprops;
227  std::vector<Real> _hprops;
228  std::vector<Real> _flowprops;
229 
231  std::string _slip_sys_file_name;
232 
235 
241 
244 
246  Real _rtol;
248  Real _abs_tol;
250  Real _gtol;
253 
255  unsigned int _maxiter;
257  unsigned int _maxiterg;
258 
261 
263  MooseEnum _tan_mod_type;
264 
266  MooseEnum _intvar_read_type;
267 
269  unsigned int _num_slip_sys_props;
270 
272 
275 
278 
280  unsigned int _rndm_seed;
281 
283  unsigned int _max_substep_iter;
284 
287 
290 
293 
295  unsigned int _lsrch_max_iter;
296 
297  // Line search method
298  MooseEnum _lsrch_method;
299 
300  MaterialProperty<RankTwoTensor> & _fp;
301  const MaterialProperty<RankTwoTensor> & _fp_old;
302  MaterialProperty<RankTwoTensor> & _pk2;
303  const MaterialProperty<RankTwoTensor> & _pk2_old;
304  MaterialProperty<RankTwoTensor> & _lag_e;
305  const MaterialProperty<RankTwoTensor> & _lag_e_old;
306  MaterialProperty<std::vector<Real>> & _gss;
307  const MaterialProperty<std::vector<Real>> & _gss_old;
308  MaterialProperty<Real> & _acc_slip;
309  const MaterialProperty<Real> & _acc_slip_old;
310  MaterialProperty<RankTwoTensor> & _update_rot;
311 
312  const MaterialProperty<RankTwoTensor> & _deformation_gradient;
313  const MaterialProperty<RankTwoTensor> & _deformation_gradient_old;
314  const MaterialProperty<RankFourTensor> & _elasticity_tensor;
315  const MaterialProperty<RankTwoTensor> & _crysrot;
316 
317  DenseVector<Real> _mo;
318  DenseVector<Real> _no;
319 
320  DenseVector<Real> _a0;
321  DenseVector<Real> _xm;
322 
323  Real _h0;
324  Real _tau_sat;
325  Real _tau_init;
326  Real _r;
327 
328  RankTwoTensor _dfgrd_tmp;
330  DenseVector<Real> _slip_incr, _tau, _dslipdtau;
331  std::vector<RankTwoTensor> _s0;
332 
333  RankTwoTensor _pk2_tmp, _pk2_tmp_old;
335  std::vector<Real> _gss_tmp;
336  std::vector<Real> _gss_tmp_old;
337 
338  DenseVector<Real> _slip_sys_props;
339 
340  DenseMatrix<Real> _dgss_dsliprate;
341 
343 
344  bool _err_tol;
345 
347  RankTwoTensor _delta_dfgrd, _dfgrd_tmp_old;
352 };
353 
354 #endif // FINITESTRAINCRYSTALPLASTICITY_H
const MaterialProperty< RankTwoTensor > & _pk2_old
void internalVariableUpdateNRiteration()
This function updates internal variables after each NewTon Raphson iteration (_fp_inv) ...
MaterialProperty< RankTwoTensor > & _lag_e
const MaterialProperty< RankTwoTensor > & _deformation_gradient
ComputeStressBase is the base class for stress tensors.
bool _input_rndm_scale_var
Input option for scaling variable to generate random stress when convergence fails.
virtual void getHardnessParams()
This function assign flow rate parameters from .i file - see test.
bool _use_line_search
Flag to activate line serach.
const MaterialProperty< std::vector< Real > > & _gss_old
const MaterialProperty< RankFourTensor > & _elasticity_tensor
virtual void calcJacobian(RankFourTensor &)
This function calculate jacobian.
virtual void computeQpStress()
This function updates the stress at a quadrature point.
MaterialProperty< RankTwoTensor > & _pk2
virtual void initQpStatefulProperties()
This function initializes the stateful properties such as stress, plastic deformation gradient...
virtual void solveStress()
This function solves for stress, updates plastic deformation gradient.
std::string _slip_sys_flow_prop_file_name
File should contain values of the flow rate equation parameters.
virtual void solveStatevar()
This function solves internal variables.
FiniteStrainCrystalPlasticity(const InputParameters &parameters)
virtual void update_slip_system_resistance()
This function updates the slip system resistances.
unsigned int _num_slip_sys_flowrate_props
Number of slip system flow rate parameters.
virtual void calcResidual(RankTwoTensor &)
This function calculate stress residual.
const MaterialProperty< RankTwoTensor > & _fp_old
RankTwoTensor getMatRot(const RankTwoTensor &a)
This function perform RU decomposition to obtain the rotation tensor.
virtual void readFileHardnessParams()
This function read hardness parameters from file.
virtual void updateGss()
This function updates the slip system resistances.
Real _rtol
Stress residual equation relative tolerance.
std::string _slip_sys_hard_prop_file_name
The hardening parameters in this class are read from .i file. The user can override to read from file...
unsigned int _maxiter
Maximum number of iterations for stress update.
virtual void assignSlipSysRes()
This function assign initial values of slip system resistances/internal variables read from getSlipSy...
Real _dfgrd_scale_factor
Scales the substepping increment to obtain deformation gradient at a substep iteration.
const unsigned int _nss
Number of slip system resistance.
virtual RankFourTensor elastoPlasticTangentModuli()
This function calculate the exact tangent moduli for preconditioner.
RankTwoTensor _delta_dfgrd
Flag to check whether convergence is achieved.
virtual RankFourTensor elasticTangentModuli()
This function calculate the elastic tangent moduli for preconditioner.
MaterialProperty< std::vector< Real > > & _gss
virtual RankFourTensor calcTangentModuli()
This function calculate the tangent moduli for preconditioner.
bool line_search_update(const Real rnorm_prev, const RankTwoTensor)
This function performs the line search update.
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
const MaterialProperty< RankTwoTensor > & _crysrot
void calc_schmid_tensor()
This function calculate the Schmid tensor.
MaterialProperty< RankTwoTensor > & _fp
virtual void postSolveStatevar()
This function update internal variable after solve.
unsigned int _lsrch_max_iter
Line search bisection method maximum iteration number.
virtual void getFlowRateParams()
This function assign flow rate parameters - see test.
Real _gtol
Internal variable update equation tolerance.
std::string _slip_sys_res_prop_file_name
File should contain initial values of the slip system resistances.
const MaterialProperty< RankTwoTensor > & _lag_e_old
virtual void postSolveQp()
This function update stress and internal variable after solve.
virtual void preSolveStress()
This function set variables for stress solve.
InputParameters validParams< FiniteStrainCrystalPlasticity >()
virtual void postSolveStress()
This function update stress and plastic deformation gradient after solve.
virtual void calc_resid_jacob(RankTwoTensor &, RankFourTensor &)
This function calls the residual and jacobian functions used in the stress update algorithm...
MooseEnum _tan_mod_type
Type of tangent moduli calculation.
Real _lsrch_tol
Line search bisection method tolerance.
virtual void computeQpElasticityTensor()
This function updates the elasticity tensor at a quadrature point.
virtual void getSlipSystems()
This function reads slip system from file - see test.
unsigned int _num_slip_sys_props
Number of slip system specific properties provided in the file containing slip system normals and dir...
virtual void readFileFlowRateParams()
This function read flow rate parameters from file - see test.
unsigned int _max_substep_iter
Maximum number of substep iterations.
bool _first_step_iter
Flags to reset variables and reinitialize variables.
virtual void initSlipSysProps()
This function initializes slip system resistances.
virtual void getInitSlipSysRes()
This function assign slip system resistances - see test.
virtual void initAdditionalProps()
This function initializes additional parameters.
virtual void preSolveStatevar()
This function set variables for internal variable solve.
MaterialProperty< RankTwoTensor > & _update_rot
Real _abs_tol
Stress residual equation absolute tolerance.
virtual void readFileInitSlipSysRes()
This function read slip system resistances from file - see test.
unsigned int _maxiterg
Maximum number of iterations for internal variable update.
virtual void getSlipIncrements()
This function updates the slip increments.
const MaterialProperty< Real > & _acc_slip_old
RankTwoTensor get_current_rotation(const RankTwoTensor &a)
This function perform RU decomposition to obtain the rotation tensor.
virtual void solveQp()
This function solves stress and internal variables.
MooseEnum _intvar_read_type
Read from options for initial values of internal variables.
Real _min_lsrch_step
Minimum line search step size.
virtual void preSolveQp()
This function set variables for stress and internal variable solve.
std::string _slip_sys_file_name
File should contain slip plane normal and direction. See test.
Real _slip_incr_tol
Slip increment tolerance.