vector_dist.hpp 32.9 KB
Newer Older
incardon's avatar
incardon committed
1
2
3
4
5
6
7
8
9
10
/*
 * Vector.hpp
 *
 *  Created on: Mar 5, 2015
 *      Author: Pietro Incardona
 */

#ifndef VECTOR_HPP_
#define VECTOR_HPP_

incardon's avatar
incardon committed
11
#include "VCluster.hpp"
incardon's avatar
incardon committed
12
#include "Space/Shape/Point.hpp"
incardon's avatar
incardon committed
13
14
15
#include "Vector/vector_dist_iterator.hpp"
#include "Space/Shape/Box.hpp"
#include "Vector/vector_dist_key.hpp"
incardon's avatar
incardon committed
16
17
#include "memory/PreAllocHeapMemory.hpp"
#include "memory/PtrMemory.hpp"
incardon's avatar
incardon committed
18
#include "NN/CellList/CellList.hpp"
19
#include "util/common.hpp"
incardon's avatar
incardon committed
20
#include "util/object_util.hpp"
incardon's avatar
incardon committed
21
#include "memory/ExtPreAlloc.hpp"
incardon's avatar
incardon committed
22
#include "CSVWriter.hpp"
23
#include "Decomposition/common.hpp"
Pietro Incardona's avatar
Pietro Incardona committed
24
#include "Vector/vector_dist_ofb.hpp"
25
26

#define V_SUB_UNIT_FACTOR 64
incardon's avatar
incardon committed
27
28
29
30

#define NO_ID false
#define ID true

incardon's avatar
incardon committed
31
// Perform a ghost get or a ghost put
incardon's avatar
incardon committed
32
33
34
#define GET	1
#define PUT 2

incardon's avatar
incardon committed
35
#define NO_POSITION 1
36
#define WITH_POSITION 2
incardon's avatar
incardon committed
37

incardon's avatar
incardon committed
38
39
40
41
// Write the particles with ghost
#define NO_GHOST 0
#define WITH_GHOST 2

incardon's avatar
incardon committed
42
/*! \brief Distributed vector
43
 *
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
 * This class reppresent a distributed vector, the distribution of the structure
 * is based on the positional information of the elements the vector store
 *
 * ## Create a vector of random elements on each processor 2D
 * \snippet vector_dist_unit_test.hpp Create a vector of random elements on each processor 2D
 *
 * ## Create a vector of random elements on each processor 3D
 * \snippet vector_dist_unit_test.hpp Create a vector of random elements on each processor 3D
 *
 * ## Create a vector of elements distributed on a grid like way
 * \snippet vector_dist_unit_test.hpp Create a vector of elements distributed on a grid like way
 *
 * ## Redistribute the particles and sync the ghost properties
 * \snippet vector_dist_unit_test.hpp Redistribute the particles and sync the ghost properties
 *
 * \tparam dim Dimensionality of the space where the elements lives
 * \tparam St type of space float, double ...
 * \tparam prop properties the vector element store in OpenFPM data structure format
 * \tparam Decomposition Decomposition strategy to use CartDecomposition ...
 * \tparam Memory Memory pool where store the information HeapMemory ...
incardon's avatar
incardon committed
64
65
66
 *
 */

67
template<unsigned int dim, typename St, typename prop, typename Decomposition , typename Memory=HeapMemory>
incardon's avatar
incardon committed
68
69
70
71
class vector_dist
{
private:

72
	//! Ghost marker, all the particle with id > g_m are ghost all with g_m < are real particle
73
74
	size_t g_m = 0;

incardon's avatar
incardon committed
75
76
77
	//! Space Decomposition
	Decomposition dec;

78
79
	//! Particle position vector, (It has 2 elements) the first has real particles assigned to a processor
	//! the second element contain unassigned particles
80
	openfpm::vector<Point<dim,St>> v_pos;
incardon's avatar
incardon committed
81

82
83
	//! Particle properties vector, (It has 2 elements) the first has real particles assigned to a processor
	//! the second element contain unassigned particles
84
	openfpm::vector<prop> v_prp;
incardon's avatar
incardon committed
85

86
	//! Virtual cluster
incardon's avatar
incardon committed
87
88
	Vcluster & v_cl;

89
	// definition of the send vector for position
90
	typedef  openfpm::vector<Point<dim,St>,ExtPreAlloc<Memory>,openfpm::grow_policy_identity> send_pos_vector;
91
92
93
94
95
96
97
98
99

	//////////////////////////////
	// COMMUNICATION variables
	//////////////////////////////

	//! It map the processor id with the communication request into map procedure
	openfpm::vector<size_t> p_map_req;

	//! For each near processor, outgoing particle id and shift vector
100
101
102
103
	openfpm::vector<openfpm::vector<size_t>> opart;

	//! For each near processor, particle shift vector
	openfpm::vector<openfpm::vector<size_t>> oshift;
104
105
106
107
108
109
110
111
112
113
114
115
116

	//! For each adjacent processor the size of the ghost sending buffer
	openfpm::vector<size_t> ghost_prc_sz;

	//! Sending buffer for the ghost particles properties
	Memory g_prp_mem;

	//! Sending buffer for the ghost particles position
	Memory g_pos_mem;

