Clean up truss elements (#12563)

- Extracted duplicated code of member `FinalizeSolutionStep` and moved it to a new member `FinalizeMaterialResponse`. The latter receives the appropriate strain measure to carry out its task. Also its implementation has been slightly simplified.
- Eliminated duplicated code for calculating the Cauchy stress vector.
- Applied the rule of zero to class `GeoTrussElement`.
- Removed unused aliases.
- Removed an unneeded intermediate variable.
- Use alias `BaseType` when appropriate.
- Unnamed an unused function parameter.
- Removed an empty statement.
- Removed some unnecessary namespace qualifiers.
- Avoid using a local (empty) process info object. Rather use the one that is passed as an argument. This is in line with the code of the StructuralMechanicsApplication.

applications/GeoMechanicsApplication/custom_elements/geo_cable_element.hpp

@@ -48,7 +48,6 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) GeoCableElement : public GeoTrussEle
     using VectorType     = Element::VectorType;
 
     using FullDofMatrixType = typename GeoTrussElementBase<TDim, TNumNodes>::FullDofMatrixType;
-    using FullDofVectorType = typename GeoTrussElementBase<TDim, TNumNodes>::FullDofVectorType;
 
     using GeoTrussElementBase<TDim, TNumNodes>::mpConstitutiveLaw;
     using GeoTrussElement<TDim, TNumNodes>::mInternalStressesFinalizedPrevious;

applications/GeoMechanicsApplication/custom_elements/geo_linear_truss_element.cpp

