www.mooseframework.org
TensileStressUpdate.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 TENSILESTRESSUPDATE_H
8 #define TENSILESTRESSUPDATE_H
9 
12 
14 
15 template <>
16 InputParameters validParams<TensileStressUpdate>();
17 
23 {
24 public:
25  TensileStressUpdate(const InputParameters & parameters);
26 
30  bool requiresIsotropicTensor() override { return true; }
31 
32 protected:
35 
37  const bool _perfect_guess;
38 
40  RankTwoTensor _eigvecs;
41 
43  virtual Real tensile_strength(const Real internal_param) const;
44 
46  virtual Real dtensile_strength(const Real internal_param) const;
47 
48  void computeStressParams(const RankTwoTensor & stress,
49  std::vector<Real> & stress_params) const override;
50 
51  std::vector<RankTwoTensor> dstress_param_dstress(const RankTwoTensor & stress) const override;
52 
53  std::vector<RankFourTensor> d2stress_param_dstress(const RankTwoTensor & stress) const override;
54 
55  virtual void setStressAfterReturnV(const RankTwoTensor & stress_trial,
56  const std::vector<Real> & stress_params,
57  Real gaE,
58  const std::vector<Real> & intnl,
59  const yieldAndFlow & smoothed_q,
60  const RankFourTensor & Eijkl,
61  RankTwoTensor & stress) const override;
62 
63  virtual void preReturnMapV(const std::vector<Real> & trial_stress_params,
64  const RankTwoTensor & stress_trial,
65  const std::vector<Real> & intnl_old,
66  const std::vector<Real> & yf,
67  const RankFourTensor & Eijkl) override;
68 
69  void setEffectiveElasticity(const RankFourTensor & Eijkl) override;
70 
71  void yieldFunctionValuesV(const std::vector<Real> & stress_params,
72  const std::vector<Real> & intnl,
73  std::vector<Real> & yf) const override;
74 
75  void computeAllQV(const std::vector<Real> & stress_params,
76  const std::vector<Real> & intnl,
77  std::vector<yieldAndFlow> & all_q) const override;
78 
79  void initializeVarsV(const std::vector<Real> & trial_stress_params,
80  const std::vector<Real> & intnl_old,
81  std::vector<Real> & stress_params,
82  Real & gaE,
83  std::vector<Real> & intnl) const override;
84 
85  void setIntnlValuesV(const std::vector<Real> & trial_stress_params,
86  const std::vector<Real> & current_stress_params,
87  const std::vector<Real> & intnl_old,
88  std::vector<Real> & intnl) const override;
89 
90  void setIntnlDerivativesV(const std::vector<Real> & trial_stress_params,
91  const std::vector<Real> & current_stress_params,
92  const std::vector<Real> & intnl,
93  std::vector<std::vector<Real>> & dintnl) const override;
94 
95  virtual void consistentTangentOperatorV(const RankTwoTensor & stress_trial,
96  const std::vector<Real> & trial_stress_params,
97  const RankTwoTensor & stress,
98  const std::vector<Real> & stress_params,
99  Real gaE,
100  const yieldAndFlow & smoothed_q,
101  const RankFourTensor & Eijkl,
102  bool compute_full_tangent_operator,
103  const std::vector<std::vector<Real>> & dvar_dtrial,
104  RankFourTensor & cto) override;
105 };
106 
107 #endif // TENSILESTRESSUPDATE_H
TensileStressUpdate(const InputParameters &parameters)
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) override
Derived classes may employ this function to record stuff or do other computations prior to the return...
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) override
Calculates the consistent tangent operator.
std::vector< RankFourTensor > d2stress_param_dstress(const RankTwoTensor &stress) const override
d2(stress_param[i])/d(stress)/d(stress) at given stress
InputParameters validParams< TensileStressUpdate >()
Struct designed to hold info about a single yield function and its derivatives, as well as the flow d...
TensileStressUpdate implements rate-independent associative tensile failure ("Rankine" plasticity) wi...
void setEffectiveElasticity(const RankFourTensor &Eijkl) override
Sets _Eij and _En and _Cij.
RankTwoTensor _eigvecs
Eigenvectors of the trial stress as a RankTwoTensor, in order to rotate the returned stress back to s...
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 override
Sets the internal parameters based on the trial values of stress_params, their current values...
const TensorMechanicsHardeningModel & _strength
Hardening model for tensile strength.
void yieldFunctionValuesV(const std::vector< Real > &stress_params, const std::vector< Real > &intnl, std::vector< Real > &yf) const override
Computes the values of the yield functions, given stress_params and intnl parameters.
std::vector< RankTwoTensor > dstress_param_dstress(const RankTwoTensor &stress) const override
d(stress_param[i])/d(stress) at given stress
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 override
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.
void computeStressParams(const RankTwoTensor &stress, std::vector< Real > &stress_params) const override
Computes stress_params, given stress.
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 override
Sets (stress_params, intnl) at "good guesses" of the solution to the Return-Map algorithm.
virtual Real dtensile_strength(const Real internal_param) const
d(tensile strength)/d(internal_param) as a function of residual value, rate, and internal_param ...
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 override
Sets stress from the admissible parameters.
virtual Real tensile_strength(const Real internal_param) const
tensile strength as a function of residual value, rate, and internal_param
const bool _perfect_guess
Whether to provide an estimate of the returned stress, based on perfect plasticity.
MultiParameterPlasticityStressUpdate performs the return-map algorithm and associated stress updates ...
void computeAllQV(const std::vector< Real > &stress_params, const std::vector< Real > &intnl, std::vector< yieldAndFlow > &all_q) const override
Completely fills all_q with correct values.
bool requiresIsotropicTensor() override
Does the model require the elasticity tensor to be isotropic?