grid_dist_id.hpp 75.2 KB
Newer Older
incardon's avatar
incardon committed
1 2 3 4
#ifndef COM_UNIT_HPP
#define COM_UNIT_HPP

#include <vector>
5
#include <unordered_map>
incardon's avatar
incardon committed
6
#include "Grid/map_grid.hpp"
incardon's avatar
incardon committed
7
#include "VCluster/VCluster.hpp"
incardon's avatar
incardon committed
8
#include "Space/SpaceBox.hpp"
9
#include "util/mathutil.hpp"
incardon's avatar
incardon committed
10 11 12
#include "Iterators/grid_dist_id_iterator_dec.hpp"
#include "Iterators/grid_dist_id_iterator.hpp"
#include "Iterators/grid_dist_id_iterator_sub.hpp"
incardon's avatar
incardon committed
13
#include "grid_dist_key.hpp"
14
#include "NN/CellList/CellDecomposer.hpp"
15 16
#include "util/object_util.hpp"
#include "memory/ExtPreAlloc.hpp"
17
#include "VTKWriter/VTKWriter.hpp"
18 19
#include "Packer_Unpacker/Packer.hpp"
#include "Packer_Unpacker/Unpacker.hpp"
20
#include "Decomposition/CartDecomposition.hpp"
21
#include "data_type/aggregate.hpp"
22
#include "hdf5.h"
Yaroslav's avatar
Yaroslav committed
23
#include "grid_dist_id_comm.hpp"
incardon's avatar
Working  
incardon committed
24
#include "HDF5_wr/HDF5_wr.hpp"
incardon's avatar
incardon committed
25
#include "SparseGrid/SparseGrid.hpp"
26 27 28 29
#ifdef __NVCC__
#include "SparseGridGpu/SparseGridGpu.hpp"
#include "cuda/grid_dist_id_kernels.cuh"
#endif
incardon's avatar
incardon committed
30

incardon's avatar
incardon committed
31

Pietro Incardona's avatar
Pietro Incardona committed
32
//! Internal ghost box sent to construct external ghost box into the other processors
Pietro Incardona's avatar
Pietro Incardona committed
33 34 35
template<unsigned int dim>
struct Box_fix
{
Pietro Incardona's avatar
Pietro Incardona committed
36
	//! Box in global unit
Pietro Incardona's avatar
Pietro Incardona committed
37
	Box<dim,size_t> bx;
Pietro Incardona's avatar
Pietro Incardona committed
38
	//! In which sector live the box
39
	comb<dim> cmb;
Pietro Incardona's avatar
Pietro Incardona committed
40
	//! Global id of the internal ghost box
Pietro Incardona's avatar
Pietro Incardona committed
41
	size_t g_id;
Pietro Incardona's avatar
Pietro Incardona committed
42
	//! from which sub-domain this internal ghost box is generated (or with which sub-domain is overlapping)
43
	size_t r_sub;
Pietro Incardona's avatar
Pietro Incardona committed
44 45
};

incardon's avatar
incardon committed
46
#define NO_GDB_EXT_SWITCH 0x1000
incardon's avatar
incardon committed
47 48 49

/*! \brief This is a distributed grid
 *
50 51
 * Implementation of a distributed grid the decomposition is geometrical, grid
 * is splitted across several processor
incardon's avatar
incardon committed
52 53
 *
 * \param dim Dimensionality of the grid
54 55
 * \param St Type of space where the grid is living
 * \param T object the grid is storing
incardon's avatar
incardon committed
56 57 58 59
 * \param Decomposition Class that decompose the grid for example CartDecomposition
 * \param Mem Is the allocator
 * \param device type of base structure is going to store the data
 *
60
 * ### Create a distributed grid and access it
Pietro Incardona's avatar
Pietro Incardona committed
61
 * \snippet grid_dist_id_unit_test.cpp Create and access a distributed grid
62
 * ### Synchronize the ghosts and check the information
Pietro Incardona's avatar
Pietro Incardona committed
63
 * \snippet grid_dist_id_unit_test.cpp Synchronize the ghost and check the information
64
 * ### Create and access a distributed grid for complex structures
Pietro Incardona's avatar
Pietro Incardona committed
65
 * \snippet grid_dist_id_unit_test.cpp Create and access a distributed grid complex
incardon's avatar
incardon committed
66
 * ### Synchronize a distributed grid for complex structures
Pietro Incardona's avatar
Pietro Incardona committed
67
 * \snippet grid_dist_id_unit_test.cpp Synchronized distributed grid complex
incardon's avatar
incardon committed
68
 * ### Usage of a grid dist iterator sub
incardon's avatar
incardon committed
69
 * \snippet grid_dist_id_iterators_unit_tests.hpp Usage of a sub_grid iterator
incardon's avatar
incardon committed
70
 * ### Construct two grid with the same decomposition
Pietro Incardona's avatar
Pietro Incardona committed
71
 * \snippet grid_dist_id_unit_test.cpp Construct two grid with the same decomposition
72
 *
incardon's avatar
incardon committed
73
 */
incardon's avatar
incardon committed
74 75 76 77 78 79
template<unsigned int dim,
		 typename St,
		 typename T,
		 typename Decomposition = CartDecomposition<dim,St>,
		 typename Memory=HeapMemory ,
		 typename device_grid=grid_cpu<dim,T> >
