[GeoMechanicsApplication] Make stress state classes more similar

applications/GeoMechanicsApplication/custom_constitutive/p_q.cpp

@@ -7,7 +7,6 @@
 //
 //  License:         geo_mechanics_application/license.txt
 //
-//
 //  Main authors:    Anne van de Graaf
 //
 
@@ -32,4 +31,16 @@ double PQ::Q() const noexcept { return mValues[1]; }
 
 double& PQ::Q() noexcept { return mValues[1]; }
 
+PQ& PQ::operator+=(const PQ& rRhs)
+{
+    std::ranges::transform(mValues, rRhs.mValues, mValues.begin(), std::plus{});
+    return *this;
+}
+
+PQ operator+(PQ Lhs, const PQ& rRhs)
+{
+    Lhs += rRhs;
+    return Lhs;
+}
+
 } // namespace Kratos::Geo

applications/GeoMechanicsApplication/custom_constitutive/p_q.hpp

@@ -7,7 +7,6 @@
 //
 //  License:         geo_mechanics_application/license.txt
 //
-//
 //  Main authors:    Anne van de Graaf
 //
 
@@ -59,6 +58,10 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) PQ
         return result;
     }
 
+    PQ& operator+=(const PQ& rRhs);
+    KRATOS_API(GEO_MECHANICS_APPLICATION)
+    friend PQ operator+(PQ Lhs, const PQ& rRhs);
+
 private:
     InternalVectorType mValues = ZeroVector{msVectorSize};
 };

applications/GeoMechanicsApplication/custom_constitutive/principal_stresses.cpp

