grid_dist_id.hpp 76.6 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
 */
74
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
75
class grid_dist_id : public grid_dist_id_comm<dim,St,T,Decomposition,Memory,device_grid>
incardon's avatar
incardon committed
76
{
Pietro Incardona's avatar
Pietro Incardona committed
77
	//! Domain
78 79
	Box<dim,St> domain;

Pietro Incardona's avatar
Pietro Incardona committed
80
	//! Ghost expansion
81
	Ghost<dim,St> ghost;
incardon's avatar
incardon committed
82

incardon's avatar
incardon committed
83 84 85
	//! Ghost expansion
	Ghost<dim,long int> ghost_int;

incardon's avatar
incardon committed
86
	//! Local grids
Yaroslav's avatar
Yaroslav committed
87
	mutable openfpm::vector<device_grid> loc_grid;
incardon's avatar
incardon committed
88

Yaroslav's avatar
Yaroslav committed
89 90
	//! Old local grids
	mutable openfpm::vector<device_grid> loc_grid_old;
incardon's avatar
incardon committed
91 92

	//! Space Decomposition
93
	Decomposition dec;
incardon's avatar
incardon committed
94

incardon's avatar
incardon committed
95 96 97 98 99 100 101 102 103 104 105
	//! 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
106 107 108
	//! Extension of each grid: Domain and ghost + domain
	openfpm::vector<GBoxes<device_grid::dims>> gdb_ext;

Yaroslav's avatar
Yaroslav committed
109 110 111 112 113
	//! 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
114

incardon's avatar
incardon committed
115 116 117
	//! Size of the grid on each dimension
	size_t g_sz[dim];

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

121 122 123
	//! size of the insert pool on gpu
	size_t gpu_insert_pool_size;

124
	//! Communicator class
125
	Vcluster<> & v_cl;
incardon's avatar
incardon committed
126

127 128 129
	//! properties names
	openfpm::vector<std::string> prp_names;

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

incardon's avatar
incardon committed
134 135 136 137 138 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
	/*! 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
204 205 206 207
	//! 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
208
	//! Receiving size
209 210
	openfpm::vector<size_t> recv_sz;

Pietro Incardona's avatar
Pietro Incardona committed
211
	//! Receiving buffer for particles ghost get
212 213
	openfpm::vector<HeapMemory> recv_mem_gg;

Pietro Incardona's avatar
Pietro Incardona committed
214
	//! Grid informations object
incardon's avatar
incardon committed
215 216
	grid_sm<dim,T> ginfo;

Pietro Incardona's avatar
Pietro Incardona committed
217
	//! Grid informations object without type
incardon's avatar
incardon committed
218 219
	grid_sm<dim,void> ginfo_v;

incardon's avatar
incardon committed
220 221 222 223 224 225
	//! 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
226 227 228 229 230 231
	//! 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
232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255
	//! 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;

256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283
	/*! \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();
	}

284 285
	/*! \brief flip box just convert and internal ghost box into an external ghost box
	 *
Pietro Incardona's avatar
Pietro Incardona committed
286 287
	 * \param box to convert
	 * \param cmb sector position of the box
288
	 *
Pietro Incardona's avatar
Pietro Incardona committed
289 290
	 * \return the converted box
	 *
291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317
	 */
	Box<dim,long int> flip_box(const Box<dim,long int> & box, const comb<dim> & cmb)
	{
		Box<dim,long int> flp;

		for (size_t i = 0 ; i < dim; i++)
		{
			if (cmb[i] == 0)
			{
				flp.setLow(i,box.getLow(i));
				flp.setHigh(i,box.getHigh(i));
			}
			else if (cmb[i] == 1)
			{
				flp.setLow(i,box.getLow(i) + ginfo.size(i));
				flp.setHigh(i,box.getHigh(i) + ginfo.size(i));
			}
			else if (cmb[i] == -1)
			{
				flp.setLow(i,box.getLow(i) - ginfo.size(i));
				flp.setHigh(i,box.getHigh(i) - ginfo.size(i));
			}
		}

		return flp;
	}

incardon's avatar
incardon committed
318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 359 360 361 362 363 364
	/*! \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);}
		}
	}

365
	/*! \brief Create per-processor internal ghost boxes list in grid units and g_id_to_external_ghost_box
366 367 368 369
	 *
	 */
	void create_ig_box()
	{
incardon's avatar
incardon committed
370 371 372
		// temporal vector used for computation
		openfpm::vector_std<result_box<dim>> ibv;

Pietro Incardona's avatar
Pietro Incardona committed
373 374
		if (init_i_g_box == true)	return;

375 376 377
		// Get the grid info
		auto g = cd_sm.getGrid();

Pietro Incardona's avatar
Pietro Incardona committed
378
		g_id_to_internal_ghost_box.resize(dec.getNNProcessors());
379 380 381 382 383 384 385 386 387 388 389

		// 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
390
				::Box<dim,St> ib_dom = dec.getProcessorIGhostBox(i,j);
391
				::Box<dim,long int> ib = cd_sm.convertDomainSpaceIntoGridUnits(ib_dom,dec.periodicity());
incardon's avatar
incardon committed
392

incardon's avatar
incardon committed
393 394 395 396 397
				// 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
398

incardon's avatar
incardon committed
399 400 401 402 403
				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
404

incardon's avatar
incardon committed
405 406
					// 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
407

incardon's avatar
incardon committed
408 409 410
					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
411

incardon's avatar
incardon committed
412 413 414 415
					bid_t.sub = convert_to_gdb_ext(dec.getProcessorIGhostSub(i,j),
							                       ibv.get(k).id,
												   gdb_ext,
												   gdb_ext_markers);
416

incardon's avatar
incardon committed
417 418 419
					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
420

incardon's avatar
incardon committed
421 422
					g_id_to_internal_ghost_box.get(i)[bid_t.g_id | (k) << 52 ] = pib.bid.size()-1;
				}
423 424 425 426 427 428
			}
		}

		init_i_g_box = true;
	}

