//!
//! Tetmesh module. Describes tetrahedron mesh data structures and possible
//! operations on them.
//!
//! The root module defines the most basic tetmesh that other tetmeshes can extend.
//!

mod surface;
mod extended;

pub use extended::*;
pub use surface::*;

use crate::mesh::attrib::*;
use crate::mesh::topology::*;
use crate::mesh::vertex_positions::*;
use crate::prim::Tetrahedron;
use crate::Real;
use std::slice::{Iter, IterMut};
use math::{SquareMatrix, Matrix3};

#[cfg(feature = "serde")]
use serde::{Deserialize, Serialize};

/// A basic mesh composed of tetrahedra. This mesh is based on vertex positions and a list of
/// vertex indices representing tetrahedra.
#[derive(Clone, Debug, PartialEq, Attrib, Intrinsic)]
#[cfg_attr(feature = "serde", derive(Serialize, Deserialize))]
pub struct TetMesh<T: Real> {
    /// Vertex positions.
    #[intrinsic(VertexPositions)]
    pub vertex_positions: IntrinsicAttribute<[T; 3], VertexIndex>,
    /// Quadruples of indices into `vertices` representing tetrahedra.
    /// The canonical non-inverted tetrahedron is indexed as follows:
    /// ```verbatim
    ///       3
    ///      /|\
    ///     / | \
    ///    /  |  \
    ///  2/...|...\0
    ///   \   |   /
    ///    \  |  /
    ///     \ | /
    ///      \|/
    ///       1
    /// ```
    /// (the dotted line is behind the dashed).
    pub indices: IntrinsicAttribute<[usize; 4], CellIndex>,
    /// Vertex Attributes.
    pub vertex_attributes: AttribDict<VertexIndex>,
    /// Cell Attributes.
    pub cell_attributes: AttribDict<CellIndex>,
    /// Cell vertex Attributes.
    pub cell_vertex_attributes: AttribDict<CellVertexIndex>,
    /// Cell face Attributes.
    pub cell_face_attributes: AttribDict<CellFaceIndex>,
}

impl<T: Real> TetMesh<T> {
    /// This constant defines the triangle faces of each tet. The rule is that `i`th face of the
    /// tet is the one opposite to the `i`th vertex. The triangle starts with the smallest index.
    pub const TET_FACES: [[usize; 3]; 4] = [[1, 3, 2], [0, 2, 3], [0, 3, 1], [0, 1, 2]];

    pub fn new(verts: Vec<[T; 3]>, indices: Vec<usize>) -> TetMesh<T> {
        assert_eq!(indices.len() % 4, 0); // check that we have the right number of indices.

        TetMesh {
            vertex_positions: IntrinsicAttribute::from_vec(verts),
            indices: IntrinsicAttribute::from_vec(crate::into_chunked_vec4(indices)),
            vertex_attributes: AttribDict::new(),
            cell_attributes: AttribDict::new(),
            cell_vertex_attributes: AttribDict::new(),
            cell_face_attributes: AttribDict::new(),
        }
    }

    #[inline]
    pub fn cell_iter(&self) -> Iter<[usize; 4]> {
        self.indices.iter()
    }

    #[cfg(feature = "rayon")]
    #[inline]
    pub fn cell_par_iter(&self) -> rayon::slice::Iter<[usize; 4]> {
        self.indices.par_iter()
    }

    #[inline]
    pub fn cell_iter_mut(&mut self) -> IterMut<[usize; 4]> {
        self.indices.iter_mut()
    }

    #[cfg(feature = "rayon")]
    #[inline]
    pub fn cell_par_iter_mut(&mut self) -> rayon::slice::IterMut<[usize; 4]> {
        self.indices.par_iter_mut()
    }

    /// Cell accessor. These are vertex indices.
    #[inline]
    pub fn cell<CI: Into<CellIndex>>(&self, cidx: CI) -> &[usize; 4] {
        &self.indices[cidx.into()]
    }

    /// Return a slice of individual cells.
    #[inline]
    pub fn cells(&self) -> &[[usize; 4]] {
        self.indices.as_slice()
    }

