www.mooseframework.org
ComputeMultiPlasticityStress.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 COMPUTEMULTIPLASTICITYSTRESS_H
8 #define COMPUTEMULTIPLASTICITYSTRESS_H
9 
10 #include "ComputeStressBase.h"
12 
14 
15 template <>
17 
29 {
30 public:
31  ComputeMultiPlasticityStress(const InputParameters & parameters);
32 
33 protected:
34  virtual void computeQpStress();
35  virtual void initQpStatefulProperties();
36 
38  unsigned int _max_iter;
39 
42 
45 
48 
51  {
56 
58  Real _epp_tol;
59 
61  std::vector<Real> _dummy_pm;
62 
67  std::vector<Real> _cumulative_pm;
68 
69  /*
70  * Scheme by which constraints are deactivated.
71  * This enum is defined here for computational
72  * efficiency. If you add to this list you must
73  * also add to the MooseEnum in the .C file
74  */
76  {
85 
88 
90  RealVectorValue _n_input;
91 
93  RealTensorValue _rot;
94 
97 
99  MaterialProperty<RankTwoTensor> & _plastic_strain;
100 
102  const MaterialProperty<RankTwoTensor> & _plastic_strain_old;
103 
105  MaterialProperty<std::vector<Real>> & _intnl;
106 
108  const MaterialProperty<std::vector<Real>> & _intnl_old;
109 
111  MaterialProperty<std::vector<Real>> & _yf;
112 
114  MaterialProperty<Real> & _iter;
115 
117  MaterialProperty<Real> & _linesearch_needed;
118 
120  MaterialProperty<Real> & _ld_encountered;
121 
123  MaterialProperty<Real> & _constraints_added;
124 
126  MaterialProperty<RealVectorValue> & _n;
127 
129  const MaterialProperty<RealVectorValue> & _n_old;
130 
132  const MaterialProperty<RankTwoTensor> & _strain_increment;
133 
135  const MaterialProperty<RankTwoTensor> & _total_strain_old;
136 
138  const MaterialProperty<RankTwoTensor> & _rotation_increment;
139 
141  const MaterialProperty<RankTwoTensor> & _stress_old;
142 
144  const MaterialProperty<RankTwoTensor> & _elastic_strain_old;
145 
147  bool _cosserat;
148 
150  const MaterialProperty<RankTwoTensor> * _curvature;
151 
153  const MaterialProperty<RankFourTensor> * _elastic_flexural_rigidity_tensor;
154 
156  MaterialProperty<RankTwoTensor> * _couple_stress;
157 
159  const MaterialProperty<RankTwoTensor> * _couple_stress_old;
160 
162  MaterialProperty<RankFourTensor> * _Jacobian_mult_couple;
163 
165  RankFourTensor _my_elasticity_tensor;
166 
168  RankTwoTensor _my_strain_increment;
169 
172 
174  RankTwoTensor _my_curvature;
175 
179  virtual bool reinstateLinearDependentConstraints(std::vector<bool> & deactivated_due_to_ld);
180 
184  virtual unsigned int numberActive(const std::vector<bool> & active);
185 
197  virtual Real residual2(const std::vector<Real> & pm,
198  const std::vector<Real> & f,
199  const RankTwoTensor & epp,
200  const std::vector<Real> & ic,
201  const std::vector<bool> & active,
202  const std::vector<bool> & deactivated_due_to_ld);
203 
240  virtual bool returnMap(const RankTwoTensor & stress_old,
241  RankTwoTensor & stress,
242  const std::vector<Real> & intnl_old,
243  std::vector<Real> & intnl,
244  const RankTwoTensor & plastic_strain_old,
245  RankTwoTensor & plastic_strain,
246  const RankFourTensor & E_ijkl,
247  const RankTwoTensor & strain_increment,
248  std::vector<Real> & f,
249  unsigned int & iter,
250  bool can_revert_to_dumb,
251  bool & linesearch_needed,
252  bool & ld_encountered,
253  bool & constraints_added,
254  bool final_step,
255  RankFourTensor & consistent_tangent_operator,
256  std::vector<Real> & cumulative_pm);
257 
290  virtual bool lineSearch(Real & nr_res2,
291  RankTwoTensor & stress,
292  const std::vector<Real> & intnl_old,
293  std::vector<Real> & intnl,
294  std::vector<Real> & pm,
295  const RankFourTensor & E_inv,
296  RankTwoTensor & delta_dp,
297  const RankTwoTensor & dstress,
298  const std::vector<Real> & dpm,
299  const std::vector<Real> & dintnl,
300  std::vector<Real> & f,
301  RankTwoTensor & epp,
302  std::vector<Real> & ic,
303  const std::vector<bool> & active,
304  const std::vector<bool> & deactivated_due_to_ld,
305  bool & linesearch_needed);
306 
333  virtual bool singleStep(Real & nr_res2,
334  RankTwoTensor & stress,
335  const std::vector<Real> & intnl_old,
336  std::vector<Real> & intnl,
337  std::vector<Real> & pm,
338  RankTwoTensor & delta_dp,
339  const RankFourTensor & E_inv,
340  std::vector<Real> & f,
341  RankTwoTensor & epp,
342  std::vector<Real> & ic,
343  std::vector<bool> & active,
344  DeactivationSchemeEnum deactivation_scheme,
345  bool & linesearch_needed,
346  bool & ld_encountered);
347 
355  virtual bool checkAdmissible(const RankTwoTensor & stress,
356  const std::vector<Real> & intnl,
357  std::vector<Real> & all_f);
358 
368  void buildDumbOrder(const RankTwoTensor & stress,
369  const std::vector<Real> & intnl,
370  std::vector<unsigned int> & dumb_order);
371 
383  virtual void incrementDumb(int & dumb_iteration,
384  const std::vector<unsigned int> & dumb_order,
385  std::vector<bool> & act);
386 
404  virtual bool checkKuhnTucker(const std::vector<Real> & f,
405  const std::vector<Real> & pm,
406  const std::vector<bool> & active);
407 
425  virtual void applyKuhnTucker(const std::vector<Real> & f,
426  const std::vector<Real> & pm,
427  std::vector<bool> & active);
428 
429  // gets called before any return-map
430  virtual void preReturnMap();
431 
432  // gets called after return-map
433  virtual void postReturnMap();
434 
437  {
440  };
441 
480  virtual bool quickStep(const RankTwoTensor & stress_old,
481  RankTwoTensor & stress,
482  const std::vector<Real> & intnl_old,
483  std::vector<Real> & intnl,
484  std::vector<Real> & pm,
485  std::vector<Real> & cumulative_pm,
486  const RankTwoTensor & plastic_strain_old,
487  RankTwoTensor & plastic_strain,
488  const RankFourTensor & E_ijkl,
489  const RankTwoTensor & strain_increment,
490  std::vector<Real> & yf,
491  unsigned int & iterations,
492  RankFourTensor & consistent_tangent_operator,
493  const quickStep_called_from_t called_from,
494  bool final_step);
495 
522  virtual bool plasticStep(const RankTwoTensor & stress_old,
523  RankTwoTensor & stress,
524  const std::vector<Real> & intnl_old,
525  std::vector<Real> & intnl,
526  const RankTwoTensor & plastic_strain_old,
527  RankTwoTensor & plastic_strain,
528  const RankFourTensor & E_ijkl,
529  const RankTwoTensor & strain_increment,
530  std::vector<Real> & yf,
531  unsigned int & iterations,
532  bool & linesearch_needed,
533  bool & ld_encountered,
534  bool & constraints_added,
535  RankFourTensor & consistent_tangent_operator);
536 
537  // bool checkAndModifyConstraints(bool nr_exit_condition, const RankTwoTensor & stress, const
538  // std::vector<Real> & intnl, const std::vector<Real> & pm, const std::vector<bool> &
539  // initial_act, bool can_revert_to_dumb, const RankTwoTensor & initial_stress, const
540  // std::vector<Real> & intnl_old, const std::vector<Real> & f, DeactivationSchemeEnum
541  // deact_scheme, std::vector<bool> & act, int & dumb_iteration, std::vector<unsigned int>
542  // dumb_order, bool & die);
543 
544  bool canChangeScheme(DeactivationSchemeEnum current_deactivation_scheme, bool can_revert_to_dumb);
545 
546  bool canIncrementDumb(int dumb_iteration);
547 
548  void changeScheme(const std::vector<bool> & initial_act,
549  bool can_revert_to_dumb,
550  const RankTwoTensor & initial_stress,
551  const std::vector<Real> & intnl_old,
552  DeactivationSchemeEnum & current_deactivation_scheme,
553  std::vector<bool> & act,
554  int & dumb_iteration,
555  std::vector<unsigned int> & dumb_order);
556 
557  bool canAddConstraints(const std::vector<bool> & act, const std::vector<Real> & all_f);
558 
559  unsigned int activeCombinationNumber(const std::vector<bool> & act);
560 
578  RankFourTensor consistentTangentOperator(const RankTwoTensor & stress,
579  const std::vector<Real> & intnl,
580  const RankFourTensor & E_ijkl,
581  const std::vector<Real> & pm_this_step,
582  const std::vector<Real> & cumulative_pm);
583 
584 private:
585  // InitialStress Deprecation: remove _step_one parameter
591  bool & _step_one;
592 
593  RankTwoTensor rot(const RankTwoTensor & tens);
594 };
595 
596 #endif // COMPUTEMULTIPLASTICITYSTRESS_H
bool _perform_finite_strain_rotations
whether to perform the rotations necessary in finite-strain simulations
virtual bool reinstateLinearDependentConstraints(std::vector< bool > &deactivated_due_to_ld)
makes all deactivated_due_to_ld false, and if >0 of them were initially true, returns true ...
ComputeStressBase is the base class for stress tensors.
bool canAddConstraints(const std::vector< bool > &act, const std::vector< Real > &all_f)
quickStep_called_from_t
The functions from which quickStep can be called.
virtual bool singleStep(Real &nr_res2, RankTwoTensor &stress, const std::vector< Real > &intnl_old, std::vector< Real > &intnl, std::vector< Real > &pm, RankTwoTensor &delta_dp, const RankFourTensor &E_inv, std::vector< Real > &f, RankTwoTensor &epp, std::vector< Real > &ic, std::vector< bool > &active, DeactivationSchemeEnum deactivation_scheme, bool &linesearch_needed, bool &ld_encountered)
Performs a single Newton-Raphson + linesearch step Constraints are deactivated and the step is re-don...
const MaterialProperty< RankTwoTensor > & _rotation_increment
Rotation increment (coming from ComputeIncrementalSmallStrain, for example)
const MaterialProperty< RealVectorValue > & _n_old
old value of transverse direction
MaterialProperty< Real > & _constraints_added
Whether constraints were added in during the latest Newton-Raphson process (1 if true, 0 otherwise)
MaterialProperty< std::vector< Real > > & _intnl
internal parameters
bool canChangeScheme(DeactivationSchemeEnum current_deactivation_scheme, bool can_revert_to_dumb)
RankFourTensor _my_flexural_rigidity_tensor
Flexual rigidity tensor that can be rotated by this class (ie, its not const)
ComputeMultiPlasticityStress performs the return-map algorithm and associated stress updates for plas...
const MaterialProperty< RankTwoTensor > & _strain_increment
strain increment (coming from ComputeIncrementalSmallStrain, for example)
TangentOperatorEnum
The type of tangent operator to return. tangent operator = d(stress_rate)/d(strain_rate).
const MaterialProperty< RankTwoTensor > & _stress_old
Old value of stress.
unsigned int _max_iter
Maximum number of Newton-Raphson iterations allowed.
RankTwoTensor _my_strain_increment
Strain increment that can be rotated by this class, and split into multiple increments (ie...
unsigned int activeCombinationNumber(const std::vector< bool > &act)
InputParameters validParams< ComputeMultiPlasticityStress >()
virtual void applyKuhnTucker(const std::vector< Real > &f, const std::vector< Real > &pm, std::vector< bool > &active)
Checks Kuhn-Tucker conditions, and alters "active" if appropriate.
bool _ignore_failures
Even if the returnMap fails, return the best values found for stress and internal parameters...
RankFourTensor _my_elasticity_tensor
Elasticity tensor that can be rotated by this class (ie, its not const)
MaterialProperty< Real > & _iter
Number of Newton-Raphson iterations used in the return-map.
const MaterialProperty< RankTwoTensor > & _elastic_strain_old
Old value of elastic strain.
bool _n_supplied
User supplied the transverse direction vector.
MaterialProperty< Real > & _ld_encountered
Whether linear-dependence was encountered in the latest Newton-Raphson process (1 if true...
void changeScheme(const std::vector< bool > &initial_act, bool can_revert_to_dumb, const RankTwoTensor &initial_stress, const std::vector< Real > &intnl_old, DeactivationSchemeEnum &current_deactivation_scheme, std::vector< bool > &act, int &dumb_iteration, std::vector< unsigned int > &dumb_order)
virtual bool checkKuhnTucker(const std::vector< Real > &f, const std::vector< Real > &pm, const std::vector< bool > &active)
Checks Kuhn-Tucker conditions, and alters "active" if appropriate.
virtual bool lineSearch(Real &nr_res2, RankTwoTensor &stress, const std::vector< Real > &intnl_old, std::vector< Real > &intnl, std::vector< Real > &pm, const RankFourTensor &E_inv, RankTwoTensor &delta_dp, const RankTwoTensor &dstress, const std::vector< Real > &dpm, const std::vector< Real > &dintnl, std::vector< Real > &f, RankTwoTensor &epp, std::vector< Real > &ic, const std::vector< bool > &active, const std::vector< bool > &deactivated_due_to_ld, bool &linesearch_needed)
Performs a line search.
RankTwoTensor rot(const RankTwoTensor &tens)
MaterialProperty< std::vector< Real > > & _yf
yield functions
MaterialProperty< RealVectorValue > & _n
current value of transverse direction
std::vector< Real > _dummy_pm
dummy "consistency parameters" (plastic multipliers) used in quickStep when called from computeQpStre...
bool & _step_one
True if this is the first timestep (timestep < 2).
virtual void incrementDumb(int &dumb_iteration, const std::vector< unsigned int > &dumb_order, std::vector< bool > &act)
Increments "dumb_iteration" by 1, and sets "act" appropriately (act[alpha] = true iff alpha_th bit of...
MaterialProperty< Real > & _linesearch_needed
Whether a line-search was needed in the latest Newton-Raphson process (1 if true, 0 otherwise) ...
const MaterialProperty< RankTwoTensor > * _couple_stress_old
the old value of Cosserat couple-stress
RealVectorValue _n_input
the supplied transverse direction vector
Real _max_stepsize_for_dumb
"dumb" deactivation will only be used if the stepsize falls below this quantity
const MaterialProperty< std::vector< Real > > & _intnl_old
old values of internal parameters
enum ComputeMultiPlasticityStress::TangentOperatorEnum _tangent_operator_type
virtual bool checkAdmissible(const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< Real > &all_f)
Checks whether the yield functions are in the admissible region.
virtual bool returnMap(const RankTwoTensor &stress_old, RankTwoTensor &stress, const std::vector< Real > &intnl_old, std::vector< Real > &intnl, const RankTwoTensor &plastic_strain_old, RankTwoTensor &plastic_strain, const RankFourTensor &E_ijkl, const RankTwoTensor &strain_increment, std::vector< Real > &f, unsigned int &iter, bool can_revert_to_dumb, bool &linesearch_needed, bool &ld_encountered, bool &constraints_added, bool final_step, RankFourTensor &consistent_tangent_operator, std::vector< Real > &cumulative_pm)
Implements the return map.
virtual Real residual2(const std::vector< Real > &pm, const std::vector< Real > &f, const RankTwoTensor &epp, const std::vector< Real > &ic, const std::vector< bool > &active, const std::vector< bool > &deactivated_due_to_ld)
The residual-squared.
ComputeMultiPlasticityStress(const InputParameters &parameters)
MaterialProperty< RankTwoTensor > & _plastic_strain
plastic strain
Real _epp_tol
Tolerance on the plastic strain increment ("direction") constraint.
RankFourTensor consistentTangentOperator(const RankTwoTensor &stress, const std::vector< Real > &intnl, const RankFourTensor &E_ijkl, const std::vector< Real > &pm_this_step, const std::vector< Real > &cumulative_pm)
Computes the consistent tangent operator (another name for the jacobian = d(stress_rate)/d(strain_rat...
Real _min_stepsize
Minimum fraction of applied strain that may be applied during adaptive stepsizing.
const MaterialProperty< RankTwoTensor > * _curvature
The Cosserat curvature strain.
MultiPlasticityDebugger computes various finite-difference things to help developers remove bugs in t...
RankTwoTensor _my_curvature
Curvature that can be rotated by this class, and split into multiple increments (ie, its not const)
std::vector< Real > _cumulative_pm
the sum of the plastic multipliers over all the sub-steps.
const MaterialProperty< RankTwoTensor > & _plastic_strain_old
Old value of plastic strain.
MaterialProperty< RankTwoTensor > * _couple_stress
the Cosserat couple-stress
virtual unsigned int numberActive(const std::vector< bool > &active)
counts the number of active constraints
RealTensorValue _rot
rotation matrix that takes _n to (0, 0, 1)
virtual bool quickStep(const RankTwoTensor &stress_old, RankTwoTensor &stress, const std::vector< Real > &intnl_old, std::vector< Real > &intnl, std::vector< Real > &pm, std::vector< Real > &cumulative_pm, const RankTwoTensor &plastic_strain_old, RankTwoTensor &plastic_strain, const RankFourTensor &E_ijkl, const RankTwoTensor &strain_increment, std::vector< Real > &yf, unsigned int &iterations, RankFourTensor &consistent_tangent_operator, const quickStep_called_from_t called_from, bool final_step)
Attempts to find an admissible (stress, intnl) by using the customized return-map algorithms defined ...
bool _cosserat
whether Cosserat mechanics should be used
void buildDumbOrder(const RankTwoTensor &stress, const std::vector< Real > &intnl, std::vector< unsigned int > &dumb_order)
Builds the order which "dumb" activation will take.
const MaterialProperty< RankFourTensor > * _elastic_flexural_rigidity_tensor
The Cosserat elastic flexural rigidity tensor.
enum ComputeMultiPlasticityStress::DeactivationSchemeEnum _deactivation_scheme
virtual bool plasticStep(const RankTwoTensor &stress_old, RankTwoTensor &stress, const std::vector< Real > &intnl_old, std::vector< Real > &intnl, const RankTwoTensor &plastic_strain_old, RankTwoTensor &plastic_strain, const RankFourTensor &E_ijkl, const RankTwoTensor &strain_increment, std::vector< Real > &yf, unsigned int &iterations, bool &linesearch_needed, bool &ld_encountered, bool &constraints_added, RankFourTensor &consistent_tangent_operator)
performs a plastic step
MaterialProperty< RankFourTensor > * _Jacobian_mult_couple
derivative of couple-stress w.r.t. curvature
const MaterialProperty< RankTwoTensor > & _total_strain_old
Old value of total strain (coming from ComputeIncrementalSmallStrain, for example) ...