Block-Structured AMR Software Framework
AMReX_ParticleTile.H
Go to the documentation of this file.
1 #ifndef AMREX_PARTICLETILE_H_
2 #define AMREX_PARTICLETILE_H_
3 #include <AMReX_Config.H>
4 
5 #include <AMReX_Extension.H>
6 #include <AMReX_Particle.H>
7 #include <AMReX_ArrayOfStructs.H>
8 #include <AMReX_StructOfArrays.H>
9 #include <AMReX_Vector.H>
10 #include <AMReX_REAL.H>
11 #include <AMReX_RealVect.H>
12 
13 #include <array>
14 #include <string>
15 #include <type_traits>
16 #include <vector>
17 
18 
19 namespace amrex {
20 
21 // Forward Declaration
22 template <int NArrayReal, int NArrayInt>
23 struct ConstSoAParticle;
24 template <int NArrayReal, int NArrayInt>
25 struct SoAParticle;
26 
27 template <typename T_ParticleType, int NArrayReal, int NArrayInt>
28 struct ConstParticleTileData;
29 
30 template <typename T_ParticleType, int NArrayReal, int NArrayInt>
32 {
33  static constexpr int NAR = NArrayReal;
34  static constexpr int NAI = NArrayInt;
35 
36  using ParticleType = T_ParticleType;
37  using ParticleRefType = T_ParticleType&;
39 
40  static constexpr int NStructReal = ParticleType::NReal;
41  static constexpr int NStructInt = ParticleType::NInt;
42 
44 
45  static constexpr bool is_particle_tile_data = true;
46 
47  Long m_size;
48 
49  using AOS_PTR = std::conditional_t<T_ParticleType::is_soa_particle,
52 
53  uint64_t* m_idcpu;
56 
61 
63  ParticleReal& pos (const int dir, const int index) const &
64  {
65  if constexpr(!ParticleType::is_soa_particle) {
66  return this->m_aos[index].pos(dir);
67  } else {
68  return this->m_rdata[dir][index];
69  }
70  }
71 
73  decltype(auto) id (const int index) const &
74  {
75  if constexpr(!ParticleType::is_soa_particle) {
76  return this->m_aos[index].id();
77  } else {
78  return ParticleIDWrapper(this->m_idcpu[index]);
79  }
80  }
81 
83  decltype(auto) cpu (const int index) const &
84  {
85  if constexpr(!ParticleType::is_soa_particle) {
86  return this->m_aos[index].cpu();
87  } else {
88  return ParticleCPUWrapper(this->m_idcpu[index]);
89  }
90  }
91 
93  decltype(auto) idcpu (const int index) const &
94  {
95  if constexpr(ParticleType::is_soa_particle) {
96  return this->m_idcpu[index];
97  } else {
98  amrex::Abort("not implemented");
99  }
100  }
101 
103  ParticleReal * rdata (const int attribute_index) const
104  {
105  return this->m_rdata[attribute_index];
106  }
107 
109  int * idata (const int attribute_index) const
110  {
111  return this->m_idata[attribute_index];
112  }
113 
115  decltype(auto) operator[] (const int index) const
116  {
117  if constexpr (!ParticleType::is_soa_particle) {
118  return m_aos[index];
119  } else {
120  return SoAParticle<NAR, NAI>(*this, index);
121  }
122  }
123 
125  void packParticleData (char* buffer, int src_index, std::size_t dst_offset,
126  const int* comm_real, const int * comm_int) const noexcept
127  {
128  AMREX_ASSERT(src_index < m_size);
129  auto* dst = buffer + dst_offset;
130  if constexpr (!ParticleType::is_soa_particle) {
131  memcpy(dst, m_aos + src_index, sizeof(ParticleType));
132  dst += sizeof(ParticleType);
133  } else {
134  memcpy(dst, m_idcpu + src_index, sizeof(uint64_t));
135  dst += sizeof(uint64_t);
136  }
137  int array_start_index = AMREX_SPACEDIM + NStructReal;
138  for (int i = 0; i < NAR; ++i)
139  {
140  if (comm_real[array_start_index + i])
141  {
142  memcpy(dst, m_rdata[i] + src_index, sizeof(ParticleReal));
143  dst += sizeof(ParticleReal);
144  }
145  }
146  int runtime_start_index = AMREX_SPACEDIM + NStructReal + NAR;
147  for (int i = 0; i < m_num_runtime_real; ++i)
148  {
149  if (comm_real[runtime_start_index + i])
150  {
151  memcpy(dst, m_runtime_rdata[i] + src_index, sizeof(ParticleReal));
152  dst += sizeof(ParticleReal);
153  }
154  }
155  array_start_index = 2 + NStructInt;
156  for (int i = 0; i < NAI; ++i)
157  {
158  if (comm_int[array_start_index + i])
159  {
160  memcpy(dst, m_idata[i] + src_index, sizeof(int));
161  dst += sizeof(int);
162  }
163  }
164  runtime_start_index = 2 + NStructInt + NAI;
165  for (int i = 0; i < m_num_runtime_int; ++i)
166  {
167  if (comm_int[runtime_start_index + i])
168  {
169  memcpy(dst, m_runtime_idata[i] + src_index, sizeof(int));
170  dst += sizeof(int);
171  }
172  }
173  }
174 
176  void unpackParticleData (const char* buffer, Long src_offset, int dst_index,
177  const int* comm_real, const int* comm_int) const noexcept
178  {
179  AMREX_ASSERT(dst_index < m_size);
180  const auto* src = buffer + src_offset;
181  if constexpr (!ParticleType::is_soa_particle) {
182  memcpy(m_aos + dst_index, src, sizeof(ParticleType));
183  src += sizeof(ParticleType);
184  } else {
185  memcpy(m_idcpu + dst_index, src, sizeof(uint64_t));
186  src += sizeof(uint64_t);
187  }
188  int array_start_index = AMREX_SPACEDIM + NStructReal;
189  for (int i = 0; i < NAR; ++i)
190  {
191  if (comm_real[array_start_index + i])
192  {
193  memcpy(m_rdata[i] + dst_index, src, sizeof(ParticleReal));
194  src += sizeof(ParticleReal);
195  }
196  }
197  int runtime_start_index = AMREX_SPACEDIM + NStructReal + NAR;
198  for (int i = 0; i < m_num_runtime_real; ++i)
199  {
200  if (comm_real[runtime_start_index + i])
201  {
202  memcpy(m_runtime_rdata[i] + dst_index, src, sizeof(ParticleReal));
203  src += sizeof(ParticleReal);
204  }
205  }
206  array_start_index = 2 + NStructInt;
207  for (int i = 0; i < NAI; ++i)
208  {
209  if (comm_int[array_start_index + i])
210  {
211  memcpy(m_idata[i] + dst_index, src, sizeof(int));
212  src += sizeof(int);
213  }
214  }
215  runtime_start_index = 2 + NStructInt + NAI;
216  for (int i = 0; i < m_num_runtime_int; ++i)
217  {
218  if (comm_int[runtime_start_index + i])
219  {
220  memcpy(m_runtime_idata[i] + dst_index, src, sizeof(int));
221  src += sizeof(int);
222  }
223  }
224  }
225 
226  template <typename T = ParticleType, std::enable_if_t<!T::is_soa_particle, int> = 0>
228  SuperParticleType getSuperParticle (int index) const noexcept
229  {
230  AMREX_ASSERT(index < m_size);
232  for (int i = 0; i < AMREX_SPACEDIM; ++i) {
233  sp.pos(i) = m_aos[index].pos(i);
234  }
235  for (int i = 0; i < NStructReal; ++i) {
236  sp.rdata(i) = m_aos[index].rdata(i);
237  }
238  for (int i = 0; i < NAR; ++i) {
239  sp.rdata(NStructReal+i) = m_rdata[i][index];
240  }
241  sp.id() = m_aos[index].id();
242  sp.cpu() = m_aos[index].cpu();
243  for (int i = 0; i < NStructInt; ++i) {
244  sp.idata(i) = m_aos[index].idata(i);
245  }
246  for (int i = 0; i < NAI; ++i) {
247  sp.idata(NStructInt+i) = m_idata[i][index];
248  }
249  return sp;
250  }
251 
252  template <typename T = ParticleType, std::enable_if_t<T::is_soa_particle, int> = 0>
254  SuperParticleType getSuperParticle (int index) const noexcept
255  {
256  AMREX_ASSERT(index < m_size);
258  sp.m_idcpu = m_idcpu[index];
259  for (int i = 0; i < AMREX_SPACEDIM; ++i) {sp.pos(i) = m_rdata[i][index];}
260  for (int i = 0; i < NAR; ++i) {
261  sp.rdata(i) = m_rdata[i][index];
262  }
263  for (int i = 0; i < NAI; ++i) {
264  sp.idata(i) = m_idata[i][index];
265  }
266  return sp;
267  }
268 
269  template <typename T = ParticleType, std::enable_if_t<!T::is_soa_particle, int> = 0>
271  void setSuperParticle (const SuperParticleType& sp, int index) const noexcept
272  {
273  for (int i = 0; i < AMREX_SPACEDIM; ++i) {
274  m_aos[index].pos(i) = sp.pos(i);
275  }
276  for (int i = 0; i < NStructReal; ++i) {
277  m_aos[index].rdata(i) = sp.rdata(i);
278  }
279  for (int i = 0; i < NAR; ++i) {
280  m_rdata[i][index] = sp.rdata(NStructReal+i);
281  }
282  m_aos[index].id() = sp.id();
283  m_aos[index].cpu() = sp.cpu();
284  for (int i = 0; i < NStructInt; ++i) {
285  m_aos[index].idata(i) = sp.idata(i);
286  }
287  for (int i = 0; i < NAI; ++i) {
288  m_idata[i][index] = sp.idata(NStructInt+i);
289  }
290  }
291 
292  template <typename T = ParticleType, std::enable_if_t<T::is_soa_particle, int> = 0>
294  void setSuperParticle (const SuperParticleType& sp, int index) const noexcept
295  {
296  m_idcpu[index] = sp.m_idcpu;
297  for (int i = 0; i < NAR; ++i) {
298  m_rdata[i][index] = sp.rdata(i);
299  }
300  for (int i = 0; i < NAI; ++i) {
301  m_idata[i][index] = sp.idata(i);
302  }
303  }
304 };
305 
306 // SOA Particle Structure
307 template <int T_NArrayReal, int T_NArrayInt>
309 {
310  static constexpr int NArrayReal = T_NArrayReal;
311  static constexpr int NArrayInt = T_NArrayInt;
314  static constexpr bool is_soa_particle = true;
315  static constexpr bool is_constsoa_particle = true;
316 
317  using RealType = ParticleReal;
318 
320  ConstSoAParticle (ConstPTD const& ptd, long i) : // Note: should this be int instead?
322  {
323  }
324 
325  //static Long the_next_id;
326 
327  //functions to get id and cpu in the SOA data
328 
331 
334 
335  //functions to get positions of the particle in the SOA data
336 
338  RealVect pos () const & {return RealVect(AMREX_D_DECL(this->m_constparticle_tile_data.m_rdata[0][m_index], this->m_constparticle_tile_data.m_rdata[1][m_index], this->m_constparticle_tile_data.m_rdata[2][m_index]));}
339 
341  const RealType& pos (int position_index) const &
342  {
343  AMREX_ASSERT(position_index < AMREX_SPACEDIM);
344  return this->m_constparticle_tile_data.m_rdata[position_index][m_index];
345  }
346 
347  static Long NextID ();
348 
352  static Long UnprotectedNextID ();
353 
359  static void NextID (Long nextid);
360 
361  private :
362 
363  static_assert(std::is_trivially_copyable<ParticleTileData<SoAParticleBase, NArrayReal, NArrayInt>>(), "ParticleTileData is not trivially copyable");
364 
366  int m_index;
367 };
368 
369 template <int T_NArrayReal, int T_NArrayInt>
371 {
372  static constexpr int NArrayReal = T_NArrayReal;
373  static constexpr int NArrayInt = T_NArrayInt;
376  static constexpr bool is_soa_particle = true;
377  static constexpr bool is_constsoa_particle = false;
378 
380  using RealType = ParticleReal;
381 
383  SoAParticle (PTD const& ptd, long i) : // Note: should this be int instead?
385  {
386  }
387 
388  static Long the_next_id;
389 
390  //functions to get id and cpu in the SOA data
391 
394 
397 
399  uint64_t& idcpu () & { return this->m_particle_tile_data.m_idcpu[m_index]; }
400 
403 
406 
408  const uint64_t& idcpu () const & { return this->m_particle_tile_data.m_idcpu[m_index]; }
409 
410  //functions to get positions of the particle in the SOA data
411 
413  RealVect pos () const & {return RealVect(AMREX_D_DECL(this->m_particle_tile_data.m_rdata[0][m_index], this->m_particle_tile_data.m_rdata[1][m_index], this->m_particle_tile_data.m_rdata[2][m_index]));}
414 
416  RealType& pos (int position_index) &
417  {
418  AMREX_ASSERT(position_index < AMREX_SPACEDIM);
419  return this->m_particle_tile_data.m_rdata[position_index][m_index];
420  }
421 
423  RealType pos (int position_index) const &
424  {
425  AMREX_ASSERT(position_index < AMREX_SPACEDIM);
426  return this->m_particle_tile_data.m_rdata[position_index][m_index];
427  }
428 
429  static Long NextID ();
430 
434  static Long UnprotectedNextID ();
435 
441  static void NextID (Long nextid);
442 
443 private :
444 
445  static_assert(std::is_trivially_copyable<ParticleTileData<SoAParticleBase, NArrayReal, NArrayInt>>(), "ParticleTileData is not trivially copyable");
446 
448  int m_index;
449 };
450 
451 //template <int NArrayReal, int NArrayInt> Long ConstSoAParticle<NArrayReal, NArrayInt>::the_next_id = 1;
452 template <int NArrayReal, int NArrayInt> Long SoAParticle<NArrayReal, NArrayInt>::the_next_id = 1;
453 
454 template <int NArrayReal, int NArrayInt>
455 Long
457 {
458  Long next;
459 // we should be able to test on _OPENMP < 201107 for capture (version 3.1)
460 // but we must work around a bug in gcc < 4.9
461 #if defined(AMREX_USE_OMP) && defined(_OPENMP) && _OPENMP < 201307
462 #pragma omp critical (amrex_particle_nextid)
463 #elif defined(AMREX_USE_OMP)
464 #pragma omp atomic capture
465 #endif
466  next = the_next_id++;
467 
468  if (next > LongParticleIds::LastParticleID) {
469  amrex::Abort("SoAParticle<NArrayReal, NArrayInt>::NextID() -- too many particles");
470  }
471 
472  return next;
473 }
474 
475 template <int NArrayReal, int NArrayInt>
476 Long
478 {
479  Long next = the_next_id++;
480  if (next > LongParticleIds::LastParticleID) {
481  amrex::Abort("SoAParticle<NArrayReal, NArrayInt>::NextID() -- too many particles");
482  }
483  return next;
484 }
485 
486 template <int NArrayReal, int NArrayInt>
487 void
489 {
490  the_next_id = nextid;
491 }
492 
493 template <typename T_ParticleType, int NArrayReal, int NArrayInt>
495 {
496  static constexpr int NAR = NArrayReal;
497  static constexpr int NAI = NArrayInt;
498  using ParticleType = T_ParticleType;
499  using ParticleRefType = T_ParticleType const&;
500 
501  static constexpr int NStructReal = ParticleType::NReal;
502  static constexpr int NStructInt = ParticleType::NInt;
503 
505 
506  static constexpr bool is_particle_tile_data = true;
507 
508  Long m_size;
509 
510  using AOS_PTR = std::conditional_t<T_ParticleType::is_soa_particle,
511  void const * AMREX_RESTRICT, ParticleType const * AMREX_RESTRICT>;
513 
514  const uint64_t* m_idcpu;
517 
519  const ParticleReal& pos (const int dir, const int index) const &
520  {
521  if constexpr(!ParticleType::is_soa_particle) {
522  return this->m_aos[index].pos(dir);
523  } else {
524  return this->m_rdata[dir][index];
525  }
526  }
527 
529  decltype(auto) id (const int index) const &
530  {
531  if constexpr(!ParticleType::is_soa_particle) {
532  return this->m_aos[index].id();
533  } else {
534  return ConstParticleIDWrapper(this->m_idcpu[index]);
535  }
536  }
537 
539  decltype(auto) cpu (const int index) const &
540  {
541  if constexpr(!ParticleType::is_soa_particle) {
542  return this->m_aos[index].cpu();
543  } else {
544  return ConstParticleCPUWrapper(this->m_idcpu[index]);
545  }
546  }
547 
549  decltype(auto) idcpu (const int index) const &
550  {
551  if constexpr(ParticleType::is_soa_particle) {
552  return this->m_idcpu[index];
553  } else {
554  amrex::Abort("not implemented");
555  }
556  }
557 
559  const ParticleReal * rdata (const int attribute_index) const
560  {
561  return this->m_rdata[attribute_index];
562  }
563 
565  const int * idata (const int attribute_index) const
566  {
567  return this->m_idata[attribute_index];
568  }
569 
571  decltype(auto) operator[] (const int index) const
572  {
573  if constexpr (!ParticleType::is_soa_particle) {
574  return m_aos[index];
575  } else {
576  return ConstSoAParticle<NAR, NAI>(*this, index);
577  }
578  }
579 
584 
586  void packParticleData(char* buffer, int src_index, Long dst_offset,
587  const int* comm_real, const int * comm_int) const noexcept
588  {
589  AMREX_ASSERT(src_index < m_size);
590  auto* dst = buffer + dst_offset;
591  if constexpr (!ParticleType::is_soa_particle) {
592  memcpy(dst, m_aos + src_index, sizeof(ParticleType));
593  dst += sizeof(ParticleType);
594  } else {
595  memcpy(dst, m_idcpu + src_index, sizeof(uint64_t));
596  dst += sizeof(uint64_t);
597  }
598  int array_start_index = AMREX_SPACEDIM + NStructReal;
599  for (int i = 0; i < NArrayReal; ++i)
600  {
601  if (comm_real[array_start_index + i])
602  {
603  memcpy(dst, m_rdata[i] + src_index, sizeof(ParticleReal));
604  dst += sizeof(ParticleReal);
605  }
606  }
607  int runtime_start_index = AMREX_SPACEDIM + NStructReal + NArrayReal;
608  for (int i = 0; i < m_num_runtime_real; ++i)
609  {
610  if (comm_real[runtime_start_index + i])
611  {
612  memcpy(dst, m_runtime_rdata[i] + src_index, sizeof(ParticleReal));
613  dst += sizeof(ParticleReal);
614  }
615  }
616  array_start_index = 2 + NStructInt;
617  for (int i = 0; i < NArrayInt; ++i)
618  {
619  if (comm_int[array_start_index + i])
620  {
621  memcpy(dst, m_idata[i] + src_index, sizeof(int));
622  dst += sizeof(int);
623  }
624  }
625  runtime_start_index = 2 + NStructInt + NArrayInt;
626  for (int i = 0; i < m_num_runtime_int; ++i)
627  {
628  if (comm_int[runtime_start_index + i])
629  {
630  memcpy(dst, m_runtime_idata[i] + src_index, sizeof(int));
631  dst += sizeof(int);
632  }
633  }
634  }
635 
636  template <typename T = ParticleType, std::enable_if_t<!T::is_soa_particle, int> = 0>
638  SuperParticleType getSuperParticle (int index) const noexcept
639  {
640  AMREX_ASSERT(index < m_size);
642  for (int i = 0; i < AMREX_SPACEDIM; ++i) {
643  sp.pos(i) = m_aos[index].pos(i);
644  }
645  for (int i = 0; i < NStructReal; ++i) {
646  sp.rdata(i) = m_aos[index].rdata(i);
647  }
648  if constexpr(NArrayReal > 0) {
649  for (int i = 0; i < NArrayReal; ++i) {
650  sp.rdata(NStructReal+i) = m_rdata[i][index];
651  }
652  }
653  sp.id() = m_aos[index].id();
654  sp.cpu() = m_aos[index].cpu();
655  for (int i = 0; i < NStructInt; ++i) {
656  sp.idata(i) = m_aos[index].idata(i);
657  }
658  if constexpr(NArrayInt > 0) {
659  for (int i = 0; i < NArrayInt; ++i) {
660  sp.idata(NStructInt+i) = m_idata[i][index];
661  }
662  }
663  return sp;
664  }
665 
666  template <typename T = ParticleType, std::enable_if_t<T::is_soa_particle, int> = 0>
668  SuperParticleType getSuperParticle (int index) const noexcept
669  {
670  AMREX_ASSERT(index < m_size);
672  for (int i = 0; i < AMREX_SPACEDIM; ++i) {sp.pos(i) = m_rdata[i][index];}
673  sp.m_idcpu = m_idcpu[index];
674  for (int i = 0; i < NAR; ++i) {
675  sp.rdata(i) = m_rdata[i][index];
676  }
677  for (int i = 0; i < NAI; ++i) {
678  sp.idata(i) = m_idata[i][index];
679  }
680  return sp;
681  }
682 };
683 
685 
688 };
689 
690 template <typename T_ParticleType, int NArrayReal, int NArrayInt,
691  template<class> class Allocator=DefaultAllocator>
693 {
694  template <typename T>
695  using AllocatorType = Allocator<T>;
696 
697  using ParticleType = T_ParticleType;
698  static constexpr int NAR = NArrayReal;
699  static constexpr int NAI = NArrayInt;
700  using RealType = typename ParticleType::RealType;
701 
702  static constexpr int NStructReal = ParticleType::NReal;
703  static constexpr int NStructInt = ParticleType::NInt;
704 
706 
707  using AoS = std::conditional_t<
708  ParticleType::is_soa_particle,
711  //using ParticleVector = typename AoS::ParticleVector;
712 
713  using SoA = std::conditional_t<
714  ParticleType::is_soa_particle,
717  using RealVector = typename SoA::RealVector;
718  using IntVector = typename SoA::IntVector;
719  using StorageParticleType = typename ParticleType::StorageParticleType;
720 
723 
724  ParticleTile () = default;
725 
726 #ifndef _WIN32 // workaround windows compiler bug
727  ~ParticleTile () = default;
728 
729  ParticleTile (ParticleTile const&) = delete;
730  ParticleTile (ParticleTile &&) noexcept = default;
731 
732  ParticleTile& operator= (ParticleTile const&) = delete;
733  ParticleTile& operator= (ParticleTile &&) noexcept = default;
734 #endif
735 
736  void define (
737  int a_num_runtime_real,
738  int a_num_runtime_int,
739  std::vector<std::string>* soa_rdata_names=nullptr,
740  std::vector<std::string>* soa_idata_names=nullptr
741  )
742  {
743  m_defined = true;
744  GetStructOfArrays().define(a_num_runtime_real, a_num_runtime_int, soa_rdata_names, soa_idata_names);
745  m_runtime_r_ptrs.resize(a_num_runtime_real);
746  m_runtime_i_ptrs.resize(a_num_runtime_int);
747  m_runtime_r_cptrs.resize(a_num_runtime_real);
748  m_runtime_i_cptrs.resize(a_num_runtime_int);
749  }
750 
751  // Get id data
752  decltype(auto) id (int index) & {
753  if constexpr (!ParticleType::is_soa_particle) {
754  return m_aos_tile[index].id();
755  } else {
756  return ParticleIDWrapper(m_soa_tile.GetIdCPUData()[index]);
757  }
758  }
759 
760  // const
761  decltype(auto) id (int index) const & {
762  if constexpr (!ParticleType::is_soa_particle) {
763  return m_aos_tile[index].id();
764  } else {
765  return ConstParticleIDWrapper(m_soa_tile.GetIdCPUData()[index]);
766  }
767  }
768 
769  // Get cpu data
770  decltype(auto) cpu (int index) & {
771  if constexpr (!ParticleType::is_soa_particle) {
772  return m_aos_tile[index].cpu();
773  } else {
774  return ParticleCPUWrapper(m_soa_tile.GetIdCPUData()[index]);
775  }
776  }
777 
778  // const
779  decltype(auto) cpu (int index) const & {
780  if constexpr (!ParticleType::is_soa_particle) {
781  return m_aos_tile[index].cpu();
782  } else {
783  return ConstParticleCPUWrapper(m_soa_tile.GetIdCPUData()[index]);
784  }
785  }
786 
787  // Get positions data
788  RealType& pos (int index, int position_index) & {
789  if constexpr (!ParticleType::is_soa_particle) {
790  return m_aos_tile[index].pos(position_index);
791  } else {
792  static_assert(NArrayReal == ParticleType::PTD::NAR, "ParticleTile mismatch in R");
793  static_assert(NArrayInt == ParticleType::PTD::NAI, "ParticleTile mismatch in I");
794  static_assert(0 == ParticleType::StorageParticleType::NReal, "ParticleTile 2 mismatch in R");
795  static_assert(0 == ParticleType::StorageParticleType::NInt, "ParticleTile 2 mismatch in I");
796 
797  return m_soa_tile.GetRealData(position_index)[index];
798  }
799  }
800 
801  // const
802  RealType pos (int index, int position_index) const &
803  {
804  if constexpr (!ParticleType::is_soa_particle) {
805  return m_aos_tile[index].pos(position_index);
806  } else {
807  return m_soa_tile.GetRealData(position_index)[index];
808  }
809  }
810 
812  const AoS& GetArrayOfStructs () const { return m_aos_tile; }
813 
815  const SoA& GetStructOfArrays () const { return m_soa_tile; }
816 
817  bool empty () const { return size() == 0; }
818 
823  std::size_t size () const
824  {
825  if constexpr (!ParticleType::is_soa_particle) {
826  return m_aos_tile.size();
827  } else {
828  return m_soa_tile.size();
829  }
830  }
831 
836  int numParticles () const
837  {
838  if constexpr (!ParticleType::is_soa_particle) {
839  return m_aos_tile.numParticles();
840  } else {
841  return m_soa_tile.numParticles();
842  }
843  }
844 
849  int numRealParticles () const
850  {
851  if constexpr (!ParticleType::is_soa_particle) {
852  return m_aos_tile.numRealParticles();
853  } else {
854  return m_soa_tile.numRealParticles();
855  }
856  }
857 
862  int numNeighborParticles () const
863  {
864  if constexpr (!ParticleType::is_soa_particle) {
865  return m_aos_tile.numNeighborParticles();
866  } else {
867  return m_soa_tile.numNeighborParticles();
868  }
869  }
870 
875  int numTotalParticles () const
876  {
877  if constexpr (!ParticleType::is_soa_particle) {
878  return m_aos_tile.numTotalParticles();
879  } else {
880  return m_soa_tile.numTotalParticles();
881  }
882  }
883 
884  void setNumNeighbors (int num_neighbors)
885  {
886  if constexpr(!ParticleType::is_soa_particle) {
887  m_aos_tile.setNumNeighbors(num_neighbors);
888  }
889  m_soa_tile.setNumNeighbors(num_neighbors);
890  }
891 
892  int getNumNeighbors () const
893  {
894  if constexpr (!ParticleType::is_soa_particle) {
895  AMREX_ASSERT( m_soa_tile.getNumNeighbors() == m_aos_tile.getNumNeighbors() );
896  return m_aos_tile.getNumNeighbors();
897  } else {
898  return m_soa_tile.getNumNeighbors();
899  }
900  }
901 
902  void resize (std::size_t count)
903  {
904  if constexpr (!ParticleType::is_soa_particle) {
905  m_aos_tile.resize(count);
906  }
907  m_soa_tile.resize(count);
908  }
909 
913  template <typename T = ParticleType, std::enable_if_t<!T::is_soa_particle, int> = 0>
914  void push_back (const ParticleType& p) { m_aos_tile().push_back(p); }
915 
919  template < int NR = NArrayReal, int NI = NArrayInt,
920  std::enable_if_t<NR != 0 || NI != 0, int> foo = 0>
921  void push_back (const SuperParticleType& sp)
922  {
923  auto np = numParticles();
924 
925  if constexpr (!ParticleType::is_soa_particle) {
926  m_aos_tile.resize(np+1);
927  for (int i = 0; i < AMREX_SPACEDIM; ++i) {
928  m_aos_tile[np].pos(i) = sp.pos(i);
929  }
930  for (int i = 0; i < NStructReal; ++i) {
931  m_aos_tile[np].rdata(i) = sp.rdata(i);
932  }
933  m_aos_tile[np].id() = sp.id();
934  m_aos_tile[np].cpu() = sp.cpu();
935  for (int i = 0; i < NStructInt; ++i) {
936  m_aos_tile[np].idata(i) = sp.idata(i);
937  }
938  }
939 
940  m_soa_tile.resize(np+1);
941  if constexpr (ParticleType::is_soa_particle) {
942  m_soa_tile.GetIdCPUData()[np] = sp.m_idcpu;
943  }
944  auto& arr_rdata = m_soa_tile.GetRealData();
945  auto& arr_idata = m_soa_tile.GetIntData();
946  for (int i = 0; i < NArrayReal; ++i) {
947  arr_rdata[i][np] = sp.rdata(NStructReal+i);
948  }
949  for (int i = 0; i < NArrayInt; ++i) {
950  arr_idata[i][np] = sp.idata(NStructInt+i);
951  }
952  }
953 
958  void push_back_real (int comp, ParticleReal v) {
959  m_soa_tile.GetRealData(comp).push_back(v);
960  }
961 
966  void push_back_real (const std::array<ParticleReal, NArrayReal>& v) {
967  for (int i = 0; i < NArrayReal; ++i) {
968  m_soa_tile.GetRealData(i).push_back(v[i]);
969  }
970  }
971 
976  void push_back_real (int comp, const ParticleReal* beg, const ParticleReal* end) {
977  auto it = m_soa_tile.GetRealData(comp).end();
978  m_soa_tile.GetRealData(comp).insert(it, beg, end);
979  }
980 
986  push_back_real(comp, &(*beg), &(*end));
987  }
988 
993  void push_back_real (int comp, amrex::Vector<amrex::ParticleReal> const & vec) {
994  push_back_real(comp, vec.cbegin(), vec.cend());
995  }
996 
1001  void push_back_real (int comp, std::size_t npar, ParticleReal v) {
1002  auto new_size = m_soa_tile.GetRealData(comp).size() + npar;
1003  m_soa_tile.GetRealData(comp).resize(new_size, v);
1004  }
1005 
1010  void push_back_int (int comp, int v) {
1011  m_soa_tile.GetIntData(comp).push_back(v);
1012  }
1013 
1018  void push_back_int (const std::array<int, NArrayInt>& v) {
1019  for (int i = 0; i < NArrayInt; ++i) {
1020  m_soa_tile.GetIntData(i).push_back(v[i]);
1021  }
1022  }
1023 
1028  void push_back_int (int comp, const int* beg, const int* end) {
1029  auto it = m_soa_tile.GetIntData(comp).end();
1030  m_soa_tile.GetIntData(comp).insert(it, beg, end);
1031  }
1032 
1038  push_back_int(comp, &(*beg), &(*end));
1039  }
1040 
1045  void push_back_int (int comp, amrex::Vector<int> const & vec) {
1046  push_back_int(comp, vec.cbegin(), vec.cend());
1047  }
1048 
1053  void push_back_int (int comp, std::size_t npar, int v) {
1054  auto new_size = m_soa_tile.GetIntData(comp).size() + npar;
1055  m_soa_tile.GetIntData(comp).resize(new_size, v);
1056  }
1057 
1058  int NumRealComps () const noexcept { return m_soa_tile.NumRealComps(); }
1059 
1060  int NumIntComps () const noexcept { return m_soa_tile.NumIntComps(); }
1061 
1062  int NumRuntimeRealComps () const noexcept { return m_runtime_r_ptrs.size(); }
1063 
1064  int NumRuntimeIntComps () const noexcept { return m_runtime_i_ptrs.size(); }
1065 
1067  {
1068  if constexpr (ParticleType::is_soa_particle) {
1069  GetStructOfArrays().GetIdCPUData().shrink_to_fit();
1070  } else {
1071  m_aos_tile().shrink_to_fit();
1072  }
1073  for (int j = 0; j < NumRealComps(); ++j)
1074  {
1075  auto& rdata = GetStructOfArrays().GetRealData(j);
1076  rdata.shrink_to_fit();
1077  }
1078 
1079  for (int j = 0; j < NumIntComps(); ++j)
1080  {
1081  auto& idata = GetStructOfArrays().GetIntData(j);
1082  idata.shrink_to_fit();
1083  }
1084  }
1085 
1086  Long capacity () const
1087  {
1088  Long nbytes = 0;
1089  if constexpr (ParticleType::is_soa_particle) {
1090  nbytes += GetStructOfArrays().GetIdCPUData().capacity() * sizeof(uint64_t);
1091  } else {
1092  nbytes += m_aos_tile().capacity() * sizeof(ParticleType);
1093  }
1094  for (int j = 0; j < NumRealComps(); ++j)
1095  {
1096  auto& rdata = GetStructOfArrays().GetRealData(j);
1097  nbytes += rdata.capacity() * sizeof(ParticleReal);
1098  }
1099 
1100  for (int j = 0; j < NumIntComps(); ++j)
1101  {
1102  auto& idata = GetStructOfArrays().GetIntData(j);
1103  nbytes += idata.capacity()*sizeof(int);
1104  }
1105  return nbytes;
1106  }
1107 
1109  {
1110  if constexpr (ParticleType::is_soa_particle) {
1111  GetStructOfArrays().GetIdCPUData().swap(other.GetStructOfArrays().GetIdCPUData());
1112  } else {
1113  m_aos_tile().swap(other.GetArrayOfStructs()());
1114  }
1115  for (int j = 0; j < NumRealComps(); ++j)
1116  {
1117  auto& rdata = GetStructOfArrays().GetRealData(j);
1118  rdata.swap(other.GetStructOfArrays().GetRealData(j));
1119  }
1120 
1121  for (int j = 0; j < NumIntComps(); ++j)
1122  {
1123  auto& idata = GetStructOfArrays().GetIntData(j);
1124  idata.swap(other.GetStructOfArrays().GetIntData(j));
1125  }
1126  }
1127 
1129  {
1130  m_runtime_r_ptrs.resize(m_soa_tile.NumRealComps() - NArrayReal);
1131  m_runtime_i_ptrs.resize(m_soa_tile.NumIntComps() - NArrayInt);
1132 #ifdef AMREX_USE_GPU
1133  bool copy_real = false;
1134  m_h_runtime_r_ptrs.resize(m_soa_tile.NumRealComps() - NArrayReal, nullptr);
1135  for (std::size_t i = 0; i < m_h_runtime_r_ptrs.size(); ++i) {
1136  if (m_h_runtime_r_ptrs[i] != m_soa_tile.GetRealData(i + NArrayReal).dataPtr()) {
1137  m_h_runtime_r_ptrs[i] = m_soa_tile.GetRealData(i + NArrayReal).dataPtr();
1138  copy_real = true;
1139  }
1140  }
1141  if (copy_real) {
1143  m_h_runtime_r_ptrs.size()*sizeof(ParticleReal*));
1144  }
1145 
1146  bool copy_int = false;
1147  m_h_runtime_i_ptrs.resize(m_soa_tile.NumIntComps() - NArrayInt, nullptr);
1148  for (std::size_t i = 0; i < m_h_runtime_i_ptrs.size(); ++i) {
1149  if (m_h_runtime_i_ptrs[i] != m_soa_tile.GetIntData(i + NArrayInt).dataPtr()) {
1150  m_h_runtime_i_ptrs[i] = m_soa_tile.GetIntData(i + NArrayInt).dataPtr();
1151  copy_int = true;
1152  }
1153  }
1154  if (copy_int) {
1156  m_h_runtime_i_ptrs.size()*sizeof(int*));
1157  }
1158 #else
1159  for (std::size_t i = 0; i < m_runtime_r_ptrs.size(); ++i) {
1160  m_runtime_r_ptrs[i] = m_soa_tile.GetRealData(i + NArrayReal).dataPtr();
1161  }
1162 
1163  for (std::size_t i = 0; i < m_runtime_i_ptrs.size(); ++i) {
1164  m_runtime_i_ptrs[i] = m_soa_tile.GetIntData(i + NArrayInt).dataPtr();
1165  }
1166 #endif
1167 
1169  if constexpr (!ParticleType::is_soa_particle) {
1170  ptd.m_aos = m_aos_tile().dataPtr();
1171  } else {
1172  ptd.m_aos = nullptr;
1173  }
1174  if constexpr (ParticleType::is_soa_particle) {
1175  ptd.m_idcpu = m_soa_tile.GetIdCPUData().dataPtr();
1176  } else {
1177  ptd.m_idcpu = nullptr;
1178  }
1179  if constexpr(NArrayReal > 0) {
1180  for (int i = 0; i < NArrayReal; ++i) {
1181  ptd.m_rdata[i] = m_soa_tile.GetRealData(i).dataPtr();
1182  }
1183  }
1184  if constexpr(NArrayInt > 0) {
1185  for (int i = 0; i < NArrayInt; ++i) {
1186  ptd.m_idata[i] = m_soa_tile.GetIntData(i).dataPtr();
1187  }
1188  }
1189  ptd.m_size = size();
1190  ptd.m_num_runtime_real = m_runtime_r_ptrs.size();
1191  ptd.m_num_runtime_int = m_runtime_i_ptrs.size();
1192  ptd.m_runtime_rdata = m_runtime_r_ptrs.dataPtr();
1193  ptd.m_runtime_idata = m_runtime_i_ptrs.dataPtr();
1194 
1195 #ifdef AMREX_USE_GPU
1196  if (copy_real || copy_int) {
1198  }
1199 #endif
1200 
1201  return ptd;
1202  }
1203 
1205  {
1206  m_runtime_r_cptrs.resize(m_soa_tile.NumRealComps() - NArrayReal);
1207  m_runtime_i_cptrs.resize(m_soa_tile.NumIntComps() - NArrayInt);
1208 #ifdef AMREX_USE_GPU
1209  bool copy_real = false;
1210  m_h_runtime_r_cptrs.resize(m_soa_tile.NumRealComps() - NArrayReal, nullptr);
1211  for (std::size_t i = 0; i < m_h_runtime_r_cptrs.size(); ++i) {
1212  if (m_h_runtime_r_cptrs[i] != m_soa_tile.GetRealData(i + NArrayReal).dataPtr()) {
1213  m_h_runtime_r_cptrs[i] = m_soa_tile.GetRealData(i + NArrayReal).dataPtr();
1214  copy_real = true;
1215  }
1216  }
1217  if (copy_real) {
1219  m_h_runtime_r_cptrs.size()*sizeof(ParticleReal*));
1220  }
1221 
1222  bool copy_int = false;
1223  m_h_runtime_i_cptrs.resize(m_soa_tile.NumIntComps() - NArrayInt, nullptr);
1224  for (std::size_t i = 0; i < m_h_runtime_i_cptrs.size(); ++i) {
1225  if (m_h_runtime_i_cptrs[i] != m_soa_tile.GetIntData(i + NArrayInt).dataPtr()) {
1226  m_h_runtime_i_cptrs[i] = m_soa_tile.GetIntData(i + NArrayInt).dataPtr();
1227  copy_int = true;
1228  }
1229  }
1230  if (copy_int) {
1232  m_h_runtime_i_cptrs.size()*sizeof(int*));
1233  }
1234 #else
1235  for (std::size_t i = 0; i < m_runtime_r_cptrs.size(); ++i) {
1236  m_runtime_r_cptrs[i] = m_soa_tile.GetRealData(i + NArrayReal).dataPtr();
1237  }
1238 
1239  for (std::size_t i = 0; i < m_runtime_i_cptrs.size(); ++i) {
1240  m_runtime_i_cptrs[i] = m_soa_tile.GetIntData(i + NArrayInt).dataPtr();
1241  }
1242 #endif
1243 
1245  if constexpr (!ParticleType::is_soa_particle) {
1246  ptd.m_aos = m_aos_tile().dataPtr();
1247  } else {
1248  ptd.m_aos = nullptr;
1249  }
1250  if constexpr (ParticleType::is_soa_particle) {
1251  ptd.m_idcpu = m_soa_tile.GetIdCPUData().dataPtr();
1252  } else {
1253  ptd.m_idcpu = nullptr;
1254  }
1255  if constexpr(NArrayReal > 0) {
1256  for (int i = 0; i < NArrayReal; ++i) {
1257  ptd.m_rdata[i] = m_soa_tile.GetRealData(i).dataPtr();
1258  }
1259  }
1260  if constexpr(NArrayInt > 0) {
1261  for (int i = 0; i < NArrayInt; ++i) {
1262  ptd.m_idata[i] = m_soa_tile.GetIntData(i).dataPtr();
1263  }
1264  }
1265  ptd.m_size = size();
1267  ptd.m_num_runtime_int = m_runtime_i_cptrs.size();
1268  ptd.m_runtime_rdata = m_runtime_r_cptrs.dataPtr();
1269  ptd.m_runtime_idata = m_runtime_i_cptrs.dataPtr();
1270 
1271 #ifdef AMREX_USE_GPU
1272  if (copy_real || copy_int) {
1274  }
1275 #endif
1276 
1277  return ptd;
1278  }
1279 
1280 private:
1281 
1284 
1285  bool m_defined = false;
1286 
1289 
1292 
1295 
1298 };
1299 
1300 } // namespace amrex
1301 
1302 #endif // AMREX_PARTICLETILE_H_
#define AMREX_ASSERT(EX)
Definition: AMReX_BLassert.H:38
#define AMREX_FORCE_INLINE
Definition: AMReX_Extension.H:119
#define AMREX_RESTRICT
Definition: AMReX_Extension.H:37
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
#define AMREX_D_DECL(a, b, c)
Definition: AMReX_SPACE.H:104
Definition: AMReX_GpuAllocators.H:114
Definition: AMReX_ArrayOfStructs.H:13
Definition: AMReX_PODVector.H:246
A Real vector in SpaceDim-dimensional space.
Definition: AMReX_RealVect.H:32
This class is a thin wrapper around std::vector. Unlike vector, Vector::operator[] provides bound che...
Definition: AMReX_Vector.H:27
void streamSynchronize() noexcept
Definition: AMReX_GpuDevice.H:237
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void * memcpy(void *dest, const void *src, std::size_t count)
Definition: AMReX_GpuUtility.H:214
void htod_memcpy_async(void *p_d, const void *p_h, const std::size_t sz) noexcept
Definition: AMReX_GpuDevice.H:251
constexpr Long LastParticleID
Definition: AMReX_Particle.H:20
Definition: AMReX_Amr.cpp:49
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 end(BoxND< dim > const &box) noexcept
Definition: AMReX_Box.H:1890
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition: AMReX.cpp:225
const int[]
Definition: AMReX_BLProfiler.cpp:1664
Definition: AMReX_Particle.H:222
Definition: AMReX_Particle.H:188
static constexpr bool is_particle_tile_data
Definition: AMReX_ParticleTile.H:506
static constexpr int NAI
Definition: AMReX_ParticleTile.H:497
GpuArray< const int *, NArrayInt > m_idata
Definition: AMReX_ParticleTile.H:516
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE SuperParticleType getSuperParticle(int index) const noexcept
Definition: AMReX_ParticleTile.H:638
const ParticleReal *AMREX_RESTRICT *AMREX_RESTRICT m_runtime_rdata
Definition: AMReX_ParticleTile.H:582
AMREX_GPU_HOST_DEVICEdecltype(auto) AMREX_FORCE_INLINE id(const int index) const &
Definition: AMReX_ParticleTile.H:529
const int *AMREX_RESTRICT *AMREX_RESTRICT m_runtime_idata
Definition: AMReX_ParticleTile.H:583
Long m_size
Definition: AMReX_ParticleTile.H:508
T_ParticleType ParticleType
Definition: AMReX_ParticleTile.H:498
static constexpr int NStructReal
Definition: AMReX_ParticleTile.H:501
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE const ParticleReal * rdata(const int attribute_index) const
Definition: AMReX_ParticleTile.H:559
int m_num_runtime_real
Definition: AMReX_ParticleTile.H:580
static constexpr int NStructInt
Definition: AMReX_ParticleTile.H:502
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE const int * idata(const int attribute_index) const
Definition: AMReX_ParticleTile.H:565
int m_num_runtime_int
Definition: AMReX_ParticleTile.H:581
std::conditional_t< T_ParticleType::is_soa_particle, void const *AMREX_RESTRICT, ParticleType const *AMREX_RESTRICT > AOS_PTR
Definition: AMReX_ParticleTile.H:511
T_ParticleType const & ParticleRefType
Definition: AMReX_ParticleTile.H:499
static constexpr int NAR
Definition: AMReX_ParticleTile.H:496
AMREX_GPU_HOST_DEVICEdecltype(auto) AMREX_FORCE_INLINE idcpu(const int index) const &
Definition: AMReX_ParticleTile.H:549
AOS_PTR m_aos
Definition: AMReX_ParticleTile.H:512
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE const ParticleReal & pos(const int dir, const int index) const &
Definition: AMReX_ParticleTile.H:519
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void packParticleData(char *buffer, int src_index, Long dst_offset, const int *comm_real, const int *comm_int) const noexcept
Definition: AMReX_ParticleTile.H:586
const uint64_t * m_idcpu
Definition: AMReX_ParticleTile.H:514
AMREX_GPU_HOST_DEVICEdecltype(auto) AMREX_FORCE_INLINE cpu(const int index) const &
Definition: AMReX_ParticleTile.H:539
GpuArray< const ParticleReal *, NArrayReal > m_rdata
Definition: AMReX_ParticleTile.H:515
Definition: AMReX_ParticleTile.H:309
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE const RealType & pos(int position_index) const &
Definition: AMReX_ParticleTile.H:341
static constexpr int NArrayReal
Definition: AMReX_ParticleTile.H:310
ParticleReal RealType
Definition: AMReX_ParticleTile.H:317
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ConstParticleCPUWrapper cpu() const
Definition: AMReX_ParticleTile.H:330
static constexpr int NArrayInt
Definition: AMReX_ParticleTile.H:311
ConstPTD m_constparticle_tile_data
Definition: AMReX_ParticleTile.H:363
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE RealVect pos() const &
Definition: AMReX_ParticleTile.H:338
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ConstParticleIDWrapper id() const
Definition: AMReX_ParticleTile.H:333
static void NextID(Long nextid)
Reset on restart.
int m_index
Definition: AMReX_ParticleTile.H:366
static constexpr bool is_constsoa_particle
Definition: AMReX_ParticleTile.H:315
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ConstSoAParticle(ConstPTD const &ptd, long i)
Definition: AMReX_ParticleTile.H:320
static Long UnprotectedNextID()
This version can only be used inside omp critical.
uint64_t m_idcpu
Definition: AMReX_Particle.H:252
Definition: AMReX_Particle.H:140
Definition: AMReX_Particle.H:37
Definition: AMReX_ParticleTile.H:32
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void unpackParticleData(const char *buffer, Long src_offset, int dst_index, const int *comm_real, const int *comm_int) const noexcept
Definition: AMReX_ParticleTile.H:176
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void setSuperParticle(const SuperParticleType &sp, int index) const noexcept
Definition: AMReX_ParticleTile.H:271
T_ParticleType & ParticleRefType
Definition: AMReX_ParticleTile.H:37
uint64_t * m_idcpu
Definition: AMReX_ParticleTile.H:53
GpuArray< ParticleReal *, NAR > m_rdata
Definition: AMReX_ParticleTile.H:54
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE SuperParticleType getSuperParticle(int index) const noexcept
Definition: AMReX_ParticleTile.H:228
AMREX_GPU_HOST_DEVICEdecltype(auto) AMREX_FORCE_INLINE idcpu(const int index) const &
Definition: AMReX_ParticleTile.H:93
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ParticleReal & pos(const int dir, const int index) const &
Definition: AMReX_ParticleTile.H:63
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ParticleReal * rdata(const int attribute_index) const
Definition: AMReX_ParticleTile.H:103
static constexpr int NStructInt
Definition: AMReX_ParticleTile.H:41
GpuArray< int *, NAI > m_idata
Definition: AMReX_ParticleTile.H:55
T_ParticleType ParticleType
Definition: AMReX_ParticleTile.H:36
std::conditional_t< T_ParticleType::is_soa_particle, void *AMREX_RESTRICT, ParticleType *AMREX_RESTRICT > AOS_PTR
Definition: AMReX_ParticleTile.H:50
int m_num_runtime_int
Definition: AMReX_ParticleTile.H:58
Long m_size
Definition: AMReX_ParticleTile.H:47
static constexpr int NAR
Definition: AMReX_ParticleTile.H:33
ParticleReal *AMREX_RESTRICT *AMREX_RESTRICT m_runtime_rdata
Definition: AMReX_ParticleTile.H:59
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void packParticleData(char *buffer, int src_index, std::size_t dst_offset, const int *comm_real, const int *comm_int) const noexcept
Definition: AMReX_ParticleTile.H:125
int m_num_runtime_real
Definition: AMReX_ParticleTile.H:57
static constexpr bool is_particle_tile_data
Definition: AMReX_ParticleTile.H:45
AMREX_GPU_HOST_DEVICEdecltype(auto) AMREX_FORCE_INLINE cpu(const int index) const &
Definition: AMReX_ParticleTile.H:83
AOS_PTR m_aos
Definition: AMReX_ParticleTile.H:51
int *AMREX_RESTRICT *AMREX_RESTRICT m_runtime_idata
Definition: AMReX_ParticleTile.H:60
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE int * idata(const int attribute_index) const
Definition: AMReX_ParticleTile.H:109
static constexpr int NAI
Definition: AMReX_ParticleTile.H:34
static constexpr int NStructReal
Definition: AMReX_ParticleTile.H:40
AMREX_GPU_HOST_DEVICEdecltype(auto) AMREX_FORCE_INLINE id(const int index) const &
Definition: AMReX_ParticleTile.H:73
Definition: AMReX_ParticleTile.H:693
typename ParticleType::StorageParticleType StorageParticleType
Definition: AMReX_ParticleTile.H:719
int NumRealComps() const noexcept
Definition: AMReX_ParticleTile.H:1058
const AoS & GetArrayOfStructs() const
Definition: AMReX_ParticleTile.H:812
void push_back_real(int comp, amrex::Vector< amrex::ParticleReal > const &vec)
Definition: AMReX_ParticleTile.H:993
int NumIntComps() const noexcept
Definition: AMReX_ParticleTile.H:1060
typename ParticleType::RealType RealType
Definition: AMReX_ParticleTile.H:700
ParticleTile()=default
int getNumNeighbors() const
Definition: AMReX_ParticleTile.H:892
AoS & GetArrayOfStructs()
Definition: AMReX_ParticleTile.H:811
amrex::PODVector< const int *, Allocator< const int * > > m_runtime_i_cptrs
Definition: AMReX_ParticleTile.H:1291
void push_back(const ParticleType &p)
Definition: AMReX_ParticleTile.H:914
ParticleTile(ParticleTile const &)=delete
std::conditional_t< ParticleType::is_soa_particle, StructOfArrays< NArrayReal, NArrayInt, Allocator, true >, StructOfArrays< NArrayReal, NArrayInt, Allocator, false > > SoA
Definition: AMReX_ParticleTile.H:716
static constexpr int NAI
Definition: AMReX_ParticleTile.H:699
ParticleTile(ParticleTile &&) noexcept=default
void push_back_real(int comp, const ParticleReal *beg, const ParticleReal *end)
Definition: AMReX_ParticleTile.H:976
void define(int a_num_runtime_real, int a_num_runtime_int, std::vector< std::string > *soa_rdata_names=nullptr, std::vector< std::string > *soa_idata_names=nullptr)
Definition: AMReX_ParticleTile.H:736
static constexpr int NStructInt
Definition: AMReX_ParticleTile.H:703
amrex::Gpu::HostVector< const ParticleReal * > m_h_runtime_r_cptrs
Definition: AMReX_ParticleTile.H:1296
void setNumNeighbors(int num_neighbors)
Definition: AMReX_ParticleTile.H:884
amrex::Gpu::HostVector< int * > m_h_runtime_i_ptrs
Definition: AMReX_ParticleTile.H:1294
Long capacity() const
Definition: AMReX_ParticleTile.H:1086
int numTotalParticles() const
Returns the total number of particles, real and neighbor.
Definition: AMReX_ParticleTile.H:875
AoS m_aos_tile
Definition: AMReX_ParticleTile.H:1282
std::size_t size() const
Returns the total number of particles (real and neighbor)
Definition: AMReX_ParticleTile.H:823
ParticleTileDataType getParticleTileData()
Definition: AMReX_ParticleTile.H:1128
void push_back(const SuperParticleType &sp)
Definition: AMReX_ParticleTile.H:921
int numNeighborParticles() const
Returns the number of neighbor particles (excluding reals)
Definition: AMReX_ParticleTile.H:862
~ParticleTile()=default
amrex::PODVector< ParticleReal *, Allocator< ParticleReal * > > m_runtime_r_ptrs
Definition: AMReX_ParticleTile.H:1287
static constexpr int NAR
Definition: AMReX_ParticleTile.H:698
ConstParticleTileDataType getConstParticleTileData() const
Definition: AMReX_ParticleTile.H:1204
void shrink_to_fit()
Definition: AMReX_ParticleTile.H:1066
T_ParticleType ParticleType
Definition: AMReX_ParticleTile.H:697
int numParticles() const
Returns the number of real particles (excluding neighbors)
Definition: AMReX_ParticleTile.H:836
void push_back_int(int comp, amrex::Vector< int > const &vec)
Definition: AMReX_ParticleTile.H:1045
amrex::Gpu::HostVector< ParticleReal * > m_h_runtime_r_ptrs
Definition: AMReX_ParticleTile.H:1293
const SoA & GetStructOfArrays() const
Definition: AMReX_ParticleTile.H:815
void swap(ParticleTile< ParticleType, NArrayReal, NArrayInt, Allocator > &other) noexcept
Definition: AMReX_ParticleTile.H:1108
void push_back_int(int comp, amrex::Vector< int >::const_iterator beg, amrex::Vector< int >::const_iterator end)
Definition: AMReX_ParticleTile.H:1037
void push_back_real(const std::array< ParticleReal, NArrayReal > &v)
Definition: AMReX_ParticleTile.H:966
int NumRuntimeRealComps() const noexcept
Definition: AMReX_ParticleTile.H:1062
SoA & GetStructOfArrays()
Definition: AMReX_ParticleTile.H:814
amrex::Gpu::HostVector< const int * > m_h_runtime_i_cptrs
Definition: AMReX_ParticleTile.H:1297
amrex::PODVector< const ParticleReal *, Allocator< const ParticleReal * > > m_runtime_r_cptrs
Definition: AMReX_ParticleTile.H:1290
RealType & pos(int index, int position_index) &
Definition: AMReX_ParticleTile.H:788
void resize(std::size_t count)
Definition: AMReX_ParticleTile.H:902
typename SoA::IntVector IntVector
Definition: AMReX_ParticleTile.H:718
void push_back_int(const std::array< int, NArrayInt > &v)
Definition: AMReX_ParticleTile.H:1018
void push_back_real(int comp, amrex::Vector< amrex::ParticleReal >::const_iterator beg, amrex::Vector< amrex::ParticleReal >::const_iterator end)
Definition: AMReX_ParticleTile.H:985
Allocator< T > AllocatorType
Definition: AMReX_ParticleTile.H:695
int NumRuntimeIntComps() const noexcept
Definition: AMReX_ParticleTile.H:1064
bool empty() const
Definition: AMReX_ParticleTile.H:817
void push_back_int(int comp, int v)
Definition: AMReX_ParticleTile.H:1010
RealType pos(int index, int position_index) const &
Definition: AMReX_ParticleTile.H:802
decltype(auto) cpu(int index) &
Definition: AMReX_ParticleTile.H:770
SoA m_soa_tile
Definition: AMReX_ParticleTile.H:1283
void push_back_int(int comp, std::size_t npar, int v)
Definition: AMReX_ParticleTile.H:1053
bool m_defined
Definition: AMReX_ParticleTile.H:1285
void push_back_real(int comp, std::size_t npar, ParticleReal v)
Definition: AMReX_ParticleTile.H:1001
void push_back_int(int comp, const int *beg, const int *end)
Definition: AMReX_ParticleTile.H:1028
typename SoA::RealVector RealVector
Definition: AMReX_ParticleTile.H:717
amrex::PODVector< int *, Allocator< int * > > m_runtime_i_ptrs
Definition: AMReX_ParticleTile.H:1288
int numRealParticles() const
Returns the number of real particles (excluding neighbors)
Definition: AMReX_ParticleTile.H:849
std::conditional_t< ParticleType::is_soa_particle, ThisParticleTileHasNoAoS, ArrayOfStructs< ParticleType, Allocator > > AoS
Definition: AMReX_ParticleTile.H:710
decltype(auto) id(int index) &
Definition: AMReX_ParticleTile.H:752
static constexpr int NStructReal
Definition: AMReX_ParticleTile.H:702
void push_back_real(int comp, ParticleReal v)
Definition: AMReX_ParticleTile.H:958
The struct used to store particles.
Definition: AMReX_Particle.H:295
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ParticleIDWrapper id() &
Definition: AMReX_Particle.H:315
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE RealVect pos() const &
Definition: AMReX_Particle.H:338
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE RealType & rdata(int index) &
Definition: AMReX_Particle.H:356
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE int & idata(int index) &
Definition: AMReX_Particle.H:427
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ParticleCPUWrapper cpu() &
Definition: AMReX_Particle.H:312
Definition: AMReX_Particle.H:281
Definition: AMReX_ParticleTile.H:371
static constexpr int NArrayReal
Definition: AMReX_ParticleTile.H:372
static Long UnprotectedNextID()
This version can only be used inside omp critical.
Definition: AMReX_ParticleTile.H:477
PTD m_particle_tile_data
Definition: AMReX_ParticleTile.H:445
static Long the_next_id
Definition: AMReX_ParticleTile.H:388
static constexpr int NArrayInt
Definition: AMReX_ParticleTile.H:373
int m_index
Definition: AMReX_ParticleTile.H:448
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ParticleCPUWrapper cpu() &
Definition: AMReX_ParticleTile.H:393
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE const uint64_t & idcpu() const &
Definition: AMReX_ParticleTile.H:408
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ConstParticleCPUWrapper cpu() const &
Definition: AMReX_ParticleTile.H:402
ParticleReal RealType
Definition: AMReX_ParticleTile.H:380
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE RealType pos(int position_index) const &
Definition: AMReX_ParticleTile.H:423
static constexpr bool is_constsoa_particle
Definition: AMReX_ParticleTile.H:377
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ParticleIDWrapper id() &
Definition: AMReX_ParticleTile.H:396
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE SoAParticle(PTD const &ptd, long i)
Definition: AMReX_ParticleTile.H:383
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE RealVect pos() const &
Definition: AMReX_ParticleTile.H:413
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE RealType & pos(int position_index) &
Definition: AMReX_ParticleTile.H:416
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ConstParticleIDWrapper id() const &
Definition: AMReX_ParticleTile.H:405
static Long NextID()
Definition: AMReX_ParticleTile.H:456
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE uint64_t & idcpu() &
Definition: AMReX_ParticleTile.H:399
Definition: AMReX_StructOfArrays.H:20
Definition: AMReX_ParticleTile.H:686
Definition: AMReX_ParticleTile.H:684
Definition: AMReX_MakeParticle.H:11