WarpX
SpectralFieldDataRZ.H
Go to the documentation of this file.
1 /* Copyright 2019-2020 David Grote
2  *
3  * This file is part of WarpX.
4  *
5  * License: BSD-3-Clause-LBNL
6  */
7 #ifndef WARPX_SPECTRAL_FIELD_DATA_RZ_H_
8 #define WARPX_SPECTRAL_FIELD_DATA_RZ_H_
9 
10 #include "SpectralBinomialFilter.H"
11 #include "SpectralFieldData.H"
13 #include "SpectralKSpaceRZ.H"
14 
16 
17 #include <AMReX_MultiFab.H>
18 
19 /* \brief Class that stores the fields in spectral space, and performs the
20  * Fourier transforms between real space and spectral space
21  */
23 {
24 
25  public:
26 
27  // Define the FFTplans type, which holds one fft plan per box
28  // (plans are only initialized for the boxes that are owned by
29  // the local MPI rank)
31 
32  // Similarly, define the Hankel transformers and filter for each box.
34 
36 
37  SpectralFieldDataRZ (int lev,
38  const amrex::BoxArray& realspace_ba,
39  const SpectralKSpaceRZ& k_space,
41  int n_field_required,
42  int n_modes);
43  SpectralFieldDataRZ () = default; // Default constructor
45 
54 
55 
56  void ForwardTransform (int lev, const amrex::MultiFab& mf, int field_index,
57  int i_comp=0);
58  void ForwardTransform (int lev, const amrex::MultiFab& mf_r, int field_index_r,
59  const amrex::MultiFab& mf_t, int field_index_t);
60  void BackwardTransform (int lev, amrex::MultiFab& mf, int field_index,
61  int i_comp=0);
62  void BackwardTransform (int lev, amrex::MultiFab& mf_r, int field_index_r,
63  amrex::MultiFab& mf_t, int field_index_t);
64 
65  void FABZForwardTransform (amrex::MFIter const & mfi, amrex::Box const & realspace_bx,
66  amrex::MultiFab const & tempHTransformedSplit,
67  int field_index, bool is_nodal_z);
68  void FABZBackwardTransform (amrex::MFIter const & mfi, amrex::Box const & realspace_bx,
69  int field_index,
70  amrex::MultiFab & tempHTransformedSplit,
71  bool is_nodal_z);
72 
80  void CopySpectralDataComp (int src_comp, int dest_comp)
81  {
82  // In spectral space fields of each mode are grouped together, so that the index
83  // of a field for a specific mode is given by field_index + mode*n_fields.
84  // That’s why, looping over all azimuthal modes, we need to take into account corresponding
85  // shift, given by imode * m_n_fields.
86  for (int imode=0 ; imode < n_rz_azimuthal_modes ; imode++)
87  {
88  const int shift_comp = imode * m_n_fields;
89  // The last two arguments represent the number of components and
90  // the number of ghost cells to perform this operation
91  Copy(fields, fields, src_comp + shift_comp, dest_comp + shift_comp, 1, 0);
92  }
93  }
94 
100  void ZeroOutDataComp(int icomp)
101  {
102  // In spectral space fields of each mode are grouped together, so that the index
103  // of a field for a specific mode is given by field_index + mode*n_fields.
104  // That’s why, looping over all azimuthal modes, we need to take into account corresponding
105  // shift, given by imode * m_n_fields.
106  for (int imode=0 ; imode < n_rz_azimuthal_modes ; imode++)
107  {
108  const int shift_comp = imode * m_n_fields;
109  // The last argument represents the number of components to perform this operation
110  fields.setVal(0., icomp + shift_comp, 1);
111  }
112  }
113 
120  void ScaleDataComp(int icomp, amrex::Real scale_factor)
121  {
122  // In spectral space fields of each mode are grouped together, so that the index
123  // of a field for a specific mode is given by field_index + mode*n_fields.
124  // That’s why, looping over all azimuthal modes, we need to take into account corresponding
125  // shift, given by imode * m_n_fields.
126  for (int imode=0 ; imode < n_rz_azimuthal_modes ; imode++)
127  {
128  const int shift_comp = imode * m_n_fields;
129  // The last argument represents the number of components to perform this operation
130  fields.mult(scale_factor, icomp + shift_comp, 1);
131  }
132  }
133 
134  void InitFilter (amrex::IntVect const & filter_npass_each_dir, bool compensation,
135  SpectralKSpaceRZ const & k_space);
136 
137  void ApplyFilter (int lev, int field_index);
138  void ApplyFilter (int lev, int field_index1,
139  int field_index2, int field_index3);
140 
141  // Returns an array that holds the kr for all of the modes
143  return multi_spectral_hankel_transformer[mfi].getKrArray();
144  }
145 
148 
152  int m_ncomps;
153 
154  private:
155 
158 
159  // tempHTransformed and tmpSpectralField store fields
160  // right before/after the z Fourier transform
161  SpectralField tempHTransformed; // contains Complexes
162  SpectralField tmpSpectralField; // contains Complexes
164  // Correcting "shift" factors when performing FFT from/to
165  // a cell-centered grid in real space, instead of a nodal grid
169 
170 };
171 
172 #endif // WARPX_SPECTRAL_FIELD_DATA_RZ_H_
Definition: SpectralFieldDataRZ.H:23
FFTplans backward_plan
Definition: SpectralFieldDataRZ.H:163
void ZeroOutDataComp(int icomp)
Set to zero the data on component icomp of fields.
Definition: SpectralFieldDataRZ.H:100
int m_n_fields
Definition: SpectralFieldDataRZ.H:157
void CopySpectralDataComp(int src_comp, int dest_comp)
Copy spectral data from component src_comp to component dest_comp of fields.
Definition: SpectralFieldDataRZ.H:80
MultiSpectralHankelTransformer multi_spectral_hankel_transformer
Definition: SpectralFieldDataRZ.H:167
SpectralFieldDataRZ()=default
SpectralFieldIndex m_spectral_index
Definition: SpectralFieldDataRZ.H:156
void BackwardTransform(int lev, amrex::MultiFab &mf, int field_index, int i_comp=0)
Definition: SpectralFieldDataRZ.cpp:576
SpectralFieldDataRZ & operator=(SpectralFieldDataRZ const &)=delete
SpectralFieldDataRZ & operator=(SpectralFieldDataRZ &&field_data)=default
SpectralField tempHTransformed
Definition: SpectralFieldDataRZ.H:161
void InitFilter(amrex::IntVect const &filter_npass_each_dir, bool compensation, SpectralKSpaceRZ const &k_space)
Definition: SpectralFieldDataRZ.cpp:756
HankelTransform::RealVector const & getKrArray(amrex::MFIter const &mfi) const
Definition: SpectralFieldDataRZ.H:142
void FABZBackwardTransform(amrex::MFIter const &mfi, amrex::Box const &realspace_bx, int field_index, amrex::MultiFab &tempHTransformedSplit, bool is_nodal_z)
Definition: SpectralFieldDataRZ.cpp:349
SpectralFieldDataRZ(SpectralFieldDataRZ const &)=delete
void ForwardTransform(int lev, const amrex::MultiFab &mf, int field_index, int i_comp=0)
Definition: SpectralFieldDataRZ.cpp:451
FFTplans forward_plan
Definition: SpectralFieldDataRZ.H:163
void FABZForwardTransform(amrex::MFIter const &mfi, amrex::Box const &realspace_bx, amrex::MultiFab const &tempHTransformedSplit, int field_index, bool is_nodal_z)
Definition: SpectralFieldDataRZ.cpp:241
void ScaleDataComp(int icomp, amrex::Real scale_factor)
Scale the data on component icomp of fields by a given scale factor.
Definition: SpectralFieldDataRZ.H:120
int m_ncomps
Number of MultiFab components, see WarpX::ncomps.
Definition: SpectralFieldDataRZ.H:152
SpectralField fields
fields stores fields in spectral space, as multicomponent FabArray
Definition: SpectralFieldDataRZ.H:147
int n_rz_azimuthal_modes
Number of modes for the RZ multi-mode version, see WarpX::n_rz_azimuthal_modes.
Definition: SpectralFieldDataRZ.H:150
SpectralShiftFactor zshift_FFTfromCell
Definition: SpectralFieldDataRZ.H:166
SpectralFieldDataRZ(SpectralFieldDataRZ &&)=default
~SpectralFieldDataRZ()
Definition: SpectralFieldDataRZ.cpp:208
SpectralField tmpSpectralField
Definition: SpectralFieldDataRZ.H:162
SpectralShiftFactor zshift_FFTtoCell
Definition: SpectralFieldDataRZ.H:166
void ApplyFilter(int lev, int field_index)
Definition: SpectralFieldDataRZ.cpp:773
BinomialFilter binomialfilter
Definition: SpectralFieldDataRZ.H:168
Definition: SpectralFieldData.H:34
Definition: SpectralKSpaceRZ.H:21
void setVal(value_type val)
void mult(value_type val, int comp, int num_comp, int nghost=0)
void Copy(FabArray< DFAB > &dst, FabArray< SFAB > const &src, int srccomp, int dstcomp, int numcomp, int nghost)