From 77da8be43b9776f010d1954de8961e4ee6d0c849 Mon Sep 17 00:00:00 2001
From: Pietro Incardona <incardon@mpi-cbg.de>
Date: Thu, 21 Jun 2018 19:19:00 +0200
Subject: [PATCH] Fixing bug of consistency between map and ghost_get

---
 configure.ac                                  |   1 +
 images/Makefile.am                            |   2 +-
 script/solve_gpp                              |   4 +-
 src/Decomposition/CartDecomposition.hpp       | 147 +++++++++++++++---
 src/Decomposition/CartDecomposition_ext.hpp   |  11 +-
 src/Decomposition/common.hpp                  |   8 +
 src/Decomposition/ie_ghost.hpp                |   2 +-
 .../tests/CartDecomposition_unit_test.cpp     |  47 ++++++
 src/Grid/grid_dist_id_comm.hpp                |   2 +-
 src/Makefile.am                               |   2 +-
 10 files changed, 200 insertions(+), 26 deletions(-)

diff --git a/configure.ac b/configure.ac
index 56b18d19..a73caf41 100644
--- a/configure.ac
+++ b/configure.ac
@@ -232,6 +232,7 @@ AC_COMPILE_IFELSE([AC_LANG_PROGRAM([])],
 )
 AC_LANG_POP([C++])
 CXXFLAGS="$my_save_cflags"
+
 AC_SUBST([AM_CXXFLAGS])
 
 
diff --git a/images/Makefile.am b/images/Makefile.am
index 732e1211..afa9f3a5 100644
--- a/images/Makefile.am
+++ b/images/Makefile.am
@@ -1,4 +1,4 @@
-LINKLIBS = $(HDF5_LDFLAGS)  $(HDF5_LIBS) $(OPENMP_LDFLAGS) $(LIBHILBERT_LIB) $(PETSC_LIB) $(METIS_LIB) $(PARMETIS_LIB)  $(PTHREAD_LIBS) $(OPT_LIBS) $(BOOST_LDFLAGS) $(BOOST_IOSTREAMS_LIB) $(CUDA_LIBS)
+LINKLIBS = $(HDF5_LDFLAGS)  $(HDF5_LIBS) $(OPENMP_LDFLAGS) $(LIBHILBERT_LIB) $(PETSC_LIB) $(SUITESPARSE_LIBS) $(LAPACK_LIBS) $(BLAS_LIBS)  $(METIS_LIB) $(PARMETIS_LIB)  $(PTHREAD_LIBS) $(OPT_LIBS) $(BOOST_LDFLAGS) $(BOOST_IOSTREAMS_LIB) $(CUDA_LIBS)
 
 noinst_PROGRAMS = cart_dec metis_dec dom_box vector_dist
 cart_dec_SOURCES = CartDecomposition_gen_vtk.cpp ../src/lib/pdata.cpp ../openfpm_devices/src/memory/HeapMemory.cpp ../openfpm_devices/src/memory/PtrMemory.cpp ../openfpm_vcluster/src/VCluster/VCluster.cpp ../openfpm_devices/src/Memleak_check.cpp
diff --git a/script/solve_gpp b/script/solve_gpp
index 6630e0c8..72060944 100755
--- a/script/solve_gpp
+++ b/script/solve_gpp
@@ -20,8 +20,8 @@ if [ x"$1" = x"osx" ]; then
 elif [ x"$1" = x"linux"  ]; then
 
 	if [ x"$pcman" == x"pacman" ]; then
-          commands[0]="su -c \"$pcman -Sy gcc sudo make\""
-          commands[1]="sudo $pcman -Sy gcc make"
+          commands[0]="su -c \"$pcman -Sy gcc make\""
+          commands[1]="sudo $pcman -Sy gcc make mpfr"
           possible_solutions "${commands[@]}"
           compiler_gcc="gcc"
           compiler_gpp="g++"
