Check whether the matrix-vector product can be calculated

KratosMultiphysics · Sep 10, 2024 · a49b75f · a49b75f
1 parent 30cd1f5
commit a49b75f
Show file tree

Hide file tree

Showing 2 changed files with 26 additions and 3 deletions.
diff --git a/applications/GeoMechanicsApplication/custom_utilities/equation_of_motion_utilities.cpp b/applications/GeoMechanicsApplication/custom_utilities/equation_of_motion_utilities.cpp
@@ -54,13 +54,22 @@ void CheckThatStressVectorsHaveSameSize(const std::vector<Vector>& rStressVector
         << "Cannot calculate the internal force vector: stress vectors have different sizes\n";
 }
 
+void CheckThatTransposedMatrixVectorProductCanBeCalculated(const Matrix& rFirstB, const Vector& rFirstStressVector)
+{
+    // B-matrix will be transposed, so check against its number of rows rather than its number of columns!
+    KRATOS_DEBUG_ERROR_IF(rFirstB.size1() != rFirstStressVector.size())
+        << "Cannot calculate the internal force vector: matrix-vector product cannot be calculated "
+           "due to size mismatch\n";
+}
+
 void CheckInputOfCalculateInternalForceVector(const std::vector<Matrix>& rBs,
                                               const std::vector<Vector>& rStressVectors,
                                               const std::vector<double>& rIntegrationCoefficients)
 {
     CheckThatVectorsHaveSameSizeAndAreNotEmpty(rBs, rStressVectors, rIntegrationCoefficients);
     CheckThatBMatricesHaveSameSizes(rBs);
     CheckThatStressVectorsHaveSameSize(rStressVectors);
+    CheckThatTransposedMatrixVectorProductCanBeCalculated(rBs.front(), rStressVectors.front());
 }
 
 } // namespace

diff --git a/applications/GeoMechanicsApplication/tests/cpp_tests/test_equation_of_motion.cpp b/applications/GeoMechanicsApplication/tests/cpp_tests/test_equation_of_motion.cpp
@@ -244,16 +244,30 @@ KRATOS_TEST_CASE_IN_SUITE(CalculatingTheInternalForceVectorFailsWhenBMatricesHav
 KRATOS_TEST_CASE_IN_SUITE(CalculatingTheInternalForceVectorFailsWhenStressVectorsHaveDifferentSizes,
                           KratosGeoMechanicsFastSuiteWithoutKernel)
 {
-    const auto b_matrix                 = Matrix{ScalarMatrix{2, 8, 1.0}};
-    const auto b_matrices               = std::vector<Matrix>{b_matrix, b_matrix};
-    const auto stress_vectors           = std::vector<Vector>{ScalarVector{2, 1.0}, ScalarVector{3, 1.0}}; // Error: vectors have different sizes
+    const auto b_matrix   = Matrix{ScalarMatrix{2, 8, 1.0}};
+    const auto b_matrices = std::vector<Matrix>{b_matrix, b_matrix};
+    const auto stress_vectors = std::vector<Vector>{ScalarVector{2, 1.0}, ScalarVector{3, 1.0}}; // Error: vectors have different sizes
     const auto integration_coefficients = std::vector<double>{0.25, 0.4};
 
     KRATOS_EXPECT_EXCEPTION_IS_THROWN(
         GeoEquationOfMotionUtilities::CalculateInternalForceVector(b_matrices, stress_vectors, integration_coefficients),
         "Cannot calculate the internal force vector: stress vectors have different sizes")
 }
 
+KRATOS_TEST_CASE_IN_SUITE(CalculatingTheInternalForceVectorFailsWhenTheMatrixVectorProductCantBeComputed,
+                          KratosGeoMechanicsFastSuiteWithoutKernel)
+{
+    // Error: transpose of the B-matrix has more rows than the number of stress components
+    const auto b_matrix                 = Matrix{ScalarMatrix{3, 8, 1.0}};
+    const auto b_matrices               = std::vector<Matrix>{b_matrix, b_matrix};
+    const auto stress_vector            = Vector{ScalarVector{2, 1.0}};
+    const auto stress_vectors           = std::vector<Vector>{stress_vector, stress_vector};
+    const auto integration_coefficients = std::vector<double>{0.25, 0.4};
+
+    KRATOS_EXPECT_EXCEPTION_IS_THROWN(
+        GeoEquationOfMotionUtilities::CalculateInternalForceVector(b_matrices, stress_vectors, integration_coefficients), "Cannot calculate the internal force vector: matrix-vector product cannot be calculated due to size mismatch")
+}
+
 #endif
 
 } // namespace Kratos::Testing