www.mooseframework.org
FiniteStrainCrystalPlasticity.h
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #pragma once
11 
12 #include "ComputeStressBase.h"
13 
22 {
23 public:
25 
27 
28 protected:
32  virtual void computeQpStress();
33 
38  virtual void computeQpElasticityTensor();
39 
44  virtual void initQpStatefulProperties();
45 
51 
56  virtual void getSlipIncrements();
57 
58  // Override to modify slip system resistance evolution
62  virtual void update_slip_system_resistance();
63 
64  // Old function: Kept to avoid code break in computeQpStress
68  virtual void updateGss();
69 
73  virtual void getSlipSystems();
74 
79  virtual void assignSlipSysRes();
80 
84  virtual void readFileInitSlipSysRes();
85 
90  virtual void getInitSlipSysRes();
91 
95  virtual void readFileFlowRateParams();
96 
101  virtual void getFlowRateParams();
102 
106  virtual void readFileHardnessParams();
107 
111  virtual void getHardnessParams();
112 
116  virtual void initSlipSysProps();
117 
121  virtual void initAdditionalProps();
122 
126  virtual void preSolveQp();
127 
131  virtual void solveQp();
132 
136  virtual void postSolveQp();
137 
141  virtual void preSolveStatevar();
142 
146  virtual void solveStatevar();
147 
151  virtual void postSolveStatevar();
152 
156  virtual void preSolveStress();
157 
161  virtual void solveStress();
162 
166  virtual void postSolveStress();
167 
171  virtual void calcResidual(RankTwoTensor &);
172 
176  virtual void calcJacobian(RankFourTensor &);
177 
185 
190 
195 
200 
202 
206 
210  void calc_schmid_tensor();
211 
215  bool line_search_update(const Real rnorm_prev, const RankTwoTensor);
216 
221 
223  const unsigned int _nss;
224 
225  std::vector<Real> _gprops;
226  std::vector<Real> _hprops;
227  std::vector<Real> _flowprops;
228 
230  std::string _slip_sys_file_name;
231 
234 
240 
243 
252 
254  unsigned int _maxiter;
256  unsigned int _maxiterg;
257 
260 
263 
266 
268  unsigned int _num_slip_sys_props;
269 
271 
274 
277 
279  unsigned int _rndm_seed;
280 
282  unsigned int _max_substep_iter;
283 
286 
289 
292 
294  unsigned int _lsrch_max_iter;
295 
296  // Line search method
298 
310 
314  const std::string _elasticity_tensor_name;
318 
319  DenseVector<Real> _mo;
320  DenseVector<Real> _no;
321 
322  DenseVector<Real> _a0;
323  DenseVector<Real> _xm;
324 
329 
332  DenseVector<Real> _slip_incr, _tau, _dslipdtau;
333  std::vector<RankTwoTensor> _s0;
334 
337  std::vector<Real> _gss_tmp;
338  std::vector<Real> _gss_tmp_old;
339 
340  DenseVector<Real> _slip_sys_props;
341 
343 
345 
346  bool _err_tol;
347 
354 };
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 computed from MOOSE&#39;s strain calculators...
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
Elasticity tensor material property.
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.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
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...
const InputParameters & parameters() const
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.
FiniteStrainCrystalPlasticity uses the multiplicative decomposition of deformation gradient and solve...
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
const std::string _elasticity_tensor_name
Name of the elasticity tensor material property.
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.