CellListMem.hpp

/*
 * CelListMem.hpp
 *
 *  Created on: Mar 22, 2015
 *  Last modified: June 25, 2015
 *      Authors: Pietro Incardona, Yaroslav Zaluzhnyi
 */

#ifndef CELLISTMEM_HPP_
#define CELLISTMEM_HPP_

#include <unordered_map>

template<unsigned int dim, typename T, typename transform, typename base>
class CellList<dim,T,MEMORY,transform,base> : public CellDecomposer_sm<dim,T,transform>
{
    // The array contain the neighborhood of the cell-id in case of asymmetric interaction
    //
    //    * * *
    //    * x *
    //    * * *

    long int NNc_full[openfpm::math::pow(3,dim)];

    // The array contain the neighborhood of the cell-id in case of symmetric interaction
    //
    //   * * *
    //     x *
    //
    long int NNc_sym[openfpm::math::pow(3,dim)/2+1];

    // The array contain the neighborhood of the cell-id in case of symmetric interaction (Optimized)
    //
    //   * *
    //   x *
    //
    long int NNc_cr[openfpm::math::pow(2,dim)];

    // each cell has a dynamic structure
    // that store the elements in the cell
    std::unordered_map<size_t,base> cl_base;

    //Origin point
    Point<dim,T> orig;

public:

    // Object type that the structure store
    typedef T value_type;

    grid_sm<dim,void> & getGrid()
    {
        CellDecomposer_sm<dim,T>::getGrid();
    }

    void Initialize(Box<dim,T> & box, size_t (&div)[dim], Point<dim,T> & orig, const size_t pad = 1)
    {
        SpaceBox<dim,T> sbox;
        Initialize(sbox,div,orig);
    }

    void Initialize(SpaceBox<dim,T> & box, size_t (&div)[dim], Point<dim,T> & orig, const size_t pad = 1)
    {
        // Add padding
        size_t div_pad[dim];
        for (size_t i = 0 ; i < dim ; i++)
            div_pad[i] = div[i] + 2;

        CellDecomposer_sm<dim,T>::setDimensions(box,div_pad, pad);

        this->orig = orig;

        //filling a map with "base" structures
        for (int i = 0; i < this->tot_n_cell; i++)
            cl_base[i] = base();


        // Calculate the NNc-arrays (for neighborhood):

        // compile-time array {0,0,0,....} and {3,3,3,...}

        typedef typename generate_array<size_t,dim, Fill_zero>::result NNzero;
        typedef typename generate_array<size_t,dim, Fill_two>::result NNtwo;
        typedef typename generate_array<size_t,dim, Fill_one>::result NNone;

        // Generate the sub-grid iterator

        grid_key_dx_iterator_sub<dim> gr_sub3(this->gr_cell,NNzero::data,NNtwo::data);

        // Calculate the NNc array

        size_t middle = this->gr_cell.LinId(NNone::data);
        size_t i = 0;
        while (gr_sub3.isNext())
        {
            NNc_full[i] = (long int)this->gr_cell.LinId(gr_sub3.get()) - middle;

            ++gr_sub3;
            i++;
        }

        // Calculate the NNc_sym array

        i = 0;
        gr_sub3.reset();
        while (gr_sub3.isNext())
        {
            auto key = gr_sub3.get();

            size_t lin = this->gr_cell.LinId(key);

            // Only the first half is considered
            if (lin < middle)
            {
                ++gr_sub3;
                continue;
            }

            NNc_sym[i] = lin - middle;

            ++gr_sub3;
            i++;
        }

        // Calculate the NNc_cross array

        i = 0;
        grid_key_dx_iterator_sub<dim> gr_sub2(this->gr_cell,NNzero::data,NNone::data);

        while (gr_sub2.isNext())
        {
            auto key = gr_sub2.get();

            NNc_cr[i] = (long int)this->gr_cell.LinId(key);

            ++gr_sub2;
            i++;
        }
    }

    CellList()
    {
    }


    CellList(Box<dim,T> & box, size_t (&div)[dim], Point<dim,T> & orig, const size_t pad = 1)
    {
        SpaceBox<dim,T> sbox(box);
        Initialize(sbox,div,orig,pad);
    }

    CellList(SpaceBox<dim,T> & box, size_t (&div)[dim], Point<dim,T> & orig, const size_t pad = 1)
    {
        Initialize(box,div,orig,pad);
    }


    void add(const T (& pos)[dim], typename base::value_type ele)
    {
        // calculate the Cell id

        size_t cell_id = this->getCell(pos);

        // Get the number of element the cell is storing

        cl_base[cell_id].add(ele);
    }

    void add(const Point<dim,T> & pos, typename base::value_type ele)
    {
        // calculate the Cell id

        size_t cell_id = this->getCell(pos);

        // add a new element

        cl_base[cell_id].add(ele);
    }


    void remove(size_t cell, size_t ele)
    {
        cl_base[cell].remove(ele);
    }

    size_t getNelements(size_t cell_id)
    {
        return cl_base[cell_id].size();
    }

    auto get(size_t cell, size_t ele) -> decltype(cl_base[cell].get(ele))
    {
        return cl_base[cell].get(ele);
    }
    void swap(CellList<dim,T,MEMORY,base> & cl)
    {
        cl_base.swap(cl.cl_base);
    }

    CellIterator<CellList<dim,T,MEMORY,base>> getIterator(size_t cell)
    {
        return CellIterator<CellList<dim,T,MEMORY,base>>(cell,*this);
    }

    template<unsigned int impl> CellNNIterator<dim,CellList<dim,T,MEMORY,base>,FULL,impl> getNNIterator(size_t cell)
    {
        CellNNIterator<dim,CellList<dim,T,MEMORY,base>,FULL,impl> cln(cell,NNc_full,*this);

        return cln;
    }

    template<unsigned int impl> CellNNIterator<dim,CellList<dim,T,MEMORY,base>,SYM,impl> getNNIteratorSym(size_t cell)
    {
        CellNNIterator<dim,CellList<dim,T,MEMORY,base>,SYM,impl> cln(cell,NNc_sym,*this);

        return cln;
    }

    template<unsigned int impl> CellNNIterator<dim,CellList<dim,T,MEMORY,base>,CRS,impl> getNNIteratorCross(size_t cell)
    {
        CellNNIterator<dim,CellList<dim,T,MEMORY,base>,CRS,impl> cln(cell,NNc_cr,*this);

        return cln;
    }
};


#endif /* CELLISTMEM_HPP_ */