incardon's avatar
incardon committed
429

430 431 432 433 434 435 436 437 438 439
	/*! \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;

440 441 442 443 444 445 446 447 448 449 450 451 452 453 454 455 456 457 458 459 460 461 462
		// 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
463
		eg_box.resize(dec.getNNProcessors());
464

Pietro Incardona's avatar
Pietro Incardona committed
465
		for (size_t i = 0 ; i < eg_box.size() ; i++)
incardon's avatar
incardon committed
466
		{eg_box.get(i).prc = dec.IDtoProc(i);}
467 468 469

		for (size_t i = 0 ; i < box_int_recv.size() ; i++)
		{
incardon's avatar
incardon committed
470
			size_t pib_cnt = 0;
471
			size_t p_id = dec.ProctoID(prc_recv.get(i));
Pietro Incardona's avatar
Pietro Incardona committed
472
			auto&& pib = eg_box.get(p_id);
473 474
			pib.prc = prc_recv.get(i);

incardon's avatar
incardon committed
475 476 477
			eg_box.get(p_id).recv_pnt = 0;
			eg_box.get(p_id).n_r_box = box_int_recv.get(i).size();

478 479 480
			// For each received internal ghost box
			for (size_t j = 0 ; j < box_int_recv.get(i).size() ; j++)
			{
incardon's avatar
incardon committed
481 482 483
				size_t vol_recv = box_int_recv.get(i).get(j).bx.getVolumeKey();

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

incardon's avatar
incardon committed
486 487 488 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 541 542 543 544 545 546 547 548 549 550 551 552 553 554 555 556 557 558 559 560 561
				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;
						Box<dim,long int> flp_i = flip_box(box_int_recv.get(i).get(j).bx,box_int_recv.get(i).get(j).cmb);

						// 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
						Box<dim,long int> output = flip_box(box_int_recv.get(i).get(j).bx,box_int_recv.get(i).get(j).cmb);

						GBoxes<dim> tmp;
						tmp.GDbox = box_int_recv.get(i).get(j).bx;
						tmp.GDbox -= tmp.GDbox.getP1();
						tmp.origin = output.getP1();
562
						for (size_t s = 0 ; s < dim ; s++)
incardon's avatar
incardon committed
563 564
						{
							// we set an invalid box, there is no-domain
565 566
							tmp.Dbox.setLow(s,0);
							tmp.Dbox.setHigh(s,-1);
incardon's avatar
incardon committed
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 601 602 603 604 605 606 607 608 609 610 611 612 613 614 615 616 617 618 619 620 621 622 623 624 625 626 627 628 629 630 631
						}
						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();
					::Box<dim,long int> ib = flip_box(box_int_recv.get(i).get(j).bx,box_int_recv.get(i).get(j).cmb);
					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;
				}
632 633 634
			}
		}

Pietro Incardona's avatar
Pietro Incardona committed
635
		init_e_g_box = true;
636 637
	}

incardon's avatar
incardon committed
638 639 640 641 642
	/*! \brief Create local internal ghost box in grid units
	 *
	 */
	void create_local_ig_box()
	{
incardon's avatar
incardon committed
643 644
		openfpm::vector_std<result_box<dim>> ibv;

incardon's avatar
incardon committed
645 646 647 648 649
		// Get the grid info
		auto g = cd_sm.getGrid();

		if (init_local_i_g_box == true)	return;

Pietro Incardona's avatar
Pietro Incardona committed
650
		// Get the number of sub-domains
Pietro Incardona's avatar
Pietro Incardona committed
651
		for (size_t i = 0 ; i < dec.getNSubDomain() ; i++)
incardon's avatar
incardon committed
652 653 654 655 656 657
		{
			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
658 659 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 685 686 687 688 689 690 691 692 693 694 695 696 697 698 699 700 701 702 703 704 705 706 707 708 709 710 711 712 713 714 715 716 717 718 719 720 721
				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)
							continue;

						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().k.add(-1)/*dec.getLocalIGhostE(i,j)*/;
						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;

					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
722

incardon's avatar
incardon committed
723 724 725 726 727
						// 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
728

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

incardon's avatar
incardon committed
734
							pib.bid.last().box = output;
incardon's avatar
incardon committed
735

incardon's avatar
incardon committed
736 737
							pib.bid.last().sub_gdb_ext = j;
							pib.bid.last().sub = i;
incardon's avatar
incardon committed
738

incardon's avatar
incardon committed
739 740 741 742 743 744
//							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
745

incardon's avatar
incardon committed
746 747 748
						}
					}
				}
