Block-Structured AMR Software Framework
Loading...
Searching...
No Matches
AMReX_MultiFabUtil.H
Go to the documentation of this file.
1#ifndef AMREX_MultiFabUtil_H_
2#define AMREX_MultiFabUtil_H_
3#include <AMReX_Config.H>
4
5#include <AMReX_MultiFab.H>
6#include <AMReX_iMultiFab.H>
7#include <AMReX_LayoutData.H>
8#include <AMReX_MFIter.H>
9#include <AMReX_Array.H>
10#include <AMReX_Vector.H>
11#include <AMReX_MultiFabUtil_C.H>
12
13#include <AMReX_MultiFabUtilI.H>
14
15namespace amrex
16{
18 void average_node_to_cellcenter (MultiFab& cc, int dcomp,
19 const MultiFab& nd, int scomp,
20 int ncomp, int ngrow = 0);
21 void average_node_to_cellcenter (MultiFab& cc, int dcomp,
22 const MultiFab& nd, int scomp,
23 int ncomp, IntVect const& ng_vect);
24
31 void average_edge_to_cellcenter (MultiFab& cc, int dcomp,
32 const Vector<const MultiFab*>& edge,
33 int ngrow = 0);
34 void average_edge_to_cellcenter (MultiFab& cc, int dcomp,
35 const Vector<const MultiFab*>& edge,
36 IntVect const& ng_vect);
37
44 void average_face_to_cellcenter (MultiFab& cc, int dcomp,
45 const Vector<const MultiFab*>& fc,
46 int ngrow = 0);
47 void average_face_to_cellcenter (MultiFab& cc, int dcomp,
48 const Vector<const MultiFab*>& fc,
49 IntVect const& ng_vect);
50
52 template <typename CMF, typename FMF,
53 std::enable_if_t<IsFabArray_v<CMF> && IsFabArray_v<FMF>, int> = 0>
54 void average_face_to_cellcenter (CMF& cc, int dcomp,
55 const Array<const FMF*,AMREX_SPACEDIM>& fc,
56 int ngrow = 0);
57
58 template <typename CMF, typename FMF,
59 std::enable_if_t<IsFabArray_v<CMF> && IsFabArray_v<FMF>, int> = 0>
60 void average_face_to_cellcenter (CMF& cc, int dcomp,
61 const Array<const FMF*,AMREX_SPACEDIM>& fc,
62 IntVect const& ng_vect);
63
65 void average_face_to_cellcenter (MultiFab& cc,
66 const Vector<const MultiFab*>& fc,
67 const Geometry& geom);
69 void average_face_to_cellcenter (MultiFab& cc,
70 const Array<const MultiFab*,AMREX_SPACEDIM>& fc,
71 const Geometry& geom);
73 void average_cellcenter_to_face (const Vector<MultiFab*>& fc,
74 const MultiFab& cc,
75 const Geometry& geom,
76 int ncomp = 1,
77 bool use_harmonic_averaging = false);
79 void average_cellcenter_to_face (const Array<MultiFab*,AMREX_SPACEDIM>& fc,
80 const MultiFab& cc,
81 const Geometry& geom,
82 int ncomp = 1,
83 bool use_harmonic_averaging = false);
84
86 template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int> = 0>
87 void average_down_faces (const Vector<const MF*>& fine,
88 const Vector<MF*>& crse,
89 const IntVect& ratio,
90 int ngcrse = 0);
92 template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int> = 0>
93 void average_down_faces (const Vector<const MF*>& fine,
94 const Vector<MF*>& crse,
95 int ratio,
96 int ngcrse = 0);
98 template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int> = 0>
99 void average_down_faces (const Array<const MF*,AMREX_SPACEDIM>& fine,
100 const Array<MF*,AMREX_SPACEDIM>& crse,
101 const IntVect& ratio,
102 int ngcrse = 0);
104 template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int> = 0>
105 void average_down_faces (const Array<const MF*,AMREX_SPACEDIM>& fine,
106 const Array<MF*,AMREX_SPACEDIM>& crse,
107 int ratio,
108 int ngcrse = 0);
115 template <typename FAB>
116 void average_down_faces (const FabArray<FAB>& fine, FabArray<FAB>& crse,
117 const IntVect& ratio, int ngcrse=0);
118
119 // This version takes periodicity into account.
120 template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int> = 0>
121 void average_down_faces (const Array<const MF*,AMREX_SPACEDIM>& fine,
122 const Array<MF*,AMREX_SPACEDIM>& crse,
123 const IntVect& ratio, const Geometry& crse_geom);
124 // This version takes periodicity into account.
125 template <typename FAB>
126 void average_down_faces (const FabArray<FAB>& fine, FabArray<FAB>& crse,
127 const IntVect& ratio, const Geometry& crse_geom);
128
130 void average_down_edges (const Vector<const MultiFab*>& fine,
131 const Vector<MultiFab*>& crse,
132 const IntVect& ratio,
133 int ngcrse = 0);
134 void average_down_edges (const Array<const MultiFab*,AMREX_SPACEDIM>& fine,
135 const Array<MultiFab*,AMREX_SPACEDIM>& crse,
136 const IntVect& ratio,
137 int ngcrse = 0);
141 void average_down_edges (const MultiFab& fine, MultiFab& crse,
142 const IntVect& ratio, int ngcrse=0);
143
145 template <typename FAB>
146 void average_down_nodal (const FabArray<FAB>& S_fine,
147 FabArray<FAB>& S_crse,
148 const IntVect& ratio,
149 int ngcrse = 0,
150 bool mfiter_is_definitely_safe=false);
151
158 void average_down (const MultiFab& S_fine, MultiFab& S_crse,
159 const Geometry& fgeom, const Geometry& cgeom,
160 int scomp, int ncomp, const IntVect& ratio);
161 void average_down (const MultiFab& S_fine, MultiFab& S_crse,
162 const Geometry& fgeom, const Geometry& cgeom,
163 int scomp, int ncomp, int rr);
164
168 template<typename FAB>
169 void average_down (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse,
170 int scomp, int ncomp, const IntVect& ratio);
171 template<typename FAB>
172 void average_down (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse,
173 int scomp, int ncomp, int rr);
174
177 void sum_fine_to_coarse (const MultiFab& S_Fine, MultiFab& S_crse,
178 int scomp, int ncomp,
179 const IntVect& ratio,
180 const Geometry& cgeom, const Geometry& fgeom);
181
183 void print_state (const MultiFab& mf, const IntVect& cell, int n=-1,
184 const IntVect& ng = IntVect::TheZeroVector());
185
187 void writeFabs (const MultiFab& mf, const std::string& name);
188 void writeFabs (const MultiFab& mf, int comp, int ncomp, const std::string& name);
189
192 std::unique_ptr<MultiFab> get_slice_data(int dir, Real coord,
193 const MultiFab& cc,
194 const Geometry& geom, int start_comp, int ncomp,
195 bool interpolate=false,
196 RealBox const& bnd_rbx = RealBox());
197
205 template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int> FOO = 0>
206 Vector<typename MF::value_type> get_cell_data (MF const& mf, IntVect const& cell);
207
214 template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int> FOO = 0>
215 MF get_line_data (MF const& mf, int dir, IntVect const& cell, Box const& bnd_bx = Box());
216
220 template <typename FAB>
221 iMultiFab makeFineMask (const FabArray<FAB>& cmf, const BoxArray& fba, const IntVect& ratio,
222 int crse_value = 0, int fine_value = 1,
223 MFInfo const& info = MFInfo());
224 iMultiFab makeFineMask (const BoxArray& cba, const DistributionMapping& cdm,
225 const BoxArray& fba, const IntVect& ratio,
226 int crse_value = 0, int fine_value = 1,
227 MFInfo const& info = MFInfo());
228 template <typename FAB>
229 iMultiFab makeFineMask (const FabArray<FAB>& cmf, const BoxArray& fba, const IntVect& ratio,
230 Periodicity const& period, int crse_value, int fine_value,
231 MFInfo const& info = MFInfo());
232 iMultiFab makeFineMask (const BoxArray& cba, const DistributionMapping& cdm,
233 const IntVect& cnghost, const BoxArray& fba, const IntVect& ratio,
234 Periodicity const& period, int crse_value, int fine_value,
235 MFInfo const& info = MFInfo());
236 template <typename FAB>
237 iMultiFab makeFineMask (const FabArray<FAB>& cmf, const FabArray<FAB>& fmf,
238 const IntVect& cnghost, const IntVect& ratio,
239 Periodicity const& period, int crse_value, int fine_value,
240 MFInfo const& info = MFInfo());
241 template <typename FAB>
242 iMultiFab makeFineMask (const FabArray<FAB>& cmf, const FabArray<FAB>& fmf,
243 const IntVect& cnghost, const IntVect& ratio,
244 Periodicity const& period, int crse_value, int fine_value,
245 LayoutData<int>& has_cf,
246 MFInfo const& info = MFInfo());
247
248 MultiFab makeFineMask (const BoxArray& cba, const DistributionMapping& cdm,
249 const BoxArray& fba, const IntVect& ratio,
250 Real crse_value, Real fine_value,
251 MFInfo const& info = MFInfo());
252
254 void computeDivergence (MultiFab& divu, const Array<MultiFab const*,AMREX_SPACEDIM>& umac,
255 const Geometry& geom);
256
258 void computeGradient (MultiFab& grad, const Array<MultiFab const*,AMREX_SPACEDIM>& umac,
259 const Geometry& geom);
260
262 MultiFab ToMultiFab (const iMultiFab& imf);
264 FabArray<BaseFab<Long> > ToLongMultiFab (const iMultiFab& imf);
265
267 MultiFab periodicShift (MultiFab const& mf, IntVect const& offset,
268 Periodicity const& period);
269
271 template <typename T, typename U>
272 T cast (U const& mf_in)
273 {
274 T mf_out(mf_in.boxArray(), mf_in.DistributionMap(), mf_in.nComp(), mf_in.nGrowVect());
275
276#ifdef AMREX_USE_OMP
277#pragma omp parallel if (Gpu::notInLaunchRegion())
278#endif
279 for (MFIter mfi(mf_in); mfi.isValid(); ++mfi)
280 {
281 const Long n = mfi.fabbox().numPts() * mf_in.nComp();
282 auto * pdst = mf_out[mfi].dataPtr();
283 auto const* psrc = mf_in [mfi].dataPtr();
285 {
286 pdst[i] = static_cast<typename T::value_type>(psrc[i]); // NOLINT(bugprone-signed-char-misuse)
287 });
288 }
289 return mf_out;
290 }
291
352 template <typename Op, typename T, typename FAB, typename F,
353 std::enable_if_t<IsBaseFab<FAB>::value
354#ifndef AMREX_USE_CUDA
355 && IsCallableR<T,F,int,int,int,int>::value
356#endif
357 , int> FOO = 0>
358 BaseFab<T>
359 ReduceToPlane (int direction, Box const& domain, FabArray<FAB> const& mf, F const& f);
360
405 template <typename Op, typename FA, typename F,
406 std::enable_if_t<IsMultiFabLike_v<FA>
407#ifndef AMREX_USE_CUDA
408 && IsCallableR<typename FA::value_type,
409 F,int,int,int,int>::value
410#endif
411 , int> FOO = 0>
412 FA ReduceToPlaneMF (int direction, Box const& domain, FA const& mf, F const& f);
413
456 template <typename Op, typename FA, typename F,
457 std::enable_if_t<IsMultiFabLike_v<FA>
458#ifndef AMREX_USE_CUDA
459 && IsCallableR<typename FA::value_type,
460 F,int,int,int,int>::value
461#endif
462 , int> FOO = 0>
463 std::pair<FA,FA>
464 ReduceToPlaneMF2 (int direction, Box const& domain, FA const& mf, F const& f);
465
490 template <typename Op, typename FA, typename F,
491 std::enable_if_t<IsMultiFabLike_v<FA>
492#ifndef AMREX_USE_CUDA
493 && IsCallableR<typename FA::value_type,
494 F,int,int,int,int>::value
495#endif
496 , int> FOO = 0>
497 std::pair<FA,FA>
498 ReduceToPlaneMF2Patchy (int direction, Box const& domain, FA const& mf, F const& f);
499
515 Gpu::HostVector<Real> sumToLine (MultiFab const& mf, int icomp, int ncomp,
516 Box const& domain, int direction, bool local = false);
517
525 Real volumeWeightedSum (Vector<MultiFab const*> const& mf, int icomp,
526 Vector<Geometry> const& geom,
527 Vector<IntVect> const& ratio,
528 bool local = false);
529
543 void FourthOrderInterpFromFineToCoarse (MultiFab& cmf, int scomp, int ncomp,
544 MultiFab const& fmf,
545 IntVect const& ratio);
546
557 void FillRandom (MultiFab& mf, int scomp, int ncomp);
558
570 void FillRandomNormal (MultiFab& mf, int scomp, int ncomp, Real mean, Real stddev);
571
583 [[nodiscard]] Vector<MultiFab> convexify (Vector<MultiFab const*> const& mf,
584 Vector<IntVect> const& refinement_ratio);
585}
586
587namespace amrex {
588
589template <typename FAB>
590iMultiFab
591makeFineMask (const FabArray<FAB>& cmf, const BoxArray& fba, const IntVect& ratio,
592 int crse_value, int fine_value, MFInfo const& info)
593{
594 return makeFineMask(cmf.boxArray(), cmf.DistributionMap(), cmf.nGrowVect(),
595 fba, ratio, Periodicity::NonPeriodic(), crse_value, fine_value,
596 info);
597}
598
599template <typename FAB>
600iMultiFab
601makeFineMask (const FabArray<FAB>& cmf, const BoxArray& fba, const IntVect& ratio,
602 Periodicity const& period, int crse_value, int fine_value,
603 MFInfo const& info)
604{
605 return makeFineMask(cmf.boxArray(), cmf.DistributionMap(), cmf.nGrowVect(),
606 fba, ratio, period, crse_value, fine_value, info);
607}
608
609template <typename FAB>
610iMultiFab
612 const IntVect& cnghost, const IntVect& ratio,
613 Periodicity const& period, int crse_value, int fine_value,
614 MFInfo const& info)
615{
616 iMultiFab mask(cmf.boxArray(), cmf.DistributionMap(), 1, cnghost, info);
617 mask.setVal(crse_value);
618
619 iMultiFab foo(amrex::coarsen(fmf.boxArray(),ratio), fmf.DistributionMap(),
620 1, 0, MFInfo().SetAlloc(false));
621 const FabArrayBase::CPC& cpc = mask.getCPC(cnghost,foo,IntVect::TheZeroVector(),period);
622 mask.setVal(fine_value, cpc, 0, 1);
623
624 return mask;
625}
626
627template <typename FAB>
628iMultiFab
630 const IntVect& cnghost, const IntVect& ratio,
631 Periodicity const& period, int crse_value, int fine_value,
632 LayoutData<int>& has_cf, MFInfo const& info)
633{
634 iMultiFab mask(cmf.boxArray(), cmf.DistributionMap(), 1, cnghost, info);
635 mask.setVal(crse_value);
636
637 iMultiFab foo(amrex::coarsen(fmf.boxArray(),ratio), fmf.DistributionMap(),
638 1, 0, MFInfo().SetAlloc(false));
639 const FabArrayBase::CPC& cpc = mask.getCPC(cnghost,foo,IntVect::TheZeroVector(),period);
640 mask.setVal(fine_value, cpc, 0, 1);
641
642 has_cf = mask.RecvLayoutMask(cpc);
643
644 return mask;
645}
646
649template <typename FAB>
651 const IntVect& ratio, int ngcrse, bool mfiter_is_definitely_safe)
652{
653 AMREX_ASSERT(fine.is_nodal());
654 AMREX_ASSERT(crse.is_nodal());
655 AMREX_ASSERT(crse.nComp() == fine.nComp());
656
657 int ncomp = crse.nComp();
658 using value_type = typename FAB::value_type;
659
660 if (mfiter_is_definitely_safe || isMFIterSafe(fine, crse))
661 {
662#ifdef AMREX_USE_OMP
663#pragma omp parallel if (Gpu::notInLaunchRegion())
664#endif
665 for (MFIter mfi(crse,TilingIfNotGPU()); mfi.isValid(); ++mfi)
666 {
667 const Box& bx = mfi.growntilebox(ngcrse);
668 Array4<value_type> const& crsearr = crse.array(mfi);
669 Array4<value_type const> const& finearr = fine.const_array(mfi);
670
672 {
673 amrex_avgdown_nodes(tbx,crsearr,finearr,0,0,ncomp,ratio);
674 });
675 }
676 }
677 else
678 {
679 FabArray<FAB> ctmp(amrex::coarsen(fine.boxArray(),ratio), fine.DistributionMap(),
680 ncomp, ngcrse);
681 average_down_nodal(fine, ctmp, ratio, ngcrse);
682 crse.ParallelCopy(ctmp,0,0,ncomp,ngcrse,ngcrse);
683 }
684}
685
686// *************************************************************************************************************
687
688// Average fine cell-based MultiFab onto crse cell-centered MultiFab.
689// We do NOT assume that the coarse layout is a coarsened version of the fine layout.
690// This version does NOT use volume-weighting
691template<typename FAB>
692void average_down (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse, int scomp, int ncomp, int rr)
693{
694 average_down(S_fine,S_crse,scomp,ncomp,rr*IntVect::TheUnitVector());
695}
696
697template<typename FAB>
698void average_down (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse,
699 int scomp, int ncomp, const IntVect& ratio)
700{
701 BL_PROFILE("amrex::average_down");
702 AMREX_ASSERT(S_crse.nComp() == S_fine.nComp());
703 AMREX_ASSERT((S_crse.is_cell_centered() && S_fine.is_cell_centered()) ||
704 (S_crse.is_nodal() && S_fine.is_nodal()));
705
706 using value_type = typename FAB::value_type;
707
708 bool is_cell_centered = S_crse.is_cell_centered();
709
710 //
711 // Coarsen() the fine stuff on processors owning the fine data.
712 //
713 BoxArray crse_S_fine_BA = S_fine.boxArray(); crse_S_fine_BA.coarsen(ratio);
714
715 if (crse_S_fine_BA == S_crse.boxArray() && S_fine.DistributionMap() == S_crse.DistributionMap())
716 {
717#ifdef AMREX_USE_GPU
718 if (Gpu::inLaunchRegion() && S_crse.isFusingCandidate()) {
719 auto const& crsema = S_crse.arrays();
720 auto const& finema = S_fine.const_arrays();
721 if (is_cell_centered) {
722 ParallelFor(S_crse, IntVect(0), ncomp,
723 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k, int n) noexcept
724 {
725 amrex_avgdown(i,j,k,n,crsema[box_no],finema[box_no],scomp,scomp,ratio);
726 });
727 } else {
728 ParallelFor(S_crse, IntVect(0), ncomp,
729 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k, int n) noexcept
730 {
731 amrex_avgdown_nodes(i,j,k,n,crsema[box_no],finema[box_no],scomp,scomp,ratio);
732 });
733 }
734 if (!Gpu::inNoSyncRegion()) {
736 }
737 } else
738#endif
739 {
740#ifdef AMREX_USE_OMP
741#pragma omp parallel if (Gpu::notInLaunchRegion())
742#endif
743 for (MFIter mfi(S_crse,TilingIfNotGPU()); mfi.isValid(); ++mfi)
744 {
745 // NOTE: The tilebox is defined at the coarse level.
746 const Box& bx = mfi.tilebox();
747 Array4<value_type> const& crsearr = S_crse.array(mfi);
748 Array4<value_type const> const& finearr = S_fine.const_array(mfi);
749
750 if (is_cell_centered) {
751 AMREX_HOST_DEVICE_PARALLEL_FOR_4D(bx, ncomp, i, j, k, n,
752 {
753 amrex_avgdown(i,j,k,n,crsearr,finearr,scomp,scomp,ratio);
754 });
755 } else {
756 AMREX_HOST_DEVICE_PARALLEL_FOR_4D(bx, ncomp, i, j, k, n,
757 {
758 amrex_avgdown_nodes(i,j,k,n,crsearr,finearr,scomp,scomp,ratio);
759 });
760 }
761 }
762 }
763 }
764 else
765 {
766 FabArray<FAB> crse_S_fine(crse_S_fine_BA, S_fine.DistributionMap(), ncomp, 0, MFInfo(),DefaultFabFactory<FAB>());
767
768#ifdef AMREX_USE_GPU
769 if (Gpu::inLaunchRegion() && crse_S_fine.isFusingCandidate()) {
770 auto const& crsema = crse_S_fine.arrays();
771 auto const& finema = S_fine.const_arrays();
772 if (is_cell_centered) {
773 ParallelFor(crse_S_fine, IntVect(0), ncomp,
774 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k, int n) noexcept
775 {
776 amrex_avgdown(i,j,k,n,crsema[box_no],finema[box_no],0,scomp,ratio);
777 });
778 } else {
779 ParallelFor(crse_S_fine, IntVect(0), ncomp,
780 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k, int n) noexcept
781 {
782 amrex_avgdown_nodes(i,j,k,n,crsema[box_no],finema[box_no],0,scomp,ratio);
783 });
784 }
785 if (!Gpu::inNoSyncRegion()) {
787 }
788 } else
789#endif
790 {
791#ifdef AMREX_USE_OMP
792#pragma omp parallel if (Gpu::notInLaunchRegion())
793#endif
794 for (MFIter mfi(crse_S_fine,TilingIfNotGPU()); mfi.isValid(); ++mfi)
795 {
796 // NOTE: The tilebox is defined at the coarse level.
797 const Box& bx = mfi.tilebox();
798 Array4<value_type> const& crsearr = crse_S_fine.array(mfi);
799 Array4<value_type const> const& finearr = S_fine.const_array(mfi);
800
801 // NOTE: We copy from component scomp of the fine fab into component 0 of the crse fab
802 // because the crse fab is a temporary which was made starting at comp 0, it is
803 // not part of the actual crse multifab which came in.
804
805 if (is_cell_centered) {
806 AMREX_HOST_DEVICE_PARALLEL_FOR_4D(bx, ncomp, i, j, k, n,
807 {
808 amrex_avgdown(i,j,k,n,crsearr,finearr,0,scomp,ratio);
809 });
810 } else {
811 AMREX_HOST_DEVICE_PARALLEL_FOR_4D(bx, ncomp, i, j, k, n,
812 {
813 amrex_avgdown_nodes(i,j,k,n,crsearr,finearr,0,scomp,ratio);
814 });
815 }
816 }
817 }
818
819 S_crse.ParallelCopy(crse_S_fine,0,scomp,ncomp);
820 }
821}
822
823
824
825
826
834template <typename F>
835Real
836NormHelper (const MultiFab& x, int xcomp,
837 const MultiFab& y, int ycomp,
838 F const& f,
839 int numcomp, IntVect nghost, bool local)
840{
841 BL_ASSERT(x.boxArray() == y.boxArray());
842 BL_ASSERT(x.DistributionMap() == y.DistributionMap());
843 BL_ASSERT(x.nGrowVect().allGE(nghost) && y.nGrowVect().allGE(nghost));
844
845 Real sm = Real(0.0);
846#ifdef AMREX_USE_GPU
847 if (Gpu::inLaunchRegion()) {
848 auto const& xma = x.const_arrays();
849 auto const& yma = y.const_arrays();
851 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k) noexcept -> GpuTuple<Real>
852 {
853 Real t = Real(0.0);
854 auto const& xfab = xma[box_no];
855 auto const& yfab = yma[box_no];
856 for (int n = 0; n < numcomp; ++n) {
857 t += f(xfab(i,j,k,xcomp+n) , yfab(i,j,k,ycomp+n));
858 }
859 return t;
860 });
861 } else
862#endif
863 {
864#ifdef AMREX_USE_OMP
865#pragma omp parallel if (!system::regtest_reduction) reduction(+:sm)
866#endif
867 for (MFIter mfi(x,true); mfi.isValid(); ++mfi)
868 {
869 Box const& bx = mfi.growntilebox(nghost);
870 Array4<Real const> const& xfab = x.const_array(mfi);
871 Array4<Real const> const& yfab = y.const_array(mfi);
872 AMREX_LOOP_4D(bx, numcomp, i, j, k, n,
873 {
874 sm += f(xfab(i,j,k,xcomp+n) , yfab(i,j,k,ycomp+n));
875 });
876 }
877 }
878
879 if (!local) {
881 }
882
883 return sm;
884}
885
894template <typename MMF, typename Pred, typename F>
895Real
896NormHelper (const MMF& mask,
897 const MultiFab& x, int xcomp,
898 const MultiFab& y, int ycomp,
899 Pred const& pf,
900 F const& f,
901 int numcomp, IntVect nghost, bool local)
902{
903 BL_ASSERT(x.boxArray() == y.boxArray());
904 BL_ASSERT(x.boxArray() == mask.boxArray());
905 BL_ASSERT(x.DistributionMap() == y.DistributionMap());
906 BL_ASSERT(x.DistributionMap() == mask.DistributionMap());
907 BL_ASSERT(x.nGrowVect().allGE(nghost) && y.nGrowVect().allGE(nghost));
908 BL_ASSERT(mask.nGrowVect().allGE(nghost));
909
910 Real sm = Real(0.0);
911#ifdef AMREX_USE_GPU
912 if (Gpu::inLaunchRegion()) {
913 auto const& xma = x.const_arrays();
914 auto const& yma = y.const_arrays();
915 auto const& mma = mask.const_arrays();
917 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k) noexcept -> GpuTuple<Real>
918 {
919 Real t = Real(0.0);
920 if (pf(mma[box_no](i,j,k))) {
921 auto const& xfab = xma[box_no];
922 auto const& yfab = yma[box_no];
923 for (int n = 0; n < numcomp; ++n) {
924 t += f(xfab(i,j,k,xcomp+n) , yfab(i,j,k,ycomp+n));
925 }
926 }
927 return t;
928 });
929 } else
930#endif
931 {
932#ifdef AMREX_USE_OMP
933#pragma omp parallel if (!system::regtest_reduction) reduction(+:sm)
934#endif
935 for (MFIter mfi(x,true); mfi.isValid(); ++mfi)
936 {
937 Box const& bx = mfi.growntilebox(nghost);
938 Array4<Real const> const& xfab = x.const_array(mfi);
939 Array4<Real const> const& yfab = y.const_array(mfi);
940 auto const& mfab = mask.const_array(mfi);
941 AMREX_LOOP_4D(bx, numcomp, i, j, k, n,
942 {
943 if (pf(mfab(i,j,k))) {
944 sm += f(xfab(i,j,k,xcomp+n) , yfab(i,j,k,ycomp+n));
945 }
946 });
947 }
948 }
949
950 if (!local) {
952 }
953
954 return sm;
955}
956
957template <typename CMF, typename FMF,
958 std::enable_if_t<IsFabArray_v<CMF> && IsFabArray_v<FMF>, int> FOO>
959void average_face_to_cellcenter (CMF& cc, int dcomp,
961 int ngrow)
962{
963 IntVect ng_vect(ngrow);
964 average_face_to_cellcenter(cc, dcomp, fc, ng_vect);
965}
966
967template <typename CMF, typename FMF,
968 std::enable_if_t<IsFabArray_v<CMF> && IsFabArray_v<FMF>, int> FOO>
969void average_face_to_cellcenter (CMF& cc, int dcomp,
971 IntVect const& ng_vect)
972{
973 AMREX_ASSERT(cc.nComp() >= dcomp + AMREX_SPACEDIM);
974 AMREX_ASSERT(fc[0]->nComp() == 1);
975
976#ifdef AMREX_USE_GPU
977 if (Gpu::inLaunchRegion() && cc.isFusingCandidate()) {
978 auto const& ccma = cc.arrays();
979 AMREX_D_TERM(auto const& fxma = fc[0]->const_arrays();,
980 auto const& fyma = fc[1]->const_arrays();,
981 auto const& fzma = fc[2]->const_arrays(););
982 ParallelFor(cc, ng_vect,
983 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k) noexcept
984 {
985#if (AMREX_SPACEDIM == 1)
986 GeometryData gd{};
987 gd.coord = 0;
988#endif
989 amrex_avg_fc_to_cc(i,j,k, ccma[box_no], AMREX_D_DECL(fxma[box_no],
990 fyma[box_no],
991 fzma[box_no]),
992 dcomp
993#if (AMREX_SPACEDIM == 1)
994 , gd
995#endif
996 );
997 });
998 if (!Gpu::inNoSyncRegion()) {
1000 }
1001 } else
1002#endif
1003 {
1004#ifdef AMREX_USE_OMP
1005#pragma omp parallel if (Gpu::notInLaunchRegion())
1006#endif
1007 for (MFIter mfi(cc,TilingIfNotGPU()); mfi.isValid(); ++mfi)
1008 {
1009 const Box bx = mfi.growntilebox(ng_vect);
1010 auto const& ccarr = cc.array(mfi);
1011 AMREX_D_TERM(auto const& fxarr = fc[0]->const_array(mfi);,
1012 auto const& fyarr = fc[1]->const_array(mfi);,
1013 auto const& fzarr = fc[2]->const_array(mfi););
1014
1015#if (AMREX_SPACEDIM == 1)
1017 {
1018 GeometryData gd;
1019 gd.coord = 0;
1020 amrex_avg_fc_to_cc(i,j,k, ccarr, fxarr, dcomp, gd);
1021 });
1022#else
1024 {
1025 amrex_avg_fc_to_cc(i,j,k, ccarr, AMREX_D_DECL(fxarr,fyarr,fzarr), dcomp);
1026 });
1027#endif
1028 }
1029 }
1030}
1031
1032template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int>>
1034 const Vector<MF*>& crse,
1035 const IntVect& ratio, int ngcrse)
1036{
1037 AMREX_ASSERT(fine.size() == AMREX_SPACEDIM && crse.size() == AMREX_SPACEDIM);
1039 {{AMREX_D_DECL(fine[0],fine[1],fine[2])}},
1041 {{AMREX_D_DECL(crse[0],crse[1],crse[2])}},
1042 ratio, ngcrse);
1043}
1044
1045template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int>>
1047 const Vector<MF*>& crse, int ratio, int ngcrse)
1048{
1049 average_down_faces(fine,crse,IntVect{ratio},ngcrse);
1050}
1051
1052template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int>>
1055 int ratio, int ngcrse)
1056{
1057 average_down_faces(fine,crse,IntVect{ratio},ngcrse);
1058}
1059
1060template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int>>
1063 const IntVect& ratio, int ngcrse)
1064{
1065 for (int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
1066 average_down_faces(*fine[idim], *crse[idim], ratio, ngcrse);
1067 }
1068}
1069
1070template <typename FAB>
1072 const IntVect& ratio, int ngcrse)
1073{
1074 BL_PROFILE("average_down_faces");
1075
1076 AMREX_ASSERT(crse.nComp() == fine.nComp());
1077 AMREX_ASSERT(fine.ixType() == crse.ixType());
1078 const auto type = fine.ixType();
1079 int dir;
1080 for (dir = 0; dir < AMREX_SPACEDIM; ++dir) {
1081 if (type.nodeCentered(dir)) { break; }
1082 }
1083 auto tmptype = type;
1084 tmptype.unset(dir);
1085 if (dir >= AMREX_SPACEDIM || !tmptype.cellCentered()) {
1086 amrex::Abort("average_down_faces: not face index type");
1087 }
1088 const int ncomp = crse.nComp();
1089 if (isMFIterSafe(fine, crse))
1090 {
1091#ifdef AMREX_USE_GPU
1092 if (Gpu::inLaunchRegion() && crse.isFusingCandidate()) {
1093 auto const& crsema = crse.arrays();
1094 auto const& finema = fine.const_arrays();
1095 ParallelFor(crse, IntVect(ngcrse), ncomp,
1096 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k, int n) noexcept
1097 {
1098 amrex_avgdown_faces(i,j,k,n, crsema[box_no], finema[box_no], 0, 0, ratio, dir);
1099 });
1100 if (!Gpu::inNoSyncRegion()) {
1102 }
1103 } else
1104#endif
1105 {
1106#ifdef AMREX_USE_OMP
1107#pragma omp parallel if (Gpu::notInLaunchRegion())
1108#endif
1109 for (MFIter mfi(crse,TilingIfNotGPU()); mfi.isValid(); ++mfi)
1110 {
1111 const Box& bx = mfi.growntilebox(ngcrse);
1112 auto const& crsearr = crse.array(mfi);
1113 auto const& finearr = fine.const_array(mfi);
1114 AMREX_HOST_DEVICE_PARALLEL_FOR_4D(bx, ncomp, i, j, k, n,
1115 {
1116 amrex_avgdown_faces(i,j,k,n, crsearr, finearr, 0, 0, ratio, dir);
1117 });
1118 }
1119 }
1120 }
1121 else
1122 {
1123 FabArray<FAB> ctmp(amrex::coarsen(fine.boxArray(),ratio), fine.DistributionMap(),
1124 ncomp, ngcrse, MFInfo(), DefaultFabFactory<FAB>());
1125 average_down_faces(fine, ctmp, ratio, ngcrse);
1126 crse.ParallelCopy(ctmp,0,0,ncomp,ngcrse,ngcrse);
1127 }
1128}
1129
1130template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int>>
1133 const IntVect& ratio, const Geometry& crse_geom)
1134{
1135 for (int idim = 0; idim < AMREX_SPACEDIM; ++idim) {
1136 average_down_faces(*fine[idim], *crse[idim], ratio, crse_geom);
1137 }
1138}
1139
1140template <typename FAB>
1142 const IntVect& ratio, const Geometry& crse_geom)
1143{
1144 FabArray<FAB> ctmp(amrex::coarsen(fine.boxArray(),ratio), fine.DistributionMap(),
1145 crse.nComp(), 0);
1146 average_down_faces(fine, ctmp, ratio, 0);
1147 crse.ParallelCopy(ctmp,0,0,crse.nComp(),0,0,crse_geom.periodicity());
1148}
1149
1150template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int> FOO>
1152{
1153 using T = typename MF::value_type;
1154 const int ncomp = mf.nComp();
1155 Gpu::DeviceVector<T> dv(ncomp);
1156 auto* dp = dv.data();
1157 bool found = false;
1158 auto loc = cell.dim3();
1159 for (MFIter mfi(mf); mfi.isValid() && !found; ++mfi)
1160 {
1161 Box const& box = mfi.validbox();
1162 if (box.contains(cell)) {
1163 found = true;
1164 auto const& fab = mf.const_array(mfi);
1165 amrex::ParallelFor(1, [=] AMREX_GPU_DEVICE (int) noexcept
1166 {
1167 for (int n = 0; n < ncomp; ++n) {
1168 dp[n] = fab(loc.x,loc.y,loc.z,n);
1169 }
1170 });
1171 }
1172 }
1173 Vector<T> hv;
1174 if (found) {
1175 hv.resize(ncomp);
1176 Gpu::copy(Gpu::deviceToHost, dv.begin(), dv.end(), hv.begin());
1177 }
1178 return hv;
1179}
1180
1181template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int> FOO>
1182MF get_line_data (MF const& mf, int dir, IntVect const& cell, Box const& bnd_bx)
1183{
1184 bool do_bnd = (!bnd_bx.isEmpty());
1185
1186 BoxArray const& ba = mf.boxArray();
1187 DistributionMapping const& dm = mf.DistributionMap();
1188 const auto nboxes = static_cast<int>(ba.size());
1189
1190 BoxList bl(ba.ixType());
1191 Vector<int> procmap;
1192 Vector<int> index_map;
1193 if (!do_bnd) {
1194 for (int i = 0; i < nboxes; ++i) {
1195 Box const& b = ba[i];
1196 IntVect lo = cell;
1197 lo[dir] = b.smallEnd(dir);
1198 if (b.contains(lo)) {
1199 IntVect hi = lo;
1200 hi[dir] = b.bigEnd(dir);
1201 Box b1d(lo,hi,b.ixType());
1202 bl.push_back(b1d);
1203 procmap.push_back(dm[i]);
1204 index_map.push_back(i);
1205 }
1206 }
1207 } else {
1208 for (int i = 0; i < nboxes; ++i) {
1209 Box const& b = ba[i];
1210 Box const& b1d = bnd_bx & b;
1211 if (b1d.ok()) {
1212 bl.push_back(b1d);
1213 procmap.push_back(dm[i]);
1214 index_map.push_back(i);
1215 }
1216 }
1217 }
1218
1219 if (bl.isEmpty()) {
1220 return MF();
1221 } else {
1222 BoxArray rba(std::move(bl));
1223 DistributionMapping rdm(std::move(procmap));
1224 MF rmf(rba, rdm, mf.nComp(), IntVect(0),
1225 MFInfo().SetArena(mf.arena()));
1226#ifdef AMREX_USE_OMP
1227#pragma omp parallel if (Gpu::notInLaunchRegion())
1228#endif
1229 for (MFIter mfi(rmf); mfi.isValid(); ++mfi) {
1230 Box const& b = mfi.validbox();
1231 auto const& dfab = rmf.array(mfi);
1232 auto const& sfab = mf.const_array(index_map[mfi.index()]);
1233 amrex::ParallelFor(b, mf.nComp(),
1234 [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) noexcept
1235 {
1236 dfab(i,j,k,n) = sfab(i,j,k,n);
1237 });
1238 }
1239 return rmf;
1240 }
1241}
1242
1244namespace detail {
1245template <typename Op, typename T, typename F>
1246void reduce_to_plane (Array4<T> const& ar, int direction, Box const& bx, int box_no,
1247 F const& f)
1248{
1249#if defined(AMREX_USE_GPU)
1250 Box b2d = bx;
1251 b2d.setRange(direction,0);
1252 const auto blo = amrex::lbound(bx);
1253 const auto len = amrex::length(bx);
1254 constexpr int nthreads = 128;
1255 auto nblocks = static_cast<int>(b2d.numPts());
1256#ifdef AMREX_USE_SYCL
1257 constexpr std::size_t shared_mem_bytes = sizeof(T)*Gpu::Device::warp_size;
1258 amrex::launch<nthreads>(nblocks, shared_mem_bytes, Gpu::gpuStream(),
1259 [=] AMREX_GPU_DEVICE (Gpu::Handler const& h)
1260 {
1261 int bid = h.blockIdx();
1262 int tid = h.threadIdx();
1263#else
1264 amrex::launch<nthreads>(nblocks, Gpu::gpuStream(),
1265 [=] AMREX_GPU_DEVICE ()
1266 {
1267 int bid = blockIdx.x;
1268 int tid = threadIdx.x;
1269#endif
1270 T tmp;
1271 Op().init(tmp);
1272 T* p;
1273 if (direction == 0) {
1274 int k = bid / len.y;
1275 int j = bid - k*len.y;
1276 k += blo.z;
1277 j += blo.y;
1278 for (int i = blo.x + tid; i < blo.x+len.x; i += nthreads) {
1279 Op().local_update(tmp, f(box_no,i,j,k));
1280 }
1281 p = ar.ptr(0,j,k);
1282 } else if (direction == 1) {
1283 int k = bid / len.x;
1284 int i = bid - k*len.x;
1285 k += blo.z;
1286 i += blo.x;
1287 for (int j = blo.y + tid; j < blo.y+len.y; j += nthreads) {
1288 Op().local_update(tmp, f(box_no,i,j,k));
1289 }
1290 p = ar.ptr(i,0,k);
1291 } else {
1292 int j = bid / len.x;
1293 int i = bid - j*len.x;
1294 j += blo.y;
1295 i += blo.x;
1296 for (int k = blo.z + tid; k < blo.z+len.z; k += nthreads) {
1297 Op().local_update(tmp, f(box_no,i,j,k));
1298 }
1299 p = ar.ptr(i,j,0);
1300 }
1301#ifdef AMREX_USE_SYCL
1302 Op().template parallel_update<T>(*p, tmp, h);
1303#else
1304 Op().template parallel_update<T,nthreads>(*p, tmp);
1305#endif
1306 });
1307#else
1308 // CPU
1309 if (direction == 0) {
1310 AMREX_LOOP_3D(bx, i, j, k,
1311 {
1312 Op().local_update(ar(0,j,k), f(box_no,i,j,k));
1313 });
1314 } else if (direction == 1) {
1315 AMREX_LOOP_3D(bx, i, j, k,
1316 {
1317 Op().local_update(ar(i,0,k), f(box_no,i,j,k));
1318 });
1319 } else {
1320 AMREX_LOOP_3D(bx, i, j, k,
1321 {
1322 Op().local_update(ar(i,j,0), f(box_no,i,j,k));
1323 });
1324 }
1325#endif
1326}
1327}
1329
1330template <typename Op, typename T, typename FAB, typename F,
1331 std::enable_if_t<IsBaseFab<FAB>::value
1332#ifndef AMREX_USE_CUDA
1333 && IsCallableR<T,F,int,int,int,int>::value
1334#endif
1335 , int> FOO>
1336BaseFab<T>
1337ReduceToPlane (int direction, Box const& a_domain, FabArray<FAB> const& mf, F const& f)
1338{
1339 Box const domain = amrex::convert(a_domain, mf.ixType());
1340
1341 Box domain2d = domain;
1342 domain2d.setRange(direction, 0);
1343
1344 T initval;
1345 Op().init(initval);
1346
1347 BaseFab<T> r(domain2d);
1348 r.template setVal<RunOn::Device>(initval);
1349 auto const& ar = r.array();
1350
1351 for (MFIter mfi(mf,MFItInfo().UseDefaultStream().DisableDeviceSync());
1352 mfi.isValid(); ++mfi)
1353 {
1354 Box bx = mfi.validbox() & domain;
1355 if (bx.ok()) {
1356 int box_no = mfi.LocalIndex();
1357 detail::reduce_to_plane<Op, T>(ar, direction, bx, box_no, f);
1358 }
1359 }
1360 if (!Gpu::inNoSyncRegion()) {
1362 }
1363
1364 return r;
1365}
1366
1368namespace detail {
1369template <typename Op, typename FA, typename F>
1370FA reduce_to_plane (int direction, Box const& domain, FA const& mf, F const& f)
1371{
1372 using T = typename FA::value_type;
1373
1374 Box const ndomain = amrex::convert(domain, mf.ixType());
1375
1376 auto npts = amrex::convert(mf.boxArray(),IntVect(0)).numPts();
1377 if (npts != amrex::convert(domain, amrex::IntVect(0)).numPts()) {
1378 amrex::Abort("ReduceToPlaneMF: mf's BoxArray must have a rectangular domain.");
1379 }
1380
1381 T initval;
1382 Op().init(initval);
1383
1384 BoxList bl = mf.boxArray().boxList();
1385 for (auto& b : bl) {
1386 b.setRange(direction, 0);
1387 }
1388 BoxArray ba(std::move(bl));
1389 FA tmpfa(ba, mf.DistributionMap(), 1, 0);
1390 tmpfa.setVal(initval);
1391
1392 for (MFIter mfi(mf); mfi.isValid(); ++mfi)
1393 {
1394 Box bx = mfi.validbox() & ndomain;
1395 if (bx.ok()) {
1396 int box_no = mfi.LocalIndex();
1397 detail::reduce_to_plane<Op, T>(tmpfa.array(mfi), direction, bx, box_no, f);
1398 }
1399 }
1400
1401 return tmpfa;
1402}
1403}
1405
1406template <typename Op, typename FA, typename F,
1407 std::enable_if_t<IsMultiFabLike_v<FA>
1408#ifndef AMREX_USE_CUDA
1409 && IsCallableR<typename FA::value_type,
1410 F,int,int,int,int>::value
1411#endif
1412 , int> FOO>
1413FA ReduceToPlaneMF (int direction, Box const& domain, FA const& mf, F const& f)
1414{
1415 auto [fa3, fa2] = ReduceToPlaneMF2<Op>(direction, domain, mf, f);
1416 fa3.ParallelCopy(fa2);
1417 return std::move(fa3);
1418}
1419
1420template <typename Op, typename FA, typename F,
1421 std::enable_if_t<IsMultiFabLike_v<FA>
1422#ifndef AMREX_USE_CUDA
1423 && IsCallableR<typename FA::value_type,
1424 F,int,int,int,int>::value
1425#endif
1426 , int> FOO>
1427std::pair<FA,FA>
1428ReduceToPlaneMF2 (int direction, Box const& domain, FA const& mf, F const& f)
1429{
1430 using T = typename FA::value_type;
1431
1432 T initval;
1433 Op().init(initval);
1434
1435 auto tmpmf = detail::reduce_to_plane<Op>(direction, domain, mf, f);
1436
1437 BoxList bl2d(mf.ixType());
1438 Vector<int> procmap2d;
1439 auto const& ba3d = mf.boxArray();
1440 auto const& dm3d = mf.DistributionMap();
1441 int dlo = domain.smallEnd(direction);
1442 for (int i = 0, N = mf.size(); i < N; ++i) {
1443 Box b = ba3d[i];
1444 if (b.smallEnd(direction) <= dlo && dlo <= b.bigEnd(direction)) {
1445 b.setRange(direction, 0);
1446 bl2d.push_back(b);
1447 procmap2d.push_back(dm3d[i]);
1448 }
1449 }
1450
1451 BoxArray ba2d(std::move(bl2d));
1452 DistributionMapping dm2d(std::move(procmap2d));
1453
1454 FA mf2d(ba2d, dm2d, 1, 0);
1455 mf2d.setVal(initval);
1456
1457 static_assert(std::is_same_v<Op, ReduceOpSum>, "Currently only ReduceOpSum is supported.");
1458 mf2d.ParallelAdd(tmpmf);
1459
1460 return std::make_pair(std::move(tmpmf), std::move(mf2d));
1461}
1462
1463template <typename Op, typename FA, typename F,
1464 std::enable_if_t<IsMultiFabLike_v<FA>
1465#ifndef AMREX_USE_CUDA
1466 && IsCallableR<typename FA::value_type,
1467 F,int,int,int,int>::value
1468#endif
1469 , int> FOO>
1470std::pair<FA,FA>
1471ReduceToPlaneMF2Patchy (int direction, Box const& domain, FA const& mf, F const& f)
1472{
1473 AMREX_ALWAYS_ASSERT_WITH_MESSAGE(mf.is_cell_centered(), "ReduceToPlaneMF2Patchy: only cell-centered data is supported.");
1474
1475 Box const ndomain = amrex::convert(domain, mf.ixType());
1476 AMREX_ALWAYS_ASSERT_WITH_MESSAGE(ndomain.contains(mf.boxArray().minimalBox()),
1477 "ReduceToPlaneMF2Patchy: subdomains are not supported; domain must cover the full extent of mf.boxArray().");
1478
1479 using T = typename FA::value_type;
1480 T initval;
1481 Op().init(initval);
1482
1483 BoxList bl_patch2d = mf.boxArray().boxList();
1484 for (auto& b : bl_patch2d) {
1485 b.setRange(direction, 0);
1486 }
1487 BoxArray ba_patch2d(std::move(bl_patch2d));
1488 FA tmpmf(ba_patch2d, mf.DistributionMap(), 1, 0);
1489 tmpmf.setVal(initval);
1490
1491 for (MFIter mfi(mf); mfi.isValid(); ++mfi)
1492 {
1493 Box bx = mfi.validbox() & ndomain;
1494 if (bx.ok()) {
1495 int box_no = mfi.LocalIndex();
1496 detail::reduce_to_plane<Op, T>(tmpmf.array(mfi), direction, bx, box_no, f);
1497 }
1498 }
1499
1500 BoxList bl_unique = ba_patch2d.boxList();
1501 bl_unique.simplify();
1502 BoxArray ba2d(std::move(bl_unique));
1503 ba2d.removeOverlap();
1504 DistributionMapping dm2d(ba2d);
1505
1506 FA mf2d(ba2d, dm2d, 1, 0);
1507 mf2d.setVal(initval);
1508
1509 static_assert(std::is_same_v<Op, ReduceOpSum>, "Currently only ReduceOpSum is supported.");
1510 mf2d.ParallelAdd(tmpmf);
1511
1512 return std::make_pair(std::move(tmpmf), std::move(mf2d));
1513}
1514
1515}
1516
1517#endif
#define BL_PROFILE(a)
Definition AMReX_BLProfiler.H:551
#define AMREX_ALWAYS_ASSERT_WITH_MESSAGE(EX, MSG)
Definition AMReX_BLassert.H:49
#define BL_ASSERT(EX)
Definition AMReX_BLassert.H:39
#define AMREX_ASSERT(EX)
Definition AMReX_BLassert.H:38
#define AMREX_HOST_DEVICE_PARALLEL_FOR_1D(...)
Definition AMReX_GpuLaunchMacrosC.nolint.H:109
#define AMREX_HOST_DEVICE_PARALLEL_FOR_3D(...)
Definition AMReX_GpuLaunchMacrosC.nolint.H:110
#define AMREX_HOST_DEVICE_PARALLEL_FOR_4D(...)
Definition AMReX_GpuLaunchMacrosC.nolint.H:111
#define AMREX_LAUNCH_HOST_DEVICE_LAMBDA(...)
Definition AMReX_GpuLaunch.nolint.H:16
#define AMREX_GPU_DEVICE
Definition AMReX_GpuQualifiers.H:18
Array4< int const > offset
Definition AMReX_HypreMLABecLap.cpp:1139
Real * pdst
Definition AMReX_HypreMLABecLap.cpp:1140
Array4< Real > fine
Definition AMReX_InterpFaceRegister.cpp:90
Array4< int const > mask
Definition AMReX_InterpFaceRegister.cpp:93
Array4< Real const > crse
Definition AMReX_InterpFaceRegister.cpp:92
#define AMREX_LOOP_3D(bx, i, j, k, block)
Definition AMReX_Loop.nolint.H:4
#define AMREX_LOOP_4D(bx, ncomp, i, j, k, n, block)
Definition AMReX_Loop.nolint.H:16
#define AMREX_D_TERM(a, b, c)
Definition AMReX_SPACE.H:172
#define AMREX_D_DECL(a, b, c)
Definition AMReX_SPACE.H:171
A FortranArrayBox(FAB)-like object.
Definition AMReX_BaseFab.H:190
Array4< T const > array() const noexcept
Definition AMReX_BaseFab.H:382
A collection of Boxes stored in an Array.
Definition AMReX_BoxArray.H:568
IndexType ixType() const noexcept
Return index type of this BoxArray.
Definition AMReX_BoxArray.H:858
BoxList boxList() const
Create a BoxList from this BoxArray.
Definition AMReX_BoxArray.cpp:949
void removeOverlap(bool simplify=true)
Change the BoxArray to one with no overlap and then simplify it (see the simplify function in BoxList...
Definition AMReX_BoxArray.cpp:1456
BoxArray & coarsen(int refinement_ratio)
Coarsen each Box in the BoxArray to the specified ratio.
Definition AMReX_BoxArray.cpp:672
Long size() const noexcept
Return the number of boxes in the BoxArray.
Definition AMReX_BoxArray.H:615
A class for managing a List of Boxes that share a common IndexType. This class implements operations ...
Definition AMReX_BoxList.H:52
bool isEmpty() const noexcept
Is this BoxList empty?
Definition AMReX_BoxList.H:135
int simplify(bool best=false)
Merge adjacent Boxes in this BoxList. Return the number of Boxes merged. If "best" is specified we do...
Definition AMReX_BoxList.cpp:654
void push_back(const Box &bn)
Append a Box to this BoxList.
Definition AMReX_BoxList.H:93
__host__ __device__ const IntVectND< dim > & bigEnd() const &noexcept
Return the inclusive upper bound of the box.
Definition AMReX_Box.H:123
__host__ __device__ bool isEmpty() const noexcept
Checks if it is an empty BoxND.
Definition AMReX_Box.H:204
__host__ __device__ bool contains(const IntVectND< dim > &p) const noexcept
Return true if argument is contained within BoxND.
Definition AMReX_Box.H:212
__host__ __device__ IndexTypeND< dim > ixType() const noexcept
Return the indexing type.
Definition AMReX_Box.H:135
__host__ __device__ BoxND & setRange(int dir, int sm_index, int n_cells=1) noexcept
Set the entire range in a given direction, starting at sm_index with length n_cells....
Definition AMReX_Box.H:1108
__host__ __device__ bool ok() const noexcept
Checks if it is a proper BoxND (including a valid type).
Definition AMReX_Box.H:208
__host__ __device__ const IntVectND< dim > & smallEnd() const &noexcept
Return the inclusive lower bound of the box.
Definition AMReX_Box.H:111
Definition AMReX_FabFactory.H:76
Calculates the distribution of FABs to MPI processes.
Definition AMReX_DistributionMapping.H:43
IntVect nGrowVect() const noexcept
Definition AMReX_FabArrayBase.H:80
bool isFusingCandidate() const noexcept
Is this a good candidate for kernel fusing?
Definition AMReX_FabArrayBase.cpp:2704
bool is_cell_centered() const noexcept
This tests on whether the FabArray is cell-centered.
Definition AMReX_FabArrayBase.cpp:2698
bool is_nodal() const noexcept
This tests on whether the FabArray is fully nodal.
Definition AMReX_FabArrayBase.cpp:2686
IndexType ixType() const noexcept
Return index type.
Definition AMReX_FabArrayBase.H:86
const DistributionMapping & DistributionMap() const noexcept
Return constant reference to associated DistributionMapping.
Definition AMReX_FabArrayBase.H:131
int nComp() const noexcept
Return number of variables (aka components) associated with each point.
Definition AMReX_FabArrayBase.H:83
const BoxArray & boxArray() const noexcept
Return a constant reference to the BoxArray that defines the valid region associated with this FabArr...
Definition AMReX_FabArrayBase.H:95
An Array of FortranArrayBox(FAB)-like Objects.
Definition AMReX_FabArray.H:349
Array4< typename FabArray< FAB >::value_type const > const_array(const MFIter &mfi) const noexcept
Definition AMReX_FabArray.H:589
void ParallelCopy(const FabArray< FAB > &src, const Periodicity &period=Periodicity::NonPeriodic(), CpOp op=FabArrayBase::COPY)
Definition AMReX_FabArray.H:849
Array4< typename FabArray< FAB >::value_type const > array(const MFIter &mfi) const noexcept
Definition AMReX_FabArray.H:565
MultiArray4< typename FabArray< FAB >::value_type > arrays() noexcept
Definition AMReX_FabArray.H:637
MultiArray4< typename FabArray< FAB >::value_type const > const_arrays() const noexcept
Definition AMReX_FabArray.H:651
Rectangular problem domain geometry.
Definition AMReX_Geometry.H:74
Periodicity periodicity() const noexcept
Definition AMReX_Geometry.H:356
GPU-compatible tuple.
Definition AMReX_Tuple.H:98
static constexpr int warp_size
Definition AMReX_GpuDevice.H:236
__host__ static __device__ constexpr IntVectND< dim > TheUnitVector() noexcept
This static member function returns a reference to a constant IntVectND object, all of whose dim argu...
Definition AMReX_IntVect.H:689
__host__ static __device__ constexpr IntVectND< dim > TheZeroVector() noexcept
This static member function returns a reference to a constant IntVectND object, all of whose dim argu...
Definition AMReX_IntVect.H:679
__host__ __device__ constexpr Dim3 dim3() const noexcept
Definition AMReX_IntVect.H:173
a one-thingy-per-box distributed object
Definition AMReX_LayoutData.H:13
Iterator for looping ever tiles and boxes of amrex::FabArray based containers.
Definition AMReX_MFIter.H:85
bool isValid() const noexcept
Is the iterator valid i.e. is it associated with a FAB?
Definition AMReX_MFIter.H:169
A collection (stored as an array) of FArrayBox objects.
Definition AMReX_MultiFab.H:40
Dynamically allocated vector for trivially copyable data.
Definition AMReX_PODVector.H:308
iterator begin() noexcept
Definition AMReX_PODVector.H:674
iterator end() noexcept
Definition AMReX_PODVector.H:678
T * data() noexcept
Definition AMReX_PODVector.H:666
This provides length of period for periodic domains. 0 means it is not periodic in that direction....
Definition AMReX_Periodicity.H:17
static const Periodicity & NonPeriodic() noexcept
Definition AMReX_Periodicity.cpp:52
This class is a thin wrapper around std::vector. Unlike vector, Vector::operator[] provides bound che...
Definition AMReX_Vector.H:28
A Collection of IArrayBoxes.
Definition AMReX_iMultiFab.H:34
amrex_real Real
Floating Point Type for Fields.
Definition AMReX_REAL.H:79
amrex_long Long
Definition AMReX_INT.H:30
__host__ __device__ Dim3 length(Array4< T > const &a) noexcept
Return the spatial extents of an Array4 in Dim3 form.
Definition AMReX_Array4.H:1345
__host__ __device__ Dim3 lbound(Array4< T > const &a) noexcept
Return the inclusive lower bounds of an Array4 in Dim3 form.
Definition AMReX_Array4.H:1317
__host__ __device__ BoxND< dim > coarsen(const BoxND< dim > &b, int ref_ratio) noexcept
Coarsen BoxND by given (positive) coarsening ratio.
Definition AMReX_Box.H:1409
std::array< T, N > Array
Definition AMReX_Array.H:26
void Sum(T &v, MPI_Comm comm)
Definition AMReX_ParallelReduce.H:221
void copy(HostToDevice, InIter begin, InIter end, OutIter result) noexcept
A host-to-device copy routine. Note this is just a wrapper around memcpy, so it assumes contiguous st...
Definition AMReX_GpuContainers.H:128
static constexpr DeviceToHost deviceToHost
Definition AMReX_GpuContainers.H:106
void streamSynchronize() noexcept
Definition AMReX_GpuDevice.H:310
bool inLaunchRegion() noexcept
Definition AMReX_GpuControl.H:92
bool inNoSyncRegion() noexcept
Definition AMReX_GpuControl.H:152
gpuStream_t gpuStream() noexcept
Definition AMReX_GpuDevice.H:291
MPI_Comm CommunicatorSub() noexcept
sub-communicator for current frame
Definition AMReX_ParallelContext.H:70
Definition AMReX_Amr.cpp:49
__host__ __device__ BoxND< dim > convert(const BoxND< dim > &b, const IntVectND< dim > &typ) noexcept
Return a BoxND with different type.
Definition AMReX_Box.H:1558
void FillRandomNormal(MultiFab &mf, int scomp, int ncomp, Real mean, Real stddev)
Fill MultiFab with random numbers from normal distribution.
Definition AMReX_MultiFabUtil.cpp:1238
int nComp(FabArrayBase const &fa)
Definition AMReX_FabArrayBase.cpp:2851
void FillRandom(MultiFab &mf, int scomp, int ncomp)
Fill MultiFab with random numbers from uniform distribution.
Definition AMReX_MultiFabUtil.cpp:1225
std::enable_if_t< std::is_integral_v< T > > ParallelFor(TypeList< CTOs... > ctos, std::array< int, sizeof...(CTOs)> const &runtime_options, T N, F &&f)
Definition AMReX_CTOParallelForImpl.H:193
void average_down(const MultiFab &S_fine, MultiFab &S_crse, const Geometry &fgeom, const Geometry &cgeom, int scomp, int ncomp, int rr)
Definition AMReX_MultiFabUtil.cpp:336
__host__ __device__ void cast(BaseFab< Tto > &tofab, BaseFab< Tfrom > const &fromfab, Box const &bx, SrcComp scomp, DestComp dcomp, NumComps ncomp) noexcept
Definition AMReX_BaseFabUtility.H:13
std::unique_ptr< MultiFab > get_slice_data(int dir, Real coord, const MultiFab &cc, const Geometry &geom, int start_comp, int ncomp, bool interpolate, RealBox const &bnd_rbx)
Definition AMReX_MultiFabUtil.cpp:574
Vector< typename MF::value_type > get_cell_data(MF const &mf, IntVect const &cell)
Get data in a cell of MultiFab/FabArray.
Definition AMReX_MultiFabUtil.H:1151
ReduceData< Ts... >::Type ParReduce(TypeList< Ops... > operation_list, TypeList< Ts... > type_list, FabArray< FAB > const &fa, IntVect const &nghost, F &&f)
Parallel reduce for MultiFab/FabArray. The reduce result is local and it's the user's responsibility ...
Definition AMReX_ParReduce.H:48
void average_face_to_cellcenter(MultiFab &cc, int dcomp, const Vector< const MultiFab * > &fc, IntVect const &ng_vect)
Definition AMReX_MultiFabUtil.cpp:155
iMultiFab makeFineMask(const BoxArray &cba, const DistributionMapping &cdm, const BoxArray &fba, const IntVect &ratio, int crse_value, int fine_value, MFInfo const &info)
Definition AMReX_MultiFabUtil.cpp:629
BaseFab< T > ReduceToPlane(int direction, Box const &domain, FabArray< FAB > const &mf, F const &f)
Reduce FabArray/MultiFab data to a plane Fab.
Definition AMReX_MultiFabUtil.H:1337
std::pair< FA, FA > ReduceToPlaneMF2(int direction, Box const &domain, FA const &mf, F const &f)
Reduce FabArray/MultiFab data to plane FabArray.
Definition AMReX_MultiFabUtil.H:1428
Vector< MultiFab > convexify(Vector< MultiFab const * > const &mf, Vector< IntVect > const &refinement_ratio)
Convexify AMR data.
Definition AMReX_MultiFabUtil.cpp:1251
FA ReduceToPlaneMF(int direction, Box const &domain, FA const &mf, F const &f)
Reduce FabArray/MultiFab data to plane FabArray.
Definition AMReX_MultiFabUtil.H:1413
std::pair< FA, FA > ReduceToPlaneMF2Patchy(int direction, Box const &domain, FA const &mf, F const &f)
Reduce FabArray/MultiFab data to plane FabArray for patchy BoxArrays.
Definition AMReX_MultiFabUtil.H:1471
void average_down_faces(const Vector< const MF * > &fine, const Vector< MF * > &crse, const IntVect &ratio, int ngcrse=0)
Average fine face-based FabArray onto crse face-based FabArray.
Definition AMReX_MultiFabUtil.H:1033
MultiFab periodicShift(MultiFab const &mf, IntVect const &offset, Periodicity const &period)
Periodic shift MultiFab.
Definition AMReX_MultiFabUtil.cpp:818
Real volumeWeightedSum(Vector< MultiFab const * > const &mf, int icomp, Vector< Geometry > const &geom, Vector< IntVect > const &ratio, bool local)
Volume weighted sum for a vector of MultiFabs.
Definition AMReX_MultiFabUtil.cpp:984
BoxND< 3 > Box
Box is an alias for amrex::BoxND instantiated with AMREX_SPACEDIM.
Definition AMReX_BaseFwd.H:30
bool isMFIterSafe(const FabArrayBase &x, const FabArrayBase &y)
Definition AMReX_MFIter.H:249
MultiFab ToMultiFab(const iMultiFab &imf)
Convert iMultiFab to MultiFab.
Definition AMReX_MultiFabUtil.cpp:564
void writeFabs(const MultiFab &mf, const std::string &name)
Write each fab individually.
Definition AMReX_MultiFabUtil.cpp:551
Real NormHelper(const MultiFab &x, int xcomp, const MultiFab &y, int ycomp, F const &f, int numcomp, IntVect nghost, bool local)
Returns part of a norm based on two MultiFabs.
Definition AMReX_MultiFabUtil.H:836
void computeDivergence(MultiFab &divu, const Array< MultiFab const *, 3 > &umac, const Geometry &geom)
Computes divergence of face-data stored in the umac MultiFab.
Definition AMReX_MultiFabUtil.cpp:717
void average_down_edges(const Vector< const MultiFab * > &fine, const Vector< MultiFab * > &crse, const IntVect &ratio, int ngcrse)
Average fine edge-based MultiFab onto crse edge-based MultiFab.
Definition AMReX_MultiFabUtil.cpp:469
Gpu::HostVector< Real > sumToLine(MultiFab const &mf, int icomp, int ncomp, Box const &domain, int direction, bool local)
Sum MultiFab data to line.
Definition AMReX_MultiFabUtil.cpp:842
IntVectND< 3 > IntVect
IntVect is an alias for amrex::IntVectND instantiated with AMREX_SPACEDIM.
Definition AMReX_BaseFwd.H:33
void print_state(const MultiFab &mf, const IntVect &cell, const int n, const IntVect &ng)
Output state data for a single zone.
Definition AMReX_MultiFabUtil.cpp:546
bool TilingIfNotGPU() noexcept
Definition AMReX_MFIter.H:12
void sum_fine_to_coarse(const MultiFab &S_fine, MultiFab &S_crse, int scomp, int ncomp, const IntVect &ratio, const Geometry &cgeom, const Geometry &)
Definition AMReX_MultiFabUtil.cpp:415
void average_cellcenter_to_face(const Vector< MultiFab * > &fc, const MultiFab &cc, const Geometry &geom, int ncomp, bool use_harmonic_averaging)
Average cell-centered MultiFab onto face-based MultiFab with geometric weighting.
Definition AMReX_MultiFabUtil.cpp:240
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition AMReX.cpp:240
const int[]
Definition AMReX_BLProfiler.cpp:1664
void average_node_to_cellcenter(MultiFab &cc, int dcomp, const MultiFab &nd, int scomp, int ncomp, int ngrow)
Average nodal-based MultiFab onto cell-centered MultiFab.
Definition AMReX_MultiFabUtil.cpp:62
void FourthOrderInterpFromFineToCoarse(MultiFab &cmf, int scomp, int ncomp, MultiFab const &fmf, IntVect const &ratio)
Fourth-order interpolation from fine to coarse level.
Definition AMReX_MultiFabUtil.cpp:1138
FabArray< BaseFab< Long > > ToLongMultiFab(const iMultiFab &imf)
Convert iMultiFab to Long.
Definition AMReX_MultiFabUtil.cpp:569
MF get_line_data(MF const &mf, int dir, IntVect const &cell, Box const &bnd_bx=Box())
Get data in a line of MultiFab/FabArray.
Definition AMReX_MultiFabUtil.H:1182
void computeGradient(MultiFab &grad, const Array< MultiFab const *, 3 > &umac, const Geometry &geom)
Computes gradient of face-data stored in the umac MultiFab.
Definition AMReX_MultiFabUtil.cpp:769
void average_edge_to_cellcenter(MultiFab &cc, int dcomp, const Vector< const MultiFab * > &edge, int ngrow)
Average edge-based MultiFab onto cell-centered MultiFab.
Definition AMReX_MultiFabUtil.cpp:104
void average_down_nodal(const FabArray< FAB > &S_fine, FabArray< FAB > &S_crse, const IntVect &ratio, int ngcrse=0, bool mfiter_is_definitely_safe=false)
Average fine node-based MultiFab onto crse node-centered MultiFab.
Definition AMReX_MultiFabUtil.H:650
A multidimensional array accessor.
Definition AMReX_Array4.H:283
parallel copy or add
Definition AMReX_FabArrayBase.H:538
Definition AMReX_Geometry.H:25
int coord
Definition AMReX_Geometry.H:59
FabArray memory allocation information.
Definition AMReX_FabArray.H:66
Definition AMReX_MFIter.H:20
Struct for holding types.
Definition AMReX_TypeList.H:13