Dec optimizer work with disconnected sub-domains

src/SubdomainGraphNodes.hpp

@@ -34,7 +34,7 @@
 struct nm_v
 {
 	//! The node contain 3 unsigned long integer for communication computation memory and id
-	typedef boost::fusion::vector<float[3], size_t, size_t, size_t, size_t, size_t, size_t> type;
+	typedef boost::fusion::vector<float[3], size_t, size_t, size_t, size_t, long int, size_t> type;
 
 	//! type of the positional field
 	typedef float s_type;

src/dec_optimizer.hpp

@@ -207,7 +207,7 @@ private:
 	 *
 	 */
 
-	template<unsigned int p_sub> void fill_domain(Graph & graph,const Box<dim,size_t> & box, size_t ids)
+	template<unsigned int p_sub> void fill_domain(Graph & graph,const Box<dim,size_t> & box, long int ids)
 	{
 		// Create a subgrid iterator
 		grid_key_dx_iterator_sub<dim,do_not_print_warning_on_adjustment<dim>> g_sub(gh,box.getKP1(),box.getKP2());
@@ -327,8 +327,8 @@ private:
 	 * starting from one initial sub-domain find the biggest hyper-cube
 	 * output the box, and fill a list of neighborhood processor
 	 *
-	 * \tparam j id of the property storing the sub-decomposition
-	 * \tparam i id of the property containing the decomposition
+	 * \tparam p_sub id of the property storing the sub-decomposition
+	 * \tparam p_id id of the property containing the decomposition
 	 *
 	 * \param start_p initial domain
 	 * \param graph representing the grid of sub-sub-domain
@@ -500,7 +500,7 @@ private:
 	 *
 	 */
 
-	template<unsigned int p_id> grid_key_dx<dim> search_first_seed(Graph & graph, long int id)
+	template<unsigned int p_id, unsigned int p_sub> grid_key_dx<dim> search_seed(Graph & graph, long int id)
 	{
 		// if no processor is selected return the first point
 		if (id < -1)
@@ -522,7 +522,7 @@ private:
 			const grid_key_dx<dim> & gk = g_sub.get();
 
 			// if the subdomain has the id we are searching stop
-			if ((long int)graph.vertex(gh.LinId(gk)).template get<p_id>() == id)
+			if ((long int)graph.vertex(gh.LinId(gk)).template get<p_id>() == id && graph.vertex(gh.LinId(gk)).template get<p_sub>() == -1)
 			{
 				return gk;
 			}
@@ -537,69 +537,6 @@ private:
 		return key;
 	}
 
-public:
-
-	/*! \brief Constructor
-	 *
-	 * \param g Graph to simplify
-	 * \param sz size of the grid on each dimension
-	 *
-	 */
-
-	dec_optimizer(Graph & g, const size_t (& sz)[dim])
-	:gh(sz)
-	{
-		// The graph g is suppose to represent a cartesian grid
-		// No check is performed on g
-	}
-
-	/*! \brief optimize the graph
-	 *
-	 * Starting from a sub-sub-domain, it create wavefronts at the boundary and expand
-	 * the boundary until the wavefronts cannot expand any more, creating a sub-domain covering more sub-sub-domain.
-	 * This procedure continue until all the domain is covered by a sub-domains
-	 *
-	 * \tparam j property containing the processor decomposition
-	 * \tparam i property to fill with the sub-domain-decomposition id
-	 *
-	 * \param start_p seed point
-	 * \param graph we are processing
-	 *
-	 */
-	template <unsigned int p_sub, unsigned int p_id> void optimize(grid_key_dx<dim> & start_p, Graph & graph, const size_t (& bc)[dim])
-	{
-		// temporal vector
-		openfpm::vector<Box<dim,size_t>> tmp;
-
-		// temporal vector
-		openfpm::vector< openfpm::vector<size_t> > box_nn_processor;
-
-		// optimize
-		optimize<p_sub,p_id>(start_p,graph,-1,tmp, box_nn_processor,bc);
-	}
-
-	/*! \brief optimize the graph
-	 *
-	 * Starting from a sub-sub-domain, it create wavefronts at the boundary and expand
-	 * the boundary until the wavefronts cannot expand any more, creating a sub-domain covering more sub-sub-domain.
-	 * This procedure continue until all the sub-domain of the processor p_id are covered by a sub-domains
-	 *
-	 * \tparam j property containing the decomposition
-	 * \tparam i property to fill with the sub-domain-decomposition id
-	 *
-	 * \param graph we are processing
-	 * \param p_id Processor id (if p_id == -1 the optimization is done for all the processors)
-	 * \param list of sub-domain boxes
-	 *
-	 */
-	template <unsigned int p_sub, unsigned int p_id> void optimize(Graph & graph, long int pr_id, openfpm::vector<Box<dim,size_t>> & lb, openfpm::vector< openfpm::vector<size_t> > & box_nn_processor, const size_t (& bc)[dim])
-	{
-		// search for the first seed
-		grid_key_dx<dim> key_seed = search_first_seed<p_id>(graph,pr_id);
-
-		// optimize
-		optimize<p_sub,p_id>(key_seed,graph,pr_id,lb,box_nn_processor,bc);
-	}
 
 	/*! \brief optimize the graph
 	 *
@@ -617,13 +554,14 @@ public:
 	 * \param list of sub-domain boxes produced by the algorithm
 	 * \param box_nn_processor for each box it list all the neighborhood processor
 	 * \param bc Boundary condition
+	 * \param init_sub_id when true p_sub property is initial set to -1 [default true]
+	 * \param sub_id starting sub_id enumeration [default 0]
+	 *
+	 * \return last assigned sub-id
 	 *
 	 */
-	template <unsigned int p_sub, unsigned int p_id> void optimize(grid_key_dx<dim> & start_p, Graph & graph, long int pr_id, openfpm::vector<Box<dim,size_t>> & lb, openfpm::vector< openfpm::vector<size_t> > & box_nn_processor, const size_t (& bc)[dim])
+	template <unsigned int p_sub, unsigned int p_id> size_t optimize(grid_key_dx<dim> & start_p, Graph & graph, long int pr_id, openfpm::vector<Box<dim,size_t>> & lb, openfpm::vector< openfpm::vector<size_t> > & box_nn_processor, const size_t (& bc)[dim], bool init_sub_id = true, size_t sub_id = 0)
 	{
-		// sub-domain id
-		size_t sub_id =  0;
-
 		// queue
 		openfpm::vector<size_t> v_q;
 
@@ -637,7 +575,9 @@ public:
 		openfpm::vector<wavefront<dim>> v_w(w_comb.size());
 
 		// fill the sub decomposition with negative number
-		fill_domain<p_sub>(graph,gh.getBox(),-1);
+
+		if (init_sub_id == true)
+			fill_domain<p_sub>(graph,gh.getBox(),-1);
 
 		// push the first domain
 		v_q.add(gh.LinId(start_p));
@@ -671,6 +611,84 @@ public:
 			// increment the sub_id
 			sub_id++;
 		}
+
+		return sub_id;
+	}
+
+public:
+
+	/*! \brief Constructor
+	 *
+	 * \param g Graph to simplify
+	 * \param sz size of the grid on each dimension
+	 *
+	 */
+
+	dec_optimizer(Graph & g, const size_t (& sz)[dim])
+	:gh(sz)
+	{
+		// The graph g is suppose to represent a cartesian grid
+		// No check is performed on g
+	}
+
+	/*! \brief optimize the graph
+	 *
+	 * Starting from a sub-sub-domain, it create wavefronts at the boundary and expand
+	 * the boundary until the wavefronts cannot expand any more, creating a sub-domain covering more sub-sub-domain.
+	 * This procedure continue until all the domain is covered by a sub-domains
+	 *
+	 * \tparam j property containing the processor decomposition
+	 * \tparam i property to fill with the sub-domain-decomposition id
+	 *
+	 * \param start_p seed point
+	 * \param graph we are processing
+	 *
+	 */
+	template <unsigned int p_sub, unsigned int p_id> void optimize(grid_key_dx<dim> & start_p, Graph & graph, const size_t (& bc)[dim])
+	{
+		// temporal vector
+		openfpm::vector<Box<dim,size_t>> tmp;
+
+		// temporal vector
+		openfpm::vector< openfpm::vector<size_t> > box_nn_processor;
+
+		// optimize
+		optimize<p_sub,p_id>(start_p,graph,-1,tmp, box_nn_processor,bc);
+	}
+
+	/*! \brief optimize the graph
+	 *
+	 * Starting from a sub-sub-domain, it create wavefronts at the boundary and expand
+	 * the boundary until the wavefronts cannot expand any more, creating a sub-domain covering more sub-sub-domain.
+	 * This procedure continue until all the sub-domain of the processor p_id are covered by a sub-domains
+	 *
+	 * \tparam j property containing the decomposition
+	 * \tparam i property to fill with the sub-domain-decomposition id
+	 *
+	 * \param graph we are processing
+	 * \param p_id Processor id (if p_id == -1 the optimization is done for all the processors)
+	 * \param list of sub-domain boxes
+	 *
+	 */
+	template <unsigned int p_sub, unsigned int p_id> void optimize(Graph & graph, long int pr_id, openfpm::vector<Box<dim,size_t>> & lb, openfpm::vector< openfpm::vector<size_t> > & box_nn_processor, const size_t (& bc)[dim])
+	{
+		size_t sub_id = 0;
+
+		// fill the sub decomposition with negative number
+		fill_domain<p_sub>(graph,gh.getBox(),-1);
+
+		grid_key_dx<dim> key_seed = search_seed<p_id,p_sub>(graph,pr_id);
+
+		while (key_seed.isValid())
+		{
+			// optimize
+			sub_id = optimize<p_sub,p_id>(key_seed,graph,pr_id,lb,box_nn_processor,bc,false,sub_id);
+
+			// new seed
+			key_seed = search_seed<p_id,p_sub>(graph,pr_id);
+
+			std::cerr << "Key seed " << key_seed.to_string() << "\n";
+		}
 	}
 };
 