incardon's avatar
incardon committed
749
			}
incardon's avatar
incardon committed
750

incardon's avatar
incardon committed
751 752
		}

incardon's avatar
incardon committed
753

incardon's avatar
incardon committed
754 755 756
		init_local_i_g_box = true;
	}

757
	/*! \brief Create per-processor external ghost boxes list in grid units
incardon's avatar
incardon committed
758 759 760 761 762 763 764 765 766
	 *
	 */
	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
767
		loc_eg_box.resize(dec.getNSubDomain());
768 769 770 771 772 773

		// 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
774 775 776 777 778 779 780 781 782 783 784 785 786 787 788 789 790 791 792 793 794 795 796 797 798 799 800 801 802 803 804 805 806 807 808 809 810 811 812 813 814 815 816 817 818 819 820 821 822 823 824 825 826 827 828 829 830 831 832 833 834 835 836 837 838 839 840 841 842 843 844 845 846 847 848 849 850 851 852 853 854 855 856 857 858 859 860 861 862 863 864 865 866 867 868 869 870 871 872 873 874 875 876 877 878 879 880 881 882 883 884 885 886 887 888
				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;
						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);

						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();

							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);
							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,
										   flp_i,
										   cmb);


							volume_linked += pib.bid.last().box.getVolumeKey();
						}
					}