@@ -279,13 +279,7 @@ void GeoLinearTrussElement<TDim, TNumNodes>::FinalizeSolutionStep(const ProcessI
 {
     KRATOS_TRY
 
-    ConstitutiveLaw::Parameters Values(this->GetGeometry(), this->GetProperties(), rCurrentProcessInfo);
-    Vector temp_strain = ZeroVector(1);
-    Vector temp_stress = ZeroVector(1);
-    temp_strain[0]     = this->CalculateLinearStrain();
-    Values.SetStrainVector(temp_strain);
-    Values.SetStressVector(temp_stress);
-    mpConstitutiveLaw->FinalizeMaterialResponse(Values, ConstitutiveLaw::StressMeasure_PK2);
+    this->FinalizeMaterialResponse(this->CalculateLinearStrain(), rCurrentProcessInfo);
 
     mInternalStressesFinalized = mInternalStresses + mInternalStressesFinalizedPrevious;
 

applications/GeoMechanicsApplication/custom_elements/geo_truss_element.cpp

@@ -44,8 +44,7 @@ Element::Pointer GeoTrussElement<TDim, TNumNodes>::Create(IndexType
                                                           NodesArrayType const&   rThisNodes,
                                                           PropertiesType::Pointer pProperties) const
 {
-    const GeometryType& rGeom = this->GetGeometry();
-    return Kratos::make_intrusive<GeoTrussElement>(NewId, rGeom.Create(rThisNodes), pProperties);
+    return make_intrusive<GeoTrussElement>(NewId, this->GetGeometry().Create(rThisNodes), pProperties);
 }
 
 //----------------------------------------------------------------------------------------
@@ -54,13 +53,7 @@ Element::Pointer GeoTrussElement<TDim, TNumNodes>::Create(IndexType
                                                           GeometryType::Pointer   pGeom,
                                                           PropertiesType::Pointer pProperties) const
 {
-    return Kratos::make_intrusive<GeoTrussElement>(NewId, pGeom, pProperties);
-}
-
-//----------------------------------------------------------------------------------------
-template <unsigned int TDim, unsigned int TNumNodes>
-GeoTrussElement<TDim, TNumNodes>::~GeoTrussElement()
-{
+    return make_intrusive<GeoTrussElement>(NewId, pGeom, pProperties);
 }
 
 //----------------------------------------------------------------------------------------
@@ -82,7 +75,7 @@ void GeoTrussElement<TDim, TNumNodes>::Initialize(const ProcessInfo& rCurrentPro
 {
     KRATOS_TRY
 
-    GeoTrussElementBase<TDim, TNumNodes>::Initialize(rCurrentProcessInfo);
+    BaseType::Initialize(rCurrentProcessInfo);
 
     if (rCurrentProcessInfo.Has(RESET_DISPLACEMENTS)) {
         bool ResetDisplacement = rCurrentProcessInfo[RESET_DISPLACEMENTS];
@@ -115,8 +108,7 @@ void GeoTrussElement<TDim, TNumNodes>::CalculateOnIntegrationPoints(const Variab
         Vector strain = ZeroVector(TDim);
         strain[0]     = this->CalculateGreenLagrangeStrain();
         rOutput[0]    = strain;
-    }
-    if (rVariable == PK2_STRESS_VECTOR) {
+    } else if (rVariable == PK2_STRESS_VECTOR) {
         ConstitutiveLaw::Parameters Values(this->GetGeometry(), this->GetProperties(), rCurrentProcessInfo);
         Vector temp_strain = ZeroVector(1);
         temp_strain[0]     = this->CalculateGreenLagrangeStrain();
@@ -129,24 +121,12 @@ void GeoTrussElement<TDim, TNumNodes>::CalculateOnIntegrationPoints(const Variab
             temp_internal_stresses[i] += mInternalStressesFinalizedPrevious[i];
 
         rOutput[0] = temp_internal_stresses;
-    }
-    if (rVariable == CAUCHY_STRESS_VECTOR) {
-        ProcessInfo temp_process_information;
-
-        ConstitutiveLaw::Parameters Values(this->GetGeometry(), this->GetProperties(), temp_process_information);
-        Vector temp_strain = ZeroVector(1);
-        temp_strain[0]     = this->CalculateGreenLagrangeStrain();
-        Values.SetStrainVector(temp_strain);
-        array_1d<double, 3> temp_internal_stresses = ZeroVector(3);
-        mpConstitutiveLaw->CalculateValue(Values, FORCE, temp_internal_stresses);
-
-        for (unsigned int i = 0; i < TDim; ++i)
-            temp_internal_stresses[i] += mInternalStressesFinalizedPrevious[i];
+    } else if (rVariable == CAUCHY_STRESS_VECTOR) {
+        CalculateOnIntegrationPoints(PK2_STRESS_VECTOR, rOutput, rCurrentProcessInfo);
 
-        const double l  = StructuralMechanicsElementUtilities::CalculateCurrentLength3D2N(*this);
-        const double L0 = StructuralMechanicsElementUtilities::CalculateReferenceLength3D2N(*this);
-
-        rOutput[0] = temp_internal_stresses * l / L0;
+        const auto l  = StructuralMechanicsElementUtilities::CalculateCurrentLength3D2N(*this);
+        const auto L0 = StructuralMechanicsElementUtilities::CalculateReferenceLength3D2N(*this);
+        rOutput[0] *= (l / L0);
     }
 
     KRATOS_CATCH("")
@@ -156,7 +136,7 @@ void GeoTrussElement<TDim, TNumNodes>::CalculateOnIntegrationPoints(const Variab
 template <unsigned int TDim, unsigned int TNumNodes>
 void GeoTrussElement<TDim, TNumNodes>::CalculateOnIntegrationPoints(const Variable<array_1d<double, 3>>& rVariable,
                                                                     std::vector<array_1d<double, 3>>& rOutput,
-                                                                    const ProcessInfo& rCurrentProcessInfo)
+                                                                    const ProcessInfo&)
 {
     KRATOS_TRY
 
@@ -245,9 +225,10 @@ void GeoTrussElement<TDim, TNumNodes>::UpdateInternalForces(BoundedVector<double
 template <unsigned int TDim, unsigned int TNumNodes>
 void GeoTrussElement<TDim, TNumNodes>::FinalizeSolutionStep(const ProcessInfo& rCurrentProcessInfo)
 {
-    KRATOS_TRY;
+    KRATOS_TRY
+
+    this->FinalizeMaterialResponse(this->CalculateGreenLagrangeStrain(), rCurrentProcessInfo);
 
-    GeoTrussElementBase<TDim, TNumNodes>::FinalizeSolutionStep(rCurrentProcessInfo);
     mInternalStressesFinalized = mInternalStresses + mInternalStressesFinalizedPrevious;
 
     KRATOS_CATCH("")

applications/GeoMechanicsApplication/custom_elements/geo_truss_element.hpp

@@ -49,8 +49,6 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) GeoTrussElement : public GeoTrussEle
     using SizeType          = Element::SizeType;
     using MatrixType        = Element::MatrixType;
     using VectorType        = Element::VectorType;
-    using FullDofMatrixType = typename GeoTrussElementBase<TDim, TNumNodes>::FullDofMatrixType;
-    using FullDofVectorType = typename GeoTrussElementBase<TDim, TNumNodes>::FullDofVectorType;
 
     using GeoTrussElementBase<TDim, TNumNodes>::mpConstitutiveLaw;
 
@@ -60,8 +58,6 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) GeoTrussElement : public GeoTrussEle
     GeoTrussElement(IndexType NewId, GeometryType::Pointer pGeometry);
     GeoTrussElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties);
 
-    ~GeoTrussElement() override;
-
     /**
      * @brief Creates a new element
      * @param NewId The Id of the new created element

applications/GeoMechanicsApplication/custom_elements/geo_truss_element_base.cpp

@@ -1051,21 +1051,15 @@ void GeoTrussElementBase<2, 2>::CalculateElasticStiffnessMatrix(MatrixType& rEla
     KRATOS_CATCH("")
 }
 
-//----------------------------------------------------------------------------------------
 template <unsigned int TDim, unsigned int TNumNodes>
-void GeoTrussElementBase<TDim, TNumNodes>::FinalizeSolutionStep(const ProcessInfo& rCurrentProcessInfo)
+void GeoTrussElementBase<TDim, TNumNodes>::FinalizeMaterialResponse(double Strain, const ProcessInfo& rCurrentProcessInfo)
 {
-    KRATOS_TRY;
-
-    ConstitutiveLaw::Parameters Values(GetGeometry(), GetProperties(), rCurrentProcessInfo);
-    Vector                      temp_strain = ZeroVector(1);
-    Vector                      temp_stress = ZeroVector(1);
-    temp_strain[0]                          = CalculateGreenLagrangeStrain();
-    Values.SetStrainVector(temp_strain);
-    Values.SetStressVector(temp_stress);
-    mpConstitutiveLaw->FinalizeMaterialResponse(Values, ConstitutiveLaw::StressMeasure_PK2);
-
-    KRATOS_CATCH("")
+    auto values = ConstitutiveLaw::Parameters(this->GetGeometry(), this->GetProperties(), rCurrentProcessInfo);
+    auto temp_strain = Vector{ScalarVector(1, Strain)};
+    auto temp_stress = Vector{ZeroVector(1)};
+    values.SetStrainVector(temp_strain);
+    values.SetStressVector(temp_stress);
+    mpConstitutiveLaw->FinalizeMaterialResponse(values, ConstitutiveLaw::StressMeasure_PK2);
 }
 
 //----------------------------------------------------------------------------------------

applications/GeoMechanicsApplication/custom_elements/geo_truss_element_base.hpp

@@ -206,7 +206,8 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) GeoTrussElementBase : public Element
 
     double ReturnTangentModulus1D(const ProcessInfo& rCurrentProcessInfo);
 
-    void FinalizeSolutionStep(const ProcessInfo& rCurrentProcessInfo) override;
+protected:
+    void FinalizeMaterialResponse(double Strain, const ProcessInfo& rCurrentProcessInfo);
 
 private:
     /**