www.mooseframework.org
FiniteStrainUObasedCP.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 
18 
34 {
35 public:
37 
39 
40 protected:
44  virtual void computeQpStress();
45 
50  virtual void initQpStatefulProperties();
51 
56  virtual void calcResidJacob();
57 
63 
67  virtual void preSolveQp();
68 
72  virtual void solveQp();
73 
77  virtual void postSolveQp();
78 
82  virtual void preSolveStatevar();
83 
87  virtual void solveStatevar();
88 
92  virtual void postSolveStatevar();
93 
97  virtual void preSolveStress();
98 
102  virtual void solveStress();
103 
107  virtual void postSolveStress();
108 
112  virtual void calcResidual();
113 
117  virtual void calcJacobian();
118 
122  virtual void getSlipRates();
123 
130  virtual void calcTangentModuli();
131 
135  virtual void elasticTangentModuli();
136 
140  virtual void elastoPlasticTangentModuli();
141 
145  bool lineSearchUpdate(const Real rnorm_prev, const RankTwoTensor);
146 
150  virtual bool isStateVariablesConverged();
151 
153  std::vector<const CrystalPlasticitySlipRate *> _uo_slip_rates;
154 
156  std::vector<const CrystalPlasticitySlipResistance *> _uo_slip_resistances;
157 
159  std::vector<const CrystalPlasticityStateVariable *> _uo_state_vars;
160 
162  std::vector<const CrystalPlasticityStateVarRateComponent *> _uo_state_var_evol_rate_comps;
163 
165  std::vector<MaterialProperty<std::vector<Real>> *> _mat_prop_slip_rates;
166 
168  std::vector<MaterialProperty<std::vector<Real>> *> _mat_prop_slip_resistances;
169 
171  std::vector<MaterialProperty<std::vector<Real>> *> _mat_prop_state_vars;
172 
174  std::vector<const MaterialProperty<std::vector<Real>> *> _mat_prop_state_vars_old;
175 
177  std::vector<MaterialProperty<std::vector<Real>> *> _mat_prop_state_var_evol_rate_comps;
178 
180  unsigned int _num_uo_slip_rates;
181 
184 
186  unsigned int _num_uo_state_vars;
187 
190 
192  std::vector<std::vector<Real>> _state_vars_old;
193 
195  std::vector<std::vector<Real>> _state_vars_old_stored;
196 
198  std::vector<std::vector<Real>> _state_vars_prev;
199 
208 
211 
214 
216  unsigned int _maxiter;
218  unsigned int _maxiterg;
219 
222 
224  unsigned int _max_substep_iter;
225 
228 
231 
234 
236  unsigned int _lsrch_max_iter;
237 
240 
248 
250  const std::string _elasticity_tensor_name;
255 
258 
261  DenseVector<Real> _tau;
262  std::vector<MaterialProperty<std::vector<RankTwoTensor>> *> _flow_direction;
263 
265  bool _err_tol;
266 
273 };
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 computed from MOOSE&#39;s strain calculators...
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...
std::vector< std::vector< Real > > _state_vars_old_stored
Local stored state variable (for sub-stepping)
RankTwoTensor _delta_dfgrd
Used for substepping; Uniformly divides the increment in deformation gradient.
static InputParameters validParams()
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.
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.
Real _substep_dt
Current substep size.
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.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
MaterialProperty< RankTwoTensor > & _update_rot
const std::string _elasticity_tensor_name
Name of the elasticity tensor material property.
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.
const InputParameters & parameters() const
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.
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor 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.