diff --git a/src/language/modules/MathFunctionRegisterForVh.cpp b/src/language/modules/MathFunctionRegisterForVh.cpp
index 38dc08b6e4bec82be8cab2666080563bd7cb9736..b31fe47d0eba87ff2e7f0823a530bf170e410d6c 100644
--- a/src/language/modules/MathFunctionRegisterForVh.cpp
+++ b/src/language/modules/MathFunctionRegisterForVh.cpp
@@ -212,6 +212,76 @@ MathFunctionRegisterForVh::MathFunctionRegisterForVh(SchemeModule& scheme_module
));
+ scheme_module._addBuiltinFunction("tensorProduct", std::function(
+
+ [](std::shared_ptr<const DiscreteFunctionVariant> a,
+ std::shared_ptr<const DiscreteFunctionVariant> b)
+ -> std::shared_ptr<const DiscreteFunctionVariant> {
+ return tensorProduct(a, b);
+ }
+
+ ));
+
+ scheme_module._addBuiltinFunction("tensorProduct",
+ std::function(
+
+ [](std::shared_ptr<const DiscreteFunctionVariant> a,
+ const TinyVector<1> b) -> std::shared_ptr<const DiscreteFunctionVariant> {
+ return tensorProduct(a, b);
+ }
+
+ ));
+
+ scheme_module._addBuiltinFunction("tensorProduct",
+ std::function(
+
+ [](std::shared_ptr<const DiscreteFunctionVariant> a,
+ const TinyVector<2> b) -> std::shared_ptr<const DiscreteFunctionVariant> {
+ return tensorProduct(a, b);
+ }
+
+ ));
+
+ scheme_module._addBuiltinFunction("tensorProduct",
+ std::function(
+
+ [](std::shared_ptr<const DiscreteFunctionVariant> a,
+ const TinyVector<3>& b) -> std::shared_ptr<const DiscreteFunctionVariant> {
+ return tensorProduct(a, b);
+ }
+
+ ));
+
+ scheme_module._addBuiltinFunction("tensorProduct",
+ std::function(
+
+ [](const TinyVector<1> a, std::shared_ptr<const DiscreteFunctionVariant> b)
+ -> std::shared_ptr<const DiscreteFunctionVariant> {
+ return tensorProduct(a, b);
+ }
+
+ ));
+
+ scheme_module._addBuiltinFunction("tensorProduct",
+ std::function(
+
+ [](const TinyVector<2> a, std::shared_ptr<const DiscreteFunctionVariant> b)
+ -> std::shared_ptr<const DiscreteFunctionVariant> {
+ return tensorProduct(a, b);
+ }
+
+ ));
+
+ scheme_module._addBuiltinFunction("tensorProduct",
+ std::function(
+
+ [](const TinyVector<3>& a, std::shared_ptr<const DiscreteFunctionVariant> b)
+ -> std::shared_ptr<const DiscreteFunctionVariant> {
+ return tensorProduct(a, b);
+ }
+
+ ));
+
scheme_module._addBuiltinFunction("det", std::function(
[](std::shared_ptr<const DiscreteFunctionVariant> A)
diff --git a/src/language/modules/MathModule.cpp b/src/language/modules/MathModule.cpp
index 4a1e7b35b76a4d9b8b05ebb1355f3cdca4f751d9..af1dfeafb4ac2bdd22894cb3ad1fde26b1ed0d8c 100644
--- a/src/language/modules/MathModule.cpp
+++ b/src/language/modules/MathModule.cpp
@@ -71,6 +71,21 @@ MathModule::MathModule()
return dot(x, y);
}));
+ this->_addBuiltinFunction("tensorProduct",
+ std::function([](const TinyVector<1> x, const TinyVector<1> y) -> TinyMatrix<1> {
+ return tensorProduct(x, y);
+ }));
+
+ this->_addBuiltinFunction("tensorProduct",
+ std::function([](const TinyVector<2> x, const TinyVector<2> y) -> TinyMatrix<2> {
+ return tensorProduct(x, y);
+ }));
+
+ this->_addBuiltinFunction("tensorProduct",
+ std::function([](const TinyVector<3>& x, const TinyVector<3>& y) -> TinyMatrix<3> {
+ return tensorProduct(x, y);
+ }));
+
this->_addBuiltinFunction("det", std::function([](const TinyMatrix<1> A) -> double { return det(A); }));
this->_addBuiltinFunction("det", std::function([](const TinyMatrix<2>& A) -> double { return det(A); }));
diff --git a/src/language/utils/EmbeddedDiscreteFunctionMathFunctions.cpp b/src/language/utils/EmbeddedDiscreteFunctionMathFunctions.cpp
index 715de1be2019070e517cb30639c199f721258cd2..961cc0e4092ba57c001866a52d7aed113e6b0f5d 100644
--- a/src/language/utils/EmbeddedDiscreteFunctionMathFunctions.cpp
+++ b/src/language/utils/EmbeddedDiscreteFunctionMathFunctions.cpp
@@ -315,6 +315,110 @@ template std::shared_ptr<const DiscreteFunctionVariant> dot(const TinyVector<2>&
template std::shared_ptr<const DiscreteFunctionVariant> dot(const TinyVector<3>&,
const std::shared_ptr<const DiscreteFunctionVariant>&);
+std::shared_ptr<const DiscreteFunctionVariant>
+tensorProduct(const std::shared_ptr<const DiscreteFunctionVariant>& f_v,
+ const std::shared_ptr<const DiscreteFunctionVariant>& g_v)
+{
+ if (not hasSameMesh({f_v, g_v})) {
+ throw NormalError("operands are defined on different meshes");
+ }
+
+ return std::visit(
+ [&](auto&& f, auto&& g) -> std::shared_ptr<DiscreteFunctionVariant> {
+ using TypeOfF = std::decay_t<decltype(f)>;
+ using TypeOfG = std::decay_t<decltype(g)>;
+ if constexpr (not std::is_same_v<TypeOfF, TypeOfG>) {
+ throw NormalError(EmbeddedDiscreteFunctionUtils::incompatibleOperandTypes(f, g));
+ } else {
+ using DataType = std::decay_t<typename TypeOfF::data_type>;
+ if constexpr (is_discrete_function_P0_v<TypeOfF>) {
+ if constexpr (is_tiny_vector_v<DataType>) {
+ return std::make_shared<DiscreteFunctionVariant>(tensorProduct(f, g));
+ } else {
+ throw NormalError(EmbeddedDiscreteFunctionUtils::invalidOperandType(f));
+ }
+ } else {
+ throw NormalError(EmbeddedDiscreteFunctionUtils::invalidOperandType(f));
+ }
+ }
+ },
+ f_v->discreteFunction(), g_v->discreteFunction());
+}
+
+template <size_t VectorDimension>
+std::shared_ptr<const DiscreteFunctionVariant>
+tensorProduct(const std::shared_ptr<const DiscreteFunctionVariant>& f, const TinyVector<VectorDimension>& a)
+{
+ return std::visit(
+ [&](auto&& discrete_function0) -> std::shared_ptr<DiscreteFunctionVariant> {
+ using DiscreteFunction0Type = std::decay_t<decltype(discrete_function0)>;
+ if constexpr (is_discrete_function_P0_v<DiscreteFunction0Type>) {
+ using DataType = std::decay_t<typename DiscreteFunction0Type::data_type>;
+ if constexpr (is_tiny_vector_v<DataType>) {
+ if constexpr (std::is_same_v<DataType, TinyVector<VectorDimension>>) {
+ return std::make_shared<DiscreteFunctionVariant>(tensorProduct<DataType>(discrete_function0, a));
+ } else {
+ throw NormalError(EmbeddedDiscreteFunctionUtils::incompatibleOperandTypes(f, a));
+ }
+ } else {
+ throw NormalError(EmbeddedDiscreteFunctionUtils::invalidOperandType(f));
+ }
+ } else {
+ throw NormalError(EmbeddedDiscreteFunctionUtils::incompatibleOperandTypes(f, a));
+ }
+ },
+ f->discreteFunction());
+}
+
+template <size_t VectorDimension>
+std::shared_ptr<const DiscreteFunctionVariant>
+tensorProduct(const TinyVector<VectorDimension>& a, const std::shared_ptr<const DiscreteFunctionVariant>& f)
+{
+ return std::visit(
+ [&](auto&& discrete_function0) -> std::shared_ptr<DiscreteFunctionVariant> {
+ using DiscreteFunction0Type = std::decay_t<decltype(discrete_function0)>;
+ if constexpr (is_discrete_function_P0_v<DiscreteFunction0Type>) {
+ using DataType = std::decay_t<typename DiscreteFunction0Type::data_type>;
+ if constexpr (is_tiny_vector_v<DataType>) {
+ if constexpr (std::is_same_v<DataType, TinyVector<VectorDimension>>) {
+ return std::make_shared<DiscreteFunctionVariant>(tensorProduct<DataType>(a, discrete_function0));
+ } else {
+ throw NormalError(EmbeddedDiscreteFunctionUtils::incompatibleOperandTypes(a, f));
+ }
+ } else {
+ throw NormalError(EmbeddedDiscreteFunctionUtils::invalidOperandType(f));
+ }
+ } else {
+ throw NormalError(EmbeddedDiscreteFunctionUtils::incompatibleOperandTypes(a, f));
+ }
+ },
+ f->discreteFunction());
+}
+
+template std::shared_ptr<const DiscreteFunctionVariant> tensorProduct(
+ const std::shared_ptr<const DiscreteFunctionVariant>&,
+ const TinyVector<1>&);
+
+template std::shared_ptr<const DiscreteFunctionVariant> tensorProduct(
+ const std::shared_ptr<const DiscreteFunctionVariant>&,
+ const TinyVector<2>&);
+
+template std::shared_ptr<const DiscreteFunctionVariant> tensorProduct(
+ const std::shared_ptr<const DiscreteFunctionVariant>&,
+ const TinyVector<3>&);
+
+template std::shared_ptr<const DiscreteFunctionVariant> tensorProduct(
+ const TinyVector<1>&,
+ const std::shared_ptr<const DiscreteFunctionVariant>&);
+
+template std::shared_ptr<const DiscreteFunctionVariant> tensorProduct(
+ const TinyVector<2>&,
+ const std::shared_ptr<const DiscreteFunctionVariant>&);
+
+template std::shared_ptr<const DiscreteFunctionVariant> tensorProduct(
+ const TinyVector<3>&,
+ const std::shared_ptr<const DiscreteFunctionVariant>&);
+
std::shared_ptr<const DiscreteFunctionVariant>
det(const std::shared_ptr<const DiscreteFunctionVariant>& A)
{
diff --git a/src/language/utils/EmbeddedDiscreteFunctionMathFunctions.hpp b/src/language/utils/EmbeddedDiscreteFunctionMathFunctions.hpp
index 66b784921da372cc7dfa2d616b94d6fbfe403ee3..bcda26113a4fac6d72e7b54edfe0f0e8d8d74b2d 100644
--- a/src/language/utils/EmbeddedDiscreteFunctionMathFunctions.hpp
+++ b/src/language/utils/EmbeddedDiscreteFunctionMathFunctions.hpp
@@ -69,6 +69,17 @@ template <size_t VectorDimension>
std::shared_ptr<const DiscreteFunctionVariant> dot(const TinyVector<VectorDimension>&,
const std::shared_ptr<const DiscreteFunctionVariant>&);
+std::shared_ptr<const DiscreteFunctionVariant> tensorProduct(const std::shared_ptr<const DiscreteFunctionVariant>&,
+ const std::shared_ptr<const DiscreteFunctionVariant>&);
+
+template <size_t VectorDimension>
+std::shared_ptr<const DiscreteFunctionVariant> tensorProduct(const std::shared_ptr<const DiscreteFunctionVariant>&,
+ const TinyVector<VectorDimension>&);
+
+template <size_t VectorDimension>
+std::shared_ptr<const DiscreteFunctionVariant> tensorProduct(const TinyVector<VectorDimension>&,
+ const std::shared_ptr<const DiscreteFunctionVariant>&);
+
std::shared_ptr<const DiscreteFunctionVariant> det(const std::shared_ptr<const DiscreteFunctionVariant>&);
std::shared_ptr<const DiscreteFunctionVariant> trace(const std::shared_ptr<const DiscreteFunctionVariant>&);
diff --git a/src/scheme/DiscreteFunctionP0.hpp b/src/scheme/DiscreteFunctionP0.hpp
index df43f52905837036d0ef843fbc5635f0d1274c6d..ea0ec56aeeebb28e4af381ae352a20c9abd43bd4 100644
--- a/src/scheme/DiscreteFunctionP0.hpp
+++ b/src/scheme/DiscreteFunctionP0.hpp
@@ -606,6 +606,51 @@ class DiscreteFunctionP0
return result;
}
+ template <typename LHSDataType>
+ PUGS_INLINE friend DiscreteFunctionP0<decltype(tensorProduct(LHSDataType{}, DataType{}))>
+ tensorProduct(const LHSDataType& u, const DiscreteFunctionP0& v)
+ {
+ static_assert(is_tiny_vector_v<std::decay_t<DataType>>);
+ static_assert(is_tiny_vector_v<std::decay_t<LHSDataType>>);
+ Assert(v.m_cell_values.isBuilt());
+ DiscreteFunctionP0<decltype(tensorProduct(LHSDataType{}, DataType{}))> result{v.m_mesh_v};
+ parallel_for(
+ v.m_mesh_v->numberOfCells(), PUGS_LAMBDA(CellId cell_id) { result[cell_id] = tensorProduct(u, v[cell_id]); });
+
+ return result;
+ }
+
+ template <typename RHSDataType>
+ PUGS_INLINE friend DiscreteFunctionP0<decltype(tensorProduct(DataType{}, RHSDataType{}))>
+ tensorProduct(const DiscreteFunctionP0& u, const RHSDataType& v)
+ {
+ static_assert(is_tiny_vector_v<std::decay_t<DataType>>);
+ static_assert(is_tiny_vector_v<std::decay_t<RHSDataType>>);
+ Assert(u.m_cell_values.isBuilt());
+ DiscreteFunctionP0<decltype(tensorProduct(DataType{}, RHSDataType{}))> result{u.m_mesh_v};
+ parallel_for(
+ u.m_mesh_v->numberOfCells(), PUGS_LAMBDA(CellId cell_id) { result[cell_id] = tensorProduct(u[cell_id], v); });
+
+ return result;
+ }
+
+ template <typename DataType2>
+ PUGS_INLINE friend DiscreteFunctionP0<decltype(tensorProduct(DataType{}, DataType2{}))>
+ tensorProduct(const DiscreteFunctionP0& u, const DiscreteFunctionP0<DataType2>& v)
+ {
+ static_assert(is_tiny_vector_v<std::decay_t<DataType>>);
+ static_assert(is_tiny_vector_v<std::decay_t<DataType2>>);
+ Assert(u.m_cell_values.isBuilt());
+ Assert(v.m_cell_values.isBuilt());
+ Assert(u.m_mesh_v->id() == v.m_mesh_v->id());
+ DiscreteFunctionP0<decltype(tensorProduct(DataType{}, DataType2{}))> result{u.m_mesh_v};
+ parallel_for(
+ u.m_mesh_v->numberOfCells(),
+ PUGS_LAMBDA(CellId cell_id) { result[cell_id] = tensorProduct(u[cell_id], v[cell_id]); });
+
+ return result;
+ }
+
PUGS_INLINE friend double
min(const DiscreteFunctionP0& f)
{
diff --git a/tests/test_BuiltinFunctionProcessor.cpp b/tests/test_BuiltinFunctionProcessor.cpp
index b678612fe2c338b641b1e66c6ccce9ae65454b8f..f600bfa4408c81cf447a8665880fb3040ea81e75 100644
--- a/tests/test_BuiltinFunctionProcessor.cpp
+++ b/tests/test_BuiltinFunctionProcessor.cpp
@@ -362,6 +362,40 @@ let s:R, s = dot(x,y);
CHECK_BUILTIN_FUNCTION_EVALUATION_RESULT(data, "s", dot(TinyVector<3>{-2, 3, 4}, TinyVector<3>{4, 3, 5}));
}
+ { // dot
+ tested_function_set.insert("tensorProduct:R^1*R^1");
+ std::string_view data = R"(
+import math;
+let x:R^1, x = -2;
+let y:R^1, y = 4;
+let s:R^1x1, s = tensorProduct(x,y);
+)";
+ CHECK_BUILTIN_FUNCTION_EVALUATION_RESULT(data, "s", TinyMatrix<1>{-2 * 4});
+ }
+
+ { // dot
+ tested_function_set.insert("tensorProduct:R^2*R^2");
+ std::string_view data = R"(
+import math;
+let x:R^2, x = [-2, 3];
+let y:R^2, y = [4, 3];
+let s:R^2x2, s = tensorProduct(x,y);
+)";
+ CHECK_BUILTIN_FUNCTION_EVALUATION_RESULT(data, "s", tensorProduct(TinyVector<2>{-2, 3}, TinyVector<2>{4, 3}));
+ }
+
+ { // dot
+ tested_function_set.insert("tensorProduct:R^3*R^3");
+ std::string_view data = R"(
+import math;
+let x:R^3, x = [-2, 3, 4];
+let y:R^3, y = [4, 3, 5];
+let s:R^3x3, s = tensorProduct(x,y);
+)";
+ CHECK_BUILTIN_FUNCTION_EVALUATION_RESULT(data, "s",
+ tensorProduct(TinyVector<3>{-2, 3, 4}, TinyVector<3>{4, 3, 5}));
+ }
+
{ // det
tested_function_set.insert("det:R^1x1");
std::string_view data = R"(
diff --git a/tests/test_DiscreteFunctionP0.cpp b/tests/test_DiscreteFunctionP0.cpp
index 5b1133b4e6c60290b06aa90e5c68b64a1c1b4973..14c96b0db712daa220951f407e66eb1f774f126e 100644
--- a/tests/test_DiscreteFunctionP0.cpp
+++ b/tests/test_DiscreteFunctionP0.cpp
@@ -2729,6 +2729,53 @@ TEST_CASE("DiscreteFunctionP0", "[scheme]")
CHECK_STD_BINARY_MATH_FUNCTION_WITH_LHS_VALUE(u, vh, dot);
}
+ SECTION("tensorProduct(uh,hv)")
+ {
+ DiscreteFunctionP0<TinyVector<2>> uh{mesh};
+ parallel_for(
+ mesh->numberOfCells(), PUGS_LAMBDA(const CellId cell_id) {
+ const double x = xj[cell_id][0];
+ uh[cell_id] = TinyVector<2>{x + 1, 2 * x - 3};
+ });
+
+ DiscreteFunctionP0<TinyVector<2>> vh{mesh};
+ parallel_for(
+ mesh->numberOfCells(), PUGS_LAMBDA(const CellId cell_id) {
+ const double x = xj[cell_id][0];
+ vh[cell_id] = TinyVector<2>{2.3 * x, 1 - x};
+ });
+
+ CHECK_STD_BINARY_MATH_FUNCTION(uh, vh, tensorProduct);
+ }
+
+ SECTION("tensorProduct(uh,v)")
+ {
+ DiscreteFunctionP0<TinyVector<2>> uh{mesh};
+ parallel_for(
+ mesh->numberOfCells(), PUGS_LAMBDA(const CellId cell_id) {
+ const double x = xj[cell_id][0];
+ uh[cell_id] = TinyVector<2>{x + 1, 2 * x - 3};
+ });
+
+ const TinyVector<2> v{1, 2};
+
+ CHECK_STD_BINARY_MATH_FUNCTION_WITH_RHS_VALUE(uh, v, tensorProduct);
+ }
+
+ SECTION("tensorProduct(u,hv)")
+ {
+ const TinyVector<2> u{3, -2};
+
+ DiscreteFunctionP0<TinyVector<2>> vh{mesh};
+ parallel_for(
+ mesh->numberOfCells(), PUGS_LAMBDA(const CellId cell_id) {
+ const double x = xj[cell_id][0];
+ vh[cell_id] = TinyVector<2>{2.3 * x, 1 - x};
+ });
+
+ CHECK_STD_BINARY_MATH_FUNCTION_WITH_LHS_VALUE(u, vh, tensorProduct);
+ }
+
SECTION("scalar sum")
{
const CellValue<const double> cell_value = positive_function.cellValues();
diff --git a/tests/test_EmbeddedDiscreteFunctionMathFunctions1D.cpp b/tests/test_EmbeddedDiscreteFunctionMathFunctions1D.cpp
index 5843a8ae766549e0828551bae8d33ff8c595935d..41f3c0d2233b1884051f0eafd5df637ce37ed7e0 100644
--- a/tests/test_EmbeddedDiscreteFunctionMathFunctions1D.cpp
+++ b/tests/test_EmbeddedDiscreteFunctionMathFunctions1D.cpp
@@ -502,6 +502,29 @@ TEST_CASE("EmbeddedDiscreteFunctionVariantMathFunctions1D", "[scheme]")
REQUIRE_THROWS_WITH(dot(p_Vector3_u, p_Vector2_w), "error: operands have different dimension");
}
+ SECTION("tensorProduct Vh*Vh -> Vh")
+ {
+ CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_R1_u, p_R1_v, tensorProduct, //
+ DiscreteFunctionR1, DiscreteFunctionR1, DiscreteFunctionR1x1);
+ CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_R2_u, p_R2_v, tensorProduct, //
+ DiscreteFunctionR2, DiscreteFunctionR2, DiscreteFunctionR2x2);
+ CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_R3_u, p_R3_v, tensorProduct, //
+ DiscreteFunctionR3, DiscreteFunctionR3, DiscreteFunctionR3x3);
+
+ REQUIRE_THROWS_WITH(tensorProduct(p_R1_u, p_other_mesh_R1_u),
+ "error: operands are defined on different meshes");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R2_u, p_other_mesh_R2_u),
+ "error: operands are defined on different meshes");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R3_u, p_other_mesh_R3_u),
+ "error: operands are defined on different meshes");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R1_u, p_R3_u),
+ "error: incompatible operand types Vh(P0:R^1) and Vh(P0:R^3)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_u, p_u), "error: invalid operand type Vh(P0:R)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R1x1_u, p_R1x1_u), "error: invalid operand type Vh(P0:R^1x1)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R2x2_u, p_R2x2_u), "error: invalid operand type Vh(P0:R^2x2)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R3x3_u, p_R3x3_u), "error: invalid operand type Vh(P0:R^3x3)");
+ }
+
SECTION("det Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_R1x1_u, det, //
@@ -675,6 +698,48 @@ TEST_CASE("EmbeddedDiscreteFunctionVariantMathFunctions1D", "[scheme]")
REQUIRE_THROWS_WITH(dot((TinyVector<1>{-1}), p_R3_u), "error: incompatible operand types R^1 and Vh(P0:R^3)");
}
+ SECTION("tensorProduct Vh*Rd -> Vh")
+ {
+ CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_R1_u, (TinyVector<1>{3}), tensorProduct, //
+ DiscreteFunctionR1, DiscreteFunctionR1x1);
+ CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_R2_u, (TinyVector<2>{-6, 2}), tensorProduct, //
+ DiscreteFunctionR2, DiscreteFunctionR2x2);
+ CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_R3_u, (TinyVector<3>{-1, 5, 2}), tensorProduct, //
+ DiscreteFunctionR3, DiscreteFunctionR3x3);
+
+ REQUIRE_THROWS_WITH(tensorProduct(p_u, TinyVector<2>{1, 2}), "error: invalid operand type Vh(P0:R)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_Vector3_u, TinyVector<2>{1, 2}),
+ "error: incompatible operand types Vh(P0Vector:R) and R^2");
+
+ REQUIRE_THROWS_WITH(tensorProduct(p_R1_u, (TinyVector<2>{-6, 2})),
+ "error: incompatible operand types Vh(P0:R^1) and R^2");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R2_u, (TinyVector<3>{-1, 5, 2})),
+ "error: incompatible operand types Vh(P0:R^2) and R^3");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R3_u, (TinyVector<1>{-1})),
+ "error: incompatible operand types Vh(P0:R^3) and R^1");
+ }
+
+ SECTION("tensorProduct Rd*Vh -> Vh")
+ {
+ CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION((TinyVector<1>{3}), p_R1_u, tensorProduct, //
+ DiscreteFunctionR1, DiscreteFunctionR1x1);
+ CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION((TinyVector<2>{-6, 2}), p_R2_u, tensorProduct, //
+ DiscreteFunctionR2, DiscreteFunctionR2x2);
+ CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION((TinyVector<3>{-1, 5, 2}), p_R3_u, tensorProduct, //
+ DiscreteFunctionR3, DiscreteFunctionR3x3);
+
+ REQUIRE_THROWS_WITH(tensorProduct(TinyVector<2>{1, 2}, p_u), "error: invalid operand type Vh(P0:R)");
+ REQUIRE_THROWS_WITH(tensorProduct(TinyVector<2>{1, 2}, p_Vector3_u),
+ "error: incompatible operand types R^2 and Vh(P0Vector:R)");
+
+ REQUIRE_THROWS_WITH(tensorProduct((TinyVector<2>{-6, 2}), p_R1_u),
+ "error: incompatible operand types R^2 and Vh(P0:R^1)");
+ REQUIRE_THROWS_WITH(tensorProduct((TinyVector<3>{-1, 5, 2}), p_R2_u),
+ "error: incompatible operand types R^3 and Vh(P0:R^2)");
+ REQUIRE_THROWS_WITH(tensorProduct((TinyVector<1>{-1}), p_R3_u),
+ "error: incompatible operand types R^1 and Vh(P0:R^3)");
+ }
+
SECTION("sum_of_R* Vh -> R*")
{
REQUIRE(sum_of<double>(p_u) == sum(p_u->get<DiscreteFunctionR>().cellValues()));
diff --git a/tests/test_EmbeddedDiscreteFunctionMathFunctions2D.cpp b/tests/test_EmbeddedDiscreteFunctionMathFunctions2D.cpp
index 11582b9820d9e2f6d910773b2084517a47069b37..919d4f8280e764d263129f70a849064867bb6f2a 100644
--- a/tests/test_EmbeddedDiscreteFunctionMathFunctions2D.cpp
+++ b/tests/test_EmbeddedDiscreteFunctionMathFunctions2D.cpp
@@ -504,6 +504,30 @@ TEST_CASE("EmbeddedDiscreteFunctionVariantMathFunctions2D", "[scheme]")
REQUIRE_THROWS_WITH(dot(p_Vector3_u, p_Vector2_w), "error: operands have different dimension");
}
+ SECTION("tensorProduct Vh*Vh -> Vh")
+ {
+ CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_R1_u, p_R1_v, tensorProduct, //
+ DiscreteFunctionR1, DiscreteFunctionR1, DiscreteFunctionR1x1);
+ CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_R2_u, p_R2_v, tensorProduct, //
+ DiscreteFunctionR2, DiscreteFunctionR2, DiscreteFunctionR2x2);
+ CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_R3_u, p_R3_v, tensorProduct, //
+ DiscreteFunctionR3, DiscreteFunctionR3, DiscreteFunctionR3x3);
+
+ REQUIRE_THROWS_WITH(tensorProduct(p_R1_u, p_other_mesh_R1_u),
+ "error: operands are defined on different meshes");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R2_u, p_other_mesh_R2_u),
+ "error: operands are defined on different meshes");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R3_u, p_other_mesh_R3_u),
+ "error: operands are defined on different meshes");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R1_u, p_R3_u),
+ "error: incompatible operand types Vh(P0:R^1) and Vh(P0:R^3)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_u, p_u), "error: invalid operand type Vh(P0:R)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R1x1_u, p_R1x1_u), "error: invalid operand type Vh(P0:R^1x1)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R2x2_u, p_R2x2_u), "error: invalid operand type Vh(P0:R^2x2)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R3x3_u, p_R3x3_u), "error: invalid operand type Vh(P0:R^3x3)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_Vector3_u, p_Vector2_w), "error: invalid operand type Vh(P0Vector:R)");
+ }
+
SECTION("det Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_R1x1_u, det, //
@@ -677,6 +701,48 @@ TEST_CASE("EmbeddedDiscreteFunctionVariantMathFunctions2D", "[scheme]")
REQUIRE_THROWS_WITH(dot((TinyVector<1>{-1}), p_R3_u), "error: incompatible operand types R^1 and Vh(P0:R^3)");
}
+ SECTION("tensorProduct Vh*Rd -> Vh")
+ {
+ CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_R1_u, (TinyVector<1>{3}), tensorProduct, //
+ DiscreteFunctionR1, DiscreteFunctionR1x1);
+ CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_R2_u, (TinyVector<2>{-6, 2}), tensorProduct, //
+ DiscreteFunctionR2, DiscreteFunctionR2x2);
+ CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_R3_u, (TinyVector<3>{-1, 5, 2}), tensorProduct, //
+ DiscreteFunctionR3, DiscreteFunctionR3x3);
+
+ REQUIRE_THROWS_WITH(tensorProduct(p_u, TinyVector<2>{1, 2}), "error: invalid operand type Vh(P0:R)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_Vector3_u, TinyVector<2>{1, 2}),
+ "error: incompatible operand types Vh(P0Vector:R) and R^2");
+
+ REQUIRE_THROWS_WITH(tensorProduct(p_R1_u, (TinyVector<2>{-6, 2})),
+ "error: incompatible operand types Vh(P0:R^1) and R^2");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R2_u, (TinyVector<3>{-1, 5, 2})),
+ "error: incompatible operand types Vh(P0:R^2) and R^3");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R3_u, (TinyVector<1>{-1})),
+ "error: incompatible operand types Vh(P0:R^3) and R^1");
+ }
+
+ SECTION("tensorProduct Rd*Vh -> Vh")
+ {
+ CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION((TinyVector<1>{3}), p_R1_u, tensorProduct, //
+ DiscreteFunctionR1, DiscreteFunctionR1x1);
+ CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION((TinyVector<2>{-6, 2}), p_R2_u, tensorProduct, //
+ DiscreteFunctionR2, DiscreteFunctionR2x2);
+ CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION((TinyVector<3>{-1, 5, 2}), p_R3_u, tensorProduct, //
+ DiscreteFunctionR3, DiscreteFunctionR3x3);
+
+ REQUIRE_THROWS_WITH(tensorProduct(TinyVector<2>{1, 2}, p_u), "error: invalid operand type Vh(P0:R)");
+ REQUIRE_THROWS_WITH(tensorProduct(TinyVector<2>{1, 2}, p_Vector3_u),
+ "error: incompatible operand types R^2 and Vh(P0Vector:R)");
+
+ REQUIRE_THROWS_WITH(tensorProduct((TinyVector<2>{-6, 2}), p_R1_u),
+ "error: incompatible operand types R^2 and Vh(P0:R^1)");
+ REQUIRE_THROWS_WITH(tensorProduct((TinyVector<3>{-1, 5, 2}), p_R2_u),
+ "error: incompatible operand types R^3 and Vh(P0:R^2)");
+ REQUIRE_THROWS_WITH(tensorProduct((TinyVector<1>{-1}), p_R3_u),
+ "error: incompatible operand types R^1 and Vh(P0:R^3)");
+ }
+
SECTION("sum_of_R* Vh -> R*")
{
REQUIRE(sum_of<double>(p_u) == sum(p_u->get<DiscreteFunctionR>().cellValues()));
diff --git a/tests/test_EmbeddedDiscreteFunctionMathFunctions3D.cpp b/tests/test_EmbeddedDiscreteFunctionMathFunctions3D.cpp
index edca0a4eb9c9e015a7e0ed31fd74ca068d4258c5..99994e9d0a902add718dee6149c3c040765a7106 100644
--- a/tests/test_EmbeddedDiscreteFunctionMathFunctions3D.cpp
+++ b/tests/test_EmbeddedDiscreteFunctionMathFunctions3D.cpp
@@ -502,6 +502,30 @@ TEST_CASE("EmbeddedDiscreteFunctionVariantMathFunctions3D", "[scheme]")
REQUIRE_THROWS_WITH(dot(p_Vector3_u, p_Vector2_w), "error: operands have different dimension");
}
+ SECTION("tensorProduct Vh*Vh -> Vh")
+ {
+ CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_R1_u, p_R1_v, tensorProduct, //
+ DiscreteFunctionR1, DiscreteFunctionR1, DiscreteFunctionR1x1);
+ CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_R2_u, p_R2_v, tensorProduct, //
+ DiscreteFunctionR2, DiscreteFunctionR2, DiscreteFunctionR2x2);
+ CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_R3_u, p_R3_v, tensorProduct, //
+ DiscreteFunctionR3, DiscreteFunctionR3, DiscreteFunctionR3x3);
+
+ REQUIRE_THROWS_WITH(tensorProduct(p_R1_u, p_other_mesh_R1_u),
+ "error: operands are defined on different meshes");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R2_u, p_other_mesh_R2_u),
+ "error: operands are defined on different meshes");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R3_u, p_other_mesh_R3_u),
+ "error: operands are defined on different meshes");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R1_u, p_R3_u),
+ "error: incompatible operand types Vh(P0:R^1) and Vh(P0:R^3)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_u, p_u), "error: invalid operand type Vh(P0:R)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R1x1_u, p_R1x1_u), "error: invalid operand type Vh(P0:R^1x1)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R2x2_u, p_R2x2_u), "error: invalid operand type Vh(P0:R^2x2)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R3x3_u, p_R3x3_u), "error: invalid operand type Vh(P0:R^3x3)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_Vector3_u, p_Vector2_w), "error: invalid operand type Vh(P0Vector:R)");
+ }
+
SECTION("det Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_R1x1_u, det, //
@@ -675,6 +699,48 @@ TEST_CASE("EmbeddedDiscreteFunctionVariantMathFunctions3D", "[scheme]")
REQUIRE_THROWS_WITH(dot((TinyVector<1>{-1}), p_R3_u), "error: incompatible operand types R^1 and Vh(P0:R^3)");
}
+ SECTION("tensorProduct Vh*Rd -> Vh")
+ {
+ CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_R1_u, (TinyVector<1>{3}), tensorProduct, //
+ DiscreteFunctionR1, DiscreteFunctionR1x1);
+ CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_R2_u, (TinyVector<2>{-6, 2}), tensorProduct, //
+ DiscreteFunctionR2, DiscreteFunctionR2x2);
+ CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_R3_u, (TinyVector<3>{-1, 5, 2}), tensorProduct, //
+ DiscreteFunctionR3, DiscreteFunctionR3x3);
+
+ REQUIRE_THROWS_WITH(tensorProduct(p_u, TinyVector<2>{1, 2}), "error: invalid operand type Vh(P0:R)");
+ REQUIRE_THROWS_WITH(tensorProduct(p_Vector3_u, TinyVector<2>{1, 2}),
+ "error: incompatible operand types Vh(P0Vector:R) and R^2");
+
+ REQUIRE_THROWS_WITH(tensorProduct(p_R1_u, (TinyVector<2>{-6, 2})),
+ "error: incompatible operand types Vh(P0:R^1) and R^2");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R2_u, (TinyVector<3>{-1, 5, 2})),
+ "error: incompatible operand types Vh(P0:R^2) and R^3");
+ REQUIRE_THROWS_WITH(tensorProduct(p_R3_u, (TinyVector<1>{-1})),
+ "error: incompatible operand types Vh(P0:R^3) and R^1");
+ }
+
+ SECTION("tensorProduct Rd*Vh -> Vh")
+ {
+ CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION((TinyVector<1>{3}), p_R1_u, tensorProduct, //
+ DiscreteFunctionR1, DiscreteFunctionR1x1);
+ CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION((TinyVector<2>{-6, 2}), p_R2_u, tensorProduct, //
+ DiscreteFunctionR2, DiscreteFunctionR2x2);
+ CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION((TinyVector<3>{-1, 5, 2}), p_R3_u, tensorProduct, //
+ DiscreteFunctionR3, DiscreteFunctionR3x3);
+
+ REQUIRE_THROWS_WITH(tensorProduct(TinyVector<2>{1, 2}, p_u), "error: invalid operand type Vh(P0:R)");
+ REQUIRE_THROWS_WITH(tensorProduct(TinyVector<2>{1, 2}, p_Vector3_u),
+ "error: incompatible operand types R^2 and Vh(P0Vector:R)");
+
+ REQUIRE_THROWS_WITH(tensorProduct((TinyVector<2>{-6, 2}), p_R1_u),
+ "error: incompatible operand types R^2 and Vh(P0:R^1)");
+ REQUIRE_THROWS_WITH(tensorProduct((TinyVector<3>{-1, 5, 2}), p_R2_u),
+ "error: incompatible operand types R^3 and Vh(P0:R^2)");
+ REQUIRE_THROWS_WITH(tensorProduct((TinyVector<1>{-1}), p_R3_u),
+ "error: incompatible operand types R^1 and Vh(P0:R^3)");
+ }
+
SECTION("sum_of_R* Vh -> R*")
{
REQUIRE(sum_of<double>(p_u) == sum(p_u->get<DiscreteFunctionR>().cellValues()));
diff --git a/tests/test_MathModule.cpp b/tests/test_MathModule.cpp
index f299d1e3678217d4fc6bcb4c08e417264b22a252..89a480e1ac4993b3cdbbd6702f5e2130a1872ec0 100644
--- a/tests/test_MathModule.cpp
+++ b/tests/test_MathModule.cpp
@@ -13,7 +13,7 @@ TEST_CASE("MathModule", "[language]")
MathModule math_module;
const auto& name_builtin_function = math_module.getNameBuiltinFunctionMap();
- REQUIRE(name_builtin_function.size() == 42);
+ REQUIRE(name_builtin_function.size() == 45);
SECTION("Z -> N")
{
@@ -458,6 +458,63 @@ TEST_CASE("MathModule", "[language]")
}
}
+ SECTION("(R^d, R^d) -> R^dxd")
+ {
+ SECTION("tensorProduct:R^1*R^1")
+ {
+ TinyVector<1> arg0{3};
+ TinyVector<1> arg1{2};
+
+ DataVariant arg0_variant = arg0;
+ DataVariant arg1_variant = arg1;
+
+ auto i_function = name_builtin_function.find("tensorProduct:R^1*R^1");
+ REQUIRE(i_function != name_builtin_function.end());
+
+ IBuiltinFunctionEmbedder& function_embedder = *i_function->second;
+ DataVariant result_variant = function_embedder.apply({arg0_variant, arg1_variant});
+
+ const TinyMatrix<1> result = tensorProduct(arg0, arg1);
+ REQUIRE(std::get<TinyMatrix<1>>(result_variant) == result);
+ }
+
+ SECTION("tensorProduct:R^2*R^2")
+ {
+ TinyVector<2> arg0{3, 2};
+ TinyVector<2> arg1{-2, 5};
+
+ DataVariant arg0_variant = arg0;
+ DataVariant arg1_variant = arg1;
+
+ auto i_function = name_builtin_function.find("tensorProduct:R^2*R^2");
+ REQUIRE(i_function != name_builtin_function.end());
+
+ IBuiltinFunctionEmbedder& function_embedder = *i_function->second;
+ DataVariant result_variant = function_embedder.apply({arg0_variant, arg1_variant});
+
+ const TinyMatrix<2> result = tensorProduct(arg0, arg1);
+ REQUIRE(std::get<TinyMatrix<2>>(result_variant) == result);
+ }
+
+ SECTION("tensorProduct:R^3*R^3")
+ {
+ TinyVector<3> arg0{3, 2, 4};
+ TinyVector<3> arg1{-2, 5, 2};
+
+ DataVariant arg0_variant = arg0;
+ DataVariant arg1_variant = arg1;
+
+ auto i_function = name_builtin_function.find("tensorProduct:R^3*R^3");
+ REQUIRE(i_function != name_builtin_function.end());
+
+ IBuiltinFunctionEmbedder& function_embedder = *i_function->second;
+ DataVariant result_variant = function_embedder.apply({arg0_variant, arg1_variant});
+
+ const TinyMatrix<3> result = tensorProduct(arg0, arg1);
+ REQUIRE(std::get<TinyMatrix<3>>(result_variant) == result);
+ }
+ }
+
SECTION("R^dxd -> double")
{
SECTION("det:R^1x1")