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.
viscous_generalized_maxwell.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 & Lucia Barbu
12 // Collaborator: Vicente Mataix Ferrandiz
13 //
14 
15 #pragma once
16 
17 // System includes
18 
19 // External includes
20 
21 // Project includes
23 #include "custom_constitutive/elastic_isotropic_3d.h"
24 
25 namespace Kratos
26 {
29 
33 
37 
41 
45 
59 template<class TElasticBehaviourLaw>
60 class KRATOS_API(CONSTITUTIVE_LAWS_APPLICATION) ViscousGeneralizedMaxwell
61  : public TElasticBehaviourLaw
62 {
63  public:
66 
69 
71  typedef TElasticBehaviourLaw BaseType;
72 
74  typedef std::size_t IndexType;
75 
77  typedef std::size_t SizeType;
78 
80  static constexpr SizeType Dimension = TElasticBehaviourLaw::Dimension;
81 
83  static constexpr SizeType VoigtSize = TElasticBehaviourLaw::VoigtSize;
84 
87 
89  typedef Node NodeType;
90 
93 
95  static constexpr double tolerance = std::numeric_limits<double>::epsilon();
96 
100 
105 
109  ConstitutiveLaw::Pointer Clone() const override;
110 
115 
119  ~ViscousGeneralizedMaxwell() override;
120 
124 
128 
133  void CalculateMaterialResponsePK1(ConstitutiveLaw::Parameters& rValues) override;
134 
139  void CalculateMaterialResponsePK2(ConstitutiveLaw::Parameters& rValues) override;
140 
145  void CalculateMaterialResponseKirchhoff(ConstitutiveLaw::Parameters& rValues) override;
146 
151  void CalculateMaterialResponseCauchy(ConstitutiveLaw::Parameters& rValues) override;
152 
157  void FinalizeMaterialResponsePK1(ConstitutiveLaw::Parameters& rValues) override;
158 
163  void FinalizeMaterialResponsePK2(ConstitutiveLaw::Parameters& rValues) override;
164 
169  void FinalizeMaterialResponseKirchhoff(ConstitutiveLaw::Parameters& rValues) override;
170 
175  void FinalizeMaterialResponseCauchy(ConstitutiveLaw::Parameters& rValues) override;
176 
184  Matrix& CalculateValue(
185  ConstitutiveLaw::Parameters& rParameterValues,
186  const Variable<Matrix>& rThisVariable,
187  Matrix& rValue
188  ) override;
189 
199  int Check(
200  const Properties& rMaterialProperties,
201  const GeometryType& rElementGeometry,
202  const ProcessInfo& rCurrentProcessInfo
203  ) const override;
204 
209  {
210  return false;
211  }
212 
217  {
218  return true;
219  }
220 
224 
228 
232 
236 
238 
239  protected:
242 
246 
250 
254 
258 
262 
266 
268  private:
271 
275 
276  // Converged values
277  Vector mPrevStressVector = ZeroVector(VoigtSize);
278  Vector mPrevStrainVector = ZeroVector(VoigtSize);
279 
283 
287 
288  Vector& GetPreviousStressVector() { return mPrevStressVector; }
289  void SetPreviousStressVector(const Vector& PrevStressVector) { mPrevStressVector = PrevStressVector; }
290 
291  Vector& GetPreviousStrainVector() { return mPrevStrainVector; }
292  void SetPreviousStrainVector(const Vector& PrevStrainVector) { mPrevStrainVector = PrevStrainVector; }
293 
297  void ComputeViscoElasticity(ConstitutiveLaw::Parameters& rValues);
298 
302 
306 
310 
311  // Serialization
312 
313  friend class Serializer;
314 
315  void save(Serializer &rSerializer) const override
316  {
318  rSerializer.save("PrevStressVector", mPrevStressVector);
319  rSerializer.save("PrevStrainVector", mPrevStrainVector);
320  }
321 
322  void load(Serializer &rSerializer) override
323  {
325  rSerializer.load("PrevStressVector", mPrevStressVector);
326  rSerializer.load("PrevStrainVector", mPrevStrainVector);
327  }
328 
330 
331 }; // Class GenericYieldSurface
332 
333 } // namespace Kratos
Definition: constitutive_law.h:47
Geometry base class.
Definition: geometry.h:71
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
Variable class contains all information needed to store and retrive data from a data container.
Definition: variable.h:63
This is a viscous law using Maxwell formulation.
Definition: viscous_generalized_maxwell.h:62
Geometry< NodeType > GeometryType
The geometry definition.
Definition: viscous_generalized_maxwell.h:92
Node NodeType
The node definition.
Definition: viscous_generalized_maxwell.h:89
KRATOS_CLASS_POINTER_DEFINITION(ViscousGeneralizedMaxwell)
Counted pointer of GenericYieldSurface.
bool RequiresFinalizeMaterialResponse() override
If the CL requires to initialize the material response, called by the element in InitializeSolutionSt...
Definition: viscous_generalized_maxwell.h:216
std::size_t IndexType
The index definition.
Definition: viscous_generalized_maxwell.h:74
ConstitutiveLaw CLBaseType
Definition of the base CL class.
Definition: viscous_generalized_maxwell.h:68
TElasticBehaviourLaw BaseType
Definition of the base class.
Definition: viscous_generalized_maxwell.h:71
std::size_t SizeType
The size definition.
Definition: viscous_generalized_maxwell.h:77
bool RequiresInitializeMaterialResponse() override
If the CL requires to initialize the material response, called by the element in InitializeSolutionSt...
Definition: viscous_generalized_maxwell.h:208
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
KratosZeroVector< double > ZeroVector
Definition: amatrix_interface.h:561
Internals::Matrix< double, AMatrix::dynamic, 1 > Vector
Definition: amatrix_interface.h:472
def load(f)
Definition: ode_solve.py:307
Definition: constitutive_law.h:189