@@ -14,14 +14,19 @@
 
 namespace Kratos::Geo
 {
-PrincipalStresses::PrincipalStresses(const std::initializer_list<double>& rValues)
-    : PrincipalStresses{std::begin(rValues), std::end(rValues)}
+PrincipalStresses::PrincipalStresses(double Sigma1, double Sigma2, double Sigma3)
 {
+    mValues[0] = Sigma1;
+    mValues[1] = Sigma2;
+    mValues[2] = Sigma3;
 }
 
-const PrincipalStresses::InternalVectorType& PrincipalStresses::Values() const { return mValues; }
+const PrincipalStresses::InternalVectorType& PrincipalStresses::Values() const noexcept
+{
+    return mValues;
+}
 
-PrincipalStresses::InternalVectorType& PrincipalStresses::Values() { return mValues; }
+PrincipalStresses::InternalVectorType& PrincipalStresses::Values() noexcept { return mValues; }
 
 PrincipalStresses& PrincipalStresses::operator+=(const PrincipalStresses& rRhs)
 {

applications/GeoMechanicsApplication/custom_constitutive/principal_stresses.hpp

@@ -17,8 +17,6 @@
 #include "includes/ublas_interface.h"
 
 #include <algorithm>
-#include <initializer_list>
-#include <iterator>
 
 namespace Kratos::Geo
 {
@@ -35,36 +33,33 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) PrincipalStresses
     using InternalVectorType                  = BoundedVector<double, msVectorSize>;
 
     PrincipalStresses() = default;
+    PrincipalStresses(double Sigma1, double Sigma2, double Sigma3);
 
     template <typename VectorType>
-    explicit PrincipalStresses(const VectorType& rStressVector)
-        : PrincipalStresses{std::begin(rStressVector), std::end(rStressVector)}
+    explicit PrincipalStresses(const VectorType& rValues)
     {
-    }
-
-    template <std::forward_iterator InputIt>
-    PrincipalStresses(InputIt First, InputIt Last)
-    {
-        KRATOS_DEBUG_ERROR_IF(std::distance(First, Last) != msVectorSize)
+        // For some reason, the `std::ranges` versions of the below algorithms don't play nicely
+        // with UBlas vector expressions. Therefore, we're using the iterator-style algorithms.
+        auto first = std::begin(rValues);
+        auto last  = std::end(rValues);
+        KRATOS_DEBUG_ERROR_IF(std::distance(first, last) != msVectorSize)
             << "Cannot construct a PrincipalStresses instance: expected " << msVectorSize
-            << " values, but got " << std::distance(First, Last) << " value(s)\n";
+            << " values, but got " << std::distance(first, last) << " value(s)\n";
 
-        std::copy(First, Last, mValues.begin());
+        std::copy(first, last, mValues.begin());
     }
 
-    explicit PrincipalStresses(const std::initializer_list<double>& rValues);
+    [[nodiscard]] const InternalVectorType& Values() const noexcept;
+    [[nodiscard]] InternalVectorType&       Values() noexcept;
 
     template <typename VectorType>
-    VectorType CopyTo()
+    VectorType CopyTo() const
     {
         VectorType result(msVectorSize);
         std::ranges::copy(mValues, result.begin());
         return result;
     }
 
-    [[nodiscard]] const InternalVectorType& Values() const;
-    [[nodiscard]] InternalVectorType&       Values();
-
     PrincipalStresses& operator+=(const PrincipalStresses& rRhs);
     KRATOS_API(GEO_MECHANICS_APPLICATION)
     friend PrincipalStresses operator+(PrincipalStresses Lhs, const PrincipalStresses& rRhs);

applications/GeoMechanicsApplication/custom_constitutive/sigma_tau.cpp

@@ -7,7 +7,6 @@
 //
 //  License:         geo_mechanics_application/license.txt
 //
-//
 //  Main authors:    Anne van de Graaf
 //
 
@@ -16,31 +15,32 @@
 namespace Kratos::Geo
 {
 
-SigmaTau::SigmaTau(const std::initializer_list<double>& rValues)
-    : SigmaTau{rValues.begin(), rValues.end()}
+SigmaTau::SigmaTau(double Sigma, double Tau)
 {
+    mValues[0] = Sigma;
+    mValues[1] = Tau;
 }
 
-const SigmaTau::InternalVectorType& SigmaTau::Values() const { return mValues; }
+const SigmaTau::InternalVectorType& SigmaTau::Values() const noexcept { return mValues; }
 
-double SigmaTau::Sigma() const { return mValues[0]; }
+double SigmaTau::Sigma() const noexcept { return mValues[0]; }
 
-double& SigmaTau::Sigma() { return mValues[0]; }
+double& SigmaTau::Sigma() noexcept { return mValues[0]; }
 
-double SigmaTau::Tau() const { return mValues[1]; }
+double SigmaTau::Tau() const noexcept { return mValues[1]; }
 
-double& SigmaTau::Tau() { return mValues[1]; }
+double& SigmaTau::Tau() noexcept { return mValues[1]; }
 
-SigmaTau& SigmaTau::operator+=(const SigmaTau& rRhsTraction)
+SigmaTau& SigmaTau::operator+=(const SigmaTau& rRhs)
 {
-    std::ranges::transform(mValues, rRhsTraction.mValues, mValues.begin(), std::plus{});
+    std::ranges::transform(mValues, rRhs.mValues, mValues.begin(), std::plus{});
     return *this;
 }
 
-SigmaTau operator+(SigmaTau LhsTraction, const SigmaTau& rRhsTraction)
+SigmaTau operator+(SigmaTau Lhs, const SigmaTau& rRhs)
 {
-    LhsTraction += rRhsTraction;
-    return LhsTraction;
+    Lhs += rRhs;
+    return Lhs;
 }
 
 } // namespace Kratos::Geo

applications/GeoMechanicsApplication/custom_constitutive/sigma_tau.hpp

@@ -7,7 +7,6 @@
 //
 //  License:         geo_mechanics_application/license.txt
 //
-//
 //  Main authors:    Anne van de Graaf
 //
 
@@ -18,8 +17,6 @@
 #include "includes/ublas_interface.h"
 
 #include <algorithm>
-#include <initializer_list>
-#include <iterator>
 
 namespace Kratos::Geo
 {
@@ -31,19 +28,27 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) SigmaTau
     using InternalVectorType                  = BoundedVector<double, msVectorSize>;
 
     SigmaTau() = default;
+    SigmaTau(double Sigma, double Tau);
 
     template <typename VectorType>
-    explicit SigmaTau(const VectorType& rValues) : SigmaTau{std::begin(rValues), std::end(rValues)}
+    explicit SigmaTau(const VectorType& rValues)
     {
-    }
+        // For some reason, the `std::ranges` versions of the below algorithms don't play nicely
+        // with UBlas vector expressions. Therefore, we're using the iterator-style algorithms.
+        auto first = std::begin(rValues);
+        auto last  = std::end(rValues);
+        KRATOS_DEBUG_ERROR_IF(std::distance(first, last) != msVectorSize)
+            << "Cannot construct a SigmaTau instance: expected " << msVectorSize
+            << " values, but got " << std::distance(first, last) << " value(s)\n";
 
-    explicit SigmaTau(const std::initializer_list<double>& rValues);
+        std::copy(first, last, mValues.begin());
+    }
 
-    [[nodiscard]] const InternalVectorType& Values() const;
-    [[nodiscard]] double                    Sigma() const;
-    double&                                 Sigma();
-    [[nodiscard]] double                    Tau() const;
-    double&                                 Tau();
+    [[nodiscard]] const InternalVectorType& Values() const noexcept;
+    [[nodiscard]] double                    Sigma() const noexcept;
+    double&                                 Sigma() noexcept;
+    [[nodiscard]] double                    Tau() const noexcept;
+    double&                                 Tau() noexcept;
 
     template <typename VectorType>
     VectorType CopyTo() const
@@ -53,21 +58,11 @@ class KRATOS_API(GEO_MECHANICS_APPLICATION) SigmaTau
         return result;
     }
 
-    SigmaTau& operator+=(const SigmaTau& rRhsTraction);
+    SigmaTau& operator+=(const SigmaTau& rRhs);
     KRATOS_API(GEO_MECHANICS_APPLICATION)