diff --git a/src/Decomposition/CartDecomposition.hpp b/src/Decomposition/CartDecomposition.hpp
index 6edb0bc8..ae045e6e 100755
--- a/src/Decomposition/CartDecomposition.hpp
+++ b/src/Decomposition/CartDecomposition.hpp
@@ -131,7 +131,6 @@ template<unsigned int dim> static void nsub_to_div(size_t (& div)[dim], size_t n
 template<unsigned int dim, typename T, typename Memory, typename Distribution>
 class CartDecomposition: public ie_loc_ghost<dim, T>, public nn_prcs<dim, T>, public ie_ghost<dim, T>, public domain_nn_calculator_cart<dim>
 {
-
 public:
 
 	//! Type of the domain we are going to decompose
@@ -163,15 +162,15 @@ protected:
 	//! the set of all local sub-domain as vector
 	openfpm::vector<SpaceBox<dim, T>> sub_domains;
 
-	//! the global set of all sub-domains as vector of 'sub_domains' vectors
-	mutable openfpm::vector<openfpm::vector<SpaceBox<dim, T>>> sub_domains_global;
+	//! the remote set of all sub-domains as vector of 'sub_domains' vectors
+	mutable openfpm::vector<Box_map<dim, T>> sub_domains_global;
 
 	//! for each sub-domain, contain the list of the neighborhood processors
 	openfpm::vector<openfpm::vector<long unsigned int> > box_nn_processor;
 
 	//! Structure that contain for each sub-sub-domain box the processor id
 	//! exist for efficient global communication
-	openfpm::vector<size_t> fine_s;
+	CellList<dim,T,Mem_fast<>,shift<dim,T>> fine_s;
 
 	//! Structure that store the cartesian grid information
 	grid_sm<dim, void> gr;
@@ -269,9 +268,79 @@ protected:
 		return sub_d;
 	}
 
+	void collect_all_sub_domains(openfpm::vector<Box_map<dim,T>> & sub_domains_global)
+	{
+#ifdef SE_CLASS2
+		check_valid(this,8);
+#endif
+
+		sub_domains_global.clear();
+		openfpm::vector<Box_map<dim,T>> bm;
+
+		for (size_t i = 0 ; i < sub_domains.size() ; i++)
+		{
+			Box_map<dim,T> tmp;
+			tmp.box = ::SpaceBox<dim,T>(sub_domains.get(i));
+			tmp.prc = v_cl.rank();
+
+			bm.add(tmp);
+
+		}
+
+		v_cl.SGather(bm,sub_domains_global,0);
+
+		size_t size = sub_domains_global.size();
+
+		v_cl.max(size);
+		v_cl.execute();
+
+		sub_domains_global.resize(size);
+
+		v_cl.Bcast(sub_domains_global,0);
+		v_cl.execute();
+	}
 
 public:
 
+	void initialize_fine_s(const ::Box<dim,T> & domain)
+	{
+		fine_s.clear();
+		size_t div_g[dim];
+
+		// We reduce the size of the cells by a factor 8 in 3d 4 in 2d
+		for (size_t i = 0 ; i < dim ; i++)
+		{div_g[i] = gr.size(i)/2;}
+
+		fine_s.Initialize(domain,div_g);
+	}
+
+	void construct_fine_s()
+	{
+		collect_all_sub_domains(sub_domains_global);
+
+		// now draw all sub-domains in fine-s
+
+		for (size_t i = 0 ; i < sub_domains_global.size() ; i++)
+		{
+
+			// get the cells this box span
+			const grid_key_dx<dim> p1 = fine_s.getCellGrid(sub_domains_global.get(i).box.getP1());
+			const grid_key_dx<dim> p2 = fine_s.getCellGrid(sub_domains_global.get(i).box.getP2());
+
+			// Get the grid and the sub-iterator
+			auto & gi = fine_s.getGrid();
+			grid_key_dx_iterator_sub<dim> g_sub(gi,p1,p2);
+
+			// add the box-id to the cell list
+			while (g_sub.isNext())
+			{
+				auto key = g_sub.get();
+				fine_s.addCell(gi.LinId(key),i);
+				++g_sub;
+			}
+		}
+	}
+
 	/*! \brief Constructor, it decompose and distribute the sub-domains across the processors
 	 *
 	 * \param v_cl Virtual cluster, used internally for communications
@@ -295,7 +364,7 @@ public:
 		}
 
 		// fill the structure that store the processor id for each sub-domain
-		fine_s.resize(gr.size());
+		initialize_fine_s(domain);
 
 		// Optimize the decomposition creating bigger spaces
 		// And reducing Ghost over-stress
@@ -355,8 +424,13 @@ public:
 		// with sub-sub-domain we mean the sub-domain decomposition before
 		// running dec_optimizer (before merging sub-domains)
 
+		///////////////////////////////// TODO //////////////////////////////////////////
+
+		construct_fine_s();
 
-		grid_key_dx_iterator<dim> git(gr);
+		/////////////////////////////////////////////////////////////////////////////////
+
+/*		grid_key_dx_iterator<dim> git(gr);
 
 		while (git.isNext())
 		{
@@ -364,15 +438,17 @@ public:
 			grid_key_dx<dim> key2;
 
 			for (size_t i = 0 ; i < dim ; i++)
-				key2.set_d(i,key.get(i) / magn[i]);
+			{key2.set_d(i,key.get(i) / magn[i]);}
 
 			size_t lin = gr_dist.LinId(key2);
 			size_t lin2 = gr.LinId(key);
 
+			// Here we draw the fine_s in the cell-list
+
 			fine_s.get(lin2) = dist.getGraph().template vertex_p<nm_v::proc_id>(lin);
 
 			++git;
-		}
+		}*/
 
 		Initialize_geo_cell_lists();
 	}
