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.
two_step_updated_lagrangian_V_P_implicit_solid_element.h
Go to the documentation of this file.
1 //
2 // Project Name: KratosFluidDynamicsApplication $
3 // Last modified by: $Author: AFranci $
4 // Date: $Date: February 2016 $
5 // Revision: $Revision: 0.0 $
6 //
7 //
8 
9 #if !defined(KRATOS_TWO_STEP_UPDATED_LAGRANGIAN_V_P_IMPLICIT_SOLID_ELEMENT_H_INCLUDED)
10 #define KRATOS_TWO_STEP_UPDATED_LAGRANGIAN_V_P_IMPLICIT_SOLID_ELEMENT_H_INCLUDED
11 
12 // System includes
13 #include <string>
14 #include <iostream>
15 
16 // External includes
17 
19 
20 namespace Kratos
21 {
22 
25 
28 
32 
36 
40 
44 
46 
48  /* class TwoStepUpdatedLagrangianVPImplicitSolidElement : public Element */
49  template <unsigned int TDim>
51  {
52 
53  public:
56 
59 
62 
64  typedef Node NodeType;
65 
68 
71 
73  typedef Vector VectorType;
74 
76  typedef Matrix MatrixType;
77 
78  typedef std::size_t IndexType;
79 
80  typedef std::size_t SizeType;
81 
82  typedef std::vector<std::size_t> EquationIdVectorType;
83 
84  typedef std::vector<Dof<double>::Pointer> DofsVectorType;
85 
87 
90 
93 
96 
98 
99  typedef typename BaseType::PropertiesType::Pointer pPropertiesType;
100 
102 
105 
107  typedef ConstitutiveLawType::Pointer ConstitutiveLawPointerType;
108 
112 
113  //Constructors.
114 
116 
120  {
121  }
122 
124 
129  {
130  }
131 
133 
137  TwoStepUpdatedLagrangianVPImplicitSolidElement(IndexType NewId, GeometryType::Pointer pGeometry) : BaseType(NewId, pGeometry)
138  {
139  }
140 
142 
147  TwoStepUpdatedLagrangianVPImplicitSolidElement(IndexType NewId, GeometryType::Pointer pGeometry, pPropertiesType pProperties) : BaseType(NewId, pGeometry, pProperties)
148  {
149  }
150 
152 
154  {
155  }
156 
159  {
160  }
161 
165 
169 
171 
178  Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes,
179  pPropertiesType pProperties) const override
180  {
181  return Element::Pointer(new TwoStepUpdatedLagrangianVPImplicitSolidElement(NewId, BaseType::GetGeometry().Create(ThisNodes), pProperties));
182  }
183 
184  Element::Pointer Clone(IndexType NewId, NodesArrayType const &ThisNodes) const override;
185 
186  void Initialize(const ProcessInfo &rCurrentProcessInfo) override;
187 
188  // /// Initializes the element and all geometric information required for the problem.
189  // void InitializeSolutionStep(const ProcessInfo &rCurrentProcessInfo) override;
190 
191  void InitializeNonLinearIteration(const ProcessInfo &rCurrentProcessInfo) override;
192 
193  void CalculateLeftHandSide(MatrixType &rLeftHandSideMatrix,
194  const ProcessInfo &rCurrentProcessInfo) override
195  {
196  KRATOS_TRY;
197  KRATOS_THROW_ERROR(std::logic_error, "TwoStepUpdatedLagrangianVPImplicitSolidElement::CalculateLeftHandSide not implemented", "");
198  KRATOS_CATCH("");
199  }
200 
201  void CalculateRightHandSide(VectorType &rRightHandSideVector,
202  const ProcessInfo &rCurrentProcessInfo) override
203  {
204  KRATOS_TRY;
205  KRATOS_THROW_ERROR(std::logic_error, "TwoStepUpdatedLagrangianVPImplicitSolidElement::CalculateRightHandSide not implemented", "");
206  KRATOS_CATCH("");
207  }
208 
209  // The following methods have different implementations depending on TDim
211 
219 
224  void UpdateCauchyStress(unsigned int g, const ProcessInfo &rCurrentProcessInfo) override;
225 
226  void InitializeElementalVariables(ElementalVariables &rElementalVariables) override;
227 
228  /* virtual void CalculateDeltaPosition (Matrix & rDeltaPosition); */
229 
233 
237 
239 
247  int Check(const ProcessInfo &rCurrentProcessInfo) const override;
248 
252 
256 
258  std::string Info() const override
259  {
260  std::stringstream buffer;
261  buffer << "TwoStepUpdatedLagrangianVPImplicitSolidElement #" << BaseType::Id();
262  return buffer.str();
263  }
264 
266  void PrintInfo(std::ostream &rOStream) const override
267  {
268  rOStream << "TwoStepUpdatedLagrangianVPImplicitSolidElement" << TDim << "D";
269  }
270 
271  // /// Print object's data.
272  // virtual void PrintData(std::ostream& rOStream) const;
273 
277 
279 
280  protected:
283 
287 
288  std::vector<Vector> mCurrentTotalCauchyStress;
289  std::vector<Vector> mCurrentDeviatoricCauchyStress;
290  std::vector<Vector> mUpdatedTotalCauchyStress;
291  std::vector<Vector> mUpdatedDeviatoricCauchyStress;
292  ConstitutiveLaw::Pointer mpConstitutiveLaw = nullptr;
293 
297 
301 
303 
310  /* void ComputeLumpedMassMatrix(Matrix& rMassMatrix, */
311  /* const double Weight, */
312  /* double& MeanValue); */
313 
315  double &MeanValue,
316  const ShapeFunctionDerivativesType &rShapeDeriv,
317  const double secondLame,
318  double &bulkModulus,
319  const double Weight,
320  double &MeanValueMass,
321  const double TimeStep){};
322 
324  MatrixType StiffnessMatrix,
325  double &meanValueStiff,
326  double &bulkCoefficient,
327  double timeStep) override{};
328 
330  const double Weight);
331 
332  void ComputeBulkMatrixForPressureVel(MatrixType &BulkVelMatrix,
333  const ShapeFunctionsType &rN,
334  const double Weight);
335 
336  void ComputeBoundLHSMatrix(MatrixType &BoundLHSMatrix,
337  const ShapeFunctionsType &rN,
338  const double Weight) override{};
339 
340  void ComputeBoundRHSVector(VectorType &BoundRHSVector,
341  const ShapeFunctionsType &rN,
342  const double TimeStep,
343  const double BoundRHSCoeffAcc,
344  const double BoundRHSCoeffDev) override{};
345 
347  ElementalVariables &rElementalVariables,
348  const unsigned int g,
349  const Vector& rN,
350  const ProcessInfo &rCurrentProcessInfo,
351  double &Density,
352  double &DeviatoricCoeff,
353  double &VolumetricCoeff) override;
354 
355  void CalculateLocalContinuityEqForPressure(MatrixType &rLeftHandSideMatrix,
356  VectorType &rRightHandSideVector,
357  const ProcessInfo &rCurrentProcessInfo) override;
358 
359  double GetThetaMomentum() override { return 1.0; };
360 
361  double GetThetaContinuity() override { return 1.0; };
362 
363  void UpdateStressTensor(ElementalVariables &rElementalVariables);
364 
368 
372 
376 
378 
379  private:
382 
386 
390 
391  friend class Serializer;
392 
393  void save(Serializer &rSerializer) const override
394  {
396  }
397 
398  void load(Serializer &rSerializer) override
399  {
401  }
402 
406 
410 
414 
418 
422 
425 
426  /* /// Copy constructor. */
427  /* TwoStepUpdatedLagrangianVPImplicitSolidElement(TwoStepUpdatedLagrangianVPImplicitSolidElement const& rOther); */
428 
430 
431  }; // Class TwoStepUpdatedLagrangianVPImplicitSolidElement
432 
434 
437 
441 
443  template <unsigned int TDim>
444  inline std::istream &operator>>(std::istream &rIStream,
446  {
447  return rIStream;
448  }
449 
451  template <unsigned int TDim>
452  inline std::ostream &operator<<(std::ostream &rOStream,
454  {
455  rThis.PrintInfo(rOStream);
456  rOStream << std::endl;
457  rThis.PrintData(rOStream);
458 
459  return rOStream;
460  }
461 
462 } // namespace Kratos.
463 
464 #endif // KRATOS_TWO_STEP_UPDATED_LAGRANGIAN_V_P_SOLID_ELEMENT defined
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
Definition: constitutive_law.h:47
Base class for all Elements.
Definition: element.h:60
virtual void MassMatrix(MatrixType &rMassMatrix, const ProcessInfo &rCurrentProcessInfo)
Definition: element.h:926
std::size_t IndexType
Definition: flags.h:74
This defines the geometrical object, base definition of the element and condition entities.
Definition: geometrical_object.h:58
GeometryType & GetGeometry()
Returns the reference of the geometry.
Definition: geometrical_object.h:158
Geometry base class.
Definition: geometry.h:71
This object defines an indexed object.
Definition: indexed_object.h:54
IndexType Id() const
Definition: indexed_object.h:107
This class defines the node.
Definition: node.h:65
PointerVector is a container like stl vector but using a vector to store pointers to its data.
Definition: pointer_vector.h:72
A sorted associative container similar to an STL set, but uses a vector to store pointers to its data...
Definition: pointer_vector_set.h:72
ProcessInfo holds the current value of different solution parameters.
Definition: process_info.h:59
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
A stabilized element for the incompressible Navier-Stokes equations.
Definition: two_step_updated_lagrangian_V_P_implicit_element.h:63
A stabilized element for the incompressible Navier-Stokes equations.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:51
#define KRATOS_THROW_ERROR(ExceptionType, ErrorMessage, MoreInfo)
Definition: define.h:77
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_CATCH(MoreInfo)
Definition: define.h:110
#define KRATOS_TRY
Definition: define.h:109
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
std::vector< std::size_t > EquationIdVectorType
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:82
std::vector< Vector > mUpdatedDeviatoricCauchyStress
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:291
GeometryType::ShapeFunctionsGradientsType ShapeFunctionDerivativesArrayType
Type for an array of shape function gradient matrices.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:95
BaseType::PropertiesType::Pointer pPropertiesType
Definition: two_step_updated_lagrangian_element.h:116
std::size_t SizeType
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:80
Kratos::Vector ShapeFunctionsType
Type for shape function values container.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:89
int Check(const ProcessInfo &rCurrentProcessInfo) const override
Checks the input and that all required Kratos variables have been registered.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.cpp:132
Vector VectorType
Vector type for local contributions to the linear system.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:73
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:86
TwoStepUpdatedLagrangianVPImplicitSolidElement(IndexType NewId=0)
Default constuctor.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:119
TwoStepUpdatedLagrangianVPImplicitSolidElement(IndexType NewId, GeometryType::Pointer pGeometry, pPropertiesType pProperties)
Constuctor using geometry and properties.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:147
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:84
ConstitutiveLaw::Pointer mpConstitutiveLaw
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:292
std::vector< Vector > mUpdatedTotalCauchyStress
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:290
BaseType::ElementalVariables ElementalVariables
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:101
Element::Pointer Clone(IndexType NewId, NodesArrayType const &ThisNodes) const override
It creates a new element pointer and clones the previous element data.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.cpp:25
void InitializeNonLinearIteration(const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.cpp:115
std::string Info() const override
Turn back information as a string.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:258
void CalculateRightHandSide(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:201
void ComputeMeanValueMaterialTangentMatrix(ElementalVariables &rElementalVariables, double &MeanValue, const ShapeFunctionDerivativesType &rShapeDeriv, const double secondLame, double &bulkModulus, const double Weight, double &MeanValueMass, const double TimeStep)
Add integration point contribution to the mass matrix.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:314
BaseType::PropertiesType::Pointer pPropertiesType
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:99
void UpdateStressTensor(ElementalVariables &rElementalVariables)
std::size_t IndexType
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:78
TwoStepUpdatedLagrangianVPImplicitElement< TDim > BaseType
base type:
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:61
TwoStepUpdatedLagrangianVPImplicitSolidElement(TwoStepUpdatedLagrangianVPImplicitSolidElement const &rOther)
copy constructor
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:153
void CalcElasticPlasticCauchySplitted(ElementalVariables &rElementalVariables, const unsigned int g, const Vector &rN, const ProcessInfo &rCurrentProcessInfo, double &Density, double &DeviatoricCoeff, double &VolumetricCoeff) override
void CalculateLocalContinuityEqForPressure(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.cpp:506
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(TwoStepUpdatedLagrangianVPImplicitSolidElement)
Pointer definition of TwoStepUpdatedLagrangianVPImplicitSolidElement.
virtual ~TwoStepUpdatedLagrangianVPImplicitSolidElement()
Destructor.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:158
void Initialize(const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.cpp:69
std::vector< Vector > mCurrentTotalCauchyStress
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:288
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:70
void CalculateLeftHandSide(MatrixType &rLeftHandSideMatrix, const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:193
Node NodeType
Node type (default is: Node)
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:64
BaseType::PropertiesType PropertiesType
Definition: two_step_updated_lagrangian_V_P_implicit_element.h:111
TwoStepUpdatedLagrangianVPImplicitSolidElement(IndexType NewId, const NodesArrayType &ThisNodes)
Constructor using an array of nodes.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:128
Matrix MatrixType
Matrix type for local contributions to the linear system.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:76
void UpdateCauchyStress(unsigned int g, const ProcessInfo &rCurrentProcessInfo) override
Provides the global indices for each one of this element's local rows.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.cpp:477
Kratos::Matrix ShapeFunctionDerivativesType
Type for a matrix containing the shape function gradients.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:92
void ComputeBulkMatrixForPressureVelLump(MatrixType &BulkVelMatrix, const double Weight)
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.cpp:243
void ComputeBulkReductionCoefficient(MatrixType MassMatrix, MatrixType StiffnessMatrix, double &meanValueStiff, double &bulkCoefficient, double timeStep) override
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:323
double GetThetaMomentum() override
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:359
std::vector< Vector > mCurrentDeviatoricCauchyStress
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:289
ConstitutiveLawType::Pointer ConstitutiveLawPointerType
Pointer type for constitutive laws.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:107
BaseType::ElementalVariables ElementalVariables
Definition: two_step_updated_lagrangian_V_P_implicit_element.h:115
BaseType::PropertiesType PropertiesType
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:97
void ComputeBulkMatrixForPressureVel(MatrixType &BulkVelMatrix, const ShapeFunctionsType &rN, const double Weight)
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.cpp:283
double GetThetaContinuity() override
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:361
Geometry< NodeType > GeometryType
Geometry type (using with given NodeType)
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:67
void ComputeBoundRHSVector(VectorType &BoundRHSVector, const ShapeFunctionsType &rN, const double TimeStep, const double BoundRHSCoeffAcc, const double BoundRHSCoeffDev) override
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:340
ConstitutiveLaw ConstitutiveLawType
Reference type definition for constitutive laws.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:104
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:266
void ComputeBoundLHSMatrix(MatrixType &BoundLHSMatrix, const ShapeFunctionsType &rN, const double Weight) override
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:336
void InitializeElementalVariables(ElementalVariables &rElementalVariables) override
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.cpp:301
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: updated_lagrangian_element.hpp:177
TwoStepUpdatedLagrangianVPImplicitSolidElement(IndexType NewId, GeometryType::Pointer pGeometry)
Constructor using a geometry object.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:137
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, pPropertiesType pProperties) const override
Create a new element of this type.
Definition: two_step_updated_lagrangian_V_P_implicit_solid_element.h:178
std::size_t IndexType
The definition of the index type.
Definition: key_hash.h:35
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
ModelPart::NodesContainerType NodesArrayType
Definition: gid_gauss_point_container.h:42
std::istream & operator>>(std::istream &rIStream, LinearMasterSlaveConstraint &rThis)
input stream function
std::ostream & operator<<(std::ostream &rOStream, const LinearMasterSlaveConstraint &rThis)
output stream function
Definition: linear_master_slave_constraint.h:432
def load(f)
Definition: ode_solve.py:307
Definition: updated_lagrangian_element.h:74