	//! For each adjacent processor it store the size of the receiving message in byte
	openfpm::vector<size_t> recv_sz;

117
	//! For each adjacent processor it store the received message for ghost get
118
119
	openfpm::vector<HeapMemory> recv_mem_gg;

120
121
	//! For each processor it store the received message for global map
	openfpm::vector<HeapMemory> recv_mem_gm;
122

123
124
125
126
127
128
129
130
131
132
133
	/*! \brief It store for each processor the position and properties vector of the particles
	 *
	 *
	 */
	struct pos_prop
	{
		//! position vector
		openfpm::vector<Point<dim,St>,PreAllocHeapMemory<2>,openfpm::grow_policy_identity> pos;
		//! properties vector
		openfpm::vector<prop,PreAllocHeapMemory<2>,openfpm::grow_policy_identity> prp;
	};
134

135
136
137
	/*! \brief Label particles for mappings
	 *
	 * \param lbl_p Particle labeled
138
	 * \param prc_sz For each processor the number of particles to send
139
140
141
	 * \param opart id of the particles to send
	 *
	 */
Pietro Incardona's avatar
Pietro Incardona committed
142
	template<typename obp> void labelParticleProcessor(openfpm::vector<openfpm::vector<size_t>> & lbl_p, openfpm::vector<size_t> & prc_sz, openfpm::vector<size_t> & opart)
143
	{
144
145
146
147
148
		// reset lbl_p
		lbl_p.resize(v_cl.getProcessingUnits());
		for (size_t i = 0 ; i < lbl_p.size() ; i++)
			lbl_p.get(i).clear();

149
150
151
152
153
154
155
156
157
158
159
160
161
		// resize the label buffer
		prc_sz.resize(v_cl.getProcessingUnits());

		auto it = v_pos.getIterator();

		// Label all the particles with the processor id where they should go
		while (it.isNext())
		{
			auto key = it.get();

			// Apply the boundary conditions
			dec.applyPointBC(v_pos.get(key));

Pietro Incardona's avatar
Pietro Incardona committed
162
163
164
			size_t p_id = 0;

			// Check if the particle is inside the domain
165
			if (dec.getDomain().isInside(v_pos.get(key)) == true)
Pietro Incardona's avatar
Pietro Incardona committed
166
				p_id = dec.processorIDBC(v_pos.get(key));
167
168
			else
				p_id = obp::out(key,v_cl.getProcessUnitID());
Pietro Incardona's avatar
Pietro Incardona committed
169

170

171
			// Particle to move
172
173
			if (p_id != v_cl.getProcessUnitID())
			{
Pietro Incardona's avatar
Pietro Incardona committed
174
				if ((long int)p_id != -1)
175
176
177
178
179
180
181
182
183
				{
					prc_sz.get(p_id)++;
					lbl_p.get(p_id).add(key);
					opart.add(key);
				}
				else
				{
					opart.add(key);
				}
184
185
186
187
188
189
190
191
			}

			// Add processors and add size

			++it;
		}
	}

192
193
194
195
	/*! \brief Label the particles
	 *
	 * It count the number of particle to send to each processors and save its ids
	 *
196
197
198
199
200
201
	 * \param prc_sz For each processor the number of particles to send
	 * \param opart id if of the particles to send
	 * \param shift_id shift correction id
	 *
	 * \see nn_prcs::getShiftvectors()
	 *
202
	 */
203
	void labelParticlesGhost()
204
205
206
207
208
209
210
211
212
	{
		// Buffer that contain the number of elements to send for each processor
		ghost_prc_sz.clear();
		ghost_prc_sz.resize(dec.getNNProcessors());

		// Buffer that contain for each processor the id of the particle to send
		opart.clear();
		opart.resize(dec.getNNProcessors());

213
214
215
216
		// Buffer that contain for each processor the id of the shift vector
		oshift.clear();
		oshift.resize(dec.getNNProcessors());

217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
		// Iterate over all particles
		auto it = v_pos.getIterator();
		while (it.isNext())
		{
			auto key = it.get();

			// Given a particle, it return which processor require it (first id) and shift id, second id
			// For an explanation about shifts vectors please consult getShiftVector in ie_ghost
			const openfpm::vector<std::pair<size_t,size_t>> & vp_id = dec.template ghost_processorID_pair<typename Decomposition::lc_processor_id, typename Decomposition::shift_id>(v_pos.get(key),UNIQUE);

			for (size_t i = 0 ; i < vp_id.size() ; i++)
			{
				// processor id
				size_t p_id = vp_id.get(i).first;

				// add particle to communicate
				ghost_prc_sz.get(p_id)++;

235
236
				opart.get(p_id).add(key);
				oshift.get(p_id).add(vp_id.get(i).second);
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
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
284
285
286
287
288
289
290
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
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
			}

			++it;
		}
	}

	/*! \brief Add local particles based on the boundary conditions
	 *
	 * In order to understand what this function use the following
	 *
		\verbatim

															[1,1]
			+---------+------------------------+---------+
			| (1,-1)  |                        | (1,1)   |
			|   |     |    (1,0) --> 7         |   |     |
			|   v     |                        |   v     |
			|   6     |                        |   8     |
			+--------------------------------------------+
			|         |                        |         |
			|         |                        |         |
			|         |                        |         |
			| (-1,0)  |                        | (1,0)   |
			|    |    |                        |   |     |
			|    v    |      (0,0) --> 4       |   v     |
			|    3    |                        |   5     |
			|         |                        |         |
		 B	|         |                        |     A   |
		*	|         |                        |    *    |
			|         |                        |         |
			|         |                        |         |
			|         |                        |         |
			+--------------------------------------------+
			| (-1,-1) |                        | (-1,1)  |
			|    |    |   (-1,0) --> 1         |    |    |
			|    v    |                        |    v    |
			|    0    |                        |    2    |
			+---------+------------------------+---------+


		\endverbatim

	 *
	 *  The box is the domain, while all boxes at the border (so not (0,0) ) are the
	 *  ghost part at the border of the domain. If a particle A is in the position in figure
	 *  a particle B must be created. This function duplicate the particle A, if A and B are
	 *  local
	 *
	 */
	void add_loc_particles_bc()
	{
		// get the shift vectors
		const openfpm::vector<Point<dim,St>> & shifts = dec.getShiftVectors();

		// Add local particles coming from periodic boundary, the only boxes that count are the one
		// at the border, filter them

		openfpm::vector_std<Box<dim,St>> box_f;
		openfpm::vector_std<comb<dim>> box_cmb;

		for (size_t i = 0 ; i < dec.getNLocalSub() ; i++)
		{
			size_t Nl = dec.getLocalNIGhost(i);

			for (size_t j = 0 ; j < Nl ; j++)
			{
				// If they are not in the border the combination is all zero
				if (dec.getLocalIGhostPos(i,j).n_zero() == dim)
					continue;

				box_f.add(dec.getLocalIGhostBox(i,j));
				box_cmb.add(dec.getLocalIGhostPos(i,j));
			}
		}


		if (box_f.size() == 0)
			return;
		else
		{
			// Label the internal (assigned) particles
			auto it = v_pos.getIterator();

			while (it.isNext())
			{
				auto key = it.get();

				// If particles are inside these boxes
				for (size_t i = 0 ; i < box_f.size() ; i++)
				{
					if (box_f.get(i).isInside(v_pos.get(key)) == true)
					{
						Point<dim,St> p = v_pos.get(key);
						// shift
						p -= shifts.get(box_cmb.get(i).lin());

						// add this particle shifting its position
						v_pos.add(p);
						v_prp.add();
						v_prp.last() = v_prp.get(key);
					}
				}

				++it;
			}
		}
	}

