Tensor migration by tuchpaul · Pull Request #2084 · 4C-multiphysics/4C

tuchpaul · 2026-06-16T14:04:37Z

Migrate some more fourth-order tensor products to the tensor API + fix two bugs

Part of the ongoing FourTensorOperations → tensor-API migration.

Added

kronecker_product(A) (symmetric box self product) and kronecker_product(A, B) (full)
adbc_tensor_product(A, B) = A_AD B_BC, non_symmetric_product(A, B) = A_AC B_DB
Conversion helpers: make_6x9 / make_9x6 / make_6x6_voigt_matrix_from_tensor

Fixed (`add_kronecker_tensor_product`, 6x9 overload)

Duplicate C(5,5) assignment (applied scalar_this twice)
C(5,8) accumulated from C(5,7) instead of C(5,8)

Refactored / removed

elast_hyper_get_derivs_of_elastic_right_cg_tensor now exemplarly uses the tensor API instead of the FourTensor/9x9-Voigt round-trip
contract_matrix_matrix → Core::LinAlg::ddot; removed unused multiply_four_tensor_matrix

Tests

New unit tests for all added ops: element-wise formula checks plus cross-checks against the legacy Voigt/9x9 implementations.

rjoussen

Thanks for tackling this!

rjoussen · 2026-06-17T10:16:05Z

+template <typename T>
+Core::LinAlg::SymmetricTensor<T, 3, 3, 3, 3> Core::LinAlg::FourTensorOperations::kronecker_product(
+    const Core::LinAlg::SymmetricTensor<T, 3, 3>& A)
+{
+  return 0.5 * (einsum_sym<"ik", "jl">(A, A) + einsum_sym<"il", "jk">(A, A));
+}


This operation returns a symmetric tensor regardless of whether A is symmetric. We should probably overload this or use SquareTensor concept directly here. This would also allow sizes different than 3.

rjoussen · 2026-06-17T10:27:43Z

+template <typename T>
+Core::LinAlg::Tensor<T, 3, 3, 3, 3> Core::LinAlg::FourTensorOperations::adbc_tensor_product(
+    const Core::LinAlg::Tensor<T, 3, 3>& A, const Core::LinAlg::Tensor<T, 3, 3>& B)
+{
+  return einsum<"ad", "bc">(A, B);
+}
+


Suggested change

template <typename T>

Core::LinAlg::Tensor<T, 3, 3, 3, 3> Core::LinAlg::FourTensorOperations::adbc_tensor_product(

const Core::LinAlg::Tensor<T, 3, 3>& A, const Core::LinAlg::Tensor<T, 3, 3>& B)

{

return einsum<"ad", "bc">(A, B);

}

Since this purely wraps einsum, which is clearer and easier to understand in my opinion, I would not implement this. Users can just use einsum directly

rjoussen · 2026-06-17T10:28:40Z

