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_j2_plasticity_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 // Manuel Caicedo
13 // Alfredo Huespe
14 // Collaborator: Vicente Mataix Ferrandiz
15 
16 #pragma once
17 
18 // System includes
19 
20 // External includes
21 
22 // Project includes
24 
25 namespace Kratos
26 {
29 
33 
37 
41 
45 
63 class KRATOS_API(CONSTITUTIVE_LAWS_APPLICATION) SmallStrainJ2Plasticity3D
64  : public ConstitutiveLaw
65 {
66 public:
67 
70 
73  typedef std::size_t SizeType;
74 
75  // Counted pointer of SmallStrainJ2Plasticity3D
77 
81 
86 
91 
95  ~SmallStrainJ2Plasticity3D() override;
96 
101  ConstitutiveLaw::Pointer Clone() const override;
102 
106 
110 
115  void GetLawFeatures(Features& rFeatures) override;
116 
121  {
122  return 3;
123  };
124 
128  SizeType GetStrainSize() const override
129  {
130  return 6;
131  };
132 
138  bool Has(const Variable<double>& rThisVariable) override;
139 
146  void SetValue(
147  const Variable<double>& rThisVariable,
148  const double& rValue,
149  const ProcessInfo& rCurrentProcessInfo
150  ) override;
151 
158  double& GetValue(
159  const Variable<double>& rThisVariable,
160  double& rValue
161  ) override;
162 
170  void InitializeMaterial(const Properties& rMaterialProperties,
171  const GeometryType& rElementGeometry,
172  const Vector& rShapeFunctionsValues) override;
173 
179  void CalculateMaterialResponsePK1(ConstitutiveLaw::Parameters& rValues) override;
180 
186  void CalculateMaterialResponsePK2(ConstitutiveLaw::Parameters& rValues) override;
187 
193  void CalculateMaterialResponseKirchhoff(ConstitutiveLaw::Parameters& rValues) override;
194 
200  void CalculateMaterialResponseCauchy(ConstitutiveLaw::Parameters& rValues) override;
201 
207  void InitializeMaterialResponsePK1(ConstitutiveLaw::Parameters& rValues) override;
208 
214  void InitializeMaterialResponsePK2(ConstitutiveLaw::Parameters& rValues) override;
215 
221  void InitializeMaterialResponseKirchhoff(ConstitutiveLaw::Parameters& rValues) override;
222 
228  void InitializeMaterialResponseCauchy(ConstitutiveLaw::Parameters& rValues) override;
229 
235  void FinalizeMaterialResponsePK1(ConstitutiveLaw::Parameters& rValues) override;
236 
242  void FinalizeMaterialResponsePK2(ConstitutiveLaw::Parameters& rValues) override;
243 
249  void FinalizeMaterialResponseKirchhoff(Parameters& rValues) override;
250 
256  void FinalizeMaterialResponseCauchy(ConstitutiveLaw::Parameters& rValues) override;
257 
265  double& CalculateValue(ConstitutiveLaw::Parameters& rParameterValues,
266  const Variable<double>& rThisVariable,
267  double& rValue) override;
268 
276  Vector& CalculateValue(ConstitutiveLaw::Parameters& rParameterValues,
277  const Variable<Vector>& rThisVariable,
278  Vector& rValue) override;
279 
288  int Check(
289  const Properties& rMaterialProperties,
290  const GeometryType& rElementGeometry,
291  const ProcessInfo& rCurrentProcessInfo
292  ) const override;
293 
300 
302  void PrintData(std::ostream& rOStream) const override {
303  rOStream << "Linear J2 Plasticity 3D constitutive law\n";
304  };
305 
306 protected:
307 
310 
314 
317 
321 
325 
333  virtual void CalculateStressResponse(ConstitutiveLaw::Parameters& rValues,
334  Vector& rPlasticStrain,
335  double& rAccumulatedPlasticStrain );
336 
343  double YieldFunction(
344  const double NormDeviationStress,
345  const Properties& rMaterialProperties,
346  const double AccumulatedPlasticStrain
347  );
348 
355  double GetAccumPlasticStrainRate(const double NormStressTrial, const Properties &rMaterialProperties,
356  const double AccumulatedPlasticStrainOld);
357 
363  double GetSaturationHardening(const Properties& rMaterialProperties, const double);
364 
370  double GetPlasticPotential(const Properties& rMaterialProperties,
371  const double accumulated_plastic_strain);
372 
381  virtual void CalculateTangentMatrix(const double DeltaGamma, const double NormStressTrial,
382  const Vector &rYFNormalVector,
383  const Properties &rMaterialProperties,
384  const double AccumulatedPlasticStrain, Matrix &rTMatrix);
385 
391  virtual void CalculateElasticMatrix(const Properties &rMaterialProperties, Matrix &rElasticMatrix);
392 
396 
397 
401 
402 
407 
408 private:
409 
412 
416 
420 
424 
429 
433 
434  friend class Serializer;
435 
436  void save(Serializer& rSerializer) const override;
437 
438  void load(Serializer& rSerializer) override;
439 
440 }; // Class SmallStrainJ2Plasticity3D
441 } // namespace Kratos.
Definition: constitutive_law.h:47
virtual void CalculateStressResponse(Parameters &rValues, Vector &rInternalVariables)
Definition: constitutive_law.cpp:698
std::size_t SizeType
Definition: constitutive_law.h:82
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
Definition: small_strain_j2_plasticity_3d.h:65
std::size_t SizeType
Definition: small_strain_j2_plasticity_3d.h:73
ConstitutiveLaw BaseType
Definition: small_strain_j2_plasticity_3d.h:72
Vector mPlasticStrain
Definition: small_strain_j2_plasticity_3d.h:315
ProcessInfo ProcessInfoType
Definition: small_strain_j2_plasticity_3d.h:71
SizeType GetStrainSize() const override
Voigt tensor size:
Definition: small_strain_j2_plasticity_3d.h:128
double mAccumulatedPlasticStrain
The previous plastic strain (one for each of the strain components)
Definition: small_strain_j2_plasticity_3d.h:316
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: small_strain_j2_plasticity_3d.h:302
KRATOS_CLASS_POINTER_DEFINITION(SmallStrainJ2Plasticity3D)
SizeType WorkingSpaceDimension() override
dimension of the constitutive law
Definition: small_strain_j2_plasticity_3d.h:120
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