    /// Tetrahedron iterator.
    #[inline]
    pub fn tet_iter<'a>(&'a self) -> impl Iterator<Item = Tetrahedron<T>> + 'a {
        self.cell_iter().map(move |tet| self.tet_from_indices(tet))
    }

    /// Get a tetrahedron primitive corresponding to the given vertex indices.
    #[inline]
    pub fn tet_from_indices(&self, indices: &[usize; 4]) -> Tetrahedron<T> {
        Tetrahedron::from_indexed_slice(indices, self.vertex_positions.as_slice())
    }

    /// Get a tetrahedron primitive corresponding to the given cell index.
    #[inline]
    pub fn tet<CI: Into<CellIndex>>(&self, cidx: CI) -> Tetrahedron<T> {
        self.tet_from_indices(self.cell(cidx))
    }

    /// Consumes the current mesh to produce a mesh with inverted tetrahedra.
    #[inline]
    pub fn inverted(mut self) -> TetMesh<T> {
        self.invert();
        self
    }

    /// A helper function to invert a single tet of the tetmesh. This keeps the inversion
    /// consistent among all methods that use it (e.g. `invert` and `canonicalize`).
    #[inline]
    fn invert_tet_cell(cell: &mut [usize; 4]) {
        cell.swap(2, 3);
    }

    /// Non consuming verion of the `inverted` function which simply modifies the given mesh.
    #[inline]
    pub fn invert(&mut self) {
        for cell in self.indices.iter_mut() {
            Self::invert_tet_cell(cell);
        }
    }

    /// Convert this mesh into canonical form. This function inverts any inverted tetrahedron such
    /// that all tetrahedra are in canonical (non-inverted) form. The canonical form is determined
    /// by the shape matrix determinant of each tetrahedron. Canonical tetrahedra have a positive
    /// shape matrix determinant (see the `gut::ops::ShapeMatrix` trait and the
    /// `gut::prim::tetrahedron` module).
    #[inline]
    pub fn canonicalized(mut self) -> TetMesh<T> {
        self.canonicalize();
        self
    }

    /// Convert this mesh into canonical form. This function inverts any inverted tetrahedron such
    /// that all tetrahedra are in canonical (non-inverted) form. This is a non-consuming version
    /// of `canonicalized`.
    #[inline]
    pub fn canonicalize(&mut self) {
        use crate::ops::ShapeMatrix;

        let TetMesh {
            ref vertex_positions,
            ref mut indices,
            ..
        } = *self;

        for cell in indices.iter_mut() {
            let tet = Tetrahedron::from_indexed_slice(cell, vertex_positions.as_slice());
            if Matrix3::from(tet.shape_matrix()).determinant() < T::zero() {
                Self::invert_tet_cell(cell);
            }
        }
    }
}

impl<T: Real> Default for TetMesh<T> {
    /// Produce an empty `TetMesh`. This is not particularly useful on its own, however it can be
    /// used as a null case for various mesh algorithms.
    fn default() -> Self {
        TetMesh::new(vec![], vec![])
    }
}

/**
 * Define `TetMesh` topology
 */

impl<T: Real> NumVertices for TetMesh<T> {
    #[inline]
    fn num_vertices(&self) -> usize {
        self.vertex_positions.len()
    }
}

impl<T: Real> NumCells for TetMesh<T> {
    #[inline]
    fn num_cells(&self) -> usize {
        self.indices.len()
    }
}

impl<T: Real> CellVertex for TetMesh<T> {
    #[inline]
    fn cell_to_vertex<CI>(&self, cidx: CI, which: usize) -> Option<VertexIndex>
    where
        CI: Copy + Into<CellIndex>,
    {
        if which >= 4 {
            None
        } else {
            Some(self.indices[cidx.into()][which].into())
        }
    }

    #[inline]
    fn cell_vertex<CI>(&self, cidx: CI, which: usize) -> Option<CellVertexIndex>
    where
        CI: Copy + Into<CellIndex>,
    {
        if which >= 4 {
            None
        } else {
            Some((4 * usize::from(cidx.into()) + which).into())
        }
    }

    #[inline]
    fn num_cell_vertices(&self) -> usize {
        self.indices.len() * 4
    }

    #[inline]
    fn num_vertices_at_cell<CI>(&self, _: CI) -> usize
    where
        CI: Copy + Into<CellIndex>,
    {
        4
    }
}

impl<T: Real> CellFace for TetMesh<T> {
    #[inline]
    fn cell_to_face<CI>(&self, cidx: CI, which: usize) -> Option<FaceIndex>
    where
        CI: Copy + Into<CellIndex>,
    {
        // Faces are indexed to be opposite to the corresponding vertices.
        if which >= 4 {
            None
        } else {
            Some(self.indices[cidx.into()][which].into())
        }
    }