	/*! \brief This function fill the send buffer for the particle position after the particles has been label with labelParticles
	 *
	 * \param g_pos_send Send buffer to fill
	 * \param prAlloc_pos Memory object for the send buffer
	 *
	 */
351
	void fill_send_ghost_pos_buf(openfpm::vector<send_pos_vector> & g_pos_send, ExtPreAlloc<Memory> * prAlloc_pos)
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
	{
		// get the shift vectors
		const openfpm::vector<Point<dim,St>> & shifts = dec.getShiftVectors();

		// create a number of send buffers equal to the near processors
		g_pos_send.resize(ghost_prc_sz.size());
		for (size_t i = 0 ; i < g_pos_send.size() ; i++)
		{
			// set the preallocated memory to ensure contiguity
			g_pos_send.get(i).setMemory(*prAlloc_pos);

			// resize the sending vector (No allocation is produced)
			g_pos_send.get(i).resize(ghost_prc_sz.get(i));
		}

		// Fill the send buffer
		for ( size_t i = 0 ; i < opart.size() ; i++ )
		{
			for (size_t j = 0 ; j < opart.get(i).size() ; j++)
			{
372
373
				Point<dim,St> s = v_pos.get(opart.get(i).get(j));
				s -= shifts.get(oshift.get(i).get(j));
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
				g_pos_send.get(i).set(j,s);
			}
		}
	}

	/*! \brief This function fill the send buffer for properties after the particles has been label with labelParticles
	 *
	 * \tparam send_vector type used to send data
	 * \tparam prp_object object containing only the properties to send
	 * \tparam prp set of properties to send
	 *
	 * \param g_send_prp Send buffer to fill
	 * \param prAlloc_prp Memory object for the send buffer
	 *
	 */
389
	template<typename send_vector,typename prp_object, int... prp> void fill_send_ghost_prp_buf(openfpm::vector<send_vector> & g_send_prp, ExtPreAlloc<Memory> * prAlloc_prp)
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
	{
		// create a number of send buffers equal to the near processors
		g_send_prp.resize(ghost_prc_sz.size());
		for (size_t i = 0 ; i < g_send_prp.size() ; i++)
		{
			// set the preallocated memory to ensure contiguity
			g_send_prp.get(i).setMemory(*prAlloc_prp);

			// resize the sending vector (No allocation is produced)
			g_send_prp.get(i).resize(ghost_prc_sz.get(i));
		}

		// Fill the send buffer
		for ( size_t i = 0 ; i < opart.size() ; i++ )
		{
			for (size_t j = 0 ; j < opart.get(i).size() ; j++)
			{
				// source object type
				typedef encapc<1,prop,typename openfpm::vector<prop>::memory_conf> encap_src;
				// destination object type
				typedef encapc<1,prp_object,typename openfpm::vector<prp_object>::memory_conf> encap_dst;

				// Copy only the selected properties
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
				object_si_d<encap_src,encap_dst,OBJ_ENCAP,prp...>(v_prp.get(opart.get(i).get(j)),g_send_prp.get(i).get(j));
			}
		}
	}

	/*! \brief allocate and fill the send buffer for the map function
	 *
	 * \param prc_r List of processor rank involved in the send
	 * \param prc_r_sz For each processor in the list the size of the message to send
	 * \param pb send buffer
	 *
	 */
	void fill_send_map_buf(openfpm::vector<size_t> & prc_r, openfpm::vector<size_t> & prc_sz_r, openfpm::vector<pos_prop> & pb)
	{
		pb.resize(prc_r.size());

		for (size_t i = 0 ;  i < prc_r.size() ; i++)
		{
			// Create the size required to store the particles position and properties to communicate
			size_t s1 = openfpm::vector<Point<dim,St>,HeapMemory,openfpm::grow_policy_identity>::calculateMem(prc_sz_r.get(i),0);
			size_t s2 = openfpm::vector<prop,HeapMemory,openfpm::grow_policy_identity>::calculateMem(prc_sz_r.get(i),0);

			// Preallocate the memory
			size_t sz[2] = {s1,s2};
			PreAllocHeapMemory<2> * mem = new PreAllocHeapMemory<2>(sz);

			// Set the memory allocator
			pb.get(i).pos.setMemory(*mem);
			pb.get(i).prp.setMemory(*mem);

			// set the size and allocate, using mem warant that pos and prp is contiguous
			pb.get(i).pos.resize(prc_sz_r.get(i));
			pb.get(i).prp.resize(prc_sz_r.get(i));
		}


		// Run through all the particles and fill the sending buffer

		for (size_t i = 0 ; i < opart.size() ; i++)
		{
			auto it = opart.get(i).getIterator();
			size_t lbl = p_map_req.get(i);

			while (it.isNext())
			{
				size_t key = it.get();
				size_t id = opart.get(i).get(key);

				pb.get(lbl).pos.set(key,v_pos.get(id));
				pb.get(lbl).prp.set(key,v_prp.get(id));

				++it;
465
466
467
468
469
470
471
472
473
474
475
			}
		}
	}