Yaroslav's avatar
Yaroslav committed
80
class grid_dist_id : public grid_dist_id_comm<dim,St,T,Decomposition,Memory,device_grid>
incardon's avatar
incardon committed
81
{
Pietro Incardona's avatar
Pietro Incardona committed
82
	//! Domain
83 84
	Box<dim,St> domain;

Pietro Incardona's avatar
Pietro Incardona committed
85
	//! Ghost expansion
86
	Ghost<dim,St> ghost;
incardon's avatar
incardon committed
87

incardon's avatar
incardon committed
88 89 90
	//! Ghost expansion
	Ghost<dim,long int> ghost_int;

incardon's avatar
incardon committed
91
	//! Local grids
Yaroslav's avatar
Yaroslav committed
92
	mutable openfpm::vector<device_grid> loc_grid;
incardon's avatar
incardon committed
93

Yaroslav's avatar
Yaroslav committed
94 95
	//! Old local grids
	mutable openfpm::vector<device_grid> loc_grid_old;
incardon's avatar
incardon committed
96 97

	//! Space Decomposition
98
	Decomposition dec;
incardon's avatar
incardon committed
99

incardon's avatar
incardon committed
100 101 102 103 104 105 106 107 108 109 110
	//! gdb_ext markers
	//! In the case where the grid is defined everywhere
	//! gdb_ext_marker is useless and so is empty
	//! in the case we have a grid defined on a smaller set
	//! of boxes gbd_ext_markers indicate the division across
	//! subdomains. For example Sub-domain 0 produce 2 grid
	//! Sub-domain 1 produce 3 grid Sub-domain 2 produce 2 grid
	//! Sub-domain 3 produce 1 grid
	//! gdb_ext_markers contain 0,2,5,7,8
	openfpm::vector<size_t> gdb_ext_markers;

incardon's avatar
incardon committed
111 112 113
	//! Extension of each grid: Domain and ghost + domain
	openfpm::vector<GBoxes<device_grid::dims>> gdb_ext;

Yaroslav's avatar
Yaroslav committed
114 115 116 117 118
	//! Global gdb_ext
	mutable openfpm::vector<GBoxes<device_grid::dims>> gdb_ext_global;

	//! Extension of each old grid (old): Domain and ghost + domain
	openfpm::vector<GBoxes<device_grid::dims>> gdb_ext_old;
Yaroslav's avatar
Yaroslav committed
119

incardon's avatar
incardon committed
120 121 122
	//! Size of the grid on each dimension
	size_t g_sz[dim];

123
	//! Structure that divide the space into cells
Pietro Incardona's avatar
Pietro Incardona committed
124
	CellDecomposer_sm<dim,St,shift<dim,St>> cd_sm;
incardon's avatar
incardon committed
125

126 127 128
	//! size of the insert pool on gpu
	size_t gpu_insert_pool_size;

129
	//! Communicator class
130
	Vcluster<> & v_cl;
incardon's avatar
incardon committed
131

132 133 134
	//! properties names
	openfpm::vector<std::string> prp_names;

135
	//! It map a global ghost id (g_id) to the external ghost box information
Pietro Incardona's avatar
Pietro Incardona committed
136
	//! It is unique across all the near processor
137 138
	std::unordered_map<size_t,size_t> g_id_to_external_ghost_box;

incardon's avatar
incardon committed
139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208
	/*! Link a received external ghost box to the linked eg_box.
	 * When the grid is defined everywhere for each received external ghost box
	 * exists one eg_box linked to it that contain the information
	 * on how to transfer the information to the associated sub-domain grid.
	 * Unfortunately when we specify where the grid is defined, a received external
	 * ghost box can be linked to multiple sub-domain grids (one sub-domain can have multiple
	 * sub grids).
	 * So in standard situation (grid defined everywhere)
	 * a received external ghost box is linked to a single
	 * eg_box entry and eb_gid_list play mainly no role.
	 * (play no role but must be filled ghost_get expect this structure to be
	 * filled consistently, it will be clear later how to do it in this case).
	 * When the grid is not defined everywhere a received ghost box can be linked
	 * to multiple external ghost boxes. (Like in figure)
	 *
	 * \verbatim
  +--------------------------------------------+------
  |        Sub-domain                          | Another sub-domain
  |          +------+           +---------+    |
  |          |      |           |         |    |
  |          |      |           |         |    |
  |          |      |           |         |    |
  |          |  3   |           |    4    |    |
  |   empty  |      |    empty  |         |    |
  |          |      |           |         |    |
  |          |      |           |         |    |
  |          |      |           |         |    |                         1
  |          |      |           |         |    |
+-+-----+----+------+-----------+---------+----+-----+-----   Processor bound
        |***##########*********#############****|****|
        |                                            |                   0
        |                                            |
        |                                            |
        |                  9                         |
        |                                            |
        |                                            |
        |                                            |
        +--------------------------------------------+

	 * \endverbatim
	 *
	 * As we can see here the grid number 9 on processo 0 has an internal ghost box
	 * The internal ghost-box is sent to processor 1 and is a received external
	 * ghost box. This external ghost box is partially shared in two separated grids.
	 * It is important to note that 3 and 4 are grid defined externally and are
	 * not defined by the sub-domain border. It is important also to note that the sub-domain
	 * granularity in processor 1 define the granularity of the internal ghost box in
	 * processor 0 and consequently every external ghost box in processor 1 is linked
	 *  uniquely with one internal ghost box in processor 0. On the other hand if we have
	 *  a secondary granularity define by external boxes like 3 and 4 this is not anymore
	 *  true and one internal ghost box in 0 can be linked with multiple grids.
	 *  The granularity of the space division is different from the granularity of where
	 *  the grid is defined. Space decomposition exist independently from the data-structure
	 *  and can be shared across multiple data-structure this mean that cannot be redefined
	 *  based on where is the grid definitions.
	 * The internal ghost box could be redefined in order to respect the granularity.
	 * We do not do this for 3 main reason.
	 *
	 * 1) The definition box must be communicated across processors.
	 * 2) An interprocessor global-id link must be established with lower sub-domain
	 *    granularty
	 * 3) Despite the points * are not linked, but must be anyway sent
	 *    to processor 1, this mean that make not too much sense to increase
	 *    the granularity in advance on processor 0, but it is better receive
	 *    the information an than solve the lower granularity locally
	 *    on processor 1
	 *
	 */
	openfpm::vector<e_box_multi<dim>> eb_gid_list;

Pietro Incardona's avatar
Pietro Incardona committed
209 210 211 212
	//! It map a global ghost id (g_id) to the internal ghost box information
	//! (is unique for processor), it is not unique across all the near processor
	openfpm::vector<std::unordered_map<size_t,size_t>> g_id_to_internal_ghost_box;

Pietro Incardona's avatar
Pietro Incardona committed
213
	//! Receiving size
214 215
	openfpm::vector<size_t> recv_sz;

Pietro Incardona's avatar
Pietro Incardona committed
216
	//! Receiving buffer for particles ghost get
217 218
	openfpm::vector<HeapMemory> recv_mem_gg;

Pietro Incardona's avatar
Pietro Incardona committed
219
	//! Grid informations object
incardon's avatar
incardon committed
220 221
	grid_sm<dim,T> ginfo;

Pietro Incardona's avatar
Pietro Incardona committed
222
	//! Grid informations object without type
incardon's avatar
incardon committed
223 224
	grid_sm<dim,void> ginfo_v;

incardon's avatar
incardon committed
225 226 227 228 229 230
	//! Set of boxes that define where the grid is defined
	openfpm::vector<Box<dim,long int>> bx_def;

	//! Indicate if we have to use bx_def to define the grid
	bool use_bx_def = false;

Pietro Incardona's avatar
Pietro Incardona committed
231 232 233 234 235 236
	//! Indicate if the local internal ghost box has been initialized
	bool init_local_i_g_box = false;

	//! Indicate if the local external ghost box has been initialized
	bool init_local_e_g_box = false;

incardon's avatar
incardon committed
237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260
	//! Flag that indicate if the external ghost box has been initialized
	bool init_e_g_box = false;

	//! Flag that indicate if the internal ghost box has been initialized
	bool init_i_g_box = false;

	//! Flag that indicate if the internal and external ghost box has been fixed
	bool init_fix_ie_g_box = false;

	//! Internal ghost boxes in grid units
	openfpm::vector<ip_box_grid<dim>> ig_box;

	//! External ghost boxes in grid units
	openfpm::vector<ep_box_grid<dim>> eg_box;

	//! Local internal ghost boxes in grid units
	openfpm::vector<i_lbox_grid<dim>> loc_ig_box;

	//! Local external ghost boxes in grid units
	openfpm::vector<e_lbox_grid<dim>> loc_eg_box;

	//! Number of sub-sub-domain for each processor
	size_t v_sub_unit_factor = 64;

261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288
	/*! \brief Call-back to allocate buffer to receive incoming objects (external ghost boxes)
	 *
	 * \param msg_i message size required to receive from i
	 * \param total_msg message size to receive from all the processors
	 * \param total_p the total number of processor want to communicate with you
	 * \param i processor id
	 * \param ri request id (it is an id that goes from 0 to total_p, and is unique
	 *           every time message_alloc is called)
	 * \param ptr void pointer parameter for additional data to pass to the call-back
	 *
	 */
	static void * msg_alloc_external_box(size_t msg_i ,size_t total_msg, size_t total_p, size_t i, size_t ri, void * ptr)
	{
		grid_dist_id<dim,St,T,Decomposition,Memory,device_grid> * g = static_cast<grid_dist_id<dim,St,T,Decomposition,Memory,device_grid> *>(ptr);

		g->recv_sz.resize(g->dec.getNNProcessors());
		g->recv_mem_gg.resize(g->dec.getNNProcessors());

		// Get the local processor id
		size_t lc_id = g->dec.ProctoID(i);

		// resize the receive buffer
		g->recv_mem_gg.get(lc_id).resize(msg_i);
		g->recv_sz.get(lc_id) = msg_i;

		return g->recv_mem_gg.get(lc_id).getPointer();
	}

incardon's avatar
incardon committed
289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306
	/*! \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)
	{
incardon's avatar
incardon committed
307
		if (g.isInvalidGhost() == true || use_bx_def == true)
incardon's avatar
incardon committed
308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335
		{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);}
		}
	}

incardon's avatar
incardon committed
336 337


338
	/*! \brief Create per-processor internal ghost boxes list in grid units and g_id_to_external_ghost_box
339 340 341 342
	 *
	 */
	void create_ig_box()
	{
incardon's avatar
incardon committed
343 344 345
		// temporal vector used for computation
		openfpm::vector_std<result_box<dim>> ibv;

Pietro Incardona's avatar
Pietro Incardona committed
346 347
		if (init_i_g_box == true)	return;

348 349 350
		// Get the grid info
		auto g = cd_sm.getGrid();

Pietro Incardona's avatar
Pietro Incardona committed
351
		g_id_to_internal_ghost_box.resize(dec.getNNProcessors());
352 353 354 355 356 357 358 359 360 361 362

		// Get the number of near processors
		for (size_t i = 0 ; i < dec.getNNProcessors() ; i++)
		{
			ig_box.add();
			auto&& pib = ig_box.last();

			pib.prc = dec.IDtoProc(i);
			for (size_t j = 0 ; j < dec.getProcessorNIGhost(i) ; j++)
			{
				// Get the internal ghost boxes and transform into grid units
incardon's avatar
incardon committed
363
				::Box<dim,St> ib_dom = dec.getProcessorIGhostBox(i,j);
364
				::Box<dim,long int> ib = cd_sm.convertDomainSpaceIntoGridUnits(ib_dom,dec.periodicity());
incardon's avatar
incardon committed
365

incardon's avatar
incardon committed
366 367 368 369 370 371 372 373 374 375 376 377
				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);

incardon's avatar
incardon committed
378 379 380 381 382
				// Here we intersect the internal ghost box with the definition boxes
				// this operation make sense when the grid is not defined in the full
				// domain and we have to intersect the internal ghost box with all the box
				// that define where the grid is defined
				bx_intersect<dim>(bx_def,use_bx_def,ib,ibv);
incardon's avatar
incardon committed
383

incardon's avatar
incardon committed
384 385 386 387 388
				for (size_t k = 0 ; k < ibv.size() ; k++)
				{
					// Check if ib is valid if not it mean that the internal ghost does not contain information so skip it
					if (ibv.get(k).bx.isValid() == false)
					{continue;}
incardon's avatar
incardon committed
389

incardon's avatar
incardon committed
390 391
					// save the box and the sub-domain id (it is calculated as the linearization of P1)
					::Box<dim,size_t> cvt = ibv.get(k).bx;
incardon's avatar
incardon committed
392

incardon's avatar
incardon committed
393 394 395
					i_box_id<dim> bid_t;
					bid_t.box = cvt;
					bid_t.g_id = dec.getProcessorIGhostId(i,j) | (k) << 52;
incardon's avatar
incardon committed
396

incardon's avatar
incardon committed
397 398 399 400
					bid_t.sub = convert_to_gdb_ext(dec.getProcessorIGhostSub(i,j),
							                       ibv.get(k).id,
												   gdb_ext,
												   gdb_ext_markers);
401

incardon's avatar
incardon committed
402 403 404
					bid_t.cmb = dec.getProcessorIGhostPos(i,j);
					bid_t.r_sub = dec.getProcessorIGhostSSub(i,j);
					pib.bid.add(bid_t);
Pietro Incardona's avatar
Pietro Incardona committed
405

incardon's avatar
incardon committed
406 407
					g_id_to_internal_ghost_box.get(i)[bid_t.g_id | (k) << 52 ] = pib.bid.size()-1;
				}
408 409 410 411 412 413
			}
		}

		init_i_g_box = true;
	}

