Select Git revision
Connectivity.cpp
Connectivity.hpp 17.54 KiB
#ifndef CONNECTIVITY_HPP
#define CONNECTIVITY_HPP
#include <PastisMacros.hpp>
#include <PastisAssert.hpp>
#include <PastisOStream.hpp>
#include <PastisUtils.hpp>
#include <PastisTraits.hpp>
#include <TinyVector.hpp>
#include <ItemValue.hpp>
#include <IConnectivity.hpp>
#include <ConnectivityMatrix.hpp>
#include <ItemToItemMatrix.hpp>
#include <SubItemValuePerItem.hpp>
#include <ConnectivityComputer.hpp>
#include <vector>
#include <algorithm>
#include <CellType.hpp>
#include <CSRGraph.hpp>
#include <RefId.hpp>
#include <ItemType.hpp>
#include <RefNodeList.hpp>
#include <RefFaceList.hpp>
#include <SynchronizerManager.hpp>
#include <tuple>
#include <algorithm>
#include <set>
class ConnectivityDescriptor
{
private:
std::vector<RefFaceList> m_ref_face_list;
public:
std::vector<std::vector<unsigned int>> cell_to_node_vector;
#warning should be set in 2d
std::vector<std::vector<unsigned int>> cell_to_face_vector;
// std::vector<std::vector<unsigned int>> face_to_cell_vector;
std::vector<std::vector<unsigned int>> face_to_node_vector;
template <typename ItemOfItemT>
auto& itemOfItemVector()
{
if constexpr (std::is_same_v<ItemOfItemT,NodeOfCell>) {
return cell_to_node_vector;
} else if constexpr (std::is_same_v<ItemOfItemT,FaceOfCell>) {
return cell_to_face_vector;
} else if constexpr (std::is_same_v<ItemOfItemT,NodeOfCell>) {
return face_to_node_vector;
} else {
static_assert(is_false_v<ItemOfItemT>, "Unexpected item of item type");
}
}
std::vector<Array<bool>> cell_face_is_reversed_vector;
std::vector<CellType> cell_type_vector;
std::vector<int> cell_number_vector;
std::vector<int> face_number_vector;
std::vector<int> node_number_vector;
std::vector<int> cell_owner_vector;
std::vector<int> face_owner_vector;
std::vector<int> node_owner_vector;
const std::vector<RefFaceList>& refFaceList() const
{
return m_ref_face_list;
}
void addRefFaceList(const RefFaceList& ref_face_list)
{
m_ref_face_list.push_back(ref_face_list);
}
ConnectivityDescriptor& operator=(const ConnectivityDescriptor&) = delete;
ConnectivityDescriptor& operator=(ConnectivityDescriptor&&) = delete;
ConnectivityDescriptor() = default;
ConnectivityDescriptor(const ConnectivityDescriptor&) = default;
ConnectivityDescriptor(ConnectivityDescriptor&&) = delete;
~ConnectivityDescriptor() = default;
};
template <size_t Dim>
class Connectivity final
: public IConnectivity
{
public:
PASTIS_INLINE
static std::shared_ptr<Connectivity<Dim>>
build(const ConnectivityDescriptor&);
private:
constexpr static auto& itemTId = ItemTypeId<Dim>::itemTId;
public:
static constexpr size_t Dimension = Dim;
PASTIS_INLINE
size_t dimension() const final
{
return Dimension;
}
private:
ConnectivityMatrix m_item_to_item_matrix[Dimension+1][Dimension+1];
WeakCellValue<const CellType> m_cell_type;
#warning is m_cell_global_index really necessary? should it be computed on demand instead?
WeakCellValue<const int> m_cell_global_index;
WeakCellValue<const int> m_cell_number;
WeakFaceValue<const int> m_face_number;
#warning check that m_edge_number is filled correctly
WeakEdgeValue<const int> m_edge_number;
WeakNodeValue<const int> m_node_number;
WeakCellValue<const int> m_cell_owner;
WeakCellValue<const bool> m_cell_is_owned;
WeakFaceValue<const int> m_face_owner;
WeakFaceValue<const bool> m_face_is_owned;
WeakEdgeValue<const int> m_edge_owner;
WeakEdgeValue<const bool> m_edge_is_owned;
WeakNodeValue<const int> m_node_owner;
WeakNodeValue<const bool> m_node_is_owned;
WeakFaceValuePerCell<const bool> m_cell_face_is_reversed;
WeakNodeValuePerCell<const unsigned short> m_cell_local_numbers_in_their_nodes;
WeakEdgeValuePerCell<const unsigned short> m_cell_local_numbers_in_their_edges;
WeakFaceValuePerCell<const unsigned short> m_cell_local_numbers_in_their_faces;
WeakCellValuePerFace<const unsigned short> m_face_local_numbers_in_their_cells;
WeakEdgeValuePerFace<const unsigned short> m_face_local_numbers_in_their_edges;
WeakNodeValuePerFace<const unsigned short> m_face_local_numbers_in_their_nodes;
WeakCellValuePerEdge<const unsigned short> m_edge_local_numbers_in_their_cells;
WeakFaceValuePerEdge<const unsigned short> m_edge_local_numbers_in_their_faces;
WeakNodeValuePerEdge<const unsigned short> m_edge_local_numbers_in_their_nodes;
WeakCellValuePerNode<const unsigned short> m_node_local_numbers_in_their_cells;
WeakEdgeValuePerNode<const unsigned short> m_node_local_numbers_in_their_edges;
WeakFaceValuePerNode<const unsigned short> m_node_local_numbers_in_their_faces;
ConnectivityComputer m_connectivity_computer;
std::vector<RefFaceList> m_ref_face_list;
std::vector<RefNodeList> m_ref_node_list;
WeakCellValue<const double> m_inv_cell_nb_nodes;
void _computeCellFaceAndFaceNodeConnectivities();
template <typename SubItemValuePerItemType>
PASTIS_INLINE
const SubItemValuePerItemType&
_lazzyBuildItemNumberInTheirChild(const SubItemValuePerItemType& sub_item_value_per_item) const
{
using ReversedItemOfItem = typename SubItemValuePerItemType::ItemOfItemType::Reversed;
if (not sub_item_value_per_item.isBuilt()) {
const_cast<SubItemValuePerItemType&>(sub_item_value_per_item)
= m_connectivity_computer
. computeLocalItemNumberInChildItem<ReversedItemOfItem>(*this);
}
return sub_item_value_per_item;
}
friend class ConnectivityComputer;
PASTIS_INLINE
const ConnectivityMatrix& _getMatrix(const ItemType& item_type_0,
const ItemType& item_type_1) const
{
const ConnectivityMatrix& connectivity_matrix
= m_item_to_item_matrix[itemTId(item_type_0)][itemTId(item_type_1)];
if (not connectivity_matrix.isBuilt()) {
const_cast<ConnectivityMatrix&>(connectivity_matrix)
= m_connectivity_computer
. computeConnectivityMatrix(*this, item_type_0, item_type_1);
}
return connectivity_matrix;
}
public:
PASTIS_INLINE
const auto& cellType() const
{
return m_cell_type;
}
PASTIS_INLINE const auto& cellNumber() const
{
return m_cell_number;
}
PASTIS_INLINE
const auto& faceNumber() const
{
return m_face_number;
}
PASTIS_INLINE
const auto& edgeNumber() const
{
return m_edge_number;
}
PASTIS_INLINE
const auto& nodeNumber() const
{
return m_node_number;
}
template <ItemType item_type>
PASTIS_INLINE
const auto& number() const
{
if constexpr(item_type == ItemType::cell) {
return m_cell_number;
} else if constexpr(item_type == ItemType::face) {
return m_face_number;
} else if constexpr(item_type == ItemType::edge) {
return m_edge_number;
} else if constexpr(item_type == ItemType::node) {
return m_node_number;
} else {
static_assert(is_false_item_type_v<item_type>, "unknown ItemType");
}
}
PASTIS_INLINE
const auto& cellOwner() const
{
return m_cell_owner;
}
PASTIS_INLINE
const auto& faceOwner() const
{
return m_face_owner;
}
PASTIS_INLINE
const auto& edgeOwner() const
{
perr() << __FILE__ << ':' << __LINE__ << ": edge owner not built\n";
std::terminate();
return m_edge_owner;
}
PASTIS_INLINE
const auto& nodeOwner() const
{
return m_node_owner;
}
template <ItemType item_type>
PASTIS_INLINE
const auto& owner() const
{ if constexpr(item_type == ItemType::cell) {
return m_cell_owner;
} else if constexpr(item_type == ItemType::face) {
return m_face_owner;
} else if constexpr(item_type == ItemType::edge) {
return m_edge_owner;
} else if constexpr(item_type == ItemType::node) {
return m_node_owner;
} else {
static_assert(is_false_item_type_v<item_type>, "unknown ItemType");
}
}
PASTIS_INLINE
const auto& cellIsOwned() const
{
return m_cell_is_owned;
}
PASTIS_INLINE
const auto& faceIsOwned() const
{
return m_face_is_owned;
}
PASTIS_INLINE
const auto& edgeIsOwned() const
{
perr() << __FILE__ << ':' << __LINE__ << ": edge is owned not built\n";
std::terminate();
return m_edge_is_owned;
}
PASTIS_INLINE
const auto& nodeIsOwned() const
{
return m_node_is_owned;
}
template <ItemType item_type>
PASTIS_INLINE
const auto& isOwned() const
{
if constexpr(item_type == ItemType::cell) {
return m_cell_is_owned;
} else if constexpr(item_type == ItemType::face) {
return m_face_is_owned;
} else if constexpr(item_type == ItemType::edge) {
return m_edge_is_owned;
} else if constexpr(item_type == ItemType::node) {
return m_node_is_owned;
} else {
static_assert(is_false_item_type_v<item_type>, "unknown ItemType");
}
}
PASTIS_INLINE
const bool& isConnectivityMatrixBuilt(const ItemType& item_type_0,
const ItemType& item_type_1) const
{
const ConnectivityMatrix& connectivity_matrix
= m_item_to_item_matrix[itemTId(item_type_0)][itemTId(item_type_1)];
return connectivity_matrix.isBuilt();
}
template <ItemType source_item_type,
ItemType target_item_type>
PASTIS_INLINE
auto getItemToItemMatrix() const
{
return ItemToItemMatrix<source_item_type, target_item_type>(_getMatrix(source_item_type,
target_item_type));
}
PASTIS_INLINE
auto cellToFaceMatrix() const
{
return this->template getItemToItemMatrix<ItemType::cell, ItemType::face>();
}
PASTIS_INLINE
auto cellToEdgeMatrix() const
{
return this->template getItemToItemMatrix<ItemType::cell, ItemType::edge>();
}
PASTIS_INLINE
auto cellToNodeMatrix() const
{
return this->template getItemToItemMatrix<ItemType::cell, ItemType::node>();
}
PASTIS_INLINE
auto faceToCellMatrix() const
{
return this->template getItemToItemMatrix<ItemType::face, ItemType::cell>();
}
PASTIS_INLINE
auto faceToEdgeMatrix() const
{
return this->template getItemToItemMatrix<ItemType::face, ItemType::edge>();
}
PASTIS_INLINE
auto faceToNodeMatrix() const
{
return this->template getItemToItemMatrix<ItemType::face, ItemType::node>();
}
PASTIS_INLINE
auto edgeToCellMatrix() const
{
return this->template getItemToItemMatrix<ItemType::edge, ItemType::cell>();
}
PASTIS_INLINE
auto edgeToFaceMatrix() const
{
return this->template getItemToItemMatrix<ItemType::edge, ItemType::face>();
}
PASTIS_INLINE
auto edgeToNodeMatrix() const
{
return this->template getItemToItemMatrix<ItemType::edge, ItemType::node>();
}
PASTIS_INLINE
auto nodeToCellMatrix() const
{
return this->template getItemToItemMatrix<ItemType::node, ItemType::cell>();
}
PASTIS_INLINE
auto nodeToFaceMatrix() const
{
return this->template getItemToItemMatrix<ItemType::node, ItemType::face>();
}
PASTIS_INLINE
auto nodeToEdgeMatrix() const
{
return this->template getItemToItemMatrix<ItemType::node, ItemType::edge>();
}
PASTIS_INLINE
const auto& cellFaceIsReversed() const
{
static_assert(Dimension>1, "reversed faces makes no sense in dimension 1");
return m_cell_face_is_reversed;
}
PASTIS_INLINE
const auto& cellLocalNumbersInTheirNodes() const
{
return _lazzyBuildItemNumberInTheirChild(m_cell_local_numbers_in_their_nodes);
}
PASTIS_INLINE
const auto& cellLocalNumbersInTheirEdges() const
{
if constexpr (Dimension>2) {
return _lazzyBuildItemNumberInTheirChild(m_cell_local_numbers_in_their_edges);
} else {
return cellLocalNumbersInTheirFaces();
}
}
PASTIS_INLINE
const auto& cellLocalNumbersInTheirFaces() const
{
if constexpr (Dimension>1) {
return _lazzyBuildItemNumberInTheirChild(m_cell_local_numbers_in_their_faces);
} else {
return cellLocalNumbersInTheirNodes();
}
}
PASTIS_INLINE
const auto& faceLocalNumbersInTheirCells() const
{
if constexpr(Dimension>1) {
return _lazzyBuildItemNumberInTheirChild(m_face_local_numbers_in_their_cells);
} else {
return nodeLocalNumbersInTheirCells();
}
}
PASTIS_INLINE
const auto& faceLocalNumbersInTheirEdges() const
{
static_assert(Dimension>2,"this function has no meaning in 1d or 2d");
return _lazzyBuildItemNumberInTheirChild(m_face_local_numbers_in_their_edges);
}
PASTIS_INLINE
const auto& faceLocalNumbersInTheirNodes() const
{
static_assert(Dimension>1,"this function has no meaning in 1d");
return _lazzyBuildItemNumberInTheirChild(m_face_local_numbers_in_their_nodes);
}
PASTIS_INLINE
const auto& edgeLocalNumbersInTheirCells() const
{
if constexpr (Dimension>2) {
return _lazzyBuildItemNumberInTheirChild(m_edge_local_numbers_in_their_cells);
} else {
return faceLocalNumbersInTheirCells(); }
}
PASTIS_INLINE
const auto& edgeLocalNumbersInTheirFaces() const
{
static_assert(Dimension>2, "this function has no meaning in 1d or 2d");
return _lazzyBuildItemNumberInTheirChild(m_edge_local_numbers_in_their_faces);
}
PASTIS_INLINE
const auto& edgeLocalNumbersInTheirNodes() const
{
static_assert(Dimension>2, "this function has no meaning in 1d or 2d");
return _lazzyBuildItemNumberInTheirChild(m_edge_local_numbers_in_their_nodes);
}
PASTIS_INLINE
const auto& nodeLocalNumbersInTheirCells() const
{
return _lazzyBuildItemNumberInTheirChild(m_node_local_numbers_in_their_cells);
}
PASTIS_INLINE
const auto& nodeLocalNumbersInTheirEdges() const
{
static_assert(Dimension>2, "this function has no meaning in 1d or 2d");
return _lazzyBuildItemNumberInTheirChild(m_node_local_numbers_in_their_edges);
}
PASTIS_INLINE
const auto& nodeLocalNumbersInTheirFaces() const
{
static_assert(Dimension>1,"this function has no meaning in 1d");
return _lazzyBuildItemNumberInTheirChild(m_node_local_numbers_in_their_faces);
}
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];
}
PASTIS_INLINE
CSRGraph cellToCellGraph() const
{
std::vector<std::set<int>> cell_cells(this->numberOfCells());
if constexpr (true) {
const auto& face_to_cell_matrix
= this->faceToCellMatrix();
for (FaceId l=0; l<this->numberOfFaces(); ++l) {
const auto& face_to_cell = face_to_cell_matrix[l];
if (face_to_cell.size() > 1) {
const CellId cell_0 = face_to_cell[0];
const CellId cell_1 = face_to_cell[1];
cell_cells[cell_0].insert(cell_1);
cell_cells[cell_1].insert(cell_0);
}
}
} else {
const auto& node_to_cell_matrix
= this->nodeToCellMatrix();
for (NodeId l=0; l<this->numberOfNodes(); ++l) {
const auto& node_to_cell = node_to_cell_matrix[l];
for (size_t i_cell=0; i_cell<node_to_cell.size(); ++i_cell) {
const CellId cell_0 = node_to_cell[i_cell];
for (size_t j_cell=0; j_cell<i_cell; ++j_cell) {
const CellId cell_1 = node_to_cell[j_cell];
cell_cells[cell_0].insert(cell_1);
cell_cells[cell_1].insert(cell_0);
}
}
}
}
Array<int> entries(this->numberOfCells()+1);
entries[0]=0;
for (size_t j=0; j<this->numberOfCells(); ++j) {
entries[j+1] = entries[j]+cell_cells[j].size();
}
Array<int> neighbors(entries[this->numberOfCells()]);
{
size_t k=0;
for (size_t j=0; j<this->numberOfCells(); ++j) {
for (CellId cell_id : cell_cells[j]) {
neighbors[k] = m_cell_global_index[cell_id];
++k;
}
}
}
return CSRGraph(entries, neighbors);
}
PASTIS_INLINE
size_t numberOfNodes() const final
{
const auto& node_to_cell_matrix
= this->_getMatrix(ItemType::node,ItemType::cell);
return node_to_cell_matrix.numRows();
}
PASTIS_INLINE
size_t numberOfEdges() const final
{
const auto& edge_to_node_matrix
= this->_getMatrix(ItemType::edge,ItemType::node);
return edge_to_node_matrix.numRows();
}
PASTIS_INLINE
size_t numberOfFaces() const final
{
const auto& face_to_node_matrix
= this->_getMatrix(ItemType::face,ItemType::cell);
return face_to_node_matrix.numRows();
}
PASTIS_INLINE size_t numberOfCells() const final
{
const auto& cell_to_node_matrix
= this->_getMatrix(ItemType::cell,ItemType::node);
return cell_to_node_matrix.numRows();
}
CellValue<const double> invCellNbNodes() const
{
#warning add calculation on demand when variables will be defined
return m_inv_cell_nb_nodes;
}
Connectivity(const Connectivity&) = delete;
private:
Connectivity();
void _buildFrom(const ConnectivityDescriptor& descriptor);
public:
~Connectivity()
{
auto& manager = SynchronizerManager::instance();
manager.deleteConnectivitySynchronizer(this);
}
};
template <size_t Dim>
PASTIS_INLINE
std::shared_ptr<Connectivity<Dim>>
Connectivity<Dim>::build(const ConnectivityDescriptor& descriptor)
{
std::shared_ptr<Connectivity<Dim>> connectivity_ptr(new Connectivity<Dim>);
connectivity_ptr->_buildFrom(descriptor);
return connectivity_ptr;
}
using Connectivity3D = Connectivity<3>;
using Connectivity2D = Connectivity<2>;
using Connectivity1D = Connectivity<1>;
#endif // CONNECTIVITY_HPP