@@ -590,6 +666,38 @@ public:
 		ie_loc_ghost<dim,T>::create(sub_domains,domain,ghost,bc);
 	}
 
+	template<typename T2> inline size_t processorID_impl(T2 & p) const
+	{
+		// Get the number of elements in the cell
+
+		size_t e;
+		size_t cl = fine_s.getCell(p);
+		size_t n_ele = fine_s.getNelements(cl);
+
+		for (size_t i = 0 ; i < n_ele ; i++)
+		{
+			e = fine_s.get(cl,i);
+
+			if (sub_domains_global.get(e).box.isInsideNP(p) == true)
+			{
+				break;
+			}
+		}
+
+#ifdef SE_CLASS1
+
+		if (n_ele == 0)
+		{
+			std::cout << __FILE__ << ":" << __LINE__ << " I cannot detect in which processor this particle go" << std::endl;
+			return -1;
+		}
+
+#endif
+
+		return sub_domains_global.get(e).prc;
+	}
+
+
 public:
 
 	//! Space dimensions
@@ -962,7 +1070,7 @@ public:
 	 */
 	template<typename Mem> size_t inline processorID(const encapc<1, Point<dim,T>, Mem> & p) const
 	{
-		return fine_s.get(cd.template getCell(p));
+		return processorID_impl(p);
 	}
 
 	/*! \brief Given a point return in which processor the particle should go
@@ -974,7 +1082,7 @@ public:
 	 */
 	size_t inline processorID(const Point<dim,T> &p) const
 	{
-		return fine_s.get(cd.getCell(p));
+		return processorID_impl(p);
 	}
 
 	/*! \brief Given a point return in which processor the particle should go
@@ -986,7 +1094,7 @@ public:
 	 */
 	size_t inline processorID(const T (&p)[dim]) const
 	{
-		return fine_s.get(cd.getCell(p));
+		return processorID_impl(p);
 	}
 
 	/*! \brief Given a point return in which processor the point/particle should go
@@ -1003,7 +1111,8 @@ public:
 		Point<dim,T> pt = p;
 		applyPointBC(pt);
 
-		return fine_s.get(cd.getCell(pt));
+
+		return processorID_impl(p);
 	}
 
 	/*! \brief Given a point return in which processor the particle should go
@@ -1015,12 +1124,14 @@ public:
 	 * \return processorID
 	 *
 	 */
-	template<typename ofb> size_t inline processorIDBC(const Point<dim,T> &p) const
+	size_t inline processorIDBC(const Point<dim,T> &p) const
 	{
 		Point<dim,T> pt = p;
 		applyPointBC(pt);
 
-		return fine_s.get(cd.getCell(p));
+		// Get the number of elements in the cell
+
+		return processorID_impl(p);
 	}
 
 	/*! \brief Given a point return in which processor the particle should go
@@ -1032,12 +1143,12 @@ public:
 	 * \return processorID
 	 *
 	 */
-	template<typename ofb> size_t inline processorIDBC(const T (&p)[dim]) const
+	size_t inline processorIDBC(const T (&p)[dim]) const
 	{
 		Point<dim,T> pt = p;
 		applyPointBC(pt);
 
-		return fine_s.get(cd.getCell(p));
+		return processorID_impl(p);
 	}
 
 	/*! \brief Get the periodicity on i dimension
@@ -1452,10 +1563,10 @@ public:
 		return sub_domains;
 	}
 
-	openfpm::vector<openfpm::vector<SpaceBox<dim, T>>> & getSubDomainsGlobal()
+/*	openfpm::vector<openfpm::vector<SpaceBox<dim, T>>> & getSubDomainsGlobal()
 	{
 		return sub_domains_global;
-	}
+	}*/
 
 	/*! \brief Check if the particle is local
 	 *
diff --git a/src/Decomposition/CartDecomposition_ext.hpp b/src/Decomposition/CartDecomposition_ext.hpp
index 7af3df06..4d3012b5 100644
--- a/src/Decomposition/CartDecomposition_ext.hpp
+++ b/src/Decomposition/CartDecomposition_ext.hpp
@@ -90,7 +90,7 @@ private:
 	 * \param dec Non-extended decomposition
 	 *
 	 */
