Commit bd5d64d2 authored by incardon's avatar incardon

CRS without debugging garbage

parent 7190e0c9
openfpm_data @ 155cc353
Subproject commit cb2d2030ab4d250e3e22cb27628582719b7a2740 Subproject commit 155cc3537485efb4993a667330a01a7d2644449a
...@@ -142,9 +142,6 @@ protected: ...@@ -142,9 +142,6 @@ protected:
//! Create distribution //! Create distribution
Distribution dist; Distribution dist;
//! Smallest subdivision on each direction
::Box<dim,T> ss_box;
//! Processor bounding box //! Processor bounding box
::Box<dim,T> bbox; ::Box<dim,T> bbox;
...@@ -261,10 +258,6 @@ public: ...@@ -261,10 +258,6 @@ public:
// optimize the decomposition // optimize the decomposition
d_o.template optimize<nm_v::sub_id, nm_v::proc_id>(dist.getGraph(), p_id, loc_box, box_nn_processor,ghe,bc); d_o.template optimize<nm_v::sub_id, nm_v::proc_id>(dist.getGraph(), p_id, loc_box, box_nn_processor,ghe,bc);
// reset ss_box
ss_box = domain;
ss_box -= ss_box.getP1();
// Initialize ss_box and bbox // Initialize ss_box and bbox
if (loc_box.size() >= 0) if (loc_box.size() >= 0)
{ {
...@@ -283,13 +276,9 @@ public: ...@@ -283,13 +276,9 @@ public:
// Calculate the bound box // Calculate the bound box
bbox.enclose(sub_d); bbox.enclose(sub_d);
proc_box.enclose(loc_box.get(s)); proc_box.enclose(loc_box.get(s));
// Create the smallest box contained in all sub-domain
ss_box.contained(sub_d);
} }
nn_prcs<dim,T>::create(box_nn_processor, sub_domains); nn_prcs<dim,T>::create(box_nn_processor, sub_domains);
nn_prcs<dim,T>::refine_ss_box(ss_box);
nn_prcs<dim,T>::applyBC(domain,ghost,bc); nn_prcs<dim,T>::applyBC(domain,ghost,bc);
// fill fine_s structure // fill fine_s structure
...@@ -318,8 +307,6 @@ public: ...@@ -318,8 +307,6 @@ public:
*/ */
void Initialize_geo_cell_lists() void Initialize_geo_cell_lists()
{ {
// Get the smallest sub-division on each direction
::Box<dim, T> unit = getSmallestSubdivision();
// Get the processor bounding Box // Get the processor bounding Box
::Box<dim,T> bound = getProcessorBounds(); ::Box<dim,T> bound = getProcessorBounds();
// Not necessary, but I prefer // Not necessary, but I prefer
...@@ -328,7 +315,7 @@ public: ...@@ -328,7 +315,7 @@ public:
// calculate the sub-divisions // calculate the sub-divisions
size_t div[dim]; size_t div[dim];
for (size_t i = 0; i < dim; i++) for (size_t i = 0; i < dim; i++)
div[i] = (size_t) ((bound.getHigh(i) - bound.getLow(i)) / unit.getHigh(i)); div[i] = (size_t) ((bound.getHigh(i) - bound.getLow(i)) / cd.getCellBox().getP2()[i]);
// Initialize the geo_cell structure // Initialize the geo_cell structure
ie_ghost<dim,T>::Initialize_geo_cell(bound,div); ie_ghost<dim,T>::Initialize_geo_cell(bound,div);
...@@ -732,7 +719,6 @@ public: ...@@ -732,7 +719,6 @@ public:
std::copy(spacing,spacing+3,cart.spacing); std::copy(spacing,spacing+3,cart.spacing);
cart.bbox = bbox; cart.bbox = bbox;
cart.ss_box = ss_box;
cart.ghost = g; cart.ghost = g;
cart.dist = dist; cart.dist = dist;
...@@ -773,7 +759,6 @@ public: ...@@ -773,7 +759,6 @@ public:
cart.ghost = ghost; cart.ghost = ghost;
cart.bbox = bbox; cart.bbox = bbox;
cart.ss_box = ss_box;
for (size_t i = 0 ; i < dim ; i++) for (size_t i = 0 ; i < dim ; i++)
cart.bc[i] = this->bc[i]; cart.bc[i] = this->bc[i];
...@@ -805,7 +790,6 @@ public: ...@@ -805,7 +790,6 @@ public:
ghost = cart.ghost; ghost = cart.ghost;
bbox = cart.bbox; bbox = cart.bbox;
ss_box = cart.ss_box;
for (size_t i = 0 ; i < dim ; i++) for (size_t i = 0 ; i < dim ; i++)
bc[i] = cart.bc[i]; bc[i] = cart.bc[i];
...@@ -837,7 +821,6 @@ public: ...@@ -837,7 +821,6 @@ public:
ghost = cart.ghost; ghost = cart.ghost;
bbox = cart.bbox; bbox = cart.bbox;
ss_box = cart.ss_box;
for (size_t i = 0 ; i < dim ; i++) for (size_t i = 0 ; i < dim ; i++)
bc[i] = cart.bc[i]; bc[i] = cart.bc[i];
...@@ -952,16 +935,6 @@ public: ...@@ -952,16 +935,6 @@ public:
return fine_s.get(cd.getCell(p)); return fine_s.get(cd.getCell(p));
} }
/*! \brief Get the smallest subdivision of the domain on each direction
*
* \return a box p1 is set to zero
*
*/
const ::Box<dim,T> & getSmallestSubdivision()
{
return ss_box;
}
/*! \brief Get the periodicity on i dimension /*! \brief Get the periodicity on i dimension
* *
* \param i dimension * \param i dimension
...@@ -1287,11 +1260,18 @@ public: ...@@ -1287,11 +1260,18 @@ public:
return domain_nn_calculator_cart<dim>::getDomainCells(shift,cell_shift,gs,proc_box,loc_box); return domain_nn_calculator_cart<dim>::getDomainCells(shift,cell_shift,gs,proc_box,loc_box);
} }
/*! \brief Get the domain anomalous cells /*! \brief Get the anomalous cells
*
* This function include also a linearization of the indexes
* *
* \param shift Shifting point * \param shift Shifting point
* \param cell_shift where the processor cell-list start (In case of symmetric)
* exist one global cell-list, but each processor span only one
* part of it
* \param gs grid extension * \param gs grid extension
* *
* \return the anomalous cells with neighborhood
*
*/ */
openfpm::vector<subsub_lin<dim>> & getAnomDomainCells(grid_key_dx<dim> & shift, grid_key_dx<dim> & cell_shift, grid_sm<dim,void> & gs) openfpm::vector<subsub_lin<dim>> & getAnomDomainCells(grid_key_dx<dim> & shift, grid_key_dx<dim> & cell_shift, grid_sm<dim,void> & gs)
{ {
......
...@@ -42,8 +42,8 @@ private: ...@@ -42,8 +42,8 @@ private:
* *
* \see duplicate (in case of extended domain) * \see duplicate (in case of extended domain)
* *
* \param cart Cartesian decomposition object * \param dec Cartesian decomposition object
* \param box Extended domain * \param ext_dom Extended domain
* *
*/ */
void extend_subdomains(const CartDecomposition<dim,T,Memory,Distribution> & dec, const ::Box<dim,T> & ext_dom) void extend_subdomains(const CartDecomposition<dim,T,Memory,Distribution> & dec, const ::Box<dim,T> & ext_dom)
...@@ -53,12 +53,6 @@ private: ...@@ -53,12 +53,6 @@ private:
this->bbox.zero(); this->bbox.zero();
for (size_t i = 0 ; i < dim ; i++)
{
this->ss_box.setLow(i,0.0);
this->ss_box.setHigh(i,ext_dom.getHigh(i) - ext_dom.getLow(i));
}
// Extend sub-domains // Extend sub-domains
for (size_t i = 0 ; i < dec.sub_domains.size() ; i++) for (size_t i = 0 ; i < dec.sub_domains.size() ; i++)
{ {
...@@ -83,16 +77,12 @@ private: ...@@ -83,16 +77,12 @@ private:
// Calculate the bound box // Calculate the bound box
this->bbox.enclose(box); this->bbox.enclose(box);
// Create the smallest box contained in all sub-domain
this->ss_box.contained(box);
} }
} }
/*! \brief Extend the fines for the new Cartesian decomposition /*! \brief Extend the fines for the new Cartesian decomposition
* *
* \param new_fines extended fine_s * \param dec Non-extended decomposition
* \param old_fines old fine_s
* *
*/ */
void extend_fines(const CartDecomposition<dim,T> & dec) void extend_fines(const CartDecomposition<dim,T> & dec)
...@@ -171,6 +161,7 @@ public: ...@@ -171,6 +161,7 @@ public:
{ {
} }
//! The non-extended decomposition base class
typedef CartDecomposition<dim,T,Memory,Distribution> base_type; typedef CartDecomposition<dim,T,Memory,Distribution> base_type;
/*! \brief It create another object that contain the same decomposition information but with different ghost boxes and an extended domain /*! \brief It create another object that contain the same decomposition information but with different ghost boxes and an extended domain
...@@ -179,36 +170,37 @@ public: ...@@ -179,36 +170,37 @@ public:
* *
* \verbatim * \verbatim
* *
+--------------^--------^----------^----------+ +--------------^--------^----------^----------+
| | | | | | | | | |
| A | E | F | N | | A | E | F | N |
| +-----------------------------------+----> | +-----------------------------------+---->
| | | | | | | | | | | | | |
| A | A | | F | | | | A | A | | F | | |
| | | | | | | | | | | | | |
| | | E +----------+ N | N | | | | E +----------+ N | N |
<--------------+ | | | | <--------------+ | | | |
| | | | | | | | | | | | | |
| | | | G | | | | | | | G | | |
| | | | +----------> | | | | +---------->
| B | B | +----------+ | | | B | B | +----------+ | |
| | +--------+ | M | M | | | +--------+ | M | M |
| | | | H | | | | | | | H | | |
| | | +-----+----+----------> | | | +-----+----+---------->
<--------------+ D | | | | <--------------+ D | | | |
| | | | I | L | L | | | | | I | L | L |
| C | C | | | | | | C | C | | | | |
| | | | | | | | | | | | | |
| +-----------------------------------+ | | +-----------------------------------+ |
| | | | | | | | | |
| C | D | I | L | | C | D | I | L |
+--------------v--------v-----v---------------+ +--------------v--------v-----v---------------+
* *
* \endverbatim * \endverbatim
* *
* \param dec Decomposition
* \param g ghost * \param g ghost
* \param domain extended domain (MUST be extended) * \param ext_domain extended domain (MUST be extended)
* *
* \return a duplicated decomposition with different ghost boxes and an extended domain * \return a duplicated decomposition with different ghost boxes and an extended domain
* *
......
...@@ -338,23 +338,6 @@ public: ...@@ -338,23 +338,6 @@ public:
return *this; return *this;
} }
/*! \brief Refine the ss_box to have the smallest size on each direction of the local sub-domains and adjacent (from other processor) one
*
* \param ss_box box that store the smallest size of the sub-domain
*
*/
void refine_ss_box(Box<dim,T> & ss_box)
{
for (size_t p = 0 ; p < getNNProcessors() ; p++)
{
const openfpm::vector< ::Box<dim,T> > & list_p_box = getNearSubdomains(IDtoProc(p));
// Create the smallest box contained in all sub-domain
for (size_t b = 0 ; b < list_p_box.size() ; b++)
ss_box.contained(list_p_box.get(b));
}
}
/*! \brief Create the list of adjacent processors and the list of adjacent sub-domains /*! \brief Create the list of adjacent processors and the list of adjacent sub-domains
* *
* \param box_nn_processors * \param box_nn_processors
......
...@@ -25,6 +25,7 @@ ...@@ -25,6 +25,7 @@
template<unsigned int dim, typename G_v, int prp> template<unsigned int dim, typename G_v, int prp>
struct fill_id struct fill_id
{ {
//! function that fill with linearization indexes
static inline void fill(G_v & g_v, const grid_key_dx<dim> & gk, const grid_sm<dim, void> & gs) static inline void fill(G_v & g_v, const grid_key_dx<dim> & gk, const grid_sm<dim, void> & gs)
{ {
g_v.template get<prp>() = gs.LinId(gk); g_v.template get<prp>() = gs.LinId(gk);
...@@ -38,6 +39,7 @@ struct fill_id ...@@ -38,6 +39,7 @@ struct fill_id
template<unsigned int dim, typename G_v> template<unsigned int dim, typename G_v>
struct fill_id<dim, G_v, NO_VERTEX_ID> struct fill_id<dim, G_v, NO_VERTEX_ID>
{ {
//! function that fill with linearization indexes
static inline void fill(G_v & g_v, const grid_key_dx<dim> & gk, const grid_sm<dim, void> & gs) static inline void fill(G_v & g_v, const grid_key_dx<dim> & gk, const grid_sm<dim, void> & gs)
{ {
} }
...@@ -219,7 +221,7 @@ public: ...@@ -219,7 +221,7 @@ public:
/*! \brief Operator for vector and scalar property /*! \brief Operator for vector and scalar property
* *
* \tparam i Size of the property * \tparam i Size of the property
* \tparam p Type of the property * \tparam p Type of the property boost mpl
* \tparam Graph Graph * \tparam Graph Graph
* \tparam pos Array of properties * \tparam pos Array of properties
*/ */
...@@ -227,7 +229,10 @@ template<int i, typename p, typename Graph, int ... pos> ...@@ -227,7 +229,10 @@ template<int i, typename p, typename Graph, int ... pos>
struct fill_prop_by_type struct fill_prop_by_type
{ {
//! Get the element 0
typedef typename boost::mpl::at<p, boost::mpl::int_<0>>::type v_element; typedef typename boost::mpl::at<p, boost::mpl::int_<0>>::type v_element;
//! Get the property v_element (v_element is a number)
typedef typename boost::mpl::at<typename Graph::V_type::type, v_element>::type pos_prop_type; typedef typename boost::mpl::at<typename Graph::V_type::type, v_element>::type pos_prop_type;
enum enum
...@@ -266,7 +271,16 @@ template<unsigned int dim, int lin_id, typename Graph, int se, typename T, unsig ...@@ -266,7 +271,16 @@ template<unsigned int dim, int lin_id, typename Graph, int se, typename T, unsig
class Graph_constructor_impl class Graph_constructor_impl
{ {
public: public:
//! Construct cartesian graph
/*! \brief Construct a cartesian graph
*
* \param sz size of the partesian graph
* \param dom domain where this cartesian graph is defined (used to fill the coordinates)
* \param bc boundary conditions (torus or cube)
*
* \return the constructed graph
*
*/
static Graph construct(const size_t (& sz)[dim], Box<dim,T> dom, const size_t(& bc)[dim]) static Graph construct(const size_t (& sz)[dim], Box<dim,T> dom, const size_t(& bc)[dim])
{ {
// Calculate the size of the hyper-cubes on each dimension // Calculate the size of the hyper-cubes on each dimension
...@@ -378,7 +392,16 @@ template<unsigned int dim, int lin_id, typename Graph, typename T, unsigned int ...@@ -378,7 +392,16 @@ template<unsigned int dim, int lin_id, typename Graph, typename T, unsigned int
class Graph_constructor_impl<dim, lin_id, Graph, NO_EDGE, T, dim_c, pos...> class Graph_constructor_impl<dim, lin_id, Graph, NO_EDGE, T, dim_c, pos...>
{ {
public: public:
//! Construct cartesian graph
/*! \brief Construct a cartesian graph
*
* \param sz size of the partesian graph
* \param dom domain where this cartesian graph is defined (used to fill the coordinates)
* \param bc boundary conditions (torus or cube)
*
* \return the constructed graph
*
*/
static Graph construct(const size_t ( & sz)[dim], Box<dim,T> dom, const size_t(& bc)[dim]) static Graph construct(const size_t ( & sz)[dim], Box<dim,T> dom, const size_t(& bc)[dim])
{ {
// Calculate the size of the hyper-cubes on each dimension // Calculate the size of the hyper-cubes on each dimension
...@@ -492,48 +515,20 @@ public: ...@@ -492,48 +515,20 @@ public:
* dim_c. One property can be used to store the contact size or the d-dimensional * dim_c. One property can be used to store the contact size or the d-dimensional
* surface in common between two connected hyper-cube. * surface in common between two connected hyper-cube.
* *
* \param sz Vector that store the size of the grid on each dimension
* \param dom Box enclosing the physical domain
*
* \tparam se Indicate which properties fill with the contact size. The * \tparam se Indicate which properties fill with the contact size. The
* contact size is the point, line , surface, d-dimensional object size * contact size is the point, line , surface, d-dimensional object size
* in contact (in common) between two hyper-cube. NO_EDGE indicate * in contact (in common) between two hyper-cube. NO_EDGE indicate
* no property will store this information * no property will store this information
* \tparam id_prp property 'id' that stores the vertex id (with -1 it skip)
* \tparam T type of the domain like (int real complex ... ) * \tparam T type of the domain like (int real complex ... )
* \tparam dim_c Connectivity dimension * \tparam dim_c Connectivity dimension
* \tparam pos... (optional)one or more integer indicating the spatial properties * \tparam pos... (optional)one or more integer indicating the spatial properties
* *
*/ * \param sz store the size of the cartesian grid on each dimension
/* template <int se,typename T, unsigned int dim_c, int... pos>
static Graph construct(const size_t (& sz)[dim], Box<dim,T> dom )
{
return Graph_constructor_impl<dim,Graph,se,T,dim_c,pos...>::construct(sz,dom,bc);
}*/
/*!
*
* \brief Construct a cartesian graph, with V and E edge properties
*
* Construct a cartesian graph, with V and E edge properties
*
* Each vertex is a subspace (Hyper-cube) of dimension dim, each vertex is
* connected with an edge if two vertex (Hyper-cube) share a element of dimension grater than
* dim_c. One property can be used to store the contact size or the d-dimensional
* surface in common between two connected hyper-cube.
*
* \param sz Vector that store the size of the grid on each dimension
* \param dom Box enclosing the physical domain * \param dom Box enclosing the physical domain
* \param bc boundary conditions {PERIODIC and NON_PERIODIC} * \param bc boundary conditions {PERIODIC = torus and NON_PERIODIC = cube}
* *
* \tparam se Indicate which properties fill with the contact size. The * \return the constructed graph
* contact size is the point, line , surface, d-dimensional object size
* in contact (in common) between two hyper-cube. NO_EDGE indicate
* no property will store this information
* \tparam id_prp property 'id' that stores the vertex id (with -1 it skip)
* \tparam T type of the domain like (int real complex ... )
* \tparam dim_c Connectivity dimension
* \tparam pos... (optional)one or more integer indicating the spatial properties
* *
*/ */
template<int se, int id_prp, typename T, unsigned int dim_c, int ... pos> template<int se, int id_prp, typename T, unsigned int dim_c, int ... pos>
......
...@@ -552,8 +552,8 @@ BOOST_AUTO_TEST_CASE( vector_dist_symmetric_crs_cell_list ) ...@@ -552,8 +552,8 @@ BOOST_AUTO_TEST_CASE( vector_dist_symmetric_crs_cell_list )
long int big_step = k / 4; long int big_step = k / 4;
big_step = (big_step == 0)?1:big_step; big_step = (big_step == 0)?1:big_step;
print_test("Testing 3D periodic vector symmetric cell-list k=",k); print_test("Testing 3D periodic vector symmetric crs cell-list k=",k);
BOOST_TEST_CHECKPOINT( "Testing 3D periodic vector symmetric cell-list k=" << k ); BOOST_TEST_CHECKPOINT( "Testing 3D periodic vector symmetric crs cell-list k=" << k );
Box<3,float> box({-L,-L,-L},{L,L,L}); Box<3,float> box({-L,-L,-L},{L,L,L});
...@@ -655,31 +655,6 @@ BOOST_AUTO_TEST_CASE( vector_dist_symmetric_crs_cell_list ) ...@@ -655,31 +655,6 @@ BOOST_AUTO_TEST_CASE( vector_dist_symmetric_crs_cell_list )
auto NN2 = vd.getCellListSym(r_cut); auto NN2 = vd.getCellListSym(r_cut);
int debug = 0;
debug++;
auto debug_it = vd.getDomainAndGhostIterator();
while (debug_it.isNext())
{
auto key = debug_it.get();
if (vd.getProp<2>(key) == 1698)
{
int debug = 0;
debug++;
}
++debug_it;
}
vd.write("debug_decomp_crs_part");
vd.getDecomposition().write("debug_decomp_crs");
std::cout << "NN2: " << NN2.getCell(vd.getPos(4853)) << std::endl;
std::cout << "NN2 Grid: " << NN2.getCellGrid(vd.getPos(4853)).to_string() << std::endl;
int debug_cnt = 0;
// In case of CRS we have to iterate particles within some cells // In case of CRS we have to iterate particles within some cells
// here we define whichone // here we define whichone
auto p_it2 = vd.getParticleIteratorCRS(NN2); auto p_it2 = vd.getParticleIteratorCRS(NN2);
...@@ -689,16 +664,8 @@ BOOST_AUTO_TEST_CASE( vector_dist_symmetric_crs_cell_list ) ...@@ -689,16 +664,8 @@ BOOST_AUTO_TEST_CASE( vector_dist_symmetric_crs_cell_list )
{ {
auto p = p_it2.get(); auto p = p_it2.get();
debug_cnt++;
Point<3,float> xp = vd.getPos(p); Point<3,float> xp = vd.getPos(p);
if (debug_cnt == 3762 && v_cl.getProcessUnitID() == 0)
{
std::cout << " P local " << p << " " << vd.getProp<2>(p) << std::endl;
int debug = 0;
debug++;
}
auto Np = p_it2.getNNIteratorCSR(vd.getPosVector()); auto Np = p_it2.getNNIteratorCSR(vd.getPosVector());
while (Np.isNext()) while (Np.isNext())
...@@ -732,16 +699,6 @@ BOOST_AUTO_TEST_CASE( vector_dist_symmetric_crs_cell_list ) ...@@ -732,16 +699,6 @@ BOOST_AUTO_TEST_CASE( vector_dist_symmetric_crs_cell_list )
vd.getProp<4>(q).last().xq = xp; vd.getProp<4>(q).last().xq = xp;
vd.getProp<4>(p).last().id = vd.getProp<2>(q); vd.getProp<4>(p).last().id = vd.getProp<2>(q);
vd.getProp<4>(q).last().id = vd.getProp<2>(p); vd.getProp<4>(q).last().id = vd.getProp<2>(p);
if (v_cl.getProcessUnitID() == 0 && vd.getProp<2>(p) == 1698)
{
std::cerr << "NN " << vd.getProp<2>(q) << std::endl;
}
if (v_cl.getProcessUnitID() == 0 && vd.getProp<2>(q) == 1689)
{
std::cerr << "NN " << vd.getProp<2>(p) << std::endl;
}
} }
++Np; ++Np;
...@@ -767,18 +724,18 @@ BOOST_AUTO_TEST_CASE( vector_dist_symmetric_crs_cell_list ) ...@@ -767,18 +724,18 @@ BOOST_AUTO_TEST_CASE( vector_dist_symmetric_crs_cell_list )
ret &= vd.getProp<3>(p).size() == vd.getProp<4>(p).size(); ret &= vd.getProp<3>(p).size() == vd.getProp<4>(p).size();
if (v_cl.getProcessUnitID() == 0) // if (v_cl.getProcessUnitID() == 0)
{ // {
std::cerr << "Particle: " << p.getKey() << " Position: " << Point<3,float>(vd.getPos(p)).toString() << std::endl; // std::cerr << "Particle: " << p.getKey() << " Position: " << Point<3,float>(vd.getPos(p)).toString() << std::endl;
for (size_t i = 0 ; i < vd.getProp<4>(p).size() ; i++) for (size_t i = 0 ; i < vd.getProp<4>(p).size() ; i++)
{ {
std::cerr << "POSITION: " << vd.getProp<3>(p).get(i).xq.toString() << std::endl; // std::cerr << "POSITION: " << vd.getProp<3>(p).get(i).xq.toString() << std::endl;
std::cerr << "ID nn " << vd.getProp<3>(p).get(i).id << " " << vd.getProp<4>(p).get(i).id << std::endl; // std::cerr << "ID nn " << vd.getProp<3>(p).get(i).id << " " << vd.getProp<4>(p).get(i).id << std::endl;
ret &= vd.getProp<3>(p).get(i).id == vd.getProp<4>(p).get(i).id; ret &= vd.getProp<3>(p).get(i).id == vd.getProp<4>(p).get(i).id;
} }
} // }
if (ret == false) if (ret == false)
break; break;
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment