www.mooseframework.org
ComputeMultiPlasticityStress.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"
14 
26 {
27 public:
29 
31 
32 protected:
33  virtual void computeQpStress();
34  virtual void initQpStatefulProperties();
35 
37  unsigned int _max_iter;
38 
41 
44 
47 
50  {
55 
58 
60  std::vector<Real> _dummy_pm;
61 
66  std::vector<Real> _cumulative_pm;
67 
68  /*
69  * Scheme by which constraints are deactivated.
70  * This enum is defined here for computational
71  * efficiency. If you add to this list you must
72  * also add to the MooseEnum in the .C file
73  */
75  {
84 
87 
90 
93 
96 
98  const std::string _elasticity_tensor_name;
101 
104 
107 
110 
113 
116 
119 
122 
125 
128 
131 
134 
137 
140 
143 
146 
149 
151  bool _cosserat;
152 
155 
158 
161 
164 
167 
170 
173 
176 
179 
183  virtual bool reinstateLinearDependentConstraints(std::vector<bool> & deactivated_due_to_ld);
184 
188  virtual unsigned int numberActive(const std::vector<bool> & active);
189 
201  virtual Real residual2(const std::vector<Real> & pm,
202  const std::vector<Real> & f,
203  const RankTwoTensor & epp,
204  const std::vector<Real> & ic,
205  const std::vector<bool> & active,
206  const std::vector<bool> & deactivated_due_to_ld);
207 
244  virtual bool returnMap(const RankTwoTensor & stress_old,
245  RankTwoTensor & stress,
246  const std::vector<Real> & intnl_old,
247  std::vector<Real> & intnl,
248  const RankTwoTensor & plastic_strain_old,
249  RankTwoTensor & plastic_strain,
250  const RankFourTensor & E_ijkl,
251  const RankTwoTensor & strain_increment,
252  std::vector<Real> & f,
253  unsigned int & iter,
254  bool can_revert_to_dumb,
255  bool & linesearch_needed,
256  bool & ld_encountered,
257  bool & constraints_added,
258  bool final_step,
259  RankFourTensor & consistent_tangent_operator,
260  std::vector<Real> & cumulative_pm);
261 
294  virtual bool lineSearch(Real & nr_res2,
295  RankTwoTensor & stress,
296  const std::vector<Real> & intnl_old,
297  std::vector<Real> & intnl,
298  std::vector<Real> & pm,
299  const RankFourTensor & E_inv,
300  RankTwoTensor & delta_dp,
301  const RankTwoTensor & dstress,
302  const std::vector<Real> & dpm,
303  const std::vector<Real> & dintnl,
304  std::vector<Real> & f,
305  RankTwoTensor & epp,
306  std::vector<Real> & ic,
307  const std::vector<bool> & active,
308  const std::vector<bool> & deactivated_due_to_ld,
309  bool & linesearch_needed);
310 
337  virtual bool singleStep(Real & nr_res2,
338  RankTwoTensor & stress,
339  const std::vector<Real> & intnl_old,
340  std::vector<Real> & intnl,
341  std::vector<Real> & pm,
342  RankTwoTensor & delta_dp,
343  const RankFourTensor & E_inv,
344  std::vector<Real> & f,
345  RankTwoTensor & epp,
346  std::vector<Real> & ic,
347  std::vector<bool> & active,
348  DeactivationSchemeEnum deactivation_scheme,
349  bool & linesearch_needed,
350  bool & ld_encountered);
351 
359  virtual bool checkAdmissible(const RankTwoTensor & stress,
360  const std::vector<Real> & intnl,
361  std::vector<Real> & all_f);
362 
372  void buildDumbOrder(const RankTwoTensor & stress,
373  const std::vector<Real> & intnl,
374  std::vector<unsigned int> & dumb_order);
375 
387  virtual void incrementDumb(int & dumb_iteration,
388  const std::vector<unsigned int> & dumb_order,
389  std::vector<bool> & act);
390 
408  virtual bool checkKuhnTucker(const std::vector<Real> & f,
409  const std::vector<Real> & pm,
410  const std::vector<bool> & active);
411 
429  virtual void applyKuhnTucker(const std::vector<Real> & f,
430  const std::vector<Real> & pm,
431  std::vector<bool> & active);
432 
433  // gets called before any return-map
434  virtual void preReturnMap();
435 
436  // gets called after return-map
437  virtual void postReturnMap();
438 
441  {
444  };
445 
484  virtual bool quickStep(const RankTwoTensor & stress_old,
485  RankTwoTensor & stress,
486  const std::vector<Real> & intnl_old,
487  std::vector<Real> & intnl,
488  std::vector<Real> & pm,
489  std::vector<Real> & cumulative_pm,
490  const RankTwoTensor & plastic_strain_old,
491  RankTwoTensor & plastic_strain,
492  const RankFourTensor & E_ijkl,
493  const RankTwoTensor & strain_increment,
494  std::vector<Real> & yf,
495  unsigned int & iterations,
496  RankFourTensor & consistent_tangent_operator,
497  const quickStep_called_from_t called_from,
498  bool final_step);
499 
526  virtual bool plasticStep(const RankTwoTensor & stress_old,
527  RankTwoTensor & stress,
528  const std::vector<Real> & intnl_old,
529  std::vector<Real> & intnl,
530  const RankTwoTensor & plastic_strain_old,
531  RankTwoTensor & plastic_strain,
532  const RankFourTensor & E_ijkl,
533  const RankTwoTensor & strain_increment,
534  std::vector<Real> & yf,
535  unsigned int & iterations,
536  bool & linesearch_needed,
537  bool & ld_encountered,
538  bool & constraints_added,
539  RankFourTensor & consistent_tangent_operator);
540 
541  // bool checkAndModifyConstraints(bool nr_exit_condition, const RankTwoTensor & stress, const
542  // std::vector<Real> & intnl, const std::vector<Real> & pm, const std::vector<bool> &
543  // initial_act, bool can_revert_to_dumb, const RankTwoTensor & initial_stress, const
544  // std::vector<Real> & intnl_old, const std::vector<Real> & f, DeactivationSchemeEnum
545  // deact_scheme, std::vector<bool> & act, int & dumb_iteration, std::vector<unsigned int>
546  // dumb_order, bool & die);
547 
548  bool canChangeScheme(DeactivationSchemeEnum current_deactivation_scheme, bool can_revert_to_dumb);
549 
550  bool canIncrementDumb(int dumb_iteration);
551 
552  void changeScheme(const std::vector<bool> & initial_act,
553  bool can_revert_to_dumb,
554  const RankTwoTensor & initial_stress,
555  const std::vector<Real> & intnl_old,
556  DeactivationSchemeEnum & current_deactivation_scheme,
557  std::vector<bool> & act,
558  int & dumb_iteration,
559  std::vector<unsigned int> & dumb_order);
560 
561  bool canAddConstraints(const std::vector<bool> & act, const std::vector<Real> & all_f);
562 
563  unsigned int activeCombinationNumber(const std::vector<bool> & act);
564 
583  const std::vector<Real> & intnl,
584  const RankFourTensor & E_ijkl,
585  const std::vector<Real> & pm_this_step,
586  const std::vector<Real> & cumulative_pm);
587 
588 private:
589  RankTwoTensor rot(const RankTwoTensor & tens);
590 };
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 computed from MOOSE&#39;s strain calculators...
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)
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.
const MaterialProperty< RankTwoTensor > *const _curvature
The Cosserat curvature strain.
bool _ignore_failures
Even if the returnMap fails, return the best values found for stress and internal parameters...
const MaterialProperty< RankFourTensor > *const _elastic_flexural_rigidity_tensor
The Cosserat elastic flexural rigidity tensor.
RankFourTensor _my_elasticity_tensor
Elasticity tensor that can be rotated by this class (ie, its not const)
const std::string _elasticity_tensor_name
Name of the elasticity tensor material property.
MaterialProperty< Real > & _iter
Number of Newton-Raphson iterations used in the return-map.
const MaterialProperty< RankTwoTensor > & _elastic_strain_old
Old value of elastic strain.
TensorValue< Real > RealTensorValue
MaterialProperty< RankFourTensor > *const _Jacobian_mult_couple
derivative of couple-stress w.r.t. curvature
Real f(Real x)
Test function for Brents method.
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...
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< RankTwoTensor > *const _couple_stress
the Cosserat couple-stress
MaterialProperty< Real > & _linesearch_needed
Whether a line-search was needed in the latest Newton-Raphson process (1 if true, 0 otherwise) ...
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
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void computeQpStress()
Compute the stress and store it in the _stress material property for the current quadrature point...
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.
MultiPlasticityDebugger computes various finite-difference things to help developers remove bugs in t...
const InputParameters & parameters() const
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.
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.
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
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor material property.
const MaterialProperty< RankTwoTensor > *const _couple_stress_old
the old value of Cosserat couple-stress
const MaterialProperty< RankTwoTensor > & _total_strain_old
Old value of total strain (coming from ComputeIncrementalSmallStrain, for example) ...