incardon's avatar
incardon committed
414

415 416 417 418 419 420 421 422 423 424
	/*! \brief Create per-processor internal ghost box list in grid units
	 *
	 */
	void create_eg_box()
	{
		// Get the grid info
		auto g = cd_sm.getGrid();

		if (init_e_g_box == true)	return;

425 426 427 428 429 430 431 432 433 434 435 436 437 438 439 440 441 442 443 444 445 446 447
		// Here we collect all the calculated internal ghost box in the sector different from 0 that this processor has

		openfpm::vector<size_t> prc;
		openfpm::vector<size_t> prc_recv;
		openfpm::vector<size_t> sz_recv;
		openfpm::vector<openfpm::vector<Box_fix<dim>>> box_int_send(dec.getNNProcessors());
		openfpm::vector<openfpm::vector<Box_fix<dim>>> box_int_recv;

		for(size_t i = 0 ; i < dec.getNNProcessors() ; i++)
		{
			for (size_t j = 0 ; j < ig_box.get(i).bid.size() ; j++)
			{
				box_int_send.get(i).add();
				box_int_send.get(i).last().bx = ig_box.get(i).bid.get(j).box;
				box_int_send.get(i).last().g_id = ig_box.get(i).bid.get(j).g_id;
				box_int_send.get(i).last().r_sub = ig_box.get(i).bid.get(j).r_sub;
				box_int_send.get(i).last().cmb = ig_box.get(i).bid.get(j).cmb;
			}
			prc.add(dec.IDtoProc(i));
		}

		v_cl.SSendRecv(box_int_send,box_int_recv,prc,prc_recv,sz_recv);

Pietro Incardona's avatar
Pietro Incardona committed
448
		eg_box.resize(dec.getNNProcessors());
449

Pietro Incardona's avatar
Pietro Incardona committed
450
		for (size_t i = 0 ; i < eg_box.size() ; i++)
incardon's avatar
incardon committed
451
		{eg_box.get(i).prc = dec.IDtoProc(i);}
452 453 454

		for (size_t i = 0 ; i < box_int_recv.size() ; i++)
		{
incardon's avatar
incardon committed
455
			size_t pib_cnt = 0;
456
			size_t p_id = dec.ProctoID(prc_recv.get(i));
Pietro Incardona's avatar
Pietro Incardona committed
457
			auto&& pib = eg_box.get(p_id);
458 459
			pib.prc = prc_recv.get(i);

incardon's avatar
incardon committed
460 461 462
			eg_box.get(p_id).recv_pnt = 0;
			eg_box.get(p_id).n_r_box = box_int_recv.get(i).size();

463 464 465
			// For each received internal ghost box
			for (size_t j = 0 ; j < box_int_recv.get(i).size() ; j++)
			{
incardon's avatar
incardon committed
466 467 468
				size_t vol_recv = box_int_recv.get(i).get(j).bx.getVolumeKey();

				eg_box.get(p_id).recv_pnt += vol_recv;
469 470
				size_t send_list_id = box_int_recv.get(i).get(j).r_sub;

incardon's avatar
incardon committed
471 472 473 474 475 476 477 478 479 480 481 482 483 484 485 486 487
				if (use_bx_def == true)
				{
					// First we understand if the internal ghost box sent intersect
					// some local extended sub-domain.

					// eb_gid_list, for an explanation check the declaration
					eb_gid_list.add();
					eb_gid_list.last().e_id = p_id;

					// Now we have to check if a received external ghost box intersect one
					// or more sub-grids
					for (size_t k = 0 ; k < gdb_ext.size() ; k++)
					{
						Box<dim,long int> bx = gdb_ext.get(k).GDbox;
						bx += gdb_ext.get(k).origin;

						Box<dim,long int> output;
incardon's avatar
incardon committed
488
						Box<dim,long int> flp_i = flip_box(box_int_recv.get(i).get(j).bx,box_int_recv.get(i).get(j).cmb,ginfo);
incardon's avatar
incardon committed
489 490 491 492 493 494 495 496 497 498 499 500 501 502 503 504 505 506 507 508 509 510 511 512 513 514 515 516 517 518 519 520 521 522 523 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540

						// it intersect one sub-grid
						if (bx.Intersect(flp_i,output))
						{
							// link

							size_t g_id = box_int_recv.get(i).get(j).g_id;
							add_eg_box<dim>(k,box_int_recv.get(i).get(j).cmb,output,
									g_id,
									gdb_ext.get(k).origin,
									box_int_recv.get(i).get(j).bx.getP1(),
									pib.bid);

							eb_gid_list.last().eb_list.add(pib.bid.size() - 1);

							g_id_to_external_ghost_box[g_id] = eb_gid_list.size() - 1;
						}
					}

					// now we check if exist a full match across the full intersected
					// ghost parts

					bool no_match = true;
					for (size_t k = 0 ; k < eb_gid_list.last().eb_list.size() ; k++)
					{
						size_t eb_id = eb_gid_list.last().eb_list.get(k);

						if (pib.bid.get(eb_id).g_e_box == box_int_recv.get(i).get(j).bx)
						{
							// full match found

							eb_gid_list.last().full_match = eb_id;
							no_match = false;

							break;
						}
					}

					// This is the case where a full match has not been found. In this case we
					// generate an additional gdb_ext and local grid with only external ghost

					if (no_match == true)
					{
						// Create a grid with the same size of the external ghost

						size_t sz[dim];
						for (size_t s = 0 ; s < dim ; s++)
						{sz[s] = box_int_recv.get(i).get(j).bx.getHigh(s) - box_int_recv.get(i).get(j).bx.getLow(s) + 1;}

						// Add an unlinked gdb_ext
						// An unlinked gdb_ext is an empty domain with only a external ghost
						// part
incardon's avatar
incardon committed
541
						Box<dim,long int> output = flip_box(box_int_recv.get(i).get(j).bx,box_int_recv.get(i).get(j).cmb,ginfo);
incardon's avatar
incardon committed
542 543 544 545 546

						GBoxes<dim> tmp;
						tmp.GDbox = box_int_recv.get(i).get(j).bx;
						tmp.GDbox -= tmp.GDbox.getP1();
						tmp.origin = output.getP1();
547
						for (size_t s = 0 ; s < dim ; s++)
incardon's avatar
incardon committed
548 549
						{
							// we set an invalid box, there is no-domain
550 551
							tmp.Dbox.setLow(s,0);
							tmp.Dbox.setHigh(s,-1);
incardon's avatar
incardon committed
552 553 554 555 556 557 558 559 560 561 562 563 564 565 566 567 568 569 570 571 572 573 574 575 576 577 578 579 580 581 582 583 584 585 586 587 588 589 590 591 592 593 594 595 596 597 598 599 600
						}
						tmp.k = -1;
						gdb_ext.add(tmp);

						// create the local grid

						loc_grid.add();
						loc_grid.last().resize(sz);

						// Add an external ghost box

						size_t g_id = box_int_recv.get(i).get(j).g_id;
						add_eg_box<dim>(gdb_ext.size()-1,box_int_recv.get(i).get(j).cmb,output,
								g_id,
								gdb_ext.get(gdb_ext.size()-1).origin,
								box_int_recv.get(i).get(j).bx.getP1(),
								pib.bid);

						// now we map the received ghost box to the information of the
						// external ghost box created
						eb_gid_list.last().full_match = pib.bid.size() - 1;
						eb_gid_list.last().eb_list.add(pib.bid.size() - 1);
						g_id_to_external_ghost_box[g_id] = eb_gid_list.size() - 1;
					}

					// now we create lr_e_box

					size_t fm = eb_gid_list.last().full_match;
					size_t sub_id = pib.bid.get(fm).sub;

					for ( ; pib_cnt < pib.bid.size() ; pib_cnt++)
					{
						pib.bid.get(pib_cnt).lr_e_box -= gdb_ext.get(sub_id).origin;
					}

				}
				else
				{
					// Get the list of the sent sub-domains
					// and recover the id of the sub-domain from
					// the sent list
					const openfpm::vector<size_t> & s_sub = dec.getSentSubdomains(p_id);

					size_t sub_id = s_sub.get(send_list_id);

					e_box_id<dim> bid_t;
					bid_t.sub = sub_id;
					bid_t.cmb = box_int_recv.get(i).get(j).cmb;
					bid_t.cmb.sign_flip();
incardon's avatar
incardon committed
601
					::Box<dim,long int> ib = flip_box(box_int_recv.get(i).get(j).bx,box_int_recv.get(i).get(j).cmb,ginfo);
incardon's avatar
incardon committed
602 603 604 605 606 607 608 609 610 611 612 613 614 615 616
					bid_t.g_e_box = ib;
					bid_t.g_id = box_int_recv.get(i).get(j).g_id;
					// Translate in local coordinate
					Box<dim,long int> tb = ib;
					tb -= gdb_ext.get(sub_id).origin;
					bid_t.l_e_box = tb;

					pib.bid.add(bid_t);
					eb_gid_list.add();
					eb_gid_list.last().eb_list.add(pib.bid.size()-1);
					eb_gid_list.last().e_id = p_id;
					eb_gid_list.last().full_match = pib.bid.size()-1;

					g_id_to_external_ghost_box[bid_t.g_id] = eb_gid_list.size()-1;
				}
617 618 619
			}
		}

Pietro Incardona's avatar
Pietro Incardona committed
620
		init_e_g_box = true;
621 622
	}