	/*! \brief This function process the receiced data for the properties and populate the ghost
	 *
	 * \tparam send_vector type used to send data
	 * \tparam prp_object object containing only the properties to send
	 * \tparam prp set of properties to send
	 *
	 */
476
	template<typename send_vector,typename prp_object, int... prp> void process_received_ghost_prp()
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
	{
		// Mark the ghost part
		g_m = v_prp.size();

		// Process the received data (recv_mem_gg != 0 if you have data)
		for (size_t i = 0 ; i < dec.getNNProcessors() && recv_mem_gg.size() != 0 ; i++)
		{
			// calculate the number of received elements
			size_t n_ele = recv_sz.get(i) / sizeof(prp_object);

			// add the received particles to the vector
			PtrMemory * ptr1 = new PtrMemory(recv_mem_gg.get(i).getPointer(),recv_sz.get(i));

			// create vector representation to a piece of memory already allocated
			openfpm::vector<prp_object,PtrMemory,openfpm::grow_policy_identity> v2;

			v2.setMemory(*ptr1);

			// resize with the number of elements
			v2.resize(n_ele);

			// Add the ghost particle
			v_prp.template add_prp<prp_object,PtrMemory,openfpm::grow_policy_identity,prp...>(v2);
		}
	}

	/*! \brief This function process the received data for the properties and populate the ghost
	 *
	 */
506
	void process_received_ghost_pos()
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
	{
		// Process the received data (recv_mem_gg != 0 if you have data)
		for (size_t i = 0 ; i < dec.getNNProcessors() && recv_mem_gg.size() != 0 ; i++)
		{
			// calculate the number of received elements
			size_t n_ele = recv_sz.get(i) / sizeof(Point<dim,St>);

			// add the received particles to the vector
			PtrMemory * ptr1 = new PtrMemory(recv_mem_gg.get(i).getPointer(),recv_sz.get(i));

			// create vector representation to a piece of memory already allocated

			openfpm::vector<Point<dim,St>,PtrMemory,openfpm::grow_policy_identity> v2;

			v2.setMemory(*ptr1);

			// resize with the number of elements
			v2.resize(n_ele);

			// Add the ghost particle
			v_pos.template add<PtrMemory,openfpm::grow_policy_identity>(v2);
		}
	}

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
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
	/*! \brief Process the received particles
	 *
	 * \param list of the out-going particles
	 *
	 */
	void process_received_map(openfpm::vector<size_t> & out_part)
	{
		size_t o_p_id = 0;

		for (size_t i = 0 ; i < recv_mem_gm.size() ; i++)
		{
			// Get the number of elements

			size_t n_ele = recv_mem_gm.get(i).size() / (sizeof(Point<dim,St>) + sizeof(prop));

			// Pointer of the received positions for each near processor
			void * ptr_pos = (unsigned char *)recv_mem_gm.get(i).getPointer();
			// Pointer of the received properties for each near processor
			void * ptr_prp = (unsigned char *)recv_mem_gm.get(i).getPointer() + n_ele * sizeof(Point<dim,St>);

			PtrMemory * ptr1 = new PtrMemory(ptr_pos,n_ele * sizeof(Point<dim,St>));
			PtrMemory * ptr2 = new PtrMemory(ptr_prp,n_ele * sizeof(prop));

			// create vector representation to a piece of memory already allocated

			openfpm::vector<Point<dim,St>,PtrMemory,openfpm::grow_policy_identity> vpos;
			openfpm::vector<prop,PtrMemory,openfpm::grow_policy_identity> vprp;

			vpos.setMemory(*ptr1);
			vprp.setMemory(*ptr2);

			vpos.resize(n_ele);
			vprp.resize(n_ele);

			// Add the received particles to v_pos and v_prp

			size_t j = 0;
			for ( ; j < vpos.size() && o_p_id < out_part.size() ; j++, o_p_id++)
			{
				v_pos.set(out_part.get(o_p_id),vpos.get(j));
				v_prp.set(out_part.get(o_p_id),vprp.get(j));
			}

			for ( ; j < vpos.size(); j++)
			{
				v_pos.add();
				v_pos.set(v_pos.size()-1,vpos.get(j));
				v_prp.add();
				v_prp.set(v_prp.size()-1,vprp.get(j));
			}
		}

		// remove the (out-going particles) in the vector

		v_pos.remove(out_part,o_p_id);
		v_prp.remove(out_part,o_p_id);
	}

589
590
591
592
593
594
595
596
597
598
	/*! \brief Calculate send buffers total size and allocation
	 *
	 * \tparam prp_object object containing only the properties to send
	 *
	 * \param size_byte_prp total size for the property buffer
	 * \param size_byte_pos total size for the position buffer
	 * \param pap_prp allocation sequence for the property buffer
	 * \param pap_pos allocation sequence for the position buffer
	 *
	 */
599
	template<typename prp_object> void calc_send_ghost_buf(size_t & size_byte_prp, size_t & size_byte_pos, std::vector<size_t> & pap_prp, std::vector<size_t> & pap_pos)
600
601
602
603
	{
		// Calculate the total size required for the sending buffer
		for ( size_t i = 0 ; i < ghost_prc_sz.size() ; i++ )
		{
604
			size_t alloc_ele = openfpm::vector<prp_object,HeapMemory,openfpm::grow_policy_identity>::calculateMem(ghost_prc_sz.get(i),0);
605
606
607
			pap_prp.push_back(alloc_ele);
			size_byte_prp += alloc_ele;

608
			alloc_ele = openfpm::vector<Point<dim,St>,HeapMemory,openfpm::grow_policy_identity>::calculateMem(ghost_prc_sz.get(i),0);
609
610
611
612
613
			pap_pos.push_back(alloc_ele);
			size_byte_pos += alloc_ele;
		}
	}

incardon's avatar
incardon committed
614
615
616
617
public:

	/*! \brief Constructor
	 *
618
619
	 * \param np number of elements
	 * \param box domain where the vector of elements live
620
	 * \param boundary conditions
621
	 * \param g Ghost margins
incardon's avatar
incardon committed
622
623
	 *
	 */
624
	vector_dist(size_t np, Box<dim,St> box, const size_t (& bc)[dim] ,const Ghost<dim,St> & g)
625
	:dec(*global_v_cluster),v_cl(*global_v_cluster)
incardon's avatar
incardon committed
626
	{
627
628
629
630
#ifdef SE_CLASS2
		check_new(this,8,VECTOR_DIST_EVENT,4);
#endif

631
		// convert to a local number of elements
632
633
634
635
636
637
638
639
		size_t p_np = np / v_cl.getProcessingUnits();

		// Get non divisible part
		size_t r = np % v_cl.getProcessingUnits();

		// Distribute the remain particles
		if (v_cl.getProcessUnitID() < r)
			p_np++;
640

incardon's avatar
incardon committed
641
		// resize the position vector
642
		v_pos.resize(p_np);
incardon's avatar
incardon committed
643
644

		// resize the properties vector
645
646
647
		v_prp.resize(p_np);

		g_m = p_np;
incardon's avatar
incardon committed
648
649
650
651
652

		// Create a valid decomposition of the space
		// Get the number of processor and calculate the number of sub-domain
		// for decomposition
		size_t n_proc = v_cl.getProcessingUnits();
653
		size_t n_sub = n_proc * getDefaultNsubsub();
incardon's avatar
incardon committed
654
655
656

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

		// Create the sub-domains
Pietro Incardona's avatar
Pietro Incardona committed
662
		dec.setParameters(div,box,bc,g);
incardon's avatar
incardon committed
663

664
665
		// and create the ghost boxes
		dec.calculateGhostBoxes();
incardon's avatar
incardon committed
666
667
	}

668
669
670
671
672
673
674
	~vector_dist()
	{
#ifdef SE_CLASS2
		check_delete(this);
#endif
	}

675
676
677
678
679
680
681
682
683
684
	/*! \brief Get the number of minimum sub-domain
	 *
	 * \return minimum number
	 *
	 */
	static size_t getDefaultNsubsub()
	{
		return  V_SUB_UNIT_FACTOR;
	}

incardon's avatar
incardon committed
685
686
687
688
689
690
691
	/*! \brief return the local size of the vector
	 *
	 * \return local size
	 *
	 */
	size_t size_local()
	{
Pietro Incardona's avatar
Pietro Incardona committed
692
		return g_m;
incardon's avatar
incardon committed
693
694
	}

695
	/*! \brief Get the position of an element
incardon's avatar
incardon committed
696
	 *
697
698
699
700
701
	 * see the vector_dist iterator usage to get an element key
	 *
	 * \param vec_key element
	 *
	 * \return the position of the element in space
incardon's avatar
incardon committed
702
703
	 *
	 */
704
	template<unsigned int id> inline auto getPos(vect_dist_key_dx vec_key) -> decltype(v_pos.template get<id>(vec_key.getKey()))
incardon's avatar
incardon committed
705
	{
706
		return v_pos.template get<id>(vec_key.getKey());
incardon's avatar
incardon committed
707
708
	}

709
	/*! \brief Get the property of an element
incardon's avatar
incardon committed
710
	 *
711
712
713
	 * see the vector_dist iterator usage to get an element key
	 *
	 * \tparam id property id
incardon's avatar
incardon committed
714
715
	 * \param vec_key vector element
	 *
716
717
	 * \return return the selected property of the vector element
	 *
incardon's avatar
incardon committed
718
	 */
719
	template<unsigned int id> inline auto getProp(vect_dist_key_dx vec_key) -> decltype(v_prp.template get<id>(vec_key.getKey()))
incardon's avatar
incardon committed
720
	{
721
		return v_prp.template get<id>(vec_key.getKey());
incardon's avatar
incardon committed
722
723
	}

724
	/*! \brief It move all the particles that does not belong to the local processor to the respective processor
Pietro Incardona's avatar
Pietro Incardona committed
725
726
	 *
	 * \tparam out of bound policy it specify what to do when the particles are detected out of bound
727
728
729
	 *
	 * In general this function is called after moving the particles to move the
	 * elements out the local processor. Or just after initialization if each processor
730
	 * contain non local particles
incardon's avatar
incardon committed
731
732
	 *
	 */
Pietro Incardona's avatar
Pietro Incardona committed
733
	template<typename obp=KillParticle> void map()
incardon's avatar
incardon committed
734
	{
incardon's avatar
incardon committed
735
		// outgoing particles-id
736
		openfpm::vector<size_t> out_part;
incardon's avatar
incardon committed
737
738
739
740
741

		// Processor communication size
		openfpm::vector<size_t> prc_sz(v_cl.getProcessingUnits());

		// It contain the list of the processors this processor should to communicate with
incardon's avatar
incardon committed
742
743
		openfpm::vector<size_t> p_list;

744
745
746
747
748
		// map completely reset the ghost part
		v_pos.resize(g_m);
		v_prp.resize(g_m);

		// Contain the processor id of each particle (basically where they have to go)
Pietro Incardona's avatar
Pietro Incardona committed
749
		labelParticleProcessor<obp>(opart,prc_sz,out_part);
incardon's avatar
incardon committed
750

751
752
		// Calculate the sending buffer size for each processor, put this information in
		// a contiguous buffer
753
		p_map_req.resize(v_cl.getProcessingUnits());
incardon's avatar
incardon committed
754
755
756
757
758
		openfpm::vector<size_t> prc_sz_r;
		openfpm::vector<size_t> prc_r;

		for (size_t i = 0 ; i < v_cl.getProcessingUnits() ; i++)
		{
759
			if (prc_sz.get(i) != 0)
incardon's avatar
incardon committed
760
			{
761
				p_map_req.get(i) = prc_r.size();
incardon's avatar
incardon committed
762
763
764
765
766
				prc_r.add(i);
				prc_sz_r.add(prc_sz.get(i));
			}
		}

767
		// Allocate the send buffers
incardon's avatar
incardon committed
768

769
		openfpm::vector<pos_prop> pb;
incardon's avatar
incardon committed
770

771
772
		// fill the send buffers
		fill_send_map_buf(prc_r,prc_sz_r,pb);
incardon's avatar
incardon committed
773
774
775
776
777
778
779
780
781

		// Create the set of pointers
		openfpm::vector<void *> ptr(prc_r.size());
		for (size_t i = 0 ; i < prc_r.size() ; i++)
		{
			ptr.get(i) = pb.get(i).pos.getPointer();
		}

		// convert the particle number to buffer size
incardon's avatar
incardon committed
782
		for (size_t i = 0 ; i < prc_sz_r.size() ; i++)
incardon's avatar
incardon committed
783
		{
784
			prc_sz_r.get(i) = prc_sz_r.get(i)*(sizeof(prop) + sizeof(Point<dim,St>));
incardon's avatar
incardon committed
785
786
787
788
		}

		// Send and receive the particles

789
790
		recv_mem_gm.clear();
		v_cl.sendrecvMultipleMessagesNBX(prc_sz_r.size(), (size_t *)prc_sz_r.getPointer(), (size_t *)prc_r.getPointer() , (void **)ptr.getPointer() , vector_dist::message_alloc_map, this ,NEED_ALL_SIZE);
incardon's avatar
incardon committed
791

792
		// Process the incoming particles
incardon's avatar
incardon committed
793

794
		process_received_map(out_part);
incardon's avatar
incardon committed
795

796
		// mark the ghost part
incardon's avatar
incardon committed
797

798
799
		g_m = v_pos.size();
	}
incardon's avatar
incardon committed
800

801
	/*! \brief It synchronize the properties and position of the ghost particles
incardon's avatar
incardon committed
802
	 *
803
804
	 * \tparam prp list of properties to get synchronize
	 * \param opt options WITH_POSITION, it send also the positional information of the particles
incardon's avatar
incardon committed
805
806
	 *
	 */
807
	template<int... prp> void ghost_get(size_t opt = WITH_POSITION)
incardon's avatar
incardon committed
808
	{
809
810
811
812
		// Sending property object
		typedef object<typename object_creator<typename prop::type,prp...>::type> prp_object;

		// send vector for each processor
813
		typedef  openfpm::vector<prp_object,ExtPreAlloc<Memory>,openfpm::grow_policy_identity> send_vector;
814

815
816
817
		// reset the ghost part
		v_pos.resize(g_m);
		v_prp.resize(g_m);
incardon's avatar
incardon committed
818

819
		// Label all the particles
820
		labelParticlesGhost();
incardon's avatar
incardon committed
821

822
		// Calculate memory and allocation for the send buffers
incardon's avatar
incardon committed
823

824
		// Total size
incardon's avatar
incardon committed
825
826
		size_t size_byte_prp = 0;
		size_t size_byte_pos = 0;
incardon's avatar
incardon committed
827

828
		// allocation patterns for property and position send buffer
incardon's avatar
incardon committed
829
830
		std::vector<size_t> pap_prp;
		std::vector<size_t> pap_pos;
incardon's avatar
incardon committed
831

832
		calc_send_ghost_buf<prp_object>(size_byte_prp,size_byte_pos,pap_prp,pap_pos);
incardon's avatar
incardon committed
833

834
		// Create memory for the send buffer
incardon's avatar
incardon committed
835

incardon's avatar
incardon committed
836
837
		g_prp_mem.resize(size_byte_prp);
		if (opt != NO_POSITION) g_pos_mem.resize(size_byte_pos);
incardon's avatar
incardon committed
838

839
		// Create and fill send buffer for particle properties
incardon's avatar
incardon committed
840

841
		ExtPreAlloc<Memory> * prAlloc_prp = new ExtPreAlloc<Memory>(pap_prp,g_prp_mem);
incardon's avatar
incardon committed
842
		openfpm::vector<send_vector> g_send_prp;
843
		fill_send_ghost_prp_buf<send_vector,prp_object,prp...>(g_send_prp,prAlloc_prp);
incardon's avatar
incardon committed
844

845
		// Create and fill the send buffer for the particle position
incardon's avatar
incardon committed
846

847
		ExtPreAlloc<Memory> * prAlloc_pos;
incardon's avatar
incardon committed
848
		openfpm::vector<send_pos_vector> g_pos_send;
incardon's avatar
incardon committed
849
850
		if (opt != NO_POSITION)
		{
851
			prAlloc_pos = new ExtPreAlloc<Memory>(pap_pos,g_pos_mem);
852
			fill_send_ghost_pos_buf(g_pos_send,prAlloc_pos);
incardon's avatar
incardon committed
853
854
		}

855
		// Create processor list
incardon's avatar
incardon committed
856
857
858
		openfpm::vector<size_t> prc;
		for (size_t i = 0 ; i < opart.size() ; i++)
			prc.add(dec.IDtoProc(i));
incardon's avatar
incardon committed
859

860
		// Send/receive the particle properties information
861
		v_cl.sendrecvMultipleMessagesNBX(prc,g_send_prp,msg_alloc_ghost_get,this);
862
		process_received_ghost_prp<send_vector,prp_object,prp...>();
863
864
865
866


		if (opt != NO_POSITION)
		{
867
			// Send/receive the particle properties information
868
			v_cl.sendrecvMultipleMessagesNBX(prc,g_pos_send,msg_alloc_ghost_get,this);
869
			process_received_ghost_pos();
870
871
		}

872
873
		add_loc_particles_bc();
	}
874

875
	/*! \brief Call-back to allocate buffer to receive incoming elements (particles)
876
877
878
879
880
881
882
883
884
885
886
887
	 *
	 * \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_ghost_get(size_t msg_i ,size_t total_msg, size_t total_p, size_t i, size_t ri, void * ptr)
	{
888
		vector_dist<dim,St,prop,Decomposition,Memory> * v = static_cast<vector_dist<dim,St,prop,Decomposition,Memory> *>(ptr);
889
890
891
892
893
894

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

		// Get the local processor id
		size_t lc_id = v->dec.ProctoID(i);
incardon's avatar
incardon committed
895

896
897
898
		// resize the receive buffer
		v->recv_mem_gg.get(lc_id).resize(msg_i);
		v->recv_sz.get(lc_id) = msg_i;
incardon's avatar
incardon committed
899

900
		return v->recv_mem_gg.get(lc_id).getPointer();
incardon's avatar
incardon committed
901
902
	}

903
	/*! \brief Call-back to allocate buffer to receive incoming elements (particles)
incardon's avatar
incardon committed
904
	 *
905
906
907
	 * \param msg_i size required to receive the message from i
	 * \param total_msg total size to receive from all the processors
	 * \param total_p the total number of processor that want to communicate with you
incardon's avatar
incardon committed
908
	 * \param i processor id
incardon's avatar
incardon committed
909
910
	 * \param ri request id (it is an id that goes from 0 to total_p, and is unique
	 *           every time message_alloc is called)
incardon's avatar
incardon committed
911
912
	 * \param ptr a pointer to the vector_dist structure
	 *
913
	 * \return the pointer where to store the message for the processor i
incardon's avatar
incardon committed
914
915
	 *
	 */
916
	static void * message_alloc_map(size_t msg_i ,size_t total_msg, size_t total_p, size_t i, size_t ri, void * ptr)
incardon's avatar
incardon committed
917
918
	{
		// cast the pointer
919
		vector_dist<dim,St,prop,Decomposition,Memory> * vd = static_cast<vector_dist<dim,St,prop,Decomposition,Memory> *>(ptr);
incardon's avatar
incardon committed
920

921
922
		vd->recv_mem_gm.resize(vd->v_cl.getProcessingUnits());
		vd->recv_mem_gm.get(i).resize(msg_i);
incardon's avatar
incardon committed
923

924
		return vd->recv_mem_gm.get(i).getPointer();
incardon's avatar
incardon committed
925
926
	}

927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
	/*! \brief Add local particle
	 *
	 * It add a local particle, with "local" we mean in this processor
	 * the particle can be also created out of the processor domain, in this
	 * case a call to map is required. Added particles are always created at the
	 * end and can be accessed with getLastPos and getLastProp
	 *
	 */
	void add()
	{
		v_prp.insert(g_m);
		v_pos.insert(g_m);

		g_m++;
	}