-	void extend_fines(const CartDecomposition<dim,T,Memory,Distribution> & dec)
+/*	void extend_fines(const CartDecomposition<dim,T,Memory,Distribution> & dec)
 	{
 		// Extension, first we calculate the extensions of the new domain compared
 		// to the old one in cell units (each cell unit is a sub-sub-domain)
@@ -152,6 +152,12 @@ private:
 
 		// the new extended CellDecomposer must be consistent with the old cellDecomposer.
 		this->cd.setDimensions(dec.cd,ext);
+	}*/
+
+	void reconstruct_fine_s_from_extended_domain(const ::Box<dim,T> & ext_domain)
+	{
+		this->initialize_fine_s(ext_domain);
+		this->construct_fine_s();
 	}
 
 public:
@@ -219,7 +225,8 @@ public:
 
 		// Calculate fine_s structure for the extended domain
 		// update the cell decomposer and gr
-		extend_fines(dec);
+//		extend_fines(dec);
+		reconstruct_fine_s_from_extended_domain(ext_domain);
 
 		// Get the old sub-sub-domain grid extension
 
diff --git a/src/Decomposition/common.hpp b/src/Decomposition/common.hpp
index 5a920679..c3c67afc 100755
--- a/src/Decomposition/common.hpp
+++ b/src/Decomposition/common.hpp
@@ -114,6 +114,14 @@ struct Box_sub_k
 	}
 };
 
+template<unsigned int dim,typename T>
+struct Box_map
+{
+	Box<dim,T> box;
+
+	long int prc;
+};
+
 //! Case for local ghost box
 template<unsigned int dim, typename T>
 struct lBox_dom
diff --git a/src/Decomposition/ie_ghost.hpp b/src/Decomposition/ie_ghost.hpp
index 7863bdb8..eccfee09 100755
--- a/src/Decomposition/ie_ghost.hpp
+++ b/src/Decomposition/ie_ghost.hpp
@@ -790,7 +790,7 @@ public:
 		{
 			size_t bid = cell_it.get();
 
-			if (vb_int.get(bid).box.isInside(p) == true)
+			if (vb_int.get(bid).box.isInsideNP(p) == true)
 			{
 				ids_p.add(std::pair<size_t,size_t>(id1::id(vb_int.get(bid),bid),id2::id(vb_int.get(bid),bid)));
 			}
diff --git a/src/Decomposition/tests/CartDecomposition_unit_test.cpp b/src/Decomposition/tests/CartDecomposition_unit_test.cpp
index e66bbbf3..4b94973e 100755
--- a/src/Decomposition/tests/CartDecomposition_unit_test.cpp
+++ b/src/Decomposition/tests/CartDecomposition_unit_test.cpp
@@ -341,6 +341,53 @@ BOOST_AUTO_TEST_CASE( CartDecomposition_ext_non_periodic_test)
 	}
 }
 