incardon's avatar
incardon committed
623 624 625 626 627
	/*! \brief Create local internal ghost box in grid units
	 *
	 */
	void create_local_ig_box()
	{
incardon's avatar
incardon committed
628 629
		openfpm::vector_std<result_box<dim>> ibv;

incardon's avatar
incardon committed
630 631 632 633 634
		// Get the grid info
		auto g = cd_sm.getGrid();

		if (init_local_i_g_box == true)	return;

Pietro Incardona's avatar
Pietro Incardona committed
635
		// Get the number of sub-domains
Pietro Incardona's avatar
Pietro Incardona committed
636
		for (size_t i = 0 ; i < dec.getNSubDomain() ; i++)
incardon's avatar
incardon committed
637 638 639 640 641 642
		{
			loc_ig_box.add();
			auto&& pib = loc_ig_box.last();

			for (size_t j = 0 ; j < dec.getLocalNIGhost(i) ; j++)
			{
incardon's avatar
incardon committed
643 644 645 646 647 648 649 650 651 652 653 654 655 656 657 658
				if (use_bx_def == true)
				{
					// Get the internal ghost boxes and transform into grid units
					::Box<dim,St> ib_dom = dec.getLocalIGhostBox(i,j);
					::Box<dim,long int> ib = cd_sm.convertDomainSpaceIntoGridUnits(ib_dom,dec.periodicity());

					// Here we intersect the internal ghost box with the definition boxes
					// this operation make sense when the grid is not defined in the full
					// domain and we have to intersect the internal ghost box with all the box
					// that define where the grid is defined
					bx_intersect<dim>(bx_def,use_bx_def,ib,ibv);

					for (size_t k = 0 ; k < ibv.size() ; k++)
					{
						// Check if ib is valid if not it mean that the internal ghost does not contain information so skip it
						if (ibv.get(k).bx.isValid() == false)
incardon's avatar
incardon committed
659
						{continue;}
incardon's avatar
incardon committed
660 661 662 663 664 665 666 667 668 669 670 671 672 673 674 675 676 677 678 679 680 681 682 683 684

						pib.bid.add();
						pib.bid.last().box = ibv.get(k).bx;

						pib.bid.last().sub_gdb_ext = convert_to_gdb_ext(i,
																ibv.get(k).id,
																gdb_ext,
																gdb_ext_markers);

						pib.bid.last().sub = dec.getLocalIGhostSub(i,j);

						// It will be filled later
						pib.bid.last().cmb = dec.getLocalIGhostPos(i,j);
					}
				}
				else
				{
					// Get the internal ghost boxes and transform into grid units
					::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;

incardon's avatar
incardon committed
685 686 687 688 689 690 691 692 693 694 695 696 697
					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;

incardon's avatar
incardon committed
698 699 700 701 702 703 704 705 706 707 708 709 710 711 712 713 714 715 716 717 718
					pib.bid.add();
					pib.bid.last().box = ib;
					pib.bid.last().sub = dec.getLocalIGhostSub(i,j);
					pib.bid.last().k.add(dec.getLocalIGhostE(i,j));
					pib.bid.last().cmb = dec.getLocalIGhostPos(i,j);
					pib.bid.last().sub_gdb_ext = i;
				}
			}

			if (use_bx_def == true)
			{
				// unfortunately boxes that define where the grid is located can generate
				// additional internal ghost boxes

				for (size_t j = gdb_ext_markers.get(i) ; j < gdb_ext_markers.get(i+1) ; j++)
				{
					// intersect within box in the save sub-domain

					for (size_t k = gdb_ext_markers.get(i) ; k < gdb_ext_markers.get(i+1) ; k++)
					{
						if (j == k)	{continue;}
incardon's avatar
incardon committed
719

incardon's avatar
incardon committed
720 721 722 723 724
						// extend k and calculate the internal ghost box
						Box<dim,long int> bx_e =  gdb_ext.get(k).GDbox;
						bx_e += gdb_ext.get(k).origin;
						Box<dim,long int> bx = gdb_ext.get(j).Dbox;
						bx += gdb_ext.get(j).origin;
incardon's avatar
incardon committed
725

incardon's avatar
incardon committed
726 727 728 729
						Box<dim,long int> output;
						if (bx.Intersect(bx_e, output) == true)
						{
							pib.bid.add();
incardon's avatar
incardon committed
730

incardon's avatar
incardon committed
731
							pib.bid.last().box = output;
incardon's avatar
incardon committed
732

incardon's avatar
incardon committed
733 734
							pib.bid.last().sub_gdb_ext = j;
							pib.bid.last().sub = i;
incardon's avatar
incardon committed
735

incardon's avatar
incardon committed
736 737 738 739 740 741
//							if (use_bx_def == true)
//							{pib.bid.last().k = -1;}
//							else
//							{pib.bid.last().k = dec.getLocalIGhostE(i,j);}
							// these ghost always in the quadrant zero
							pib.bid.last().cmb.zero();
incardon's avatar
incardon committed
742

incardon's avatar
incardon committed
743 744 745
						}
					}
				}
incardon's avatar
incardon committed
746
			}
incardon's avatar
incardon committed
747

incardon's avatar
incardon committed
748 749
		}

incardon's avatar
incardon committed
750

incardon's avatar
incardon committed
751 752 753
		init_local_i_g_box = true;
	}

