improved performance for grid

example/Grid/3_gray_scott_3d_vectorization/main.cpp

@@ -36,7 +36,7 @@
 
 //! \cond [constants] \endcond
 
-//#define FORTRAN_UPDATE
+#define FORTRAN_UPDATE
 
 constexpr int x = 0;
 constexpr int y = 1;
@@ -269,8 +269,50 @@ int main(int argc, char* argv[])
 	// New grid with the decomposition of the old grid
     grid_dist_id<3, double, aggregate<double>> NewU(OldU.getDecomposition(),sz,g);
     grid_dist_id<3, double, aggregate<double>> NewV(OldV.getDecomposition(),sz,g);
-	
+
+	//////////// DEBUG /////////////////////////
+	//
+	//
+	//
+
+    	auto debug_it = OldU.getDomainIterator();
+
+	int count_ = 0;
+
+	while (debug_it.isNext())
+	{
+		auto key = debug_it.get();
+
+		count_++;
+
+		++debug_it;
+	}
+
+        auto debug_it2 = OldU.getDomainGhostIterator();
+
+        int count_dg = 0;
+
+        while (debug_it2.isNext())
+        {       
+                auto key = debug_it2.get();
+
+                count_dg++;
+
+                ++debug_it2;
+        }
+
+	auto & v_cl = create_vcluster();
+
+	count_dg -= count_;
+	v_cl.sum(count_dg);
+	v_cl.execute();
+
+	std::cout << "Ghost points: " << count_dg << std::endl;
+
+	//////////////////////////////////////////////
+
 	// spacing of the grid on x and y
+
 	double spacing[3] = {OldU.spacing(0),OldU.spacing(1),OldU.spacing(2)};
 
 	init(OldU,OldV,NewU,NewV,domain);
@@ -298,7 +340,7 @@ int main(int argc, char* argv[])
 						    {-1,0,0},
 						    {1,0,0}};
 
