use std::collections::VecDeque;

use itertools::Itertools;
use log::debug;
use serde::{Deserialize, Serialize};

use crate::graph::Graph;

use super::{simplices::CausalSimplex2 as Triangle, triangulation2d_edge::Triangulation2D};

/// Triangulation structure designed to store a causal "gluing" of triangles.
///
/// This can store any orientable triangulation where every edge is between exactly two triangles.
/// Internally three half-edges of each triangle are stored with a pointer to the half-edge
/// they are adjacent to. Additionally this stores the vertex labels of the three vertices of
/// each triangle opposite the half-edges.
#[derive(Debug, Default, Clone, PartialEq, Eq, Hash, Serialize, Deserialize)]
pub struct CausalTriangulation2D {
    /// Number of time-slices / spatial slices in the triangulation
    pub ntime: u16,
    /// The time labels of each vertex
    pub vertex_times: Vec<u16>,
    /// List holding triangles of the triangulation.
    ///
    /// Note that the `neighbours` field of each triangle is used to store the indices of
    /// the half-edge each half-edge is adjacent to. This encodes both the label of the
    /// neighbouring triangle and the back index (the index in the `neighbour` list back).
    ///
    /// The encoding works as follows, given the "half-edge" label `i`:
    /// - `triangle_label = i / 3`
    /// - `edge_index = i % 3`
    ///
    /// Additionally the `edge_index = 0` must correspond to the time-like neighbour.
    pub triangles: Vec<Triangle>,
}

impl CausalTriangulation2D {
    /// Creates and empty triangulation
    pub fn empty() -> Self {
        Self::default()
    }
}

impl CausalTriangulation2D {
    pub fn to_triangulation2d(&self) -> Triangulation2D {
        Triangulation2D {
            num_vertices: self.num_vertices(),
            triangles: self.triangles.iter().map(From::from).collect(),
        }
    }

    pub fn into_triangulation2d(self) -> Triangulation2D {
        Triangulation2D {
            num_vertices: self.num_vertices(),
            triangles: self.triangles.into_iter().map(From::from).collect(),
        }
    }
}

impl CausalTriangulation2D {
    /// Returns the number of time-slices/spatial slices in the triangulation
    #[inline]
    pub fn num_timeslices(&self) -> u16 {
        self.ntime
    }

    /// Returns the number of vertices in the triangulation
    #[inline]
    pub fn num_vertices(&self) -> usize {
        self.vertex_times.len()
    }

    /// Returns the number of triangles in the triangulation
    #[inline]
    pub fn num_triangles(&self) -> usize {
        // Division must be exact for a valid triangulation
        self.triangles.len()
    }

    /// Returns the number of edges in the triangulation
    ///
    /// Note: this is simply `num_triangles * 3 / 2`
    #[inline]
    pub fn num_edges(&self) -> usize {
        // Division must be exact for a valid triangulation
        self.triangles.len() * 3 / 2
    }

    /// Returns the number of half edges of the triangulation
    ///
    /// Note: this is simply `num_triangles * 3`
    #[inline]
    pub fn num_halfedges(&self) -> usize {
        self.triangles.len() * 3
    }
}

/// Methods and associated functions to access and work with half-edges
impl CausalTriangulation2D {
    /// Returns the label of the adjacent `halfedge`
    #[inline]
    pub fn adjacent_halfedge(&self, halfedge: usize) -> usize {
        self.triangles[halfedge / 3].neighbours[halfedge % 3]
    }

    /// Returns the label of the triangle the `halfedge` is part of.
    ///
    /// Note: this is simply `halfedge / 3`
    #[inline]
    pub fn triangle_label(halfedge: usize) -> usize {
        halfedge / 3
    }

    /// Returns the index of the halfedge in the corresponding triangulation.
    ///
    /// This is simply `halfedge % 3` and means the index of neighbouring triangle
    /// in the `neighbours` list of [`Triangle`] corresponding to the half-edge.
    #[inline]
    pub fn halfedge_index(halfedge: usize) -> usize {
        halfedge % 3
    }

    /// Returns the label of the next `halfedge`.
    ///
    /// This is simply `3 * (halfedge / 3) + (halfedge + 1) % 3`
    #[inline]
    pub fn next_halfedge(halfedge: usize) -> usize {
        3 * (halfedge / 3) + (halfedge + 1) % 3
    }

    /// Returns the label of the previous `halfedge`.
    ///
    /// This is simply `3 * (halfedge / 3) + (halfedge + 2) % 3`
    #[inline]
    pub fn prev_halfedge(halfedge: usize) -> usize {
        3 * (halfedge / 3) + (halfedge + 2) % 3
    }

