#ifndef CONNECTIVITY_1D_HPP
#define CONNECTIVITY_1D_HPP

#include <Kokkos_Core.hpp>
#include <PastisAssert.hpp>

#include <TinyVector.hpp>
#include <ConnectivityUtils.hpp>

#include <RefId.hpp>
#include <RefNodeList.hpp>

class Connectivity1D
{
public:
  static constexpr size_t dimension = 1;
  ConnectivityMatrix m_node_to_cell_matrix;
  ConnectivityMatrix m_node_to_cell_local_node_matrix;
private:
  std::vector<RefNodeList> m_ref_node_list;

  size_t  m_number_of_nodes;
  const size_t& m_number_of_faces = m_number_of_nodes;
  const size_t  m_number_of_cells;

  const Kokkos::View<const unsigned int**> m_cell_nodes;
  const Kokkos::View<const unsigned int**>& m_cell_faces = m_cell_nodes;

  const Kokkos::View<const unsigned short*> m_cell_nb_nodes;
  Kokkos::View<double*> m_inv_cell_nb_nodes;
  const Kokkos::View<const unsigned short*>& m_cell_nb_faces = m_cell_nb_nodes;

  Kokkos::View<const unsigned short*> m_node_nb_cells;
  const Kokkos::View<const unsigned short*>& m_face_nb_cells = m_node_nb_cells;

  Kokkos::View<const unsigned int**> m_node_cells;
  const Kokkos::View<const unsigned int**>& m_face_cells = m_node_cells;

  Kokkos::View<const unsigned short**> m_node_cell_local_node;
  const Kokkos::View<const unsigned short**>& m_face_cell_local_face = m_node_cell_local_node;

  size_t  m_max_nb_node_per_cell;

  const Kokkos::View<const unsigned int**>
  _buildCellNodes(const size_t& number_of_cells)
  {
    Kokkos::View<unsigned int*[2]> cell_nodes("cell_nodes", number_of_cells);

    Kokkos::parallel_for(m_number_of_cells, KOKKOS_LAMBDA(const size_t& j) {
        cell_nodes(j,0) = j;
        cell_nodes(j,1) = j+1;
      });

    return cell_nodes;
  }


  const Kokkos::View<const unsigned short*>
  _buildCellNbNodes(const size_t& number_of_cells)
  {
    Kokkos::View<unsigned short*>  cell_nb_nodes("cell_nb_nodes", number_of_cells);
    Kokkos::parallel_for(m_number_of_cells, KOKKOS_LAMBDA(const size_t& j) {
        cell_nb_nodes[j] = 2;
      });

    return cell_nb_nodes;
  }

public:
  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];
  }

  const size_t& numberOfNodes() const
  {
    return m_number_of_nodes;
  }

  const size_t& numberOfFaces() const
  {
    return m_number_of_faces;
  }

  const size_t& numberOfCells() const
  {
    return m_number_of_cells;
  }

  const size_t& maxNbNodePerCell() const
  {
    return m_max_nb_node_per_cell;
  }

  const Kokkos::View<const unsigned int**> cellNodes() const
  {
    return m_cell_nodes;
  }

  const Kokkos::View<const unsigned int**> cellFaces() const
  {
    return m_cell_faces;
  }

  const Kokkos::View<const unsigned short*> nodeNbCells() const
  {
    return m_node_nb_cells;
  }

  const Kokkos::View<const unsigned short*> cellNbNodes() const
  {
    return m_cell_nb_nodes;
  }

  const Kokkos::View<const double*> invCellNbNodes() const
  {
    return m_inv_cell_nb_nodes;
  }

  const Kokkos::View<const unsigned short*> cellNbFaces() const
  {
    return m_cell_nb_faces;
  }

  const Kokkos::View<const unsigned short*> faceNbCells() const
  {
    return m_face_nb_cells;
  }

  const Kokkos::View<const unsigned int**> nodeCells() const
  {
    return m_node_cells;
  }

  const Kokkos::View<const unsigned int**> faceCells() const
  {
    return m_face_cells;
  }

  const Kokkos::View<const unsigned short**> nodeCellLocalNode() const
  {
    return m_node_cell_local_node;
  }

  const Kokkos::View<const unsigned short**> faceCellLocalFace() const
  {
    return m_face_cell_local_face;
  }

  Connectivity1D(const Connectivity1D&) = delete;

  Connectivity1D(const size_t& number_of_cells)
    : m_number_of_cells (number_of_cells),
      m_cell_nodes (_buildCellNodes(number_of_cells)),
      m_cell_nb_nodes (_buildCellNbNodes(number_of_cells)),
      m_inv_cell_nb_nodes ("inv_cell_nb_nodes", m_number_of_cells)
  {
    Assert(number_of_cells>0);
    Kokkos::parallel_for(m_number_of_cells, KOKKOS_LAMBDA(const size_t& j) {
        m_inv_cell_nb_nodes[j] = 1./m_cell_nb_nodes[j];
      });

    ConnectivityUtils utils;
    utils.computeNodeCellConnectivity(m_cell_nodes,
                                      m_cell_nb_nodes,
                                      m_number_of_cells,
                                      m_max_nb_node_per_cell,
                                      m_number_of_nodes,
                                      m_node_nb_cells,
                                      m_node_cells,
                                      m_node_cell_local_node,
                                      m_node_to_cell_matrix,
                                      m_node_to_cell_local_node_matrix);
  }

  Connectivity1D(const Kokkos::View<const unsigned short*> cell_nb_nodes,
                 const Kokkos::View<const unsigned int**> cell_nodes)
    : m_number_of_cells (cell_nb_nodes.extent(0)),
      m_cell_nodes (cell_nodes),
      m_cell_nb_nodes (cell_nb_nodes),
      m_inv_cell_nb_nodes ("inv_cell_nb_nodes", m_number_of_cells)
  {
    Kokkos::parallel_for(m_number_of_cells, KOKKOS_LAMBDA(const size_t& j) {
        m_inv_cell_nb_nodes[j] = 1./m_cell_nb_nodes[j];
      });

    ConnectivityUtils utils;
    utils.computeNodeCellConnectivity(m_cell_nodes,
                                      m_cell_nb_nodes,
                                      m_number_of_cells,
                                      m_max_nb_node_per_cell,
                                      m_number_of_nodes,
                                      m_node_nb_cells,
                                      m_node_cells,
                                      m_node_cell_local_node,
                                      m_node_to_cell_matrix,
                                      m_node_to_cell_local_node_matrix);
  }

  ~Connectivity1D()
  {
    ;
  }
};

#endif // CONNECTIVITY_1D_HPP