9 #if !defined(KRATOS_THREE_STEP_UPDATED_LAGRANGIAN_ELEMENT_H_INCLUDED)
10 #define KRATOS_THREE_STEP_UPDATED_LAGRANGIAN_ELEMENT_H_INCLUDED
24 #include "utilities/geometry_utilities.h"
61 template <
unsigned int TDim>
204 const ProcessInfo &rCurrentProcessInfo)
override{};
208 const double Weight);
211 const double Density,
213 const double OldPressure,
216 const double Weight);
221 const double theta_velocity);
225 const double Weight);
229 const double TimeStep,
230 const double BoundRHSCoeffAcc,
231 const double BoundRHSCoeffDev);
234 const double TimeStep,
235 const double BoundRHSCoeffAcc,
236 const double BoundRHSCoeffDev,
241 const double Weight);
245 const double Density,
246 const double Viscosity,
255 KRATOS_THROW_ERROR(std::logic_error,
"ThreeStepUpdatedLagrangianElement::CalculateLeftHandSide not implemented",
"");
266 const ProcessInfo &rCurrentProcessInfo)
override{};
269 const ProcessInfo &rCurrentProcessInfo)
override{};
280 const ProcessInfo &rCurrentProcessInfo)
const override;
288 const ProcessInfo &rCurrentProcessInfo)
const override;
322 std::string
Info()
const override
324 std::stringstream buffer;
325 buffer <<
"ThreeStepUpdatedLagrangianElement #" << this->
Id();
330 void PrintInfo(std::ostream &rOStream)
const override {}
360 const ProcessInfo &rCurrentProcessInfo)
override{};
373 const double lagMultiplier);
377 const double Density,
378 const double Viscosity,
383 std::cout <<
"I SHOULD NOT ENTER HERE!" << std::endl;
389 std::cout <<
"I SHOULD NOT ENTER HERE!" << std::endl;
394 const ProcessInfo &rCurrentProcessInfo)
override{};
397 const double Weight);
400 const double Density,
402 const double Weight)
override;
406 const unsigned int g,
410 double &DeviatoricCoeff,
411 double &VolumetricCoeff)
override
442 void save(
Serializer &rSerializer)
const override
492 template <
unsigned int TDim>
500 template <
unsigned int TDim>
505 rOStream << std::endl;
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
Definition: constitutive_law.h:47
Base class for all Elements.
Definition: element.h:60
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: three_step_updated_lagrangian_element.h:63
A stabilized element for the incompressible Navier-Stokes equations.
Definition: updated_lagrangian_element.h:64
Variable class contains all information needed to store and retrive data from a data container.
Definition: variable.h:63
#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
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: three_step_updated_lagrangian_element.h:80
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, pPropertiesType pProperties) const override
Create a new element of this type.
Definition: three_step_updated_lagrangian_element.h:186
void CalculateLeftHandSide(MatrixType &rLeftHandSideMatrix, const ProcessInfo &rCurrentProcessInfo) override
Definition: three_step_updated_lagrangian_element.h:251
void GetDofList(DofsVectorType &rElementalDofList, const ProcessInfo &rCurrentProcessInfo) const override
Returns a list of the element's Dofs.
Definition: three_step_updated_lagrangian_element.cpp:63
UpdatedLagrangianElement< TDim > BaseType
Definition: three_step_updated_lagrangian_element.h:71
ConstitutiveLaw ConstitutiveLawType
Reference type definition for constitutive laws.
Definition: three_step_updated_lagrangian_element.h:108
virtual void Calculate(const Variable< array_1d< double, 3 >> &rVariable, array_1d< double, 3 > &rOutput, const ProcessInfo &rCurrentProcessInfo) override
Definition: three_step_updated_lagrangian_element.h:264
Kratos::Matrix ShapeFunctionDerivativesType
Type for a matrix containing the shape function gradients.
Definition: three_step_updated_lagrangian_element.h:102
Vector VectorType
Vector type for local contributions to the linear system.
Definition: three_step_updated_lagrangian_element.h:83
double ElementSize()
Definition: three_step_updated_lagrangian_element.cpp:244
virtual void CalculateLocalSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Calculate the element's local contribution to the system for the current step.
Definition: three_step_updated_lagrangian_element.h:202
void CalculateTauPSPG(double &TauOne, double ElemSize, const double Density, const double Viscosity, const ProcessInfo &rCurrentProcessInfo)
Definition: three_step_updated_lagrangian_element.cpp:981
Element::Pointer Clone(IndexType NewId, NodesArrayType const &ThisNodes) const override
It creates a new element pointer and clones the previous element data.
Definition: three_step_updated_lagrangian_element.cpp:21
void ComputeBoundRHSVectorComplete(VectorType &BoundRHSVector, const double TimeStep, const double BoundRHSCoeffAcc, const double BoundRHSCoeffDev, const VectorType SpatialDefRate)
virtual void UpdateCauchyStress(unsigned int g, const ProcessInfo &rCurrentProcessInfo) override
Definition: three_step_updated_lagrangian_element.h:290
Geometry< NodeType > GeometryType
Geometry type (using with given NodeType)
Definition: three_step_updated_lagrangian_element.h:77
std::size_t IndexType
Definition: three_step_updated_lagrangian_element.h:88
void ComputeLumpedMassMatrix(Matrix &rMassMatrix, const double Weight)
Definition: three_step_updated_lagrangian_element.cpp:269
void CalculatePSPGPressureSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo)
Definition: three_step_updated_lagrangian_element.cpp:341
std::string Info() const override
Turn back information as a string.
Definition: three_step_updated_lagrangian_element.h:322
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: three_step_updated_lagrangian_element.h:94
void CalculateTauFIC(double &TauOne, double ElemSize, const double Density, const double Viscosity, const ProcessInfo &rCurrentProcessInfo)
Definition: three_step_updated_lagrangian_element.cpp:766
BaseType::ElementalVariables ElementalVariables
Definition: three_step_updated_lagrangian_element.h:117
ThreeStepUpdatedLagrangianElement(IndexType NewId, GeometryType::Pointer pGeometry, pPropertiesType pProperties)
Constuctor using geometry and properties.
Definition: three_step_updated_lagrangian_element.h:157
void ComputeBoundLHSMatrix(MatrixType &BoundLHSMatrix, const ShapeFunctionsType &rN, const double Weight)
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: three_step_updated_lagrangian_element.h:96
virtual void CalculateLocalMomentumEquations(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: three_step_updated_lagrangian_element.h:358
ThreeStepUpdatedLagrangianElement(IndexType NewId, GeometryType::Pointer pGeometry)
Constructor using a geometry object.
Definition: three_step_updated_lagrangian_element.h:147
void CalculateMassMatrix(Matrix &rMassMatrix, const ProcessInfo &rCurrentProcessInfo) override
Add integration point contribution to the mass matrix.
Definition: three_step_updated_lagrangian_element.h:393
virtual double GetThetaMomentum() override
Definition: three_step_updated_lagrangian_element.h:381
virtual double GetThetaContinuity() override
Definition: three_step_updated_lagrangian_element.h:387
Kratos::Vector ShapeFunctionsType
Type for shape function values container.
Definition: three_step_updated_lagrangian_element.h:99
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: three_step_updated_lagrangian_element.h:330
void AddMomentumRHSTerms(Vector &rRHSVector, const double Density, const array_1d< double, 3 > &rBodyForce, const double OldPressure, const ShapeFunctionsType &rN, const ShapeFunctionDerivativesType &rDN_DX, const double Weight)
Definition: three_step_updated_lagrangian_element.cpp:94
int Check(const ProcessInfo &rCurrentProcessInfo) const override
Checks the input and that all required Kratos variables have been registered.
Definition: three_step_updated_lagrangian_element.h:311
void AddExternalForces(Vector &rRHSVector, const double Density, const ShapeFunctionsType &rN, const double Weight) override
Definition: three_step_updated_lagrangian_element.cpp:123
ConstitutiveLawType::Pointer ConstitutiveLawPointerType
Pointer type for constitutive laws.
Definition: three_step_updated_lagrangian_element.h:111
GeometryType::ShapeFunctionsGradientsType ShapeFunctionDerivativesArrayType
Type for an array of shape function gradient matrices.
Definition: three_step_updated_lagrangian_element.h:105
BaseType::PropertiesType::Pointer pPropertiesType
Definition: three_step_updated_lagrangian_element.h:115
ThreeStepUpdatedLagrangianElement(ThreeStepUpdatedLagrangianElement const &rOther)
copy constructor
Definition: three_step_updated_lagrangian_element.h:163
void Initialize(const ProcessInfo &rCurrentProcessInfo) override
Definition: three_step_updated_lagrangian_element.h:194
Node NodeType
Node type (default is: Node)
Definition: three_step_updated_lagrangian_element.h:74
void AddMomentumMassTerm(Matrix &rMassMatrix, const ShapeFunctionsType &rN, const double Weight)
Definition: three_step_updated_lagrangian_element.cpp:317
void ComputeStabLaplacianMatrix(MatrixType &StabLaplacianMatrix, const ShapeFunctionDerivativesType &rShapeDeriv, const double Weight)
Definition: three_step_updated_lagrangian_element.cpp:1021
Matrix MatrixType
Matrix type for local contributions to the linear system.
Definition: three_step_updated_lagrangian_element.h:86
void CalculateRightHandSide(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: three_step_updated_lagrangian_element.h:268
ThreeStepUpdatedLagrangianElement(IndexType NewId=0)
Default constuctor.
Definition: three_step_updated_lagrangian_element.h:129
void ComputeBoundRHSVector(VectorType &BoundRHSVector, const ShapeFunctionsType &rN, const double TimeStep, const double BoundRHSCoeffAcc, const double BoundRHSCoeffDev)
Element::PropertiesType PropertiesType
Definition: updated_lagrangian_element.h:142
void AddViscousTerm(MatrixType &rDampingMatrix, const ShapeFunctionDerivativesType &rShapeDeriv, const double Weight, const double theta_velocity)
virtual void CalcElasticPlasticCauchySplitted(ElementalVariables &rElementalVariables, const unsigned int g, const Vector &rN, const ProcessInfo &rCurrentProcessInfo, double &Density, double &DeviatoricCoeff, double &VolumetricCoeff) override
Definition: three_step_updated_lagrangian_element.h:404
std::size_t SizeType
Definition: three_step_updated_lagrangian_element.h:90
void InitializeElementalVariables(ElementalVariables &rElementalVariables) override
Definition: three_step_updated_lagrangian_element.cpp:1043
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(ThreeStepUpdatedLagrangianElement)
Pointer definition of ThreeStepUpdatedLagrangianElement.
void EquationIdVector(EquationIdVectorType &rResult, const ProcessInfo &rCurrentProcessInfo) const override
Provides the global indices for each one of this element's local rows.
Definition: three_step_updated_lagrangian_element.cpp:32
void InitializeSolutionStep(const ProcessInfo &rCurrentProcessInfo) override
Initializes the element and all geometric information required for the problem.
Definition: three_step_updated_lagrangian_element.h:197
ConstitutiveLaw::Pointer mpConstitutiveLaw
Definition: three_step_updated_lagrangian_element.h:348
std::vector< std::size_t > EquationIdVectorType
Definition: three_step_updated_lagrangian_element.h:92
BaseType::PropertiesType PropertiesType
Definition: three_step_updated_lagrangian_element.h:113
void InitializeNonLinearIteration(const ProcessInfo &rCurrentProcessInfo) override
Definition: three_step_updated_lagrangian_element.h:199
virtual ~ThreeStepUpdatedLagrangianElement()
Destructor.
Definition: three_step_updated_lagrangian_element.h:166
ThreeStepUpdatedLagrangianElement(IndexType NewId, const NodesArrayType &ThisNodes)
Constructor using an array of nodes.
Definition: three_step_updated_lagrangian_element.h:138
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: updated_lagrangian_element.hpp:177
void ComputeBoundaryTermsForPressureSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ShapeFunctionsType &rN, const double lagMultiplier)
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