    /// Returns the label of the halfedge given by a `triangle_label` and `halfedge_index`
    ///
    /// Note that `halfedge_index` must be 0, 1, or 2. If another index is used, this function
    /// will give a result, but it does not correspond to a half-edge of the given triangle.
    /// Note this function is simply `3 * triangle_label + halfedge_index`
    #[inline]
    pub fn halfedge_label(triangle_label: usize, halfedge_index: usize) -> usize {
        3 * triangle_label + halfedge_index
    }
}

impl CausalTriangulation2D {
    /// Returns a reference to the triangle list of the [`Triangulation2D`]
    #[inline]
    pub fn get_triangles(&self) -> &[Triangle] {
        &self.triangles
    }

    /// Returns a reference to time labels of the vertices of the [`Triangulation2D`]
    #[inline]
    pub fn get_vertex_times(&self) -> &[u16] {
        &self.vertex_times
    }

    /// Returns a reference to the triangle given by `label`
    ///
    /// Note that the `neighbours` field hold the labels of the adjacent half-edges, encoding
    /// the triangle label as `i / 3` and the back-index as `i % 3`. Also see [`Triangulation2D`].
    ///
    /// # Panic
    /// This method will panic if `label >= self.num_triangles()`
    #[inline]
    pub fn get_triangle(&self, label: usize) -> &Triangle {
        &self.triangles[label]
    }

    /// Returns a mutable reference to the triangle given by `label`
    ///
    /// See [`Triangulation2D::get_triangle()`] for more details.
    ///
    /// # Panic
    /// This method will panic if `label >= self.num_triangles()`
    #[inline]
    pub fn get_mut_triangle(&mut self, label: usize) -> &mut Triangle {
        &mut self.triangles[label]
    }

    /// Returns the label of the triangle neighbouring `label` on the `index` side.
    ///
    /// # Panic
    /// This method will panic if `label >= self.num_triangles()` or `index >= 3`
    #[inline]
    pub fn get_neighbour(&self, label: usize, index: usize) -> usize {
        self.get_triangle(label).get_neighbour(index)
    }

    /// Returns the label and backindex of triangle neighbouring `label` on the `index` side.
    ///
    /// # Panic
    /// This method will panic if `label >= self.num_triangles()` or `index >= 3`
    #[inline]
    pub fn get_neighbour_backindex(&self, label: usize, index: usize) -> (usize, usize) {
        self.get_triangle(label).get_neighbour_backindex(index)
    }

    /// Returns the neighbouring triangle labels of the triangle given by `label`
    ///
    /// # Panic
    /// This method will panic if `label >= self.num_triangles()`
    #[inline]
    pub fn get_neighbours(&self, label: usize) -> [usize; 3] {
        self.get_triangle(label).get_neighbours()
    }

    #[inline]
    pub fn get_neighbours_backindices(&self, label: usize) -> [(usize, usize); 3] {
        self.get_triangle(label).get_neighbours_backindices()
    }

    /// Returns the vertices of the triangle by `label`
    ///
    /// # Panic
    /// This method will panic if `label >= self.num_triangles()`
    #[inline]
    pub fn get_triangle_vertices(&self, label: usize) -> [usize; 3] {
        self.get_triangle(label).vertices
    }
}

impl CausalTriangulation2D {
    /// Returns the time-slice/spatial slab volume profile, i.e. the number of vertices
    /// in each timeslice.
    pub fn timeslice_volume_profile(&self) -> Vec<usize> {
        let mut profile = vec![0; self.ntime as usize];
        for t in self.vertex_times.iter().copied() {
            profile[t as usize] += 1;
        }
        profile
    }

    /// Returns the spatial slab volume profile, i.e. the number of triangles in each
    /// slab/thick timeslice.
    ///
    /// Note: for 2D CDT this should be equal to the moving sum of the
    /// [`timeslice_volume_profile()`](CausalTriangulation2D::timeslice_volume_profile)
    /// with a window of 2.
    pub fn timeslab_volume_profile(&self) -> Vec<usize> {
        let mut profile = vec![0; self.ntime as usize];
        for triangle in &self.triangles {
            profile[triangle.t as usize] += 1;
        }
        profile
    }
}

impl CausalTriangulation2D {
    /// Returns oriented dual graph of the [`Triangulation`]
    ///
    /// Note that this graph is oriented in the sense that the links are given in counter-clockwise
    /// (positive) order.
    pub fn dual_graph(&self) -> Graph {
        Graph::new_from_iter(
            self.triangles
                .iter()
                .map(|triangle| triangle.neighbours.iter().map(|&nbr_edge| nbr_edge / 3)),
        )
    }

