Skip to content

Commit 39a4bdf

Browse files
eebassoeebasso
committed
More changes
1 parent 0e22f78 commit 39a4bdf

File tree

2 files changed

+294
-230
lines changed

2 files changed

+294
-230
lines changed

‎Src/Base/AMReX_MultiFabUtil.H‎

Lines changed: 162 additions & 168 deletions
Original file line numberDiff line numberDiff line change
@@ -55,41 +55,33 @@ namespace amrex
5555
const Geometry& geom,
5656
int ncomp = 1,
5757
bool use_harmonic_averaging = false);
58-
/********** End average-to-different-IndexType methods **********/
5958

6059
/********** Begin average-down methods **********/
6160
/**
62-
* \brief Geometric weighed average of fine FabArray onto coarse FabArray of same IndexType.
61+
* \brief Volume-weighed average of fine cell-based MultiFab onto coarse cell-based MultiFab.
6362
*
64-
* Works for any IndexType shared by both \p crse and \p fine. This routine DOES NOT assume that
65-
* the \p crse BoxArray is a coarsened version of the \p fine BoxArray.
63+
* Both MultiFabs are assumed to be cell-centered. This routine DOES NOT assume that
64+
* the BoxArray of \p S_crse is a coarsened version of the BoxArray of \p S_fine.
6665
*/
67-
template <typename MF, std::enable_if_t<IsFabArray_v<MF>,int> = 0>
68-
void average_down (const MF& S_fine, MF& S_crse,
69-
const Geometry& fgeom, const Geometry& cgeom,
70-
int scomp, int ncomp, const IntVect& ratio,
71-
const IntVect& ngcrse = IntVect(0));
72-
//! Geometric weighed average of fine FabArray onto coarse FabArray of same IndexType.
73-
template <typename MF, std::enable_if_t<IsFabArray_v<MF>,int> = 0>
74-
void average_down (const MF& S_fine, MF& S_crse,
75-
const Geometry& fgeom, const Geometry& cgeom,
76-
int scomp, int ncomp, int rr,
77-
const IntVect& ngcrse = IntVect(0));
66+
void average_down_cells (const MultiFab& S_fine, MultiFab& S_crse,
67+
const Geometry& fgeom, const Geometry& cgeom,
68+
int scomp, int ncomp, const IntVect& ratio);
69+
void average_down_cells (const MultiFab& S_fine, MultiFab& S_crse,
70+
const Geometry& fgeom, const Geometry& cgeom,
71+
int scomp, int ncomp, int rr);
72+
7873
/**
79-
* \brief Average fine FabArray onto coarse FabArray of same IndexType without geometric weighting.
74+
* \brief Average fine cell-based FabArray onto coarse cell-based FabArray without volume weighting.
8075
*
81-
* Works for any IndexType shared by both \p crse and \p fine. This routine DOES NOT assume that
76+
* Work for both cell-centered and nodal MultiFabs. This routine DOES NOT assume that
8277
* the \p crse BoxArray is a coarsened version of the \p fine BoxArray.
8378
*/
84-
template <typename MF, std::enable_if_t<IsFabArray_v<MF>,int> = 0>
85-
void average_down (const MF& S_fine, MF& S_crse,
86-
int scomp, int ncomp, const IntVect& ratio,
87-
const IntVect& ngcrse = IntVect(0));
88-
//! Average fine FabArray onto coarse FabArray of same IndexType without geometric weighting.
89-
template <typename MF, std::enable_if_t<IsFabArray_v<MF>,int> = 0>
90-
void average_down (const MF& S_fine, MF& S_crse,
91-
int scomp, int ncomp, int rr,
92-
const IntVect& ngcrse = IntVect(0));
79+
template<typename FAB>
80+
void average_down_cells (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse,
81+
int scomp, int ncomp, const IntVect& ratio);
82+
template<typename FAB>
83+
void average_down_cells (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse,
84+
int scomp, int ncomp, int rr);
9385

