tetrahedron.hh

#ifndef VTESTBED_GEOMETRY_TETRAHEDRON_HH
#define VTESTBED_GEOMETRY_TETRAHEDRON_HH

#include <array>
#include <iosfwd>

#include <vtestbed/geometry/math.hh>
#include <vtestbed/geometry/triangle.hh>
#include <vtestbed/geometry/types.hh>

namespace vtb {

    namespace geometry {

        template <class T, int N>
        class Tetrahedron: public std::array<Vertex<T,N>,4> {

        private:
            using base_type = std::array<Vertex<T,N>,4>;

        public:
            static constexpr const int dimensions = N;
            using scalar_type = T;
            using value_type = typename base_type::value_type;
            using pointer = typename base_type::pointer;
            using const_pointer = typename base_type::const_pointer;
            using reference = typename base_type::reference;
            using const_reference = typename base_type::const_reference;
            using iterator = typename base_type::iterator;
            using const_iterator = typename base_type::const_iterator;
            using size_type = typename base_type::size_type;
            using difference_type = typename base_type::difference_type;
            using reverse_iterator = typename base_type::reverse_iterator;
            using const_reverse_iterator = typename base_type::const_reverse_iterator;

        public:
            using base_type::base_type;
            using base_type::operator=;
            using base_type::operator[];

        public:

            inline explicit
            Tetrahedron(
                const_reference p0,
                const_reference p1,
                const_reference p2,
                const_reference p3
            ): base_type{p0,p1,p2,p3} {}

            inline explicit
            Tetrahedron(const_reference origin, const Triangle<T,N>& base):
            base_type{origin,base[0],base[1],base[2]} {}

        };

        template <class T, int N>
        struct is_figure<Tetrahedron<T,N>>: public std::true_type {};

        template <class T, int N>
        inline Vertex<T,N>
        centroid(const Tetrahedron<T,N>& t) {
            return (t[0]+t[1]+t[2]+t[3])*T{0.25};
        }

        template <class T, int N>
        inline T
        signed_volume(const Tetrahedron<T,N>& t) {
            Vertex<T,N> e1 = t[1]-t[0];
            Vertex<T,N> e2 = t[2]-t[0];
            Vertex<T,N> e3 = t[3]-t[0];
            return dot(e3, cross(e1,e2)) / T{6};
        }

        template <class T, int N>
        inline T
        volume(const Tetrahedron<T,N>& t) {
            return std::abs(signed_volume(t));
        }

        template <class T>
        inline T
        solid_angle(const Tetrahedron<T,3>& t) {
            using std::atan2;
            Vertex<T,3> A = t[1]-t[0];
            Vertex<T,3> B = t[2]-t[0];
            Vertex<T,3> C = t[3]-t[0];
            auto a = length(A);
            auto b = length(B);
            auto c = length(C);
            return 2*atan2(det(A,B,C), a*b*c + c*dot(A,B) + a*dot(B,C) + b*dot(C,A));
        }

    }

}

#endif // vim:filetype=cpp