-    friend SigmaTau operator+(SigmaTau LhsTraction, const SigmaTau& rRhsTraction);
+    friend SigmaTau operator+(SigmaTau Lhs, const SigmaTau& rRhs);
 
 private:
-    template <std::forward_iterator Iter>
-    SigmaTau(Iter First, Iter Last)
-    {
-        KRATOS_DEBUG_ERROR_IF(std::distance(First, Last) != msVectorSize)
-            << "Cannot construct a SigmaTau instance: expected " << msVectorSize
-            << " values, but got " << std::distance(First, Last) << " value(s)\n";
-
-        std::copy(First, Last, mValues.begin());
-    }
-
     InternalVectorType mValues = ZeroVector{msVectorSize};
 };
 

applications/GeoMechanicsApplication/custom_utilities/stress_strain_utilities.cpp

@@ -233,15 +233,15 @@ Vector StressStrainUtilities::RotatePrincipalStresses(const Vector& rPrincipalSt
 
 Geo::SigmaTau StressStrainUtilities::TransformPrincipalStressesToSigmaTau(const Geo::PrincipalStresses& rPrincipalStresses)
 {
-    return Geo::SigmaTau{{0.5 * (rPrincipalStresses.Values()[0] + rPrincipalStresses.Values()[2]),
-                          0.5 * (rPrincipalStresses.Values()[0] - rPrincipalStresses.Values()[2])}};
+    return Geo::SigmaTau{0.5 * (rPrincipalStresses.Values()[0] + rPrincipalStresses.Values()[2]),
+                         0.5 * (rPrincipalStresses.Values()[0] - rPrincipalStresses.Values()[2])};
 }
 
 Geo::PrincipalStresses StressStrainUtilities::TransformSigmaTauToPrincipalStresses(
     const Geo::SigmaTau& rSigmaTau, const Geo::PrincipalStresses& rPrincipalStresses)
 {
-    return Geo::PrincipalStresses{{rSigmaTau.Sigma() + rSigmaTau.Tau(), rPrincipalStresses.Values()[1],
-                                   rSigmaTau.Sigma() - rSigmaTau.Tau()}};
+    return Geo::PrincipalStresses{rSigmaTau.Sigma() + rSigmaTau.Tau(), rPrincipalStresses.Values()[1],
+                                  rSigmaTau.Sigma() - rSigmaTau.Tau()};
 }
 
 Geo::PQ StressStrainUtilities::TransformPrincipalStressesToPandQ(const Geo::PrincipalStresses& rPrincipalStresses)

applications/GeoMechanicsApplication/tests/cpp_tests/custom_constitutive/test_p_q.cpp

@@ -7,7 +7,6 @@
 //
 //  License:         geo_mechanics_application/license.txt
 //
-//
 //  Main authors:    Anne van de Graaf
 //
 
@@ -34,8 +33,8 @@ class TestPQFixture : public ::testing::Test
 {
 };
 
-using TestVectorTypes = ::testing::Types<Vector, BoundedVector<double, 2>, std::vector<double>>;
-TYPED_TEST_SUITE(TestPQFixture, TestVectorTypes);
+using VectorTypesForTypedPQTest = ::testing::Types<Vector, BoundedVector<double, 2>, std::vector<double>>;
+TYPED_TEST_SUITE(TestPQFixture, VectorTypesForTypedPQTest);
 
 TYPED_TEST(TestPQFixture, PQ_CanBeConstructedFromAnyVectorWithSizeOf2)
 {
@@ -115,4 +114,25 @@ TYPED_TEST(TestPQFixture, PQ_CanBeCopiedToAnyVectorTypeWithSizeOf2)
     KRATOS_EXPECT_VECTOR_NEAR(stress_state.CopyTo<TypeParam>(), (std::vector{1.0, 2.0}), Defaults::absolute_tolerance);
 }
 
+TEST_F(KratosGeoMechanicsFastSuiteWithoutKernel, PQ_SupportsCompoundAssignment)
+{
+    // Arrange
+    auto stress_state = Geo::PQ{1.0, 2.0};
+
+    // Act
+    stress_state += Geo::PQ{3.0, 4.0};
+
+    // Assert
+    KRATOS_EXPECT_VECTOR_NEAR(stress_state.Values(), (std::vector{4.0, 6.0}), Defaults::absolute_tolerance);
+}
+
+TEST_F(KratosGeoMechanicsFastSuiteWithoutKernel, PQ_SupportsAdditionOfTwoInstances)
+{
+    // Arrange & Act
+    const auto summed_stress_state = Geo::PQ{1.0, 3.0} + Geo::PQ{2.0, 4.0};
+
+    // Assert
+    KRATOS_EXPECT_VECTOR_NEAR(summed_stress_state.Values(), (std::vector{3.0, 7.0}), Defaults::absolute_tolerance);
+}
+
 } // namespace Kratos::Testing