9486
//! Average fine face-based MultiFab onto crse face-based MultiFab.
9587
template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int> = 0>
@@ -164,34 +156,39 @@ namespace amrex
164156
bool mfiter_is_definitely_safe=false);
165157

166158
/**
167-
* \brief Volume-weighed average of fine cell-based MultiFab onto coarse cell-based MultiFab.
159+
* \brief Geometric weighed average of fine FabArray onto coarse FabArray of same IndexType.
168160
*
169-
* Both MultiFabs are assumed to be cell-centered. This routine DOES NOT assume that
170-
* the BoxArray of \p S_crse is a coarsened version of the BoxArray of \p S_fine.
161+
* Works for any IndexType shared by both \p crse and \p fine. This routine DOES NOT assume that
162+
* the \p crse BoxArray is a coarsened version of the \p fine BoxArray.
171163
*/
172-
void average_down_cells (const MultiFab& S_fine, MultiFab& S_crse,
173-
const Geometry& fgeom, const Geometry& cgeom,
174-
int scomp, int ncomp, const IntVect& ratio);
175-
void average_down_cells (const MultiFab& S_fine, MultiFab& S_crse,
176-
const Geometry& fgeom, const Geometry& cgeom,
177-
int scomp, int ncomp, int rr);
178-
164+
void average_down (const MultiFab& S_fine, MultiFab& S_crse,
165+
const Geometry& fgeom, const Geometry& cgeom,
166+
int scomp, int ncomp, const IntVect& ratio,
167+
const IntVect& ngcrse = IntVect(0));
168+
//! Geometric weighed average of fine FabArray onto coarse FabArray of same IndexType.
169+
void average_down (const MultiFab& S_fine, MultiFab& S_crse,
170+
const Geometry& fgeom, const Geometry& cgeom,
171+
int scomp, int ncomp, int rr,
172+
const IntVect& ngcrse = IntVect(0));
179173
/**
180-
* \brief Average fine cell-based FabArray onto coarse cell-based FabArray without volume weighting.
174+
* \brief Average fine FabArray onto coarse FabArray of same IndexType without geometric weighting.
181175
*
182-
* Work for both cell-centered and nodal MultiFabs. This routine DOES NOT assume that
176+
* Works for any IndexType shared by both \p crse and \p fine. This routine DOES NOT assume that
183177
* the \p crse BoxArray is a coarsened version of the \p fine BoxArray.
184178
*/
185179
template<typename FAB>
186-
void average_down_cells (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse,
187-
int scomp, int ncomp, const IntVect& ratio);
180+
void average_down (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse,
181+
int scomp, int ncomp, const IntVect& ratio,
182+
const IntVect& ngcrse = IntVect(0));
183+
//! Average fine FabArray onto coarse FabArray of same IndexType without geometric weighting.
188184
template<typename FAB>
189-
void average_down_cells (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse,
190-
int scomp, int ncomp, int rr);
191-
/* End average-down methods */
185+
void average_down (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse,
186+
int scomp, int ncomp, int rr,
187+
const IntVect& ngcrse = IntVect(0));
188+
/********** End average-down methods **********/
192189

