From 0adccd02c00ae4407ee09e65c4e3d5855a94941d Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?St=C3=A9phane=20Del=20Pino?= <stephane.delpino44@gmail.com>
Date: Mon, 8 Jul 2024 16:26:26 +0200
Subject: [PATCH] Add column weighting as a preconditioner for LS problems
Also add pseudo-options for preconditioning and row weighting
---
src/scheme/PolynomialReconstruction.cpp | 67 ++++++++++---
src/scheme/PolynomialReconstruction.hpp | 9 +-
src/scheme/test_reconstruction.cpp | 2 +-
tests/test_PolynomialReconstruction.cpp | 121 ++++++++++++------------
4 files changed, 125 insertions(+), 74 deletions(-)
diff --git a/src/scheme/PolynomialReconstruction.cpp b/src/scheme/PolynomialReconstruction.cpp
index 788c6df64..e3fe01930 100644
--- a/src/scheme/PolynomialReconstruction.cpp
+++ b/src/scheme/PolynomialReconstruction.cpp
@@ -103,6 +103,8 @@ template <MeshConcept MeshType>
[[nodiscard]] std::vector<std::shared_ptr<const DiscreteFunctionDPkVariant>>
PolynomialReconstruction::_build(
const size_t degree,
+ const bool preconditioning,
+ const bool row_weighting,
const std::shared_ptr<const MeshType>& p_mesh,
const std::vector<std::shared_ptr<const DiscreteFunctionVariant>>& discrete_function_variant_list) const
{
@@ -194,11 +196,13 @@ PolynomialReconstruction::_build(
SmallArray<SmallMatrix<double>> A_pool(Kokkos::DefaultExecutionSpace::concurrency());
SmallArray<SmallMatrix<double>> B_pool(Kokkos::DefaultExecutionSpace::concurrency());
+ SmallArray<SmallArray<double>> G_pool(Kokkos::DefaultExecutionSpace::concurrency());
SmallArray<SmallMatrix<double>> X_pool(Kokkos::DefaultExecutionSpace::concurrency());
for (size_t i = 0; i < A_pool.size(); ++i) {
A_pool[i] = SmallMatrix<double>(max_stencil_size, basis_dimension - 1);
B_pool[i] = SmallMatrix<double>(max_stencil_size, number_of_columns);
+ G_pool[i] = SmallArray<double>(basis_dimension - 1);
X_pool[i] = SmallMatrix<double>(basis_dimension - 1, number_of_columns);
}
@@ -278,19 +282,55 @@ PolynomialReconstruction::_build(
}
}
- for (size_t i = 0; i < stencil_cell_list.size(); ++i) {
- const CellId cell_i_id = stencil_cell_list[i];
- const double weight = 1. / l2Norm(xj[cell_i_id] - Xj);
- for (size_t l = 0; l < basis_dimension - 1; ++l) {
- A(i, l) *= weight;
- }
- for (size_t l = 0; l < number_of_columns; ++l) {
- B(i, l) *= weight;
+ if (row_weighting) {
+ // Add row weighting (give more importance to cells that are
+ // closer to j)
+ for (size_t i = 0; i < stencil_cell_list.size(); ++i) {
+ const CellId cell_i_id = stencil_cell_list[i];
+ const double weight = 1. / l2Norm(xj[cell_i_id] - Xj);
+ for (size_t l = 0; l < basis_dimension - 1; ++l) {
+ A(i, l) *= weight;
+ }
+ for (size_t l = 0; l < number_of_columns; ++l) {
+ B(i, l) *= weight;
+ }
}
}
const SmallMatrix<double>& X = X_pool[t];
- Givens::solveCollectionInPlace(A, X, B);
+
+ if (preconditioning) {
+ // Add column weighting preconditioning (increase the presition)
+ SmallArray<double>& G = G_pool[t];
+
+ for (size_t l = 0; l < basis_dimension - 1; ++l) {
+ double g = 0;
+ for (size_t i = 0; i < stencil_cell_list.size(); ++i) {
+ const double Ail = A(i, l);
+
+ g += Ail * Ail;
+ }
+ G[l] = std::sqrt(g);
+ }
+
+ for (size_t l = 0; l < basis_dimension - 1; ++l) {
+ const double Gl = G[l];
+ for (size_t i = 0; i < stencil_cell_list.size(); ++i) {
+ A(i, l) *= Gl;
+ }
+ }
+
+ Givens::solveCollectionInPlace(A, X, B);
+
+ for (size_t l = 0; l < basis_dimension - 1; ++l) {
+ const double Gl = G[l];
+ for (size_t i = 0; i < number_of_columns; ++i) {
+ X(l, i) *= Gl;
+ }
+ }
+ } else {
+ Givens::solveCollectionInPlace(A, X, B);
+ }
column_begin = 0;
for (size_t i_dpk_variant = 0; i_dpk_variant < mutable_discrete_function_dpk_variant_list.size();
@@ -408,6 +448,8 @@ PolynomialReconstruction::_build(
std::vector<std::shared_ptr<const DiscreteFunctionDPkVariant>>
PolynomialReconstruction::build(
const size_t degree,
+ const bool preconditioning,
+ const bool row_weighting,
const std::vector<std::shared_ptr<const DiscreteFunctionVariant>>& discrete_function_variant_list) const
{
if (not hasSameMesh(discrete_function_variant_list)) {
@@ -416,6 +458,9 @@ PolynomialReconstruction::build(
auto mesh_v = getCommonMesh(discrete_function_variant_list);
- return std::visit([&](auto&& p_mesh) { return this->_build(degree, p_mesh, discrete_function_variant_list); },
- mesh_v->variant());
+ return std::visit(
+ [&](auto&& p_mesh) {
+ return this->_build(degree, preconditioning, row_weighting, p_mesh, discrete_function_variant_list);
+ },
+ mesh_v->variant());
}
diff --git a/src/scheme/PolynomialReconstruction.hpp b/src/scheme/PolynomialReconstruction.hpp
index 78f7bc5ba..74bb82785 100644
--- a/src/scheme/PolynomialReconstruction.hpp
+++ b/src/scheme/PolynomialReconstruction.hpp
@@ -29,6 +29,8 @@ class PolynomialReconstruction
template <MeshConcept MeshType>
[[nodiscard]] std::vector<std::shared_ptr<const DiscreteFunctionDPkVariant>> _build(
const size_t degree,
+ const bool preconditioning,
+ const bool row_weighting,
const std::shared_ptr<const MeshType>& p_mesh,
const std::vector<std::shared_ptr<const DiscreteFunctionVariant>>& discrete_function_variant_list) const;
@@ -265,11 +267,14 @@ class PolynomialReconstruction
[[nodiscard]] std::vector<std::shared_ptr<const DiscreteFunctionDPkVariant>> build(
const size_t degree,
+ const bool preconditioning,
+ const bool row_weighting,
const std::vector<std::shared_ptr<const DiscreteFunctionVariant>>& discrete_function_variant_list) const;
+#warning add reconstruction descriptor/option to handle the options to be used for reconstruction
template <typename... DiscreteFunctionT>
[[nodiscard]] std::vector<std::shared_ptr<const DiscreteFunctionDPkVariant>>
- build(const size_t degree, DiscreteFunctionT... input) const
+ build(const size_t degree, const bool preconditioning, const bool row_weighting, DiscreteFunctionT... input) const
{
std::vector<std::shared_ptr<const DiscreteFunctionVariant>> variant_vector;
auto convert = [&variant_vector](auto&& df) {
@@ -289,7 +294,7 @@ class PolynomialReconstruction
};
(convert(std::forward<DiscreteFunctionT>(input)), ...);
- return this->build(degree, variant_vector);
+ return this->build(degree, preconditioning, row_weighting, variant_vector);
}
PolynomialReconstruction() {}
diff --git a/src/scheme/test_reconstruction.cpp b/src/scheme/test_reconstruction.cpp
index 87b9bcc16..1267aabec 100644
--- a/src/scheme/test_reconstruction.cpp
+++ b/src/scheme/test_reconstruction.cpp
@@ -87,7 +87,7 @@ test_reconstruction(const std::vector<std::shared_ptr<const DiscreteFunctionVari
PolynomialReconstruction reconstruction;
Timer t;
for (size_t i = 0; i < 50; ++i) {
- auto res = reconstruction.build(1, discrete_function_variant_list);
+ auto res = reconstruction.build(1, true, true, discrete_function_variant_list);
}
t.pause();
std::cout << "t = " << t << '\n';
diff --git a/tests/test_PolynomialReconstruction.cpp b/tests/test_PolynomialReconstruction.cpp
index e019cc016..ada21453a 100644
--- a/tests/test_PolynomialReconstruction.cpp
+++ b/tests/test_PolynomialReconstruction.cpp
@@ -43,7 +43,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
parallel_for(
mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id) { fh[cell_id] = R_affine(xj[cell_id]); });
- auto reconstructions = PolynomialReconstruction{}.build(1, fh);
+ auto reconstructions = PolynomialReconstruction{}.build(1, true, true, fh);
auto dpk_fh = reconstructions[0]->get<DiscreteFunctionDPk<1, const double>>();
@@ -53,7 +53,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error =
std::max(max_mean_error, std::abs(dpk_fh[cell_id](xj[cell_id]) - R_affine(xj[cell_id])));
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-14));
}
{
@@ -64,7 +64,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_slope_error = std::max(max_slope_error, std::abs(reconstructed_slope - 1.7));
}
- REQUIRE(max_slope_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_slope_error) == Catch::Approx(0).margin(1E-14));
}
}
}
@@ -92,7 +92,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
parallel_for(
mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id) { uh[cell_id] = R3_affine(xj[cell_id]); });
- auto reconstructions = PolynomialReconstruction{}.build(1, uh);
+ auto reconstructions = PolynomialReconstruction{}.build(1, true, true, uh);
auto dpk_uh = reconstructions[0]->get<DiscreteFunctionDPk<1, const R3>>();
@@ -102,7 +102,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error =
std::max(max_mean_error, l2Norm(dpk_uh[cell_id](xj[cell_id]) - R3_affine(xj[cell_id])));
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-14));
}
{
@@ -113,7 +113,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_slope_error = std::max(max_slope_error, l2Norm(reconstructed_slope - R3{1.7, -0.6, 3.1}));
}
- REQUIRE(max_slope_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_slope_error) == Catch::Approx(0).margin(1E-14));
}
}
}
@@ -143,7 +143,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
parallel_for(
mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id) { Ah[cell_id] = R3x3_affine(xj[cell_id]); });
- auto reconstructions = PolynomialReconstruction{}.build(1, Ah);
+ auto reconstructions = PolynomialReconstruction{}.build(1, true, true, Ah);
auto dpk_Ah = reconstructions[0]->get<DiscreteFunctionDPk<1, const R3x3>>();
@@ -153,7 +153,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error =
std::max(max_mean_error, frobeniusNorm(dpk_Ah[cell_id](xj[cell_id]) - R3x3_affine(xj[cell_id])));
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-13));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-13));
}
{
@@ -169,7 +169,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_slope_error = std::max(max_slope_error, //
frobeniusNorm(reconstructed_slope - slops));
}
- REQUIRE(max_slope_error == Catch::Approx(0).margin(1E-13));
+ REQUIRE(parallel::allReduceMax(max_slope_error) == Catch::Approx(0).margin(1E-13));
}
}
}
@@ -200,7 +200,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
}
});
- auto reconstructions = PolynomialReconstruction{}.build(1, Vh);
+ auto reconstructions = PolynomialReconstruction{}.build(1, true, true, Vh);
auto dpk_Vh = reconstructions[0]->get<DiscreteFunctionDPkVector<1, const double>>();
@@ -212,7 +212,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error = std::max(max_mean_error, std::abs(dpk_Vh(cell_id, i)(xj[cell_id]) - vector[i]));
}
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-13));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-13));
}
{
@@ -226,8 +226,8 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_slope_error = std::max(max_slope_error, std::abs(reconstructed_slope - slope[i]));
}
- REQUIRE(max_slope_error == Catch::Approx(0).margin(1E-13));
}
+ REQUIRE(parallel::allReduceMax(max_slope_error) == Catch::Approx(0).margin(1E-13));
}
}
}
@@ -287,9 +287,10 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
}
});
- auto reconstructions = PolynomialReconstruction{}.build(1, std::make_shared<DiscreteFunctionVariant>(fh),
- uh, std::make_shared<DiscreteFunctionP0<R3x3>>(Ah),
- DiscreteFunctionVariant(Vh));
+ auto reconstructions =
+ PolynomialReconstruction{}.build(1, true, true, std::make_shared<DiscreteFunctionVariant>(fh), uh,
+ std::make_shared<DiscreteFunctionP0<R3x3>>(Ah),
+ DiscreteFunctionVariant(Vh));
auto dpk_fh = reconstructions[0]->get<DiscreteFunctionDPk<1, const double>>();
@@ -299,7 +300,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error =
std::max(max_mean_error, std::abs(dpk_fh[cell_id](xj[cell_id]) - R_affine(xj[cell_id])));
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-14));
}
{
@@ -310,7 +311,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_slope_error = std::max(max_slope_error, std::abs(reconstructed_slope - 1.7));
}
- REQUIRE(max_slope_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_slope_error) == Catch::Approx(0).margin(1E-14));
}
auto dpk_uh = reconstructions[1]->get<DiscreteFunctionDPk<1, const R3>>();
@@ -321,7 +322,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error =
std::max(max_mean_error, l2Norm(dpk_uh[cell_id](xj[cell_id]) - R3_affine(xj[cell_id])));
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-14));
}
{
@@ -332,7 +333,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_slope_error = std::max(max_slope_error, l2Norm(reconstructed_slope - R3{1.7, -0.6, 3.1}));
}
- REQUIRE(max_slope_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_slope_error) == Catch::Approx(0).margin(1E-14));
}
auto dpk_Ah = reconstructions[2]->get<DiscreteFunctionDPk<1, const R3x3>>();
@@ -343,7 +344,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error =
std::max(max_mean_error, frobeniusNorm(dpk_Ah[cell_id](xj[cell_id]) - R3x3_affine(xj[cell_id])));
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-13));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-13));
}
{
@@ -359,7 +360,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_slope_error = std::max(max_slope_error, //
frobeniusNorm(reconstructed_slope - slops));
}
- REQUIRE(max_slope_error == Catch::Approx(0).margin(1E-13));
+ REQUIRE(parallel::allReduceMax(max_slope_error) == Catch::Approx(0).margin(1E-13));
}
auto dpk_Vh = reconstructions[3]->get<DiscreteFunctionDPkVector<1, const double>>();
@@ -372,7 +373,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error = std::max(max_mean_error, std::abs(dpk_Vh(cell_id, i)(xj[cell_id]) - vector[i]));
}
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-13));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-13));
}
{
@@ -386,8 +387,8 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_slope_error = std::max(max_slope_error, std::abs(reconstructed_slope - slope[i]));
}
- REQUIRE(max_slope_error == Catch::Approx(0).margin(1E-13));
}
+ REQUIRE(parallel::allReduceMax(max_slope_error) == Catch::Approx(0).margin(1E-13));
}
}
}
@@ -414,7 +415,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
parallel_for(
mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id) { fh[cell_id] = R_affine(xj[cell_id]); });
- auto reconstructions = PolynomialReconstruction{}.build(1, fh);
+ auto reconstructions = PolynomialReconstruction{}.build(1, true, true, fh);
auto dpk_fh = reconstructions[0]->get<DiscreteFunctionDPk<2, const double>>();
@@ -424,7 +425,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error =
std::max(max_mean_error, std::abs(dpk_fh[cell_id](xj[cell_id]) - R_affine(xj[cell_id])));
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-14));
}
{
@@ -435,7 +436,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_x_slope_error = std::max(max_x_slope_error, std::abs(reconstructed_slope - 1.7));
}
- REQUIRE(max_x_slope_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_x_slope_error) == Catch::Approx(0).margin(1E-14));
}
{
@@ -446,7 +447,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_y_slope_error = std::max(max_y_slope_error, std::abs(reconstructed_slope - (-1.3)));
}
- REQUIRE(max_y_slope_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_y_slope_error) == Catch::Approx(0).margin(1E-14));
}
}
}
@@ -474,7 +475,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
parallel_for(
mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id) { uh[cell_id] = R3_affine(xj[cell_id]); });
- auto reconstructions = PolynomialReconstruction{}.build(1, uh);
+ auto reconstructions = PolynomialReconstruction{}.build(1, true, true, uh);
auto dpk_uh = reconstructions[0]->get<DiscreteFunctionDPk<2, const R3>>();
@@ -484,7 +485,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error =
std::max(max_mean_error, l2Norm(dpk_uh[cell_id](xj[cell_id]) - R3_affine(xj[cell_id])));
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-14));
}
{
@@ -495,7 +496,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_x_slope_error = std::max(max_x_slope_error, l2Norm(reconstructed_slope - R3{1.7, -0.6, 3.1}));
}
- REQUIRE(max_x_slope_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_x_slope_error) == Catch::Approx(0).margin(1E-14));
}
{
@@ -506,7 +507,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_y_slope_error = std::max(max_y_slope_error, l2Norm(reconstructed_slope - R3{-2.2, 1.3, -1.1}));
}
- REQUIRE(max_y_slope_error == Catch::Approx(0).margin(1E-13));
+ REQUIRE(parallel::allReduceMax(max_y_slope_error) == Catch::Approx(0).margin(1E-14));
}
}
}
@@ -533,7 +534,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
parallel_for(
mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id) { Ah[cell_id] = R2x2_affine(xj[cell_id]); });
- auto reconstructions = PolynomialReconstruction{}.build(1, Ah);
+ auto reconstructions = PolynomialReconstruction{}.build(1, true, true, Ah);
auto dpk_Ah = reconstructions[0]->get<DiscreteFunctionDPk<2, const R2x2>>();
@@ -543,7 +544,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error =
std::max(max_mean_error, frobeniusNorm(dpk_Ah[cell_id](xj[cell_id]) - R2x2_affine(xj[cell_id])));
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-13));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-13));
}
{
@@ -556,7 +557,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
std::max(max_x_slope_error, frobeniusNorm(reconstructed_slope - R2x2{+1.7, +2.1, //
-0.6, -2.3}));
}
- REQUIRE(max_x_slope_error == Catch::Approx(0).margin(1E-12));
+ REQUIRE(parallel::allReduceMax(max_x_slope_error) == Catch::Approx(0).margin(1E-12));
}
{
@@ -569,7 +570,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
std::max(max_y_slope_error, frobeniusNorm(reconstructed_slope - R2x2{+1.2, -2.2, //
-2.1, +1.3}));
}
- REQUIRE(max_y_slope_error == Catch::Approx(0).margin(1E-12));
+ REQUIRE(parallel::allReduceMax(max_y_slope_error) == Catch::Approx(0).margin(1E-12));
}
}
}
@@ -601,7 +602,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
}
});
- auto reconstructions = PolynomialReconstruction{}.build(1, Vh);
+ auto reconstructions = PolynomialReconstruction{}.build(1, true, true, Vh);
auto dpk_Vh = reconstructions[0]->get<DiscreteFunctionDPkVector<2, const double>>();
@@ -613,7 +614,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error = std::max(max_mean_error, std::abs(dpk_Vh(cell_id, i)(xj[cell_id]) - vector[i]));
}
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-13));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-13));
}
{
@@ -627,7 +628,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_x_slope_error = std::max(max_x_slope_error, std::abs(reconstructed_slope - slope[i]));
}
}
- REQUIRE(max_x_slope_error == Catch::Approx(0).margin(1E-12));
+ REQUIRE(parallel::allReduceMax(max_x_slope_error) == Catch::Approx(0).margin(1E-12));
}
{
@@ -641,7 +642,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_y_slope_error = std::max(max_y_slope_error, std::abs(reconstructed_slope - slope[i]));
}
}
- REQUIRE(max_y_slope_error == Catch::Approx(0).margin(1E-12));
+ REQUIRE(parallel::allReduceMax(max_y_slope_error) == Catch::Approx(0).margin(1E-12));
}
}
}
@@ -668,7 +669,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
parallel_for(
mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id) { fh[cell_id] = R_affine(xj[cell_id]); });
- auto reconstructions = PolynomialReconstruction{}.build(1, fh);
+ auto reconstructions = PolynomialReconstruction{}.build(1, true, true, fh);
auto dpk_fh = reconstructions[0]->get<DiscreteFunctionDPk<3, const double>>();
@@ -678,7 +679,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error =
std::max(max_mean_error, std::abs(dpk_fh[cell_id](xj[cell_id]) - R_affine(xj[cell_id])));
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-14));
}
{
@@ -689,7 +690,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_x_slope_error = std::max(max_x_slope_error, std::abs(reconstructed_slope - 1.7));
}
- REQUIRE(max_x_slope_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_x_slope_error) == Catch::Approx(0).margin(1E-14));
}
{
@@ -700,7 +701,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_y_slope_error = std::max(max_y_slope_error, std::abs(reconstructed_slope - (-1.3)));
}
- REQUIRE(max_y_slope_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_y_slope_error) == Catch::Approx(0).margin(1E-14));
}
{
@@ -711,7 +712,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_z_slope_error = std::max(max_z_slope_error, std::abs(reconstructed_slope - 2.1));
}
- REQUIRE(max_z_slope_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_z_slope_error) == Catch::Approx(0).margin(1E-14));
}
}
}
@@ -737,7 +738,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
parallel_for(
mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id) { uh[cell_id] = R3_affine(xj[cell_id]); });
- auto reconstructions = PolynomialReconstruction{}.build(1, uh);
+ auto reconstructions = PolynomialReconstruction{}.build(1, true, true, uh);
auto dpk_uh = reconstructions[0]->get<DiscreteFunctionDPk<3, const R3>>();
@@ -747,7 +748,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error =
std::max(max_mean_error, l2Norm(dpk_uh[cell_id](xj[cell_id]) - R3_affine(xj[cell_id])));
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-14));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-14));
}
{
@@ -758,7 +759,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_x_slope_error = std::max(max_x_slope_error, l2Norm(reconstructed_slope - R3{1.7, -0.6, 3.1}));
}
- REQUIRE(max_x_slope_error == Catch::Approx(0).margin(1E-13));
+ REQUIRE(parallel::allReduceMax(max_x_slope_error) == Catch::Approx(0).margin(1E-13));
}
{
@@ -769,7 +770,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_y_slope_error = std::max(max_y_slope_error, l2Norm(reconstructed_slope - R3{-2.2, 1.3, -1.1}));
}
- REQUIRE(max_y_slope_error == Catch::Approx(0).margin(1E-13));
+ REQUIRE(parallel::allReduceMax(max_y_slope_error) == Catch::Approx(0).margin(1E-13));
}
{
@@ -780,7 +781,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_z_slope_error = std::max(max_z_slope_error, l2Norm(reconstructed_slope - R3{1.8, -3.7, 1.9}));
}
- REQUIRE(max_z_slope_error == Catch::Approx(0).margin(1E-13));
+ REQUIRE(parallel::allReduceMax(max_z_slope_error) == Catch::Approx(0).margin(1E-13));
}
}
}
@@ -808,7 +809,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
parallel_for(
mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id) { Ah[cell_id] = R2x2_affine(xj[cell_id]); });
- auto reconstructions = PolynomialReconstruction{}.build(1, Ah);
+ auto reconstructions = PolynomialReconstruction{}.build(1, true, true, Ah);
auto dpk_Ah = reconstructions[0]->get<DiscreteFunctionDPk<3, const R2x2>>();
@@ -818,7 +819,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error =
std::max(max_mean_error, frobeniusNorm(dpk_Ah[cell_id](xj[cell_id]) - R2x2_affine(xj[cell_id])));
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-13));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-13));
}
{
@@ -830,7 +831,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_x_slope_error = std::max(max_x_slope_error, frobeniusNorm(reconstructed_slope - R2x2{+1.7, 2.1, //
-2.3, +3.1}));
}
- REQUIRE(max_x_slope_error == Catch::Approx(0).margin(1E-12));
+ REQUIRE(parallel::allReduceMax(max_x_slope_error) == Catch::Approx(0).margin(1E-12));
}
{
@@ -843,7 +844,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
std::max(max_y_slope_error, frobeniusNorm(reconstructed_slope - R2x2{+1.2, -2.2, //
1.3, +0.8}));
}
- REQUIRE(max_y_slope_error == Catch::Approx(0).margin(1E-12));
+ REQUIRE(parallel::allReduceMax(max_y_slope_error) == Catch::Approx(0).margin(1E-12));
}
{
@@ -856,7 +857,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
std::max(max_z_slope_error, frobeniusNorm(reconstructed_slope - R2x2{-1.3, -2.4, //
+1.4, -1.8}));
}
- REQUIRE(max_z_slope_error == Catch::Approx(0).margin(1E-12));
+ REQUIRE(parallel::allReduceMax(max_z_slope_error) == Catch::Approx(0).margin(1E-12));
}
}
}
@@ -889,7 +890,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
}
});
- auto reconstructions = PolynomialReconstruction{}.build(1, Vh);
+ auto reconstructions = PolynomialReconstruction{}.build(1, true, true, Vh);
auto dpk_Vh = reconstructions[0]->get<DiscreteFunctionDPkVector<3, const double>>();
@@ -901,7 +902,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_mean_error = std::max(max_mean_error, std::abs(dpk_Vh(cell_id, i)(xj[cell_id]) - vector[i]));
}
}
- REQUIRE(max_mean_error == Catch::Approx(0).margin(1E-13));
+ REQUIRE(parallel::allReduceMax(max_mean_error) == Catch::Approx(0).margin(1E-13));
}
{
@@ -915,7 +916,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_x_slope_error = std::max(max_x_slope_error, std::abs(reconstructed_slope - slope[i]));
}
}
- REQUIRE(max_x_slope_error == Catch::Approx(0).margin(1E-12));
+ REQUIRE(parallel::allReduceMax(max_x_slope_error) == Catch::Approx(0).margin(1E-12));
}
{
@@ -929,7 +930,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_y_slope_error = std::max(max_y_slope_error, std::abs(reconstructed_slope - slope[i]));
}
}
- REQUIRE(max_y_slope_error == Catch::Approx(0).margin(1E-12));
+ REQUIRE(parallel::allReduceMax(max_y_slope_error) == Catch::Approx(0).margin(1E-12));
}
{
@@ -943,7 +944,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
max_z_slope_error = std::max(max_z_slope_error, std::abs(reconstructed_slope - slope[i]));
}
}
- REQUIRE(max_z_slope_error == Catch::Approx(0).margin(1E-12));
+ REQUIRE(parallel::allReduceMax(max_z_slope_error) == Catch::Approx(0).margin(1E-12));
}
}
}
@@ -958,7 +959,7 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
auto p_mesh2 = MeshDataBaseForTests::get().cartesian1DMesh()->get<Mesh<1>>();
DiscreteFunctionP0<double> f2{p_mesh2};
- REQUIRE_THROWS_WITH(PolynomialReconstruction{}.build(1, f1, f2),
+ REQUIRE_THROWS_WITH(PolynomialReconstruction{}.build(1, true, true, f1, f2),
"error: cannot reconstruct functions living of different meshes simultaneously");
}
}
--
GitLab