www.mooseframework.org
SolidModel.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 SOLIDMODEL_H
8 #define SOLIDMODEL_H
9 
10 #include "DerivativeMaterialInterface.h"
11 #include "SymmTensor.h"
12 #include "RankTwoTensor.h"
13 
14 // Forward declarations
15 class ConstitutiveModel;
16 class SolidModel;
18 class PiecewiseLinear;
19 class VolumetricModel;
20 namespace SolidMechanics
21 {
22 class Element;
23 }
24 
25 template <>
26 InputParameters validParams<SolidModel>();
27 
31 class SolidModel : public DerivativeMaterialInterface<Material>
32 {
33 public:
34  SolidModel(const InputParameters & parameters);
35  virtual ~SolidModel();
36 
37  virtual void initStatefulProperties(unsigned n_points);
38 
39  virtual void applyThermalStrain();
40  virtual void applyVolumetricStrain();
41 
42  static void
43  rotateSymmetricTensor(const ColumnMajorMatrix & R, const SymmTensor & T, SymmTensor & result);
44 
46  {
47  CR_ABRUPT = 0,
50  CR_UNKNOWN
51  };
52 
53  QBase * qrule() { return _qrule; }
54  const Point & q_point(unsigned i) const { return _q_point[i]; }
55  Real JxW(unsigned i) const { return _JxW[i]; }
56 
57 protected:
58  Moose::CoordinateSystemType _coord_type;
59 
60  const std::string _appended_property_name;
61 
67 
69  Real _lambda;
73 
76 
80  const Real _cracking_beta;
81  const std::string _compute_method;
82  Function * const _cracking_stress_function;
83 
85  std::vector<unsigned int> _active_crack_planes;
86  const unsigned int _max_cracks;
88  // std::map<Point, unsigned> _cracked_this_step_count;
89  // std::map<Point, unsigned> _cracked_this_step;
90 
91  const bool _has_temp;
92  const VariableValue & _temperature;
93  const VariableValue & _temperature_old;
94  const VariableGradient & _temp_grad;
95  const Real _alpha;
96  Function * _alpha_function;
101  Real _ref_temp;
102 
103  std::map<SubdomainID, std::vector<MooseSharedPointer<VolumetricModel>>> _volumetric_models;
104  std::set<std::string> _dep_matl_props;
105 
106  MaterialProperty<SymmTensor> & _stress;
107 
108 private:
109  const MaterialProperty<SymmTensor> & _stress_old_prop;
110 
111 protected:
113 
114  MaterialProperty<SymmTensor> & _total_strain;
115  const MaterialProperty<SymmTensor> & _total_strain_old;
116 
117  MaterialProperty<SymmTensor> & _elastic_strain;
118  const MaterialProperty<SymmTensor> & _elastic_strain_old;
119 
120  MaterialProperty<RealVectorValue> * _crack_flags;
121  const MaterialProperty<RealVectorValue> * _crack_flags_old;
122  RealVectorValue _crack_flags_local;
123  MaterialProperty<RealVectorValue> * _crack_count;
124  const MaterialProperty<RealVectorValue> * _crack_count_old;
125  MaterialProperty<ColumnMajorMatrix> * _crack_rotation;
126  const MaterialProperty<ColumnMajorMatrix> * _crack_rotation_old;
127  MaterialProperty<RealVectorValue> * _crack_strain;
128  const MaterialProperty<RealVectorValue> * _crack_strain_old;
129  MaterialProperty<RealVectorValue> * _crack_max_strain;
130  const MaterialProperty<RealVectorValue> * _crack_max_strain_old;
131  ColumnMajorMatrix _principal_strain;
132 
133  MaterialProperty<SymmElasticityTensor> & _elasticity_tensor;
134  MaterialProperty<SymmElasticityTensor> & _Jacobian_mult;
135 
136  // Accumulate derivatives of strain tensors with respect to Temperature into this
138 
139  // The derivative of the stress with respect to Temperature
140  MaterialProperty<SymmTensor> & _d_stress_dT;
141 
144 
145  const bool _compute_JIntegral;
148 
149  // These are used in calculation of the J integral
150  MaterialProperty<Real> * _SED;
151  const MaterialProperty<Real> * _SED_old;
152  MaterialProperty<RankTwoTensor> * _Eshelby_tensor;
153  MaterialProperty<RealVectorValue> * _J_thermal_term_vec;
154 
155  // This is used in calculation of the J Integral and Interaction Integral
157 
158  virtual void initQpStatefulProperties();
159 
160  virtual void initialSetup();
161  virtual void timestepSetup();
162  virtual void jacobianSetup();
163 
164  virtual void computeProperties();
165 
166  void computeElasticityTensor();
170  virtual bool updateElasticityTensor(SymmElasticityTensor & tensor);
171 
172  virtual void elementInit() {}
173 
175  virtual void modifyStrainIncrement();
176 
178  virtual void crackingStrainDirections();
179 
180  virtual unsigned int getNumKnownCrackDirs() const;
181 
183  virtual void computeStress()
184  {
185  mooseError("SolidModel::computeStress must be defined by the derived class");
186  }
187 
188  // Compute Eshelby tensor, used in J Integral calculation
189  virtual void computeEshelby();
190 
191  // Compute strain energy density, used in Eshelby tensor calculation
192  virtual void computeStrainEnergyDensity();
193 
194  // Compute quantity used in thermal term of J Integral
195  virtual void computeThermalJvec();
196 
197  // Compute current thermal expansion coefficient, used in J Integral and Interaction Integral
198  virtual void computeCurrentInstantaneousThermalExpansionCoefficient();
199 
200  /*
201  * Determine whether new cracks have formed.
202  * Rotate old and new stress to global, if cracking active
203  */
204  virtual void crackingStressRotation();
205 
206  virtual Real computeCrackFactor(int i, Real & sigma, Real & flagVal);
207 
209  virtual void finalizeStress();
210 
211  virtual void computePreconditioning();
212 
213  void applyCracksToTensor(SymmTensor & tensor, const RealVectorValue & sigma);
214 
215  void elasticityTensor(SymmElasticityTensor * e);
216 
217  SymmElasticityTensor * elasticityTensor() const { return _local_elasticity_tensor; }
218 
219  const SolidMechanics::Element * element() const { return _element; }
220 
221  int delta(int i, int j) const { return i == j; }
222 
223  template <typename T>
224  MaterialProperty<T> & createProperty(const std::string & prop_name)
225  {
226  std::string name(prop_name + _appended_property_name);
227  return declareProperty<T>(name);
228  }
229 
230  template <typename T>
231  const MaterialProperty<T> & getPropertyOld(const std::string & prop_name)
232  {
233  std::string name(prop_name + _appended_property_name);
234  return getMaterialPropertyOld<T>(name);
235  }
236 
237  virtual void checkElasticConstants();
238 
239  virtual void createElasticityTensor();
240 
241  std::vector<SubdomainID> _block_id;
242 
243  std::map<SubdomainID, MooseSharedPointer<ConstitutiveModel>> _constitutive_model;
244  // This set keeps track of the dynamic memory allocated in this object
245  std::set<MooseSharedPointer<ConstitutiveModel>> _models_to_free;
247 
249  virtual void computeConstitutiveModelStress();
250 
251  void createConstitutiveModel(const std::string & cm_name);
252 
254  bool & _step_zero;
255  bool & _step_one;
257 
258 private:
259  void computeCrackStrainAndOrientation(ColumnMajorMatrix & principal_strain);
260 
261  SolidMechanics::Element * createElement();
262 
264 
266 };
267 
268 #endif
This class defines a basic set of capabilities any elasticity tensor should have. ...
SymmElasticityTensor * elasticityTensor() const
Definition: SolidModel.h:217
const std::string _appended_property_name
Definition: SolidModel.h:60
std::set< std::string > _dep_matl_props
Definition: SolidModel.h:104
MaterialProperty< ColumnMajorMatrix > * _crack_rotation
Definition: SolidModel.h:125
const MaterialProperty< RealVectorValue > * _crack_count_old
Definition: SolidModel.h:124
Real _poissons_ratio
Definition: SolidModel.h:70
const MaterialProperty< Real > * _SED_old
Definition: SolidModel.h:151
MaterialProperty< Real > * _SED
Definition: SolidModel.h:150
const bool _compute_JIntegral
Definition: SolidModel.h:145
const unsigned int _max_cracks
Definition: SolidModel.h:86
const Real _cracking_residual_stress
Definition: SolidModel.h:79
const MaterialProperty< T > & getPropertyOld(const std::string &prop_name)
Definition: SolidModel.h:231
Real _cracking_stress
Definition: SolidModel.h:78
const Real _alpha
Definition: SolidModel.h:95
SymmElasticityTensor * _local_elasticity_tensor
Definition: SolidModel.h:265
bool _constitutive_active
Definition: SolidModel.h:246
SymmTensor _total_strain_increment
Definition: SolidModel.h:142
Element is the base class for all of this module&#39;s solid mechanics element formulations.
Definition: Element.h:23
MaterialProperty< SymmElasticityTensor > & _Jacobian_mult
Definition: SolidModel.h:134
MaterialProperty< RealVectorValue > * _crack_count
Definition: SolidModel.h:123
const bool _compute_InteractionIntegral
Definition: SolidModel.h:146
const VariableValue & _temperature_old
Definition: SolidModel.h:93
MaterialProperty< T > & createProperty(const std::string &prop_name)
Definition: SolidModel.h:224
const Point & q_point(unsigned i) const
Definition: SolidModel.h:54
Real _stress_free_temp
Definition: SolidModel.h:99
Real _bulk_modulus
Definition: SolidModel.h:68
bool _youngs_modulus_set
Definition: SolidModel.h:66
bool _poissons_ratio_set
Definition: SolidModel.h:64
bool _store_stress_older
Definition: SolidModel.h:147
Function * _youngs_modulus_function
Definition: SolidModel.h:74
const CRACKING_RELEASE _cracking_release
Definition: SolidModel.h:77
Function *const _cracking_stress_function
Definition: SolidModel.h:82
MaterialProperty< SymmTensor > & _total_strain
Definition: SolidModel.h:114
std::vector< SubdomainID > _block_id
Definition: SolidModel.h:241
const MaterialProperty< SymmTensor > & _total_strain_old
Definition: SolidModel.h:115
const MaterialProperty< SymmTensor > & _stress_old_prop
Definition: SolidModel.h:109
MaterialProperty< SymmTensor > & _d_stress_dT
Definition: SolidModel.h:140
SolidModel is the base class for all this module&#39;s solid mechanics material models.
Definition: SolidModel.h:31
QBase * qrule()
Definition: SolidModel.h:53
MaterialProperty< SymmElasticityTensor > & _elasticity_tensor
Definition: SolidModel.h:133
const std::string _compute_method
Definition: SolidModel.h:81
virtual void elementInit()
Definition: SolidModel.h:172
bool & _step_zero
Restartable data to check for the zeroth and first time steps for thermal calculations.
Definition: SolidModel.h:254
const VariableGradient & _temp_grad
Definition: SolidModel.h:94
MaterialProperty< RealVectorValue > * _J_thermal_term_vec
Definition: SolidModel.h:153
bool & _step_one
Definition: SolidModel.h:255
bool _mean_alpha_function
Definition: SolidModel.h:100
bool _shear_modulus_set
Definition: SolidModel.h:65
const MaterialProperty< RealVectorValue > * _crack_flags_old
Definition: SolidModel.h:121
PiecewiseLinear * _piecewise_linear_alpha_function
Definition: SolidModel.h:97
Real _cracking_alpha
Definition: SolidModel.h:84
Function * _poissons_ratio_function
Definition: SolidModel.h:75
InputParameters validParams< SolidModel >()
Definition: SolidModel.C:27
const Real _cracking_neg_fraction
Definition: SolidModel.h:87
const bool _has_temp
Definition: SolidModel.h:91
SymmTensor _d_strain_dT
Definition: SolidModel.h:137
std::set< MooseSharedPointer< ConstitutiveModel > > _models_to_free
Definition: SolidModel.h:245
bool _bulk_modulus_set
Definition: SolidModel.h:62
Moose::CoordinateSystemType _coord_type
Definition: SolidModel.h:58
std::map< SubdomainID, MooseSharedPointer< ConstitutiveModel > > _constitutive_model
Definition: SolidModel.h:243
Real _youngs_modulus
Definition: SolidModel.h:72
MaterialProperty< SymmTensor > & _elastic_strain
Definition: SolidModel.h:117
MaterialProperty< RealVectorValue > * _crack_flags
Definition: SolidModel.h:120
ColumnMajorMatrix _principal_strain
Definition: SolidModel.h:131
Real JxW(unsigned i) const
Definition: SolidModel.h:55
SymmTensor _strain_increment
Definition: SolidModel.h:143
SymmTensor _stress_old
Definition: SolidModel.h:112
const MaterialProperty< ColumnMajorMatrix > * _crack_rotation_old
Definition: SolidModel.h:126
const MaterialProperty< RealVectorValue > * _crack_max_strain_old
Definition: SolidModel.h:130
const VariableValue & _temperature
Definition: SolidModel.h:92
const MaterialProperty< SymmTensor > & _elastic_strain_old
Definition: SolidModel.h:118
MaterialProperty< RealVectorValue > * _crack_max_strain
Definition: SolidModel.h:129
const MaterialProperty< RealVectorValue > * _crack_strain_old
Definition: SolidModel.h:128
bool _has_stress_free_temp
Definition: SolidModel.h:98
SolidMechanics::Element * _element
Definition: SolidModel.h:263
MaterialProperty< RealVectorValue > * _crack_strain
Definition: SolidModel.h:127
MaterialProperty< SymmTensor > & _stress
Definition: SolidModel.h:106
MaterialProperty< Real > * _current_instantaneous_thermal_expansion_coef
Definition: SolidModel.h:156
Function * _alpha_function
Definition: SolidModel.h:96
Real _ref_temp
Definition: SolidModel.h:101
bool _lambda_set
Definition: SolidModel.h:63
Real _lambda
Definition: SolidModel.h:69
virtual void computeStress()
Compute the stress (sigma += deltaSigma)
Definition: SolidModel.h:183
std::map< SubdomainID, std::vector< MooseSharedPointer< VolumetricModel > > > _volumetric_models
Definition: SolidModel.h:103
std::vector< unsigned int > _active_crack_planes
Definition: SolidModel.h:85
Real _shear_modulus
Definition: SolidModel.h:71
const SolidMechanics::Element * element() const
Definition: SolidModel.h:219
int delta(int i, int j) const
Definition: SolidModel.h:221
MaterialProperty< RankTwoTensor > * _Eshelby_tensor
Definition: SolidModel.h:152
RealVectorValue _crack_flags_local
Definition: SolidModel.h:122
const Real _cracking_beta
Definition: SolidModel.h:80