src/dec_optimizer_unit_test.hpp

@@ -13,7 +13,7 @@
 #include "Graph/map_graph.hpp"
 #include "Decomposition/Distribution/metis_util.hpp"
 #include "dec_optimizer.hpp"
-
+#include "util/SimpleRNG.hpp"
 
 #undef GS_SIZE
 #define GS_SIZE 8
@@ -111,12 +111,8 @@ BOOST_AUTO_TEST_CASE( dec_optimizer_test_use_p)
 	// For each sub-domain check the neighborhood processors
 	openfpm::vector< openfpm::vector<size_t> > box_nn_processor;
 
-	// key
-	grid_key_dx<3> zero;
-	zero.zero();
-
 	// gp,p_id,loc_box,box_nn_processor,bc
-	d_o.optimize<nm_v::sub_id,nm_v::id>(zero,g,-1,dec_o,box_nn_processor,bc);
+	d_o.optimize<nm_v::sub_id,nm_v::id>(g,-1,dec_o,box_nn_processor,bc);
 
 	BOOST_REQUIRE_EQUAL(box_nn_processor.size(),8ul);
 
@@ -145,6 +141,75 @@ BOOST_AUTO_TEST_CASE( dec_optimizer_test_use_p)
 	// check
 }
 
+BOOST_AUTO_TEST_CASE( dec_optimizer_disconnected_subdomains_np)
+{
+	// Vcluster
+	Vcluster & vcl = *global_v_cluster;
+
+	CartesianGraphFactory<2,Graph_CSR<nm_v,nm_e>> g_factory;
+
+	// Cartesian grid
+	size_t sz[2] = {GS_SIZE,GS_SIZE};
+
+	// Box
+	Box<2,float> box({0.0,0.0},{1.0,1.0});
+
+	// Boundary conditions, non periodic
+	size_t bc[] = {NON_PERIODIC,NON_PERIODIC};
+
+	// Graph to decompose
+	Graph_CSR<nm_v,nm_e> g = g_factory.construct<nm_e::communication,NO_VERTEX_ID,float,1,0,1>(sz,box,bc);
+
+	SimpleRNG rng;
+
+	auto vit = g.getVertexIterator();
+
+	while (vit.isNext())
+	{
+		auto vk = vit.get();
+
+		g.template vertex_p<nm_v::proc_id>(vk) = rng.GetUniform() * 2.9999;
+		g.template vertex_p<nm_v::sub_id>(vk) = 100;
+
+		++vit;
+	}
+
+	// optimize
+	dec_optimizer<2,Graph_CSR<nm_v,nm_e>> d_o(g,sz);
+
+	// set of Boxes produced by the decomposition optimizer
+	openfpm::vector<::Box<2, size_t>> loc_box;
+
+	//! for each sub-domain, contain the list of the neighborhood processors
+	openfpm::vector<openfpm::vector<long unsigned int> > box_nn_processor;
+
+	d_o.optimize<nm_v::sub_id, nm_v::proc_id>(g, vcl.getProcessUnitID(), loc_box, box_nn_processor,bc);
+
+	std::stringstream str_g;
+	str_g << "dec_optimizer_disc_graph" << vcl.getProcessUnitID() << ".vtk";
+	std::stringstream str_gt;
+	str_gt << "dec_optimizer_disc_graph" << vcl.getProcessUnitID() << "_test.vtk";
+
+	std::stringstream str_s;
+	str_s << "dec_optimizer_disc_sub" << vcl.getProcessUnitID() << ".vtk";
+	std::stringstream str_st;
+	str_st << "dec_optimizer_disc_sub" << vcl.getProcessUnitID() << "_test.vtk";
+
+
+	VTKWriter<Graph_CSR<nm_v,nm_e>,VTK_GRAPH> wrt(g);
+	wrt.write("dec_optimizer_disc_graph" + std::to_string(vcl.getProcessUnitID()) + ".vtk");
+
+	VTKWriter<openfpm::vector<::Box<2, size_t>>, VECTOR_BOX> vtk_box1;
+	vtk_box1.add(loc_box);
+	vtk_box1.write("dec_optimizer_disc_sub" + std::to_string(vcl.getProcessUnitID()) + std::string(".vtk"));
+
+	bool test = compare(str_g.str(), str_gt.str());
+	BOOST_REQUIRE_EQUAL(true,test);
+
+	test = compare(str_s.str(),str_st.str());
+	BOOST_REQUIRE_EQUAL(true,test);
+}
+
 BOOST_AUTO_TEST_SUITE_END()