193-
//! Add a coarsened version of the data contained in the \p S_fine MultiFab to
194-
//! \p S_crse, including ghost cells.
190+
//! Add a coarsened version of the data contained in the S_fine MultiFab to
191+
//! S_crse, including ghost cells.
195192
void sum_fine_to_coarse (const MultiFab& S_Fine, MultiFab& S_crse,
196193
int scomp, int ncomp,
197194
const IntVect& ratio,
@@ -499,128 +496,6 @@ makeFineMask (const FabArray<FAB>& cmf, const FabArray<FAB>& fmf,
499496
}
500497

501498
/********** Begin average-down methods **********/
502-
template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int> FOO>
503-
void average_down (const MF& S_fine, MF& S_crse,
504-
int scomp, int ncomp, const IntVect& ratio,
505-
const IntVect& ngcrse)
506-
{
507-
BL_PROFILE("amrex::average_down");
508-
AMREX_ASSERT(S_crse.nComp() == S_fine.nComp());
509-
AMREX_ASSERT(S_crse.ixType() == S_fine.ixType());
510-
511-
// Coarsen() the fine stuff on processors owning the fine data.
512-
BoxArray crse_S_fine_BA = S_fine.boxArray();
513-
if (ngcrse != IntVect{0}) {
514-
crse_S_fine_BA.grow(ngcrse);
515-
}
516-
crse_S_fine_BA.coarsen(ratio);
517-
518-
if (crse_S_fine_BA == S_crse.boxArray() && S_fine.DistributionMap() == S_crse.DistributionMap())
519-
{
520-
average_down_matching(S_crse, S_fine, scomp, scomp, ncomp, ratio, ngcrse);
521-
}
522-
else
523-
{
524-
MF crse_S_fine(crse_S_fine_BA, S_fine.DistributionMap(),
525-
ncomp, ngcrse, MFInfo(), DefaultFabFactory<FAB>());
526-
average_down_matching(S_fine, crse_S_fine, 0, scomp, ncomp, ratio, ngcrse);
527-
S_crse.ParallelCopy(crse_S_fine,0,scomp,ncomp);
528-
}
529-
}
530-
531-
template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int> FOO>
532-
void average_down (const MF& S_fine, MF& S_crse,
533-
int scomp, int ncomp, int rr,
534-
const IntVect& ngcrse)
535-
{
536-
average_down(S_fine,S_crse,scomp,ncomp,rr*IntVect::TheUnitVector(),ngcrse);
537-
}
538-
539-
template <typename MF, std::enable_if_t<IsFabArray<MF>::value,int> FOO>
540-
void average_down_matching (const MF& S_fine, MF& S_crse,
541-
int fcomp, int ccomp, int ncomp, const IntVect& ratio,
542-
const IntVect& ngcrse)
543-
{
544-
using VT = typename MF::value_type;
545-
IndexType ixType = S_fine.ixType();
546-
547-
Dim3 ratioDim3, range;
548-
range_for_average_down_matching(ratioDim3, range, ratio, ixType);
549-
550-
const VT denom = VT(range.x * range.y * range.z);
551-
552-
#ifdef AMREX_USE_GPU
553-
if (Gpu::inLaunchRegion() && S_crse.isFusingCandidate()) {
554-
auto const& crsema = S_crse.arrays();
555-
auto const& finema = S_fine.const_arrays();
556-
ParallelFor(S_crse, ngcrse, ncomp,
557-
[=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k, int n) noexcept
558-
{
559-
amrex_avgdown(i,j,k,n,crsema[box_no],finema[box_no],ccomp,fcomp,ratioDim3,range,denom);
560-
});
561-
if (!Gpu::inNoSyncRegion()) {
562-
Gpu::streamSynchronize();
563-
}
564-
} else
565-
#endif
566-
{
567-
#ifdef AMREX_USE_OMP
568-
#pragma omp parallel if (Gpu::notInLaunchRegion())
569-
#endif
570-
for (MFIter mfi(S_crse,TilingIfNotGPU()); mfi.isValid(); ++mfi)
571-
{
572-
// NOTE: The tilebox is defined at the coarse level.
573-
const Box& bx = mfi.growntilebox(ngcrse);
574-
Array4<VT> const& crsearr = S_crse.array(mfi);
575-
Array4<VT const> const& finearr = S_fine.const_array(mfi);
576-
AMREX_HOST_DEVICE_PARALLEL_FOR_4D(bx, ncomp, i, j, k, n,
577-
{
578-
amrex_avgdown(i,j,k,n,crsearr,finearr,ccomp,fcomp,ratioDim3,range,denom);
579-
});
580-
}
581-
}
582-
}
583-
584-
void range_for_average_down_matching (Dim3& ratioDim3, Dim3& range, const IntVect& ratio, const IndexType ixType)
585-
{
586-
ratioDim3 = Dim3{0,0,0};
587-
range = Dim3{1,1,1};
588-
589-
ratioDim3.x = ratio[0];
590-
if (ixType.cellCentered(0)) {
591-
range.x = ratio[0];
592-
}
593-
#if AMREX_SPACEDIM >= 2
594-
ratioDim3.y = ratio[1];
595-
if (ixType.cellCentered(1)) {
596-
range.y = ratio[1];
597-
}
598-
#endif
599-
#if AMREX_SPACEDIM >= 3
600-
ratioDim3.z = ratio[2];
601-
if (ixType.cellCentered(2)) {
602-
range.z = ratio[2];
603-
}
604-
#endif
605-
}
606-
607-
template <typename T>
608-
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
609-
void amrex_avgdown (int i, int j, int k, int n,
610-
const Array4<T> &crse, const Array4<const T> &fine, int ccomp, int fcomp,
611-
const Dim3 &ratio, const Dim3 &range, const T &denom) noexcept
612-
{
613-
T c = T(0.0);
614-
for ( int kk = k*ratio.z; kk < k*ratio.z + range.z; ++kk) {
615-
for ( int jj = j*ratio.y; jj < j*ratio.y + range.y; ++jj) {
616-
for (int ii = i*ratio.x; ii < i*ratio.x + range.x; ++ii) {
617-
c += fine(ii,jj,kk,n+fcomp);
618-
}
619-
}
620-
}
621-
crse(i,j,k,n+ccomp) = c / denom;
622-
}
623-
624499
template <typename FAB>
625500
void average_down_nodal (const FabArray<FAB>& fine, FabArray<FAB>& crse,
626501
const IntVect& ratio, int ngcrse, bool mfiter_is_definitely_safe)
@@ -666,7 +541,7 @@ void average_down_cells (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse, int
666541

667542
template<typename FAB>
668543
void average_down_cells (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse,
669-
int scomp, int ncomp, const IntVect& ratio)
544+
int scomp, int ncomp, const IntVect& ratio)
670545
{
671546
BL_PROFILE("amrex::average_down_cells");
672547
AMREX_ASSERT(S_crse.nComp() == S_fine.nComp());
@@ -790,6 +665,125 @@ void average_down_cells (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse,
790665
}
791666
}
792667

668+
template <typename FAB>
669+
void average_down (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse,
670+
int scomp, int ncomp, int rr,
671+
const IntVect& ngcrse)
672+
{
673+
average_down(S_fine,S_crse,scomp,ncomp,rr*IntVect::TheUnitVector(),ngcrse);
674+
}
675+
676+
template <typename FAB>
677+
void average_down (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse,
678+
int scomp, int ncomp, const IntVect& ratio,
679+
const IntVect& ngcrse)
680+
{
681+
BL_PROFILE("amrex::average_down");
682+
AMREX_ASSERT(S_crse.nComp() == S_fine.nComp());
683+
AMREX_ASSERT(S_crse.ixType() == S_fine.ixType());
684+
685+
// Coarsen() the fine stuff on processors owning the fine data.
686+
BoxArray crse_S_fine_BA = S_fine.boxArray();
687+
if (ngcrse != IntVect{0}) {
688+
crse_S_fine_BA.grow(ngcrse);
689+
}
690+
crse_S_fine_BA.coarsen(ratio);
691+
692+
if (crse_S_fine_BA == S_crse.boxArray() && S_fine.DistributionMap() == S_crse.DistributionMap())
693+
{
694+
average_down_matching(S_crse, S_fine, scomp, scomp, ncomp, ratio, ngcrse);
695+
}
696+
else
697+
{
698+
FabArray<FAB> crse_S_fine(crse_S_fine_BA, S_fine.DistributionMap(),
699+
ncomp, ngcrse, MFInfo(), DefaultFabFactory<FAB>());
700+
average_down_matching(S_fine, crse_S_fine, 0, scomp, ncomp, ratio, ngcrse);
701+
S_crse.ParallelCopy(crse_S_fine,0,scomp,ncomp);
702+
}
703+
}
704+
705+
template <typename FAB>
706+
void average_down_matching (const FabArray<FAB>& S_fine, FabArray<FAB>& S_crse,
707+
int fcomp, int ccomp, int ncomp, const IntVect& ratio,
708+
const IntVect& ngcrse)
709+
{
710+
using VT = typename MF::value_type;
711+
IndexType ixType = S_fine.ixType();
712+
713+
Dim3 ratioDim3, range;
714+
range_for_average_down_matching(ratioDim3, range, ratio, ixType);
715+
716+
const VT denom = VT(range.x * range.y * range.z);
717+
718+
#ifdef AMREX_USE_GPU
719+
if (Gpu::inLaunchRegion() && S_crse.isFusingCandidate()) {
720+
auto const& crsema = S_crse.arrays();
721+
auto const& finema = S_fine.const_arrays();
722+
ParallelFor(S_crse, ngcrse, 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],ccomp,fcomp,ratioDim3,range,denom);
726+
});
727+
if (!Gpu::inNoSyncRegion()) {
728+
Gpu::streamSynchronize();
729+
}
730+
} else
731+
#endif
732+
{
733+
#ifdef AMREX_USE_OMP
734+
#pragma omp parallel if (Gpu::notInLaunchRegion())
735+
#endif
736+
for (MFIter mfi(S_crse,TilingIfNotGPU()); mfi.isValid(); ++mfi)
737+
{
738+
// NOTE: The tilebox is defined at the coarse level.
739+
const Box& bx = mfi.growntilebox(ngcrse);
740+
Array4<VT> const& crsearr = S_crse.array(mfi);
741+
Array4<VT const> const& finearr = S_fine.const_array(mfi);
742+
AMREX_HOST_DEVICE_PARALLEL_FOR_4D(bx, ncomp, i, j, k, n,
743+
{
744+
amrex_avgdown(i,j,k,n,crsearr,finearr,ccomp,fcomp,ratioDim3,range,denom);
745+
});
746+
}
747+
}
748+
}
749+
750+
void range_for_average_down_matching (Dim3& ratioDim3, Dim3& range,
751+
const IntVect& ratio, const IndexType ixType)
752+
{
753+
ratioDim3 = Dim3{0,0,0};
754+
range = Dim3{1,1,1};
755+
756+
AMREX_D_TERM(
757+
ratioDim3.x = ratio[0];
758+
if (ixType.cellCentered(0)) {
759+
range.x = ratio[0];
760+
},
761+
ratioDim3.y = ratio[1];
762+
if (ixType.cellCentered(1)) {
763+
range.y = ratio[1];
764+
},
765+
ratioDim3.z = ratio[2];
766+
if (ixType.cellCentered(2)) {
767+
range.z = ratio[2];
768+
});
769+
}
770+
771+
template <typename T>
772+
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
773+
void amrex_avgdown (int i, int j, int k, int n,
774+
const Array4<T> &crse, const Array4<const T> &fine, int ccomp, int fcomp,
775+
const Dim3 &ratio, const Dim3 &range, const T &denom) noexcept
776+
{
777+
T c = T(0.0);
778+
for ( int kk = k*ratio.z; kk < k*ratio.z + range.z; ++kk) {
779+
for ( int jj = j*ratio.y; jj < j*ratio.y + range.y; ++jj) {
780+
for (int ii = i*ratio.x; ii < i*ratio.x + range.x; ++ii) {
781+
c += fine(ii,jj,kk,n+fcomp);
782+
}
783+
}
784+
}
785+
crse(i,j,k,n+ccomp) = c/denom;
786+
}
793787

794788

795789

0 commit comments

Comments
 (0)