/*					if (volume_linked != loc_ig_box.get(i).bid.get(j).box.getVolumeKey())
					{
						// Create a grid with the same size of the external ghost
						// and mark all the linked points

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

						// Add an unlinked gdb_ext
						// An unlinked gdb_ext is an empty domain with only a ghost
						// part
						GBoxes<dim> tmp;
						tmp.GDbox = loc_ig_box.get(i).bid.get(j).box;
						tmp.GDbox -= tmp.GDbox.getP1();
						tmp.origin = loc_ig_box.get(i).bid.get(j).box.getP1();
						for (size_t i = 0 ; i < dim ; i++)
						{
							// we set an invalid box, there is no-domain
							tmp.Dbox.setLow(i,0);
							tmp.Dbox.setHigh(i,-1);
						}
						tmp.k = -1;
						gdb_ext.add(tmp);

						// create the local grid

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

						// Add an external ghost box

						Box<dim,long int> output = flip_box(loc_ig_box.get(i).bid.get(j).box,loc_ig_box.get(i).bid.get(j).cmb);

						// fill the link variable
						loc_ig_box.get(i).bid.get(j).k = pib.bid.size();

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


						add_loc_eg_box(le_sub,
									   dec.getLocalEGhostSub(le_sub,s),
									   j,
									   gdb_ext.size() - 1,
									   pib.bid,
									   output,
									   cmb);
					}*/
				}
				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));

					pib.bid.get(s).box = flip_box(loc_ig_box.get(i).bid.get(j).box,loc_ig_box.get(i).bid.get(j).cmb);
					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;
				}
889 890 891
			}
		}

Pietro Incardona's avatar
Pietro Incardona committed
892
		init_local_e_g_box = true;
incardon's avatar
incardon committed
893 894
	}

incardon's avatar
incardon committed
895

896
	/*! \brief Check the grid has a valid size
Pietro Incardona's avatar
Pietro Incardona committed
897 898
	 *
	 * \param g_sz size of the grid
incardon's avatar
incardon committed
899 900 901 902 903 904 905
	 *
	 */
	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
906
				std::cerr << "Error: " << __FILE__ << ":" << __LINE__ << " distributed grids with size smaller than 2 are not supported\n";
incardon's avatar
incardon committed
907 908
		}
	}
incardon's avatar
incardon committed
909

incardon's avatar
incardon committed
910
	/*! \brief Create the grids on memory
incardon's avatar
incardon committed
911 912 913
	 *
	 * \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
914 915
	 *
	 */
incardon's avatar
incardon committed
916 917 918
	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