754
	/*! \brief Create per-processor external ghost boxes list in grid units
incardon's avatar
incardon committed
755 756 757 758 759 760 761 762 763
	 *
	 */
	void create_local_eg_box()
	{
		// Get the grid info
		auto g = cd_sm.getGrid();

		if (init_local_e_g_box == true)	return;

Pietro Incardona's avatar
Pietro Incardona committed
764
		loc_eg_box.resize(dec.getNSubDomain());
765 766 767 768 769 770

		// Get the number of sub-domain
		for (size_t i = 0 ; i < dec.getNSubDomain() ; i++)
		{
			for (size_t j = 0 ; j < loc_ig_box.get(i).bid.size() ; j++)
			{
incardon's avatar
incardon committed
771 772 773 774 775 776 777 778 779 780 781 782 783 784 785 786 787 788 789
				long int volume_linked = 0;

				size_t le_sub = loc_ig_box.get(i).bid.get(j).sub;
				auto & pib = loc_eg_box.get(le_sub);

				if (use_bx_def == true)
				{

					// We check if an external local ghost box intersect one
					// or more sub-grids
					for (size_t k = 0 ; k < gdb_ext.size() ; k++)
					{
						Box<dim,long int> bx = gdb_ext.get(k).Dbox;
						bx += gdb_ext.get(k).origin;

						Box<dim,long int> gbx = gdb_ext.get(k).GDbox;
						gbx += gdb_ext.get(k).origin;

						Box<dim,long int> output;
incardon's avatar
incardon committed
790
						Box<dim,long int> flp_i = flip_box(loc_ig_box.get(i).bid.get(j).box,loc_ig_box.get(i).bid.get(j).cmb,ginfo);
incardon's avatar
incardon committed
791 792 793 794 795 796 797 798 799 800 801 802 803 804 805 806 807 808 809

						bool intersect_domain = bx.Intersect(flp_i,output);
						bool intersect_gdomain = gbx.Intersect(flp_i,output);

						// it intersect one sub-grid
						if (intersect_domain == false && intersect_gdomain == true)
						{
							// fill the link variable
							loc_ig_box.get(i).bid.get(j).k.add(pib.bid.size());
							size_t s = loc_ig_box.get(i).bid.get(j).k.last();

							comb<dim> cmb = loc_ig_box.get(i).bid.get(j).cmb;
							cmb.sign_flip();

							add_loc_eg_box(le_sub,
										   loc_ig_box.get(i).bid.get(j).sub_gdb_ext,
										   j,
										   k,
										   pib.bid,
incardon's avatar
incardon committed
810
										   output,
incardon's avatar
incardon committed
811 812 813
										   cmb);


incardon's avatar
incardon committed
814
							volume_linked += pib.bid.last().ebox.getVolumeKey();
incardon's avatar
incardon committed
815 816 817 818 819 820 821 822 823 824 825
						}
					}
				}
				else
				{
					size_t k = loc_ig_box.get(i).bid.get(j).sub;
					auto & pib = loc_eg_box.get(k);

					size_t s = loc_ig_box.get(i).bid.get(j).k.get(0);
					pib.bid.resize(dec.getLocalNEGhost(k));

incardon's avatar
incardon committed
826
					pib.bid.get(s).ebox = flip_box(loc_ig_box.get(i).bid.get(j).box,loc_ig_box.get(i).bid.get(j).cmb,ginfo);
incardon's avatar
incardon committed
827 828 829 830 831 832 833
					pib.bid.get(s).sub = dec.getLocalEGhostSub(k,s);
					pib.bid.get(s).cmb = loc_ig_box.get(i).bid.get(j).cmb;
					pib.bid.get(s).cmb.sign_flip();
					pib.bid.get(s).k = j;
					pib.bid.get(s).initialized = true;
					pib.bid.get(s).sub_gdb_ext = k;
				}
834 835 836
			}
		}

Pietro Incardona's avatar
Pietro Incardona committed
837
		init_local_e_g_box = true;
incardon's avatar
incardon committed
838 839
	}