	/*! \brief Get the position of the last element
	 *
	 * \return the position of the element in space
	 *
	 */
	template<unsigned int id> inline auto getLastPos() -> decltype(v_pos.template get<id>(0))
	{
		return v_pos.template get<id>(g_m-1);
	}

	/*! \brief Get the property of the last element
	 *
	 * see the vector_dist iterator usage to get an element key
	 *
	 * \tparam id property id
	 * \param vec_key vector element
	 *
	 * \return return the selected property of the vector element
	 *
	 */
	template<unsigned int id> inline auto getLastProp() -> decltype(v_prp.template get<id>(0))
	{
		return v_prp.template get<id>(g_m-1);
	}

Pietro Incardona's avatar
Pietro Incardona committed
968
969
970
971
	/*! \brief Construct a cell list starting from the stored particles
	 *
	 * \tparam CellL CellList type to construct
	 *
972
973
974
975
	 * \param r_cut interation radius, or size of each cell
	 *
	 * \return the Cell list
	 *
Pietro Incardona's avatar
Pietro Incardona committed
976
	 */
977
	template<typename CellL=CellList<dim,St,FAST,shift<dim,St> > > CellL getCellList(St r_cut)
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
	{
		return getCellList(r_cut,dec.getGhost());
	}

	/*! \brief Construct a cell list starting from the stored particles
	 *
	 * It differ from the get getCellList for an additional parameter, in case the
	 * domain + ghost is not big enough to contain additional padding particles, a Cell list
	 * with bigger space can be created
	 * (padding particles in general are particles added by the user out of the domains)
	 *
	 * \tparam CellL CellList type to construct
	 *
	 * \param r_cut interation radius, or size of each cell
	 * \param enlarge In case of padding particles the cell list must be enlarged, like a ghost this parameter say how much must be enlarged
	 *
	 */
	template<typename CellL=CellList<dim,St,FAST,shift<dim,St> > > CellL getCellList(St r_cut, const Ghost<dim,St> & enlarge)
Pietro Incardona's avatar
Pietro Incardona committed
996
	{
997

Pietro Incardona's avatar
Pietro Incardona committed
998
999
		CellL cell_list;

1000
1001
1002
1003
		// calculate the parameters of the cell list

		// get the processor bounding box
		Box<dim,St> pbox = dec.getProcessorBounds();
1004
		// extend by the ghost
1005
		pbox.enlarge(enlarge);
1006
1007

		Box<dim,St> cell_box;
1008
1009
1010

		size_t div[dim];

1011
		// Calculate the division array and the cell box
1012
		for (size_t i = 0 ; i < dim ; i++)
1013
		{
1014
			div[i] = (pbox.getP2().get(i) - pbox.getP1().get(i))/ r_cut;
1015
1016
1017
1018
1019
			div[i]++;

			cell_box.setLow(i,0.0);
			cell_box.setHigh(i,div[i]*r_cut);
		}
1020

1021
		cell_list.Initialize(cell_box,div,pbox.getP1());
1022

Pietro Incardona's avatar
Pietro Incardona committed
1023
1024
		// for each particle add the particle to the cell list

1025
		auto it = getIterator();
Pietro Incardona's avatar
Pietro Incardona committed
1026
1027
1028
1029
1030
1031
1032
1033
1034
1035
1036
1037

		while (it.isNext())
		{
			auto key = it.get();

			cell_list.add(this->template getPos<0>(key),key.getKey());

			++it;
		}

		return cell_list;
	}
incardon's avatar
incardon committed
1038

1039
1040
1041
1042
1043
1044
1045
1046
1047
1048
1049
1050
1051
1052
1053
1054
1055
1056
1057
1058
1059
1060
1061
1062
1063
1064
1065
1066
1067
1068
1069
1070
1071
1072
1073
	/*! \brief It return the number of particles contained by the previous processors
	 *
	 * \Warning It only work with the initial decomposition
	 *
	 * Given 1000 particles and 3 processors, you will get
	 *
	 * * Processor 0: 0
	 * * Processor 1: 334
	 * * Processor 2: 667
	 *
	 * \param np initial number of particles
	 *
	 */
	size_t init_size_accum(size_t np)
	{
		size_t accum = 0;

		// convert to a local number of elements
		size_t p_np = np / v_cl.getProcessingUnits();

		// Get non divisible part
		size_t r = np % v_cl.getProcessingUnits();

		accum = p_np * v_cl.getProcessUnitID();

		// Distribute the remain particles
		if (v_cl.getProcessUnitID() <= r)
			accum += v_cl.getProcessUnitID();
		else
			accum += r;

		return accum;
	}

