From 04468e71800a835b7daf6912e5871b6249b58b12 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Miguel=20Mas=C3=B3?= <miguel.maso@upc.edu>
Date: Tue, 14 Jul 2026 19:41:17 +0200
Subject: [PATCH] added contractions

---
 Project.toml                    |  2 +-
 src/Exports.jl                  |  4 ++-
 src/TensorAlgebra/Operations.jl | 60 ++++++++++++++++++++++++++++++++-
 3 files changed, 63 insertions(+), 3 deletions(-)

diff --git a/Project.toml b/Project.toml
index 90d81dc..f40c403 100644
--- a/Project.toml
+++ b/Project.toml
@@ -1,7 +1,7 @@
 name = "HyperFEM"
 uuid = "ddbf423a-3df4-409e-9685-3dc416cae46c"
 version = "0.0.4"
-authors = ["MultiSimo_Group"]
+authors = ["MultiSimOLab"]
 
 [deps]
 AbstractTrees = "1520ce14-60c1-5f80-bbc7-55ef81b5835c"
diff --git a/src/Exports.jl b/src/Exports.jl
index 1f74631..0b613cf 100644
--- a/src/Exports.jl
+++ b/src/Exports.jl
@@ -10,12 +10,14 @@ end
 @publish TensorAlgebra VectorValue  # Reexport from Gridap for convenience
 @publish TensorAlgebra (*)
 @publish TensorAlgebra (×ᵢ⁴)
+@publish TensorAlgebra (⊗₁²)
 @publish TensorAlgebra (⊗₁₂³)
 @publish TensorAlgebra (⊗₁₃²)
 @publish TensorAlgebra (⊗₁²³)
 @publish TensorAlgebra (⊗₁₃²⁴)
 @publish TensorAlgebra (⊗₁₂³⁴)
-@publish TensorAlgebra (⊗₁²)
+@publish TensorAlgebra (⊗₁₂₃⁴)
+@publish TensorAlgebra (⊗₁₂₃₄²⁴)
 @publish TensorAlgebra logreg
 @publish TensorAlgebra Box
 @publish TensorAlgebra Ellipsoid
diff --git a/src/TensorAlgebra/Operations.jl b/src/TensorAlgebra/Operations.jl
index a7dc491..3d15693 100644
--- a/src/TensorAlgebra/Operations.jl
+++ b/src/TensorAlgebra/Operations.jl
@@ -27,6 +27,14 @@ function Gridap.TensorValues.outer(A::VectorValue{D,Float64}, B::VectorValue{D,F
   return (A ⊗₁² B)
 end
 
+function Gridap.TensorValues.outer(A::TensorValue{4,2,Float64}, B::VectorValue{2,Float64})
+  return (A ⊗₁₂₃⁴ B)
+end
+
+function Gridap.TensorValues.outer(A::TensorValue{9,3,Float64}, B::VectorValue{3,Float64})
+  return (A ⊗₁₂₃⁴ B)
+end
+
 
 """
     ⊗₁²(A::VectorValue{D}, B::VectorValue{D})::TensorValue{D,D}
@@ -121,7 +129,7 @@ end
 
 
 """
-    ⊗₁²³(A::TensorValue{D}, B::VectorValue{D})::TensorValue{D,D*D}
+    ⊗₁₂³(A::TensorValue{D}, B::VectorValue{D})::TensorValue{D,D*D}
 
 Outer product of a second-order and first-order tensors (matrix and vector),
 returning a third-order tensor represented in a `D x D²` flattened matrix using combined indices.
@@ -156,6 +164,29 @@ returning a third-order tensor represented in a `D x D²` flattened matrix using
 end
 
 
+"""
+    ⊗₁₂₃⁴(A::TensorValue{D²,D}, B::TensorValue{D})::TensorValue{D,D*D}
+
+Outer product of a third-order and first-order tensors (tensor and vector),
+returning a fourth-order tensor represented in a `D² x D²` flattened matrix using combined indices.
+"""
+@generated function (⊗₁₂₃⁴)(A::TensorValue{D²,D}, V::VectorValue{D}) where {D,D²}
+  @assert D*D == D² "Third- and first-order tensors size mismatch with $D² × $D and $D"
+  str = ""
+  for l in 1:D
+    for k in 1:D
+      for j in 1:D
+        for i in 1:D
+          a = _flat_idx(i,j,D)
+          str *= "A[$a,$k]*V[$l],"
+        end
+      end
+    end
+  end
+  Meta.parse("TensorValue{D*D,D*D}($str)")
+end
+
+
 function (×ᵢ⁴)(A::TensorValue{3,3,Float64})
 
   TensorValue(0.0, 0.0, 0.0, 0.0, A[9], -A[8], 0.0, -A[6], A[5], 0.0, 0.0, 0.0, -A[9],
@@ -516,3 +547,30 @@ The operation follows the **index contraction pattern**, where addition is perfo
   end
   Meta.parse("TensorValue{D²}($str)")
 end
+
+
+"""
+    contraction_IJKL_JL(A::TensorValue{D*D}, H::TensorValue{D})::TensorValue{D*D}
+
+Performs a tensor contraction between a fourth-order tensor (represented as a `D² × D²` matrix in flattened index notation)
+and a second-order tensor (of size `D × D`).
+The operation follows the **index contraction pattern**, where addition is performed for repeated indices.
+"""
+@generated function contraction_IJKL_JL(H::TensorValue{D²}, A::TensorValue{D}) where {D, D²}
+  @assert D*D == D² "Fourth- and Second-order tensors size mismatch"
+  str = ""
+  for i in 1:D
+    for k in 1:D
+      for j in 1:D
+        for l in 1:D
+          a = _flat_idx(i,j,k,l,D)
+          str *= "+H[$a]*A[$j,$l]"
+        end
+      end
+      str *= ","
+    end
+  end
+  Meta.parse("TensorValue{D}($str)")
+end
+
+(⊗₁₂₃₄²⁴) = contraction_IJKL_JL