Vector_eigen.hpp

/*
 * Vector_eigen.hpp
 *
 *  Created on: Nov 27, 2015
 *      Author: i-bird
 */

#ifndef OPENFPM_NUMERICS_SRC_VECTOR_VECTOR_EIGEN_HPP_
#define OPENFPM_NUMERICS_SRC_VECTOR_VECTOR_EIGEN_HPP_

#include <type_traits>
#include "util/mul_array_extents.hpp"
#include <fstream>
#include "Grid/staggered_dist_grid.hpp"
#include "Space/Ghost.hpp"
#include "FiniteDifference/util/common.hpp"
#include <boost/mpl/vector_c.hpp>
#include <unordered_map>
#include "Vector_util.hpp"

#define EIGEN_RVAL 1

template<typename T>
class rval<T,EIGEN_RVAL>
{
    long int r;

    T val;

public:

    long int & row()
    {
        return r;
    }

    T & value()
    {
        return val;
    }

    rval()
    :r(0)
    {}

    rval(long int i, T val)
    {
        row() = i;
        value() = val;
    }

    static inline bool noPointers()
    {
        return true;
    }
};

template<typename T>
class Vector<T, EIGEN_BASE>
{
    mutable Eigen::Matrix<T, Eigen::Dynamic, 1> v;

    mutable openfpm::vector<rval<T,EIGEN_RVAL>> row_val;

    mutable openfpm::vector<rval<T,EIGEN_RVAL>> row_val_recv;

    mutable std::unordered_map<size_t,size_t> map;

    T invalid;

    mutable openfpm::vector<size_t> prc;

    //size of each chunk
    mutable openfpm::vector<size_t> sz;


    void collect() const
    {
        Vcluster & vcl = create_vcluster();

        row_val_recv.clear();

        // here we collect all the triplet in one array on the root node
        vcl.SGather(row_val,row_val_recv,prc,sz,0);

        if (vcl.getProcessUnitID() != 0)
            row_val.resize(0);
        else
            row_val.swap(row_val_recv);

        build_map();
    }

    void setEigen() const
    {
        // set the vector

        for (size_t i = 0 ; i < row_val.size() ; i++)
            v[row_val.get(i).row()] = row_val.get(i).value();
    }

    void build_map() const
    {
        map.clear();

        for (size_t i = 0 ; i < row_val.size() ; i++)
            map[row_val.get(i).row()] = i;
    }

public:

    Vector(const Vector<T> & v)
    :invalid(0)
    {
        this->operator=(v);
    }

    Vector(const Vector<T> && v)
    :invalid(0)
    {
        this->operator=(v);
    }

    Vector(size_t n)
    {
        resize(n,0);
    }

    Vector()
    {
    }

    void resize(size_t row, size_t l_row)
    {
        v.resize(row);
    }

    void insert(size_t i, T val)
    {
        row_val.add();

        // Map
        map[i] = row_val.size()-1;

        row_val.last().row() = i;
        row_val.last().value() = val;
    }

    inline T & insert(size_t i)
    {
        row_val.add();

        // Map
        map[i] = row_val.size()-1;

        row_val.last().row() = i;
        return row_val.last().value();
    }

    inline const T & insert(size_t i) const
    {
        row_val.add();

        // Map
        map[i] = row_val.size()-1;

        row_val.last().row() = i;
        return row_val.last().value();
    }

    const T & operator()(size_t i) const
    {
        // Search if exist

        std::unordered_map<size_t,size_t>::iterator it = map.find(i);

        if ( it != map.end() )
            return row_val.get(it->second).value();

        return insert(i);
    }

    T & operator()(size_t i)
    {
        // Search if exist

        std::unordered_map<size_t,size_t>::iterator it = map.find(i);

        if ( it != map.end() )
            return row_val.get(it->second).value();

        return insert(i);
    }

    const Eigen::Matrix<T, Eigen::Dynamic, 1> & getVec() const
    {
        collect();
        setEigen();

        return v;
    }

    Eigen::Matrix<T, Eigen::Dynamic, 1> & getVec()
    {
        collect();
        setEigen();

        return v;
    }

    void scatter()
    {
        row_val_recv.clear();
        Vcluster & vcl = create_vcluster();

        vcl.SScatter(row_val,row_val_recv,prc,sz,0);

        // if we do not receive anything a previous collect has not been performed
        // and so nothing is scattered
        if (row_val_recv.size() != 0)
        {
            row_val.clear();
            row_val.add(row_val_recv);
            build_map();
        }
    }

    void fromFile(std::string file)
    {
        std::ifstream inputf;
        inputf.open(file);

        for (size_t i = 0 ; i < v.size() ; i++)
            inputf >> v(i);

        inputf.close();

    }

    Vector<T> & operator=(const Vector<T> & v)
    {
        prc = v.prc;
        sz = v.sz;
        map = v.map;
        row_val = v.row_val;

        return *this;
    }

    Vector<T> & operator=(const Vector<T> && v)
    {
        prc = v.prc;
        sz = v.sz;
        map = v.map;
        row_val = v.row_val;

        return *this;
    }

    Vector<T> & operator=(Eigen::Matrix<T, Eigen::Dynamic, 1> & v)
    {
        for (size_t i = 0 ; i < row_val.size() ; i++)
            row_val.get(i).value() = v(row_val.get(i).row());

        return *this;
    }
};


#endif /* OPENFPM_NUMERICS_SRC_VECTOR_VECTOR_EIGEN_HPP_ */