diff --git a/src/mesh/MeshData.hpp b/src/mesh/MeshData.hpp
index 9977c06a18f6dee5ea6443f4bc38f00190613e9f..c8cef0b7a7af3a27c93ca2e43f27c5a30b7e3843 100644
--- a/src/mesh/MeshData.hpp
+++ b/src/mesh/MeshData.hpp
@@ -11,7 +11,7 @@
template <typename M>
class MeshData
{
-public:
+ public:
typedef M MeshType;
static constexpr size_t dimension = MeshType::dimension;
@@ -21,11 +21,11 @@ public:
static constexpr double inv_dimension = 1./dimension;
-private:
+ private:
const MeshType& m_mesh;
- SubItemValuePerItem<Rd> m_Cjr;
- SubItemValuePerItem<double> m_ljr;
- SubItemValuePerItem<Rd> m_njr;
+ NodeValuePerCell<Rd> m_Cjr;
+ NodeValuePerCell<double> m_ljr;
+ NodeValuePerCell<Rd> m_njr;
Kokkos::View<Rd*> m_xj;
Kokkos::View<double*> m_Vj;
@@ -74,7 +74,7 @@ private:
void _updateCjr() {
if constexpr (dimension == 1) {
// Cjr/njr/ljr are constant overtime
- }
+ }
else if constexpr (dimension == 2) {
const Kokkos::View<const Rd*> xr = m_mesh.xr();
@@ -88,12 +88,12 @@ private:
}
});
- const SubItemValuePerItem<Rd>& Cjr = m_Cjr;
+ const NodeValuePerCell<Rd>& Cjr = m_Cjr;
Kokkos::parallel_for(m_Cjr.numberOfValues(), KOKKOS_LAMBDA(const int& j){
m_ljr[j] = l2Norm(Cjr[j]);
});
- const SubItemValuePerItem<double>& ljr = m_ljr;
+ const NodeValuePerCell<double>& ljr = m_ljr;
Kokkos::parallel_for(m_Cjr.numberOfValues(), KOKKOS_LAMBDA(const int& j){
m_njr[j] = (1./ljr[j])*Cjr[j];
});
@@ -150,12 +150,12 @@ private:
}
});
- const SubItemValuePerItem<Rd>& Cjr = m_Cjr;
+ const NodeValuePerCell<Rd>& Cjr = m_Cjr;
Kokkos::parallel_for(m_Cjr.numberOfValues(), KOKKOS_LAMBDA(const int& j){
m_ljr[j] = l2Norm(Cjr[j]);
});
- const SubItemValuePerItem<double>& ljr = m_ljr;
+ const NodeValuePerCell<double>& ljr = m_ljr;
Kokkos::parallel_for(m_Cjr.numberOfValues(), KOKKOS_LAMBDA(const int& j){
m_njr[j] = (1./ljr[j])*Cjr[j];
});
@@ -163,23 +163,23 @@ private:
static_assert((dimension<=3), "only 1d, 2d and 3d are implemented");
}
-public:
+ public:
const MeshType& mesh() const
{
return m_mesh;
}
- const SubItemValuePerItem<Rd>& Cjr() const
+ const NodeValuePerCell<Rd>& Cjr() const
{
return m_Cjr;
}
- const SubItemValuePerItem<double>& ljr() const
+ const NodeValuePerCell<double>& ljr() const
{
return m_ljr;
}
- const SubItemValuePerItem<Rd>& njr() const
+ const NodeValuePerCell<Rd>& njr() const
{
return m_njr;
}
@@ -202,12 +202,12 @@ public:
}
MeshData(const MeshType& mesh)
- : m_mesh(mesh),
- m_Cjr(mesh.connectivity().m_node_id_per_cell_matrix),
- m_ljr(mesh.connectivity().m_node_id_per_cell_matrix),
- m_njr(mesh.connectivity().m_node_id_per_cell_matrix),
- m_xj("xj", mesh.numberOfCells()),
- m_Vj("Vj", mesh.numberOfCells())
+ : m_mesh(mesh),
+ m_Cjr(mesh.connectivity().m_node_id_per_cell_matrix),
+ m_ljr(mesh.connectivity().m_node_id_per_cell_matrix),
+ m_njr(mesh.connectivity().m_node_id_per_cell_matrix),
+ m_xj("xj", mesh.numberOfCells()),
+ m_Vj("Vj", mesh.numberOfCells())
{
if constexpr (dimension==1) {
// in 1d Cjr are computed once for all
diff --git a/src/scheme/AcousticSolver.hpp b/src/scheme/AcousticSolver.hpp
index b831ec7d0e5438dead4e9a6b7db69af3c4129074..f23dce1c6ee53ac26cb0ebf65dcb7bd37d8ca950 100644
--- a/src/scheme/AcousticSolver.hpp
+++ b/src/scheme/AcousticSolver.hpp
@@ -32,13 +32,13 @@ class AcousticSolver
typedef TinyVector<dimension> Rd;
typedef TinyMatrix<dimension> Rdd;
-private:
+ private:
struct ReduceMin
{
- private:
+ private:
const Kokkos::View<const double*> x_;
- public:
+ public:
typedef Kokkos::View<const double*>::non_const_value_type value_type;
ReduceMin(const Kokkos::View<const double*>& x) : x_ (x) {}
@@ -82,9 +82,9 @@ private:
KOKKOS_INLINE_FUNCTION
void computeAjr(const Kokkos::View<const double*>& rhocj,
- const SubItemValuePerItem<Rd>& Cjr,
- const SubItemValuePerItem<double>& ljr,
- const SubItemValuePerItem<Rd>& njr)
+ const NodeValuePerCell<Rd>& Cjr,
+ const NodeValuePerCell<double>& ljr,
+ const NodeValuePerCell<Rd>& njr)
{
Kokkos::parallel_for("new nested Ajr", m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j) {
const size_t& nb_nodes =m_Ajr.numberOfSubValues(j);
@@ -97,7 +97,7 @@ private:
KOKKOS_INLINE_FUNCTION
const Kokkos::View<const Rdd*>
- computeAr(const SubItemValuePerItem<Rdd>& Ajr) {
+ computeAr(const NodeValuePerCell<Rdd>& Ajr) {
Kokkos::parallel_for(m_mesh.numberOfNodes(), KOKKOS_LAMBDA(const int& r) {
Rdd sum = zero;
const auto& node_to_cell = m_connectivity.m_node_to_cell_matrix.rowConst(r);
@@ -115,8 +115,8 @@ private:
KOKKOS_INLINE_FUNCTION
const Kokkos::View<const Rd*>
- computeBr(const SubItemValuePerItem<Rdd>& Ajr,
- const SubItemValuePerItem<Rd>& Cjr,
+ computeBr(const NodeValuePerCell<Rdd>& Ajr,
+ const NodeValuePerCell<Rd>& Cjr,
const Kokkos::View<const Rd*>& uj,
const Kokkos::View<const double*>& pj) {
Kokkos::parallel_for(m_mesh.numberOfNodes(), KOKKOS_LAMBDA(const int& r) {
@@ -138,40 +138,40 @@ private:
{
for (const auto& handler : m_boundary_condition_list) {
switch (handler.boundaryCondition().type()) {
- case BoundaryCondition::normal_velocity: {
- std::cerr << __FILE__ << ':' << __LINE__ << ": normal_velocity BC NIY\n";
- std::exit(0);
- break;
- }
- case BoundaryCondition::velocity: {
- std::cerr << __FILE__ << ':' << __LINE__ << ": velocity BC NIY\n";
- std::exit(0);
- break;
- }
- case BoundaryCondition::pressure: {
- // const PressureBoundaryCondition& pressure_bc
- // = dynamic_cast<const PressureBoundaryCondition&>(handler.boundaryCondition());
- std::cerr << __FILE__ << ':' << __LINE__ << ": pressure BC NIY\n";
- std::exit(0);
- break;
- }
- case BoundaryCondition::symmetry: {
- const SymmetryBoundaryCondition<dimension>& symmetry_bc
- = dynamic_cast<const SymmetryBoundaryCondition<dimension>&>(handler.boundaryCondition());
- const Rd& n = symmetry_bc.outgoingNormal();
-
- const Rdd I = identity;
- const Rdd nxn = tensorProduct(n,n);
- const Rdd P = I-nxn;
-
- Kokkos::parallel_for(symmetry_bc.numberOfNodes(), KOKKOS_LAMBDA(const int& r_number) {
- const int r = symmetry_bc.nodeList()[r_number];
-
- m_Ar(r) = P*m_Ar(r)*P + nxn;
- m_br(r) = P*m_br(r);
- });
- break;
- }
+ case BoundaryCondition::normal_velocity: {
+ std::cerr << __FILE__ << ':' << __LINE__ << ": normal_velocity BC NIY\n";
+ std::exit(0);
+ break;
+ }
+ case BoundaryCondition::velocity: {
+ std::cerr << __FILE__ << ':' << __LINE__ << ": velocity BC NIY\n";
+ std::exit(0);
+ break;
+ }
+ case BoundaryCondition::pressure: {
+ // const PressureBoundaryCondition& pressure_bc
+ // = dynamic_cast<const PressureBoundaryCondition&>(handler.boundaryCondition());
+ std::cerr << __FILE__ << ':' << __LINE__ << ": pressure BC NIY\n";
+ std::exit(0);
+ break;
+ }
+ case BoundaryCondition::symmetry: {
+ const SymmetryBoundaryCondition<dimension>& symmetry_bc
+ = dynamic_cast<const SymmetryBoundaryCondition<dimension>&>(handler.boundaryCondition());
+ const Rd& n = symmetry_bc.outgoingNormal();
+
+ const Rdd I = identity;
+ const Rdd nxn = tensorProduct(n,n);
+ const Rdd P = I-nxn;
+
+ Kokkos::parallel_for(symmetry_bc.numberOfNodes(), KOKKOS_LAMBDA(const int& r_number) {
+ const int r = symmetry_bc.nodeList()[r_number];
+
+ m_Ar(r) = P*m_Ar(r)*P + nxn;
+ m_br(r) = P*m_br(r);
+ });
+ break;
+ }
}
}
}
@@ -190,9 +190,9 @@ private:
}
void
- computeFjr(const SubItemValuePerItem<Rdd>& Ajr,
+ computeFjr(const NodeValuePerCell<Rdd>& Ajr,
const Kokkos::View<const Rd*>& ur,
- const SubItemValuePerItem<Rd>& Cjr,
+ const NodeValuePerCell<Rd>& Cjr,
const Kokkos::View<const Rd*>& uj,
const Kokkos::View<const double*>& pj)
{
@@ -221,9 +221,9 @@ private:
const Kokkos::View<const double*>& pj,
const Kokkos::View<const double*>& cj,
const Kokkos::View<const double*>& Vj,
- const SubItemValuePerItem<Rd>& Cjr,
- const SubItemValuePerItem<double>& ljr,
- const SubItemValuePerItem<Rd>& njr) {
+ const NodeValuePerCell<Rd>& Cjr,
+ const NodeValuePerCell<double>& ljr,
+ const NodeValuePerCell<Rd>& njr) {
const Kokkos::View<const double*> rhocj = computeRhoCj(rhoj, cj);
computeAjr(rhocj, Cjr, ljr, njr);
@@ -238,30 +238,30 @@ private:
}
Kokkos::View<Rd*> m_br;
- SubItemValuePerItem<Rdd> m_Ajr;
+ NodeValuePerCell<Rdd> m_Ajr;
Kokkos::View<Rdd*> m_Ar;
Kokkos::View<Rdd*> m_inv_Ar;
- SubItemValuePerItem<Rd> m_Fjr;
+ NodeValuePerCell<Rd> m_Fjr;
Kokkos::View<Rd*> m_ur;
Kokkos::View<double*> m_rhocj;
Kokkos::View<double*> m_Vj_over_cj;
-public:
+ public:
AcousticSolver(MeshData& mesh_data,
UnknownsType& unknowns,
const std::vector<BoundaryConditionHandler>& bc_list)
- : m_mesh_data(mesh_data),
- m_mesh(mesh_data.mesh()),
- m_connectivity(m_mesh.connectivity()),
- m_boundary_condition_list(bc_list),
- m_br("br", m_mesh.numberOfNodes()),
- m_Ajr(m_connectivity.m_node_id_per_cell_matrix),
- m_Ar("Ar", m_mesh.numberOfNodes()),
- m_inv_Ar("inv_Ar", m_mesh.numberOfNodes()),
- m_Fjr(m_connectivity.m_node_id_per_cell_matrix),
- m_ur("ur", m_mesh.numberOfNodes()),
- m_rhocj("rho_c", m_mesh.numberOfCells()),
- m_Vj_over_cj("Vj_over_cj", m_mesh.numberOfCells())
+ : m_mesh_data(mesh_data),
+ m_mesh(mesh_data.mesh()),
+ m_connectivity(m_mesh.connectivity()),
+ m_boundary_condition_list(bc_list),
+ m_br("br", m_mesh.numberOfNodes()),
+ m_Ajr(m_connectivity.m_node_id_per_cell_matrix),
+ m_Ar("Ar", m_mesh.numberOfNodes()),
+ m_inv_Ar("inv_Ar", m_mesh.numberOfNodes()),
+ m_Fjr(m_connectivity.m_node_id_per_cell_matrix),
+ m_ur("ur", m_mesh.numberOfNodes()),
+ m_rhocj("rho_c", m_mesh.numberOfCells()),
+ m_Vj_over_cj("Vj_over_cj", m_mesh.numberOfCells())
{
;
}
@@ -270,7 +270,7 @@ public:
double acoustic_dt(const Kokkos::View<const double*>& Vj,
const Kokkos::View<const double*>& cj) const
{
- const SubItemValuePerItem<double>& ljr = m_mesh_data.ljr();
+ const NodeValuePerCell<double>& ljr = m_mesh_data.ljr();
Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j){
const auto& cell_nodes = m_mesh.connectivity().m_cell_to_node_matrix.rowConst(j);
@@ -303,14 +303,14 @@ public:
const Kokkos::View<const Rd*> xj = m_mesh_data.xj();
const Kokkos::View<const double*> Vj = m_mesh_data.Vj();
- const SubItemValuePerItem<Rd>& Cjr = m_mesh_data.Cjr();
- const SubItemValuePerItem<double>& ljr = m_mesh_data.ljr();
- const SubItemValuePerItem<Rd>& njr = m_mesh_data.njr();
+ const NodeValuePerCell<Rd>& Cjr = m_mesh_data.Cjr();
+ const NodeValuePerCell<double>& ljr = m_mesh_data.ljr();
+ const NodeValuePerCell<Rd>& njr = m_mesh_data.njr();
Kokkos::View<Rd*> xr = m_mesh.xr();
computeExplicitFluxes(xr, xj, rhoj, uj, pj, cj, Vj, Cjr, ljr, njr);
- const SubItemValuePerItem<Rd>& Fjr = m_Fjr;
+ const NodeValuePerCell<Rd>& Fjr = m_Fjr;
const Kokkos::View<const Rd*> ur = m_ur;
const Kokkos::View<const double*> inv_mj = unknowns.invMj();
diff --git a/src/scheme/SubItemValuePerItem.hpp b/src/scheme/SubItemValuePerItem.hpp
index 102f046dbe41fc7fc26a6e861207a5c5edd36e14..f185b519ae4fdd810b6ceb5c2fe9d54b59500286 100644
--- a/src/scheme/SubItemValuePerItem.hpp
+++ b/src/scheme/SubItemValuePerItem.hpp
@@ -7,12 +7,21 @@
typedef Kokkos::StaticCrsGraph<unsigned int, Kokkos::HostSpace> ConnectivityMatrix;
-template <typename DataType>
+enum class TypeOfItem {
+ node,
+ edge,
+ face,
+ cell
+};
+
+template <typename DataType,
+ TypeOfItem SubItemType,
+ TypeOfItem ItemType>
class SubItemValuePerItem
{
#warning should eventually filter const from DataType
private:
- ConnectivityMatrix m_node_id_per_cell_matrix;
+ ConnectivityMatrix m_subitem_id_per_item_matrix;
Kokkos::View<DataType*> m_values;
@@ -101,13 +110,13 @@ class SubItemValuePerItem
KOKKOS_FORCEINLINE_FUNCTION
DataType& operator()(const size_t& j, const size_t& r)
{
- return m_values[m_node_id_per_cell_matrix.row_map[j]+r];
+ return m_values[m_subitem_id_per_item_matrix.row_map[j]+r];
}
KOKKOS_FORCEINLINE_FUNCTION
const DataType& operator()(const size_t& j, const size_t& r) const
{
- return m_values[m_node_id_per_cell_matrix.row_map[j]+r];
+ return m_values[m_subitem_id_per_item_matrix.row_map[j]+r];
}
KOKKOS_INLINE_FUNCTION
@@ -129,37 +138,37 @@ class SubItemValuePerItem
}
KOKKOS_INLINE_FUNCTION
- size_t numberOfEntities() const
+ size_t numberOfItems() const
{
- return m_node_id_per_cell_matrix.numRows();
+ return m_subitem_id_per_item_matrix.numRows();
}
KOKKOS_INLINE_FUNCTION
size_t numberOfSubValues(const size_t& i_cell) const
{
- return m_node_id_per_cell_matrix.row_map[i_cell+1]-m_node_id_per_cell_matrix.row_map[i_cell];
+ return m_subitem_id_per_item_matrix.row_map[i_cell+1]-m_subitem_id_per_item_matrix.row_map[i_cell];
}
KOKKOS_INLINE_FUNCTION
- SubView cellValues(const size_t& i_cell)
+ SubView itemValues(const size_t& i_cell)
{
- const ConnectivityMatrix::size_type& cell_begin = m_node_id_per_cell_matrix.row_map[i_cell];
- const ConnectivityMatrix::size_type& cell_end = m_node_id_per_cell_matrix.row_map[i_cell+1];
+ const ConnectivityMatrix::size_type& cell_begin = m_subitem_id_per_item_matrix.row_map[i_cell];
+ const ConnectivityMatrix::size_type& cell_end = m_subitem_id_per_item_matrix.row_map[i_cell+1];
return SubView(m_values, cell_begin, cell_end);
}
KOKKOS_INLINE_FUNCTION
- SubViewConst cellValues(const size_t& i_cell) const
+ SubViewConst itemValues(const size_t& i_cell) const
{
- const ConnectivityMatrix::size_type& cell_begin = m_node_id_per_cell_matrix.row_map[i_cell];
- const ConnectivityMatrix::size_type& cell_end = m_node_id_per_cell_matrix.row_map[i_cell+1];
+ const ConnectivityMatrix::size_type& cell_begin = m_subitem_id_per_item_matrix.row_map[i_cell];
+ const ConnectivityMatrix::size_type& cell_end = m_subitem_id_per_item_matrix.row_map[i_cell+1];
return SubViewConst(m_values, cell_begin, cell_end);
}
SubItemValuePerItem() = default;
- SubItemValuePerItem(const ConnectivityMatrix& node_id_per_cell_matrix)
- : m_node_id_per_cell_matrix(node_id_per_cell_matrix),
- m_values("values", m_node_id_per_cell_matrix.entries.extent(0))
+ SubItemValuePerItem(const ConnectivityMatrix& subitem_id_per_item_matrix)
+ : m_subitem_id_per_item_matrix(subitem_id_per_item_matrix),
+ m_values("values", m_subitem_id_per_item_matrix.entries.extent(0))
{
;
}
@@ -167,4 +176,8 @@ class SubItemValuePerItem
~SubItemValuePerItem() = default;
};
+
+template <typename DataType>
+using NodeValuePerCell = SubItemValuePerItem<DataType, TypeOfItem::node, TypeOfItem::cell>;
+
#endif // SUBITEM_VALUE_PER_ITEM_HPP