Select Git revision
test_EmbeddedDiscreteFunctionMathFunctions1D.cpp 33.22 KiB
#include <catch2/catch_test_macros.hpp>
#include <catch2/matchers/catch_matchers_all.hpp>
#include <MeshDataBaseForTests.hpp>
#include <scheme/DiscreteFunctionP0.hpp>
#include <language/utils/EmbeddedIDiscreteFunctionMathFunctions.hpp>
#include <scheme/DiscreteFunctionP0Vector.hpp>
#include <scheme/DiscreteFunctionVariant.hpp>
#include <test_EmbeddedDiscreteFunctionMathFunctions.hpp>
// clazy:excludeall=non-pod-global-static
TEST_CASE("EmbeddedDiscreteFunctionVariantMathFunctions1D", "[scheme]")
{
constexpr size_t Dimension = 1;
using Rd = TinyVector<Dimension>;
std::array mesh_list = MeshDataBaseForTests::get().all1DMeshes();
using DiscreteFunctionR = DiscreteFunctionP0<Dimension, const double>;
using DiscreteFunctionR1 = DiscreteFunctionP0<Dimension, const TinyVector<1>>;
using DiscreteFunctionR2 = DiscreteFunctionP0<Dimension, const TinyVector<2>>;
using DiscreteFunctionR3 = DiscreteFunctionP0<Dimension, const TinyVector<3>>;
using DiscreteFunctionR1x1 = DiscreteFunctionP0<Dimension, const TinyMatrix<1>>;
using DiscreteFunctionR2x2 = DiscreteFunctionP0<Dimension, const TinyMatrix<2>>;
using DiscreteFunctionR3x3 = DiscreteFunctionP0<Dimension, const TinyMatrix<3>>;
using DiscreteFunctionVector = DiscreteFunctionP0Vector<Dimension, const double>;
for (const auto& named_mesh : mesh_list) {
SECTION(named_mesh.name())
{
auto mesh = named_mesh.mesh();
std::shared_ptr other_mesh =
std::make_shared<Mesh<Connectivity<Dimension>>>(mesh->shared_connectivity(), mesh->xr());
CellValue<const Rd> xj = MeshDataManager::instance().getMeshData(*mesh).xj();
CellValue<double> values = [=] {
CellValue<double> build_values{mesh->connectivity()};
parallel_for(
build_values.numberOfItems(),
PUGS_LAMBDA(const CellId cell_id) { build_values[cell_id] = 0.2 + std::cos(l2Norm(xj[cell_id])); });
return build_values;
}();
CellValue<double> positive_values = [=] {
CellValue<double> build_values{mesh->connectivity()};
parallel_for(
build_values.numberOfItems(),
PUGS_LAMBDA(const CellId cell_id) { build_values[cell_id] = 2 + std::sin(l2Norm(xj[cell_id])); });
return build_values;
}();
CellValue<double> bounded_values = [=] {
CellValue<double> build_values{mesh->connectivity()};
parallel_for(
build_values.numberOfItems(),
PUGS_LAMBDA(const CellId cell_id) { build_values[cell_id] = 0.9 * std::sin(l2Norm(xj[cell_id])); });
return build_values;
}();
std::shared_ptr p_u = std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionR(mesh, values));
std::shared_ptr p_other_mesh_u =
std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionR(other_mesh, values)); std::shared_ptr p_positive_u =
std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionR(mesh, positive_values));
std::shared_ptr p_bounded_u =
std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionR(mesh, bounded_values));
std::shared_ptr p_R1_u = [=] {
CellValue<TinyVector<1>> uj{mesh->connectivity()};
parallel_for(
uj.numberOfItems(), PUGS_LAMBDA(const CellId cell_id) { uj[cell_id][0] = 2 * xj[cell_id][0] + 1; });
return std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionR1(mesh, uj));
}();
std::shared_ptr p_R1_v = [=] {
CellValue<TinyVector<1>> vj{mesh->connectivity()};
parallel_for(
vj.numberOfItems(),
PUGS_LAMBDA(const CellId cell_id) { vj[cell_id][0] = xj[cell_id][0] * xj[cell_id][0] + 1; });
return std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionR1(mesh, vj));
}();
std::shared_ptr p_other_mesh_R1_u = std::make_shared<const DiscreteFunctionVariant>(
DiscreteFunctionR1(other_mesh, p_R1_u->get<DiscreteFunctionR1>().cellValues()));
constexpr auto to_2d = [&](const TinyVector<Dimension>& x) -> TinyVector<2> {
if constexpr (Dimension == 1) {
return TinyVector<2>{x[0], 1 + x[0] * x[0]};
} else if constexpr (Dimension == 2) {
return TinyVector<2>{x[0], x[1]};
} else if constexpr (Dimension == 3) {
return TinyVector<2>{x[0], x[1] + x[2]};
}
};
std::shared_ptr p_R2_u = [=] {
CellValue<TinyVector<2>> uj{mesh->connectivity()};
parallel_for(
uj.numberOfItems(), PUGS_LAMBDA(const CellId cell_id) {
const TinyVector<2> x = to_2d(xj[cell_id]);
uj[cell_id] = TinyVector<2>{2 * x[0] + 1, 1 - x[1]};
});
return std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionR2(mesh, uj));
}();
std::shared_ptr p_R2_v = [=] {
CellValue<TinyVector<2>> vj{mesh->connectivity()};
parallel_for(
vj.numberOfItems(), PUGS_LAMBDA(const CellId cell_id) {
const TinyVector<2> x = to_2d(xj[cell_id]);
vj[cell_id] = TinyVector<2>{x[0] * x[1] + 1, 2 * x[1]};
});
return std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionR2(mesh, vj));
}();
std::shared_ptr p_other_mesh_R2_u = std::make_shared<const DiscreteFunctionVariant>(
DiscreteFunctionR2(other_mesh, p_R2_u->get<DiscreteFunctionR2>().cellValues()));
constexpr auto to_3d = [&](const TinyVector<Dimension>& x) -> TinyVector<3> {
if constexpr (Dimension == 1) {
return TinyVector<3>{x[0], 1 + x[0] * x[0], 2 - x[0]};
} else if constexpr (Dimension == 2) {
return TinyVector<3>{x[0], x[1], x[0] + x[1]};
} else if constexpr (Dimension == 3) {
return TinyVector<3>{x[0], x[1], x[2]};
}
};
std::shared_ptr p_R3_u = [=] {
CellValue<TinyVector<3>> uj{mesh->connectivity()};
parallel_for(
uj.numberOfItems(), PUGS_LAMBDA(const CellId cell_id) {
const TinyVector<3> x = to_3d(xj[cell_id]);
uj[cell_id] = TinyVector<3>{2 * x[0] + 1, 1 - x[1] * x[2], x[0] + x[2]};
});
return std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionR3(mesh, uj));
}();
std::shared_ptr p_R3_v = [=] {
CellValue<TinyVector<3>> vj{mesh->connectivity()};
parallel_for(
vj.numberOfItems(), PUGS_LAMBDA(const CellId cell_id) {
const TinyVector<3> x = to_3d(xj[cell_id]);
vj[cell_id] = TinyVector<3>{x[0] * x[1] + 1, 2 * x[1], x[2] * x[0]};
});
return std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionR3(mesh, vj));
}();
std::shared_ptr p_other_mesh_R3_u = std::make_shared<const DiscreteFunctionVariant>(
DiscreteFunctionR3(other_mesh, p_R3_u->get<DiscreteFunctionR3>().cellValues()));
std::shared_ptr p_R1x1_u = [=] {
CellValue<TinyMatrix<1>> uj{mesh->connectivity()};
parallel_for(
uj.numberOfItems(),
PUGS_LAMBDA(const CellId cell_id) { uj[cell_id] = TinyMatrix<1>{2 * xj[cell_id][0] + 1}; });
return std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionR1x1(mesh, uj));
}();
std::shared_ptr p_R2x2_u = [=] {
CellValue<TinyMatrix<2>> uj{mesh->connectivity()};
parallel_for(
uj.numberOfItems(), PUGS_LAMBDA(const CellId cell_id) {
const TinyVector<2> x = to_2d(xj[cell_id]);
uj[cell_id] = TinyMatrix<2>{2 * x[0] + 1, 1 - x[1], //
2 * x[1], -x[0]};
});
return std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionR2x2(mesh, uj));
}();
std::shared_ptr p_R3x3_u = [=] {
CellValue<TinyMatrix<3>> uj{mesh->connectivity()};
parallel_for(
uj.numberOfItems(), PUGS_LAMBDA(const CellId cell_id) {
const TinyVector<3> x = to_3d(xj[cell_id]);
uj[cell_id] = TinyMatrix<3>{2 * x[0] + 1, 1 - x[1], 3, //
2 * x[1], -x[0], x[0] - x[1], //
3 * x[2] - x[1], x[1] - 2 * x[2], x[2] - x[0]};
});
return std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionR3x3(mesh, uj));
}();
std::shared_ptr p_Vector3_u = [=] {
CellArray<double> uj_vector{mesh->connectivity(), 3};
parallel_for(
uj_vector.numberOfItems(), PUGS_LAMBDA(const CellId cell_id) {
const TinyVector<3> x = to_3d(xj[cell_id]);
uj_vector[cell_id][0] = 2 * x[0] + 1;
uj_vector[cell_id][1] = 1 - x[1] * x[2];
uj_vector[cell_id][2] = x[0] + x[2];
});
return std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionVector(mesh, uj_vector));
}();
std::shared_ptr p_Vector3_v = [=] {
CellArray<double> vj_vector{mesh->connectivity(), 3};
parallel_for(
vj_vector.numberOfItems(), PUGS_LAMBDA(const CellId cell_id) {
const TinyVector<3> x = to_3d(xj[cell_id]);
vj_vector[cell_id][0] = x[0] * x[1] + 1;
vj_vector[cell_id][1] = 2 * x[1];
vj_vector[cell_id][2] = x[2] * x[0];
});
return std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionVector(mesh, vj_vector));
}();
std::shared_ptr p_Vector2_w = [=] {
CellArray<double> wj_vector{mesh->connectivity(), 2};
parallel_for(
wj_vector.numberOfItems(), PUGS_LAMBDA(const CellId cell_id) {
const TinyVector<3> x = to_3d(xj[cell_id]);
wj_vector[cell_id][0] = x[0] + x[1] * 2;
wj_vector[cell_id][1] = x[0] * x[1];
});
return std::make_shared<const DiscreteFunctionVariant>(DiscreteFunctionVector(mesh, wj_vector));
}();
SECTION("sqrt Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_positive_u, sqrt, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(sqrt(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("abs Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_u, abs, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(abs(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("sin Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_u, sin, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(sin(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("cos Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_u, cos, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(cos(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("tan Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_u, tan, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(tan(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("asin Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_bounded_u, asin, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(asin(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("acos Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_bounded_u, acos, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(acos(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("atan Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_bounded_u, atan, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(atan(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("sinh Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_u, sinh, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(sinh(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("cosh Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_u, cosh, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(cosh(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("tanh Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_u, tanh, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(tanh(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("asinh Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_positive_u, asinh, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(asinh(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("acosh Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_positive_u, acosh, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(acosh(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("atanh Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_bounded_u, atanh, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(atanh(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("exp Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_u, exp, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(exp(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("log Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_positive_u, log, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(log(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("atan2 Vh*Vh -> Vh")
{
CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_positive_u, p_bounded_u, atan2, //
DiscreteFunctionR, DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(atan2(p_u, p_R1_u), "error: incompatible operand types Vh(P0:R) and Vh(P0:R^1)");
REQUIRE_THROWS_WITH(atan2(p_R1_u, p_u), "error: incompatible operand types Vh(P0:R^1) and Vh(P0:R)");
REQUIRE_THROWS_WITH(atan2(p_u, p_other_mesh_u), "error: operands are defined on different meshes");
}
SECTION("atan2 Vh*R -> Vh")
{
CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_u, 3.6, atan2, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(atan2(p_R1_u, 2.1), "error: incompatible operand types Vh(P0:R^1) and R");
}
SECTION("atan2 R*Vh -> Vh")
{
CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION(2.4, p_u, atan2, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(atan2(2.1, p_R1_u), "error: incompatible operand types R and Vh(P0:R^1)");
}
SECTION("min Vh*Vh -> Vh")
{
CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_u, p_bounded_u, min, //
DiscreteFunctionR, DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(::min(p_u, p_other_mesh_u), "error: operands are defined on different meshes");
REQUIRE_THROWS_WITH(::min(p_u, p_R1_u), "error: incompatible operand types Vh(P0:R) and Vh(P0:R^1)");
REQUIRE_THROWS_WITH(::min(p_R1_u, p_u), "error: incompatible operand types Vh(P0:R^1) and Vh(P0:R)");
}
SECTION("min Vh*R -> Vh")
{
CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_u, 1.2, min, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(min(p_R1_u, 2.1), "error: incompatible operand types Vh(P0:R^1) and R");
}
SECTION("min R*Vh -> Vh")
{
CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION(0.4, p_u, min, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(min(3.1, p_R1_u), "error: incompatible operand types R and Vh(P0:R^1)");
}
SECTION("min Vh -> R")
{
REQUIRE(min(p_u) == min(p_u->get<DiscreteFunctionR>().cellValues()));
REQUIRE_THROWS_WITH(min(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("max Vh*Vh -> Vh")
{
CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_u, p_bounded_u, max, //
DiscreteFunctionR, DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(::max(p_u, p_other_mesh_u), "error: operands are defined on different meshes");
REQUIRE_THROWS_WITH(::max(p_u, p_R1_u), "error: incompatible operand types Vh(P0:R) and Vh(P0:R^1)");
REQUIRE_THROWS_WITH(::max(p_R1_u, p_u), "error: incompatible operand types Vh(P0:R^1) and Vh(P0:R)");
}
SECTION("max Vh*R -> Vh")
{
CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_u, 1.2, max, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(max(p_R1_u, 2.1), "error: incompatible operand types Vh(P0:R^1) and R");
}
SECTION("max Vh -> R") {
REQUIRE(max(p_u) == max(p_u->get<DiscreteFunctionR>().cellValues()));
REQUIRE_THROWS_WITH(max(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("max R*Vh -> Vh")
{
CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION(0.4, p_u, max, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(max(3.1, p_R1_u), "error: incompatible operand types R and Vh(P0:R^1)");
}
SECTION("pow Vh*Vh -> Vh")
{
CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_positive_u, p_bounded_u, pow, //
DiscreteFunctionR, DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(pow(p_u, p_other_mesh_u), "error: operands are defined on different meshes");
REQUIRE_THROWS_WITH(pow(p_u, p_R1_u), "error: incompatible operand types Vh(P0:R) and Vh(P0:R^1)");
REQUIRE_THROWS_WITH(pow(p_R1_u, p_u), "error: incompatible operand types Vh(P0:R^1) and Vh(P0:R)");
}
SECTION("pow Vh*R -> Vh")
{
CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_positive_u, 3.3, pow, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(pow(p_R1_u, 3.1), "error: incompatible operand types Vh(P0:R^1) and R");
}
SECTION("pow R*Vh -> Vh")
{
CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION(2.1, p_u, pow, //
DiscreteFunctionR, DiscreteFunctionR);
REQUIRE_THROWS_WITH(pow(3.1, p_R1_u), "error: incompatible operand types R and Vh(P0:R^1)");
}
SECTION("dot Vh*Vh -> Vh")
{
CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_R1_u, p_R1_v, dot, //
DiscreteFunctionR1, DiscreteFunctionR1, DiscreteFunctionR);
CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_R2_u, p_R2_v, dot, //
DiscreteFunctionR2, DiscreteFunctionR2, DiscreteFunctionR);
CHECK_EMBEDDED_VH2_TO_VH_FUNCTION_EVALUATION(p_R3_u, p_R3_v, dot, //
DiscreteFunctionR3, DiscreteFunctionR3, DiscreteFunctionR);
{
std::shared_ptr p_fuv = ::dot(p_Vector3_u, p_Vector3_v);
REQUIRE(p_fuv.use_count() > 0);
const DiscreteFunctionVector& u = p_Vector3_u->get<DiscreteFunctionVector>();
const DiscreteFunctionVector& v = p_Vector3_v->get<DiscreteFunctionVector>();
const DiscreteFunctionR& fuv = p_fuv->get<DiscreteFunctionR>();
bool is_same = true;
auto u_arrays = u.cellArrays();
auto v_arrays = v.cellArrays();
for (CellId cell_id = 0; cell_id < values.numberOfItems(); ++cell_id) {
using namespace std;
double dot_u_v = [&](auto&& a, auto&& b) {
double sum = 0;
for (size_t i = 0; i < a.size(); ++i) {
sum += a[i] * b[i];
}
return sum;
}(u_arrays[cell_id], v_arrays[cell_id]);
if (fuv[cell_id] != dot_u_v) {
is_same = false;
break;
}
}
REQUIRE(is_same);
}
REQUIRE_THROWS_WITH(dot(p_R1_u, p_other_mesh_R1_u), "error: operands are defined on different meshes");
REQUIRE_THROWS_WITH(dot(p_R2_u, p_other_mesh_R2_u), "error: operands are defined on different meshes");
REQUIRE_THROWS_WITH(dot(p_R3_u, p_other_mesh_R3_u), "error: operands are defined on different meshes");
REQUIRE_THROWS_WITH(dot(p_R1_u, p_R3_u), "error: incompatible operand types Vh(P0:R^1) and Vh(P0:R^3)");
REQUIRE_THROWS_WITH(dot(p_Vector3_u, p_Vector2_w), "error: operands have different dimension");
}
SECTION("det Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_R1x1_u, det, //
DiscreteFunctionR1x1, DiscreteFunctionR);
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_R2x2_u, det, //
DiscreteFunctionR2x2, DiscreteFunctionR);
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_R3x3_u, det, //
DiscreteFunctionR3x3, DiscreteFunctionR);
REQUIRE_THROWS_WITH(det(p_u), "error: invalid operand type Vh(P0:R)");
REQUIRE_THROWS_WITH(det(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
REQUIRE_THROWS_WITH(det(p_R2_u), "error: invalid operand type Vh(P0:R^2)");
REQUIRE_THROWS_WITH(det(p_R3_u), "error: invalid operand type Vh(P0:R^3)");
}
SECTION("trace Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_R1x1_u, trace, //
DiscreteFunctionR1x1, DiscreteFunctionR);
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_R2x2_u, trace, //
DiscreteFunctionR2x2, DiscreteFunctionR);
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_R3x3_u, trace, //
DiscreteFunctionR3x3, DiscreteFunctionR);
REQUIRE_THROWS_WITH(trace(p_u), "error: invalid operand type Vh(P0:R)");
REQUIRE_THROWS_WITH(trace(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
REQUIRE_THROWS_WITH(trace(p_R2_u), "error: invalid operand type Vh(P0:R^2)");
REQUIRE_THROWS_WITH(trace(p_R3_u), "error: invalid operand type Vh(P0:R^3)");
}
SECTION("inverse Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_R1x1_u, inverse, //
DiscreteFunctionR1x1, DiscreteFunctionR1x1);
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_R2x2_u, inverse, //
DiscreteFunctionR2x2, DiscreteFunctionR2x2);
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_R3x3_u, inverse, //
DiscreteFunctionR3x3, DiscreteFunctionR3x3);
REQUIRE_THROWS_WITH(inverse(p_u), "error: invalid operand type Vh(P0:R)");
REQUIRE_THROWS_WITH(inverse(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
REQUIRE_THROWS_WITH(inverse(p_R2_u), "error: invalid operand type Vh(P0:R^2)");
REQUIRE_THROWS_WITH(inverse(p_R3_u), "error: invalid operand type Vh(P0:R^3)");
}
SECTION("transpose Vh -> Vh")
{
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_R1x1_u, transpose, //
DiscreteFunctionR1x1, DiscreteFunctionR1x1);
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_R2x2_u, transpose, //
DiscreteFunctionR2x2, DiscreteFunctionR2x2);
CHECK_EMBEDDED_VH_TO_VH_FUNCTION_EVALUATION(p_R3x3_u, transpose, //
DiscreteFunctionR3x3, DiscreteFunctionR3x3);
REQUIRE_THROWS_WITH(transpose(p_u), "error: invalid operand type Vh(P0:R)");
REQUIRE_THROWS_WITH(transpose(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
REQUIRE_THROWS_WITH(transpose(p_R2_u), "error: invalid operand type Vh(P0:R^2)");
REQUIRE_THROWS_WITH(transpose(p_R3_u), "error: invalid operand type Vh(P0:R^3)");
}
SECTION("sum_of_Vh Vh -> Vh")
{
{
auto p_sum_components = sum_of_Vh_components(p_Vector3_u);
REQUIRE(p_sum_components.use_count() == 1);
const DiscreteFunctionR& sum_components = p_sum_components->get<DiscreteFunctionR>();
const DiscreteFunctionVector& vector3_u = p_Vector3_u->get<DiscreteFunctionVector>();
DiscreteFunctionP0<Dimension, double> direct_sum(mesh);
for (CellId cell_id = 0; cell_id < mesh->numberOfCells(); ++cell_id) {
double sum = 0;
for (size_t i = 0; i < vector3_u.size(); ++i) {
sum += vector3_u[cell_id][i];
}
direct_sum[cell_id] = sum;
}
bool is_same = true;
for (CellId cell_id = 0; cell_id < mesh->numberOfCells(); ++cell_id) {
if (sum_components[cell_id] != direct_sum[cell_id]) {
is_same = false;
break;
}
}
REQUIRE(is_same);
}
REQUIRE_THROWS_WITH(sum_of_Vh_components(p_R3_u), "error: invalid operand type Vh(P0:R^3)");
}
SECTION("vectorize (Vh) -> Vh")
{
{
std::shared_ptr p_vector3 = vectorize(std::vector{p_u, p_positive_u, p_bounded_u});
REQUIRE(p_vector3.use_count() == 1);
const DiscreteFunctionVector vector3 = p_vector3->get<DiscreteFunctionVector>();
const DiscreteFunctionR& u = p_u->get<DiscreteFunctionR>();
const DiscreteFunctionR& positive_u = p_positive_u->get<DiscreteFunctionR>();
const DiscreteFunctionR& bounded_u = p_bounded_u->get<DiscreteFunctionR>();
REQUIRE(vector3.size() == 3);
bool is_same = true;
for (CellId cell_id = 0; cell_id < mesh->numberOfCells(); ++cell_id) {
is_same &= (u[cell_id] == vector3[cell_id][0]);
is_same &= (positive_u[cell_id] == vector3[cell_id][1]);
is_same &= (bounded_u[cell_id] == vector3[cell_id][2]);
}
REQUIRE(is_same);
}
REQUIRE_THROWS_WITH(vectorize(std::vector{p_u, p_other_mesh_u}),
"error: discrete functions are not defined on the same mesh");
REQUIRE_THROWS_WITH(vectorize(std::vector{p_R1_u}), "error: invalid operand type Vh(P0:R^1)");
}
SECTION("dot Vh*Rd -> Vh") {
CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_R1_u, (TinyVector<1>{3}), dot, //
DiscreteFunctionR1, DiscreteFunctionR);
CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_R2_u, (TinyVector<2>{-6, 2}), dot, //
DiscreteFunctionR2, DiscreteFunctionR);
CHECK_EMBEDDED_VHxW_TO_VH_FUNCTION_EVALUATION(p_R3_u, (TinyVector<3>{-1, 5, 2}), dot, //
DiscreteFunctionR3, DiscreteFunctionR);
REQUIRE_THROWS_WITH(dot(p_R1_u, (TinyVector<2>{-6, 2})),
"error: incompatible operand types Vh(P0:R^1) and R^2");
REQUIRE_THROWS_WITH(dot(p_R2_u, (TinyVector<3>{-1, 5, 2})),
"error: incompatible operand types Vh(P0:R^2) and R^3");
REQUIRE_THROWS_WITH(dot(p_R3_u, (TinyVector<1>{-1})), "error: incompatible operand types Vh(P0:R^3) and R^1");
}
SECTION("dot Rd*Vh -> Vh")
{
CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION((TinyVector<1>{3}), p_R1_u, dot, //
DiscreteFunctionR1, DiscreteFunctionR);
CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION((TinyVector<2>{-6, 2}), p_R2_u, dot, //
DiscreteFunctionR2, DiscreteFunctionR);
CHECK_EMBEDDED_WxVH_TO_VH_FUNCTION_EVALUATION((TinyVector<3>{-1, 5, 2}), p_R3_u, dot, //
DiscreteFunctionR3, DiscreteFunctionR);
REQUIRE_THROWS_WITH(dot((TinyVector<2>{-6, 2}), p_R1_u),
"error: incompatible operand types R^2 and Vh(P0:R^1)");
REQUIRE_THROWS_WITH(dot((TinyVector<3>{-1, 5, 2}), p_R2_u),
"error: incompatible operand types R^3 and Vh(P0:R^2)");
REQUIRE_THROWS_WITH(dot((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()));
REQUIRE(sum_of<TinyVector<1>>(p_R1_u) == sum(p_R1_u->get<DiscreteFunctionR1>().cellValues()));
REQUIRE(sum_of<TinyVector<2>>(p_R2_u) == sum(p_R2_u->get<DiscreteFunctionR2>().cellValues()));
REQUIRE(sum_of<TinyVector<3>>(p_R3_u) == sum(p_R3_u->get<DiscreteFunctionR3>().cellValues()));
REQUIRE(sum_of<TinyMatrix<1>>(p_R1x1_u) == sum(p_R1x1_u->get<DiscreteFunctionR1x1>().cellValues()));
REQUIRE(sum_of<TinyMatrix<2>>(p_R2x2_u) == sum(p_R2x2_u->get<DiscreteFunctionR2x2>().cellValues()));
REQUIRE(sum_of<TinyMatrix<3>>(p_R3x3_u) == sum(p_R3x3_u->get<DiscreteFunctionR3x3>().cellValues()));
REQUIRE_THROWS_WITH(sum_of<TinyVector<1>>(p_u), "error: invalid operand type Vh(P0:R)");
REQUIRE_THROWS_WITH(sum_of<double>(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
REQUIRE_THROWS_WITH(sum_of<double>(p_R2_u), "error: invalid operand type Vh(P0:R^2)");
REQUIRE_THROWS_WITH(sum_of<double>(p_R3_u), "error: invalid operand type Vh(P0:R^3)");
REQUIRE_THROWS_WITH(sum_of<double>(p_R1x1_u), "error: invalid operand type Vh(P0:R^1x1)");
REQUIRE_THROWS_WITH(sum_of<double>(p_R2x2_u), "error: invalid operand type Vh(P0:R^2x2)");
REQUIRE_THROWS_WITH(sum_of<double>(p_R3x3_u), "error: invalid operand type Vh(P0:R^3x3)");
}
SECTION("integral_of_R* Vh -> R*")
{
auto integrate_locally = [&](const auto& cell_values) {
const auto& Vj = MeshDataManager::instance().getMeshData(*mesh).Vj();
using DataType = decltype(double{} * cell_values[CellId{0}]);
CellValue<DataType> local_integral{mesh->connectivity()};
parallel_for(
local_integral.numberOfItems(),
PUGS_LAMBDA(const CellId cell_id) { local_integral[cell_id] = Vj[cell_id] * cell_values[cell_id]; });
return local_integral;
};
REQUIRE(integral_of<double>(p_u) == sum(integrate_locally(p_u->get<DiscreteFunctionR>().cellValues())));
REQUIRE(integral_of<TinyVector<1>>(p_R1_u) ==
sum(integrate_locally(p_R1_u->get<DiscreteFunctionR1>().cellValues())));
REQUIRE(integral_of<TinyVector<2>>(p_R2_u) ==
sum(integrate_locally(p_R2_u->get<DiscreteFunctionR2>().cellValues())));
REQUIRE(integral_of<TinyVector<3>>(p_R3_u) ==
sum(integrate_locally(p_R3_u->get<DiscreteFunctionR3>().cellValues())));
REQUIRE(integral_of<TinyMatrix<1>>(p_R1x1_u) ==
sum(integrate_locally(p_R1x1_u->get<DiscreteFunctionR1x1>().cellValues())));
REQUIRE(integral_of<TinyMatrix<2>>(p_R2x2_u) == sum(integrate_locally(p_R2x2_u->get<DiscreteFunctionR2x2>().cellValues())));
REQUIRE(integral_of<TinyMatrix<3>>(p_R3x3_u) ==
sum(integrate_locally(p_R3x3_u->get<DiscreteFunctionR3x3>().cellValues())));
REQUIRE_THROWS_WITH(integral_of<TinyVector<1>>(p_u), "error: invalid operand type Vh(P0:R)");
REQUIRE_THROWS_WITH(integral_of<double>(p_R1_u), "error: invalid operand type Vh(P0:R^1)");
REQUIRE_THROWS_WITH(integral_of<double>(p_R2_u), "error: invalid operand type Vh(P0:R^2)");
REQUIRE_THROWS_WITH(integral_of<double>(p_R3_u), "error: invalid operand type Vh(P0:R^3)");
REQUIRE_THROWS_WITH(integral_of<double>(p_R1x1_u), "error: invalid operand type Vh(P0:R^1x1)");
REQUIRE_THROWS_WITH(integral_of<double>(p_R2x2_u), "error: invalid operand type Vh(P0:R^2x2)");
REQUIRE_THROWS_WITH(integral_of<double>(p_R3x3_u), "error: invalid operand type Vh(P0:R^3x3)");
}
}
}
}