+template <typename T>
+Core::LinAlg::Tensor<T, 3, 3, 3, 3> Core::LinAlg::FourTensorOperations::non_symmetric_product(
+    const Core::LinAlg::Tensor<T, 3, 3>& A, const Core::LinAlg::Tensor<T, 3, 3>& B)
 {
-  if (clear_result_tensor) four_tensor_result.clear();
-  for (int i = 0; i < dim; ++i)
-  {
-    for (int j = 0; j < dim; ++j)
-    {
-      for (int k = 0; k < dim; ++k)
-      {
-        for (int l = 0; l < dim; ++l)
-        {
-          for (int m = 0; m < dim; ++m)
-          {  // C^ijkl = A^ijkm * B_m^l
-            four_tensor_result(i, j, k, l) += four_tensor(i, j, k, m) * matrix(m, l);
-          }
-        }
-      }
-    }
-  }
+  return einsum<"ac", "db">(A, B);


Same here, it is much clearer to use einsum directly, I think.

rjoussen · 2026-06-17T10:33:18Z

These new conversion methods will help a lot in the future, thanks for implementing them!

rjoussen · 2026-06-17T13:14:03Z

+  template <typename T>
+  Core::LinAlg::Tensor<T, 3, 3, 3, 3> kronecker_product(
+      const Core::LinAlg::Tensor<T, 3, 3>& A, const Core::LinAlg::Tensor<T, 3, 3>& B);
+


This is really a very special kind of kronecker product. For a normal kronecker product, I would probably expect something like $C = A\otimes B$ such that $C_{ijkl}=A_{ij}B_{kl}$. In other places, we refer to this exact product as Holzapfel product, but it's already inconsistent in the existing code.

Mainly, I am not sure whether we should introduce tensor helper functions with unclear semantics (unlike dot or ddot which are easy to understand). In contrast to the old Matrix api, specialized tensor products are easy to do on the fly when you need them, and the resulting code is probably easier to understand than having to look up what exactly a Kronecker product means here.

Yes, I totally agree. I found myself quite confused with the terminology and implementation of the kronecker product (which clashes sometimes).
Therefore, I am also happy with dropping the helper functions, as the einsum is actually self-documenting. But I'm also happy to hear other opinions, from e.g. @dragos-ana and @c-p-schmidt.

I agree as well. One could easily introduce "his" / "her" own helper functions directly in the cpp where the tensor products are required.

However, I also find the idea of having some products already available quite nice. Since the naming is confusing, maybe this is more explicitly addressed by naming them based on the indices involved, e.g. like in add_adbc_tensor_product. For the composed products such as the Holzapfel product, this gets slightly uglier, e.g., add_half_adbc_half_acbd_product.

dragos-ana

Looks very nice! I added some minor suggestions.
I am not entirely sure what is the best approach now regarding the added tensor product counterparts for some of the matrix products such as the Kronecker product.
The idea of having helper functions for these is nice, but this should be done with consistent, easily understandable names.

dragos-ana · 2026-06-19T18:47:28Z

+  template <typename T>
+  Core::LinAlg::Tensor<T, 3, 3, 3, 3> kronecker_product(
+      const Core::LinAlg::Tensor<T, 3, 3>& A, const Core::LinAlg::Tensor<T, 3, 3>& B);
+


I agree as well. One could easily introduce "his" / "her" own helper functions directly in the cpp where the tensor products are required.

However, I also find the idea of having some products already available quite nice. Since the naming is confusing, maybe this is more explicitly addressed by naming them based on the indices involved, e.g. like in add_adbc_tensor_product. For the composed products such as the Holzapfel product, this gets slightly uglier, e.g., add_half_adbc_half_acbd_product.

dragos-ana · 2026-06-19T19:00:59Z

+  /*!
+   * @brief Creates a 6x9 Voigt matrix from a 4th order tensor that has the minor symmetry in its
+   * first index pair (C_ijkl = C_jikl).


Really nice, the entire tensor conversion part. The conversion to 9x9 is missing as far as I see, but this is easy to implement based on your nice implementation here, if this will be required in the future.

dragos-ana · 2026-06-19T19:02:41Z

  double Mtheta_dev_sym_contract_Mtheta_dev_sym =
-      Core::LinAlg::FourTensorOperations::contract_matrix_matrix(
-          state_quantities.curr_Mtheta_dev_sym_M, state_quantities.curr_Mtheta_dev_sym_M);
+      Core::LinAlg::ddot(Core::LinAlg::make_tensor_view(state_quantities.curr_Mtheta_dev_sym_M),


Nice start to the tensor migration within the viscoplastic material ;)

dragos-ana · 2026-06-19T19:15:51Z

+    FOUR_C_EXPECT_NEAR(new_voigt, old_voigt, 1e-12);
+  }
+
+  TEST(TensorConversionTest, MakeStressLikeVoigtMatrixFromTensorTest)


Suggested change

TEST(TensorConversionTest, MakeStressLikeVoigtMatrixFromTensorTest)

TEST(TensorConversionTest, Make6x6StressLikeVoigtMatrixFromTensorTest)

Maybe slightly more descriptive and consistent with the other added tests?

tuchpaul added 2 commits June 16, 2026 13:02

Remove unused tensoroperations

5ab51eb

Add some tensor products and fix bugs

1b33278

tuchpaul self-assigned this Jun 16, 2026

tuchpaul requested a review from dragos-ana June 16, 2026 19:39

rjoussen reviewed Jun 17, 2026

View reviewed changes

dragos-ana reviewed Jun 19, 2026

View reviewed changes

	TEST(TensorConversionTest, MakeStressLikeVoigtMatrixFromTensorTest)
	TEST(TensorConversionTest, Make6x6StressLikeVoigtMatrixFromTensorTest)

Conversation

tuchpaul commented Jun 16, 2026

Migrate some more fourth-order tensor products to the tensor API + fix two bugs

Added

Fixed (add_kronecker_tensor_product, 6x9 overload)

Refactored / removed

Tests

Uh oh!

rjoussen left a comment

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

dragos-ana left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

3 participants

Fixed (`add_kronecker_tensor_product`, 6x9 overload)