    #[inline]
    fn cell_face<CI>(&self, cidx: CI, which: usize) -> Option<CellFaceIndex>
    where
        CI: Copy + Into<CellIndex>,
    {
        // Faces are indexed to be opposite to the corresponding vertices.
        if which >= 4 {
            None
        } else {
            Some((4 * usize::from(cidx.into()) + which).into())
        }
    }

    #[inline]
    fn num_cell_faces(&self) -> usize {
        self.indices.len() * 4
    }

    #[inline]
    fn num_faces_at_cell<CI>(&self, _: CI) -> usize
    where
        CI: Copy + Into<CellIndex>,
    {
        4
    }
}


#[cfg(test)]
mod tests {
    use super::*;
    use crate::index::Index;

    fn simple_tetmesh() -> TetMesh<f64> {
        let points = vec![
            [0.0, 0.0, 0.0],
            [1.0, 0.0, 0.0],
            [0.0, 1.0, 0.0],
            [1.0, 1.0, 0.0],
            [0.0, 0.0, 1.0],
            [1.0, 0.0, 1.0],
        ];
        let indices = vec![0, 4, 2, 5, 0, 5, 2, 3, 5, 3, 0, 1];

        TetMesh::new(points, indices)
    }

    #[test]
    fn tetmesh_test() {
        let mesh = simple_tetmesh();
        assert_eq!(mesh.num_vertices(), 6);
        assert_eq!(mesh.num_cells(), 3);
        assert_eq!(mesh.num_cell_vertices(), 12);
        assert_eq!(mesh.num_cell_faces(), 12);

        assert_eq!(Index::from(mesh.cell_vertex(1, 1)), 5);
        assert_eq!(Index::from(mesh.cell_vertex(0, 2)), 2);
        assert_eq!(Index::from(mesh.cell_face(2, 3)), 11);
    }

    #[test]
    fn tet_iter_test() {
        use math::Vector3;
        let mesh = simple_tetmesh();
        let points = mesh.vertex_positions().to_vec();

        let pt = |i| Vector3::from(points[i]);

        let tets = vec![
            Tetrahedron(pt(0), pt(4), pt(2), pt(5)),
            Tetrahedron(pt(0), pt(5), pt(2), pt(3)),
            Tetrahedron(pt(5), pt(3), pt(0), pt(1)),
        ];

        for (tet, exptet) in mesh.tet_iter().zip(tets.into_iter()) {
            assert_eq!(tet, exptet);
        }
    }

    /// Verify that inverting canonical tets causes their signed volume to become negative.
    #[test]
    fn invert_test() {
        use crate::ops::Volume;
        let mut mesh = simple_tetmesh();

        // Before inversion, canonical tets should have positive volume.
        let mut vols = Vec::new();
        for ref tet in mesh.tet_iter() {
            vols.push(tet.volume());
            assert!(tet.signed_volume() > 0.0);
        }

        mesh.invert();

        // After inversion, all tets should have negative volume.
        for (tet, vol) in mesh.tet_iter().zip(vols) {
            assert_eq!(tet.signed_volume(), -vol);
        }
    }

    /// Test that canonicalizing tets fixes all inverted tets but doesn't touch tets that are
    /// already in canonical form (not inverted).
    #[test]
    fn canonicalize_test() {
        use crate::ops::Volume;

        let points = vec![
            [0.0, 0.0, 0.0],
            [1.0, 0.0, 0.0],
            [0.0, 1.0, 0.0],
            [1.0, 1.0, 0.0],
            [0.0, 0.0, 1.0],
            [1.0, 0.0, 1.0],
        ];
        let indices = vec![0, 4, 2, 5, 0, 5, 3, 2, 5, 3, 1, 0];

        let mut mesh = TetMesh::new(points.clone(), indices);

        // Two tets are inverted
        let vols: Vec<_> = mesh.tet_iter().map(|t| t.volume()).collect();
        assert!(mesh.tet(0).signed_volume() > 0.0);
        assert!(mesh.tet(1).signed_volume() < 0.0);
        assert!(mesh.tet(2).signed_volume() < 0.0);

        mesh.canonicalize();

        // Canonicalization fixes up all inverted tets
        for (tet, vol) in mesh.tet_iter().zip(vols) {
            assert_eq!(tet.signed_volume(), vol);
        }
    }
}