Block-Structured AMR Software Framework
Loading...
Searching...
No Matches
AMReX_ParticleReduce.H
Go to the documentation of this file.
1#ifndef AMREX_PARTICLEREDUCE_H_
2#define AMREX_PARTICLEREDUCE_H_
3#include <AMReX_Config.H>
4
5#include <AMReX_IntVect.H>
6#include <AMReX_Box.H>
7#include <AMReX_Gpu.H>
8#include <AMReX_Print.H>
9#include <AMReX_GpuUtility.H>
10#include <AMReX_TypeTraits.H>
11#include <AMReX_ParticleUtil.H>
12#include <AMReX_Vector.H>
13
14#include <iterator>
15#include <limits>
16#include <type_traits>
17
18namespace amrex {
19
21namespace particle_detail {
22
23template <typename F, typename T_ParticleType, int NAR, int NAI>
24AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
25auto call_f (F const& f,
26 const ConstParticleTileData<T_ParticleType, NAR, NAI>& p,
27 const int i) noexcept
28{
29 if constexpr ( ! T_ParticleType::is_soa_particle &&
30 IsCallable<F, T_ParticleType const&>::value) {
31 return f(p.m_aos[i]);
32 } else if constexpr (IsCallable<F, decltype(p.getSuperParticle(i))>::value) {
33 return f(p.getSuperParticle(i));
34 } else {
35 return f(p, i);
36 }
37}
38
39template <typename F, template<class, class> class PTDType, class RType, class IType>
40AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
41auto call_f (F const& f,
42 const PTDType<RType, IType>& p,
43 const int i) noexcept
44{
45 return f(p, i);
46}
47}
49
90template <class PC, class F>
91requires (IsParticleContainer<PC>::value)
92auto
93ReduceSum (PC const& pc, F&& f)
94 -> decltype(particle_detail::call_f(f, typename PC::ConstPTDType(), int()))
95{
96 return ReduceSum(pc, 0, pc.finestLevel(), std::forward<F>(f));
97}
98
139template <class PC, class F>
140requires (IsParticleContainer<PC>::value)
141auto
142ReduceSum (PC const& pc, int lev, F&& f)
143 -> decltype(particle_detail::call_f(f, typename PC::ConstPTDType(), int()))
144{
145 return ReduceSum(pc, lev, lev, std::forward<F>(f));
146}
147
189template <class PC, class F>
190requires (IsParticleContainer<PC>::value)
191auto
192ReduceSum (PC const& pc, int lev_min, int lev_max, F const& f)
193 -> decltype(particle_detail::call_f(f, typename PC::ConstPTDType(), int()))
194{
195 using value_type = decltype(particle_detail::call_f(f, typename PC::ConstPTDType(), int()));
196 value_type sm = 0;
197
198#ifdef AMREX_USE_GPU
199 if (Gpu::inLaunchRegion())
200 {
201 ReduceOps<ReduceOpSum> reduce_op;
202 ReduceData<value_type> reduce_data(reduce_op);
203 using ReduceTuple = typename decltype(reduce_data)::Type;
204
205 for (int lev = lev_min; lev <= lev_max; ++lev)
206 {
207 const auto& plev = pc.GetParticles(lev);
208 for (const auto& kv : plev)
209 {
210 const auto& tile = plev.at(kv.first);
211 const auto np = tile.numParticles();
212 const auto& ptd = tile.getConstParticleTileData();
213 reduce_op.eval(np, reduce_data,
214 [=] AMREX_GPU_DEVICE (const int i) -> ReduceTuple {
215 return particle_detail::call_f(f, ptd, i);
216 });
217 }
218 }
219
220 ReduceTuple hv = reduce_data.value(reduce_op);
221 sm = amrex::get<0>(hv);
222 }
223 else
224#endif
225 {
226 for (int lev = lev_min; lev <= lev_max; ++lev)
227 {
228 const auto& plev = pc.GetParticles(lev);
229 Vector<std::pair<int, int> > grid_tile_ids;
230 Vector<const typename PC::ParticleTileType*> ptile_ptrs;
231 for (auto& kv : plev)
232 {
233 grid_tile_ids.push_back(kv.first);
234 ptile_ptrs.push_back(&(kv.second));
235 }
236
237 // We avoid ParIter here to allow reductions after grids change but prior to calling Redistribute()
238#ifdef AMREX_USE_OMP
239#pragma omp parallel for if (!system::regtest_reduction) reduction(+:sm)
240#endif
241 for (int pmap_it = 0; pmap_it < std::ssize(ptile_ptrs); ++pmap_it)
242 {
243 const auto& tile = plev.at(grid_tile_ids[pmap_it]);
244 const auto np = tile.numParticles();
245 const auto& ptd = tile.getConstParticleTileData();
246 for (int i = 0; i < np; ++i) {
247 sm += particle_detail::call_f(f, ptd, i);
248 }
249 }
250 }
251 }
252
253 return sm;
254}
255
295template <class PC, class F>
296requires (IsParticleContainer<PC>::value)
297auto
298ReduceMax (PC const& pc, F&& f)
299 -> decltype(particle_detail::call_f(f, typename PC::ConstPTDType(), int()))
300{
301 return ReduceMax(pc, 0, pc.finestLevel(), std::forward<F>(f));
302}
303
345template <class PC, class F>
346requires (IsParticleContainer<PC>::value)
347auto
348ReduceMax (PC const& pc, int lev, F&& f)
349 -> decltype(particle_detail::call_f(f, typename PC::ConstPTDType(), int()))
350{
351 return ReduceMax(pc, lev, lev, std::forward<F>(f));
352}
353
395template <class PC, class F>
396requires (IsParticleContainer<PC>::value)
397auto
398ReduceMax (PC const& pc, int lev_min, int lev_max, F const& f)
399 -> decltype(particle_detail::call_f(f, typename PC::ConstPTDType(), int()))
400{
401 using value_type = decltype(particle_detail::call_f(f, typename PC::ConstPTDType(), int()));
402 constexpr value_type value_lowest = std::numeric_limits<value_type>::lowest();
403 value_type r = value_lowest;
404
405#ifdef AMREX_USE_GPU
406 if (Gpu::inLaunchRegion())
407 {
408 ReduceOps<ReduceOpMax> reduce_op;
409 ReduceData<value_type> reduce_data(reduce_op);
410 using ReduceTuple = typename decltype(reduce_data)::Type;
411
412 for (int lev = lev_min; lev <= lev_max; ++lev)
413 {
414 const auto& plev = pc.GetParticles(lev);
415 for (const auto& kv : plev)
416 {
417 const auto& tile = plev.at(kv.first);
418 const auto np = tile.numParticles();
419 const auto& ptd = tile.getConstParticleTileData();
420 reduce_op.eval(np, reduce_data,
421 [=] AMREX_GPU_DEVICE (const int i) -> ReduceTuple {
422 return particle_detail::call_f(f, ptd, i);
423 });
424 }
425 }
426
427 ReduceTuple hv = reduce_data.value(reduce_op);
428 r = amrex::get<0>(hv);
429 }
430 else
431#endif
432 {
433 for (int lev = lev_min; lev <= lev_max; ++lev)
434 {
435 const auto& plev = pc.GetParticles(lev);
436 Vector<std::pair<int, int> > grid_tile_ids;
437 Vector<const typename PC::ParticleTileType*> ptile_ptrs;
438 for (auto& kv : plev)
439 {
440 grid_tile_ids.push_back(kv.first);
441 ptile_ptrs.push_back(&(kv.second));
442 }
443
444 // We avoid ParIter here to allow reductions after grids change but prior to calling Redistribute()
445#ifdef AMREX_USE_OMP
446#pragma omp parallel for if (!system::regtest_reduction) reduction(max:r)
447#endif
448 for (int pmap_it = 0; pmap_it < std::ssize(ptile_ptrs); ++pmap_it)
449 {
450 const auto& tile = plev.at(grid_tile_ids[pmap_it]);
451 const auto np = tile.numParticles();
452 const auto& ptd = tile.getConstParticleTileData();
453 for (int i = 0; i < np; ++i) {
454 r = std::max(r, particle_detail::call_f(f, ptd, i));
455 }
456 }
457 }
458 }
459
460 return r;
461}
462
502template <class PC, class F>
503requires (IsParticleContainer<PC>::value)
504auto
505ReduceMin (PC const& pc, F&& f)
506 -> decltype(particle_detail::call_f(f, typename PC::ConstPTDType(), int()))
507{
508 return ReduceMin(pc, 0, pc.finestLevel(), std::forward<F>(f));
509}
510
551template <class PC, class F>
552requires (IsParticleContainer<PC>::value)
553auto
554ReduceMin (PC const& pc, int lev, F&& f)
555 -> decltype(particle_detail::call_f(f, typename PC::ConstPTDType(), int()))
556{
557 return ReduceMin(pc, lev, lev, std::forward<F>(f));
558}
559
601template <class PC, class F>
602requires (IsParticleContainer<PC>::value)
603auto
604ReduceMin (PC const& pc, int lev_min, int lev_max, F const& f)
605 -> decltype(particle_detail::call_f(f, typename PC::ConstPTDType(), int()))
606{
607 using value_type = decltype(particle_detail::call_f(f, typename PC::ConstPTDType(), int()));
608 constexpr value_type value_max = std::numeric_limits<value_type>::max();
609 value_type r = value_max;
610
611#ifdef AMREX_USE_GPU
612 if (Gpu::inLaunchRegion())
613 {
614 ReduceOps<ReduceOpMin> reduce_op;
615 ReduceData<value_type> reduce_data(reduce_op);
616 using ReduceTuple = typename decltype(reduce_data)::Type;
617
618 for (int lev = lev_min; lev <= lev_max; ++lev)
619 {
620 const auto& plev = pc.GetParticles(lev);
621 for (const auto& kv : plev)
622 {
623 const auto& tile = plev.at(kv.first);
624 const auto np = tile.numParticles();
625 const auto& ptd = tile.getConstParticleTileData();
626 reduce_op.eval(np, reduce_data,
627 [=] AMREX_GPU_DEVICE (const int i) -> ReduceTuple {
628 return particle_detail::call_f(f, ptd, i);
629 });
630 }
631 }
632
633 ReduceTuple hv = reduce_data.value(reduce_op);
634 r = amrex::get<0>(hv);
635 }
636 else
637#endif
638 {
639 for (int lev = lev_min; lev <= lev_max; ++lev)
640 {
641 const auto& plev = pc.GetParticles(lev);
642 Vector<std::pair<int, int> > grid_tile_ids;
643 Vector<const typename PC::ParticleTileType*> ptile_ptrs;
644 for (auto& kv : plev)
645 {
646 grid_tile_ids.push_back(kv.first);
647 ptile_ptrs.push_back(&(kv.second));
648 }
649
650 // We avoid ParIter here to allow reductions after grids change but prior to calling Redistribute()
651#ifdef AMREX_USE_OMP
652#pragma omp parallel for if (!system::regtest_reduction) reduction(min:r)
653#endif
654 for (int pmap_it = 0; pmap_it < std::ssize(ptile_ptrs); ++pmap_it)
655 {
656 const auto& tile = plev.at(grid_tile_ids[pmap_it]);
657 const auto np = tile.numParticles();
658 const auto& ptd = tile.getConstParticleTileData();
659 for (int i = 0; i < np; ++i) {
660 r = std::min(r, particle_detail::call_f(f, ptd, i));
661 }
662 }
663 }
664 }
665
666 return r;
667}
668
708template <class PC, class F>
709requires (IsParticleContainer<PC>::value)
710bool
711ReduceLogicalAnd (PC const& pc, F&& f)
712{
713 return ReduceLogicalAnd(pc, 0, pc.finestLevel(), std::forward<F>(f));
714}
715
756template <class PC, class F>
757requires (IsParticleContainer<PC>::value)
758bool
759ReduceLogicalAnd (PC const& pc, int lev, F&& f)
760{
761 return ReduceLogicalAnd(pc, lev, lev, std::forward<F>(f));
762}
763
805template <class PC, class F>
806requires (IsParticleContainer<PC>::value)
807bool
808ReduceLogicalAnd (PC const& pc, int lev_min, int lev_max, F const& f)
809{
810 int r = true;
811
812#ifdef AMREX_USE_GPU
813 if (Gpu::inLaunchRegion())
814 {
815 ReduceOps<ReduceOpLogicalAnd> reduce_op;
816 ReduceData<int> reduce_data(reduce_op);
817 using ReduceTuple = typename decltype(reduce_data)::Type;
818
819 for (int lev = lev_min; lev <= lev_max; ++lev)
820 {
821 const auto& plev = pc.GetParticles(lev);
822 for (const auto& kv : plev)
823 {
824 const auto& tile = plev.at(kv.first);
825 const auto np = tile.numParticles();
826 const auto& ptd = tile.getConstParticleTileData();
827 reduce_op.eval(np, reduce_data,
828 [=] AMREX_GPU_DEVICE (const int i) -> ReduceTuple {
829 return particle_detail::call_f(f, ptd, i);
830 });
831 }
832 }
833
834 ReduceTuple hv = reduce_data.value(reduce_op);
835 r = amrex::get<0>(hv);
836 }
837 else
838#endif
839 {
840 for (int lev = lev_min; lev <= lev_max; ++lev)
841 {
842 const auto& plev = pc.GetParticles(lev);
843 Vector<std::pair<int, int> > grid_tile_ids;
844 Vector<const typename PC::ParticleTileType*> ptile_ptrs;
845 for (auto& kv : plev)
846 {
847 grid_tile_ids.push_back(kv.first);
848 ptile_ptrs.push_back(&(kv.second));
849 }
850
851 // We avoid ParIter here to allow reductions after grids change but prior to calling Redistribute()
852#ifdef AMREX_USE_OMP
853#pragma omp parallel for if (!system::regtest_reduction) reduction(&&:r)
854#endif
855 for (int pmap_it = 0; pmap_it < std::ssize(ptile_ptrs); ++pmap_it)
856 {
857 const auto& tile = plev.at(grid_tile_ids[pmap_it]);
858 const auto np = tile.numParticles();
859 const auto& ptd = tile.getConstParticleTileData();
860 for (int i = 0; i < np; ++i) {
861 r = r && particle_detail::call_f(f, ptd, i);
862 }
863 }
864 }
865 }
866
867 return r;
868}
869
909template <class PC, class F>
910requires (IsParticleContainer<PC>::value)
911bool
912ReduceLogicalOr (PC const& pc, F&& f)
913{
914 return ReduceLogicalOr(pc, 0, pc.finestLevel(), std::forward<F>(f));
915}
916
957template <class PC, class F>
958requires (IsParticleContainer<PC>::value)
959bool
960ReduceLogicalOr (PC const& pc, int lev, F&& f)
961{
962 return ReduceLogicalOr(pc, lev, lev, std::forward<F>(f));
963}
964
1006template <class PC, class F>
1007requires (IsParticleContainer<PC>::value)
1008bool
1009ReduceLogicalOr (PC const& pc, int lev_min, int lev_max, F const& f)
1010{
1011 int r = false;
1012
1013#ifdef AMREX_USE_GPU
1014 if (Gpu::inLaunchRegion())
1015 {
1016 ReduceOps<ReduceOpLogicalOr> reduce_op;
1017 ReduceData<int> reduce_data(reduce_op);
1018 using ReduceTuple = typename decltype(reduce_data)::Type;
1019
1020 for (int lev = lev_min; lev <= lev_max; ++lev)
1021 {
1022 const auto& plev = pc.GetParticles(lev);
1023 for (const auto& kv : plev)
1024 {
1025 const auto& tile = plev.at(kv.first);
1026 const auto np = tile.numParticles();
1027 const auto& ptd = tile.getConstParticleTileData();
1028 reduce_op.eval(np, reduce_data,
1029 [=] AMREX_GPU_DEVICE (const int i) -> ReduceTuple
1030 {
1031 return particle_detail::call_f(f, ptd, i);
1032 });
1033 }
1034 }
1035
1036 ReduceTuple hv = reduce_data.value(reduce_op);
1037 r = amrex::get<0>(hv);
1038 }
1039 else
1040#endif
1041 {
1042 for (int lev = lev_min; lev <= lev_max; ++lev)
1043 {
1044 const auto& plev = pc.GetParticles(lev);
1045 Vector<std::pair<int, int> > grid_tile_ids;
1046 Vector<const typename PC::ParticleTileType*> ptile_ptrs;
1047 for (auto& kv : plev)
1048 {
1049 grid_tile_ids.push_back(kv.first);
1050 ptile_ptrs.push_back(&(kv.second));
1051 }
1052
1053 // We avoid ParIter here to allow reductions after grids change but prior to calling Redistribute()
1054#ifdef AMREX_USE_OMP
1055#pragma omp parallel for if (!system::regtest_reduction) reduction(||:r)
1056#endif
1057 for (int pmap_it = 0; pmap_it < std::ssize(ptile_ptrs); ++pmap_it)
1058 {
1059 const auto& tile = plev.at(grid_tile_ids[pmap_it]);
1060 const auto np = tile.numParticles();
1061 const auto& ptd = tile.getConstParticleTileData();
1062 for (int i = 0; i < np; ++i) {
1063 r = r || particle_detail::call_f(f, ptd, i);
1064 }
1065 }
1066 }
1067 }
1068
1069 return r;
1070}
1071
1134template <class RD, class PC, class F, class ReduceOps>
1135requires (IsParticleContainer<PC>::value)
1136typename RD::Type
1137ParticleReduce (PC const& pc, F&& f, ReduceOps& reduce_ops)
1138{
1139 return ParticleReduce<RD>(pc, 0, pc.finestLevel(), std::forward<F>(f), reduce_ops);
1140}
1141
1205template <class RD, class PC, class F, class ReduceOps>
1206requires (IsParticleContainer<PC>::value)
1207typename RD::Type
1208ParticleReduce (PC const& pc, int lev, F&& f, ReduceOps& reduce_ops)
1209{
1210 return ParticleReduce<RD>(pc, lev, lev, std::forward<F>(f), reduce_ops);
1211}
1212
1277template <class RD, class PC, class F, class ReduceOps>
1278requires (IsParticleContainer<PC>::value)
1279typename RD::Type
1280ParticleReduce (PC const& pc, int lev_min, int lev_max, F const& f, ReduceOps& reduce_ops)
1281{
1282 RD reduce_data(reduce_ops);
1283 for (int lev = lev_min; lev <= lev_max; ++lev) {
1284 const auto& plev = pc.GetParticles(lev);
1285 Vector<std::pair<int, int> > grid_tile_ids;
1286 Vector<const typename PC::ParticleTileType*> ptile_ptrs;
1287 for (auto& kv : plev)
1288 {
1289 grid_tile_ids.push_back(kv.first);
1290 ptile_ptrs.push_back(&(kv.second));
1291 }
1292
1293 // We avoid ParIter here to allow reductions after grids change but prior to calling Redistribute()
1294#if !defined(AMREX_USE_GPU) && defined(AMREX_USE_OMP)
1295#pragma omp parallel for
1296#endif
1297 for (int pmap_it = 0; pmap_it < std::ssize(ptile_ptrs); ++pmap_it)
1298 {
1299 const auto& tile = plev.at(grid_tile_ids[pmap_it]);
1300 const auto np = tile.numParticles();
1301 const auto& ptd = tile.getConstParticleTileData();
1302 reduce_ops.eval(np, reduce_data,
1303 [=] AMREX_GPU_DEVICE (const int i) noexcept
1304 {
1305 return particle_detail::call_f(f, ptd, i);
1306 });
1307 }
1308 }
1309 return reduce_data.value(reduce_ops);
1310}
1311}
1312#endif
AMReX_Box.H
AMREX_FORCE_INLINE
#define AMREX_FORCE_INLINE
Definition AMReX_Extension.H:119
AMREX_GPU_DEVICE
#define AMREX_GPU_DEVICE
Definition AMReX_GpuQualifiers.H:18
AMREX_GPU_HOST_DEVICE
#define AMREX_GPU_HOST_DEVICE
Definition AMReX_GpuQualifiers.H:20
AMReX_GpuUtility.H
AMReX_Gpu.H
AMReX_IntVect.H
AMReX_ParticleUtil.H
AMReX_Print.H
AMReX_TypeTraits.H
AMReX_Vector.H
amrex::ReduceData
Definition AMReX_Reduce.H:438
amrex::ReduceData::value
Type value()
Definition AMReX_Reduce.H:473
amrex::ReduceOps
Definition AMReX_Reduce.H:597
amrex::ReduceOps::eval
void eval(MF const &mf, IntVect const &nghost, D &reduce_data, F &&f)
Definition AMReX_Reduce.H:731
amrex::Vector
This class is a thin wrapper around std::vector. Unlike vector, Vector::operator[] provides bound che...
Definition AMReX_Vector.H:29
amrex::Gpu::inLaunchRegion
bool inLaunchRegion() noexcept
Definition AMReX_GpuControl.H:88
amrex
Definition AMReX_Amr.cpp:50
amrex::ReduceMin
FAB::value_type ReduceMin(FabArray< FAB > const &fa, int nghost, F &&f)
Definition AMReX_FabArrayUtility.H:304
amrex::ParticleReduce
RD::Type ParticleReduce(PC const &pc, F &&f, ReduceOps &reduce_ops)
A general reduction method for the particles in a ParticleContainer that can run on either CPUs or GP...
Definition AMReX_ParticleReduce.H:1137
amrex::ReduceLogicalOr
bool ReduceLogicalOr(FabArray< FAB > const &fa, int nghost, F &&f)
Definition AMReX_FabArrayUtility.H:904
amrex::ReduceLogicalAnd
bool ReduceLogicalAnd(FabArray< FAB > const &fa, int nghost, F &&f)
Definition AMReX_FabArrayUtility.H:757
amrex::ReduceMax
FAB::value_type ReduceMax(FabArray< FAB > const &fa, int nghost, F &&f)
Definition AMReX_FabArrayUtility.H:530
amrex::ReduceSum
FAB::value_type ReduceSum(FabArray< FAB > const &fa, int nghost, F &&f)
Definition AMReX_FabArrayUtility.H:16