+BOOST_AUTO_TEST_CASE( CartDecomposition_check_cross_consistency_between_proc_idbc_and_ghost )
+{
+	// Vcluster
+	Vcluster & vcl = create_vcluster();
+
+	if (vcl.size() != 3)
+	{return;}
+
+	CartDecomposition<3, double> dec(vcl);
+
+	size_t bc[3] = {PERIODIC,PERIODIC,PERIODIC};
+
+	// Physical domain
+	Box<3, double> box( { -0.01, -0.01, 0.0 }, { 0.01, 0.01, 0.003 });
+
+	Ghost<3,double> g(0.0015);
+
+	dec.setGoodParameters(box, bc, g, 512);
+
+	dec.decompose();
+
+	// Now we check the point
+
+	Point<3,double> p1({-0.0067499999999999999237,-0.0012499999999999995923,0.001250000000000000026});
+	Point<3,double> p2({-0.0067499999999999999237,-0.0012499999999999993755,0.001250000000000000026});
+
+	size_t proc1 = dec.processorIDBC(p1);
+	size_t proc2 = dec.processorIDBC(p2);
+
+	const openfpm::vector<std::pair<size_t, size_t>> & vp_id1 = dec.template ghost_processorID_pair<typename CartDecomposition<3, double>::lc_processor_id, typename CartDecomposition<3, double>::shift_id>(p1, UNIQUE);
+	const openfpm::vector<std::pair<size_t, size_t>> & vp_id2 = dec.template ghost_processorID_pair<typename CartDecomposition<3, double>::lc_processor_id, typename CartDecomposition<3, double>::shift_id>(p2, UNIQUE);
+
+	BOOST_REQUIRE_EQUAL(proc1,1);
+	BOOST_REQUIRE_EQUAL(proc2,2);
+
+	if (vcl.rank() == 2)
+	{
+		BOOST_REQUIRE(vp_id2.size() != 0);
+		BOOST_REQUIRE(vp_id1.size() == 0);
+	}
+
+	if (vcl.rank() == 1)
+	{
+		BOOST_REQUIRE(vp_id2.size() == 0 );
+		BOOST_REQUIRE(vp_id1.size() != 0 );
+	}
+}
 
 BOOST_AUTO_TEST_CASE( CartDecomposition_non_periodic_test_dist_grid)
 {
diff --git a/src/Grid/grid_dist_id_comm.hpp b/src/Grid/grid_dist_id_comm.hpp
index 2b894424..5533a3ca 100644
--- a/src/Grid/grid_dist_id_comm.hpp
+++ b/src/Grid/grid_dist_id_comm.hpp
@@ -328,7 +328,7 @@ class grid_dist_id_comm
 				// create 2 sub grid iterator
 
 				if (bx_dst.isValid() == false)
-					continue;
+				{continue;}
 
 				grid_key_dx_iterator_sub<dim> sub_src(loc_grid.get(i).getGrid(),bx_src.getKP1(),bx_src.getKP2());
 				grid_key_dx_iterator_sub<dim> sub_dst(loc_grid.get(sub_id_dst).getGrid(),bx_dst.getKP1(),bx_dst.getKP2());
diff --git a/src/Makefile.am b/src/Makefile.am
index 94fb6991..cb497146 100644
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -1,4 +1,4 @@
-LINKLIBS = $(HDF5_LDFLAGS)  $(HDF5_LIBS) $(OPENMP_LDFLAGS) $(LIBHILBERT_LIB)  $(METIS_LIB) $(PTHREAD_LIBS) $(OPT_LIBS) $(BOOST_LDFLAGS) $(BOOST_IOSTREAMS_LIB) $(CUDA_LIBS) $(PETSC_LIB) $(PARMETIS_LIB) $(BOOST_UNIT_TEST_FRAMEWORK_LIB) $(BOOST_CHRONO_LIB) $(BOOST_TIMER_LIB) $(BOOST_SYSTEM_LIB) $(LIBIFCORE)
+LINKLIBS = $(HDF5_LDFLAGS)  $(HDF5_LIBS) $(OPENMP_LDFLAGS) $(LIBHILBERT_LIB)  $(METIS_LIB) $(PTHREAD_LIBS) $(OPT_LIBS) $(BOOST_LDFLAGS) $(BOOST_IOSTREAMS_LIB) $(CUDA_LIBS) $(PETSC_LIB) $(SUITESPARSE_LIBS) $(LAPACK_LIBS) $(BLAS_LIBS)  $(PARMETIS_LIB) $(BOOST_UNIT_TEST_FRAMEWORK_LIB) $(BOOST_CHRONO_LIB) $(BOOST_TIMER_LIB) $(BOOST_SYSTEM_LIB) $(LIBIFCORE)
 
 noinst_PROGRAMS = pdata
 pdata_SOURCES = main.cpp Grid/tests/grid_dist_id_HDF5_chckpnt_restart_test.cpp Grid/tests/grid_dist_id_unit_test.cpp Grid/tests/staggered_grid_dist_unit_test.cpp Vector/tests/vector_dist_cell_list_tests.cpp Vector/tests/vector_dist_complex_prp_unit_test.cpp Vector/tests/vector_dist_HDF5_chckpnt_restart_test.cpp Vector/tests/vector_dist_MP_unit_tests.cpp Vector/tests/vector_dist_NN_tests.cpp Vector/tests/vector_dist_unit_test.cpp  pdata_performance.cpp Decomposition/tests/CartDecomposition_unit_test.cpp Decomposition/tests/shift_vect_converter_tests.cpp Vector/performance/vector_dist_performance_util.cpp  lib/pdata.cpp test_multiple_o.cpp ../openfpm_devices/src/memory/HeapMemory.cpp ../openfpm_devices/src/memory/PtrMemory.cpp ../openfpm_vcluster/src/VCluster/VCluster.cpp ../openfpm_devices/src/Memleak_check.cpp
-- 
GitLab