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.
small_strain_isotropic_damage_implex_3d.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: Marcelo Raschi
12 
13 #pragma once
14 
15 // System includes
16 
17 // External includes
18 
19 // Project includes
20 #include "custom_constitutive/elastic_isotropic_3d.h"
21 
22 namespace Kratos
23 {
26 
30 
34 
38 
42 
60 class KRATOS_API(CONSTITUTIVE_LAWS_APPLICATION) SmallStrainIsotropicDamageImplex3D
61  : public ElasticIsotropic3D
62 {
63 public:
64 
68  typedef std::size_t SizeType;
69 
70  // Counted pointer
72 
76 
81 
86 
91 
96  ConstitutiveLaw::Pointer Clone() const override;
97 
102 
105 
110  void GetLawFeatures(Features& rFeatures) override;
111 
117  bool Has(const Variable<double>& rThisVariable) override;
118 
124  bool Has(const Variable<Vector>& rThisVariable) override;
125 
132  Vector& GetValue(
133  const Variable<Vector>& rThisVariable,
134  Vector& rValue
135  ) override;
136 
143  void SetValue(
144  const Variable<Vector>& rThisVariable,
145  const Vector& rValue,
146  const ProcessInfo& rProcessInfo
147  ) override;
148 
156  void InitializeMaterial(const Properties& rMaterialProperties,
157  const GeometryType& rElementGeometry,
158  const Vector& rShapeFunctionsValues) override;
159 
165  void CalculateStressResponse(ConstitutiveLaw::Parameters& rValues,
166  Vector& rInternalVariables) override;
167 
173  void CalculateMaterialResponsePK2(Parameters& rValues) override;
174 
180  {
181  return false;
182  }
183 
189  {
190  return true;
191  }
192 
198  void FinalizeMaterialResponseCauchy(Parameters& rValues) override;
199 
207  double& CalculateValue(Parameters& rValues,
208  const Variable<double>& rThisVariable,
209  double& rValue) override;
210 
218  Vector& CalculateValue(Parameters& rValues,
219  const Variable<Vector>& rThisVariable,
220  Vector& rValue) override;
221 
230  int Check(
231  const Properties& rMaterialProperties,
232  const GeometryType& rElementGeometry,
233  const ProcessInfo& rCurrentProcessInfo
234  ) const override;
235 
242 
244  void PrintData(std::ostream& rOStream) const override {
245  rOStream << "Small Strain Isotropic Damage 3D constitutive law with Implex intergration scheme\n";
246  };
247 
248 protected:
249 
253 
259 
263 
266 
272  virtual void ComputePositiveStressVector(
273  Vector& rStressVectorPos,
274  Vector& rStressVector);
275 
281  double EvaluateHardeningModulus(
282  double StrainVariable,
283  const Properties &rMaterialProperties);
284 
290  double EvaluateHardeningLaw(
291  double StrainVariable,
292  const Properties &rMaterialProperties);
293 
295 
296 private:
297 
300 
302 
305 
307 
310 
312 
315 
317 
321 
324 
325  friend class Serializer;
326 
327  void save(Serializer& rSerializer) const override;
328 
329  void load(Serializer& rSerializer) override;
330 
331 }; // class SmallStrainIsotropicDamage3D
332 } // namespace Kratos
Definition: elastic_isotropic_3d.h:53
Geometry base class.
Definition: geometry.h:71
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
Damage with hardening constitutive law in 3D, using Implex integration scheme (see J Oliver et al,...
Definition: small_strain_isotropic_damage_implex_3d.h:62
ProcessInfo ProcessInfoType
Definition: small_strain_isotropic_damage_implex_3d.h:67
SmallStrainIsotropicDamageImplex3D(const SmallStrainIsotropicDamageImplex3D &rOther)
Copy constructor.
bool RequiresFinalizeMaterialResponse() override
Indicates if this CL requires finalization of the material response, called by the element in Finaliz...
Definition: small_strain_isotropic_damage_implex_3d.h:188
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: small_strain_isotropic_damage_implex_3d.h:244
bool RequiresInitializeMaterialResponse() override
Indicates if this CL requires initialization of the material response, called by the element in Initi...
Definition: small_strain_isotropic_damage_implex_3d.h:179
double mStrainVariablePrevious
Definition: small_strain_isotropic_damage_implex_3d.h:257
~SmallStrainIsotropicDamageImplex3D() override
Destructor.
double mStrainVariable
Definition: small_strain_isotropic_damage_implex_3d.h:256
KRATOS_CLASS_POINTER_DEFINITION(SmallStrainIsotropicDamageImplex3D)
std::size_t SizeType
Definition: small_strain_isotropic_damage_implex_3d.h:68
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
def SetValue(entity, variable, value)
Definition: coupling_interface_data.py:256
def load(f)
Definition: ode_solve.py:307
Definition: constitutive_law.h:137
Definition: constitutive_law.h:189