919
	{
920
		// create gdb
incardon's avatar
incardon committed
921 922 923
		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();
924

incardon's avatar
incardon committed
925
		// create local grids for each hyper-cube
926
		loc_grid.resize(n_grid);
incardon's avatar
incardon committed
927 928 929 930 931 932 933

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

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

			// Get the size of the local grid
937 938 939 940
			// 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
941
			for (size_t j = 0 ; j < dim ; j++)
incardon's avatar
incardon committed
942
			{l_res[j] = (sp_tg.getHigh(j) >= 0)?(sp_tg.getHigh(j)+1):0;}
incardon's avatar
incardon committed
943 944 945 946 947 948 949

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


Pietro Incardona's avatar
Pietro Incardona committed
950 951 952 953 954 955
    /*! \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
956
	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
957 958 959 960 961
	{
		// Initialize the cell decomposer
		cd_sm.setDimensions(cd_old,ext);
	}

Pietro Incardona's avatar
Pietro Incardona committed
962
    /*! \brief Initialize the Cell decomposer of the grid
incardon's avatar
incardon committed
963
	 *
Pietro Incardona's avatar
Pietro Incardona committed
964 965
	 * \param g_sz Size of the grid
	 * \param bc boundary conditions
incardon's avatar
incardon committed
966 967
	 *
	 */
968
	inline void InitializeCellDecomposer(const size_t (& g_sz)[dim], const size_t (& bc)[dim])
incardon's avatar
incardon committed
969
	{
incardon's avatar
incardon committed
970
		// check that the grid has valid size
incardon's avatar
incardon committed
971 972
		check_size(g_sz);

973
		// get the size of the cell decomposer
incardon's avatar
incardon committed
974
		size_t c_g[dim];
975
		getCellDecomposerPar<dim>(c_g,g_sz,bc);
incardon's avatar
incardon committed
976 977 978

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

incardon's avatar
incardon committed
981 982 983
	/*! \brief Initialize the grid
	 *
	 * \param g_sz Global size of the grid
Pietro Incardona's avatar
Pietro Incardona committed
984
	 * \param bc boundary conditions
incardon's avatar
incardon committed
985 986
	 *
	 */
987
	inline void InitializeDecomposition(const size_t (& g_sz)[dim], const size_t (& bc)[dim])
incardon's avatar
incardon committed
988
	{
incardon's avatar
incardon committed
989
		// fill the global size of the grid
incardon's avatar
incardon committed
990
		for (size_t i = 0 ; i < dim ; i++)	{this->g_sz[i] = g_sz[i];}
incardon's avatar
incardon committed
991 992 993 994

		// 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
995
		size_t n_sub = n_proc * v_sub_unit_factor;
incardon's avatar
incardon committed
996 997 998 999

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

		// Create the sub-domains
Pietro Incardona's avatar
Pietro Incardona committed
1004
		dec.setParameters(div,domain,bc,ghost);
Pietro Incardona's avatar
Pietro Incardona committed
1005
		dec.decompose();
incardon's avatar
incardon committed
1006 1007
	}

incardon's avatar
incardon committed
1008 1009 1010 1011 1012 1013
	/*! \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
1014 1015 1016 1017 1018 1019 1020 1021 1022 1023 1024 1025 1026 1027 1028 1029 1030 1031 1032 1033 1034
	{
		// 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
1035 1036 1037 1038 1039
	{
		// 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
1040
		Create(bx,g,use_bx_def);
incardon's avatar
incardon committed
1041 1042
	}

incardon's avatar
incardon committed
1043 1044 1045
	// Ghost as integer
	Ghost<dim,long int> gint = Ghost<dim,long int>(0);

Pietro Incardona's avatar
Pietro Incardona committed
1046
protected:
incardon's avatar
incardon committed
1047 1048

	/*! \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
1049 1050
	 *
	 * \param i sub-domain
incardon's avatar
incardon committed
1051 1052 1053 1054 1055 1056 1057 1058 1059
	 *
	 * \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
1060 1061 1062 1063 1064 1065 1066 1067 1068 1069 1070 1071 1072 1073 1074 1075 1076 1077 1078 1079 1080
	/*! \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
1081 1082
			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
1083 1084 1085 1086 1087
		}

		return gc;
	}

incardon's avatar
incardon committed
1088 1089 1090 1091 1092 1093 1094 1095 1096 1097 1098 1099 1100 1101 1102 1103 1104 1105 1106 1107 1108 1109 1110 1111 1112 1113 1114
	/*! \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
1115 1116
public:

Pietro Incardona's avatar
Pietro Incardona committed
1117
	//! Which kind of grid the structure store
incardon's avatar
incardon committed
1118 1119
	typedef device_grid d_grid;

Pietro Incardona's avatar
Pietro Incardona committed
1120
	//! Decomposition used
incardon's avatar
incardon committed
1121 1122
	typedef Decomposition decomposition;

Pietro Incardona's avatar
Pietro Incardona committed
1123
	//! value_type
Pietro Incardona's avatar
Pietro Incardona committed
1124 1125
	typedef T value_type;

Pietro Incardona's avatar
Pietro Incardona committed
1126
	//! Type of space
1127 1128
	typedef St stype;

Pietro Incardona's avatar
Pietro Incardona committed
1129
	//! Type of Memory
1130 1131
	typedef Memory memory_type;

Pietro Incardona's avatar
Pietro Incardona committed
1132
	//! Type of device grid
1133 1134
	typedef device_grid device_grid_type;

Pietro Incardona's avatar
Pietro Incardona committed
1135
	//! Number of dimensions
1136 1137 1138 1139
	static const unsigned int dims = dim;

	/*! \brief Get the domain where the grid is defined
	 *
Pietro Incardona's avatar
Pietro Incardona committed
1140
	 * \return the domain of the grid
1141 1142 1143 1144 1145 1146 1147
	 *
	 */
	inline const Box<dim,St> getDomain() const
	{
		return domain;
	}

incardon's avatar
incardon committed
1148 1149 1150 1151 1152 1153 1154 1155 1156 1157 1158 1159
	/*! \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();
	}

1160
    /*! \brief Get the spacing of the grid in direction i
Pietro Incardona's avatar
Pietro Incardona committed
1161 1162
     *
     * \param i dimension
1163 1164 1165 1166 1167 1168 1169 1170 1171
     *
     * \return the spacing
     *
     */
    inline St spacing(size_t i) const
    {
    	return cd_sm.getCellBox().getHigh(i);
    }

Pietro Incardona's avatar
Pietro Incardona committed
1172 1173 1174 1175 1176 1177 1178 1179 1180 1181
	/*! \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
1182 1183 1184 1185 1186 1187 1188 1189 1190 1191 1192 1193
	/*! \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());}
	}

1194 1195 1196 1197 1198 1199 1200 1201 1202 1203 1204 1205 1206
	/*! \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
1207 1208 1209 1210 1211 1212 1213 1214 1215 1216 1217 1218 1219 1220 1221 1222 1223 1224
	/*! \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
1225 1226 1227 1228 1229 1230 1231 1232 1233
	/*! \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
	{
1234
		return ginfo_v.size(i);
Pietro Incardona's avatar
Pietro Incardona committed
1235 1236
	}

incardon's avatar
incardon committed
1237 1238 1239 1240 1241 1242 1243 1244 1245 1246 1247 1248 1249 1250 1251 1252 1253 1254 1255 1256 1257 1258 1259 1260 1261 1262 1263 1264 1265 1266 1267 1268 1269 1270
	/*! \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),
	 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
1271 1272 1273 1274 1275 1276
	/*! \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
1277
	 * \param gh Ghost part in grid units
Pietro Incardona's avatar
Pietro Incardona committed
1278 1279 1280
	 * \param ext extension of the grid (must be positive on every direction)
	 *
	 */
incardon's avatar
incardon committed
1281 1282 1283 1284 1285
	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
1286 1287 1288 1289 1290
	{
#ifdef SE_CLASS2
		check_new(this,8,GRID_DIST_EVENT,4);
#endif

1291 1292 1293 1294 1295 1296 1297 1298 1299 1300
		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
1301 1302 1303 1304 1305

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

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

1308 1309 1310 1311 1312 1313 1314 1315 1316 1317
			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
1318 1319
		}

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

incardon's avatar
incardon committed
1322 1323 1324 1325
		// an empty
		openfpm::vector<Box<dim,long int>> empty;

		InitializeStructures(g.getGridInfoVoid().getSize(),empty,gh,false);
Pietro Incardona's avatar
Pietro Incardona committed
1326
	}
incardon's avatar
incardon committed
1327

1328 1329 1330 1331 1332 1333 1334
    /*! It constructs a grid of a specified size, defined on a specified Box space, forcing to follow a specified decomposition and with a specified ghost size
     *
     * \param dec Decomposition
     * \param g_sz grid size on each dimension
     * \param ghost Ghost part
     *
     */
incardon's avatar