    /// Returns vertex graph of the [`Triangulation`]
    ///
    /// Note that this graph is not necessarily oriented, i.e. the ordering of the neighbours
    /// of each node is arbitrary.
    pub fn vertex_graph(&self) -> Graph {
        let mut vertex_neighbours = vec![Vec::new(); self.num_vertices()];
        for triangle in self.triangles.iter() {
            let [v0, v1, v2] = triangle.vertices;

            vertex_neighbours[v0].push(v1);
            vertex_neighbours[v1].push(v2);
            vertex_neighbours[v2].push(v0);
        }

        Graph::new_from_iter(vertex_neighbours.into_iter().map(|nbrs| nbrs.into_iter()))
    }

    /// Returns vertex graph of the [`Triangulation`] explicity making all edges single edges
    pub fn vertex_graph_noduplicates(&self) -> Graph {
        let mut vertex_neighbours = vec![Vec::new(); self.num_vertices()];
        for triangle in self.triangles.iter() {
            let [v0, v1, v2] = triangle.vertices;

            vertex_neighbours[v0].push(v1);
            vertex_neighbours[v1].push(v2);
            vertex_neighbours[v2].push(v0);
        }

        Graph::new_from_iter(
            vertex_neighbours
                .into_iter()
                .map(|nbrs| nbrs.into_iter().unique()),
        )
    }
}

/// Checks to check the validity of the [`Triangulation2D`] structure
// TODO: Return not `bool`, but `Result<(), SomeProperError>`
impl CausalTriangulation2D {
    pub fn check_validity(&self) -> bool {
        self.check_biconnectivity()
            && self.check_num_vertices()
            && self.check_neighbour_ordering()
            && self.check_vertex_consistency()
            && self.check_vertex_times()
    }

    /// Check if the [`Triangulation2D`] is properly biconnective, i.e. that each triangle
    /// neighbour points back the the triangle at the correct index.
    pub fn check_biconnectivity(&self) -> bool {
        for (label, triangle) in self.triangles.iter().enumerate() {
            for (nbr_label, backindex) in triangle.get_neighbours_backindices() {
                if self.get_neighbour(nbr_label, backindex) != label {
                    debug!("Neighbour {nbr_label} is not correctly connected back to {label}.");
                    return false;
                }
            }
        }
        true
    }

    /// Count the number of vertices by looking at all the vertex labels on each triangle.
    ///
    /// Note: if the [`Triangulation2D`] is a correct structure, this should be equal to
    /// `self.num_vertices()` which is just a direct lookup. Whereas this method can be slow.
    pub fn count_vertices(&self) -> usize {
        self.triangles
            .iter()
            .flat_map(|triangle| triangle.vertices.into_iter())
            .unique()
            .count()
    }

    /// Checks if `num_vertices` is set correctly based on the vertices of each triangle
    pub fn check_num_vertices(&self) -> bool {
        self.count_vertices() == self.num_vertices()
    }

    /// Checks whether the two vertices between each triangle are also correctly ordered vertices
    /// of the corresponding neighbouring triangle.
    pub fn check_vertex_consistency(&self) -> bool {
        for triangle in self.triangles.iter() {
            for (index, (nbr_label, backindex)) in triangle
                .get_neighbours_backindices()
                .into_iter()
                .enumerate()
            {
                let vertices = triangle.vertices;
                let nbr_vertices = self.get_triangle_vertices(nbr_label);

                if vertices[(index + 1) % 3] != nbr_vertices[(backindex + 2) % 3] {
                    debug!("Triangulation does not have correctly ordered vertices");
                    return false;
                }

                if vertices[(index + 2) % 3] != nbr_vertices[(backindex + 1) % 3] {
                    debug!("Triangulation does not have correctly ordered vertices");
                    return false;
                }
            }
        }

        true
    }

    /// Returns whether the triangulation obeys the absolute ordering of triangle neighbours
    /// such that the first (index 0) neighbour is the time-like neighbour.
    ///
    /// See also notes on fields of [`CausalTriangle`](super::simplices::CausalTriangle).
    pub fn check_neighbour_ordering(&self) -> bool {
        for triangle in self.triangles.iter() {
            let t_self = triangle.t;
            let t_nbr = self.triangles[triangle.get_neighbour(0)].t;

            use super::simplices::CausalSimplex2Kind as TriangleType;
            match triangle.kind {
                TriangleType::S21 => {
                    if (t_nbr + 1) % self.ntime != t_self {
                        return false;
                    }
                }
                TriangleType::S12 => {
                    if t_nbr != (t_self + 1) % self.ntime {
                        return false;
                    }
                }
            }
        }
        true
    }