incardon's avatar
incardon committed
840

841
	/*! \brief Check the grid has a valid size
Pietro Incardona's avatar
Pietro Incardona committed
842 843
	 *
	 * \param g_sz size of the grid
incardon's avatar
incardon committed
844 845 846 847 848 849 850
	 *
	 */
	inline void check_size(const size_t (& g_sz)[dim])
	{
		for (size_t i = 0 ; i < dim ; i++)
		{
			if (g_sz[i] < 2)
incardon's avatar
incardon committed
851
				std::cerr << "Error: " << __FILE__ << ":" << __LINE__ << " distributed grids with size smaller than 2 are not supported\n";
incardon's avatar
incardon committed
852 853
		}
	}
incardon's avatar
incardon committed
854

incardon's avatar
incardon committed
855
	/*! \brief Create the grids on memory
incardon's avatar
incardon committed
856 857 858
	 *
	 * \param bx_def Where the grid is defined
	 * \param use_bx_def use the array that define where the grid is defined
incardon's avatar
incardon committed
859 860
	 *
	 */
incardon's avatar
incardon committed
861 862 863
	void Create(openfpm::vector<Box<dim,long int>> & bx_def,
			    const Ghost<dim,long int> & g,
			    bool use_bx_def)
incardon's avatar
incardon committed
864
	{
865
		// create gdb
incardon's avatar
incardon committed
866 867 868
		create_gdb_ext<dim,Decomposition>(gdb_ext,gdb_ext_markers,dec,cd_sm,bx_def,g,use_bx_def);

		size_t n_grid = gdb_ext.size();
869

incardon's avatar
incardon committed
870
		// create local grids for each hyper-cube
871
		loc_grid.resize(n_grid);
incardon's avatar
incardon committed
872 873 874 875 876 877 878

		// Size of the grid on each dimension
		size_t l_res[dim];

		// Allocate the grids
		for (size_t i = 0 ; i < n_grid ; i++)
		{
879
			SpaceBox<dim,long int> sp_tg = gdb_ext.get(i).GDbox;
incardon's avatar
incardon committed
880 881

			// Get the size of the local grid
882 883 884 885
			// The boxes indicate the extension of the index the size
			// is this extension +1
			// for example a 1D box (interval) from 0 to 3 in one dimension have
			// the points 0,1,2,3 = so a total of 4 points
incardon's avatar
incardon committed
886
			for (size_t j = 0 ; j < dim ; j++)
incardon's avatar
incardon committed
887
			{l_res[j] = (sp_tg.getHigh(j) >= 0)?(sp_tg.getHigh(j)+1):0;}
incardon's avatar
incardon committed
888 889 890 891 892 893 894

			// Set the dimensions of the local grid
			loc_grid.get(i).resize(l_res);
		}
	}


Pietro Incardona's avatar
Pietro Incardona committed
895 896 897 898 899 900
    /*! \brief Initialize the Cell decomposer of the grid enforcing perfect overlap of the cells
	 *
	 * \param cd_old the CellDecomposer we are trying to mach
	 * \param ext extension of the domain
	 *
	 */
Pietro Incardona's avatar
Pietro Incardona committed
901
	inline void InitializeCellDecomposer(const CellDecomposer_sm<dim,St,shift<dim,St>> & cd_old, const Box<dim,size_t> & ext)
Pietro Incardona's avatar
Pietro Incardona committed
902 903 904 905 906
	{
		// Initialize the cell decomposer
		cd_sm.setDimensions(cd_old,ext);
	}

Pietro Incardona's avatar
Pietro Incardona committed
907
    /*! \brief Initialize the Cell decomposer of the grid
incardon's avatar
incardon committed
908
	 *
Pietro Incardona's avatar
Pietro Incardona committed
909 910
	 * \param g_sz Size of the grid
	 * \param bc boundary conditions
incardon's avatar
incardon committed
911 912
	 *
	 */
913
	inline void InitializeCellDecomposer(const size_t (& g_sz)[dim], const size_t (& bc)[dim])
incardon's avatar
incardon committed
914
	{
incardon's avatar
incardon committed
915
		// check that the grid has valid size
incardon's avatar
incardon committed
916 917
		check_size(g_sz);

918
		// get the size of the cell decomposer
incardon's avatar
incardon committed
919
		size_t c_g[dim];
920
		getCellDecomposerPar<dim>(c_g,g_sz,bc);
incardon's avatar
incardon committed
921 922 923

		// Initialize the cell decomposer
		cd_sm.setDimensions(domain,c_g,0);
incardon's avatar
incardon committed
924
	}
incardon's avatar
incardon committed
925

incardon's avatar
incardon committed
926 927 928
	/*! \brief Initialize the grid
	 *
	 * \param g_sz Global size of the grid
Pietro Incardona's avatar
Pietro Incardona committed
929
	 * \param bc boundary conditions
incardon's avatar
incardon committed
930 931
	 *
	 */
932
	inline void InitializeDecomposition(const size_t (& g_sz)[dim], const size_t (& bc)[dim])
