14 #if !defined(KRATOS_BASE_SOLID_ELEMENT_H_INCLUDED )
15 #define KRATOS_BASE_SOLID_ELEMENT_H_INCLUDED
99 B =
ZeroMatrix(StrainSize, Dimension * NumberOfNodes);
104 Displacements =
ZeroVector(Dimension * NumberOfNodes);
173 ,mThisIntegrationMethod(rOther.mThisIntegrationMethod)
192 void Initialize(
const ProcessInfo& rCurrentProcessInfo)
override;
197 void ResetConstitutiveLaw()
override;
209 void InitializeNonLinearIteration(
const ProcessInfo& rCurrentProcessInfo)
override;
215 void FinalizeNonLinearIteration(
const ProcessInfo& rCurrentProcessInfo)
override;
221 void FinalizeSolutionStep(
const ProcessInfo& rCurrentProcessInfo)
override;
230 Element::Pointer Clone (
240 void EquationIdVector(
241 EquationIdVectorType& rResult,
251 DofsVectorType& rElementalDofList,
261 return mThisIntegrationMethod;
269 void GetValuesVector(
279 void GetFirstDerivativesVector(
289 void GetSecondDerivativesVector(
301 void CalculateLocalSystem(
312 void CalculateLeftHandSide(
322 void CalculateRightHandSide(
336 void AddExplicitContribution(
363 void CalculateMassMatrix(
373 void CalculateDampingMatrix(
386 std::vector<bool>& rOutput,
398 std::vector<int>& rOutput,
410 std::vector<double>& rOutput,
446 std::vector<Vector>& rOutput,
458 std::vector<Matrix>& rOutput,
470 const std::vector<bool>& rValues,
482 const std::vector<int>& rValues,
494 const std::vector<double>& rValues,
506 const std::vector<Vector>& rValues,
518 const std::vector<Matrix>& rValues,
530 const std::vector<ConstitutiveLaw::Pointer>& rValues,
542 std::vector<ConstitutiveLaw::Pointer>& rValues,
583 int Check(
const ProcessInfo& rCurrentProcessInfo )
const override;
589 const std::size_t MatrixSize);
606 std::string
Info()
const override
608 std::stringstream buffer;
609 buffer <<
"Base Solid Element #" << Id() <<
"\nConstitutive law: " << mConstitutiveLawVector[0]->Info();
616 rOStream <<
"Base Solid Element #" << Id() <<
"\nConstitutive law: " << mConstitutiveLawVector[0]->Info();
622 pGetGeometry()->PrintData(rOStream);
656 mThisIntegrationMethod = ThisIntegrationMethod;
665 mConstitutiveLawVector = ThisConstitutiveLawVector;
671 virtual void InitializeMaterial();
681 virtual bool UseElementProvidedStrain()
const;
691 virtual void CalculateAll(
695 const bool CalculateStiffnessMatrixFlag,
696 const bool CalculateResidualVectorFlag
705 virtual void CalculateKinematicVariables(
706 KinematicVariables& rThisKinematicVariables,
719 virtual void SetConstitutiveVariables(
720 KinematicVariables& rThisKinematicVariables,
721 ConstitutiveVariables& rThisConstitutiveVariables,
736 virtual void CalculateConstitutiveVariables(
737 KinematicVariables& rThisKinematicVariables,
738 ConstitutiveVariables& rThisConstitutiveVariables,
750 Matrix& CalculateDeltaDisplacement(
Matrix& DeltaDisplacement)
const;
761 virtual double CalculateDerivativesOnReferenceConfiguration(
778 double CalculateDerivativesOnCurrentConfiguration(
804 virtual void CalculateAndAddKm(
808 const double IntegrationWeight
818 void CalculateAndAddKg(
821 const Vector& rStressVector,
822 const double IntegrationWeight
834 virtual void CalculateAndAddResidualVector(
836 const KinematicVariables& rThisKinematicVariables,
839 const Vector& rStressVector,
840 const double IntegrationWeight
851 void CalculateAndAddExtForceContribution(
856 const double IntegrationWeight
865 virtual double GetIntegrationWeight(
912 void CalculateLumpedMassVector(
914 const ProcessInfo& rCurrentProcessInfo)
const override;
921 void CalculateDampingMatrixWithLumpedMass(
933 template<
class TType>
934 void GetValueOnConstitutiveLaw(
936 std::vector<TType>& rOutput
941 for (
IndexType point_number = 0; point_number <integration_points.size(); ++point_number ) {
942 mConstitutiveLawVector[point_number]->GetValue( rVariable,rOutput[point_number]);
954 template<
class TType>
955 void CalculateOnConstitutiveLaw(
956 const Variable<TType>& rVariable,
957 std::vector<TType>& rOutput,
958 const ProcessInfo& rCurrentProcessInfo
963 const SizeType number_of_nodes = GetGeometry().size();
968 ConstitutiveVariables this_constitutive_variables(
strain_size);
971 ConstitutiveLaw::Parameters Values(GetGeometry(),GetProperties(),rCurrentProcessInfo);
974 Flags& ConstitutiveLawOptions=Values.GetOptions();
975 ConstitutiveLawOptions.Set(ConstitutiveLaw::USE_ELEMENT_PROVIDED_STRAIN, UseElementProvidedStrain());
976 ConstitutiveLawOptions.Set(ConstitutiveLaw::COMPUTE_STRESS,
true);
977 ConstitutiveLawOptions.Set(ConstitutiveLaw::COMPUTE_CONSTITUTIVE_TENSOR,
false);
979 Values.SetStrainVector(this_constitutive_variables.StrainVector);
981 for (
IndexType point_number = 0; point_number < integration_points.size(); ++point_number) {
983 this->CalculateKinematicVariables(this_kinematic_variables, point_number, this->GetIntegrationMethod());
986 this->SetConstitutiveVariables(this_kinematic_variables, this_constitutive_variables, Values, point_number, integration_points);
988 rOutput[point_number] = mConstitutiveLawVector[point_number]->CalculateValue( Values, rVariable, rOutput[point_number] );
1006 void save(
Serializer& rSerializer )
const override;
This is base class used to define the solid elements.
Definition: base_solid_element.h:67
ConstitutiveLaw ConstitutiveLawType
Reference type definition for constitutive laws.
Definition: base_solid_element.h:134
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(BaseSolidElement)
BaseSolidElement(BaseSolidElement const &rOther)
Definition: base_solid_element.h:171
ConstitutiveLawType::StressMeasure StressMeasureType
StressMeasure from constitutive laws.
Definition: base_solid_element.h:140
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: base_solid_element.h:614
std::string Info() const override
Turn back information as a string.
Definition: base_solid_element.h:606
~BaseSolidElement() override
Definition: base_solid_element.h:177
Element BaseType
The base element type.
Definition: base_solid_element.h:149
Node NodeType
This is the definition of the node.
Definition: base_solid_element.h:146
BaseSolidElement(IndexType NewId, GeometryType::Pointer pGeometry)
Definition: base_solid_element.h:163
BaseSolidElement()
Definition: base_solid_element.h:159
ConstitutiveLawType::Pointer ConstitutiveLawPointerType
Pointer type for constitutive laws.
Definition: base_solid_element.h:137
void AddExplicitContribution(const VectorType &rRHSVector, const Variable< VectorType > &rRHSVariable, const Variable< array_1d< double, 3 > > &rDestinationVariable, const ProcessInfo &rCurrentProcessInfo) override
This function is designed to make the element to assemble an rRHS vector identified by a variable rRH...
BaseSolidElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
Definition: base_solid_element.h:167
GeometryData::IntegrationMethod IntegrationMethod
Type definition for integration methods.
Definition: base_solid_element.h:143
void SetConstitutiveLawVector(const std::vector< ConstitutiveLaw::Pointer > &ThisConstitutiveLawVector)
Sets the used constitutive laws.
Definition: base_solid_element.h:663
IntegrationMethod GetIntegrationMethod() const override
Returns the used integration method.
Definition: base_solid_element.h:259
std::vector< ConstitutiveLaw::Pointer > mConstitutiveLawVector
Currently selected integration methods.
Definition: base_solid_element.h:640
IntegrationMethod mThisIntegrationMethod
Definition: base_solid_element.h:638
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: base_solid_element.h:620
void SetIntegrationMethod(const IntegrationMethod &ThisIntegrationMethod)
Sets the used integration method.
Definition: base_solid_element.h:654
Definition: constitutive_law.h:47
StressMeasure
Definition: constitutive_law.h:69
@ StressMeasure_PK2
Definition: constitutive_law.h:71
Base class for all Elements.
Definition: element.h:60
std::size_t SizeType
Definition: element.h:94
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
IntegrationMethod
Definition: geometry_data.h:76
std::vector< IntegrationPointType > IntegrationPointsArrayType
Definition: geometry.h:161
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
Variable class contains all information needed to store and retrive data from a data container.
Definition: variable.h:63
std::size_t IndexType
The definition of the index type.
Definition: key_hash.h:35
void SetValuesOnIntegrationPoints(TObject &dummy, const Variable< TDataType > &rVariable, const std::vector< TDataType > &values, const ProcessInfo &rCurrentProcessInfo)
Definition: add_mesh_to_python.cpp:185
pybind11::list CalculateOnIntegrationPoints(TObject &dummy, const Variable< TDataType > &rVariable, const ProcessInfo &rProcessInfo)
Definition: add_mesh_to_python.cpp:142
void InitializeSolutionStep(ConstructionUtility &rThisUtil, std::string ThermalSubModelPartName, std::string MechanicalSubModelPartName, std::string HeatFluxSubModelPartName, std::string HydraulicPressureSubModelPartName, bool thermal_conditions, bool mechanical_conditions, int phase)
Definition: add_custom_utilities_to_python.cpp:45
void CalculateRayleighDampingMatrix(Element &rElement, Element::MatrixType &rDampingMatrix, const ProcessInfo &rCurrentProcessInfo, const std::size_t MatrixSize)
Method to calculate the rayleigh damping-matrix.
Definition: structural_mechanics_element_utilities.cpp:193
double GetRayleighAlpha(const Properties &rProperties, const ProcessInfo &rCurrentProcessInfo)
Method to get the rayleigh-alpha parameter.
Definition: structural_mechanics_element_utilities.cpp:136
double GetRayleighBeta(const Properties &rProperties, const ProcessInfo &rCurrentProcessInfo)
Method to get the rayleigh-beta parameter.
Definition: structural_mechanics_element_utilities.cpp:154
array_1d< double, 3 > GetBodyForce(const Element &rElement, const GeometryType::IntegrationPointsArrayType &rIntegrationPoints, const IndexType PointNumber)
This method returns the computed the computed body force.
Definition: structural_mechanics_element_utilities.cpp:71
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
KratosZeroVector< double > ZeroVector
Definition: amatrix_interface.h:561
AMatrix::IdentityMatrix< double > IdentityMatrix
Definition: amatrix_interface.h:564
KratosZeroMatrix< double > ZeroMatrix
Definition: amatrix_interface.h:559
std::size_t SizeType
The definition of the size type.
Definition: mortar_classes.h:43
ModelPart::NodesContainerType NodesArrayType
Definition: gid_gauss_point_container.h:42
int strain_size
Definition: generate_hyper_elastic_simo_taylor_neo_hookean.py:16
F
Definition: hinsberg_optimization.py:144
int dimension
Definition: isotropic_damage_automatic_differentiation.py:123
def load(f)
Definition: ode_solve.py:307
N
Definition: sensitivityMatrix.py:29
B
Definition: sensitivityMatrix.py:76
Definition: base_solid_element.h:112
Vector StressVector
Definition: base_solid_element.h:114
Matrix D
Definition: base_solid_element.h:115
Vector StrainVector
Definition: base_solid_element.h:113
ConstitutiveVariables(const SizeType StrainSize)
Definition: base_solid_element.h:121
Definition: base_solid_element.h:73
Matrix InvJ0
Definition: base_solid_element.h:80
Vector Displacements
Definition: base_solid_element.h:82
Matrix DN_DX
Definition: base_solid_element.h:81
Vector N
Definition: base_solid_element.h:74
double detF
Definition: base_solid_element.h:76
Matrix J0
Definition: base_solid_element.h:79
double detJ0
Definition: base_solid_element.h:78
Matrix B
Definition: base_solid_element.h:75
KinematicVariables(const SizeType StrainSize, const SizeType Dimension, const SizeType NumberOfNodes)
Definition: base_solid_element.h:90
Matrix F
Definition: base_solid_element.h:77
Definition: constitutive_law.h:189
Definition: geometrical_sensitivity_utility.h:33