KratosMultiphysics
KRATOS Multiphysics (Kratos) is a framework for building parallel, multi-disciplinary simulation software, aiming at modularity, extensibility, and high performance. Kratos is written in C++, and counts with an extensive Python interface.
generic_small_strain_d_plus_d_minus_damage.h
Go to the documentation of this file.
1 // KRATOS ___ _ _ _ _ _ __ _
2 // / __\___ _ __ ___| |_(_) |_ _ _| |_(_)_ _____ / / __ ___ _____ /_\ _ __ _ __
3 // / / / _ \| '_ \/ __| __| | __| | | | __| \ \ / / _ \/ / / _` \ \ /\ / / __| //_\\| '_ \| '_ |
4 // / /__| (_) | | | \__ \ |_| | |_| |_| | |_| |\ V / __/ /__| (_| |\ V V /\__ \/ _ \ |_) | |_) |
5 // \____/\___/|_| |_|___/\__|_|\__|\__,_|\__|_| \_/ \___\____/\__,_| \_/\_/ |___/\_/ \_/ .__/| .__/
6 // |_| |_|
7 //
8 // License: BSD License
9 // license: structural_mechanics_application/license.txt
10 //
11 // Main authors: Alejandro Cornejo
12 //
13 //
14 
15 #pragma once
16 
17 
18 // System includes
19 
20 // External includes
21 
22 // Project includes
23 #include "custom_constitutive/linear_plane_strain.h"
24 
25 namespace Kratos
26 {
29 
33 
34  // The size type definition
35  typedef std::size_t SizeType;
36 
40 
44 
48 
56 template <class TConstLawIntegratorTensionType, class TConstLawIntegratorCompressionType>
57 class KRATOS_API(CONSTITUTIVE_LAWS_APPLICATION) GenericSmallStrainDplusDminusDamage
58  : public std::conditional<TConstLawIntegratorTensionType::VoigtSize == 6, ElasticIsotropic3D, LinearPlaneStrain >::type
59 {
60 public:
63 
65  static constexpr SizeType Dimension = TConstLawIntegratorTensionType::Dimension;
66 
68  static constexpr SizeType VoigtSize = TConstLawIntegratorTensionType::VoigtSize;
69 
72 
75 
77  typedef Node NodeType;
78 
81 
83  static constexpr double tolerance = std::numeric_limits<double>::epsilon();
84 
86  double DamageTension = 0.0;
87  double DamageCompression = 0.0;
88  double ThresholdTension = 0.0;
89  double ThresholdCompression = 0.0;
92  double UniaxialTensionStress = 0.0;
93  double UniaxialCompressionStress = 0.0;
94  };
98 
103  {
104  }
105 
109  ConstitutiveLaw::Pointer Clone() const override
110  {
111  return Kratos::make_shared<GenericSmallStrainDplusDminusDamage<TConstLawIntegratorTensionType, TConstLawIntegratorCompressionType>>(*this);
112  }
113 
118  : BaseType(rOther),
119  mTensionDamage(rOther.mTensionDamage),
120  mTensionThreshold(rOther.mTensionThreshold),
121  mNonConvTensionDamage(rOther.mNonConvTensionDamage),
122  mNonConvTensionThreshold(rOther.mNonConvTensionThreshold),
123  mCompressionDamage(rOther.mCompressionDamage),
124  mCompressionThreshold(rOther.mCompressionThreshold),
125  mNonConvCompressionDamage(rOther.mNonConvCompressionDamage),
126  mNonConvCompressionThreshold(rOther.mNonConvCompressionThreshold)
127  {
128  }
129 
134  {
135  }
136 
140 
144 
149  void CalculateMaterialResponsePK1(ConstitutiveLaw::Parameters &rValues) override;
150 
155  void CalculateMaterialResponsePK2(ConstitutiveLaw::Parameters &rValues) override;
156 
161  void CalculateMaterialResponseKirchhoff(ConstitutiveLaw::Parameters &rValues) override;
162 
167  void CalculateMaterialResponseCauchy(ConstitutiveLaw::Parameters &rValues) override;
168 
173  bool IntegrateStressTensionIfNecessary(
174  const double F_tension,
175  DamageParameters& Parameters,
176  array_1d<double, VoigtSize>& IntegratedStressVectorTension,
178  );
179 
184  bool IntegrateStressCompressionIfNecessary(
185  const double F_compression,
186  DamageParameters& Parameters,
187  array_1d<double, VoigtSize>& IntegratedStressVectorCompression,
189  );
190 
194  void CalculateIntegratedStressVector(
195  Vector& IntegratedStressVectorTension,
196  const DamageParameters& Parameters,
198  );
199 
207  void InitializeMaterial(
208  const Properties& rMaterialProperties,
209  const GeometryType& rElementGeometry,
210  const Vector& rShapeFunctionsValues
211  ) override;
212 
221  void FinalizeSolutionStep(
222  const Properties &rMaterialProperties,
223  const GeometryType &rElementGeometry,
224  const Vector& rShapeFunctionsValues,
225  const ProcessInfo& rCurrentProcessInfo
226  ) override;
227 
232  void FinalizeMaterialResponsePK1(ConstitutiveLaw::Parameters &rValues) override;
233 
238  void FinalizeMaterialResponsePK2(ConstitutiveLaw::Parameters &rValues) override;
239 
244  void FinalizeMaterialResponseKirchhoff(ConstitutiveLaw::Parameters &rValues) override;
249  void FinalizeMaterialResponseCauchy(ConstitutiveLaw::Parameters &rValues) override;
250 
256  bool Has(const Variable<double> &rThisVariable) override;
257 
263  bool Has(const Variable<Vector> &rThisVariable) override;
264 
270  bool Has(const Variable<Matrix> &rThisVariable) override;
271 
278  void SetValue(
279  const Variable<double> &rThisVariable,
280  const double& rValue,
281  const ProcessInfo& rCurrentProcessInfo
282  ) override;
283 
290  double& GetValue(
291  const Variable<double> &rThisVariable,
292  double& rValue
293  ) override;
294 
301  Vector& GetValue(
302  const Variable<Vector> &rThisVariable,
303  Vector& rValue
304  ) override;
305 
312  Matrix& GetValue(
313  const Variable<Matrix>& rThisVariable,
314  Matrix& rValue
315  ) override;
316 
321  {
322  return true;
323  }
324 
329  {
330  return false;
331  }
332 
340  double& CalculateValue(
341  ConstitutiveLaw::Parameters& rParameterValues,
342  const Variable<double>& rThisVariable,
343  double& rValue) override;
344 
352  Vector& CalculateValue(
353  ConstitutiveLaw::Parameters& rParameterValues,
354  const Variable<Vector>& rThisVariable,
355  Vector& rValue
356  ) override;
357 
365  Matrix& CalculateValue(
366  ConstitutiveLaw::Parameters& rParameterValues,
367  const Variable<Matrix>& rThisVariable,
368  Matrix& rValue
369  ) override;
370 
380  int Check(
381  const Properties& rMaterialProperties,
382  const GeometryType& rElementGeometry,
383  const ProcessInfo& rCurrentProcessInfo
384  ) const override;
385 
389 
393 
397 
401 
403 
404 protected:
407 
411 
415 
419  // Tension values
420  double& GetTensionThreshold() { return mTensionThreshold; }
421  double& GetTensionDamage() { return mTensionDamage; }
422  double& GetNonConvTensionThreshold() { return mNonConvTensionThreshold; }
423  double& GetNonConvTensionDamage() { return mNonConvTensionDamage; }
424 
425  void SetTensionThreshold(const double toThreshold) { mTensionThreshold = toThreshold; }
426  void SetTensionDamage(const double toDamage) { mTensionDamage = toDamage; }
427  void SetNonConvTensionThreshold(const double toThreshold) { mNonConvTensionThreshold = toThreshold; }
428  void SetNonConvTensionDamage(const double toDamage) { mNonConvTensionDamage = toDamage; }
429 
430  // Compression values
431  double& GetCompressionThreshold() { return mCompressionThreshold; }
432  double& GetCompressionDamage() { return mCompressionDamage; }
433  double& GetNonConvCompressionThreshold() { return mNonConvCompressionThreshold; }
434  double& GetNonConvCompressionDamage() { return mNonConvCompressionDamage; }
435 
436  void SetCompressionThreshold(const double toThreshold) { mCompressionThreshold = toThreshold; }
437  void SetCompressionDamage(const double toDamage) { mCompressionDamage = toDamage; }
438  void SetNonConvCompressionThreshold(const double toThreshold) { mNonConvCompressionThreshold = toThreshold; }
439  void SetNonConvCompressionDamage(const double toDamage) { mNonConvCompressionDamage = toDamage; }
440 
441  void SetTensionStress(const double toS){mTensionUniaxialStress = toS;}
442  void SetCompressionStress(const double toS){mCompressionUniaxialStress = toS;}
443 
447 
451 
455 
457 private:
460 
464 
465  // Converged values
466  double mTensionDamage = 0.0;
467  double mTensionThreshold = 0.0;
468 
469  // Non Converged values
470  double mNonConvTensionDamage = 0.0;
471  double mNonConvTensionThreshold = 0.0;
472 
473  double mCompressionDamage = 0.0;
474  double mCompressionThreshold = 0.0;
475  // double mUniaxialStress = 0.0;
476 
477  // Non Converged values
478  double mNonConvCompressionDamage = 0.0;
479  double mNonConvCompressionThreshold = 0.0;
480 
481  double mTensionUniaxialStress = 0.0;
482  double mCompressionUniaxialStress = 0.0;
486 
490 
495  void CalculateTangentTensor(ConstitutiveLaw::Parameters &rValues);
496 
501  void CalculateSecantTensor(ConstitutiveLaw::Parameters& rValues, Matrix& rSecantTensor);
502 
506 
510 
514 
515  // Serialization
516 
517  friend class Serializer;
518 
519  void save(Serializer &rSerializer) const override
520  {
522  rSerializer.save("TensionDamage", mTensionDamage);
523  rSerializer.save("TensionThreshold", mTensionThreshold);
524  rSerializer.save("NonConvTensionDamage", mNonConvTensionDamage);
525  rSerializer.save("NonConvTensionThreshold", mNonConvTensionThreshold);
526  rSerializer.save("CompressionDamage", mCompressionDamage);
527  rSerializer.save("CompressionThreshold", mCompressionThreshold);
528  rSerializer.save("NonConvCompressionnDamage", mNonConvCompressionDamage);
529  rSerializer.save("NonConvCompressionThreshold", mNonConvCompressionThreshold);
530  }
531 
532  void load(Serializer &rSerializer) override
533  {
535  rSerializer.load("TensionDamage", mTensionDamage);
536  rSerializer.load("TensionThreshold", mTensionThreshold);
537  rSerializer.load("NonConvTensionDamage", mNonConvTensionDamage);
538  rSerializer.load("NonConvTensionThreshold", mNonConvTensionThreshold);
539  rSerializer.load("CompressionDamage", mCompressionDamage);
540  rSerializer.load("CompressionThreshold", mCompressionThreshold);
541  rSerializer.load("NonConvCompressionnDamage", mNonConvCompressionDamage);
542  rSerializer.load("NonConvCompressionThreshold", mNonConvCompressionThreshold);
543  }
544 
546 
547 }; // Class
548 
549 } // namespace Kratos
Definition: constitutive_law.h:47
This class is the base class which define all the constitutive laws for damage in small deformation.
Definition: generic_small_strain_d_plus_d_minus_damage.h:59
ConstitutiveLaw::Pointer Clone() const override
Definition: generic_small_strain_d_plus_d_minus_damage.h:109
double & GetCompressionThreshold()
Definition: generic_small_strain_d_plus_d_minus_damage.h:431
std::conditional< VoigtSize==6, ElasticIsotropic3D, LinearPlaneStrain >::type BaseType
Definition of the base class.
Definition: generic_small_strain_d_plus_d_minus_damage.h:71
Geometry< NodeType > GeometryType
The geometry definition.
Definition: generic_small_strain_d_plus_d_minus_damage.h:80
void SetCompressionDamage(const double toDamage)
Definition: generic_small_strain_d_plus_d_minus_damage.h:437
double & GetNonConvCompressionDamage()
Definition: generic_small_strain_d_plus_d_minus_damage.h:434
GenericSmallStrainDplusDminusDamage()
Definition: generic_small_strain_d_plus_d_minus_damage.h:102
void SetNonConvTensionThreshold(const double toThreshold)
Definition: generic_small_strain_d_plus_d_minus_damage.h:427
double & GetNonConvTensionDamage()
Definition: generic_small_strain_d_plus_d_minus_damage.h:423
double & GetTensionThreshold()
Definition: generic_small_strain_d_plus_d_minus_damage.h:420
bool RequiresFinalizeMaterialResponse() override
If the CL requires to initialize the material response, called by the element in InitializeSolutionSt...
Definition: generic_small_strain_d_plus_d_minus_damage.h:320
bool RequiresInitializeMaterialResponse() override
If the CL requires to initialize the material response, called by the element in InitializeSolutionSt...
Definition: generic_small_strain_d_plus_d_minus_damage.h:328
double & GetTensionDamage()
Definition: generic_small_strain_d_plus_d_minus_damage.h:421
void SetNonConvCompressionDamage(const double toDamage)
Definition: generic_small_strain_d_plus_d_minus_damage.h:439
void SetTensionThreshold(const double toThreshold)
Definition: generic_small_strain_d_plus_d_minus_damage.h:425
void SetNonConvCompressionThreshold(const double toThreshold)
Definition: generic_small_strain_d_plus_d_minus_damage.h:438
void SetCompressionThreshold(const double toThreshold)
Definition: generic_small_strain_d_plus_d_minus_damage.h:436
Node NodeType
The node definition.
Definition: generic_small_strain_d_plus_d_minus_damage.h:77
void SetCompressionStress(const double toS)
Definition: generic_small_strain_d_plus_d_minus_damage.h:442
KRATOS_CLASS_POINTER_DEFINITION(GenericSmallStrainDplusDminusDamage)
Counted pointer of GenericYieldSurface.
GenericSmallStrainDplusDminusDamage(const GenericSmallStrainDplusDminusDamage &rOther)
Definition: generic_small_strain_d_plus_d_minus_damage.h:117
~GenericSmallStrainDplusDminusDamage() override
Definition: generic_small_strain_d_plus_d_minus_damage.h:133
void SetTensionStress(const double toS)
Definition: generic_small_strain_d_plus_d_minus_damage.h:441
void SetTensionDamage(const double toDamage)
Definition: generic_small_strain_d_plus_d_minus_damage.h:426
void SetNonConvTensionDamage(const double toDamage)
Definition: generic_small_strain_d_plus_d_minus_damage.h:428
double & GetCompressionDamage()
Definition: generic_small_strain_d_plus_d_minus_damage.h:432
double & GetNonConvTensionThreshold()
Definition: generic_small_strain_d_plus_d_minus_damage.h:422
double & GetNonConvCompressionThreshold()
Definition: generic_small_strain_d_plus_d_minus_damage.h:433
Geometry base class.
Definition: geometry.h:71
This class defines the node.
Definition: node.h:65
This class provides to Kratos a data structure for I/O based on the standard of JSON.
Definition: kratos_parameters.h:59
ProcessInfo holds the current value of different solution parameters.
Definition: process_info.h:59
Properties encapsulates data shared by different Elements or Conditions. It can store any type of dat...
Definition: properties.h:69
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
void load(std::string const &rTag, TDataType &rObject)
Definition: serializer.h:207
void save(std::string const &rTag, std::array< TDataType, TDataSize > const &rObject)
Definition: serializer.h:545
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
bool Has(const std::string &ModelerName)
Checks if the modeler is registered.
Definition: modeler_factory.cpp:24
Parameters GetValue(Parameters &rParameters, const std::string &rEntry)
Definition: add_kratos_parameters_to_python.cpp:53
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
std::size_t SizeType
The definition of the size type.
Definition: mortar_classes.h:43
def SetValue(entity, variable, value)
Definition: coupling_interface_data.py:256
type
Definition: generate_gid_list_file.py:35
def load(f)
Definition: ode_solve.py:307
Definition: constitutive_law.h:189
Definition: generic_small_strain_d_plus_d_minus_damage.h:85
array_1d< double, VoigtSize > CompressionStressVector
Definition: generic_small_strain_d_plus_d_minus_damage.h:91
array_1d< double, VoigtSize > TensionStressVector
Definition: generic_small_strain_d_plus_d_minus_damage.h:90