Block-Structured AMR Software Framework
AMReX_FFT_R2C.H
Go to the documentation of this file.
1 #ifndef AMREX_FFT_R2C_H_
2 #define AMREX_FFT_R2C_H_
3 #include <AMReX_Config.H>
4 
5 #include <AMReX_MultiFab.H>
6 #include <AMReX_FFT_Helper.H>
7 #include <algorithm>
8 #include <numeric>
9 #include <tuple>
10 
11 namespace amrex::FFT
12 {
13 
14 template <typename T> class OpenBCSolver;
15 
30 template <typename T = Real, FFT::Direction D = FFT::Direction::both,
32  // Don't change the default. Otherwise OpenBCSolver might break.
33 class R2C
34 {
35 public:
36  using MF = std::conditional_t<std::is_same_v<T,Real>,
39 
40  template <typename U> friend class OpenBCSolver;
41 
48  explicit R2C (Box const& domain, Info const& info = Info{});
49 
50  ~R2C ();
51 
52  R2C (R2C const&) = delete;
53  R2C (R2C &&) = delete;
54  R2C& operator= (R2C const&) = delete;
55  R2C& operator= (R2C &&) = delete;
56 
76  template <typename F, Direction DIR=D,
77  std::enable_if_t<DIR == Direction::both, int> = 0>
78  void forwardThenBackward (MF const& inmf, MF& outmf, F const& post_forward)
79  {
80  BL_PROFILE("FFT::R2C::forwardbackward");
81  this->forward(inmf);
82  this->post_forward_doit(post_forward);
83  this->backward(outmf);
84  }
85 
95  template <Direction DIR=D, std::enable_if_t<DIR == Direction::forward ||
96  DIR == Direction::both, int> = 0>
97  void forward (MF const& inmf);
98 
108  template <Direction DIR=D, std::enable_if_t<DIR == Direction::forward ||
109  DIR == Direction::both, int> = 0>
110  void forward (MF const& inmf, cMF& outmf);
111 
120  template <Direction DIR=D, std::enable_if_t<DIR == Direction::both, int> = 0>
121  void backward (MF& outmf);
122 
132  template <Direction DIR=D, std::enable_if_t<DIR == Direction::backward ||
133  DIR == Direction::both, int> = 0>
134  void backward (cMF const& inmf, MF& outmf);
135 
139  [[nodiscard]] T scalingFactor () const;
140 
150  template <Direction DIR=D, std::enable_if_t<DIR == Direction::forward ||
151  DIR == Direction::both, int> = 0>
152  std::pair<cMF*,IntVect> getSpectralData ();
153 
161  [[nodiscard]] std::pair<BoxArray,DistributionMapping> getSpectralDataLayout () const;
162 
163  // public for cuda
164  template <typename F>
165  void post_forward_doit (F const& post_forward);
166 
167  void prepare_openbc ();
168 
169  void backward_doit (MF& outmf, IntVect const& ngout = IntVect(0));
170 
171 private:
172 
173  static std::pair<Plan<T>,Plan<T>> make_c2c_plans (cMF& inout);
174 
183 
184  // Comm meta-data. In the forward phase, we start with (x,y,z),
185  // transpose to (y,x,z) and then (z,x,y). In the backward phase, we
186  // perform inverse transpose.
187  std::unique_ptr<MultiBlockCommMetaData> m_cmd_x2y; // (x,y,z) -> (y,x,z)
188  std::unique_ptr<MultiBlockCommMetaData> m_cmd_y2x; // (y,x,z) -> (x,y,z)
189  std::unique_ptr<MultiBlockCommMetaData> m_cmd_y2z; // (y,x,z) -> (z,x,y)
190  std::unique_ptr<MultiBlockCommMetaData> m_cmd_z2y; // (z,x,y) -> (y,x,z)
191  std::unique_ptr<MultiBlockCommMetaData> m_cmd_x2z; // (x,y,z) -> (z,x,y)
192  std::unique_ptr<MultiBlockCommMetaData> m_cmd_z2x; // (z,x,y) -> (x,y,z)
193  std::unique_ptr<MultiBlockCommMetaData> m_cmd_x2z_half; // for openbc
194  std::unique_ptr<MultiBlockCommMetaData> m_cmd_z2x_half; // for openbc
201 
206 
207  std::unique_ptr<char,DataDeleter> m_data_1;
208  std::unique_ptr<char,DataDeleter> m_data_2;
209 
214 
216 
217  bool m_do_alld_fft = false;
218  bool m_slab_decomp = false;
219  bool m_openbc_half = false;
220 };
221 
222 template <typename T, Direction D, DomainStrategy S>
223 R2C<T,D,S>::R2C (Box const& domain, Info const& info)
224  : m_real_domain(domain),
225  m_spectral_domain_x(IntVect(0), IntVect(AMREX_D_DECL(domain.length(0)/2,
226  domain.length(1)-1,
227  domain.length(2)-1)),
228  domain.ixType()),
229 #if (AMREX_SPACEDIM >= 2)
230  m_spectral_domain_y(IntVect(0), IntVect(AMREX_D_DECL(domain.length(1)-1,
231  domain.length(0)/2,
232  domain.length(2)-1)),
233  domain.ixType()),
234 #if (AMREX_SPACEDIM == 3)
235  m_spectral_domain_z(IntVect(0), IntVect(AMREX_D_DECL(domain.length(2)-1,
236  domain.length(0)/2,
237  domain.length(1)-1)),
238  domain.ixType()),
239 #endif
240 #endif
241  m_info(info)
242 {
243  BL_PROFILE("FFT::R2C");
244 
245  static_assert(std::is_same_v<float,T> || std::is_same_v<double,T>);
247 #if (AMREX_SPACEDIM == 3)
250 #else
252 #endif
253 
254  int myproc = ParallelContext::MyProcSub();
256 
257 #if (AMREX_SPACEDIM == 3)
258  if (S == DomainStrategy::slab && (m_real_domain.length(1) > 1)) {
259  m_slab_decomp = true;
260  }
261 #endif
262 
263  auto bax = amrex::decompose(m_real_domain, nprocs,
264  {AMREX_D_DECL(false,!m_slab_decomp,true)}, true);
266  m_rx.define(bax, dmx, 1, 0, MFInfo().SetAlloc(false));
267 
268  {
269  BoxList bl = bax.boxList();
270  for (auto & b : bl) {
271  b.shift(-m_real_domain.smallEnd());
272  b.setBig(0, m_spectral_domain_x.bigEnd(0));
273  }
274  BoxArray cbax(std::move(bl));
275  m_cx.define(cbax, dmx, 1, 0, MFInfo().SetAlloc(false));
276  }
277 
279 
280  if (!m_do_alld_fft) // do a series of 1d or 2d ffts
281  {
282  //
283  // make data containers
284  //
285 
286 #if (AMREX_SPACEDIM >= 2)
287  DistributionMapping cdmy;
288  if ((m_real_domain.length(1) > 1) && !m_slab_decomp) {
289  auto cbay = amrex::decompose(m_spectral_domain_y, nprocs,
290  {AMREX_D_DECL(false,true,true)}, true);
291  if (cbay.size() == dmx.size()) {
292  cdmy = dmx;
293  } else {
294  cdmy = detail::make_iota_distromap(cbay.size());
295  }
296  m_cy.define(cbay, cdmy, 1, 0, MFInfo().SetAlloc(false));
297  }
298 #endif
299 
300 #if (AMREX_SPACEDIM == 3)
301  if (m_real_domain.length(1) > 1 &&
302  (! m_info.batch_mode && m_real_domain.length(2) > 1))
303  {
304  auto cbaz = amrex::decompose(m_spectral_domain_z, nprocs,
305  {false,true,true}, true);
306  DistributionMapping cdmz;
307  if (cbaz.size() == dmx.size()) {
308  cdmz = dmx;
309  } else if (cbaz.size() == cdmy.size()) {
310  cdmz = cdmy;
311  } else {
312  cdmz = detail::make_iota_distromap(cbaz.size());
313  }
314  m_cz.define(cbaz, cdmz, 1, 0, MFInfo().SetAlloc(false));
315  }
316 #endif
317 
318  if (m_slab_decomp) {
321  } else {
324  }
325 
326  //
327  // make copiers
328  //
329 
330 #if (AMREX_SPACEDIM >= 2)
331  if (! m_cy.empty()) {
332  // comm meta-data between x and y phases
333  m_cmd_x2y = std::make_unique<MultiBlockCommMetaData>
335  m_cmd_y2x = std::make_unique<MultiBlockCommMetaData>
337  }
338 #endif
339 #if (AMREX_SPACEDIM == 3)
340  if (! m_cz.empty() ) {
341  if (m_slab_decomp) {
342  // comm meta-data between xy and z phases
343  m_cmd_x2z = std::make_unique<MultiBlockCommMetaData>
345  m_cmd_z2x = std::make_unique<MultiBlockCommMetaData>
347  } else {
348  // comm meta-data between y and z phases
349  m_cmd_y2z = std::make_unique<MultiBlockCommMetaData>
351  m_cmd_z2y = std::make_unique<MultiBlockCommMetaData>
353  }
354  }
355 #endif
356 
357  //
358  // make plans
359  //
360 
361  if (myproc < m_rx.size())
362  {
363  Box const& box = m_rx.box(myproc);
364  auto* pr = m_rx[myproc].dataPtr();
365  auto* pc = (typename Plan<T>::VendorComplex *)m_cx[myproc].dataPtr();
366 #ifdef AMREX_USE_SYCL
367  m_fft_fwd_x.template init_r2c<Direction::forward>(box, pr, pc, m_slab_decomp);
369 #else
370  if constexpr (D == Direction::both || D == Direction::forward) {
371  m_fft_fwd_x.template init_r2c<Direction::forward>(box, pr, pc, m_slab_decomp);
372  }
373  if constexpr (D == Direction::both || D == Direction::backward) {
374  m_fft_bwd_x.template init_r2c<Direction::backward>(box, pr, pc, m_slab_decomp);
375  }
376 #endif
377  }
378 
379 #if (AMREX_SPACEDIM >= 2)
380  if (! m_cy.empty()) {
382  }
383 #endif
384 #if (AMREX_SPACEDIM == 3)
385  if (! m_cz.empty()) {
387  }
388 #endif
389  }
390  else // do fft in all dimensions at the same time
391  {
394 
395  auto const& len = m_real_domain.length();
396  auto* pr = (void*)m_rx[0].dataPtr();
397  auto* pc = (void*)m_cx[0].dataPtr();
398 #ifdef AMREX_USE_SYCL
399  m_fft_fwd_x.template init_r2c<Direction::forward>(len, pr, pc, false);
401 #else
402  if constexpr (D == Direction::both || D == Direction::forward) {
403  m_fft_fwd_x.template init_r2c<Direction::forward>(len, pr, pc, false);
404  }
405  if constexpr (D == Direction::both || D == Direction::backward) {
406  m_fft_bwd_x.template init_r2c<Direction::backward>(len, pr, pc, false);
407  }
408 #endif
409  }
410 }
411 
412 template <typename T, Direction D, DomainStrategy S>
414 {
415  if (m_fft_bwd_x.plan != m_fft_fwd_x.plan) {
416  m_fft_bwd_x.destroy();
417  }
418  if (m_fft_bwd_y.plan != m_fft_fwd_y.plan) {
419  m_fft_bwd_y.destroy();
420  }
421  if (m_fft_bwd_z.plan != m_fft_fwd_z.plan) {
422  m_fft_bwd_z.destroy();
423  }
424  m_fft_fwd_x.destroy();
425  m_fft_fwd_y.destroy();
426  m_fft_fwd_z.destroy();
427  if (m_fft_bwd_x_half.plan != m_fft_fwd_x_half.plan) {
428  m_fft_bwd_x_half.destroy();
429  }
430  m_fft_fwd_x_half.destroy();
431 }
432 
433 template <typename T, Direction D, DomainStrategy S>
435 {
436 #if (AMREX_SPACEDIM == 3)
437  if (m_do_alld_fft) { return; }
438 
439  if (m_slab_decomp) {
440  auto* fab = detail::get_fab(m_rx);
441  if (fab) {
442  Box bottom_half = m_real_domain;
443  bottom_half.growHi(2,-m_real_domain.length(2)/2);
444  Box box = fab->box() & bottom_half;
445  if (box.ok()) {
446  auto* pr = fab->dataPtr();
447  auto* pc = (typename Plan<T>::VendorComplex *)
448  detail::get_fab(m_cx)->dataPtr();
449 #ifdef AMREX_USE_SYCL
450  m_fft_fwd_x_half.template init_r2c<Direction::forward>
451  (box, pr, pc, m_slab_decomp);
452  m_fft_bwd_x_half = m_fft_fwd_x_half;
453 #else
454  if constexpr (D == Direction::both || D == Direction::forward) {
455  m_fft_fwd_x_half.template init_r2c<Direction::forward>
456  (box, pr, pc, m_slab_decomp);
457  }
458  if constexpr (D == Direction::both || D == Direction::backward) {
459  m_fft_bwd_x_half.template init_r2c<Direction::backward>
460  (box, pr, pc, m_slab_decomp);
461  }
462 #endif
463  }
464  }
465  } // else todo
466 
467  if (m_cmd_x2z && ! m_cmd_x2z_half) {
468  Box bottom_half = m_spectral_domain_z;
469  // Note that z-direction's index is 0 because we z is the
470  // unit-stride direction here.
471  bottom_half.growHi(0,-m_spectral_domain_z.length(0)/2);
472  m_cmd_x2z_half = std::make_unique<MultiBlockCommMetaData>
473  (m_cz, bottom_half, m_cx, IntVect(0), m_dtos_x2z);
474  }
475 
476  if (m_cmd_z2x && ! m_cmd_z2x_half) {
477  Box bottom_half = m_spectral_domain_x;
478  bottom_half.growHi(2,-m_spectral_domain_x.length(2)/2);
479  m_cmd_z2x_half = std::make_unique<MultiBlockCommMetaData>
480  (m_cx, bottom_half, m_cz, IntVect(0), m_dtos_z2x);
481  }
482 #endif
483 }
484 
485 template <typename T, Direction D, DomainStrategy S>
486 template <Direction DIR, std::enable_if_t<DIR == Direction::forward ||
487  DIR == Direction::both, int> >
488 void R2C<T,D,S>::forward (MF const& inmf)
489 {
490  BL_PROFILE("FFT::R2C::forward(in)");
491 
492  if (&m_rx != &inmf) {
493  m_rx.ParallelCopy(inmf, 0, 0, 1);
494  }
495 
496  if (m_do_alld_fft) {
497  m_fft_fwd_x.template compute_r2c<Direction::forward>();
498  return;
499  }
500 
501  auto& fft_x = m_openbc_half ? m_fft_fwd_x_half : m_fft_fwd_x;
502  fft_x.template compute_r2c<Direction::forward>();
503 
504  if ( m_cmd_x2y) {
505  ParallelCopy(m_cy, m_cx, *m_cmd_x2y, 0, 0, 1, m_dtos_x2y);
506  }
507  m_fft_fwd_y.template compute_c2c<Direction::forward>();
508 
509  if ( m_cmd_y2z) {
510  ParallelCopy(m_cz, m_cy, *m_cmd_y2z, 0, 0, 1, m_dtos_y2z);
511  }
512 #if (AMREX_SPACEDIM == 3)
513  else if ( m_cmd_x2z) {
514  if (m_openbc_half) {
516  {NonLocalBC::PackComponents{}, m_dtos_x2z};
517  auto handler = ParallelCopy_nowait(m_cz, m_cx, *m_cmd_x2z_half, packing);
518 
519  Box upper_half = m_spectral_domain_z;
520  // Note that z-direction's index is 0 because we z is the
521  // unit-stride direction here.
522  upper_half.growLo (0,-m_spectral_domain_z.length(0)/2);
523  m_cz.setVal(0, upper_half, 0, 1);
524 
525  ParallelCopy_finish(m_cz, std::move(handler), *m_cmd_x2z_half, packing);
526  } else {
527  ParallelCopy(m_cz, m_cx, *m_cmd_x2z, 0, 0, 1, m_dtos_x2z);
528  }
529  }
530 #endif
531  m_fft_fwd_z.template compute_c2c<Direction::forward>();
532 }
533 
534 template <typename T, Direction D, DomainStrategy S>
535 template <Direction DIR, std::enable_if_t<DIR == Direction::both, int> >
537 {
538  backward_doit(outmf);
539 }
540 
541 template <typename T, Direction D, DomainStrategy S>
542 void R2C<T,D,S>::backward_doit (MF& outmf, IntVect const& ngout)
543 {
544  BL_PROFILE("FFT::R2C::backward(out)");
545 
546  if (m_do_alld_fft) {
547  m_fft_bwd_x.template compute_r2c<Direction::backward>();
548  outmf.ParallelCopy(m_rx, 0, 0, 1, IntVect(0), ngout);
549  return;
550  }
551 
552  m_fft_bwd_z.template compute_c2c<Direction::backward>();
553  if ( m_cmd_z2y) {
554  ParallelCopy(m_cy, m_cz, *m_cmd_z2y, 0, 0, 1, m_dtos_z2y);
555  }
556 #if (AMREX_SPACEDIM == 3)
557  else if ( m_cmd_z2x) {
558  auto const& cmd = m_openbc_half ? m_cmd_z2x_half : m_cmd_z2x;
559  ParallelCopy(m_cx, m_cz, *cmd, 0, 0, 1, m_dtos_z2x);
560  }
561 #endif
562 
563  m_fft_bwd_y.template compute_c2c<Direction::backward>();
564  if ( m_cmd_y2x) {
565  ParallelCopy(m_cx, m_cy, *m_cmd_y2x, 0, 0, 1, m_dtos_y2x);
566  }
567 
568  auto& fft_x = m_openbc_half ? m_fft_bwd_x_half : m_fft_bwd_x;
569  fft_x.template compute_r2c<Direction::backward>();
570  outmf.ParallelCopy(m_rx, 0, 0, 1, IntVect(0), ngout);
571 }
572 
573 template <typename T, Direction D, DomainStrategy S>
574 std::pair<Plan<T>, Plan<T>>
576 {
577  Plan<T> fwd;
578  Plan<T> bwd;
579 
580  auto* fab = detail::get_fab(inout);
581  if (!fab) { return {fwd, bwd};}
582 
583  Box const& box = fab->box();
584  auto* pio = (typename Plan<T>::VendorComplex *)fab->dataPtr();
585 
586 #ifdef AMREX_USE_SYCL
587  fwd.template init_c2c<Direction::forward>(box, pio);
588  bwd = fwd;
589 #else
590  if constexpr (D == Direction::both || D == Direction::forward) {
591  fwd.template init_c2c<Direction::forward>(box, pio);
592  }
593  if constexpr (D == Direction::both || D == Direction::backward) {
594  bwd.template init_c2c<Direction::backward>(box, pio);
595  }
596 #endif
597 
598  return {fwd, bwd};
599 }
600 
601 template <typename T, Direction D, DomainStrategy S>
602 template <typename F>
603 void R2C<T,D,S>::post_forward_doit (F const& post_forward)
604 {
605  if (m_info.batch_mode) {
606  amrex::Abort("xxxxx todo: post_forward");
607  } else {
608  if ( ! m_cz.empty()) {
609  auto* spectral_fab = detail::get_fab(m_cz);
610  if (spectral_fab) {
611  auto const& a = spectral_fab->array(); // m_cz's ordering is z,x,y
612  ParallelFor(spectral_fab->box(),
613  [=] AMREX_GPU_DEVICE (int iz, int jx, int ky)
614  {
615  post_forward(jx,ky,iz,a(iz,jx,ky));
616  });
617  }
618  } else if ( ! m_cy.empty()) {
619  auto* spectral_fab = detail::get_fab(m_cy);
620  if (spectral_fab) {
621  auto const& a = spectral_fab->array(); // m_cy's ordering is y,x,z
622  ParallelFor(spectral_fab->box(),
623  [=] AMREX_GPU_DEVICE (int iy, int jx, int k)
624  {
625  post_forward(jx,iy,k,a(iy,jx,k));
626  });
627  }
628  } else {
629  auto* spectral_fab = detail::get_fab(m_cx);
630  if (spectral_fab) {
631  auto const& a = spectral_fab->array();
632  ParallelFor(spectral_fab->box(),
633  [=] AMREX_GPU_DEVICE (int i, int j, int k)
634  {
635  post_forward(i,j,k,a(i,j,k));
636  });
637  }
638  }
639  }
640 }
641 
642 template <typename T, Direction D, DomainStrategy S>
644 {
645 #if (AMREX_SPACEDIM == 3)
646  if (m_info.batch_mode) {
647  return T(1)/T(Long(m_real_domain.length(0)) *
648  Long(m_real_domain.length(1)));
649  } else
650 #endif
651  {
652  return T(1)/T(m_real_domain.numPts());
653  }
654 }
655 
656 template <typename T, Direction D, DomainStrategy S>
657 template <Direction DIR, std::enable_if_t<DIR == Direction::forward ||
658  DIR == Direction::both, int> >
659 std::pair<typename R2C<T,D,S>::cMF *, IntVect>
661 {
662  if (!m_cz.empty()) {
663  return std::make_pair(&m_cz, IntVect{AMREX_D_DECL(2,0,1)});
664  } else if (!m_cy.empty()) {
665  return std::make_pair(&m_cy, IntVect{AMREX_D_DECL(1,0,2)});
666  } else {
667  return std::make_pair(&m_cx, IntVect{AMREX_D_DECL(0,1,2)});
668  }
669 }
670 
671 template <typename T, Direction D, DomainStrategy S>
672 template <Direction DIR, std::enable_if_t<DIR == Direction::forward ||
673  DIR == Direction::both, int> >
674 void R2C<T,D,S>::forward (MF const& inmf, cMF& outmf)
675 {
676  BL_PROFILE("FFT::R2C::forward(inout)");
677 
678  forward(inmf);
679  if (!m_cz.empty()) { // m_cz's order (z,x,y) -> (x,y,z)
680  RotateBwd dtos{};
682  (outmf, m_spectral_domain_x, m_cz, IntVect(0), dtos);
683  ParallelCopy(outmf, m_cz, cmd, 0, 0, 1, dtos);
684  } else if (!m_cy.empty()) { // m_cy's order (y,x,z) -> (x,y,z)
686  (outmf, m_spectral_domain_x, m_cy, IntVect(0), m_dtos_y2x);
687  ParallelCopy(outmf, m_cy, cmd, 0, 0, 1, m_dtos_y2x);
688  } else {
689  outmf.ParallelCopy(m_cx, 0, 0, 1);
690  }
691 }
692 
693 template <typename T, Direction D, DomainStrategy S>
694 template <Direction DIR, std::enable_if_t<DIR == Direction::backward ||
695  DIR == Direction::both, int> >
696 void R2C<T,D,S>::backward (cMF const& inmf, MF& outmf)
697 {
698  BL_PROFILE("FFT::R2C::backward(inout)");
699 
700  if (!m_cz.empty()) { // (x,y,z) -> m_cz's order (z,x,y)
701  RotateFwd dtos{};
703  (m_cz, m_spectral_domain_z, inmf, IntVect(0), dtos);
704  ParallelCopy(m_cz, inmf, cmd, 0, 0, 1, dtos);
705  } else if (!m_cy.empty()) { // (x,y,z) -> m_cy's ordering (y,x,z)
707  (m_cy, m_spectral_domain_y, inmf, IntVect(0), m_dtos_x2y);
708  ParallelCopy(m_cy, inmf, cmd, 0, 0, 1, m_dtos_x2y);
709  } else {
710  m_cx.ParallelCopy(inmf, 0, 0, 1);
711  }
712  backward_doit(outmf);
713 }
714 
715 template <typename T, Direction D, DomainStrategy S>
716 std::pair<BoxArray,DistributionMapping>
718 {
719 #if (AMREX_SPACEDIM == 3)
720  if (!m_cz.empty()) {
721  BoxList bl = m_cz.boxArray().boxList();
722  for (auto& b : bl) {
723  auto lo = b.smallEnd();
724  auto hi = b.bigEnd();
725  std::swap(lo[0], lo[1]);
726  std::swap(lo[1], lo[2]);
727  std::swap(hi[0], hi[1]);
728  std::swap(hi[1], hi[2]);
729  b.setSmall(lo);
730  b.setBig(hi);
731  }
732  return std::make_pair(BoxArray(std::move(bl)), m_cz.DistributionMap());
733  } else
734 #endif
735 #if (AMREX_SPACEDIM >= 2)
736  if (!m_cy.empty()) {
737  BoxList bl = m_cy.boxArray().boxList();
738  for (auto& b : bl) {
739  auto lo = b.smallEnd();
740  auto hi = b.bigEnd();
741  std::swap(lo[0], lo[1]);
742  std::swap(hi[0], hi[1]);
743  b.setSmall(lo);
744  b.setBig(hi);
745  }
746  return std::make_pair(BoxArray(std::move(bl)), m_cy.DistributionMap());
747  } else
748 #endif
749  {
750  return std::make_pair(m_cx.boxArray(), m_cx.DistributionMap());
751  }
752 }
753 
754 }
755 
756 #endif
#define BL_PROFILE(a)
Definition: AMReX_BLProfiler.H:551
#define AMREX_ALWAYS_ASSERT(EX)
Definition: AMReX_BLassert.H:50
#define AMREX_GPU_DEVICE
Definition: AMReX_GpuQualifiers.H:18
#define AMREX_D_DECL(a, b, c)
Definition: AMReX_SPACE.H:104
if(!(yy_init))
Definition: amrex_iparser.lex.nolint.H:935
A collection of Boxes stored in an Array.
Definition: AMReX_BoxArray.H:550
A class for managing a List of Boxes that share a common IndexType. This class implements operations ...
Definition: AMReX_BoxList.H:52
AMREX_GPU_HOST_DEVICE const IntVectND< dim > & smallEnd() const &noexcept
Get the smallend of the BoxND.
Definition: AMReX_Box.H:105
AMREX_GPU_HOST_DEVICE const IntVectND< dim > & bigEnd() const &noexcept
Get the bigend.
Definition: AMReX_Box.H:116
AMREX_GPU_HOST_DEVICE IntVectND< dim > length() const noexcept
Return the length of the BoxND.
Definition: AMReX_Box.H:146
AMREX_GPU_HOST_DEVICE BoxND & growHi(int idir, int n_cell=1) noexcept
Grow the BoxND on the high end by n_cell cells in direction idir. NOTE: n_cell negative shrinks the B...
Definition: AMReX_Box.H:659
AMREX_GPU_HOST_DEVICE BoxND & growLo(int idir, int n_cell=1) noexcept
Grow the BoxND on the low end by n_cell cells in direction idir. NOTE: n_cell negative shrinks the Bo...
Definition: AMReX_Box.H:648
AMREX_GPU_HOST_DEVICE bool ok() const noexcept
Checks if it is a proper BoxND (including a valid type).
Definition: AMReX_Box.H:200
Calculates the distribution of FABs to MPI processes.
Definition: AMReX_DistributionMapping.H:41
Long size() const noexcept
Length of the underlying processor map.
Definition: AMReX_DistributionMapping.H:127
Definition: AMReX_FFT_OpenBCSolver.H:11
Parallel Discrete Fourier Transform.
Definition: AMReX_FFT_R2C.H:34
R2C(R2C const &)=delete
T scalingFactor() const
Definition: AMReX_FFT_R2C.H:643
void prepare_openbc()
Definition: AMReX_FFT_R2C.H:434
std::pair< cMF *, IntVect > getSpectralData()
Get the internal spectral data.
RotateFwd m_dtos_x2z
Definition: AMReX_FFT_R2C.H:199
std::unique_ptr< MultiBlockCommMetaData > m_cmd_x2z_half
Definition: AMReX_FFT_R2C.H:193
Swap01 m_dtos_x2y
Definition: AMReX_FFT_R2C.H:195
Plan< T > m_fft_bwd_y
Definition: AMReX_FFT_R2C.H:178
static std::pair< Plan< T >, Plan< T > > make_c2c_plans(cMF &inout)
Definition: AMReX_FFT_R2C.H:575
Box m_real_domain
Definition: AMReX_FFT_R2C.H:210
Plan< T > m_fft_fwd_x_half
Definition: AMReX_FFT_R2C.H:181
Plan< T > m_fft_fwd_x
Definition: AMReX_FFT_R2C.H:175
bool m_do_alld_fft
Definition: AMReX_FFT_R2C.H:217
Swap02 m_dtos_y2z
Definition: AMReX_FFT_R2C.H:197
std::unique_ptr< MultiBlockCommMetaData > m_cmd_x2y
Definition: AMReX_FFT_R2C.H:187
std::unique_ptr< MultiBlockCommMetaData > m_cmd_z2y
Definition: AMReX_FFT_R2C.H:190
Plan< T > m_fft_bwd_x
Definition: AMReX_FFT_R2C.H:176
R2C(R2C &&)=delete
void forward(MF const &inmf, cMF &outmf)
Forward transform.
Definition: AMReX_FFT_R2C.H:674
cMF m_cz
Definition: AMReX_FFT_R2C.H:205
std::unique_ptr< MultiBlockCommMetaData > m_cmd_z2x_half
Definition: AMReX_FFT_R2C.H:194
Plan< T > m_fft_bwd_z
Definition: AMReX_FFT_R2C.H:180
bool m_slab_decomp
Definition: AMReX_FFT_R2C.H:218
Plan< T > m_fft_bwd_x_half
Definition: AMReX_FFT_R2C.H:182
cMF m_cx
Definition: AMReX_FFT_R2C.H:203
void backward(MF &outmf)
Backward transform.
Definition: AMReX_FFT_R2C.H:536
Box m_spectral_domain_y
Definition: AMReX_FFT_R2C.H:212
RotateBwd m_dtos_z2x
Definition: AMReX_FFT_R2C.H:200
std::pair< BoxArray, DistributionMapping > getSpectralDataLayout() const
Get BoxArray and DistributionMapping for spectral data.
Definition: AMReX_FFT_R2C.H:717
Info m_info
Definition: AMReX_FFT_R2C.H:215
std::unique_ptr< MultiBlockCommMetaData > m_cmd_y2z
Definition: AMReX_FFT_R2C.H:189
Box m_spectral_domain_z
Definition: AMReX_FFT_R2C.H:213
std::unique_ptr< MultiBlockCommMetaData > m_cmd_z2x
Definition: AMReX_FFT_R2C.H:192
std::conditional_t< std::is_same_v< T, Real >, MultiFab, FabArray< BaseFab< T > > > MF
Definition: AMReX_FFT_R2C.H:37
Plan< T > m_fft_fwd_y
Definition: AMReX_FFT_R2C.H:177
R2C(Box const &domain, Info const &info=Info{})
Constructor.
Definition: AMReX_FFT_R2C.H:223
std::unique_ptr< MultiBlockCommMetaData > m_cmd_y2x
Definition: AMReX_FFT_R2C.H:188
bool m_openbc_half
Definition: AMReX_FFT_R2C.H:219
void post_forward_doit(F const &post_forward)
Definition: AMReX_FFT_R2C.H:603
Plan< T > m_fft_fwd_z
Definition: AMReX_FFT_R2C.H:179
std::unique_ptr< char, DataDeleter > m_data_2
Definition: AMReX_FFT_R2C.H:208
Swap02 m_dtos_z2y
Definition: AMReX_FFT_R2C.H:198
R2C & operator=(R2C const &)=delete
MF m_rx
Definition: AMReX_FFT_R2C.H:202
void backward_doit(MF &outmf, IntVect const &ngout=IntVect(0))
Definition: AMReX_FFT_R2C.H:542
Box m_spectral_domain_x
Definition: AMReX_FFT_R2C.H:211
Swap01 m_dtos_y2x
Definition: AMReX_FFT_R2C.H:196
std::unique_ptr< char, DataDeleter > m_data_1
Definition: AMReX_FFT_R2C.H:207
cMF m_cy
Definition: AMReX_FFT_R2C.H:204
~R2C()
Definition: AMReX_FFT_R2C.H:413
void backward(cMF const &inmf, MF &outmf)
Backward transform.
Definition: AMReX_FFT_R2C.H:696
void forward(MF const &inmf)
Forward transform.
Definition: AMReX_FFT_R2C.H:488
std::unique_ptr< MultiBlockCommMetaData > m_cmd_x2z
Definition: AMReX_FFT_R2C.H:191
void forwardThenBackward(MF const &inmf, MF &outmf, F const &post_forward)
Forward and then backward transform.
Definition: AMReX_FFT_R2C.H:78
bool empty() const noexcept
Definition: AMReX_FabArrayBase.H:88
An Array of FortranArrayBox(FAB)-like Objects.
Definition: AMReX_FabArray.H:344
void ParallelCopy(const FabArray< FAB > &src, const Periodicity &period=Periodicity::NonPeriodic(), CpOp op=FabArrayBase::COPY)
Definition: AMReX_FabArray.H:778
void define(const BoxArray &bxs, const DistributionMapping &dm, int nvar, int ngrow, const MFInfo &info=MFInfo(), const FabFactory< FAB > &factory=DefaultFabFactory< FAB >())
Define this FabArray identically to that performed by the constructor having an analogous function si...
Definition: AMReX_FabArray.H:2027
A collection (stored as an array) of FArrayBox objects.
Definition: AMReX_MultiFab.H:38
FA::FABType::value_type * get_fab(FA &fa)
Definition: AMReX_FFT_Helper.H:1303
std::unique_ptr< char, DataDeleter > make_mfs_share(FA1 &fa1, FA2 &fa2)
Definition: AMReX_FFT_Helper.H:1314
DistributionMapping make_iota_distromap(Long n)
Definition: AMReX_FFT.cpp:88
Definition: AMReX_FFT.cpp:7
Direction
Definition: AMReX_FFT_Helper.H:46
DomainStrategy
Definition: AMReX_FFT_Helper.H:48
std::enable_if_t< IsBaseFab< FAB >) &&IsDataPacking< DataPacking, FAB >)> ParallelCopy_finish(FabArray< FAB > &dest, CommHandler handler, const FabArrayBase::CommMetaData &cmd, const DataPacking &data_packing)
Definition: AMReX_NonLocalBC.H:793
AMREX_NODISCARD CommHandler ParallelCopy_nowait(NoLocalCopy, FabArray< FAB > &dest, const FabArray< FAB > &src, const FabArrayBase::CommMetaData &cmd, const DataPacking &data_packing)
Definition: AMReX_NonLocalBC.H:701
int MyProcSub() noexcept
my sub-rank in current frame
Definition: AMReX_ParallelContext.H:76
int NProcsSub() noexcept
number of ranks in current frame
Definition: AMReX_ParallelContext.H:74
int NProcs() noexcept
return the number of MPI ranks local to the current Parallel Context
Definition: AMReX_ParallelDescriptor.H:243
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void swap(T &a, T &b) noexcept
Definition: AMReX_algoim_K.H:113
@ min
Definition: AMReX_ParallelReduce.H:18
constexpr int iz
Definition: AMReX_Interp_3D_C.H:37
constexpr int iy
Definition: AMReX_Interp_2D_C.H:33
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:200
IntVectND< AMREX_SPACEDIM > IntVect
Definition: AMReX_BaseFwd.H:30
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 length(Array4< T > const &a) noexcept
Definition: AMReX_Array4.H:322
BoxArray decompose(Box const &domain, int nboxes, Array< bool, AMREX_SPACEDIM > const &decomp={AMREX_D_DECL(true, true, true)}, bool no_overlap=false)
Decompose domain box into BoxArray.
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition: AMReX.cpp:225
void ParallelCopy(MF &dst, MF const &src, int scomp, int dcomp, int ncomp, IntVect const &ng_src=IntVect(0), IntVect const &ng_dst=IntVect(0), Periodicity const &period=Periodicity::NonPeriodic())
dst = src w/ MPI communication
Definition: AMReX_FabArrayUtility.H:1672
Definition: AMReX_FFT_Helper.H:56
bool batch_mode
Definition: AMReX_FFT_Helper.H:60
int nprocs
Max number of processes to use.
Definition: AMReX_FFT_Helper.H:63
Definition: AMReX_FFT_Helper.H:111
std::conditional_t< std::is_same_v< float, T >, cuComplex, cuDoubleComplex > VendorComplex
Definition: AMReX_FFT_Helper.H:115
Definition: AMReX_FFT_Helper.H:1426
Definition: AMReX_FFT_Helper.H:1401
Definition: AMReX_FFT_Helper.H:1355
Definition: AMReX_FFT_Helper.H:1378
FabArray memory allocation information.
Definition: AMReX_FabArray.H:66
This class specializes behaviour on local copies and unpacking receive buffers.
Definition: AMReX_NonLocalBC.H:615
This is the index mapping based on the DTOS MultiBlockDestToSrc.
Definition: AMReX_NonLocalBC.H:210
Contains information about which components take part of the data transaction.
Definition: AMReX_NonLocalBC.H:528