Select Git revision
BiCGStab.hpp
-
Stéphane Del Pino authoredStéphane Del Pino authored
Connectivity.hpp 9.46 KiB
#ifndef CONNECTIVITY_HPP
#define CONNECTIVITY_HPP
#include <PastisAssert.hpp>
#include <TinyVector.hpp>
#include <Kokkos_Core.hpp>
#include <IConnectivity.hpp>
#include <ConnectivityMatrix.hpp>
#include <SubItemValuePerItem.hpp>
#include <ConnectivityComputer.hpp>
#include <vector>
#include <unordered_map>
#include <algorithm>
#include <RefId.hpp>
#include <ItemType.hpp>
#include <RefNodeList.hpp>
#include <RefFaceList.hpp>
#include <tuple>
#include <algorithm>
template <size_t Dimension>
class Connectivity;
template <size_t Dimension>
class ConnectivityFace;
template<>
class ConnectivityFace<1>
{
public:
friend struct Hash;
struct Hash
{
size_t operator()(const ConnectivityFace& f) const;
};
};
template<>
class ConnectivityFace<2>
{
public:
friend struct Hash;
struct Hash
{
size_t operator()(const ConnectivityFace& f) const {
size_t hash = 0;
hash ^= std::hash<unsigned int>()(f.m_node0_id);
hash ^= std::hash<unsigned int>()(f.m_node1_id) >> 1;
return hash;
}
};
unsigned int m_node0_id;
unsigned int m_node1_id;
friend std::ostream& operator<<(std::ostream& os, const ConnectivityFace& f)
{
os << f.m_node0_id << ' ' << f.m_node1_id << ' ';
return os;
}
KOKKOS_INLINE_FUNCTION
bool operator==(const ConnectivityFace& f) const
{
return ((m_node0_id == f.m_node0_id) and
(m_node1_id == f.m_node1_id));
}
KOKKOS_INLINE_FUNCTION
bool operator<(const ConnectivityFace& f) const
{
return ((m_node0_id<f.m_node0_id) or
((m_node0_id == f.m_node0_id) and
(m_node1_id<f.m_node1_id)));
}
KOKKOS_INLINE_FUNCTION
ConnectivityFace& operator=(const ConnectivityFace&) = default;
KOKKOS_INLINE_FUNCTION
ConnectivityFace& operator=(ConnectivityFace&&) = default;
KOKKOS_INLINE_FUNCTION
ConnectivityFace(const std::vector<unsigned int>& given_node_id_list)
{
Assert(given_node_id_list.size()==2);
#warning rework this dirty constructor
const auto& [min, max] = std::minmax(given_node_id_list[0], given_node_id_list[1]);
m_node0_id = min;
m_node1_id = max;
}
KOKKOS_INLINE_FUNCTION
ConnectivityFace(const ConnectivityFace&) = default;
KOKKOS_INLINE_FUNCTION
ConnectivityFace(ConnectivityFace&&) = default;
KOKKOS_INLINE_FUNCTION
~ConnectivityFace() = default;
};
template <>
class ConnectivityFace<3>
{
private:
friend class Connectivity<3>;
friend struct Hash;
struct Hash
{
size_t operator()(const ConnectivityFace& f) const {
size_t hash = 0;
for (size_t i=0; i<f.m_node_id_list.size(); ++i) {
hash ^= std::hash<unsigned int>()(f.m_node_id_list[i]) >> i;
}
return hash;
}
};
bool m_reversed;
std::vector<unsigned int> m_node_id_list;
friend std::ostream& operator<<(std::ostream& os, const ConnectivityFace& f)
{
for (auto id : f.m_node_id_list) {
std::cout << id << ' ';
}
return os;
}
KOKKOS_INLINE_FUNCTION
const bool& reversed() const
{
return m_reversed;
}
KOKKOS_INLINE_FUNCTION
const std::vector<unsigned int>& nodeIdList() const
{
return m_node_id_list;
}
KOKKOS_INLINE_FUNCTION
std::vector<unsigned int> _sort(const std::vector<unsigned int>& node_list)
{
const auto min_id = std::min_element(node_list.begin(), node_list.end());
const int shift = std::distance(node_list.begin(), min_id);
std::vector<unsigned int> rotated_node_list(node_list.size());
if (node_list[(shift+1)%node_list.size()] > node_list[(shift+node_list.size()-1)%node_list.size()]) {
for (size_t i=0; i<node_list.size(); ++i) {
rotated_node_list[i] = node_list[(shift+node_list.size()-i)%node_list.size()];
m_reversed = true;
}
} else {
for (size_t i=0; i<node_list.size(); ++i) {
rotated_node_list[i] = node_list[(shift+i)%node_list.size()];
}
}
return rotated_node_list;
}
KOKKOS_INLINE_FUNCTION
ConnectivityFace(const std::vector<unsigned int>& given_node_id_list)
: m_reversed(false),
m_node_id_list(_sort(given_node_id_list))
{
;
}
public:
bool operator==(const ConnectivityFace& f) const
{
if (m_node_id_list.size() == f.nodeIdList().size()) {
for (size_t j=0; j<m_node_id_list.size(); ++j) {
if (m_node_id_list[j] != f.nodeIdList()[j]) {
return false;
}
}
return true;
}
return false;
}
KOKKOS_INLINE_FUNCTION
bool operator<(const ConnectivityFace& f) const
{
const size_t min_nb_nodes = std::min(f.m_node_id_list.size(), m_node_id_list.size());
for (size_t i=0; i<min_nb_nodes; ++i) {
if (m_node_id_list[i] < f.m_node_id_list[i]) return true;
if (m_node_id_list[i] != f.m_node_id_list[i]) return false;
}
return m_node_id_list.size() < f.m_node_id_list.size();
}
KOKKOS_INLINE_FUNCTION
ConnectivityFace& operator=(const ConnectivityFace&) = default;
KOKKOS_INLINE_FUNCTION
ConnectivityFace& operator=(ConnectivityFace&&) = default;
KOKKOS_INLINE_FUNCTION
ConnectivityFace(const ConnectivityFace&) = default;
KOKKOS_INLINE_FUNCTION
ConnectivityFace(ConnectivityFace&&) = default;
KOKKOS_INLINE_FUNCTION
ConnectivityFace() = delete;
KOKKOS_INLINE_FUNCTION
~ConnectivityFace() = default;
};
template <size_t Dimension>
class Connectivity final
: public IConnectivity
{
private:
constexpr static auto& itemId = ItemId<Dimension>::itemId;
public:
static constexpr size_t dimension = Dimension;
template <ItemType item_type_0, ItemType item_type_1>
const ConnectivityMatrix& itemToItemMatrix() const
{
return m_item_to_item_matrix[itemId(item_type_0)][itemId(item_type_1)];
};
const ConnectivityMatrix& getMatrix(const ItemType& item_type_0,
const ItemType& item_type_1) const
{
return m_item_to_item_matrix[itemId(item_type_0)][itemId(item_type_1)];
}
const auto& cellFaceIsReversed() const
{
static_assert(dimension>1, "reversed faces makes no sense in dimension 1");
return m_cell_face_is_reversed;
}
const auto& cellLocalNumbersInTheirNodes() const
{
return m_cell_local_numbers_in_their_nodes;
}
const auto& nodeLocalNumbersInTheirCells() const
{
return m_node_local_numbers_in_their_cells;
}
const auto& faceLocalNumbersInTheirCells() const
{
if constexpr(dimension>1) {
return m_face_local_numbers_in_their_cells;
} else {
return m_node_local_numbers_in_their_cells;
}
}
const auto& nodeToFaceLocalNode() const
{
static_assert(Dimension==1,"this function has no meaning in 1d");
return m_node_local_numbers_in_their_faces;
}
private:
ConnectivityMatrix m_item_to_item_matrix[Dimension+1][Dimension+1];
FaceValuePerCell<const bool> m_cell_face_is_reversed;
NodeValuePerCell<const unsigned short> m_cell_local_numbers_in_their_nodes;
CellValuePerFace<const unsigned short> m_face_local_numbers_in_their_cells;
CellValuePerNode<const unsigned short> m_node_local_numbers_in_their_cells;
// not plugged ...
#warning remaining def
NodeValuePerFace<unsigned short> m_node_local_numbers_in_their_faces;
FaceValuePerNode<unsigned short> m_face_local_numbers_in_their_nodes;
// ... not plugged
ConnectivityComputer m_connectivity_computer;
std::vector<RefFaceList> m_ref_face_list;
std::vector<RefNodeList> m_ref_node_list;
Kokkos::View<double*> m_inv_cell_nb_nodes;
using Face = ConnectivityFace<Dimension>;
std::unordered_map<Face, unsigned int, typename Face::Hash> m_face_number_map;
void _computeFaceCellConnectivities();
public:
void addRefFaceList(const RefFaceList& ref_face_list)
{
m_ref_face_list.push_back(ref_face_list);
}
size_t numberOfRefFaceList() const
{
return m_ref_face_list.size();
}
const RefFaceList& refFaceList(const size_t& i) const
{
return m_ref_face_list[i];
}
void addRefNodeList(const RefNodeList& ref_node_list)
{
m_ref_node_list.push_back(ref_node_list);
}
size_t numberOfRefNodeList() const
{
return m_ref_node_list.size();
}
const RefNodeList& refNodeList(const size_t& i) const
{
return m_ref_node_list[i];
}
KOKKOS_INLINE_FUNCTION
size_t numberOfNodes() const
{
const auto& node_to_cell_matrix
= m_item_to_item_matrix[itemId(ItemType::node)][itemId(ItemType::cell)];
return node_to_cell_matrix.numRows();
}
KOKKOS_INLINE_FUNCTION
size_t numberOfFaces() const
{
const auto& face_to_node_matrix
= m_item_to_item_matrix[itemId(ItemType::face)][itemId(ItemType::node)];
return face_to_node_matrix.numRows();
}
KOKKOS_INLINE_FUNCTION
size_t numberOfCells() const
{
const auto& cell_to_node_matrix
= m_item_to_item_matrix[itemId(ItemType::cell)][itemId(ItemType::node)];
return cell_to_node_matrix.numRows();
}
const Kokkos::View<const double*> invCellNbNodes() const
{
return m_inv_cell_nb_nodes;
}
unsigned int getFaceNumber(const std::vector<unsigned int>& face_nodes) const
{
const Face face(face_nodes);
auto i_face = m_face_number_map.find(face);
if (i_face == m_face_number_map.end()) {
std::cerr << "Face " << face << " not found!\n";
throw std::exception();
std::exit(0);
}
return i_face->second;
}
Connectivity(const Connectivity&) = delete;
Connectivity(const std::vector<std::vector<unsigned int>>& cell_by_node_vector);
~Connectivity()
{
;
}
};
using Connectivity3D = Connectivity<3>;
using Connectivity2D = Connectivity<2>;
using Connectivity1D = Connectivity<1>;
#endif // CONNECTIVITY_HPP