incardon's avatar
incardon committed
933
	{
incardon's avatar
incardon committed
934
		// fill the global size of the grid
incardon's avatar
incardon committed
935
		for (size_t i = 0 ; i < dim ; i++)	{this->g_sz[i] = g_sz[i];}
incardon's avatar
incardon committed
936 937 938 939

		// Get the number of processor and calculate the number of sub-domain
		// for decomposition
		size_t n_proc = v_cl.getProcessingUnits();
incardon's avatar
incardon committed
940
		size_t n_sub = n_proc * v_sub_unit_factor;
incardon's avatar
incardon committed
941 942 943 944

		// Calculate the maximum number (before merging) of sub-domain on
		// each dimension
		size_t div[dim];
incardon's avatar
incardon committed
945
		for (size_t i = 0 ; i < dim ; i++)
946
		{div[i] = openfpm::math::round_big_2(pow(n_sub,1.0/dim));}
incardon's avatar
incardon committed
947 948

		// Create the sub-domains
Pietro Incardona's avatar
Pietro Incardona committed
949
		dec.setParameters(div,domain,bc,ghost);
Pietro Incardona's avatar
Pietro Incardona committed
950
		dec.decompose();
incardon's avatar
incardon committed
951 952
	}

incardon's avatar
incardon committed
953 954 955 956 957 958
	/*! \brief Initialize the grid
	 *
	 * \param g_sz Global size of the grid
	 *
	 */
	inline void InitializeStructures(const size_t (& g_sz)[dim])
incardon's avatar
incardon committed
959 960 961 962 963 964 965 966 967 968 969 970 971 972 973 974 975 976 977 978 979
	{
		// an empty
		openfpm::vector<Box<dim,long int>> empty;

		// Ghost zero
		Ghost<dim,long int> zero;

		InitializeStructures(g_sz,empty,zero,false);
	}

	/*! \brief Initialize the grid
	 *
	 * \param g_sz Global size of the grid
	 * \param g ghost extension of the grid in integer unit
	 * \param bx set of boxes that define where is defined the grid
	 *
	 */
	inline void InitializeStructures(const size_t (& g_sz)[dim],
			                         openfpm::vector<Box<dim,long int>> & bx,
									 const Ghost<dim,long int> & g,
									 bool use_bx_def)
incardon's avatar
incardon committed
980 981 982 983 984
	{
		// fill the global size of the grid
		for (size_t i = 0 ; i < dim ; i++)	{this->g_sz[i] = g_sz[i];}

		// Create local grid
incardon's avatar
incardon committed
985
		Create(bx,g,use_bx_def);
incardon's avatar
incardon committed
986 987
	}

incardon's avatar
incardon committed
988 989 990
	// Ghost as integer
	Ghost<dim,long int> gint = Ghost<dim,long int>(0);

Pietro Incardona's avatar
Pietro Incardona committed
991
protected:
incardon's avatar
incardon committed
992 993

	/*! \brief Given a local sub-domain i with a local grid Domain + ghost return the part of the local grid that is domain
Pietro Incardona's avatar
Pietro Incardona committed
994 995
	 *
	 * \param i sub-domain
incardon's avatar
incardon committed
996 997 998 999 1000 1001 1002 1003 1004
	 *
	 * \return the Box defining the domain in the local grid
	 *
	 */
	Box<dim,size_t> getDomain(size_t i)
	{
		return gdb_ext.get(i).Dbox;
	}

Pietro Incardona's avatar
Pietro Incardona committed
1005 1006 1007 1008 1009 1010 1011 1012 1013 1014 1015 1016 1017 1018 1019 1020 1021 1022 1023 1024 1025
	/*! \brief Convert a ghost from grid point units into continus space
	 *
	 * \param gd Ghost in continuous space
	 * \param cd_sm CellDecomposer of the grid
	 *
	 * \return the ghost in continuous unit
	 *
	 */
	static inline Ghost<dim,float> convert_ghost(const Ghost<dim,long int> & gd, const CellDecomposer_sm<dim,St,shift<dim,St>> & cd_sm)
	{
		Ghost<dim,float> gc;

		// get the grid spacing
		Box<dim,St> sp = cd_sm.getCellBox();

		// enlarge 0.001 of the spacing
		sp.magnify_fix_P1(1.1);

		// set the ghost
		for (size_t i = 0 ; i < dim ; i++)
		{
incardon's avatar
incardon committed
1026 1027
			gc.setLow(i,gd.getLow(i)*(sp.getHigh(i)));
			gc.setHigh(i,gd.getHigh(i)*(sp.getHigh(i)));
Pietro Incardona's avatar
Pietro Incardona committed
1028 1029 1030 1031 1032
		}

		return gc;
	}

incardon's avatar
incardon committed
1033 1034 1035 1036 1037 1038 1039 1040 1041 1042 1043 1044 1045 1046 1047 1048 1049 1050 1051 1052 1053 1054 1055 1056 1057 1058 1059
	/*! \brief Set the minimum number of sub-domain per processor
	 *
	 * \param n_sub
	 *
	 */
	void setDecompositionGranularity(size_t n_sub)
	{
		this->v_sub_unit_factor = n_sub;
	}

	void reset_ghost_structures()
	{
		g_id_to_internal_ghost_box.clear();
		ig_box.clear();
		init_i_g_box = false;

		eg_box.clear();
		eb_gid_list.clear();
		init_e_g_box = false;

		init_local_i_g_box = false;
		loc_ig_box.clear();

		init_local_e_g_box = false;
		loc_eg_box.clear();
	}

incardon's avatar
incardon committed
1060 1061
public:

Pietro Incardona's avatar
Pietro Incardona committed
1062
	//! Which kind of grid the structure store
incardon's avatar
incardon committed
1063 1064
	typedef device_grid d_grid;

Pietro Incardona's avatar
Pietro Incardona committed
1065
	//! Decomposition used
incardon's avatar
incardon committed
1066 1067
	typedef Decomposition decomposition;

Pietro Incardona's avatar
Pietro Incardona committed
1068
	//! value_type
Pietro Incardona's avatar
Pietro Incardona committed
1069 1070
	typedef T value_type;

Pietro Incardona's avatar
Pietro Incardona committed
1071
	//! Type of space
1072 1073
	typedef St stype;

Pietro Incardona's avatar
Pietro Incardona committed
1074
	//! Type of Memory
1075 1076
	typedef Memory memory_type;

Pietro Incardona's avatar
Pietro Incardona committed
1077
	//! Type of device grid
1078 1079
	typedef device_grid device_grid_type;

Pietro Incardona's avatar
Pietro Incardona committed
1080
	//! Number of dimensions
1081 1082 1083 1084
	static const unsigned int dims = dim;

	/*! \brief Get the domain where the grid is defined
	 *
Pietro Incardona's avatar
Pietro Incardona committed
1085
	 * \return the domain of the grid
1086 1087 1088 1089 1090 1091 1092
	 *
	 */
	inline const Box<dim,St> getDomain() const
	{
		return domain;
	}

incardon's avatar
incardon committed
1093 1094 1095 1096 1097 1098 1099 1100 1101 1102 1103 1104
	/*! \brief Get the point where it start the origin of the grid of the sub-domain i
	 *
	 * \param i sub-domain
	 *
	 * \return the point
	 *
	 */
	Point<dim,St> getOffset(size_t i)
	{
		return pmul(Point<dim,St>(gdb_ext.get(i).origin), cd_sm.getCellBox().getP2()) + getDomain().getP1();
	}

1105
    /*! \brief Get the spacing of the grid in direction i
Pietro Incardona's avatar
Pietro Incardona committed
1106 1107
     *
     * \param i dimension
1108 1109 1110 1111 1112 1113 1114 1115 1116
     *
     * \return the spacing
     *
     */
    inline St spacing(size_t i) const
    {
    	return cd_sm.getCellBox().getHigh(i);
    }

