Skip to content
Snippets Groups Projects
Commit b4c6ee70 authored by Stéphane Del Pino's avatar Stéphane Del Pino
Browse files

Begin code cleaning and interface improvements

parent a03881be
No related branches found
No related tags found
1 merge request!205High-order polynomial reconstruction
Show whitespace changes
Inline Side-by-side
src/scheme/DiscreteFunctionDPk.hpp
+ 8
− 8
View file @ b4c6ee70
......@@ -135,7 +135,7 @@ class PolynomialCenteredCanonicalBasisView
if constexpr (Dimension == 1) {
const double x_xj = (x - m_xj)[0];
DataType result = m_coefficients[m_degree];
std::remove_const_t<DataType> result = m_coefficients[m_degree];
for (ssize_t i_coeffiencient = m_degree - 1; i_coeffiencient >= 0; --i_coeffiencient) {
result = x_xj * result + m_coefficients[i_coeffiencient];
}
......@@ -148,9 +148,9 @@ class PolynomialCenteredCanonicalBasisView
size_t i = m_coefficients.size() - 1;
DataType result = m_coefficients[i--];
std::remove_const_t<DataType> result = m_coefficients[i--];
for (ssize_t i_y = m_degree - 1; i_y >= 0; --i_y) {
DataType x_result = m_coefficients[i--];
std::remove_const_t<DataType> x_result = m_coefficients[i--];
for (ssize_t i_x = m_degree - i_y - 1; i_x >= 0; --i_x) {
x_result = x_xj * x_result + m_coefficients[i--];
}
......@@ -166,11 +166,11 @@ class PolynomialCenteredCanonicalBasisView
size_t i = m_coefficients.size() - 1;
DataType result = m_coefficients[i--];
std::remove_const_t<DataType> result = m_coefficients[i--];
for (ssize_t i_z = m_degree - 1; i_z >= 0; --i_z) {
DataType y_result = m_coefficients[i--];
std::remove_const_t<DataType> y_result = m_coefficients[i--];
for (ssize_t i_y = m_degree - i_z - 1; i_y >= 0; --i_y) {
DataType x_result = m_coefficients[i--];
std::remove_const_t<DataType> x_result = m_coefficients[i--];
for (ssize_t i_x = m_degree - i_z - i_y - 1; i_x >= 0; --i_x) {
x_result = x_xj * x_result + m_coefficients[i--];
}
......@@ -297,7 +297,7 @@ class DiscreteFunctionDPk
m_xj{MeshDataManager::instance().getMeshData(*m_mesh_v->get<Mesh<Dimension>>()).xj()}
{}
DiscreteFunctionDPk(const std::shared_ptr<const MeshVariant>& mesh_v, const CellValue<DataType>& cell_array)
DiscreteFunctionDPk(const std::shared_ptr<const MeshVariant>& mesh_v, const CellArray<DataType>& cell_array)
: m_mesh_v{mesh_v},
m_degree{BasisView::degreeFromDimension(cell_array.sizeOfArrays())},
m_by_cell_coefficients{cell_array},
......@@ -312,7 +312,7 @@ class DiscreteFunctionDPk
{}
template <MeshConcept MeshType>
DiscreteFunctionDPk(const std::shared_ptr<const MeshType>& p_mesh, const CellValue<DataType>& cell_array)
DiscreteFunctionDPk(const std::shared_ptr<const MeshType>& p_mesh, const CellArray<DataType>& cell_array)
: DiscreteFunctionDPk{p_mesh->meshVariant(), cell_array}
{
Assert(m_mesh_v->connectivity().id() == cell_array.connectivity_ptr()->id());
......
src/scheme/PolynomialReconstruction.cpp
+ 9
− 133
View file @ b4c6ee70
......@@ -305,141 +305,17 @@ PolynomialReconstruction::_build(
}
});
return {};
// parallel_for(
// mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_j_id) {
// if (cell_is_owned[cell_j_id]) {
// const int32_t t = tokens.acquire();
// auto stencil_cell_list = stencil[cell_j_id];
// ShrinkMatrixView B(B_pool[t], stencil_cell_list.size());
// size_t column_begin = 0;
// for (size_t i_discrete_function_variant = 0;
// i_discrete_function_variant < discrete_function_variant_list.size(); ++i_discrete_function_variant) {
// // const auto& discrete_function_variant = discrete_function_variant_list[i_discrete_function_variant];
// // std::visit(
// // [&](auto&& discrete_function) {
// // using DiscreteFunctionT = std::decay_t<decltype(discrete_function)>;
// // if constexpr (is_discrete_function_P0_v<DiscreteFunctionT>) {
// // using DataType = std::decay_t<typename DiscreteFunctionT::data_type>;
// // const DataType& pj = discrete_function[cell_j_id];
// // for (size_t i = 0; i < stencil_cell_list.size(); ++i) {
// // const CellId cell_i_id = stencil_cell_list[i];
// // const DataType& pi_pj = discrete_function[cell_i_id] - pj;
// // if constexpr (std::is_arithmetic_v<DataType>) {
// // B(i, column_begin) = pi_pj;
// // } else if constexpr (is_tiny_vector_v<DataType>) {
// // for (size_t kB = column_begin, k = 0; k < DataType::Dimension; ++k, ++kB) {
// // B(i, kB) = pi_pj[k];
// // }
// // } else if constexpr (is_tiny_matrix_v<DataType>) {
// // for (size_t k = 0; k < DataType::NumberOfRows; ++k) {
// // for (size_t l = 0; l < DataType::NumberOfColumns; ++l) {
// // B(i, column_begin + k * DataType::NumberOfColumns + l) = pi_pj(k, l);
// // }
// // }
// // }
// // }
// // if constexpr (std::is_arithmetic_v<DataType>) {
// // ++column_begin;
// // } else if constexpr (is_tiny_vector_v<DataType> or is_tiny_matrix_v<DataType>) {
// // column_begin += DataType::Dimension;
// // }
// // }
// // },
// // discrete_function_variant->discreteFunction());
// using DataType = double;
// const DataType& pj = my_discrete_function[cell_j_id];
// for (size_t i = 0; i < stencil_cell_list.size(); ++i) {
// const CellId cell_i_id = stencil_cell_list[i];
// const DataType& pi_pj = my_discrete_function[cell_i_id] - pj;
// B(i, column_begin) = pi_pj;
// }
// ++column_begin;
// }
// ShrinkMatrixView A(A_pool[t], stencil_cell_list.size());
// const Rd& Xj = xj[cell_j_id];
// for (size_t i = 0; i < stencil_cell_list.size(); ++i) {
// const CellId cell_i_id = stencil_cell_list[i];
// const Rd Xi_Xj = xj[cell_i_id] - Xj;
// for (size_t l = 0; l < MeshType::Dimension; ++l) {
// A(i, l) = Xi_Xj[l];
// }
// }
// const SmallMatrix<double>& X = X_pool[t];
// Givens::solveCollectionInPlace(A, X, B);
// column_begin = 0;
// for (size_t i_discrete_function_variant = 0;
// i_discrete_function_variant < discrete_function_variant_list.size(); ++i_discrete_function_variant) {
// auto discrete_function_variant = discrete_function_variant_list[i_discrete_function_variant];
// std::visit(
// [&](auto&& discrete_function) {
// using DiscreteFunctionT = std::decay_t<decltype(discrete_function)>;
// if constexpr (is_discrete_function_P0_v<DiscreteFunctionT>) {
// using DataType = std::decay_t<typename DiscreteFunctionT::data_type>;
// const DataType& pj = discrete_function[cell_j_id];
// auto discrete_function_dpk = mutable_discrete_function_dpk_variant_list[i_discrete_function_variant]
// .get<DiscreteFunctionDPk<MeshType::Dimension, DataType>>();
// auto dpk_j = discrete_function_dpk.coefficients(cell_j_id);
// dpk_j[0] = pj;
// if constexpr (std::is_arithmetic_v<DataType>) {
// for (size_t i = 0; i < MeshType::Dimension; ++i) {
// auto& dpk_j_ip1 = dpk_j[i + 1];
// dpk_j_ip1 = X(i, column_begin);
// }
// ++column_begin;
// } else if constexpr (is_tiny_vector_v<DataType>) {
// for (size_t i = 0; i < MeshType::Dimension; ++i) {
// auto& dpk_j_ip1 = dpk_j[i + 1];
// for (size_t k = 0; k < DataType::Dimension; ++k) {
// dpk_j_ip1[k] = X(i, column_begin + k);
// }
// }
// column_begin += DataType::Dimension;
// } else if constexpr (is_tiny_matrix_v<DataType>) {
// for (size_t i = 0; i < MeshType::Dimension; ++i) {
// auto& dpk_j_ip1 = dpk_j[i + 1];
// for (size_t k = 0; k < DataType::NumberOfRows; ++k) {
// for (size_t l = 0; l < DataType::NumberOfColumns; ++l) {
// dpk_j_ip1(k, l) = X(i, column_begin + k * DataType::NumberOfColumns + l);
// }
// }
// }
// column_begin += DataType::Dimension;
// }
// } else {
// throw NotImplementedError("discrete function type");
// }
// },
// discrete_function_variant->discreteFunction());
// }
// tokens.release(t);
// }
// });
std::vector<std::shared_ptr<const DiscreteFunctionDPkVariant>> discrete_function_dpk_variant_list;
// for (auto discrete_function_dpk_variant_p : mutable_discrete_function_dpk_variant_list) {
// discrete_function_dpk_variant_list.push_back(
// std::make_shared<const DiscreteFunctionDPkVariant>(discrete_function_dpk_variant_p));
// }
for (auto discrete_function_dpk_variant_p : mutable_discrete_function_dpk_variant_list) {
std::visit(
[&](auto&& mutable_function_dpk) {
synchronize(mutable_function_dpk.cellArrays());
discrete_function_dpk_variant_list.push_back(
std::make_shared<DiscreteFunctionDPkVariant>(mutable_function_dpk));
},
discrete_function_dpk_variant_p.mutableDiscreteFunctionDPk());
}
return discrete_function_dpk_variant_list;
}
......
tests/test_PolynomialReconstruction.cpp
+ 14
− 2
View file @ b4c6ee70
......@@ -12,10 +12,21 @@
#include <mesh/Mesh.hpp>
#include <mesh/MeshDataManager.hpp>
#include <scheme/DiscreteFunctionP0.hpp>
#include <scheme/DiscreteFunctionVariant.hpp>
#include <scheme/PolynomialReconstruction.hpp>
#include <MeshDataBaseForTests.hpp>
template <typename... DataType>
std::vector<std::shared_ptr<const DiscreteFunctionVariant>>
build_list(DiscreteFunctionP0<DataType>... input)
{
std::vector<std::shared_ptr<const DiscreteFunctionVariant>> variant_vector;
variant_vector.push_back(std::make_shared<DiscreteFunctionVariant>(input...));
return variant_vector;
}
// clazy:excludeall=non-pod-global-static
TEST_CASE("PolynomialReconstruction", "[scheme]")
......@@ -40,8 +51,9 @@ TEST_CASE("PolynomialReconstruction", "[scheme]")
parallel_for(
mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id) { fh[cell_id] = affine(xj[cell_id]); });
PolynomialReconstruction reconstruction;
auto dpk = reconstruction.build(mesh, fh);
auto reconstructions = PolynomialReconstruction{}.build(build_list(fh));
auto dpk = reconstructions[0]->get<DiscreteFunctionDPk<1, const double>>();
{
double max_mean_error = 0;
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please to comment