-	for (size_t i = 0; i < timeSteps; ++i)
+	for (size_t i = 0; i < 1 /*timeSteps*/; ++i)
 	{
 		if (i % 300 == 0)
 			std::cout << "STEP: " << i << std::endl;
@@ -329,9 +371,11 @@ int main(int argc, char* argv[])
 		{
 			step(OldU,OldV,NewU,NewV,star_stencil_3D,uFactor,vFactor,deltaT,F,K);
 
+			for (size_t i = 0 ; i < 10000 ; i++)
+			{
 			NewU.ghost_get<0>();
 			NewV.ghost_get<0>();
-
+			}
 		}
 		else
 		{
@@ -355,9 +399,9 @@ int main(int argc, char* argv[])
 		//! \cond [save hdf5] \endcond
 
 		// Every 500 time step we output the configuration on hdf5
-		if (i % 500 == 0)
+		if (i % 2000 == 0)
 		{
-			OldU.save("output_u_" + std::to_string(count));
+//			OldU.save("output_u_" + std::to_string(count));
 			OldV.save("output_v_" + std::to_string(count));
 			count++;
 		}
@@ -368,6 +412,9 @@ int main(int argc, char* argv[])
 	tot_sim.stop();
 	std::cout << "Total simulation: " << tot_sim.getwct() << std::endl;
 
+	// We frite the final configuration
+	OldV.write("final");
+
 	//! \cond [time stepping] \endcond
 
 	/*!

example/Numerics/Vortex_in_cell/Makefile

-include ../../../example.mk
+include ../../example.mk
 
 CC=mpic++
 

example/VCluster/0_simple/main.cpp

@@ -7,6 +7,7 @@
  *
  * \subpage VCluster_0_simple
  * \subpage VCluster_1_semantic
+ * \subpage VCluster_2_serial_and_parallel
  *
  */
 

images/vector.cpp

@@ -34,6 +34,11 @@ public:
 
 #endif
 
+	static inline bool noPointers()
+	{
+		return true;
+	}
+
 };
 
 int main(int argc, char* argv[])

openfpm_data @ 8a861515

-Subproject commit 4fecaeb91f2651821aa8ec5ebbd9932b937d6bc8
+Subproject commit 8a861515b74d2f9ca063e045818c153615eb07fd

src/Decomposition/ORB.hpp

@@ -8,7 +8,6 @@
 #ifndef ORB_HPP_
 #define ORB_HPP_
 
-#include "data_type/scalar.hpp"
 #include "util/mathutil.hpp"
 
 /*! \brief this class is a functor for "for_each" algorithm
@@ -77,7 +76,7 @@ struct do_when_dim_gr_i<dim,i,ORB,typename boost::enable_if< boost::mpl::bool_<(
  *
  */
 
-template<typename T> class ORB_node : public scalar<T>
+template<typename T> class ORB_node : public aggregate<T>
 {
 public:
 

src/Grid/Iterators/grid_dist_id_iterators_unit_tests.hpp

@@ -46,7 +46,7 @@ void Test2D_sub(const Box<2,float> & domain, long int k)
 		Ghost<2,float> g(0.01 / factor);
 
 		// Distributed grid with id decomposition
-		grid_dist_id<2, float, scalar<float>> g_dist(sz,domain,g);
+		grid_dist_id<2, float, aggregate<float>> g_dist(sz,domain,g);
 
 		// check the consistency of the decomposition
 		bool val = g_dist.getDecomposition().check_consistency();

src/Grid/grid_dist_id.hpp

@@ -74,6 +74,9 @@ class grid_dist_id : public grid_dist_id_comm<dim,St,T,Decomposition,Memory,devi
 	//! Ghost expansion
 	Ghost<dim,St> ghost;
 
+	//! Ghost expansion
+	Ghost<dim,long int> ghost_int;
+
 	//! Local grids
 	mutable openfpm::vector<device_grid> loc_grid;
 
@@ -192,6 +195,53 @@ class grid_dist_id : public grid_dist_id_comm<dim,St,T,Decomposition,Memory,devi
 		return flp;
 	}
 
+	/*! \brief this function is for optimization of the ghost size
+	 *
+	 * Because the decomposition work in continuum and discrete ghost is
+	 *  converted in continuum, in some case continuum ghost because of
+	 *  rounding-off error can produce ghost bigger than the discrete selected
+	 *   one. This function adjust for this round-off error
+	 *
+	 * \param sub_domain the sub-domain
+	 * \param sub_domain_other the other sub-domain
+	 * \param ib internal ghost box to adjust
+	 *
+	 */
+	void set_for_adjustment(const Box<dim,long int> & sub_domain,
+							const Box<dim,St> & sub_domain_other,
+							const comb<dim> & cmb,
+							Box<dim,long int> & ib,
+							Ghost<dim,long int> & g)
+	{
+		if (g.isInvalidGhost() == true)
+		{return;}
+
+		// Convert from SpaceBox<dim,St> to SpaceBox<dim,long int>
+		Box<dim,long int> sub_domain_other_exp = cd_sm.convertDomainSpaceIntoGridUnits(sub_domain_other,dec.periodicity());
+
+		// translate sub_domain_other based on cmb
+		for (size_t i = 0 ; i < dim ; i++)
+		{
+			if (cmb.c[i] == 1)
+			{
+				sub_domain_other_exp.setLow(i,sub_domain_other_exp.getLow(i) - ginfo.size(i));
+				sub_domain_other_exp.setHigh(i,sub_domain_other_exp.getHigh(i) - ginfo.size(i));
+			}
+			else if (cmb.c[i] == -1)
+			{
+				sub_domain_other_exp.setLow(i,sub_domain_other_exp.getLow(i) + ginfo.size(i));
+				sub_domain_other_exp.setHigh(i,sub_domain_other_exp.getHigh(i) + ginfo.size(i));
+			}
+		}
+
+		sub_domain_other_exp.enlarge(g);
+		if (sub_domain_other_exp.Intersect(sub_domain,ib) == false)
+		{
+			for (size_t i = 0 ; i < dim ; i++)
+			{ib.setHigh(i,ib.getLow(i) - 1);}
+		}
+	}
+
 	/*! \brief Create per-processor internal ghost boxes list in grid units and g_id_to_external_ghost_box
 	 *
 	 */
@@ -221,6 +271,21 @@ class grid_dist_id : public grid_dist_id_comm<dim,St,T,Decomposition,Memory,devi
 				if (ib.isValid() == false)
 					continue;
 
+				size_t sub_id = dec.getProcessorIGhostSub(i,j);
+				size_t r_sub = dec.getProcessorIGhostSSub(i,j);
+
+				auto & n_box = dec.getNearSubdomains(dec.IDtoProc(i));
+
+				Box<dim,long int> sub = gdb_ext.get(sub_id).Dbox;
+				sub += gdb_ext.get(sub_id).origin;
+
+				set_for_adjustment(sub,
+						           n_box.get(r_sub),dec.getProcessorIGhostPos(i,j),
+								   ib,ghost_int);
+
+				if (ib.isValid() == false)
+					continue;
+
 				// save the box and the sub-domain id (it is calculated as the linearization of P1)
 				::Box<dim,size_t> cvt = ib;
 
@@ -337,6 +402,19 @@ class grid_dist_id : public grid_dist_id_comm<dim,St,T,Decomposition,Memory,devi
 				::Box<dim,St> ib_dom = dec.getLocalIGhostBox(i,j);
 				::Box<dim,long int> ib = cd_sm.convertDomainSpaceIntoGridUnits(ib_dom,dec.periodicity());
 
+				// Check if ib is valid if not it mean that the internal ghost does not contain information so skip it
+				if (ib.isValid() == false)
+					continue;
+
+				size_t sub_id = i;
+				size_t r_sub = dec.getLocalIGhostSub(i,j);
+
+				Box<dim,long int> sub = gdb_ext.get(sub_id).Dbox;
+				sub += gdb_ext.get(sub_id).origin;
+
+				set_for_adjustment(sub,dec.getSubDomain(r_sub),
+						           dec.getLocalIGhostPos(i,j),ib,ghost_int);
+
 				// Check if ib is valid if not it mean that the internal ghost does not contain information so skip it
 				if (ib.isValid() == false)
 					continue;
@@ -653,8 +731,11 @@ public:
 	 * \param ext extension of the grid (must be positive on every direction)
 	 *
 	 */
-	template<typename H> grid_dist_id(const grid_dist_id<dim,St,H,typename Decomposition::base_type,Memory,grid_cpu<dim,H>> & g, const Ghost<dim,long int> & gh, Box<dim,size_t> ext)
-	:dec(create_vcluster()),v_cl(create_vcluster())
+	template<typename H>
+	grid_dist_id(const grid_dist_id<dim,St,H,typename Decomposition::base_type,Memory,grid_cpu<dim,H>> & g,
+			     const Ghost<dim,long int> & gh,
+				 Box<dim,size_t> ext)
+	:ghost_int(gh),dec(create_vcluster()),v_cl(create_vcluster())
 	{
 #ifdef SE_CLASS2
 		check_new(this,8,GRID_DIST_EVENT,4);
@@ -701,8 +782,11 @@ public:
      * \param ghost Ghost part
      *
      */
-    grid_dist_id(const Decomposition & dec, const size_t (& g_sz)[dim], const Ghost<dim,St> & ghost)
-    :domain(dec.getDomain()),ghost(ghost),dec(dec),v_cl(create_vcluster()),ginfo(g_sz),ginfo_v(g_sz)
+    grid_dist_id(const Decomposition & dec,
+    		     const size_t (& g_sz)[dim],
+				 const Ghost<dim,St> & ghost)
+    :domain(dec.getDomain()),ghost(ghost),ghost_int(INVALID_GHOST),dec(dec),v_cl(create_vcluster()),
+	 ginfo(g_sz),ginfo_v(g_sz)
 	{
 #ifdef SE_CLASS2
 		check_new(this,8,GRID_DIST_EVENT,4);
@@ -719,8 +803,10 @@ public:
      * \param ghost Ghost part
      *
      */
-    grid_dist_id(Decomposition && dec, const size_t (& g_sz)[dim], const Ghost<dim,St> & ghost)
-    :domain(dec.getDomain()),ghost(ghost),dec(dec),ginfo(g_sz),ginfo_v(g_sz),v_cl(create_vcluster())
+    grid_dist_id(Decomposition && dec, const size_t (& g_sz)[dim],
+    		     const Ghost<dim,St> & ghost)
+    :domain(dec.getDomain()),ghost(ghost),dec(dec),ginfo(g_sz),
+	 ginfo_v(g_sz),v_cl(create_vcluster()),ghost_int(INVALID_GHOST)
 	{
 #ifdef SE_CLASS2
 		check_new(this,8,GRID_DIST_EVENT,4);
@@ -739,8 +825,10 @@ public:
      * \warning In very rare case the ghost part can be one point bigger than the one specified
      *
      */
-	grid_dist_id(const Decomposition & dec, const size_t (& g_sz)[dim], const Ghost<dim,long int> & g)
-	:domain(dec.getDomain()),dec(create_vcluster()),v_cl(create_vcluster()),ginfo(g_sz),ginfo_v(g_sz)
+	grid_dist_id(const Decomposition & dec, const size_t (& g_sz)[dim],
+			     const Ghost<dim,long int> & g)
+	:domain(dec.getDomain()),ghost_int(g),dec(create_vcluster()),v_cl(create_vcluster()),
+	 ginfo(g_sz),ginfo_v(g_sz)
 	{
 #ifdef SE_CLASS2
 		check_new(this,8,GRID_DIST_EVENT,4);
@@ -764,8 +852,10 @@ public:
      * \warning In very rare case the ghost part can be one point bigger than the one specified
      *
      */
-	grid_dist_id(Decomposition && dec, const size_t (& g_sz)[dim], const Ghost<dim,long int> & g)
-	:domain(dec.getDomain()),dec(dec),v_cl(create_vcluster()),ginfo(g_sz),ginfo_v(g_sz)
+	grid_dist_id(Decomposition && dec, const size_t (& g_sz)[dim],
+			     const Ghost<dim,long int> & g)
+	:domain(dec.getDomain()),dec(dec),v_cl(create_vcluster()),ginfo(g_sz),
+	 ginfo_v(g_sz),ghost_int(g)
 	{
 #ifdef SE_CLASS2
 		check_new(this,8,GRID_DIST_EVENT,4);
@@ -787,7 +877,8 @@ public:
      * \warning In very rare case the ghost part can be one point bigger than the one specified
      *
      */
-	grid_dist_id(const size_t (& g_sz)[dim],const Box<dim,St> & domain, const Ghost<dim,St> & g)
+	grid_dist_id(const size_t (& g_sz)[dim],const Box<dim,St> & domain,
+			     const Ghost<dim,St> & g)
 	:grid_dist_id(g_sz,domain,g,create_non_periodic<dim>())
 	{
 	}
@@ -816,8 +907,10 @@ public:
      * \warning In very rare case the ghost part can be one point bigger than the one specified
      *
      */
-	grid_dist_id(const size_t (& g_sz)[dim],const Box<dim,St> & domain, const Ghost<dim,St> & g, const periodicity<dim> & p)
-	:domain(domain),ghost(g),dec(create_vcluster()),v_cl(create_vcluster()),ginfo(g_sz),ginfo_v(g_sz)
+	grid_dist_id(const size_t (& g_sz)[dim],const Box<dim,St> & domain,
+			     const Ghost<dim,St> & g, const periodicity<dim> & p)
+	:domain(domain),ghost(g),ghost_int(INVALID_GHOST),dec(create_vcluster()),v_cl(create_vcluster()),
+	 ginfo(g_sz),ginfo_v(g_sz)
 	{
 #ifdef SE_CLASS2
 		check_new(this,8,GRID_DIST_EVENT,4);
@@ -838,8 +931,10 @@ public:
      * \warning In very rare case the ghost part can be one point bigger than the one specified
      *
      */
-	grid_dist_id(const size_t (& g_sz)[dim],const Box<dim,St> & domain, const Ghost<dim,long int> & g, const periodicity<dim> & p)
-	:domain(domain),dec(create_vcluster()),v_cl(create_vcluster()),ginfo(g_sz),ginfo_v(g_sz)
+	grid_dist_id(const size_t (& g_sz)[dim],const Box<dim,St> & domain,
+			     const Ghost<dim,long int> & g, const periodicity<dim> & p)
+	:domain(domain),ghost_int(g),dec(create_vcluster()),v_cl(create_vcluster()),ginfo(g_sz),
+	 ginfo_v(g_sz)
 	{
 #ifdef SE_CLASS2
 		check_new(this,8,GRID_DIST_EVENT,4);
@@ -1398,11 +1493,6 @@ public:
 																						  	  	  	 g_id_to_internal_ghost_box);
 	}
 
-	// copy bench test
-	double mem_mem_time = 0.0;
-	double mem_ite_time = 0.0;
-
-	int mem_select = 0;
 
 	/*! \brief Copy the give grid into this grid
 	 *
@@ -1713,6 +1803,26 @@ public:
 		map();
 	}
 
+	/*! \brief Get the internal local ghost box
+	 *
+	 * \return the internal local ghost box
+	 *
+	 */
+	const openfpm::vector<i_lbox_grid<dim>> & get_loc_ig_box()
+	{
+		return this->loc_ig_box;
+	}
+
+	/*! \brief Get the internal ghost box
+	 *
+	 * \return the internal local ghost box
+	 *
+	 */
+	const openfpm::vector<i_lbox_grid<dim>> & get_ig_box()
+	{
+		return this->ig_box;
+	}
+
 	//! Define friend classes
 	//\cond
 	friend grid_dist_id<dim,St,T,typename Decomposition::extended_type,Memory,device_grid>;

src/Grid/grid_dist_id_HDF5_chckpnt_restart_test.hpp

@@ -38,7 +38,7 @@ BOOST_AUTO_TEST_CASE( grid_dist_id_hdf5_save_test )
 	Ghost<2,float> g(ghost_part);
 
 	// Distributed grid with id decomposition
-	grid_dist_id<2, float, scalar<float>, CartDecomposition<2,float>> g_dist(sz,domain,g);
+	grid_dist_id<2, float, aggregate<float>, CartDecomposition<2,float>> g_dist(sz,domain,g);
 
 	// get the decomposition
 	auto & dec = g_dist.getDecomposition();
@@ -106,7 +106,7 @@ BOOST_AUTO_TEST_CASE( grid_dist_id_hdf5_load_test )
 	Ghost<2,float> g(ghost_part);
 
 	// Distributed grid with id decomposition
-	grid_dist_id<2, float, scalar<float>, CartDecomposition<2,float>> g_dist(sz,domain,g);
+	grid_dist_id<2, float, aggregate<float>, CartDecomposition<2,float>> g_dist(sz,domain,g);
 
 	g_dist.getDecomposition().write("Before_load_grid_decomposition");
 	g_dist.write("Before_Loaded_grid");

src/Grid/grid_dist_id_unit_test.cpp

@@ -6,7 +6,6 @@
 
 #include "Point_test.hpp"
 #include "grid_dist_id.hpp"
-#include "data_type/scalar.hpp"
 #include "data_type/aggregate.hpp"
 
 
@@ -56,7 +55,7 @@ BOOST_AUTO_TEST_CASE( grid_dist_id_domain_grid_unit_converter3D_test)
 		Ghost<3,float> g(0.01);
 
 		// Distributed grid with id decomposition
-		grid_dist_id<3, float, scalar<float>, CartDecomposition<3,float>> g_dist(sz,domain,g);
+		grid_dist_id<3, float, aggregate<float>, CartDecomposition<3,float>> g_dist(sz,domain,g);
 
 		// get the decomposition
 		auto & dec = g_dist.getDecomposition();
@@ -127,7 +126,7 @@ BOOST_AUTO_TEST_CASE( grid_dist_id_domain_grid_unit_converter_test)
 		Ghost<2,float> g(0.01);
 
 		// Distributed grid with id decomposition
-		grid_dist_id<2, float, scalar<float>, CartDecomposition<2,float>> g_dist(sz,domain,g);
+		grid_dist_id<2, float, aggregate<float>, CartDecomposition<2,float>> g_dist(sz,domain,g);
 
 		// get the decomposition
 		auto & dec = g_dist.getDecomposition();
@@ -187,7 +186,7 @@ void Test2D(const Box<2,float> & domain, long int k)
 		Ghost<2,float> g(0.01 / factor);
 
 		// Distributed grid with id decomposition
-		grid_dist_id<2, float, scalar<float>> g_dist(sz,domain,g);
+		grid_dist_id<2, float, aggregate<float>> g_dist(sz,domain,g);
 
 		// check the consistency of the decomposition
 		bool val = g_dist.getDecomposition().check_consistency();
@@ -307,7 +306,7 @@ void Test1D(const Box<1,float> & domain, long int k)
 		Ghost<1,float> g(0.01 / factor);
 
 		// Distributed grid with id decomposition
-		grid_dist_id<1, float, scalar<float>> g_dist(sz,domain,g);
+		grid_dist_id<1, float, aggregate<float>> g_dist(sz,domain,g);
 
 		// check the consistency of the decomposition
 		bool val = g_dist.getDecomposition().check_consistency();
@@ -425,7 +424,7 @@ void Test3D_sub(const Box<3,float> & domain, long int k)
 		Ghost<3,float> g(0.01 / factor);
 
 		// Distributed grid with id decomposition
-		grid_dist_id<3, float, scalar<float>, CartDecomposition<3,float>> g_dist(sz,domain,g);
+		grid_dist_id<3, float, aggregate<float>, CartDecomposition<3,float>> g_dist(sz,domain,g);
 
 		// check the consistency of the decomposition
 		bool val = g_dist.getDecomposition().check_consistency();
@@ -534,7 +533,7 @@ void Test3D(const Box<3,float> & domain, long int k)
 		Ghost<3,float> g(0.01 / factor);
 
 		// Distributed grid with id decomposition
-		grid_dist_id<3, float, scalar<float>, CartDecomposition<3,float>> g_dist(sz,domain,g);
+		grid_dist_id<3, float, aggregate<float>, CartDecomposition<3,float>> g_dist(sz,domain,g);
 
 		// check the consistency of the decomposition
 		bool val = g_dist.getDecomposition().check_consistency();
@@ -644,7 +643,7 @@ void Test3D_gg(const Box<3,float> & domain, long int k, long int gk)
 		Ghost<3,long int> g(gk);
 
 		// Distributed grid with id decomposition
-		grid_dist_id<3, float, scalar<float>, CartDecomposition<3,float>> g_dist(sz,domain,g);
+		grid_dist_id<3, float, aggregate<float>, CartDecomposition<3,float>> g_dist(sz,domain,g);
 
 		// check the consistency of the decomposition
 		bool val = g_dist.getDecomposition().check_consistency();
@@ -1590,6 +1589,89 @@ void Test_grid_copy(const Box<3,float> & domain, long int k)
 	}
 }
 
+void Test_ghost_correction(Box<3,double> & domain, long int k, long int g_)
+{
+	size_t sz[3] = {(size_t)k,(size_t)k,(size_t)k};
+	periodicity<3> bc = {PERIODIC,PERIODIC,PERIODIC};
+
+	Ghost<3,long int> g(g_);
+
+	grid_dist_id<3, double, aggregate<double>> grid(sz,domain,g,bc);
+
+    auto itg = grid.getDomainGhostIterator();
+
+    while (itg.isNext())
+    {
+    	auto key = itg.get();
+
+    	grid.template get<0>(key) = 0.0;
+
+    	++itg;
+    }
+
+	// Fill everything with 5
+
+    auto it = grid.getDomainIterator();
+
+    while (it.isNext())
+    {
+    	auto key = it.get();
+    	auto gkey = it.getGKey(key);
+
+    	if (gkey.get(0) == -4 && gkey.get(1) == 20 && gkey.get(2) == -4)
+    	{
+    		grid.template get<0>(key) = 20.0;
+    	}
+    	else
+    	{
+    		grid.template get<0>(key) = 5.0;
+    	}
+
+    	++it;
+    }
+
+    grid.ghost_get<0>();
+    auto it2 = grid.getDomainGhostIterator();
+
+    bool is_inside = true;
+
+    while (it2.isNext())
+    {
+    	auto key = it2.get();
+    	auto gkey = it2.getGKey(key);
+
+    	if (grid.template get<0>(key) == 5.0)
+    	{
+    		// Here we check that the point is with in one stencil point
+    		// from one sub-domain
+
+    		bool is_inside_point = false;
+    		for (size_t i = 0 ; i < grid.getN_loc_grid() ; i++)
+    		{
+    			Box<3,long int> bx = grid.getLocalGridsInfo().get(i).Dbox;
+    			bx += grid.getLocalGridsInfo().get(i).origin;
+
+    			bx.enlarge(g);
+
+    			if (bx.isInside(gkey.toPoint()) == true)
+    			{
+    				is_inside_point |= true;
+    			}
+    		}
+
+    		is_inside &= is_inside_point;
+    	}
+
+        ++it2;
+    }
+
+
+    grid.getDecomposition().write("dec_set_for_adj");
+    grid.write("dec_for_adj");
+
+    BOOST_REQUIRE_EQUAL(is_inside,true);
+}
+
 #include "grid_dist_id_unit_test_ext_dom.hpp"
 #include "grid_dist_id_unit_test_unb_ghost.hpp"
 
@@ -1756,6 +1838,39 @@ BOOST_AUTO_TEST_CASE( grid_dist_id_periodic_put_test )
 	Test3D_periodic_put(domain3,k);
 }
 
+BOOST_AUTO_TEST_CASE ( grid_ghost_correction )
+{
+	Box<3,double> domain({0.0,0.0,0.0},{2.5,2.5,2.5});
+
+	long int k = 128;
+
+	Test_ghost_correction(domain,k,1);
+	Test_ghost_correction(domain,k,2);
+	Test_ghost_correction(domain,k,3);
+	Test_ghost_correction(domain,k,4);
+
+	k = 64;
+
+	Test_ghost_correction(domain,k,1);
+	Test_ghost_correction(domain,k,2);