www.mooseframework.org
MultiParameterPlasticityStressUpdate.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 MULTIPARAMETERPLASTICITYSTRESSUPDATE_H
8 #define MULTIPARAMETERPLASTICITYSTRESSUPDATE_H
9 
10 #include "StressUpdateBase.h"
11 
12 #include <array>
13 
15 
16 template <>
18 
98 {
99 public:
100  MultiParameterPlasticityStressUpdate(const InputParameters & parameters,
101  unsigned num_sp,
102  unsigned num_yf,
103  unsigned num_intnl);
104 
105 protected:
106  virtual void initQpStatefulProperties() override;
107  virtual void updateState(RankTwoTensor & strain_increment,
108  RankTwoTensor & inelastic_strain_increment,
109  const RankTwoTensor & rotation_increment,
110  RankTwoTensor & stress_new,
111  const RankTwoTensor & stress_old,
112  const RankFourTensor & elasticity_tensor,
113  const RankTwoTensor & elastic_strain_old,
114  bool compute_full_tangent_operator,
115  RankFourTensor & tangent_operator) override;
116 
117  virtual void propagateQpStatefulProperties() override;
118 
120  constexpr static unsigned _tensor_dimensionality = 3;
121 
123  const unsigned _num_sp;
124 
126  const std::vector<Real> _definitely_ok_sp;
127 
129  std::vector<std::vector<Real>> _Eij;
130 
132  Real _En;
133 
135  std::vector<std::vector<Real>> _Cij;
136 
138  const unsigned _num_yf;
139 
141  const unsigned _num_intnl;
142 
144  const std::string _base_name;
145 
147  const unsigned _max_nr_its;
148 
151 
153  const Real _smoothing_tol;
154 
156  const Real _f_tol;
157 
159  const Real _f_tol2;
160 
166  const Real _min_step_size;
167 
169  bool _step_one;
170 
173 
175  MaterialProperty<RankTwoTensor> & _plastic_strain;
176 
178  const MaterialProperty<RankTwoTensor> & _plastic_strain_old;
179 
181  MaterialProperty<std::vector<Real>> & _intnl;
182 
184  const MaterialProperty<std::vector<Real>> & _intnl_old;
185 
187  MaterialProperty<std::vector<Real>> & _yf;
188 
190  MaterialProperty<Real> & _iter;
191 
193  MaterialProperty<Real> & _linesearch_needed;
194 
200  {
201  Real f; // yield function value
202  std::vector<Real> df; // df/d(stress_param[i])
203  std::vector<Real> df_di; // df/d(intnl_variable[a])
204  std::vector<Real> dg; // d(flow)/d(stress_param[i])
205  std::vector<std::vector<Real>> d2g; // d^2(flow)/d(sp[i])/d(sp[j])
206  std::vector<std::vector<Real>> d2g_di; // d^2(flow)/d(sp[i])/d(intnl[a])
207 
209 
210  yieldAndFlow(unsigned num_var, unsigned num_intnl)
211  : f(0.0),
212  df(num_var),
213  df_di(num_intnl),
214  dg(num_var),
215  d2g(num_var, std::vector<Real>(num_var, 0.0)),
216  d2g_di(num_var, std::vector<Real>(num_intnl, 0.0))
217  {
218  }
219 
220  // this may be involved in the smoothing of a group of yield functions
221  bool operator<(const yieldAndFlow & fd) const { return f < fd.f; }
222  };
223 
231  Real yieldF(const std::vector<Real> & stress_params, const std::vector<Real> & intnl) const;
232 
238  Real yieldF(const std::vector<Real> & yfs) const;
239 
251  Real ismoother(Real f_diff) const;
252 
256  Real smoother(Real f_diff) const;
257 
261  Real dsmoother(Real f_diff) const;
262 
270  yieldAndFlow smoothAllQuantities(const std::vector<Real> & stress_params,
271  const std::vector<Real> & intnl) const;
272 
312  int lineSearch(Real & res2,
313  std::vector<Real> & stress_params,
314  Real & gaE,
315  const std::vector<Real> & trial_stress_params,
316  yieldAndFlow & smoothed_q,
317  const std::vector<Real> & intnl_ok,
318  std::vector<Real> & intnl,
319  std::vector<Real> & rhs,
320  Real & linesearch_needed) const;
321 
335  int nrStep(const yieldAndFlow & smoothed_q,
336  const std::vector<Real> & trial_stress_params,
337  const std::vector<Real> & stress_params,
338  const std::vector<Real> & intnl,
339  Real gaE,
340  std::vector<Real> & rhs) const;
341 
347  Real calculateRes2(const std::vector<Real> & rhs) const;
348 
365  void calculateRHS(const std::vector<Real> & trial_stress_params,
366  const std::vector<Real> & stress_params,
367  Real gaE,
368  const yieldAndFlow & smoothed_q,
369  std::vector<Real> & rhs) const;
370 
384  void dnRHSdVar(const yieldAndFlow & smoothed_q,
385  const std::vector<std::vector<Real>> & dintnl,
386  const std::vector<Real> & stress_params,
387  Real gaE,
388  std::vector<double> & jac) const;
389 
394  virtual void errorHandler(const std::string & message) const;
395 
404  virtual void yieldFunctionValuesV(const std::vector<Real> & stress_params,
405  const std::vector<Real> & intnl,
406  std::vector<Real> & yf) const = 0;
407 
420  virtual void computeAllQV(const std::vector<Real> & stress_params,
421  const std::vector<Real> & intnl,
422  std::vector<yieldAndFlow> & all_q) const = 0;
423 
435  virtual void preReturnMapV(const std::vector<Real> & trial_stress_params,
436  const RankTwoTensor & stress_trial,
437  const std::vector<Real> & intnl_old,
438  const std::vector<Real> & yf,
439  const RankFourTensor & Eijkl);
440 
456  virtual void initializeVarsV(const std::vector<Real> & trial_stress_params,
457  const std::vector<Real> & intnl_old,
458  std::vector<Real> & stress_params,
459  Real & gaE,
460  std::vector<Real> & intnl) const;
461 
472  virtual void setIntnlValuesV(const std::vector<Real> & trial_stress_params,
473  const std::vector<Real> & current_stress_params,
474  const std::vector<Real> & intnl_old,
475  std::vector<Real> & intnl) const = 0;
476 
487  virtual void setIntnlDerivativesV(const std::vector<Real> & trial_stress_params,
488  const std::vector<Real> & current_stress_params,
489  const std::vector<Real> & intnl,
490  std::vector<std::vector<Real>> & dintnl) const = 0;
491 
498  virtual void computeStressParams(const RankTwoTensor & stress,
499  std::vector<Real> & stress_params) const = 0;
500 
508  virtual void initializeReturnProcess();
509 
516  virtual void finalizeReturnProcess(const RankTwoTensor & rotation_increment);
517 
535  virtual void setStressAfterReturnV(const RankTwoTensor & stress_trial,
536  const std::vector<Real> & stress_params,
537  Real gaE,
538  const std::vector<Real> & intnl,
539  const yieldAndFlow & smoothed_q,
540  const RankFourTensor & Eijkl,
541  RankTwoTensor & stress) const = 0;
542 
558  virtual void
559  setInelasticStrainIncrementAfterReturn(const RankTwoTensor & stress_trial,
560  Real gaE,
561  const yieldAndFlow & smoothed_q,
562  const RankFourTensor & elasticity_tensor,
563  const RankTwoTensor & returned_stress,
564  RankTwoTensor & inelastic_strain_increment) const;
565 
587  virtual void consistentTangentOperatorV(const RankTwoTensor & stress_trial,
588  const std::vector<Real> & trial_stress_params,
589  const RankTwoTensor & stress,
590  const std::vector<Real> & stress_params,
591  Real gaE,
592  const yieldAndFlow & smoothed_q,
593  const RankFourTensor & Eijkl,
594  bool compute_full_tangent_operator,
595  const std::vector<std::vector<Real>> & dvar_dtrial,
596  RankFourTensor & cto);
597 
603  virtual std::vector<RankTwoTensor> dstress_param_dstress(const RankTwoTensor & stress) const = 0;
604 
610  virtual std::vector<RankFourTensor>
611  d2stress_param_dstress(const RankTwoTensor & stress) const = 0;
612 
614  virtual void setEffectiveElasticity(const RankFourTensor & Eijkl) = 0;
615 
637  void dVardTrial(bool elastic_only,
638  const std::vector<Real> & trial_stress_params,
639  const std::vector<Real> & stress_params,
640  Real gaE,
641  const std::vector<Real> & intnl,
642  const yieldAndFlow & smoothed_q,
643  Real step_size,
644  bool compute_full_tangent_operator,
645  std::vector<std::vector<Real>> & dvar_dtrial) const;
646 
654  bool precisionLoss(const std::vector<Real> & solution,
655  const std::vector<Real> & stress_params,
656  Real gaE) const;
657 
658 private:
666  std::vector<Real> _trial_sp;
667 
673  RankTwoTensor _stress_trial;
674 
683  std::vector<Real> _rhs;
684 
688  std::vector<std::vector<Real>> _dvar_dtrial;
689 
698  std::vector<Real> _ok_sp;
699 
703  std::vector<Real> _ok_intnl;
704 
712  std::vector<Real> _del_stress_params;
713 
717  std::vector<Real> _current_sp;
718 
722  std::vector<Real> _current_intnl;
723 };
724 
725 #endif // MULTIPARAMETERPLASTICITYSTRESSUPDATE_H
int nrStep(const yieldAndFlow &smoothed_q, const std::vector< Real > &trial_stress_params, const std::vector< Real > &stress_params, const std::vector< Real > &intnl, Real gaE, std::vector< Real > &rhs) const
Performs a Newton-Raphson step to attempt to zero rhs Upon return, rhs will contain the solution...
virtual void initializeVarsV(const std::vector< Real > &trial_stress_params, const std::vector< Real > &intnl_old, std::vector< Real > &stress_params, Real &gaE, std::vector< Real > &intnl) const
Sets (stress_params, intnl) at "good guesses" of the solution to the Return-Map algorithm.
MaterialProperty< std::vector< Real > > & _yf
yield functions
std::vector< Real > _trial_sp
"Trial" value of stress_params that initializes the return-map process This is derived from stress = ...
virtual void computeAllQV(const std::vector< Real > &stress_params, const std::vector< Real > &intnl, std::vector< yieldAndFlow > &all_q) const =0
Completely fills all_q with correct values.
MultiParameterPlasticityStressUpdate(const InputParameters &parameters, unsigned num_sp, unsigned num_yf, unsigned num_intnl)
virtual void setEffectiveElasticity(const RankFourTensor &Eijkl)=0
Sets _Eij and _En and _Cij.
const bool _warn_about_precision_loss
Output a warning message if precision loss is encountered during the return-map process.
InputParameters validParams< MultiParameterPlasticityStressUpdate >()
Real smoother(Real f_diff) const
Derivative of ismoother.
StressUpdateBase is a material that is not called by MOOSE because of the compute=false flag set in t...
const unsigned _num_intnl
Number of internal parameters.
Struct designed to hold info about a single yield function and its derivatives, as well as the flow d...
void calculateRHS(const std::vector< Real > &trial_stress_params, const std::vector< Real > &stress_params, Real gaE, const yieldAndFlow &smoothed_q, std::vector< Real > &rhs) const
Calculates the RHS in the following 0 = rhs[0] = S[0] - S[0]^trial + ga * E[0, j] * dg/dS[j] 0 = rhs[...
std::vector< std::vector< Real > > _Cij
_Cij[i, j] * _Eij[j, k] = 1 iff j == k
virtual void consistentTangentOperatorV(const RankTwoTensor &stress_trial, const std::vector< Real > &trial_stress_params, const RankTwoTensor &stress, const std::vector< Real > &stress_params, Real gaE, const yieldAndFlow &smoothed_q, const RankFourTensor &Eijkl, bool compute_full_tangent_operator, const std::vector< std::vector< Real >> &dvar_dtrial, RankFourTensor &cto)
Calculates the consistent tangent operator.
virtual void setStressAfterReturnV(const RankTwoTensor &stress_trial, const std::vector< Real > &stress_params, Real gaE, const std::vector< Real > &intnl, const yieldAndFlow &smoothed_q, const RankFourTensor &Eijkl, RankTwoTensor &stress) const =0
Sets stress from the admissible parameters.
virtual void yieldFunctionValuesV(const std::vector< Real > &stress_params, const std::vector< Real > &intnl, std::vector< Real > &yf) const =0
Computes the values of the yield functions, given stress_params and intnl parameters.
std::vector< Real > _ok_intnl
The state (ok_sp, ok_intnl) is known to be admissible.
const Real _f_tol2
Square of the yield-function tolerance.
const Real _min_step_size
In order to help the Newton-Raphson procedure, the applied strain increment may be applied in sub-inc...
virtual void updateState(RankTwoTensor &strain_increment, RankTwoTensor &inelastic_strain_increment, const RankTwoTensor &rotation_increment, RankTwoTensor &stress_new, const RankTwoTensor &stress_old, const RankFourTensor &elasticity_tensor, const RankTwoTensor &elastic_strain_old, bool compute_full_tangent_operator, RankFourTensor &tangent_operator) override
Given a strain increment that results in a trial stress, perform some procedure (such as an iterative...
MaterialProperty< std::vector< Real > > & _intnl
internal parameters
std::vector< std::vector< Real > > _Eij
E[i, j] in the system of equations to be solved.
const Real _f_tol
The yield-function tolerance.
const unsigned _num_yf
Number of yield functions.
int lineSearch(Real &res2, std::vector< Real > &stress_params, Real &gaE, const std::vector< Real > &trial_stress_params, yieldAndFlow &smoothed_q, const std::vector< Real > &intnl_ok, std::vector< Real > &intnl, std::vector< Real > &rhs, Real &linesearch_needed) const
Performs a line-search to find stress_params and gaE Upon entry:
Real ismoother(Real f_diff) const
Smooths yield functions.
const std::vector< Real > _definitely_ok_sp
An admissible value of stress_params at the initial time.
const MaterialProperty< std::vector< Real > > & _intnl_old
old values of internal parameters
RankTwoTensor _stress_trial
"Trial" value of stress that is set at the beginning of the return-map process.
virtual void propagateQpStatefulProperties() override
If updateState is not called during a timestep, this will be.
virtual void setInelasticStrainIncrementAfterReturn(const RankTwoTensor &stress_trial, Real gaE, const yieldAndFlow &smoothed_q, const RankFourTensor &elasticity_tensor, const RankTwoTensor &returned_stress, RankTwoTensor &inelastic_strain_increment) const
Sets inelastic strain increment from the returned configuration This is called after the return-map p...
const std::string _base_name
String prepended to various MaterialProperties that are defined by this class.
yieldAndFlow smoothAllQuantities(const std::vector< Real > &stress_params, const std::vector< Real > &intnl) const
Calculates all yield functions and derivatives, and then performs the smoothing scheme.
virtual void computeStressParams(const RankTwoTensor &stress, std::vector< Real > &stress_params) const =0
Computes stress_params, given stress.
MaterialProperty< Real > & _iter
Number of Newton-Raphson iterations used in the return-map.
MaterialProperty< RankTwoTensor > & _plastic_strain
plastic strain
const Real _smoothing_tol
Smoothing tolerance: edges of the yield surface get smoothed by this amount.
Real dsmoother(Real f_diff) const
Derivative of smoother.
virtual std::vector< RankTwoTensor > dstress_param_dstress(const RankTwoTensor &stress) const =0
d(stress_param[i])/d(stress) at given stress
std::vector< Real > _current_intnl
The current values of the internal params during the Newton-Raphson.
virtual void finalizeReturnProcess(const RankTwoTensor &rotation_increment)
Derived classes may use this to perform calculations after the return-map process has completed succe...
const bool _perform_finite_strain_rotations
Whether to perform finite-strain rotations.
Real yieldF(const std::vector< Real > &stress_params, const std::vector< Real > &intnl) const
Computes the smoothed yield function.
virtual std::vector< RankFourTensor > d2stress_param_dstress(const RankTwoTensor &stress) const =0
d2(stress_param[i])/d(stress)/d(stress) at given stress
virtual void setIntnlDerivativesV(const std::vector< Real > &trial_stress_params, const std::vector< Real > &current_stress_params, const std::vector< Real > &intnl, std::vector< std::vector< Real >> &dintnl) const =0
Sets the derivatives of internal parameters, based on the trial values of stress_params, their current values, and the current values of the internal parameters.
bool precisionLoss(const std::vector< Real > &solution, const std::vector< Real > &stress_params, Real gaE) const
Check whether precision loss has occurred.
MaterialProperty< Real > & _linesearch_needed
Whether a line-search was needed in the latest Newton-Raphson process (1 if true, 0 otherwise) ...
MultiParameterPlasticityStressUpdate performs the return-map algorithm and associated stress updates ...
std::vector< std::vector< Real > > _dvar_dtrial
d({stress_param[i], gaE})/d(trial_stress_param[j])
Real calculateRes2(const std::vector< Real > &rhs) const
Calculates the residual-squared for the Newton-Raphson + line-search.
virtual void setIntnlValuesV(const std::vector< Real > &trial_stress_params, const std::vector< Real > &current_stress_params, const std::vector< Real > &intnl_old, std::vector< Real > &intnl) const =0
Sets the internal parameters based on the trial values of stress_params, their current values...
std::vector< Real > _current_sp
The current values of the stress params during the Newton-Raphson.
bool _step_one
handles case of initial_stress that is inadmissible being supplied
virtual void errorHandler(const std::string &message) const
Performs any necessary cleaning-up, then throw MooseException(message)
void dnRHSdVar(const yieldAndFlow &smoothed_q, const std::vector< std::vector< Real >> &dintnl, const std::vector< Real > &stress_params, Real gaE, std::vector< double > &jac) const
Derivative of -RHS with respect to the stress_params and gaE, placed into an array ready for solving ...
virtual void preReturnMapV(const std::vector< Real > &trial_stress_params, const RankTwoTensor &stress_trial, const std::vector< Real > &intnl_old, const std::vector< Real > &yf, const RankFourTensor &Eijkl)
Derived classes may employ this function to record stuff or do other computations prior to the return...
std::vector< Real > _del_stress_params
_del_stress_params = trial_stress_params - ok_sp This is fixed at the beginning of the return-map pro...
const unsigned _max_nr_its
Maximum number of Newton-Raphson iterations allowed in the return-map process.
virtual void initializeReturnProcess()
Derived classes may use this to perform calculations before any return-map process is performed...
const MaterialProperty< RankTwoTensor > & _plastic_strain_old
Old value of plastic strain.
static constexpr unsigned _tensor_dimensionality
Internal dimensionality of tensors (currently this is 3 throughout tensor_mechanics) ...
std::vector< Real > _rhs
0 = rhs[0] = S[0] - S[0]^trial + ga * E[0, i] * dg/dS[i] 0 = rhs[1] = S[1] - S[1]^trial + ga * E[1...
std::vector< Real > _ok_sp
The state (ok_sp, ok_intnl) is known to be admissible, so ok_sp are stress_params that are "OK"...
void dVardTrial(bool elastic_only, const std::vector< Real > &trial_stress_params, const std::vector< Real > &stress_params, Real gaE, const std::vector< Real > &intnl, const yieldAndFlow &smoothed_q, Real step_size, bool compute_full_tangent_operator, std::vector< std::vector< Real >> &dvar_dtrial) const
Calculates derivatives of the stress_params and gaE with repect to the trial values of the stress_par...
const unsigned _num_sp
Number of stress parameters.