Point.hpp

#ifndef POINT_HPP
#define POINT_HPP

#include <boost/fusion/sequence/intrinsic/at_c.hpp>
#include <boost/fusion/include/at_c.hpp>
#include <boost/fusion/container/vector.hpp>
#include <boost/fusion/include/vector.hpp>
#include <boost/fusion/container/vector/vector_fwd.hpp>
#include <boost/fusion/include/vector_fwd.hpp>
#include "boost/multi_array.hpp"
#include "Grid/Encap.hpp"
#include "Point_operators.hpp"


template<unsigned int dim ,typename T> class Point
{
    public:

    typedef int is_vtk_writable;

    typedef T coord_type;

    typedef boost::fusion::vector<T[dim]> type;

    type data;

    static const unsigned int x = 0;


    template<typename orig, typename exp1, typename exp2, unsigned int op> Point(const point_expression_op<orig,exp1,exp2,op> & p_exp)
    {
        this->operator=(p_exp);
    }

    inline Point(const Point<dim,T> && p)
    {
        for(size_t i = 0; i < dim ; i++)
        {get(i) = p.get(i);}
    }

    inline Point(const Point<dim,T> & p)
    {
        for(size_t i = 0; i < dim ; i++)
        {get(i) = p.get(i);}
    }

    inline Point(const T (&p)[dim])
    {
        for(size_t i = 0; i < dim ; i++)
        {get(i) = p[i];}
    }

    inline Point(T d)
    {
        this->operator=(d);
    }

    template <typename S> inline Point(const Point<dim,S> & p)
    {
        for (size_t i = 0 ; i < dim ; i++)
            get(i) = static_cast<S>(p.get(i));
    }

    template <unsigned int d, typename M> inline Point(const encapc<d,Point<dim,T>,M> & p)
    {
        for (size_t i = 0 ; i < dim ; i++)
            get(i) = p.template get<0>()[i];
    }

    inline Point(std::initializer_list<T> p1)
    {
        size_t i = 0;
        for(T x : p1)
        {get(i) = x;i++;}
    }

    inline Point()
    {}

    inline const T & get(size_t i) const
    {
        return boost::fusion::at_c<x>(data)[i];
    }

    inline T get_vtk(size_t i) const
    {
        return boost::fusion::at_c<x>(data)[i];
    }

    inline T& get(size_t i)
    {
        return boost::fusion::at_c<x>(data)[i];
    }

    inline T& operator[](size_t i)
    {
        return get(i);
    }

    inline const T& operator[](size_t i) const
    {
        return get(i);
    }

    T norm()
    {
        T n = 0.0;

        for (size_t i = 0 ; i < dim ; i++)
            n+=get(i) * get(i);

        return sqrt(n);
    }

    T distance(const Point<dim,T> & q)
    {
        T tot = 0.0;

        for (size_t i = 0 ; i < dim ; i++)
            tot += (this->get(i)  - q.get(i)) * (this->get(i)  - q.get(i));

        return sqrt(tot);
    }

    T distance2(const Point<dim,T> & q)
    {
        T tot = 0.0;

        for (size_t i = 0 ; i < dim ; i++)
            tot += (this->get(i)  - q.get(i)) * (this->get(i)  - q.get(i));

        return tot;
    }


    inline void zero()
    {
        for (size_t i = 0 ; i < dim ; i++)
        {
            get(i) = 0;
        }
    }

    inline void one()
    {
        for (size_t i = 0 ; i < dim ; i++)
        {
            get(i) = 1;
        }
    }

    inline static Point<dim,T> zero_p()
    {
        Point<dim,T> p;

        for (size_t i = 0 ; i < dim ; i++)
        {
            p.get(i) = 0;
        }

        return p;
    }

    std::string toPointString() const
    {
        std::stringstream ps;

        for (size_t i = 0 ; i < dim ; i++)
            ps << "x[" << i << "]=" << get(i) << " ";

        ps << "\n";

        return ps.str();
    }

    void swap(Point<dim,T> & p)
    {
        for (size_t i = 0 ; i < dim ; i++)
        {
            T tmp = get(i);
            get(i) = p.get(i);
            p.get(i) = tmp;
        }
    }

    inline bool operator==(const Point<dim,T> & p)
    {
        for (size_t i = 0 ; i < dim ; i++)
        {
            if (p.get(i) != get(i))
                return false;
        }

        return true;
    }

    inline bool operator!=(const Point<dim,T> & p)
    {
        return !this->operator==(p);

        return true;
    }

    std::string to_string() const
    {
        return toString();
    }

    std::string toString() const
    {
        std::string str;

        for (size_t i = 0 ; i < dim - 1 ; i++)
        {
            str += std::to_string(static_cast<double>(get(i))) + " ";
        }
        str += std::to_string(static_cast<double>(get(dim-1)));

        return str;
    }

    T & value(size_t i)
    {
        return get(i);
    }

    inline T value(size_t i) const
    {
        return get(i);
    }

    T (&asArray())[dim]
    {
        return boost::fusion::at_c<x>(data);
    }

    template<typename A> Point<dim,A> convertPoint() const
    {
        Point<dim,A> p;

        for (size_t i = 0; i < dim ; i++)
            p.get(i) = static_cast<A>(get(i));

        return p;
    }


    static bool noPointers()
    {
        return true;
    }


    template<typename orig, typename exp1, typename exp2, unsigned int op> Point<dim,T> & operator=(const point_expression_op<orig,exp1,exp2,op> & p_exp)
    {
        p_exp.init();

        for (size_t i = 0; i < dim ; i++)
            get(i) = p_exp.value(i);

        return *this;
    }

    Point<dim,T> & operator=(const point_expression<T[dim]> & p_exp)
    {
        p_exp.init();

        for (size_t i = 0; i < dim ; i++)
            get(i) = p_exp.value(i);

        return *this;
    }

    Point<dim,T> & operator=(const T (& p)[dim])
    {
        for (size_t i = 0; i < dim ; i++)
            get(i) = p[i];

        return *this;
    }

    template<typename T1, typename check = typename std::enable_if<std::is_const<T1>::value == false>::type> Point<dim,T> & operator=(const point_expression<const T1[dim]> & p_exp)
    {
        p_exp.init();

        for (size_t i = 0; i < dim ; i++)
            get(i) = p_exp.value(i);

        return *this;
    }

    template<typename aT> inline Point<dim,T> operator/(const aT (&ar)[dim])
    {
        Point<dim,T> result;

        for (size_t i = 0 ; i < dim ; i++)
            result.get(i) = get(i) / ar[i];

        return result;
    }

    template<typename aT> inline Point<dim,T> operator/=(const aT c)
    {
        for (size_t i = 0 ; i < dim ; i++)
            get(i) = get(i) / c;

        return *this;
    }


    Point<dim,T> & operator=(T d)
    {
        for (size_t i = 0; i < dim ; i++)
            get(i) = d;

        return *this;
    }

    inline Point<dim,T> & operator=(const Point<dim,T> & p)
    {
        for (size_t i = 0 ; i < dim ; i++)
            get(i) = p.get(i);

        return *this;
    }

    inline Point<dim,T> & operator-=(const Point<dim,T> & p)
    {
        for (size_t i = 0 ; i < dim ; i++)
            get(i) -= p.get(i);

        return *this;
    }

    inline Point<dim,T> & operator+=(const Point<dim,T> & p)
    {
        for (size_t i = 0 ; i < dim ; i++)
            get(i) += p.get(i);

        return *this;
    }

    inline void init() const
    {}


    static const unsigned int max_prop = 1;
    static const unsigned int max_prop_real = 1;

    static const unsigned int dims = dim;

    static const unsigned int nvals = dim;
};

template <unsigned int N, typename T> std::string toPointString(const T (&p)[N] )
{
    std::stringstream ps;

    for (size_t i = 0 ; i < N ; i++)
        ps << "x[" << i << "]=" << p[i] << " ";

    ps << "\n";

    return ps.str();
}

template <unsigned int N, typename T, typename Mem> std::string toPointString(const encapc<1,Point<N,T>,Mem> & p )
{
    std::stringstream ps;

    for (size_t i = 0 ; i < N ; i++)
        ps << "x[" << i << "]=" << p.template get<Point<N,T>::x>()[i] << " ";

    ps << "\n";

    return ps.str();
}



template<unsigned int dim, typename T>  using VectorS = Point<dim,T>;

#endif