	/*! \brief Get an iterator that traverse domain and ghost particles
incardon's avatar
incardon committed
1074
1075
1076
1077
	 *
	 * \return an iterator
	 *
	 */
1078
	vector_dist_iterator getIterator()
incardon's avatar
incardon committed
1079
	{
1080
		return vector_dist_iterator(0,v_pos.size());
incardon's avatar
incardon committed
1081
1082
	}

1083
	/*! \brief Get an iterator that traverse the ghost particles
1084
1085
1086
1087
	 *
	 * \return an iterator
	 *
	 */
1088
	vector_dist_iterator getGhostIterator()
1089
	{
1090
		return vector_dist_iterator(g_m,v_pos.size());
1091
1092
	}

incardon's avatar
incardon committed
1093

1094
	/*! \brief Get an iterator that traverse the particles in the domain
1095
1096
1097
1098
	 *
	 * \return an iterator
	 *
	 */
1099
	vector_dist_iterator getDomainIterator()
1100
	{
1101
		return vector_dist_iterator(0,g_m);
1102
1103
	}

incardon's avatar
incardon committed
1104
1105
1106
1107
1108
	/*! \brief Get the decomposition
	 *
	 * \return
	 *
	 */
1109
	const Decomposition & getDecomposition()
incardon's avatar
incardon committed
1110
1111
1112
	{
		return dec;
	}
incardon's avatar
incardon committed
1113
1114
1115