    /// Returns whether the triangulation obeys the absolute ordering of triangle neighbours
    /// such that the first (index 0) neighbour is the time-like neighbour.
    ///
    /// See also notes on fields of [`CausalTriangle`](super::simplices::CausalTriangle).
    pub fn check_vertex_times(&self) -> bool {
        for triangle in self.triangles.iter() {
            let t_self = triangle.t;
            let t_next = (t_self + 1) % self.ntime;

            use super::simplices::CausalSimplex2Kind as TriangleType;
            match triangle.kind {
                TriangleType::S21 => {
                    if self.vertex_times[triangle.vertices[0]] != t_next {
                        return false;
                    }
                    if self.vertex_times[triangle.vertices[1]] != t_self {
                        return false;
                    }
                    if self.vertex_times[triangle.vertices[2]] != t_self {
                        return false;
                    }
                }
                TriangleType::S12 => {
                    if self.vertex_times[triangle.vertices[0]] != t_self {
                        return false;
                    }
                    if self.vertex_times[triangle.vertices[1]] != t_next {
                        return false;
                    }
                    if self.vertex_times[triangle.vertices[2]] != t_next {
                        return false;
                    }
                }
            }
        }
        true
    }
}

impl CausalTriangulation2D {
    pub fn get_future_neighbour(&self, label: usize) -> Option<usize> {
        let triangle = self.get_triangle(label);
        match triangle.kind {
            super::simplices::CausalSimplex2Kind::S21 => None,
            super::simplices::CausalSimplex2Kind::S12 => Some(triangle.get_neighbour(0)),
        }
    }

    pub fn get_past_neighbour(&self, label: usize) -> Option<usize> {
        let triangle = self.get_triangle(label);
        match triangle.kind {
            super::simplices::CausalSimplex2Kind::S21 => Some(triangle.get_neighbour(0)),
            super::simplices::CausalSimplex2Kind::S12 => None,
        }
    }

    pub fn get_timelike_neighbour(&self, label: usize) -> usize {
        self.get_triangle(label).get_neighbour(0)
    }

    pub fn get_left_neighbour(&self, label: usize) -> usize {
        let triangle = self.get_triangle(label);
        match triangle.kind {
            super::simplices::CausalSimplex2Kind::S21 => triangle.get_neighbour(2),
            super::simplices::CausalSimplex2Kind::S12 => triangle.get_neighbour(1),
        }
    }

    pub fn get_right_neighbour(&self, label: usize) -> usize {
        let triangle = self.get_triangle(label);
        match triangle.kind {
            super::simplices::CausalSimplex2Kind::S21 => triangle.get_neighbour(1),
            super::simplices::CausalSimplex2Kind::S12 => triangle.get_neighbour(2),
        }
    }

    pub fn get_spacelike_neighbours(&self, label: usize) -> impl Iterator<Item = usize> {
        self.get_triangle(label)
            .get_neighbours()
            .into_iter()
            .skip(1)
    }
}

impl CausalTriangulation2D {
    /// Walk along triangles until the original triangle is reached
    ///
    /// Restrict to stay between the origin timeslice and `ntime` further timeslices, find
    /// the shortest loop under this restriction.
    ///
    /// Note: the returned path is in the past direction, ending not stating with the `origin`
    pub fn identify_timelike_loop(&self, origin: usize) -> Vec<usize> {
        let origin_triangle = self.get_triangle(origin);
        // Make sure the starting triangle is of 21-type, get neighbour if not.
        let (origin, origin_triangle) = match origin_triangle.kind {
            super::simplices::CausalSimplex2Kind::S21 => (origin, origin_triangle),
            super::simplices::CausalSimplex2Kind::S12 => {
                let upper = origin_triangle.get_neighbour(0);
                (upper, self.get_triangle(upper))
            }
        };
        let t_origin = origin_triangle.t;

        let mut visited = vec![None; self.num_triangles()];
        visited[origin] = Some(usize::MAX);

        let mut queue = VecDeque::with_capacity(2 * self.num_triangles() / (self.ntime as usize));
        queue.push_back(origin);

        'outer: loop {
            let label = queue
                .pop_front()
                .expect("Loop should always be found in proper CDT");
            if let Some(nbr) = self.get_future_neighbour(label) {
                if nbr == origin {
                    // Went back to original triangle
                    visited[nbr] = Some(label);
                    break 'outer;
                }

                if self.triangles[nbr].t != t_origin && visited[nbr].is_none() {
                    visited[nbr] = Some(label);
                    {
                        queue.push_back(nbr);
                    }
                }
            }
            for nbr in self.get_spacelike_neighbours(label) {
                // Only consider unvisited nodes
                if visited[nbr].is_some() {
                    continue;
                }
                visited[nbr] = Some(label);
                queue.push_back(nbr);
            }
        }

        // Rewalk path
        let mut path = Vec::with_capacity((self.ntime * 2).into());
        let mut current_triangle = origin;
        loop {
            current_triangle = visited[current_triangle]
                .expect("The path should have been fully explored by construction");
            path.push(current_triangle);
            if current_triangle == origin {
                break;
            }
        }

        path
    }
}
