-
Pietro Incardona authoredPietro Incardona authored
grid_dist_amr.hpp 27.13 KiB
/*
* grid_amr_dist.hpp
*
* Created on: Sep 21, 2017
* Author: i-bird
*/
#ifndef AMR_GRID_AMR_DIST_HPP_
#define AMR_GRID_AMR_DIST_HPP_
#define OPENFPM_DATA_ENABLE_IO_MODULE
#include "Grid/grid_dist_id.hpp"
#include "Amr/grid_dist_amr_key_iterator.hpp"
#ifdef __NVCC__
#include "SparseGridGpu/SparseGridGpu.hpp"
#endif
#define AMR_IMPL_TRIVIAL 1
#define AMR_IMPL_PATCHES 2
#define AMR_IMPL_OPENVDB 3
template<typename Decomposition, typename garray>
class Decomposition_encap
{
Decomposition & dec;
garray & gd_array;
public:
Decomposition_encap(Decomposition & dec, garray & gd_array)
:dec(dec),gd_array(gd_array)
{}
Decomposition & internal_dec() const
{
return dec;
}
/*! \brief Start decomposition
*
*/
void decompose()
{
dec.decompose();
for(size_t i = 0 ; i < gd_array.size() ; i++)
{
Ghost<Decomposition::dims,typename Decomposition::stype> gold = gd_array.get(i).getDecomposition().getGhost();
gd_array.get(i).getDecomposition() = dec.duplicate(gold);
}
}
/*! \brief Refine the decomposition, available only for ParMetis distribution, for Metis it is a null call
*
* \param ts number of time step from the previous load balancing
*
*/
void refine(size_t ts)
{
dec.refine();
for(size_t i = 0 ; i < gd_array.size() ; i++)
{
Ghost<Decomposition::dims,typename Decomposition::stype> gold = gd_array.get(i).getDecomposition().getGhost();
gd_array.get(i).getDecomposition() = dec.duplicate(gold);
}
}
/*! \brief Refine the decomposition, available only for ParMetis distribution, for Metis it is a null call
*
* \param ts number of time step from the previous load balancing
*
*/
void redecompose(size_t ts)
{
dec.redecompose();
for(size_t i = 0 ; i < gd_array.size() ; i++)
{
Ghost<Decomposition::dims,typename Decomposition::stype> gold = gd_array.get(i).getDecomposition().getGhost();
gd_array.get(i).getDecomposition() = dec.duplicate(gold);
}
}
auto getDistribution() -> decltype(dec.getDistribution())
{
return dec.getDistribution();
}
Decomposition_encap<Decomposition,garray> operator=(const Decomposition_encap<Decomposition,garray> & de) const
{
for(size_t i = 0 ; i < gd_array.size() ; i++)
{gd_array.get(i).getDecomposition() = de.gd_array.get(i).getDecomposition();}
return *this;
}
bool write(std::string output) const
{
return dec.write(output);
}
};
template<unsigned int dim,
typename St,
typename T,
unsigned int impl=AMR_IMPL_TRIVIAL ,
typename Decomposition = CartDecomposition<dim,St>,
typename Memory=HeapMemory ,
typename device_grid=grid_cpu<dim,T> >
class grid_dist_amr
{
};
/*! \brief AMR Adaptive Multi Resolution Grid
*
* \tparam dim Dimensionality
* \tparam St type of space
* \tparam T what each point of the grid store
* \tparam Decomposition type of decomposition
*
*/
template<unsigned int dim,
typename St,
typename T,
typename Decomposition,
typename Memory,
typename device_grid >
class grid_dist_amr<dim,St,T,AMR_IMPL_TRIVIAL,Decomposition,Memory,device_grid>
{
//! Simulation domain
Box<dim,St> domain;
//! Ghost integer
Ghost<dim,long int> g_int;
//! Boundary conditions of the structure
periodicity<dim> bc;
//! array of grids
//
openfpm::vector<grid_dist_id<dim,St,T,Decomposition,Memory,device_grid>,
HeapMemory,
memory_traits_lin,
openfpm::grow_policy_identity,STD_VECTOR> gd_array;
//! Type of structure sub-grid iterator
typedef decltype(device_grid::type_of_subiterator()) device_sub_it;
//! Type of structure for the grid iterator
typedef decltype(device_grid::type_of_iterator()) device_it;
//! Domain iterator for each distributed grid
openfpm::vector<grid_dist_iterator<dim,device_grid,device_sub_it,FREE>> git;
//! Domain and ghost iterator for each distributed grid
openfpm::vector<grid_dist_iterator<dim,device_grid,device_it,FIXED>> git_g;
//! Iterator for each distributed grid
openfpm::vector<grid_dist_iterator_sub<dim,device_grid>> git_sub;
//! Moving offsets
openfpm::vector<openfpm::vector<offset_mv<dim>>> mv_off;
//! background level
T bck;
/*! \brief Initialize the others levels
*
* \param n_grid_dist_id<dim,St,T,Decomposition,Memory,device_grid>lvl number of levels
* \param g_sz_lvl grid size on each level
*
*/
void initialize_other(size_t n_lvl, size_t (& g_sz_lvl)[dim])
{
for (size_t i = 0; i < n_lvl - 1 ; i++)
{
for (size_t j = 0 ; j < dim ; j++)
{
if (bc.bc[j] == NON_PERIODIC)
{g_sz_lvl[j] = (g_sz_lvl[j]-1)*2 + 1;}
else
{g_sz_lvl[j] = g_sz_lvl[j]*2;}
}
gd_array.add(grid_dist_id<dim,St,T,Decomposition,Memory,device_grid>(gd_array.get(0).getDecomposition(),g_sz_lvl,g_int));
gd_array.last().setBackgroundValue(bck);
gd_array.last().getDecomposition().free_geo_cell();
gd_array.last().getDecomposition().getDistribution().destroy_internal_graph();
gd_array.last().getDecomposition().free_fines();
}
recalculate_mvoff();
}
public:
/*! \brief Constructor
*
* \param domain Simulation domain
* \param g ghost extension
*
*/
grid_dist_amr(const Box<dim,St> & domain, const Ghost<dim,long int> & g) :domain(domain),g_int(g)
{
// set boundary consitions to non periodic
for (size_t i = 0; i < dim ; i++)
{bc.bc[i] = NON_PERIODIC;}
}
/*! \brief Constructor
*
* \param domain Simulation domain
* \param g ghost extension
* \param bc boundary conditions
*
*/
grid_dist_amr(const Box<dim,St> & domain, const Ghost<dim,long int> & g, periodicity<dim> & bc)
:domain(domain),g_int(g),bc(bc)
{
}
/*! \brief Initialize the amr grid
*
* \param dec Decomposition (this parameter is useful in case we want to constrain the AMR to an external decomposition)
* \param n_lvl maximum number of levels (0 mean no additional levels)
* \param g_sz coarsest grid size on each direction
*
*/
void initLevels(const Decomposition & dec, size_t n_lvl,const size_t (& g_sz)[dim])
{
size_t g_sz_lvl[dim];
for (size_t i = 0; i < dim ; i++)
{g_sz_lvl[i] = g_sz[i];}
// Add the coarse level
gd_array.add(grid_dist_id<dim,St,T,Decomposition,Memory,device_grid>(dec,g_sz,g_int));
gd_array.last().setBackgroundValue(bck);
initialize_other(n_lvl,g_sz_lvl);
}
/*! \brief Initialize the amr grid
*
* \param dec Decomposition (this parameter is useful in case we want to constrain the AMR to an external decomposition)
* \param n_lvl maximum number of levels (0 mean no additional levels)
* \param g_sz coarsest grid size on each direction
*
*/
template<typename TT> void initLevels(const Decomposition_encap<Decomposition,TT> & dec, size_t n_lvl,const size_t (& g_sz)[dim])
{
initLevels(dec.internal_dec(),n_lvl,g_sz);
}
/*! \brief Recalculate the offset array for the moveLvlUp and moveLvlDw
*
*
*
*/
void recalculate_mvoff()
{
// Here we calculate the offset to move one level up and one level down
// in global coordinated moving one level up is multiply the coordinates by 2
// and moving one level down is dividing by 2. In local coordinates is the same
// with the exception that because of the decomposition you can have an offset
// look at the picture below
//
// (-1) (0)
// * | * * coarse level
// * |* * * * finer level
// |(0)(1) //
// Line of the decomposition
//
// The coarse level point 0 in local coordinates converted to the finer level is not
// just 2*0 = 0 but is 2*(0) + 1 so a formula like 2*x+offset is required. here we calculate
// these offset. In the case of moving from finer to coarse is the same the formula is
// Integer_round(x+1)/2 - 1
//
mv_off.resize(gd_array.size());
for (size_t i = 1 ; i < gd_array.size() ; i++)
{
auto & g_box_c = gd_array.get(i-1).getLocalGridsInfo();
auto & g_box_f = gd_array.get(i).getLocalGridsInfo();
#ifdef SE_CLASS1
if (g_box_c.size() != g_box_f.size())
{
std::cerr << __FILE__ << ":" << __LINE__ << " error it seem that the AMR construction between level " <<
i << " and " << i-1 << " is inconsistent" << std::endl;
}
#endif
mv_off.get(i-1).resize(g_box_f.size());
mv_off.get(i).resize(g_box_f.size());
for (size_t j = 0 ; j < g_box_f.size() ; j++)
{
for (size_t s = 0 ; s < dim ; s++)
{
size_t d_orig_c = g_box_c.get(j).origin.get(s);
size_t d_orig_f = g_box_f.get(j).origin.get(s);
mv_off.get(i-1).get(j).dw.get(s) = d_orig_c*2 - d_orig_f;
mv_off.get(i).get(j).up.get(s) = d_orig_c*2 - d_orig_f;
}
}
}
}
/*! \brief Initialize the amr grid
*
* \param n_lvl maximum number of levels (0 mean no additional levels)
* \param g_sz coarsest grid size on each direction
* \param opt options
*
*/
void initLevels(size_t n_lvl,const size_t (& g_sz)[dim], size_t opt = 0)
{
size_t g_sz_lvl[dim];
for (size_t i = 0; i < dim ; i++)
{g_sz_lvl[i] = g_sz[i];}
// Add the coarse level
gd_array.add(grid_dist_id<dim,St,T,Decomposition,Memory,device_grid>(g_sz,domain,g_int,bc,opt));
initialize_other(n_lvl,g_sz_lvl);
}
/*! \brief Add the computation cost on the decomposition using a resolution function
*
*
* \param md Model to use
* \param ts It is an optional parameter approximately should be the number of ghost get between two
* rebalancing at first decomposition this number can be ignored (default = 1) because not used
*
*/ template <typename Model>inline void addComputationCosts(Model md=Model(), size_t ts = 1)
{
gd_array.get(0).addComputationCosts(md,ts);
}
/*! \brief Get the object that store the information about the decomposition
*
* \return the decomposition object
*
*/
Decomposition_encap<Decomposition,decltype(gd_array)> getDecomposition()
{
Decomposition_encap<Decomposition,decltype(gd_array)> tmp(gd_array.get(0).getDecomposition(),gd_array);
return tmp;
}
/*! \brief Get the underlying grid level
*
* \param lvl level
*
* \return the grid level
*
*/
grid_dist_id<dim,St,T,Decomposition,Memory,device_grid> & getLevel(size_t lvl)
{
return gd_array.get(lvl);
}
grid_dist_amr_key_iterator<dim,device_grid,
decltype(grid_dist_id<dim,St,T,Decomposition,Memory,device_grid>::type_of_subiterator()),
decltype(grid_dist_id<dim,St,T,Decomposition,Memory,device_grid>::type_of_subiterator()) >
getDomainIteratorCells()
{
git_sub.clear();
for (size_t i = 0 ; i < gd_array.size() ; i++)
{
grid_key_dx<dim> start;
grid_key_dx<dim> stop;
for (size_t j = 0 ; j < dim ; j++)
{
start.set_d(j,0);
if (bc.bc[j] == NON_PERIODIC)
{stop.set_d(j,getGridInfoVoid(i).size(j) - 2);}
else
{stop.set_d(j,getGridInfoVoid(i).size(j) - 1);}
}
git_sub.add(gd_array.get(i).getSubDomainIterator(start,stop));
}
return grid_dist_amr_key_iterator<dim,device_grid,
decltype(grid_dist_id<dim,St,T,Decomposition,Memory,device_grid>::type_of_subiterator()),
decltype(grid_dist_id<dim,St,T,Decomposition,Memory,device_grid>::type_of_subiterator())>(git_sub);
}
grid_dist_iterator_sub<dim,device_grid> getDomainIteratorCells(size_t lvl)
{
grid_key_dx<dim> start;
grid_key_dx<dim> stop;
for (size_t j = 0 ; j < dim ; j++)
{
start.set_d(j,0);
if (bc.bc[j] == NON_PERIODIC)
{stop.set_d(j,getGridInfoVoid(lvl).size(j) - 2);}
else {stop.set_d(j,getGridInfoVoid(lvl).size(j) - 1);}
}
return gd_array.get(lvl).getSubDomainIterator(start,stop);
}
/*! \brief Get an iterator to the grid
*
* \return an iterator to the grid
*
*/
auto getGridGhostIterator(size_t lvl) -> decltype(gd_array.get(lvl).getGridGhostIterator(grid_key_dx<dim>(),grid_key_dx<dim>()))
{
grid_key_dx<dim> key_start;
grid_key_dx<dim> key_stop;
for (size_t i = 0 ; i < dim ; i++)
{
key_start.set_d(i,g_int.getLow(i));
key_stop.set_d(i,g_int.getHigh(i) + getGridInfoVoid(lvl).size(i) -1);
}
return gd_array.get(lvl).getGridGhostIterator(key_start,key_stop);
}
/*! \brief Get an iterator to the grid
*
* \return an iterator to the grid
*
*/
auto getGridIterator(size_t lvl) -> decltype(gd_array.get(lvl).getGridIterator())
{
return gd_array.get(lvl).getGridIterator();
}
/*! \brief Get an iterator to the grid
*
* \return an iterator to the grid
*
*/
auto getGridIterator(size_t lvl, grid_key_dx<dim> & start, grid_key_dx<dim> & stop) -> decltype(gd_array.get(lvl).getGridIterator(start,stop))
{
return gd_array.get(lvl).getGridIterator(start,stop);
}
#ifdef __NVCC__
/*! \brief Get an iterator to the grid
*
* \return an iterator to the grid
*
*/
auto getGridIteratorGPU(size_t lvl) -> decltype(gd_array.get(lvl).getGridIteratorGPU())
{
return gd_array.get(lvl).getGridIteratorGPU();
}
#endif
/*! \brief Get an iterator to the grid
*
* \return an iterator to the grid
*
*/
auto getGridIteratorCells(size_t lvl) -> decltype(gd_array.get(lvl).getGridIterator())
{
grid_key_dx<dim> start;
grid_key_dx<dim> stop;
for (size_t j = 0 ; j < dim ; j++) {
start.set_d(j,0);
if (bc.bc[j] == NON_PERIODIC)
{stop.set_d(j,getGridInfoVoid(lvl).size(j) - 2);}
else
{stop.set_d(j,getGridInfoVoid(lvl).size(j) - 1);}
}
return gd_array.get(lvl).getGridIterator(start,stop);
}
/*! \brief return an iterator over the level lvl
*
* \param lvl level
*
* \return an iterator over the level lvl selected
*
*/
grid_dist_iterator<dim,device_grid,decltype(device_grid::type_of_subiterator()),FREE>
getDomainIterator(size_t lvl) const
{
return gd_array.get(lvl).getDomainIterator();
}
/*! \brief return an iterator over the level lvl
*
* \param lvl level
*
* \return an iterator over the level lvl selected
*
*/
grid_dist_iterator<dim,device_grid,
decltype(device_grid::type_of_iterator()),
FIXED>
getDomainGhostIterator(size_t lvl) const
{
return gd_array.get(lvl).getDomainGhostIterator();
}
/*! \brief Get domain iterator
*
* \return an iterator over all the grid levels
*
*/
grid_dist_amr_key_iterator<dim,device_grid, decltype(device_grid::type_of_subiterator())>
getDomainIterator()
{
git.clear();
for (size_t i = 0 ; i < gd_array.size() ; i++)
{
git.add(gd_array.get(i).getDomainIterator());
}
return grid_dist_amr_key_iterator<dim,device_grid,decltype(device_grid::type_of_subiterator())>(git);
}
/*! \brief Get domain iterator
*
* \return an iterator over all the grid levels
*
*/
grid_dist_amr_key_iterator<dim,device_grid, decltype(device_grid::type_of_iterator()),
grid_dist_iterator<dim,device_grid,decltype(device_grid::type_of_iterator()),FIXED>>
getDomainGhostIterator()
{
git_g.clear();
for (size_t i = 0 ; i < gd_array.size() ; i++) {
git_g.add(gd_array.get(i).getDomainGhostIterator());
}
return grid_dist_amr_key_iterator<dim,device_grid,decltype(device_grid::type_of_iterator()),
grid_dist_iterator<dim,device_grid,decltype(device_grid::type_of_iterator()),FIXED>>(git_g);
}
/*! \brief Get the reference of the selected element
*
* \tparam p property to get (is an integer)
* \param v1 grid_key that identify the element in the grid
*
* \return the selected element
*
*/
template <unsigned int p>inline auto get(const grid_dist_amr_key<dim> & v1) const -> decltype(gd_array.get(v1.getLvl()).template get<p>(v1.getKey()))
{
#ifdef SE_CLASS2
check_valid(this,8);
#endif
return gd_array.get(v1.getLvl()).template get<p>(v1.getKey());
}
/*! \brief Get the reference of the selected element
*
* \tparam p property to get (is an integer)
* \param v1 grid_key that identify the element in the grid
*
* \return the selected element
*
*/
template <unsigned int p>inline auto get(const grid_dist_amr_key<dim> & v1) -> decltype(gd_array.get(v1.getLvl()).template get<p>(v1.getKey()))
{
#ifdef SE_CLASS2
check_valid(this,8);
#endif
return gd_array.get(v1.getLvl()).template get<p>(v1.getKey());
}
/*! \brief Get the reference of the selected element
*
* \tparam p property to get (is an integer)
* \param v1 grid_key that identify the element in the grid
*
* \return the selected element
*
*/
template <unsigned int p>inline auto get(size_t lvl, const grid_dist_key_dx<dim> & v1) const -> decltype(gd_array.get(lvl).template get<p>(v1))
{
#ifdef SE_CLASS2
check_valid(this,8);
#endif
return gd_array.get(lvl).template get<p>(v1);
}
/*! \brief Get the reference of the selected element
*
* \tparam p property to get (is an integer)
* \param v1 grid_key that identify the element in the grid
*
* \return the selected element
*
*/
template <unsigned int p>inline auto get(size_t lvl, const grid_dist_key_dx<dim> & v1) -> decltype(gd_array.get(lvl).template get<p>(v1))
{
#ifdef SE_CLASS2
check_valid(this,8);
#endif return gd_array.get(lvl).template get<p>(v1);
}
//////////////////// Insert functions
/*! \brief Get the reference of the selected element
*
* \tparam p property to get (is an integer)
* \param v1 grid_key that identify the element in the grid
*
* \return the selected element
*
*/
template <unsigned int p>
inline auto insert(const grid_dist_amr_key<dim> & v1)
-> decltype(gd_array.get(v1.getLvl()).template insert<p>(v1.getKey()))
{
#ifdef SE_CLASS2
check_valid(this,8);
#endif
return gd_array.get(v1.getLvl()).template insert<p>(v1.getKey());
}
/*! \brief Get the reference of the selected element
*
* \tparam p property to get (is an integer)
* \param v1 grid_key that identify the element in the grid
*
* \return the selected element
*
*/
template <unsigned int p>inline auto insert(size_t lvl, const grid_dist_key_dx<dim> & v1)
-> decltype(gd_array.get(lvl).template insert<p>(v1))
{
#ifdef SE_CLASS2
check_valid(this,8);
#endif
return gd_array.get(lvl).template insert<p>(v1);
}
//////////////////////////////////////
/*! \brief Get the internal distributed grid
*
* \param lvl level
*
* \return the internal distributed grid
*
*/
grid_dist_id<dim,St,T,Decomposition,Memory,device_grid> & getDistGrid(size_t lvl)
{
return gd_array.get(lvl);
}
//////////////////// Remove functions
/*! \brief Remove a grid point (this function make sense only in case of
* sparse grid)
*
* \param v1 grid_key that identify the element in the AMR grid to eleminate
*
*/
inline void remove(const grid_dist_amr_key<dim> & v1)
{
#ifdef SE_CLASS2
check_valid(this,8);
#endif return gd_array.get(v1.getLvl()).remove(v1.getKey());
}
/*! \brief Remove a grid point (this function make sense only in case of
* sparse grid)
*
* \param v1 grid_key that identify the element in the AMR grid to eleminate
*
*/
void remove(size_t lvl, const grid_dist_key_dx<dim> & v1)
{
#ifdef SE_CLASS2
check_valid(this,8);
#endif
return gd_array.get(lvl).remove(v1);
}
/*! \brief construct level connections for padding particles
*
*
*/
void construct_level_connections()
{
for (int lvl = 0 ; lvl < gd_array.size() ; lvl++)
{
if (lvl == 0)
{
gd_array.get(lvl).construct_link_dw(gd_array.get(lvl+1),mv_off.get(lvl));
}
else if (lvl == gd_array.size() - 1)
{gd_array.get(lvl).construct_link_up(gd_array.get(lvl-1),mv_off.get(lvl));}
else
{
gd_array.get(lvl).construct_link_dw(gd_array.get(lvl+1),mv_off.get(lvl));
gd_array.get(lvl).construct_link_up(gd_array.get(lvl-1),mv_off.get(lvl));
}
}
}
/*! \brief construct level connections for padding particles
*
* \tparam stencil_type type of stencil
*
*/
template<typename stencil_type>
void tagBoundaries()
{
for (int lvl = 0 ; lvl < gd_array.size() ; lvl++)
{
gd_array.get(lvl).template tagBoundaries<stencil_type>();
}
}
//////////////////////////////////////
/*! \brief It synchronize the ghost parts
*
* \tparam prp... Properties to synchronize
*
*/
template<int... prp> void ghost_get(size_t opt = 0)
{
for (size_t i = 0 ; i < gd_array.size() ; i++)
{
gd_array.get(i).template ghost_get<prp...>(opt);
}
}
/*! \brief It move all the grid parts that do not belong to the local processor to the respective processor
* */
void map(size_t opt = 0)
{
for (size_t i = 0 ; i < gd_array.size() ; i++)
{
gd_array.get(i).map();
}
recalculate_mvoff();
}
/*! \brief Apply the ghost put
*
* \tparam prp... Properties to apply ghost put
*
*/
template<template<typename,typename> class op,int... prp> void ghost_put()
{
for (size_t i = 0 ; i < gd_array.size() ; i++)
{
gd_array.get(i).template ghost_put<op,prp...>();
}
}
/*! \brief Return the number of inserted points on a particular level
*
* \return the number of inserted points
*
*/
size_t size_inserted(size_t lvl)
{
return gd_array.get(lvl).size_local_inserted();
}
/*! \brief set the background value
*
* You can use this function make sense in case of sparse in case of dense
* it does nothing
*
*/
void setBackgroundValue(T & bv)
{
for (size_t i = 0 ; i < getNLvl() ; i++)
{gd_array.get(i).setBackgroundValue(bv);}
meta_copy<T>::meta_copy_(bv,bck);
}
/*! \brief delete all the points in the grid
*
* In case of sparse grid in delete all the inserted points, in case
* of dense it does nothing
*
*/
void clear()
{
for (size_t i = 0 ; i < getNLvl() ; i++)
{gd_array.get(i).clear();}
}
/*! \brief Get an object containing the grid informations for a specific level
*
* \param lvl level
*
* \return an information object about this grid
*
*/
const grid_sm<dim,void> & getGridInfoVoid(size_t lvl) const
{
return gd_array.get(lvl).getGridInfoVoid(); }
/*! \brief Return the maximum number of levels in the AMR struct
*
* \return the number of levels
*
*/
size_t getNLvl()
{
return gd_array.size();
}
/*! \brief Move down (to finer level) the key
*
* \param key multi-resolution AMR key
*
*/
void moveLvlDw(grid_dist_amr_key<dim> & key)
{
#ifdef SE_CLASS1
if (key.getLvl() >= getNLvl() - 1)
{std::cerr << __FILE__ << ":" << __LINE__ << " error: we are already at the last level, we cannot go one level down" << std::endl;}
#endif
auto & key_ref = key.getKeyRef().getKeyRef();
size_t lvl = key.getLvl();
for (size_t i = 0 ; i < dim ; i++)
{
key_ref.set_d(i,(key_ref.get(i) << 1) + mv_off.get(key.getLvl()).get(key.getKeyRef().getSub()).dw.get(i) );
}
key.setLvl(lvl+1);
}
/*! \brief Move down (to finer level) the key
*
* \param lvl level
* \param key multi-resolution AMR key
*
*/
grid_dist_key_dx<dim> moveDw(int lvl, const grid_dist_key_dx<dim> & key)
{
#ifdef SE_CLASS1
if (lvl >= getNLvl() - 1)
{std::cerr << __FILE__ << ":" << __LINE__ << " error: we are already at the last level, we cannot go one level down" << std::endl;}
#endif
grid_dist_key_dx<dim> out;
for (size_t i = 0 ; i < dim ; i++)
{
out.getKeyRef().set_d(i,(key.getKeyRef().get(i) << 1) + mv_off.get(lvl).get(key.getSub()).dw.get(i) );
}
out.setSub(key.getSub());
return out;
}
/*! \brief From a distributed key it return a AMR key that contain also the grid level
*
* \param lvl level
* \param key distributed key
*
*/ inline grid_dist_amr_key<dim> getAMRKey(size_t lvl, grid_dist_key_dx<dim> key)
{
return grid_dist_amr_key<dim>(lvl,key);
}
/*! \brief Move up (to coarser level) the key
*
* \param key multi-resolution AMR key
*
*/
void moveLvlUp(grid_dist_amr_key<dim> & key)
{
#ifdef SE_CLASS1
if (key.getLvl() == 0)
{std::cerr << __FILE__ << ":" << __LINE__ << " error: we are already at the top level, we cannot go one level up" << std::endl;}
#endif
auto & key_ref = key.getKeyRef().getKeyRef();
size_t lvl = key.getLvl();
for (size_t i = 0 ; i < dim ; i++)
{
key_ref.set_d(i,(key_ref.get(i) - mv_off.get(key.getLvl()).get(key.getKeyRef().getSub()).up.get(i)) >> 1);
}
key.setLvl(lvl-1);
}
/*! \brief Move up (to coarser level) the key
*
* \param lvl level
* \param key multi-resolution AMR key
*
*/
grid_dist_key_dx<dim> moveUp(int lvl, const grid_dist_key_dx<dim> & key)
{
#ifdef SE_CLASS1
if (lvl == 0)
{std::cerr << __FILE__ << ":" << __LINE__ << " error: we are already at the top level, we cannot go one level up" << std::endl;}
#endif
grid_dist_key_dx<dim> out;
for (size_t i = 0 ; i < dim ; i++)
{
out.getKeyRef().set_d(i,(key.getKeyRef().get(i) - mv_off.get(lvl).get(key.getSub()).up.get(i)) >> 1);
}
out.setSub(key.getSub());
return out;
}
/*! \brief Get the position on the grid in global coordinates
*
* \param v1 amr key
*
* \return the position in global coordinates
*
*/
grid_key_dx<dim> getGKey(const grid_dist_amr_key<dim> & v1)
{
return gd_array.get(v1.getLvl()).getGKey(v1.getKey());
}
/*! \brief Get the position on the grid in global coordinates *
* \param lvl level
* \param v1 point in the level
*
* \return the position in global coordinates
*
*/
grid_key_dx<dim> getGKey(int lvl, const grid_dist_key_dx<dim> & v1)
{
return gd_array.get(lvl).getGKey(v1);
}
/*! \brief Get the the position of the point
*
* \param v1 AMR point
*
* \return the position in global coordinates
*
*/
Point<dim,St> getPos(const grid_dist_amr_key<dim> & v1)
{
return gd_array.get(v1.getLvl()).getPos(v1.getKey());
}
/*! \brief Get the the position of the point
*
* \param v1 AMR point
*
* \return the position in global coordinates
*
*/
Point<dim,St> getPos(int lvl, const grid_dist_key_dx<dim> & v1)
{
return gd_array.get(lvl).getPos(v1);
}
/*! \brief return the spacing for the grid in the level lvl
*
* \param lvl level
*
* \return return the spacing
*
*/
Point<dim,St> getSpacing(size_t lvl)
{
return gd_array.get(lvl).getSpacing();
}
/* \brief Check if a point exist
*
* \param v1 point to checl
*
*/
bool existPoint(const grid_dist_amr_key<dim> & v1)
{
return gd_array.get(v1.getLvl()).existPoint(v1.getKey());
}
/* \brief Check if a point exist
*
* \param v1 point to checl
*
*/
bool existPoint(int lvl, const grid_dist_key_dx<dim> & v1)
{
return gd_array.get(lvl).existPoint(v1);
}
/*! \brief Write on vtk file *
* \param output filename output
*
*/
bool write(std::string output, size_t opt = VTK_WRITER | FORMAT_ASCII )
{
bool ret = true;
for (size_t i = 0 ; i < gd_array.size() ; i++)
{
ret &= gd_array.get(i).write(output + "_" + std::to_string(i),opt);
}
return ret;
}
#ifdef __NVCC__
/*! \brief Move the memory from the device to host memory
*
*/
template<unsigned int ... prp> void deviceToHost()
{
for (size_t i = 0 ; i < gd_array.size() ; i++)
{
gd_array.get(i).template deviceToHost<prp ...>();
}
}
/*! \brief Move the memory from the device to host memory
*
*/
template<unsigned int ... prp> void hostToDevice()
{
for (size_t i = 0 ; i < gd_array.size() ; i++)
{
gd_array.get(i).template hostToDevice<prp ...>();
}
}
#endif
};
template<unsigned int dim, typename St, typename T>
using sgrid_dist_amr = grid_dist_amr<dim,St,T,AMR_IMPL_TRIVIAL,CartDecomposition<dim,St>,HeapMemory,sgrid_cpu<dim,T,HeapMemory>>;
#ifdef __NVCC__
template<unsigned int dim, typename St, typename T, unsigned int blockEdgeSize = 8>
using sgrid_dist_amr_gpu = grid_dist_amr<dim,St,T,AMR_IMPL_TRIVIAL,CartDecomposition<dim,St,CudaMemory,memory_traits_inte>,CudaMemory,SparseGridGpu<dim,T,blockEdgeSize,IntPow<blockEdgeSize,dim>::value >>;
#endif
#endif /* AMR_GRID_AMR_DIST_HPP_ */