WarpX
ComputeDiagFunctor.H
Go to the documentation of this file.
1 #ifndef WARPX_COMPUTEDIAGFUNCTOR_H_
2 #define WARPX_COMPUTEDIAGFUNCTOR_H_
3 
5 #include "Utils/TextMsg.H"
6 
8 
9 #include <AMReX.H>
10 #include <AMReX_MultiFab.H>
11 
17 {
18 public:
19  ComputeDiagFunctor( int ncomp, amrex::IntVect crse_ratio) :
20  m_ncomp(ncomp), m_crse_ratio(crse_ratio) {}
21  //** Virtual Destructor to handle clean destruction of derived classes */
22  virtual ~ComputeDiagFunctor() = default;
23 
24  // Default move and copy operations
29 
37  virtual void operator() (amrex::MultiFab& mf_dst, int dcomp, int i_buffer = 0) const = 0;
40  [[nodiscard]] int nComp () const { return m_ncomp; }
41 
58  virtual void PrepareFunctorData ( int i_buffer, bool z_slice_in_domain,
59  amrex::Real current_z_boost,
60  amrex::Box buffer_box, const int k_index_zlab,
61  const int snapshot_full) {
62  amrex::ignore_unused(i_buffer,
63  z_slice_in_domain,
64  current_z_boost, buffer_box,
65  k_index_zlab, snapshot_full);
66  }
67  virtual void InitData() {}
68 
70  amrex::MultiFab& mf_dst, const amrex::MultiFab& mf_src, int dcomp,
71  const amrex::DistributionMapping& dm, bool convertRZmodes2cartesian ) const
72  {
73 #ifdef WARPX_DIM_RZ
74  if (convertRZmodes2cartesian) {
75  // In cylindrical geometry, sum real part of all modes of mf_src in
76  // temporary MultiFab mf_dst_stag, and cell-center it to mf_dst
78  nComp()==1,
79  "The RZ averaging over modes must write into one single component");
80  amrex::MultiFab mf_dst_stag( mf_src.boxArray(), dm, 1, mf_src.nGrowVect() );
81  // Mode 0
82  amrex::MultiFab::Copy( mf_dst_stag, mf_src, 0, 0, 1, mf_src.nGrowVect() );
83  for (int ic=1 ; ic < mf_src.nComp() ; ic += 2) {
84  // Real part of all modes > 0
85  amrex::MultiFab::Add( mf_dst_stag, mf_src, ic, 0, 1, mf_src.nGrowVect() );
86  }
87  ablastr::coarsen::sample::Coarsen( mf_dst, mf_dst_stag, dcomp, 0, nComp(), 0, m_crse_ratio );
88  } else {
89  ablastr::coarsen::sample::Coarsen( mf_dst, mf_src, dcomp, 0, nComp(), 0, m_crse_ratio );
90  }
91 #else
92  // In Cartesian geometry, coarsen and interpolate from temporary MultiFab mf_src
93  // to output diagnostic MultiFab mf_dst
94  ablastr::coarsen::sample::Coarsen(mf_dst, mf_src, dcomp, 0, nComp(), mf_dst.nGrowVect(), m_crse_ratio );
95  amrex::ignore_unused(convertRZmodes2cartesian, dm);
96 #endif
97  }
98 
99 private:
101  int m_ncomp;
102 
103 protected:
106 };
107 
108 #endif // WARPX_COMPUTEDIAGFUNCTOR_H_
#define WARPX_ALWAYS_ASSERT_WITH_MESSAGE(EX, MSG)
Definition: TextMsg.H:13
Functor to compute a diagnostic and store the result in existing MultiFab.
Definition: ComputeDiagFunctor.H:17
ComputeDiagFunctor(ComputeDiagFunctor &&)=default
virtual void operator()(amrex::MultiFab &mf_dst, int dcomp, int i_buffer=0) const =0
ComputeDiagFunctor & operator=(ComputeDiagFunctor &&)=default
int nComp() const
Definition: ComputeDiagFunctor.H:40
virtual void PrepareFunctorData(int i_buffer, bool z_slice_in_domain, amrex::Real current_z_boost, amrex::Box buffer_box, const int k_index_zlab, const int snapshot_full)
Prepare data required to process fields in the operator() Note that this function has parameters that...
Definition: ComputeDiagFunctor.H:58
ComputeDiagFunctor & operator=(const ComputeDiagFunctor &)=default
void InterpolateMFForDiag(amrex::MultiFab &mf_dst, const amrex::MultiFab &mf_src, int dcomp, const amrex::DistributionMapping &dm, bool convertRZmodes2cartesian) const
Definition: ComputeDiagFunctor.H:69
amrex::IntVect m_crse_ratio
Definition: ComputeDiagFunctor.H:105
ComputeDiagFunctor(const ComputeDiagFunctor &)=default
virtual void InitData()
Definition: ComputeDiagFunctor.H:67
ComputeDiagFunctor(int ncomp, amrex::IntVect crse_ratio)
Definition: ComputeDiagFunctor.H:19
virtual ~ComputeDiagFunctor()=default
int m_ncomp
Definition: ComputeDiagFunctor.H:101
IntVect nGrowVect() const noexcept
const BoxArray & boxArray() const noexcept
int nComp() const noexcept
static void Add(MultiFab &dst, const MultiFab &src, int srccomp, int dstcomp, int numcomp, int nghost)
static void Copy(MultiFab &dst, const MultiFab &src, int srccomp, int dstcomp, int numcomp, int nghost)
void Coarsen(amrex::MultiFab &mf_dst, const amrex::MultiFab &mf_src, const int dcomp, const int scomp, const int ncomp, const int ngrow, const amrex::IntVect crse_ratio)
Stores in the coarsened MultiFab mf_dst the values obtained by interpolating the data contained in th...
Definition: sample.cpp:83
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void ignore_unused(const Ts &...)