Pietro Incardona's avatar
Pietro Incardona committed
1117 1118 1119 1120 1121 1122 1123 1124 1125 1126
	/*! \brief Return the total number of points in the grid
	 *
	 * \return number of points
	 *
	 */
	size_t size() const
	{
		return ginfo_v.size();
	}

incardon's avatar
incardon committed
1127 1128 1129 1130 1131 1132 1133 1134 1135 1136 1137 1138
	/*! \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 < loc_grid.size() ; i++)
		{meta_copy<T>::meta_copy_(bv,loc_grid.get(i).getBackgroundValue());}
	}

1139 1140 1141 1142 1143 1144 1145 1146 1147 1148 1149 1150 1151
	/*! \brief set the background value
	 *
	 * You can use this function make sense in case of sparse in case of dense
	 * it does nothing
	 *
	 */
	template<unsigned int p>
	void setBackgroundValue(const typename boost::mpl::at<typename T::type,boost::mpl::int_<p>>::type & bv)
	{
		for (size_t i = 0 ; i < loc_grid.size() ; i++)
		{loc_grid.get(i).template setBackgroundValue<p>(bv);}
	}

incardon's avatar
incardon committed
1152 1153 1154 1155 1156 1157 1158 1159 1160 1161 1162 1163 1164 1165 1166 1167 1168 1169
	/*! \brief Return the local total number of points inserted in the grid
	 *
	 * in case of dense grid it return the number of local points, in case of
	 * sparse it return the number of inserted points
	 *
	 * \return number of points
	 *
	 */
	size_t size_local_inserted() const
	{
		size_t lins = 0;

		for (size_t i = 0 ; i < loc_grid.size() ; i++)
		{lins += loc_grid.get(i).size_inserted();}

		return lins;
	}

Pietro Incardona's avatar
Pietro Incardona committed
1170 1171 1172 1173 1174 1175 1176 1177 1178
	/*! \brief Return the total number of points in the grid
	 *
	 * \param i direction
	 *
	 * \return number of points on direction i
	 *
	 */
	size_t size(size_t i) const
	{
1179
		return ginfo_v.size(i);
Pietro Incardona's avatar
Pietro Incardona committed
1180 1181
	}

incardon's avatar
incardon committed
1182 1183 1184 1185 1186 1187 1188 1189 1190
	/*! \brief Default Copy constructor on this class make no sense and is unsafe, this definition disable it
	 *
	 * \param g grid to copy
	 *
	 */
	grid_dist_id(const grid_dist_id<dim,St,T,Decomposition,Memory,device_grid> & g)
	:grid_dist_id_comm<dim,St,T,Decomposition,Memory,device_grid>(g),
	 domain(g.domain),
	 ghost(g.ghost),
incardon's avatar
incardon committed
1191
	 ghost_int(g.ghost_int),
incardon's avatar
incardon committed
1192 1193 1194 1195 1196 1197 1198 1199 1200 1201 1202 1203 1204 1205 1206 1207 1208 1209 1210 1211 1212 1213 1214 1215 1216
	 loc_grid(g.loc_grid),
	 loc_grid_old(g.loc_grid_old),
	 dec(g.dec),
	 gdb_ext_markers(g.gdb_ext_markers),
	 gdb_ext(g.gdb_ext),
	 gdb_ext_global(g.gdb_ext_global),
	 gdb_ext_old(g.gdb_ext_old),
	 cd_sm(g.cd_sm),
	 v_cl(g.v_cl),
	 prp_names(g.prp_names),
	 g_id_to_external_ghost_box(g.g_id_to_external_ghost_box),
	 g_id_to_internal_ghost_box(g.g_id_to_internal_ghost_box),
	 ginfo(g.ginfo),
	 ginfo_v(g.ginfo_v),
	 init_local_i_g_box(g.init_local_i_g_box),
	 init_local_e_g_box(g.init_local_e_g_box)
	{
#ifdef SE_CLASS2
		check_new(this,8,GRID_DIST_EVENT,4);
#endif

		for (size_t i = 0 ; i < dim ; i++)
		{g_sz[i] = g.g_sz[i];}
	}

Pietro Incardona's avatar
Pietro Incardona committed
1217 1218 1219 1220 1221 1222
	/*! \brief This constructor is special, it construct an expanded grid that perfectly overlap with the previous
	 *
	 * The key-word here is "perfectly overlap". Using the default constructor you could create
	 * something similar, but because of rounding-off error it can happen that it is not perfectly overlapping
	 *
	 * \param g previous grid
Pietro Incardona's avatar
Pietro Incardona committed
1223
	 * \param gh Ghost part in grid units
Pietro Incardona's avatar
Pietro Incardona committed
1224 1225 1226
	 * \param ext extension of the grid (must be positive on every direction)
	 *
	 */
incardon's avatar
incardon committed
1227 1228 1229 1230 1231
	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())
Pietro Incardona's avatar
Pietro Incardona committed
1232 1233 1234 1235 1236
	{
#ifdef SE_CLASS2
		check_new(this,8,GRID_DIST_EVENT,4);
#endif

1237 1238 1239 1240 1241 1242 1243 1244 1245 1246
		size_t ext_dim[dim];
		for (size_t i = 0 ; i < dim ; i++) {ext_dim[i] = g.getGridInfoVoid().size(i) + ext.getKP1().get(i) + ext.getKP2().get(i);}

		// Set the grid info of the extended grid
		ginfo.setDimensions(ext_dim);
		ginfo_v.setDimensions(ext_dim);

		InitializeCellDecomposer(g.getCellDecomposer(),ext);

		ghost = convert_ghost(gh,cd_sm);
Pietro Incardona's avatar
Pietro Incardona committed
1247 1248 1249 1250 1251

		// Extend the grid by the extension part and calculate the domain

		for (size_t i = 0 ; i < dim ; i++)
		{
1252
			g_sz[i] = g.size(i) + ext.getLow(i) + ext.getHigh(i);
Pietro Incardona's avatar
Pietro Incardona committed
1253

1254 1255 1256 1257 1258 1259 1260 1261 1262 1263
			if (g.getDecomposition().periodicity(i) == NON_PERIODIC)
			{
				this->domain.setLow(i,g.getDomain().getLow(i) - ext.getLow(i) * g.spacing(i) - g.spacing(i) / 2.0);
				this->domain.setHigh(i,g.getDomain().getHigh(i) + ext.getHigh(i) * g.spacing(i) + g.spacing(i) / 2.0);
			}
			else
			{
				this->domain.setLow(i,g.getDomain().getLow(i) - ext.getLow(i) * g.spacing(i));
				this->domain.setHigh(i,g.getDomain().getHigh(i) + ext.getHigh(i) * g.spacing(i));
			}
Pietro Incardona's avatar
Pietro Incardona committed
1264 1265
		}

1266
		dec.setParameters(g.getDecomposition(),ghost,this->domain);
Pietro Incardona's avatar
Pietro Incardona committed
1267

incardon's avatar
incardon committed
1268 1269 1270 1271
		// an empty
		openfpm::vector<Box<dim,long int>> empty;

		InitializeStructures(g.getGridInfoVoid().getSize(),empty,gh,false);
Pietro Incardona's avatar