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_kinematic_plasticity.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 // Collaborator:
13 //
14 
15 #pragma once
16 
17 // System includes
18 
19 // External includes
20 
21 // Project includes
22 #include "custom_constitutive/elastic_isotropic_3d.h"
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 
57 template <class TConstLawIntegratorType>
58 class KRATOS_API(CONSTITUTIVE_LAWS_APPLICATION) GenericSmallStrainKinematicPlasticity
59  : public std::conditional<TConstLawIntegratorType::VoigtSize == 6, ElasticIsotropic3D, LinearPlaneStrain >::type
60 {
61 public:
64 
66  static constexpr SizeType Dimension = TConstLawIntegratorType::Dimension;
67 
69  static constexpr SizeType VoigtSize = TConstLawIntegratorType::VoigtSize;
70 
73 
76 
79 
82 
84  typedef Node NodeType;
85 
88 
92 
97  {
98  }
99 
103  ConstitutiveLaw::Pointer Clone() const override
104  {
105  return Kratos::make_shared<GenericSmallStrainKinematicPlasticity<TConstLawIntegratorType>>(*this);
106  }
107 
112  : BaseType(rOther),
113  mPlasticDissipation(rOther.mPlasticDissipation),
114  mThreshold(rOther.mThreshold),
115  mPlasticStrain(rOther.mPlasticStrain),
116  mPreviousStressVector(rOther.mPreviousStressVector),
117  mBackStressVector(rOther.mBackStressVector)
118  {
119  }
120 
125  {
126  }
127 
131 
135 
140  void CalculateMaterialResponsePK1(ConstitutiveLaw::Parameters &rValues) override;
141 
146  void CalculateMaterialResponsePK2(ConstitutiveLaw::Parameters &rValues) override;
147 
152  void CalculateMaterialResponseKirchhoff(ConstitutiveLaw::Parameters &rValues) override;
153 
158  void CalculateMaterialResponseCauchy(ConstitutiveLaw::Parameters &rValues) override;
159 
167  void InitializeMaterial(
168  const Properties& rMaterialProperties,
169  const GeometryType& rElementGeometry,
170  const Vector& rShapeFunctionsValues
171  ) override;
172 
177  void FinalizeMaterialResponsePK1(ConstitutiveLaw::Parameters &rValues) override;
178 
183  void FinalizeMaterialResponsePK2(ConstitutiveLaw::Parameters &rValues) override;
184 
189  void FinalizeMaterialResponseKirchhoff(ConstitutiveLaw::Parameters &rValues) override;
194  void FinalizeMaterialResponseCauchy(ConstitutiveLaw::Parameters &rValues) override;
195 
201  bool Has(const Variable<double> &rThisVariable) override;
202 
208  bool Has(const Variable<Vector> &rThisVariable) override;
209 
215  bool Has(const Variable<Matrix> &rThisVariable) override;
216 
223  void SetValue(
224  const Variable<double> &rThisVariable,
225  const double& rValue,
226  const ProcessInfo& rCurrentProcessInfo
227  ) override;
228 
235  void SetValue(
236  const Variable<Vector> &rThisVariable,
237  const Vector& rValue,
238  const ProcessInfo& rCurrentProcessInfo
239  ) override;
240 
247  double& GetValue(
248  const Variable<double> &rThisVariable,
249  double& rValue
250  ) override;
251 
258  Vector& GetValue(
259  const Variable<Vector> &rThisVariable,
260  Vector& rValue
261  ) override;
262 
269  Matrix& GetValue(
270  const Variable<Matrix>& rThisVariable,
271  Matrix& rValue
272  ) override;
273 
278  {
279  return true;
280  }
281 
286  {
287  return false;
288  }
289 
297  double& CalculateValue(
298  ConstitutiveLaw::Parameters& rParameterValues,
299  const Variable<double>& rThisVariable,
300  double& rValue) override;
301 
309  Vector& CalculateValue(
310  ConstitutiveLaw::Parameters& rParameterValues,
311  const Variable<Vector>& rThisVariable,
312  Vector& rValue
313  ) override;
314 
322  Matrix& CalculateValue(
323  ConstitutiveLaw::Parameters& rParameterValues,
324  const Variable<Matrix>& rThisVariable,
325  Matrix& rValue
326  ) override;
327 
337  int Check(
338  const Properties& rMaterialProperties,
339  const GeometryType& rElementGeometry,
340  const ProcessInfo& rCurrentProcessInfo
341  ) const override;
342 
346 
350 
354 
358 
360 
361 protected:
364 
368 
372 
373  double& GetThreshold() { return mThreshold; }
374  double& GetPlasticDissipation() { return mPlasticDissipation; }
375  Vector& GetPlasticStrain() { return mPlasticStrain; }
376 
377  void SetThreshold(const double Threshold) { mThreshold = Threshold; }
378  void SetPlasticDissipation(const double PlasticDissipation) { mPlasticDissipation = PlasticDissipation; }
379  void SetPlasticStrain(const array_1d<double, VoigtSize>& rPlasticStrain) { mPlasticStrain = rPlasticStrain; }
380 
381  void SetPreviousStressVector(const Vector& toBS) {mPreviousStressVector = toBS; }
382  Vector& GetPreviousStressVector() { return mPreviousStressVector;}
383 
384  void SetBackStressVector(const Vector& toBS) {mBackStressVector = toBS; }
385  Vector& GetBackStressVector() { return mBackStressVector;}
386 
390 
394 
398 
402 
404  private:
407 
411 
412  // Converged values
413  double mPlasticDissipation = 0.0;
414  double mThreshold = 0.0;
415  Vector mPlasticStrain = ZeroVector(VoigtSize);
416 
417  // Kinematic variables
418  Vector mPreviousStressVector = ZeroVector(VoigtSize);
419  Vector mBackStressVector = ZeroVector(VoigtSize);
420 
424 
428 
433  void CalculateTangentTensor(ConstitutiveLaw::Parameters &rValues);
434 
438 
442 
446 
447  // Serialization
448 
449  friend class Serializer;
450 
451  void save(Serializer &rSerializer) const override
452  {
454  rSerializer.save("PlasticDissipation", mPlasticDissipation);
455  rSerializer.save("Threshold", mThreshold);
456  rSerializer.save("PlasticStrain", mPlasticStrain);
457  rSerializer.save("PreviousStressVector", mPreviousStressVector);
458  rSerializer.save("BackStressVector", mBackStressVector);
459  }
460 
461  void load(Serializer &rSerializer) override
462  {
464  rSerializer.load("PlasticDissipation", mPlasticDissipation);
465  rSerializer.load("Threshold", mThreshold);
466  rSerializer.load("PlasticStrain", mPlasticStrain);
467  rSerializer.load("PreviousStressVector", mPreviousStressVector);
468  rSerializer.load("BackStressVector", mBackStressVector);
469  }
470 
472 
473 }; // Class GenericYieldSurface
474 
475 } // namespace Kratos
Definition: constitutive_law.h:47
This class is the base class which define all the constitutive laws for kinematic plasticity in small...
Definition: generic_small_strain_kinematic_plasticity.h:60
Geometry< NodeType > GeometryType
The geometry definition.
Definition: generic_small_strain_kinematic_plasticity.h:87
Vector & GetPreviousStressVector()
Definition: generic_small_strain_kinematic_plasticity.h:382
void SetPlasticDissipation(const double PlasticDissipation)
Definition: generic_small_strain_kinematic_plasticity.h:378
std::conditional< VoigtSize==6, ElasticIsotropic3D, LinearPlaneStrain >::type BaseType
Definition of the base class.
Definition: generic_small_strain_kinematic_plasticity.h:72
void SetPlasticStrain(const array_1d< double, VoigtSize > &rPlasticStrain)
Definition: generic_small_strain_kinematic_plasticity.h:379
void SetThreshold(const double Threshold)
Definition: generic_small_strain_kinematic_plasticity.h:377
void SetPreviousStressVector(const Vector &toBS)
Definition: generic_small_strain_kinematic_plasticity.h:381
double & GetThreshold()
Definition: generic_small_strain_kinematic_plasticity.h:373
Node NodeType
The node definition.
Definition: generic_small_strain_kinematic_plasticity.h:84
~GenericSmallStrainKinematicPlasticity() override
Definition: generic_small_strain_kinematic_plasticity.h:124
double & GetPlasticDissipation()
Definition: generic_small_strain_kinematic_plasticity.h:374
bool RequiresFinalizeMaterialResponse() override
If the CL requires to initialize the material response, called by the element in InitializeSolutionSt...
Definition: generic_small_strain_kinematic_plasticity.h:277
void SetBackStressVector(const Vector &toBS)
Definition: generic_small_strain_kinematic_plasticity.h:384
GenericSmallStrainKinematicPlasticity()
Definition: generic_small_strain_kinematic_plasticity.h:96
BoundedMatrix< double, Dimension, Dimension > BoundedMatrixType
The definition of the bounded matrix type.
Definition: generic_small_strain_kinematic_plasticity.h:78
GenericSmallStrainKinematicPlasticity(const GenericSmallStrainKinematicPlasticity &rOther)
Definition: generic_small_strain_kinematic_plasticity.h:111
bool RequiresInitializeMaterialResponse() override
If the CL requires to initialize the material response, called by the element in InitializeSolutionSt...
Definition: generic_small_strain_kinematic_plasticity.h:285
Vector & GetBackStressVector()
Definition: generic_small_strain_kinematic_plasticity.h:385
ConstitutiveLaw::Pointer Clone() const override
Definition: generic_small_strain_kinematic_plasticity.h:103
KRATOS_CLASS_POINTER_DEFINITION(GenericSmallStrainKinematicPlasticity)
Counted pointer of GenericSmallStrainKinematicPlasticity.
array_1d< double, VoigtSize > BoundedArrayType
The definition of the Voigt array type.
Definition: generic_small_strain_kinematic_plasticity.h:75
Vector & GetPlasticStrain()
Definition: generic_small_strain_kinematic_plasticity.h:375
Geometry base class.
Definition: geometry.h:71
Definition: amatrix_interface.h:41
This class defines the node.
Definition: node.h:65
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
KratosZeroVector< double > ZeroVector
Definition: amatrix_interface.h:561
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