Block-Structured AMR Software Framework
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
Loading...
Searching...
No Matches
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>
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
19namespace amrex {
20
21// Forward Declaration
22template <int NArrayReal, int NArrayInt>
23struct ConstSoAParticle;
24template <int NArrayReal, int NArrayInt>
25struct SoAParticle;
26
27template <typename T_ParticleType, int NArrayReal, int NArrayInt>
28struct ConstParticleTileData;
29
30template <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 = 0;
138 if constexpr (!ParticleType::is_soa_particle) {
139 array_start_index = AMREX_SPACEDIM + NStructReal;
140 }
141 for (int i = 0; i < NAR; ++i)
142 {
143 if (comm_real[array_start_index + i])
144 {
145 memcpy(dst, m_rdata[i] + src_index, sizeof(ParticleReal));
146 dst += sizeof(ParticleReal);
147 }
148 }
149 int runtime_start_index = array_start_index + NAR;
150 for (int i = 0; i < m_num_runtime_real; ++i)
151 {
152 if (comm_real[runtime_start_index + i])
153 {
154 memcpy(dst, m_runtime_rdata[i] + src_index, sizeof(ParticleReal));
155 dst += sizeof(ParticleReal);
156 }
157 }
158 array_start_index = 2 + NStructInt;
159 for (int i = 0; i < NAI; ++i)
160 {
161 if (comm_int[array_start_index + i])
162 {
163 memcpy(dst, m_idata[i] + src_index, sizeof(int));
164 dst += sizeof(int);
165 }
166 }
167 runtime_start_index = 2 + NStructInt + NAI;
168 for (int i = 0; i < m_num_runtime_int; ++i)
169 {
170 if (comm_int[runtime_start_index + i])
171 {
172 memcpy(dst, m_runtime_idata[i] + src_index, sizeof(int));
173 dst += sizeof(int);
174 }
175 }
176 }
177
179 void unpackParticleData (const char* buffer, Long src_offset, int dst_index,
180 const int* comm_real, const int* comm_int) const noexcept
181 {
182 AMREX_ASSERT(dst_index < m_size);
183 const auto* src = buffer + src_offset;
184 if constexpr (!ParticleType::is_soa_particle) {
185 memcpy(m_aos + dst_index, src, sizeof(ParticleType));
186 src += sizeof(ParticleType);
187 } else {
188 memcpy(m_idcpu + dst_index, src, sizeof(uint64_t));
189 src += sizeof(uint64_t);
190 }
191 int array_start_index = 0;
192 if constexpr (!ParticleType::is_soa_particle) {
193 array_start_index = AMREX_SPACEDIM + NStructReal;
194 }
195 for (int i = 0; i < NAR; ++i)
196 {
197 if (comm_real[array_start_index + i])
198 {
199 memcpy(m_rdata[i] + dst_index, src, sizeof(ParticleReal));
200 src += sizeof(ParticleReal);
201 }
202 }
203 int runtime_start_index = array_start_index + NAR;
204 for (int i = 0; i < m_num_runtime_real; ++i)
205 {
206 if (comm_real[runtime_start_index + i])
207 {
208 memcpy(m_runtime_rdata[i] + dst_index, src, sizeof(ParticleReal));
209 src += sizeof(ParticleReal);
210 }
211 }
212 array_start_index = 2 + NStructInt;
213 for (int i = 0; i < NAI; ++i)
214 {
215 if (comm_int[array_start_index + i])
216 {
217 memcpy(m_idata[i] + dst_index, src, sizeof(int));
218 src += sizeof(int);
219 }
220 }
221 runtime_start_index = 2 + NStructInt + NAI;
222 for (int i = 0; i < m_num_runtime_int; ++i)
223 {
224 if (comm_int[runtime_start_index + i])
225 {
226 memcpy(m_runtime_idata[i] + dst_index, src, sizeof(int));
227 src += sizeof(int);
228 }
229 }
230 }
231
232 template <typename T = ParticleType, std::enable_if_t<!T::is_soa_particle, int> = 0>
234 SuperParticleType getSuperParticle (int index) const noexcept
235 {
236 AMREX_ASSERT(index < m_size);
238 for (int i = 0; i < AMREX_SPACEDIM; ++i) {
239 sp.pos(i) = m_aos[index].pos(i);
240 }
241 for (int i = 0; i < NStructReal; ++i) {
242 sp.rdata(i) = m_aos[index].rdata(i);
243 }
244 for (int i = 0; i < NAR; ++i) {
245 sp.rdata(NStructReal+i) = m_rdata[i][index];
246 }
247 sp.id() = m_aos[index].id();
248 sp.cpu() = m_aos[index].cpu();
249 for (int i = 0; i < NStructInt; ++i) {
250 sp.idata(i) = m_aos[index].idata(i);
251 }
252 for (int i = 0; i < NAI; ++i) {
253 sp.idata(NStructInt+i) = m_idata[i][index];
254 }
255 return sp;
256 }
257
258 template <typename T = ParticleType, std::enable_if_t<T::is_soa_particle, int> = 0>
260 SuperParticleType getSuperParticle (int index) const noexcept
261 {
262 AMREX_ASSERT(index < m_size);
264 sp.m_idcpu = m_idcpu[index];
265 for (int i = 0; i < AMREX_SPACEDIM; ++i) {sp.pos(i) = m_rdata[i][index];}
266 for (int i = 0; i < NAR; ++i) {
267 sp.rdata(i) = m_rdata[i][index];
268 }
269 for (int i = 0; i < NAI; ++i) {
270 sp.idata(i) = m_idata[i][index];
271 }
272 return sp;
273 }
274
275 template <typename T = ParticleType, std::enable_if_t<!T::is_soa_particle, int> = 0>
277 void setSuperParticle (const SuperParticleType& sp, int index) const noexcept
278 {
279 for (int i = 0; i < AMREX_SPACEDIM; ++i) {
280 m_aos[index].pos(i) = sp.pos(i);
281 }
282 for (int i = 0; i < NStructReal; ++i) {
283 m_aos[index].rdata(i) = sp.rdata(i);
284 }
285 for (int i = 0; i < NAR; ++i) {
286 m_rdata[i][index] = sp.rdata(NStructReal+i);
287 }
288 m_aos[index].id() = sp.id();
289 m_aos[index].cpu() = sp.cpu();
290 for (int i = 0; i < NStructInt; ++i) {
291 m_aos[index].idata(i) = sp.idata(i);
292 }
293 for (int i = 0; i < NAI; ++i) {
294 m_idata[i][index] = sp.idata(NStructInt+i);
295 }
296 }
297
298 template <typename T = ParticleType, std::enable_if_t<T::is_soa_particle, int> = 0>
300 void setSuperParticle (const SuperParticleType& sp, int index) const noexcept
301 {
302 m_idcpu[index] = sp.m_idcpu;
303 for (int i = 0; i < NAR; ++i) {
304 m_rdata[i][index] = sp.rdata(i);
305 }
306 for (int i = 0; i < NAI; ++i) {
307 m_idata[i][index] = sp.idata(i);
308 }
309 }
310};
311
312// SOA Particle Structure
313template <int T_NArrayReal, int T_NArrayInt>
315{
316 static constexpr int NArrayReal = T_NArrayReal;
317 static constexpr int NArrayInt = T_NArrayInt;
320 static constexpr bool is_soa_particle = true;
321 static constexpr bool is_constsoa_particle = true;
322
323 using RealType = ParticleReal;
324
326 ConstSoAParticle (ConstPTD const& ptd, long i) : // Note: should this be int instead?
328 {
329 }
330
331 //static Long the_next_id;
332
333 //functions to get id and cpu in the SOA data
334
337
340
341 //functions to get positions of the particle in the SOA data
342
344 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]));}
345
347 const RealType& pos (int position_index) const &
348 {
349 AMREX_ASSERT(position_index < AMREX_SPACEDIM);
350 return this->m_constparticle_tile_data.m_rdata[position_index][m_index];
351 }
352
353 static Long NextID ();
354
358 static Long UnprotectedNextID ();
359
365 static void NextID (Long nextid);
366
367 private :
368
369 static_assert(std::is_trivially_copyable<ParticleTileData<SoAParticleBase, NArrayReal, NArrayInt>>(), "ParticleTileData is not trivially copyable");
370
373};
374
375template <int T_NArrayReal, int T_NArrayInt>
377{
378 static constexpr int NArrayReal = T_NArrayReal;
379 static constexpr int NArrayInt = T_NArrayInt;
382 static constexpr bool is_soa_particle = true;
383 static constexpr bool is_constsoa_particle = false;
384
386 using RealType = ParticleReal;
387
389 SoAParticle (PTD const& ptd, long i) : // Note: should this be int instead?
391 {
392 }
393
394 static Long the_next_id;
395
396 //functions to get id and cpu in the SOA data
397
400
403
405 uint64_t& idcpu () & { return this->m_particle_tile_data.m_idcpu[m_index]; }
406
409
412
414 const uint64_t& idcpu () const & { return this->m_particle_tile_data.m_idcpu[m_index]; }
415
416 //functions to get positions of the particle in the SOA data
417
419 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]));}
420
422 RealType& pos (int position_index) &
423 {
424 AMREX_ASSERT(position_index < AMREX_SPACEDIM);
425 return this->m_particle_tile_data.m_rdata[position_index][m_index];
426 }
427
429 RealType pos (int position_index) const &
430 {
431 AMREX_ASSERT(position_index < AMREX_SPACEDIM);
432 return this->m_particle_tile_data.m_rdata[position_index][m_index];
433 }
434
435 static Long NextID ();
436
440 static Long UnprotectedNextID ();
441
447 static void NextID (Long nextid);
448
449private :
450
451 static_assert(std::is_trivially_copyable<ParticleTileData<SoAParticleBase, NArrayReal, NArrayInt>>(), "ParticleTileData is not trivially copyable");
452
455};
456
457//template <int NArrayReal, int NArrayInt> Long ConstSoAParticle<NArrayReal, NArrayInt>::the_next_id = 1;
458template <int NArrayReal, int NArrayInt> Long SoAParticle<NArrayReal, NArrayInt>::the_next_id = 1;
459
460template <int NArrayReal, int NArrayInt>
461Long
463{
464 Long next;
465// we should be able to test on _OPENMP < 201107 for capture (version 3.1)
466// but we must work around a bug in gcc < 4.9
467#if defined(AMREX_USE_OMP) && defined(_OPENMP) && _OPENMP < 201307
468#pragma omp critical (amrex_particle_nextid)
469#elif defined(AMREX_USE_OMP)
470#pragma omp atomic capture
471#endif
472 next = the_next_id++;
473
475 amrex::Abort("SoAParticle<NArrayReal, NArrayInt>::NextID() -- too many particles");
476 }
477
478 return next;
479}
480
481template <int NArrayReal, int NArrayInt>
482Long
484{
485 Long next = the_next_id++;
487 amrex::Abort("SoAParticle<NArrayReal, NArrayInt>::NextID() -- too many particles");
488 }
489 return next;
490}
491
492template <int NArrayReal, int NArrayInt>
493void
495{
496 the_next_id = nextid;
497}
498
499template <typename T_ParticleType, int NArrayReal, int NArrayInt>
501{
502 static constexpr int NAR = NArrayReal;
503 static constexpr int NAI = NArrayInt;
504 using ParticleType = T_ParticleType;
505 using ParticleRefType = T_ParticleType const&;
506
507 static constexpr int NStructReal = ParticleType::NReal;
508 static constexpr int NStructInt = ParticleType::NInt;
509
511
512 static constexpr bool is_particle_tile_data = true;
513
514 Long m_size;
515
516 using AOS_PTR = std::conditional_t<T_ParticleType::is_soa_particle,
517 void const * AMREX_RESTRICT, ParticleType const * AMREX_RESTRICT>;
519
520 const uint64_t* m_idcpu;
523
525 const ParticleReal& pos (const int dir, const int index) const &
526 {
527 if constexpr(!ParticleType::is_soa_particle) {
528 return this->m_aos[index].pos(dir);
529 } else {
530 return this->m_rdata[dir][index];
531 }
532 }
533
535 decltype(auto) id (const int index) const &
536 {
537 if constexpr(!ParticleType::is_soa_particle) {
538 return this->m_aos[index].id();
539 } else {
540 return ConstParticleIDWrapper(this->m_idcpu[index]);
541 }
542 }
543
545 decltype(auto) cpu (const int index) const &
546 {
547 if constexpr(!ParticleType::is_soa_particle) {
548 return this->m_aos[index].cpu();
549 } else {
550 return ConstParticleCPUWrapper(this->m_idcpu[index]);
551 }
552 }
553
555 decltype(auto) idcpu (const int index) const &
556 {
557 if constexpr(ParticleType::is_soa_particle) {
558 return this->m_idcpu[index];
559 } else {
560 amrex::Abort("not implemented");
561 }
562 }
563
565 const ParticleReal * rdata (const int attribute_index) const
566 {
567 return this->m_rdata[attribute_index];
568 }
569
571 const int * idata (const int attribute_index) const
572 {
573 return this->m_idata[attribute_index];
574 }
575
577 decltype(auto) operator[] (const int index) const
578 {
579 if constexpr (!ParticleType::is_soa_particle) {
580 return m_aos[index];
581 } else {
582 return ConstSoAParticle<NAR, NAI>(*this, index);
583 }
584 }
585
590
592 void packParticleData(char* buffer, int src_index, Long dst_offset,
593 const int* comm_real, const int * comm_int) const noexcept
594 {
595 AMREX_ASSERT(src_index < m_size);
596 auto* dst = buffer + dst_offset;
597 if constexpr (!ParticleType::is_soa_particle) {
598 memcpy(dst, m_aos + src_index, sizeof(ParticleType));
599 dst += sizeof(ParticleType);
600 } else {
601 memcpy(dst, m_idcpu + src_index, sizeof(uint64_t));
602 dst += sizeof(uint64_t);
603 }
604 int array_start_index = 0;
605 if constexpr (!ParticleType::is_soa_particle) {
606 array_start_index = AMREX_SPACEDIM + NStructReal;
607 }
608 for (int i = 0; i < NArrayReal; ++i)
609 {
610 if (comm_real[array_start_index + i])
611 {
612 memcpy(dst, m_rdata[i] + src_index, sizeof(ParticleReal));
613 dst += sizeof(ParticleReal);
614 }
615 }
616 int runtime_start_index = array_start_index + NArrayReal;
617 for (int i = 0; i < m_num_runtime_real; ++i)
618 {
619 if (comm_real[runtime_start_index + i])
620 {
621 memcpy(dst, m_runtime_rdata[i] + src_index, sizeof(ParticleReal));
622 dst += sizeof(ParticleReal);
623 }
624 }
625 array_start_index = 2 + NStructInt;
626 for (int i = 0; i < NArrayInt; ++i)
627 {
628 if (comm_int[array_start_index + i])
629 {
630 memcpy(dst, m_idata[i] + src_index, sizeof(int));
631 dst += sizeof(int);
632 }
633 }
634 runtime_start_index = 2 + NStructInt + NArrayInt;
635 for (int i = 0; i < m_num_runtime_int; ++i)
636 {
637 if (comm_int[runtime_start_index + i])
638 {
639 memcpy(dst, m_runtime_idata[i] + src_index, sizeof(int));
640 dst += sizeof(int);
641 }
642 }
643 }
644
645 template <typename T = ParticleType, std::enable_if_t<!T::is_soa_particle, int> = 0>
647 SuperParticleType getSuperParticle (int index) const noexcept
648 {
649 AMREX_ASSERT(index < m_size);
651 for (int i = 0; i < AMREX_SPACEDIM; ++i) {
652 sp.pos(i) = m_aos[index].pos(i);
653 }
654 for (int i = 0; i < NStructReal; ++i) {
655 sp.rdata(i) = m_aos[index].rdata(i);
656 }
657 if constexpr(NArrayReal > 0) {
658 for (int i = 0; i < NArrayReal; ++i) {
659 sp.rdata(NStructReal+i) = m_rdata[i][index];
660 }
661 }
662 sp.id() = m_aos[index].id();
663 sp.cpu() = m_aos[index].cpu();
664 for (int i = 0; i < NStructInt; ++i) {
665 sp.idata(i) = m_aos[index].idata(i);
666 }
667 if constexpr(NArrayInt > 0) {
668 for (int i = 0; i < NArrayInt; ++i) {
669 sp.idata(NStructInt+i) = m_idata[i][index];
670 }
671 }
672 return sp;
673 }
674
675 template <typename T = ParticleType, std::enable_if_t<T::is_soa_particle, int> = 0>
677 SuperParticleType getSuperParticle (int index) const noexcept
678 {
679 AMREX_ASSERT(index < m_size);
681 for (int i = 0; i < AMREX_SPACEDIM; ++i) {sp.pos(i) = m_rdata[i][index];}
682 sp.m_idcpu = m_idcpu[index];
683 for (int i = 0; i < NAR; ++i) {
684 sp.rdata(i) = m_rdata[i][index];
685 }
686 for (int i = 0; i < NAI; ++i) {
687 sp.idata(i) = m_idata[i][index];
688 }
689 return sp;
690 }
691};
692
694
698
699template <typename T_ParticleType, int NArrayReal, int NArrayInt,
700 template<class> class Allocator=DefaultAllocator>
702{
703 template <typename T>
704 using AllocatorType = Allocator<T>;
705
706 using ParticleType = T_ParticleType;
707 static constexpr int NAR = NArrayReal;
708 static constexpr int NAI = NArrayInt;
709 using RealType = typename ParticleType::RealType;
710
711 static constexpr int NStructReal = ParticleType::NReal;
712 static constexpr int NStructInt = ParticleType::NInt;
713
715
716 using AoS = std::conditional_t<
717 ParticleType::is_soa_particle,
720 //using ParticleVector = typename AoS::ParticleVector;
721
722 using SoA = std::conditional_t<
723 ParticleType::is_soa_particle,
726 using RealVector = typename SoA::RealVector;
727 using IntVector = typename SoA::IntVector;
728 using StorageParticleType = typename ParticleType::StorageParticleType;
729
732
733 ParticleTile () = default;
734
735#ifndef _WIN32 // workaround windows compiler bug
736 ~ParticleTile () = default;
737
738 ParticleTile (ParticleTile const&) = delete;
739 ParticleTile (ParticleTile &&) noexcept = default;
740
741 ParticleTile& operator= (ParticleTile const&) = delete;
742 ParticleTile& operator= (ParticleTile &&) noexcept = default;
743#endif
744
745 void define (
746 int a_num_runtime_real,
747 int a_num_runtime_int,
748 std::vector<std::string>* soa_rdata_names=nullptr,
749 std::vector<std::string>* soa_idata_names=nullptr
750 )
751 {
752 m_defined = true;
753 GetStructOfArrays().define(a_num_runtime_real, a_num_runtime_int, soa_rdata_names, soa_idata_names);
754 m_runtime_r_ptrs.resize(a_num_runtime_real);
755 m_runtime_i_ptrs.resize(a_num_runtime_int);
756 m_runtime_r_cptrs.resize(a_num_runtime_real);
757 m_runtime_i_cptrs.resize(a_num_runtime_int);
758 }
759
760 // Get id data
761 decltype(auto) id (int index) & {
762 if constexpr (!ParticleType::is_soa_particle) {
763 return m_aos_tile[index].id();
764 } else {
765 return ParticleIDWrapper(m_soa_tile.GetIdCPUData()[index]);
766 }
767 }
768
769 // const
770 decltype(auto) id (int index) const & {
771 if constexpr (!ParticleType::is_soa_particle) {
772 return m_aos_tile[index].id();
773 } else {
774 return ConstParticleIDWrapper(m_soa_tile.GetIdCPUData()[index]);
775 }
776 }
777
778 // Get cpu data
779 decltype(auto) cpu (int index) & {
780 if constexpr (!ParticleType::is_soa_particle) {
781 return m_aos_tile[index].cpu();
782 } else {
783 return ParticleCPUWrapper(m_soa_tile.GetIdCPUData()[index]);
784 }
785 }
786
787 // const
788 decltype(auto) cpu (int index) const & {
789 if constexpr (!ParticleType::is_soa_particle) {
790 return m_aos_tile[index].cpu();
791 } else {
792 return ConstParticleCPUWrapper(m_soa_tile.GetIdCPUData()[index]);
793 }
794 }
795
796 // Get positions data
797 RealType& pos (int index, int position_index) & {
798 if constexpr (!ParticleType::is_soa_particle) {
799 return m_aos_tile[index].pos(position_index);
800 } else {
801 static_assert(NArrayReal == ParticleType::PTD::NAR, "ParticleTile mismatch in R");
802 static_assert(NArrayInt == ParticleType::PTD::NAI, "ParticleTile mismatch in I");
803 static_assert(0 == ParticleType::StorageParticleType::NReal, "ParticleTile 2 mismatch in R");
804 static_assert(0 == ParticleType::StorageParticleType::NInt, "ParticleTile 2 mismatch in I");
805
806 return m_soa_tile.GetRealData(position_index)[index];
807 }
808 }
809
810 // const
811 RealType pos (int index, int position_index) const &
812 {
813 if constexpr (!ParticleType::is_soa_particle) {
814 return m_aos_tile[index].pos(position_index);
815 } else {
816 return m_soa_tile.GetRealData(position_index)[index];
817 }
818 }
819
821 const AoS& GetArrayOfStructs () const { return m_aos_tile; }
822
824 const SoA& GetStructOfArrays () const { return m_soa_tile; }
825
826 bool empty () const { return size() == 0; }
827
832 std::size_t size () const
833 {
834 if constexpr (!ParticleType::is_soa_particle) {
835 return m_aos_tile.size();
836 } else {
837 return m_soa_tile.size();
838 }
839 }
840
845 int numParticles () const
846 {
847 if constexpr (!ParticleType::is_soa_particle) {
848 return m_aos_tile.numParticles();
849 } else {
850 return m_soa_tile.numParticles();
851 }
852 }
853
858 int numRealParticles () const
859 {
860 if constexpr (!ParticleType::is_soa_particle) {
861 return m_aos_tile.numRealParticles();
862 } else {
863 return m_soa_tile.numRealParticles();
864 }
865 }
866
872 {
873 if constexpr (!ParticleType::is_soa_particle) {
874 return m_aos_tile.numNeighborParticles();
875 } else {
876 return m_soa_tile.numNeighborParticles();
877 }
878 }
879
884 int numTotalParticles () const
885 {
886 if constexpr (!ParticleType::is_soa_particle) {
887 return m_aos_tile.numTotalParticles();
888 } else {
889 return m_soa_tile.numTotalParticles();
890 }
891 }
892
893 void setNumNeighbors (int num_neighbors)
894 {
895 if constexpr(!ParticleType::is_soa_particle) {
896 m_aos_tile.setNumNeighbors(num_neighbors);
897 }
898 m_soa_tile.setNumNeighbors(num_neighbors);
899 }
900
901 int getNumNeighbors () const
902 {
903 if constexpr (!ParticleType::is_soa_particle) {
904 AMREX_ASSERT( m_soa_tile.getNumNeighbors() == m_aos_tile.getNumNeighbors() );
905 return m_aos_tile.getNumNeighbors();
906 } else {
907 return m_soa_tile.getNumNeighbors();
908 }
909 }
910
911 void resize (std::size_t count)
912 {
913 if constexpr (!ParticleType::is_soa_particle) {
914 m_aos_tile.resize(count);
915 }
916 m_soa_tile.resize(count);
917 }
918
922 template <typename T = ParticleType, std::enable_if_t<!T::is_soa_particle, int> = 0>
923 void push_back (const ParticleType& p) { m_aos_tile().push_back(p); }
924
928 template < int NR = NArrayReal, int NI = NArrayInt,
929 std::enable_if_t<NR != 0 || NI != 0, int> foo = 0>
931 {
932 auto np = numParticles();
933
934 if constexpr (!ParticleType::is_soa_particle) {
935 m_aos_tile.resize(np+1);
936 for (int i = 0; i < AMREX_SPACEDIM; ++i) {
937 m_aos_tile[np].pos(i) = sp.pos(i);
938 }
939 for (int i = 0; i < NStructReal; ++i) {
940 m_aos_tile[np].rdata(i) = sp.rdata(i);
941 }
942 m_aos_tile[np].id() = sp.id();
943 m_aos_tile[np].cpu() = sp.cpu();
944 for (int i = 0; i < NStructInt; ++i) {
945 m_aos_tile[np].idata(i) = sp.idata(i);
946 }
947 }
948
949 m_soa_tile.resize(np+1);
950 if constexpr (ParticleType::is_soa_particle) {
951 m_soa_tile.GetIdCPUData()[np] = sp.m_idcpu;
952 }
953 auto& arr_rdata = m_soa_tile.GetRealData();
954 auto& arr_idata = m_soa_tile.GetIntData();
955 for (int i = 0; i < NArrayReal; ++i) {
956 arr_rdata[i][np] = sp.rdata(NStructReal+i);
957 }
958 for (int i = 0; i < NArrayInt; ++i) {
959 arr_idata[i][np] = sp.idata(NStructInt+i);
960 }
961 }
962
967 void push_back_real (int comp, ParticleReal v) {
968 m_soa_tile.GetRealData(comp).push_back(v);
969 }
970
975 void push_back_real (const std::array<ParticleReal, NArrayReal>& v) {
976 for (int i = 0; i < NArrayReal; ++i) {
977 m_soa_tile.GetRealData(i).push_back(v[i]);
978 }
979 }
980
985 void push_back_real (int comp, const ParticleReal* beg, const ParticleReal* end) {
986 auto it = m_soa_tile.GetRealData(comp).end();
987 m_soa_tile.GetRealData(comp).insert(it, beg, end);
988 }
989
997
1003 push_back_real(comp, vec.cbegin(), vec.cend());
1004 }
1005
1010 void push_back_real (int comp, std::size_t npar, ParticleReal v) {
1011 auto new_size = m_soa_tile.GetRealData(comp).size() + npar;
1012 m_soa_tile.GetRealData(comp).resize(new_size, v);
1013 }
1014
1019 void push_back_int (int comp, int v) {
1020 m_soa_tile.GetIntData(comp).push_back(v);
1021 }
1022
1027 void push_back_int (const std::array<int, NArrayInt>& v) {
1028 for (int i = 0; i < NArrayInt; ++i) {
1029 m_soa_tile.GetIntData(i).push_back(v[i]);
1030 }
1031 }
1032
1037 void push_back_int (int comp, const int* beg, const int* end) {
1038 auto it = m_soa_tile.GetIntData(comp).end();
1039 m_soa_tile.GetIntData(comp).insert(it, beg, end);
1040 }
1041
1049
1054 void push_back_int (int comp, amrex::Vector<int> const & vec) {
1055 push_back_int(comp, vec.cbegin(), vec.cend());
1056 }
1057
1062 void push_back_int (int comp, std::size_t npar, int v) {
1063 auto new_size = m_soa_tile.GetIntData(comp).size() + npar;
1064 m_soa_tile.GetIntData(comp).resize(new_size, v);
1065 }
1066
1067 int NumRealComps () const noexcept { return m_soa_tile.NumRealComps(); }
1068
1069 int NumIntComps () const noexcept { return m_soa_tile.NumIntComps(); }
1070
1071 int NumRuntimeRealComps () const noexcept { return m_runtime_r_ptrs.size(); }
1072
1073 int NumRuntimeIntComps () const noexcept { return m_runtime_i_ptrs.size(); }
1074
1076 {
1077 if constexpr (ParticleType::is_soa_particle) {
1078 GetStructOfArrays().GetIdCPUData().shrink_to_fit();
1079 } else {
1080 m_aos_tile().shrink_to_fit();
1081 }
1082 for (int j = 0; j < NumRealComps(); ++j)
1083 {
1084 auto& rdata = GetStructOfArrays().GetRealData(j);
1085 rdata.shrink_to_fit();
1086 }
1087
1088 for (int j = 0; j < NumIntComps(); ++j)
1089 {
1090 auto& idata = GetStructOfArrays().GetIntData(j);
1091 idata.shrink_to_fit();
1092 }
1093 }
1094
1095 Long capacity () const
1096 {
1097 Long nbytes = 0;
1098 if constexpr (ParticleType::is_soa_particle) {
1099 nbytes += GetStructOfArrays().GetIdCPUData().capacity() * sizeof(uint64_t);
1100 } else {
1101 nbytes += m_aos_tile().capacity() * sizeof(ParticleType);
1102 }
1103 for (int j = 0; j < NumRealComps(); ++j)
1104 {
1105 auto& rdata = GetStructOfArrays().GetRealData(j);
1106 nbytes += rdata.capacity() * sizeof(ParticleReal);
1107 }
1108
1109 for (int j = 0; j < NumIntComps(); ++j)
1110 {
1111 auto& idata = GetStructOfArrays().GetIntData(j);
1112 nbytes += idata.capacity()*sizeof(int);
1113 }
1114 return nbytes;
1115 }
1116
1118 {
1119 if constexpr (ParticleType::is_soa_particle) {
1120 GetStructOfArrays().GetIdCPUData().swap(other.GetStructOfArrays().GetIdCPUData());
1121 } else {
1122 m_aos_tile().swap(other.GetArrayOfStructs()());
1123 }
1124 for (int j = 0; j < NumRealComps(); ++j)
1125 {
1126 auto& rdata = GetStructOfArrays().GetRealData(j);
1127 rdata.swap(other.GetStructOfArrays().GetRealData(j));
1128 }
1129
1130 for (int j = 0; j < NumIntComps(); ++j)
1131 {
1132 auto& idata = GetStructOfArrays().GetIntData(j);
1133 idata.swap(other.GetStructOfArrays().GetIntData(j));
1134 }
1135 }
1136
1138 {
1139 m_runtime_r_ptrs.resize(m_soa_tile.NumRealComps() - NArrayReal);
1140 m_runtime_i_ptrs.resize(m_soa_tile.NumIntComps() - NArrayInt);
1141#ifdef AMREX_USE_GPU
1142 bool copy_real = false;
1143 m_h_runtime_r_ptrs.resize(m_soa_tile.NumRealComps() - NArrayReal, nullptr);
1144 for (std::size_t i = 0; i < m_h_runtime_r_ptrs.size(); ++i) {
1145 if (m_h_runtime_r_ptrs[i] != m_soa_tile.GetRealData(i + NArrayReal).dataPtr()) {
1146 m_h_runtime_r_ptrs[i] = m_soa_tile.GetRealData(i + NArrayReal).dataPtr();
1147 copy_real = true;
1148 }
1149 }
1150 if (copy_real) {
1152 m_h_runtime_r_ptrs.size()*sizeof(ParticleReal*));
1153 }
1154
1155 bool copy_int = false;
1156 m_h_runtime_i_ptrs.resize(m_soa_tile.NumIntComps() - NArrayInt, nullptr);
1157 for (std::size_t i = 0; i < m_h_runtime_i_ptrs.size(); ++i) {
1158 if (m_h_runtime_i_ptrs[i] != m_soa_tile.GetIntData(i + NArrayInt).dataPtr()) {
1159 m_h_runtime_i_ptrs[i] = m_soa_tile.GetIntData(i + NArrayInt).dataPtr();
1160 copy_int = true;
1161 }
1162 }
1163 if (copy_int) {
1165 m_h_runtime_i_ptrs.size()*sizeof(int*));
1166 }
1167#else
1168 for (std::size_t i = 0; i < m_runtime_r_ptrs.size(); ++i) {
1169 m_runtime_r_ptrs[i] = m_soa_tile.GetRealData(i + NArrayReal).dataPtr();
1170 }
1171
1172 for (std::size_t i = 0; i < m_runtime_i_ptrs.size(); ++i) {
1173 m_runtime_i_ptrs[i] = m_soa_tile.GetIntData(i + NArrayInt).dataPtr();
1174 }
1175#endif
1176
1178 if constexpr (!ParticleType::is_soa_particle) {
1179 ptd.m_aos = m_aos_tile().dataPtr();
1180 } else {
1181 ptd.m_aos = nullptr;
1182 }
1183 if constexpr (ParticleType::is_soa_particle) {
1184 ptd.m_idcpu = m_soa_tile.GetIdCPUData().dataPtr();
1185 } else {
1186 ptd.m_idcpu = nullptr;
1187 }
1188 if constexpr(NArrayReal > 0) {
1189 for (int i = 0; i < NArrayReal; ++i) {
1190 ptd.m_rdata[i] = m_soa_tile.GetRealData(i).dataPtr();
1191 }
1192 }
1193 if constexpr(NArrayInt > 0) {
1194 for (int i = 0; i < NArrayInt; ++i) {
1195 ptd.m_idata[i] = m_soa_tile.GetIntData(i).dataPtr();
1196 }
1197 }
1198 ptd.m_size = size();
1201 ptd.m_runtime_rdata = m_runtime_r_ptrs.dataPtr();
1202 ptd.m_runtime_idata = m_runtime_i_ptrs.dataPtr();
1203
1204#ifdef AMREX_USE_GPU
1205 if (copy_real || copy_int) {
1207 }
1208#endif
1209
1210 return ptd;
1211 }
1212
1214 {
1215 m_runtime_r_cptrs.resize(m_soa_tile.NumRealComps() - NArrayReal);
1216 m_runtime_i_cptrs.resize(m_soa_tile.NumIntComps() - NArrayInt);
1217#ifdef AMREX_USE_GPU
1218 bool copy_real = false;
1219 m_h_runtime_r_cptrs.resize(m_soa_tile.NumRealComps() - NArrayReal, nullptr);
1220 for (std::size_t i = 0; i < m_h_runtime_r_cptrs.size(); ++i) {
1221 if (m_h_runtime_r_cptrs[i] != m_soa_tile.GetRealData(i + NArrayReal).dataPtr()) {
1222 m_h_runtime_r_cptrs[i] = m_soa_tile.GetRealData(i + NArrayReal).dataPtr();
1223 copy_real = true;
1224 }
1225 }
1226 if (copy_real) {
1228 m_h_runtime_r_cptrs.size()*sizeof(ParticleReal*));
1229 }
1230
1231 bool copy_int = false;
1232 m_h_runtime_i_cptrs.resize(m_soa_tile.NumIntComps() - NArrayInt, nullptr);
1233 for (std::size_t i = 0; i < m_h_runtime_i_cptrs.size(); ++i) {
1234 if (m_h_runtime_i_cptrs[i] != m_soa_tile.GetIntData(i + NArrayInt).dataPtr()) {
1235 m_h_runtime_i_cptrs[i] = m_soa_tile.GetIntData(i + NArrayInt).dataPtr();
1236 copy_int = true;
1237 }
1238 }
1239 if (copy_int) {
1241 m_h_runtime_i_cptrs.size()*sizeof(int*));
1242 }
1243#else
1244 for (std::size_t i = 0; i < m_runtime_r_cptrs.size(); ++i) {
1245 m_runtime_r_cptrs[i] = m_soa_tile.GetRealData(i + NArrayReal).dataPtr();
1246 }
1247
1248 for (std::size_t i = 0; i < m_runtime_i_cptrs.size(); ++i) {
1249 m_runtime_i_cptrs[i] = m_soa_tile.GetIntData(i + NArrayInt).dataPtr();
1250 }
1251#endif
1252
1254 if constexpr (!ParticleType::is_soa_particle) {
1255 ptd.m_aos = m_aos_tile().dataPtr();
1256 } else {
1257 ptd.m_aos = nullptr;
1258 }
1259 if constexpr (ParticleType::is_soa_particle) {
1260 ptd.m_idcpu = m_soa_tile.GetIdCPUData().dataPtr();
1261 } else {
1262 ptd.m_idcpu = nullptr;
1263 }
1264 if constexpr(NArrayReal > 0) {
1265 for (int i = 0; i < NArrayReal; ++i) {
1266 ptd.m_rdata[i] = m_soa_tile.GetRealData(i).dataPtr();
1267 }
1268 }
1269 if constexpr(NArrayInt > 0) {
1270 for (int i = 0; i < NArrayInt; ++i) {
1271 ptd.m_idata[i] = m_soa_tile.GetIntData(i).dataPtr();
1272 }
1273 }
1274 ptd.m_size = size();
1277 ptd.m_runtime_rdata = m_runtime_r_cptrs.dataPtr();
1278 ptd.m_runtime_idata = m_runtime_i_cptrs.dataPtr();
1279
1280#ifdef AMREX_USE_GPU
1281 if (copy_real || copy_int) {
1283 }
1284#endif
1285
1286 return ptd;
1287 }
1288
1289private:
1290
1293
1294 bool m_defined = false;
1295
1298
1301
1304
1307};
1308
1309} // namespace amrex
1310
1311#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
Definition AMReX_GpuAllocators.H:114
Definition AMReX_ArrayOfStructs.H:13
Definition AMReX_PODVector.H:262
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
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:230
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:512
static constexpr int NAI
Definition AMReX_ParticleTile.H:503
GpuArray< const int *, NArrayInt > m_idata
Definition AMReX_ParticleTile.H:522
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE SuperParticleType getSuperParticle(int index) const noexcept
Definition AMReX_ParticleTile.H:647
const ParticleReal *AMREX_RESTRICT *AMREX_RESTRICT m_runtime_rdata
Definition AMReX_ParticleTile.H:588
const int *AMREX_RESTRICT *AMREX_RESTRICT m_runtime_idata
Definition AMReX_ParticleTile.H:589
Long m_size
Definition AMReX_ParticleTile.H:514
T_ParticleType ParticleType
Definition AMReX_ParticleTile.H:504
static constexpr int NStructReal
Definition AMReX_ParticleTile.H:507
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE const ParticleReal & pos(const int dir, const int index) const &
Definition AMReX_ParticleTile.H:525
int m_num_runtime_real
Definition AMReX_ParticleTile.H:586
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE const int * idata(const int attribute_index) const
Definition AMReX_ParticleTile.H:571
static constexpr int NStructInt
Definition AMReX_ParticleTile.H:508
int m_num_runtime_int
Definition AMReX_ParticleTile.H:587
std::conditional_t< T_ParticleType::is_soa_particle, void const *AMREX_RESTRICT, ParticleType const *AMREX_RESTRICT > AOS_PTR
Definition AMReX_ParticleTile.H:517
T_ParticleType const & ParticleRefType
Definition AMReX_ParticleTile.H:505
static constexpr int NAR
Definition AMReX_ParticleTile.H:502
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE decltype(auto) idcpu(const int index) const &
Definition AMReX_ParticleTile.H:555
AOS_PTR m_aos
Definition AMReX_ParticleTile.H:518
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:592
const uint64_t * m_idcpu
Definition AMReX_ParticleTile.H:520
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE decltype(auto) cpu(const int index) const &
Definition AMReX_ParticleTile.H:545
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE const ParticleReal * rdata(const int attribute_index) const
Definition AMReX_ParticleTile.H:565
GpuArray< const ParticleReal *, NArrayReal > m_rdata
Definition AMReX_ParticleTile.H:521
Definition AMReX_ParticleTile.H:315
static constexpr int NArrayReal
Definition AMReX_ParticleTile.H:316
ParticleReal RealType
Definition AMReX_ParticleTile.H:323
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ConstParticleCPUWrapper cpu() const
Definition AMReX_ParticleTile.H:336
static constexpr int NArrayInt
Definition AMReX_ParticleTile.H:317
ConstPTD m_constparticle_tile_data
Definition AMReX_ParticleTile.H:371
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE RealVect pos() const &
Definition AMReX_ParticleTile.H:344
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ConstParticleIDWrapper id() const
Definition AMReX_ParticleTile.H:339
static void NextID(Long nextid)
Reset on restart.
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE const RealType & pos(int position_index) const &
Definition AMReX_ParticleTile.H:347
int m_index
Definition AMReX_ParticleTile.H:372
static constexpr bool is_constsoa_particle
Definition AMReX_ParticleTile.H:321
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ConstSoAParticle(ConstPTD const &ptd, long i)
Definition AMReX_ParticleTile.H:326
static Long UnprotectedNextID()
This version can only be used inside omp critical.
Definition AMReX_Array.H:34
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:179
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void setSuperParticle(const SuperParticleType &sp, int index) const noexcept
Definition AMReX_ParticleTile.H:277
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 decltype(auto) cpu(const int index) const &
Definition AMReX_ParticleTile.H:83
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE SuperParticleType getSuperParticle(int index) const noexcept
Definition AMReX_ParticleTile.H:234
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE decltype(auto) idcpu(const int index) const &
Definition AMReX_ParticleTile.H:93
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_DEVICE AMREX_FORCE_INLINE int * idata(const int attribute_index) const
Definition AMReX_ParticleTile.H:109
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ParticleReal & pos(const int dir, const int index) const &
Definition AMReX_ParticleTile.H:63
AOS_PTR m_aos
Definition AMReX_ParticleTile.H:51
int *AMREX_RESTRICT *AMREX_RESTRICT m_runtime_idata
Definition AMReX_ParticleTile.H:60
static constexpr int NAI
Definition AMReX_ParticleTile.H:34
static constexpr int NStructReal
Definition AMReX_ParticleTile.H:40
Definition AMReX_ParticleTile.H:702
typename ParticleType::StorageParticleType StorageParticleType
Definition AMReX_ParticleTile.H:728
int NumRealComps() const noexcept
Definition AMReX_ParticleTile.H:1067
void push_back_real(int comp, amrex::Vector< amrex::ParticleReal > const &vec)
Definition AMReX_ParticleTile.H:1002
int NumIntComps() const noexcept
Definition AMReX_ParticleTile.H:1069
typename ParticleType::RealType RealType
Definition AMReX_ParticleTile.H:709
ParticleTile()=default
int getNumNeighbors() const
Definition AMReX_ParticleTile.H:901
amrex::PODVector< const int *, Allocator< const int * > > m_runtime_i_cptrs
Definition AMReX_ParticleTile.H:1300
void push_back(const ParticleType &p)
Definition AMReX_ParticleTile.H:923
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:725
const SoA & GetStructOfArrays() const
Definition AMReX_ParticleTile.H:824
static constexpr int NAI
Definition AMReX_ParticleTile.H:708
ParticleTile(ParticleTile &&) noexcept=default
void push_back_real(int comp, const ParticleReal *beg, const ParticleReal *end)
Definition AMReX_ParticleTile.H:985
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:745
static constexpr int NStructInt
Definition AMReX_ParticleTile.H:712
amrex::Gpu::HostVector< const ParticleReal * > m_h_runtime_r_cptrs
Definition AMReX_ParticleTile.H:1305
void setNumNeighbors(int num_neighbors)
Definition AMReX_ParticleTile.H:893
amrex::Gpu::HostVector< int * > m_h_runtime_i_ptrs
Definition AMReX_ParticleTile.H:1303
Long capacity() const
Definition AMReX_ParticleTile.H:1095
int numTotalParticles() const
Returns the total number of particles, real and neighbor.
Definition AMReX_ParticleTile.H:884
AoS m_aos_tile
Definition AMReX_ParticleTile.H:1291
std::size_t size() const
Returns the total number of particles (real and neighbor)
Definition AMReX_ParticleTile.H:832
ParticleTileDataType getParticleTileData()
Definition AMReX_ParticleTile.H:1137
void push_back(const SuperParticleType &sp)
Definition AMReX_ParticleTile.H:930
int numNeighborParticles() const
Returns the number of neighbor particles (excluding reals)
Definition AMReX_ParticleTile.H:871
~ParticleTile()=default
amrex::PODVector< ParticleReal *, Allocator< ParticleReal * > > m_runtime_r_ptrs
Definition AMReX_ParticleTile.H:1296
const AoS & GetArrayOfStructs() const
Definition AMReX_ParticleTile.H:821
static constexpr int NAR
Definition AMReX_ParticleTile.H:707
ConstParticleTileDataType getConstParticleTileData() const
Definition AMReX_ParticleTile.H:1213
void shrink_to_fit()
Definition AMReX_ParticleTile.H:1075
T_ParticleType ParticleType
Definition AMReX_ParticleTile.H:706
int numParticles() const
Returns the number of real particles (excluding neighbors)
Definition AMReX_ParticleTile.H:845
void push_back_int(int comp, amrex::Vector< int > const &vec)
Definition AMReX_ParticleTile.H:1054
amrex::Gpu::HostVector< ParticleReal * > m_h_runtime_r_ptrs
Definition AMReX_ParticleTile.H:1302
void swap(ParticleTile< ParticleType, NArrayReal, NArrayInt, Allocator > &other) noexcept
Definition AMReX_ParticleTile.H:1117
void push_back_int(int comp, amrex::Vector< int >::const_iterator beg, amrex::Vector< int >::const_iterator end)
Definition AMReX_ParticleTile.H:1046
void push_back_real(const std::array< ParticleReal, NArrayReal > &v)
Definition AMReX_ParticleTile.H:975
int NumRuntimeRealComps() const noexcept
Definition AMReX_ParticleTile.H:1071
amrex::Gpu::HostVector< const int * > m_h_runtime_i_cptrs
Definition AMReX_ParticleTile.H:1306
amrex::PODVector< const ParticleReal *, Allocator< const ParticleReal * > > m_runtime_r_cptrs
Definition AMReX_ParticleTile.H:1299
AoS & GetArrayOfStructs()
Definition AMReX_ParticleTile.H:820
RealType & pos(int index, int position_index) &
Definition AMReX_ParticleTile.H:797
void resize(std::size_t count)
Definition AMReX_ParticleTile.H:911
typename SoA::IntVector IntVector
Definition AMReX_ParticleTile.H:727
void push_back_int(const std::array< int, NArrayInt > &v)
Definition AMReX_ParticleTile.H:1027
void push_back_real(int comp, amrex::Vector< amrex::ParticleReal >::const_iterator beg, amrex::Vector< amrex::ParticleReal >::const_iterator end)
Definition AMReX_ParticleTile.H:994
Allocator< T > AllocatorType
Definition AMReX_ParticleTile.H:704
int NumRuntimeIntComps() const noexcept
Definition AMReX_ParticleTile.H:1073
bool empty() const
Definition AMReX_ParticleTile.H:826
void push_back_int(int comp, int v)
Definition AMReX_ParticleTile.H:1019
RealType pos(int index, int position_index) const &
Definition AMReX_ParticleTile.H:811
decltype(auto) cpu(int index) &
Definition AMReX_ParticleTile.H:779
SoA m_soa_tile
Definition AMReX_ParticleTile.H:1292
SoA & GetStructOfArrays()
Definition AMReX_ParticleTile.H:823
void push_back_int(int comp, std::size_t npar, int v)
Definition AMReX_ParticleTile.H:1062
bool m_defined
Definition AMReX_ParticleTile.H:1294
void push_back_real(int comp, std::size_t npar, ParticleReal v)
Definition AMReX_ParticleTile.H:1010
void push_back_int(int comp, const int *beg, const int *end)
Definition AMReX_ParticleTile.H:1037
typename SoA::RealVector RealVector
Definition AMReX_ParticleTile.H:726
amrex::PODVector< int *, Allocator< int * > > m_runtime_i_ptrs
Definition AMReX_ParticleTile.H:1297
int numRealParticles() const
Returns the number of real particles (excluding neighbors)
Definition AMReX_ParticleTile.H:858
std::conditional_t< ParticleType::is_soa_particle, ThisParticleTileHasNoAoS, ArrayOfStructs< ParticleType, Allocator > > AoS
Definition AMReX_ParticleTile.H:719
static constexpr int NStructReal
Definition AMReX_ParticleTile.H:711
void push_back_real(int comp, ParticleReal v)
Definition AMReX_ParticleTile.H:967
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 ParticleCPUWrapper cpu() &
Definition AMReX_Particle.H:312
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE int & idata(int index) &
Definition AMReX_Particle.H:427
Definition AMReX_Particle.H:281
Definition AMReX_ParticleTile.H:377
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE const uint64_t & idcpu() const &
Definition AMReX_ParticleTile.H:414
static constexpr int NArrayReal
Definition AMReX_ParticleTile.H:378
static Long UnprotectedNextID()
This version can only be used inside omp critical.
Definition AMReX_ParticleTile.H:483
PTD m_particle_tile_data
Definition AMReX_ParticleTile.H:453
static Long the_next_id
Definition AMReX_ParticleTile.H:394
static constexpr int NArrayInt
Definition AMReX_ParticleTile.H:379
int m_index
Definition AMReX_ParticleTile.H:454
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ParticleCPUWrapper cpu() &
Definition AMReX_ParticleTile.H:399
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ConstParticleCPUWrapper cpu() const &
Definition AMReX_ParticleTile.H:408
ParticleReal RealType
Definition AMReX_ParticleTile.H:386
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE RealType pos(int position_index) const &
Definition AMReX_ParticleTile.H:429
static constexpr bool is_constsoa_particle
Definition AMReX_ParticleTile.H:383
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ParticleIDWrapper id() &
Definition AMReX_ParticleTile.H:402
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE SoAParticle(PTD const &ptd, long i)
Definition AMReX_ParticleTile.H:389
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE RealVect pos() const &
Definition AMReX_ParticleTile.H:419
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE uint64_t & idcpu() &
Definition AMReX_ParticleTile.H:405
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ConstParticleIDWrapper id() const &
Definition AMReX_ParticleTile.H:411
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE RealType & pos(int position_index) &
Definition AMReX_ParticleTile.H:422
static Long NextID()
Definition AMReX_ParticleTile.H:462
Definition AMReX_StructOfArrays.H:20
Definition AMReX_ParticleTile.H:695
Definition AMReX_ParticleTile.H:693
Definition AMReX_MakeParticle.H:13