	/*! \brief Output particle position and properties
	 *
1116
	 * \param out output
incardon's avatar
incardon committed
1117
1118
1119
1120
1121
1122
1123
	 * \param opt NO_GHOST or WITH_GHOST
	 *
	 * \return if the file has been written correctly
	 *
	 */
	inline bool write(std::string out, int opt = NO_GHOST)
	{
1124
1125
		// CSVWriter test
		CSVWriter<openfpm::vector<Point<dim,St>>, openfpm::vector<prop> > csv_writer;
incardon's avatar
incardon committed
1126

1127
		std::string output = std::to_string(out + std::to_string(v_cl.getProcessUnitID()) + std::to_string(".csv"));
incardon's avatar
incardon committed
1128

1129
		// Write the CSV
1130
		return csv_writer.write(output,v_pos,v_prp);
incardon's avatar
incardon committed
1131
	}
1132
1133
1134
1135
1136
1137
1138
1139
1140
1141

	/* \brief It return the id of structure in the allocation list
	 *
	 * \see print_alloc and SE_CLASS2
	 *
	 */
	long int who()
	{
#ifdef SE_CLASS2
		return check_whoami(this,8);
Pietro Incardona's avatar
Pietro Incardona committed
1142
1143
#else
		return -1;
1144
1145
#endif
	}
1146
1147
1148
1149
1150
1151
1152
1153
1154
1155
1156
1157
1158
1159

	/*! \brief Get the Virtual Cluster machine
	 *
	 * \return the Virtual cluster machine
	 *
	 */

	Vcluster & getVC()
	{
#ifdef SE_CLASS2
		check_valid(this,8);
#endif
		return v_cl;
	}
incardon's avatar
incardon committed
1160
1161
1162
1163
};


#endif /* VECTOR_HPP_ */