Block-Structured AMR Software Framework
AMReX_NonLocalBCImpl.H
Go to the documentation of this file.
1 #ifndef AMREX_NONLOCAL_BC_H_
2 #include "AMReX_NonLocalBC.H"
3 #endif
4 
5 #ifndef AMREX_NONLOCAL_BC_IMPL_H_
6 #define AMREX_NONLOCAL_BC_IMPL_H_
7 #include <AMReX_Config.H>
8 #include <AMReX_TypeTraits.H>
9 
10 namespace amrex::NonLocalBC {
11 struct Rotate90ClockWise {
13  IntVect operator() (IntVect const& iv) const noexcept {
14  return IntVect{AMREX_D_DECL(iv[1], -1-iv[0], iv[2])};
15  }
16 
18  Dim3 operator() (Dim3 const& a) const noexcept {
19  return Dim3{a.y, -1-a.x, a.z};
20  }
21 
22  Box operator() (Box const& box) const noexcept {
23  return Box(operator()(IntVect{AMREX_D_DECL(box.bigEnd (0),
24  box.smallEnd(1),
25  box.smallEnd(2))}),
26  operator()(IntVect{AMREX_D_DECL(box.smallEnd(0),
27  box.bigEnd (1),
28  box.bigEnd (2))}));
29  }
30 };
31 
32 struct Rotate90CounterClockWise {
34  IntVect operator() (IntVect const& iv) const noexcept {
35  return IntVect{AMREX_D_DECL(-1-iv[1], iv[0], iv[2])};
36  }
37 
39  Dim3 operator() (Dim3 const& a) const noexcept {
40  return Dim3{-1-a.y, a.x, a.z};
41  }
42 
43  Box operator() (Box const& box) const noexcept {
44  return Box(operator()(IntVect{AMREX_D_DECL(box.smallEnd(0),
45  box.bigEnd (1),
46  box.smallEnd(2))}),
47  operator()(IntVect{AMREX_D_DECL(box.bigEnd (0),
48  box.smallEnd(1),
49  box.bigEnd (2))}));
50  }
51 };
52 
53 struct Rotate90DstToSrc
54 {
56  Dim3 operator() (Dim3 const& a) const noexcept {
57  if (a.x < 0) {
58  return Rotate90ClockWise()(a);
59  } else {
60  return Rotate90CounterClockWise()(a);
61  }
62  }
63 };
64 
65 struct Rotate180Fn {
66  int Ly;
67 
69  IntVect operator() (IntVect const& iv) const noexcept {
70  return IntVect{AMREX_D_DECL(-1-iv[0], Ly-1-iv[1], iv[2])};
71  }
72 
74  Dim3 operator() (Dim3 const& a) const noexcept {
75  return Dim3{-1-a.x, Ly-1-a.y, a.z};
76  }
77 
78  Box operator() (Box const& box) const noexcept {
79  return Box(operator()(IntVect{AMREX_D_DECL(box.bigEnd (0),
80  box.bigEnd (1),
81  box.smallEnd(2))}),
82  operator()(IntVect{AMREX_D_DECL(box.smallEnd(0),
83  box.smallEnd(1),
84  box.bigEnd (2))}));
85  }
86 };
87 
88 struct PolarFn {
89  int Lx, Ly;
90 
91  [[nodiscard]] AMREX_GPU_HOST_DEVICE
92  int i_index (int i) const noexcept {
93  return (i < Lx/2) ? -1-i : 2*Lx-1-i;
94  }
95 
96  [[nodiscard]] AMREX_GPU_HOST_DEVICE
97  int j_index (int j) const noexcept {
98  return (j < Ly/2) ? j+Ly/2 : j-Ly/2;
99  }
100 
101  [[nodiscard]] AMREX_GPU_HOST_DEVICE
102  IntVect operator() (IntVect const& iv) const noexcept {
103  return IntVect{AMREX_D_DECL(i_index(iv[0]), j_index(iv[1]), iv[2])};
104  }
105 
106  [[nodiscard]] AMREX_GPU_HOST_DEVICE
107  Dim3 operator() (Dim3 const& a) const noexcept {
108  return Dim3{i_index(a.x), j_index(a.y), a.z};
109  }
110 
111  [[nodiscard]] Box operator() (Box const& box) const noexcept {
112  return Box(operator()(IntVect{AMREX_D_DECL(box.bigEnd (0),
113  box.smallEnd(1),
114  box.smallEnd(2))}),
115  operator()(IntVect{AMREX_D_DECL(box.smallEnd(0),
116  box.bigEnd (1),
117  box.bigEnd (2))}));
118  }
119 };
120 
121 struct PolarFn2 { // for the x-y corners
122  int Lx, Ly;
123 
124  [[nodiscard]] AMREX_GPU_HOST_DEVICE
125  int i_index (int i) const noexcept {
126  return (i < Lx/2) ? -1-i : 2*Lx-1-i;
127  }
128 
129  [[nodiscard]] AMREX_GPU_HOST_DEVICE
130  int j_index (int j) const noexcept {
131  if (j < 0) { // NOLINT
132  return j+Ly/2;
133  } else if (j >= Ly) { // NOLINT
134  return j-Ly/2;
135  } else if (j < Ly/2) {
136  return j-Ly/2;
137  } else {
138  return j+Ly/2;
139  }
140  }
141 
142  [[nodiscard]] AMREX_GPU_HOST_DEVICE
143  IntVect operator() (IntVect const& iv) const noexcept {
144  return IntVect{AMREX_D_DECL(i_index(iv[0]), j_index(iv[1]), iv[2])};
145  }
146 
147  [[nodiscard]] AMREX_GPU_HOST_DEVICE
148  Dim3 operator() (Dim3 const& a) const noexcept {
149  return Dim3{i_index(a.x), j_index(a.y), a.z};
150  }
151 
152  [[nodiscard]] Box operator() (Box const& box) const noexcept {
153  return Box(operator()(IntVect{AMREX_D_DECL(box.bigEnd (0),
154  box.smallEnd(1),
155  box.smallEnd(2))}),
156  operator()(IntVect{AMREX_D_DECL(box.smallEnd(0),
157  box.bigEnd (1),
158  box.bigEnd (2))}));
159  }
160 };
161 
186 template <class FAB, class DTOS, class Proj>
187 std::enable_if_t<IsBaseFab<FAB>() && IsCallableR<Dim3, DTOS, Dim3>() && IsFabProjection<Proj, FAB>()>
188 local_copy_cpu (FabArray<FAB>& dest, const FabArray<FAB>& src, int dcomp, int scomp, int ncomp,
189  FabArrayBase::CopyComTagsContainer const& local_tags, DTOS const& dtos, Proj const& proj) noexcept {
190  const auto N_locs = static_cast<int>(local_tags.size());
191  if (N_locs == 0) { return; }
192 #ifdef AMREX_USE_OMP
193 #pragma omp parallel for
194 #endif
195  for (int itag = 0; itag < N_locs; ++itag) {
196  const auto& tag = local_tags[itag];
197  auto const& sfab = src.const_array(tag.srcIndex);
198  auto const& dfab = dest.array (tag.dstIndex);
199  amrex::LoopConcurrentOnCpu(tag.dbox, ncomp, [=] (int i, int j, int k, int n) noexcept
200  {
201  auto const si = dtos(Dim3{i,j,k});
202  dfab(i,j,k,dcomp+n) = proj(sfab,si,scomp+n);
203  });
204  }
205 }
206 
212 template <class FAB, class DTOS, class Proj>
213 std::enable_if_t<IsBaseFab<FAB>() && IsCallableR<Dim3, DTOS, Dim3>() && IsFabProjection<Proj, FAB>()>
214 unpack_recv_buffer_cpu (FabArray<FAB>& mf, int dcomp, int ncomp, Vector<char*> const& recv_data,
215  Vector<std::size_t> const& recv_size,
216  Vector<FabArrayBase::CopyComTagsContainer const*> const& recv_cctc,
217  DTOS const& dtos, Proj const& proj) noexcept {
218  amrex::ignore_unused(recv_size);
219 
220  const auto N_rcvs = static_cast<int>(recv_cctc.size());
221  if (N_rcvs == 0) { return; }
222 
223  using T = typename FAB::value_type;
224 #ifdef AMREX_USE_OMP
225 #pragma omp parallel for
226 #endif
227  for (int ircv = 0; ircv < N_rcvs; ++ircv) {
228  const char* dptr = recv_data[ircv];
229  auto const& cctc = *recv_cctc[ircv];
230  for (auto const& tag : cctc) {
231  auto const& dfab = mf.array(tag.dstIndex);
232  auto const& sfab = amrex::makeArray4((T const*)(dptr), tag.sbox, ncomp);
233  amrex::LoopConcurrentOnCpu(tag.dbox, ncomp, [=](int i, int j, int k, int n) noexcept {
234  auto const si = dtos(Dim3{i, j, k});
235  dfab(i, j, k, dcomp + n) = proj(sfab, si, n);
236  });
237  dptr += tag.sbox.numPts() * ncomp * sizeof(T);
238  AMREX_ASSERT(dptr <= recv_data[ircv] + recv_size[ircv]);
239  }
240  }
241 }
242 
243 #ifdef AMREX_USE_GPU
244 template <class T>
245 struct Array4Array4Box {
246  Array4<T > dfab;
247  Array4<T const> sfab;
248  Box dbox;
249 
251  Box const& box () const noexcept { return dbox; }
252 };
253 
254 template <class FAB, class DTOS, class Proj>
255 std::enable_if_t<IsBaseFab<FAB>() && IsCallableR<Dim3, DTOS, Dim3>() && IsFabProjection<Proj, FAB>()>
256 local_copy_gpu (FabArray<FAB>& dest, const FabArray<FAB>& src, int dcomp, int scomp, int ncomp,
257  FabArrayBase::CopyComTagsContainer const& local_tags, DTOS const& dtos, Proj const& proj) noexcept {
258  int N_locs = local_tags.size();
259  if (N_locs == 0) { return; }
260 
261  using T = typename FAB::value_type;
262  Vector<Array4Array4Box<T> > loc_copy_tags;
263  loc_copy_tags.reserve(N_locs);
264  for (auto const& tag : local_tags) {
265  loc_copy_tags.push_back({dest.array(tag.dstIndex), src.const_array(tag.srcIndex), tag.dbox});
266  }
267 
268  ParallelFor(loc_copy_tags, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n,
269  Array4Array4Box<T> const& tag) noexcept
270  {
271  auto const si = dtos(Dim3{i,j,k});
272  tag.dfab(i,j,k,dcomp+n) = proj(tag.sfab, si, scomp+n);
273  });
274 }
275 
276 template <class FAB, class DTOS, class Proj>
277 std::enable_if_t<IsBaseFab<FAB>() && IsCallableR<Dim3, DTOS, Dim3>() && IsFabProjection<Proj, FAB>()>
278 unpack_recv_buffer_gpu (FabArray<FAB>& mf, int scomp, int ncomp,
279  Vector<char*> const& recv_data,
280  Vector<std::size_t> const& recv_size,
281  Vector<FabArrayBase::CopyComTagsContainer const*> const& recv_cctc,
282  DTOS const& dtos, Proj const& proj)
283 {
284  amrex::ignore_unused(recv_size);
285 
286  const int N_rcvs = recv_cctc.size();
287  if (N_rcvs == 0) { return; }
288 
289  char* pbuffer = recv_data[0];
290 #if 0
291  std::size_t szbuffer = 0;
292  // For linear solver test on summit, this is slower than writing to
293  // pinned memory directly on device.
295  // Memory in recv_data is pinned.
296  szbuffer = (recv_data[N_rcvs-1]-recv_data[0]) + recv_size[N_rcvs-1];
297  pbuffer = (char*)The_Arena()->alloc(szbuffer);
298  Gpu::copyAsync(Gpu::hostToDevice,recv_data[0],recv_data[0]+szbuffer,pbuffer);
300  }
301 #endif
302 
303  using T = typename FAB::value_type;
304  using TagType = Array4Array4Box<T>;
305  Vector<TagType> tags;
306  tags.reserve(N_rcvs);
307 
308  for (int k = 0; k < N_rcvs; ++k)
309  {
310  std::size_t offset = recv_data[k]-recv_data[0];
311  const char* dptr = pbuffer + offset;
312  auto const& cctc = *recv_cctc[k];
313  for (auto const& tag : cctc)
314  {
315  tags.emplace_back(TagType{mf.array(tag.dstIndex),
316  amrex::makeArray4((T const*)dptr, tag.sbox, ncomp),
317  tag.dbox});
318  dptr += tag.dbox.numPts() * ncomp * sizeof(T);
319  BL_ASSERT(dptr <= pbuffer + offset + recv_size[k]);
320  }
321  }
322 
323  ParallelFor(tags, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n,
324  Array4Array4Box<T> const& tag) noexcept
325  {
326  auto const si = dtos(Dim3{i,j,k});
327  tag.dfab(i,j,k,scomp+n) = proj(tag.sfab, si ,n);
328  });
329 
330  // There is Gpu::streamSynchronize in ParallelFor above
331 
332  if (pbuffer != recv_data[0]) {
333  The_Arena()->free(pbuffer);
334  }
335 }
336 #endif
337 
338 template <typename DTOS, typename>
339 MultiBlockCommMetaData::MultiBlockCommMetaData (const FabArrayBase& dst, const Box& dstbox, const FabArrayBase& src,
340  const IntVect& ngrow, DTOS const& dtos)
341  : MultiBlockCommMetaData(dst.boxArray(), dst.DistributionMap(), dstbox, src.boxArray(),
342  src.DistributionMap(), ngrow, dtos) {}
343 
344 template <typename DTOS, typename>
345 MultiBlockCommMetaData::MultiBlockCommMetaData (const BoxArray& dstba, const DistributionMapping& dstdm,
346  const Box& dstbox, const BoxArray& srcba,
347  const DistributionMapping& srcdm, const IntVect& ngrow, DTOS const& dtos) {
348  define(dstba, dstdm, dstbox, srcba, srcdm, ngrow, dtos);
349 }
350 
351 template <typename DTOS>
352 std::enable_if_t<IsIndexMapping<DTOS>::value>
353 MultiBlockCommMetaData::define (const BoxArray& dstba, const DistributionMapping& dstdm, const Box& dstbox,
354  const BoxArray& srcba, const DistributionMapping& srcdm, const IntVect& ngrow,
355  DTOS const& dtos) {
356  m_LocTags = std::make_unique<FabArrayBase::CopyComTagsContainer>();
357  m_SndTags = std::make_unique<FabArrayBase::MapOfCopyComTagContainers>();
358  m_RcvTags = std::make_unique<FabArrayBase::MapOfCopyComTagContainers>();
359  const int myproc = ParallelDescriptor::MyProc();
360  for (int i = 0, N = static_cast<int>(dstba.size()); i < N; ++i) {
361  const int dest_owner = dstdm[i];
362  const Box partial_dstbox = grow(dstba[i], ngrow) & dstbox;
363  if (partial_dstbox.isEmpty()) {
364  continue;
365  }
366  const Box partial_dstbox_mapped_in_src = Image(dtos, partial_dstbox).setType(srcba.ixType());
367  enum { not_first_only = 0, first_only = 1 };
368  std::vector<std::pair<int, Box>> boxes_from_src =
369  srcba.intersections(partial_dstbox_mapped_in_src, not_first_only, ngrow);
370  for (std::pair<int, Box> counted_box : boxes_from_src) {
371  const int k = counted_box.first;
372  const Box src_box = counted_box.second;
373  AMREX_ASSERT(k < srcba.size());
374  const int src_owner = srcdm[k];
375  if (dest_owner == myproc || src_owner == myproc) {
376  if (src_owner == dest_owner) {
377  const BoxList tilelist(src_box, FabArrayBase::comm_tile_size);
378  for (const Box& tilebox : tilelist) {
379  const Box inverse_image = InverseImage(dtos, tilebox).setType(dstba.ixType());
380  if ((inverse_image & partial_dstbox).ok()) {
381  m_LocTags->emplace_back(inverse_image, tilebox, i, k);
382  }
383  }
384  } else {
385  const Box inverse_image = InverseImage(dtos, src_box).setType(dstba.ixType());
386  if ((inverse_image & partial_dstbox).ok()) {
387  FabArrayBase::CopyComTagsContainer& copy_tags =
388  (src_owner == myproc) ? (*m_SndTags)[dest_owner]
389  : (*m_RcvTags)[src_owner];
390  copy_tags.emplace_back(inverse_image, src_box, i, k);
391  }
392  }
393  }
394  }
395  }
396 }
397 
398 template <class FAB, class DTOS, class Proj>
399 #ifdef AMREX_USE_MPI
401 #endif
402 CommHandler
403 Comm_nowait (FabArray<FAB>& mf, int scomp, int ncomp, FabArrayBase::CommMetaData const& cmd,
404  DTOS const& dtos, Proj const& proj)
405 {
406 #ifdef AMREX_USE_MPI
407  if (ParallelContext::NProcsSub() == 1)
408 #endif
409  {
410  if (cmd.m_LocTags->empty()) { return CommHandler{}; }
411 #ifdef AMREX_USE_GPU
412  if (Gpu::inLaunchRegion()) {
413  local_copy_gpu(mf, mf, scomp, scomp, ncomp, *cmd.m_LocTags, dtos, proj);
414  } else
415 #endif
416  {
417  local_copy_cpu(mf, mf, scomp, scomp, ncomp, *cmd.m_LocTags, dtos, proj);
418  }
419  return CommHandler{};
420  }
421 
422 #ifdef AMREX_USE_MPI
423  //
424  // Do this before prematurely exiting if running in parallel.
425  // Otherwise sequence numbers will not match across MPI processes.
426  //
428 
429  const auto N_locs = cmd.m_LocTags->size();
430  const auto N_rcvs = cmd.m_RcvTags->size();
431  const auto N_snds = cmd.m_SndTags->size();
432 
433  if (N_locs == 0 && N_rcvs == 0 && N_snds == 0) {
434  // No work to do.
435  return CommHandler{};
436  }
437 
438  CommHandler handler{};
439  handler.mpi_tag = SeqNum;
440 
441  if (N_rcvs > 0) {
442  handler.recv.the_data = FabArray<FAB>::PostRcvs(*cmd.m_RcvTags, handler.recv.data, handler.recv.size,
443  handler.recv.rank, handler.recv.request, ncomp, SeqNum);
444  }
445 
446  if (N_snds > 0) {
447  handler.send.the_data =
448  FabArray<FAB>::PrepareSendBuffers(*cmd.m_SndTags, handler.send.data, handler.send.size,
449  handler.send.rank, handler.send.request, handler.send.cctc, ncomp);
450 
451 #ifdef AMREX_USE_GPU
452  if (Gpu::inLaunchRegion()) {
453  FabArray<FAB>::pack_send_buffer_gpu(mf, scomp, ncomp, handler.send.data,
454  handler.send.size, handler.send.cctc);
455  } else
456 #endif
457  {
458  FabArray<FAB>::pack_send_buffer_cpu(mf, scomp, ncomp, handler.send.data,
459  handler.send.size, handler.send.cctc);
460  }
461 
462  FabArray<FAB>::PostSnds(handler.send.data, handler.send.size, handler.send.rank, handler.send.request, SeqNum);
463  }
464 
465  if (N_locs > 0)
466  {
467 #ifdef AMREX_USE_GPU
468  if (Gpu::inLaunchRegion()) {
469  local_copy_gpu(mf, mf, scomp, scomp, ncomp, *cmd.m_LocTags, dtos, proj);
470  } else
471 #endif
472  {
473  local_copy_cpu(mf, mf, scomp, scomp, ncomp, *cmd.m_LocTags, dtos, proj);
474  }
475  }
476 
477  return handler;
478 #endif
479 }
480 
481 #ifdef AMREX_USE_MPI
482 template <class FAB, class DTOS, class Proj>
483 void
484 Comm_finish (FabArray<FAB>& mf, int scomp, int ncomp, FabArrayBase::CommMetaData const& cmd,
485  CommHandler handler, DTOS const& dtos, Proj const& proj)
486 {
487  if (ParallelContext::NProcsSub() == 1) { return; }
488 
489  const auto N_rcvs = static_cast<int>(cmd.m_RcvTags->size());
490  if (N_rcvs > 0)
491  {
492  handler.recv.cctc.resize(N_rcvs, nullptr);
493  for (int k = 0; k < N_rcvs; ++k) {
494  auto const& cctc = cmd.m_RcvTags->at(handler.recv.rank[k]);
495  handler.recv.cctc[k] = &cctc;
496  }
497  handler.recv.stats.resize(handler.recv.request.size());
498  ParallelDescriptor::Waitall(handler.recv.request, handler.recv.stats);
499 #ifdef AMREX_DEBUG
500  if (!CheckRcvStats(handler.recv.stats, handler.recv.size, handler.mpi_tag)) {
501  amrex::Abort("NonLocalBC::Comm_finish failed with wrong message size");
502  }
503 #endif
504 
505 #ifdef AMREX_USE_GPU
506  if (Gpu::inLaunchRegion())
507  {
508  unpack_recv_buffer_gpu(mf, scomp, ncomp, handler.recv.data,
509  handler.recv.size, handler.recv.cctc, dtos, proj);
510  } else
511 #endif
512  {
513  unpack_recv_buffer_cpu(mf, scomp, ncomp, handler.recv.data,
514  handler.recv.size, handler.recv.cctc, dtos, proj);
515  }
516  }
517 
518  if ( ! cmd.m_SndTags->empty() ) {
519  handler.send.stats.resize(handler.send.request.size());
520  ParallelDescriptor::Waitall(handler.send.request, handler.send.stats);
521  }
522 }
523 #endif
524 
525 template <class FAB>
526 std::enable_if_t<IsBaseFab<FAB>::value>
527 Rotate90 (FabArray<FAB>& mf, int scomp, int ncomp, IntVect const& nghost, Box const& domain)
528 {
529  BL_PROFILE("Rotate90");
530 
531  AMREX_ASSERT(domain.cellCentered());
532  AMREX_ASSERT(domain.smallEnd() == 0);
533  AMREX_ASSERT(domain.length(0) == domain.length(1));
534  AMREX_ASSERT(mf.is_cell_centered());
535  AMREX_ASSERT(scomp < mf.nComp() && scomp+ncomp <= mf.nComp());
536  AMREX_ASSERT(nghost.allLE(mf.nGrowVect()) && nghost[0] == nghost[1]);
537 
538  if (nghost[0] <= 0) { return; }
539 
540  const FabArrayBase::RB90& TheRB90 = mf.getRB90(nghost, domain);
541 
542  auto handler = Comm_nowait(mf, scomp, ncomp, TheRB90,Rotate90DstToSrc{},
543  Identity{});
544 
545  Box corner(-nghost, IntVect{AMREX_D_DECL(-1,-1,domain.bigEnd(2)+nghost[2])});
546 #ifdef AMREX_USE_OMP
547 #pragma omp parallel if (Gpu::notInLaunchRegion())
548 #endif
549  for (MFIter mfi(mf); mfi.isValid(); ++mfi) {
550  Box const& bx = corner & mfi.fabbox();
551  if (bx.ok()) {
552  auto const& fab = mf.array(mfi);
553  AMREX_HOST_DEVICE_PARALLEL_FOR_4D(bx,ncomp,i,j,k,n,
554  {
555  fab(i,j,k,n) = fab(-i-1,-j-1,k,n);
556  });
557  }
558  }
559 
560 #ifdef AMREX_USE_MPI
561  Comm_finish(mf, scomp, ncomp, TheRB90, std::move(handler), Rotate90DstToSrc{},
562  Identity{});
563 #else
564  amrex::ignore_unused(handler);
565 #endif
566 }
567 
568 template <class FAB>
569 std::enable_if_t<IsBaseFab<FAB>::value>
570 Rotate90 (FabArray<FAB>& mf, Box const& domain)
571 {
572  Rotate90(mf, 0, mf.nComp(), mf.nGrowVect(), domain);
573 }
574 
575 template <class FAB>
576 std::enable_if_t<IsBaseFab<FAB>::value>
577 Rotate180 (FabArray<FAB>& mf, int scomp, int ncomp, IntVect const& nghost, Box const& domain)
578 {
579  BL_PROFILE("Rotate180");
580 
581  AMREX_ASSERT(domain.cellCentered());
582  AMREX_ASSERT(domain.smallEnd() == 0);
583  AMREX_ASSERT(domain.length(1) % 2 == 0);
584  AMREX_ASSERT(mf.is_cell_centered());
585  AMREX_ASSERT(scomp < mf.nComp() && scomp+ncomp <= mf.nComp());
586  AMREX_ASSERT(nghost.allLE(mf.nGrowVect()));
587 
588  if (nghost[0] <= 0) { return; }
589 
590  const FabArrayBase::RB180& TheRB180 = mf.getRB180(nghost, domain);
591 
592  auto handler = Comm_nowait(mf, scomp, ncomp, TheRB180,
593  Rotate180Fn{domain.length(1)}, Identity{});
594 
595 #ifdef AMREX_USE_MPI
596  Comm_finish(mf, scomp, ncomp, TheRB180, std::move(handler),
597  Rotate180Fn{domain.length(1)}, Identity{});
598 #else
599  amrex::ignore_unused(handler);
600 #endif
601 }
602 
603 template <class FAB>
604 std::enable_if_t<IsBaseFab<FAB>::value>
605 Rotate180 (FabArray<FAB>& mf, Box const& domain)
606 {
607  Rotate180(mf, 0, mf.nComp(), mf.nGrowVect(), domain);
608 }
609 
610 template <class FAB>
611 std::enable_if_t<IsBaseFab<FAB>::value>
612 FillPolar (FabArray<FAB>& mf, int scomp, int ncomp, IntVect const& nghost, Box const& domain)
613 {
614  BL_PROFILE("FillPolar");
615 
616  AMREX_ASSERT(domain.cellCentered());
617  AMREX_ASSERT(domain.smallEnd() == 0);
618  AMREX_ASSERT(domain.length(1) % 2 == 0);
619  AMREX_ASSERT(mf.is_cell_centered());
620  AMREX_ASSERT(scomp < mf.nComp() && scomp+ncomp <= mf.nComp());
621  AMREX_ASSERT(nghost.allLE(mf.nGrowVect()));
622 
623  if (nghost[0] <= 0) { return; }
624 
625  const FabArrayBase::PolarB& ThePolarB = mf.getPolarB(nghost, domain);
626 
627  auto handler = Comm_nowait(mf, scomp, ncomp, ThePolarB,
628  PolarFn{domain.length(0), domain.length(1)},
629  Identity{});
630 
631 #ifdef AMREX_USE_MPI
632  Comm_finish(mf, scomp, ncomp, ThePolarB, std::move(handler),
633  PolarFn{domain.length(0), domain.length(1)}, Identity{});
634 #else
635  amrex::ignore_unused(handler);
636 #endif
637 }
638 
639 template <class FAB>
640 std::enable_if_t<IsBaseFab<FAB>::value>
641 FillPolar (FabArray<FAB>& mf, Box const& domain)
642 {
643  FillPolar(mf, 0, mf.nComp(), mf.nGrowVect(), domain);
644 }
645 
646 template <typename FAB, typename DTOS, typename Proj>
647 std::enable_if_t<IsBaseFab<FAB>() &&
648  IsCallableR<Dim3,DTOS,Dim3>() &&
649  IsFabProjection<Proj,FAB>(),
650  CommHandler>
651 FillBoundary_nowait (FabArray<FAB>& mf, const FabArrayBase::CommMetaData& cmd,
652  int scomp, int ncomp, DTOS const& dtos, Proj const& proj)
653 {
654  BL_PROFILE("FillBoundary_nowait(cmd)");
655  AMREX_ASSERT(scomp < mf.nComp() && scomp+ncomp <= mf.nComp());
656  return Comm_nowait(mf, scomp, ncomp, cmd, dtos, proj);
657 }
658 
659 template <typename FAB, typename DTOS, typename Proj>
660 std::enable_if_t<IsBaseFab<FAB>() &&
661  IsCallableR<Dim3,DTOS,Dim3>() &&
662  IsFabProjection<Proj,FAB>()>
663 FillBoundary_finish (CommHandler handler,
664  FabArray<FAB>& mf, const FabArrayBase::CommMetaData& cmd,
665  int scomp, int ncomp, DTOS const& dtos, Proj const& proj)
666 {
667 #ifdef AMREX_USE_MPI
668  BL_PROFILE("FillBoundary_finish(cmd)");
669  Comm_finish(mf, scomp, ncomp, cmd, std::move(handler), dtos, proj);
670 #else
671  amrex::ignore_unused(handler,mf,cmd,scomp,ncomp,dtos,proj);
672 #endif
673 }
674 
675 template <typename DTOS>
676 Vector<std::pair<Box,Box>>
677 get_src_dst_boxes (DTOS const& dtos, Box const& dstbox, Box const& domain)
678 {
679  Vector<std::pair<Box,Box>> r;
680  IntVect mapped_smallend(dtos(amrex::lbound(dstbox)));
681  IntVect mapped_bigend (dtos(amrex::ubound(dstbox)));
682  if (!domain.contains(mapped_smallend) || !domain.contains(mapped_bigend)) {
683  return r;
684  }
685 
686  auto sign = dtos.sign(amrex::lbound(dstbox));
687  auto perm = dtos.permutation(amrex::lbound(dstbox));
688  auto dtype = dstbox.type();
689  IntVect stype{AMREX_D_DECL(dtype[perm[0]],
690  dtype[perm[1]],
691  dtype[perm[2]])};
692  Array<Array<std::pair<int,int>,2>,AMREX_SPACEDIM> ends;
693  Array<Array<std::pair<int,int>,2>,AMREX_SPACEDIM> dst_ends;
694  Array<int,AMREX_SPACEDIM> nboxes;
695  for (int ddim = 0; ddim < AMREX_SPACEDIM; ++ddim) {
696  int sdim = perm[ddim];
697  auto mm = std::minmax(mapped_smallend[sdim],mapped_bigend[sdim]);
698  if (((sign[ddim] > 0) && (mapped_smallend[sdim] <= mapped_bigend[sdim])) ||
699  ((sign[ddim] < 0) && (mapped_bigend[sdim] <= mapped_smallend[sdim])))
700  {
701  nboxes[sdim] = 1;
702  ends[sdim][0] = mm;
703  dst_ends[ddim][0] = std::make_pair(dstbox.smallEnd(ddim),
704  dstbox.bigEnd(ddim));
705  } else {
706  nboxes[sdim] = 2;
707  ends[sdim][0].first = domain.smallEnd(sdim);
708  ends[sdim][0].second = mm.first;
709  ends[sdim][1].first = mm.second;
710  ends[sdim][1].second = domain.bigEnd(sdim);
711  int n0 = ends[sdim][0].second - ends[sdim][0].first;
712  int n1 = ends[sdim][1].second - ends[sdim][1].first;
713  if (mm.first == mapped_smallend[sdim]) {
714  dst_ends[ddim][0] = std::make_pair(dstbox.smallEnd(ddim),
715  dstbox.smallEnd(ddim)+n0);
716  dst_ends[ddim][1] = std::make_pair(dstbox.bigEnd(ddim)-n1,
717  dstbox.bigEnd(ddim));
718  } else {
719  dst_ends[ddim][0] = std::make_pair(dstbox.bigEnd(ddim)-n0,
720  dstbox.bigEnd(ddim));
721  dst_ends[ddim][1] = std::make_pair(dstbox.smallEnd(ddim),
722  dstbox.smallEnd(ddim)+n1);
723  }
724  }
725  }
726 
727  r.reserve(AMREX_D_TERM(nboxes[0],*nboxes[1],*nboxes[2]));
728 
729 #if (AMREX_SPACEDIM == 3)
730  for (int kbox = 0; kbox < nboxes[2]; ++kbox) {
731 #endif
732 #if (AMREX_SPACEDIM >=2 )
733  for (int jbox = 0; jbox < nboxes[1]; ++jbox) {
734 #endif
735  for (int ibox = 0; ibox < nboxes[0]; ++ibox)
736  {
737  IntVect siv(AMREX_D_DECL(ibox,jbox,kbox));
738  IntVect div(AMREX_D_DECL(siv[perm[0]],siv[perm[1]],siv[perm[2]]));
739  r.emplace_back(Box(IntVect(AMREX_D_DECL(ends[0][ibox].first,
740  ends[1][jbox].first,
741  ends[2][kbox].first)),
742  IntVect(AMREX_D_DECL(ends[0][ibox].second,
743  ends[1][jbox].second,
744  ends[2][kbox].second)),
745  stype),
746  Box(IntVect(AMREX_D_DECL(dst_ends[0][div[0]].first,
747  dst_ends[1][div[1]].first,
748  dst_ends[2][div[2]].first)),
749  IntVect(AMREX_D_DECL(dst_ends[0][div[0]].second,
750  dst_ends[1][div[1]].second,
751  dst_ends[2][div[2]].second)),
752  dtype));
753  AMREX_D_TERM(},},})
754 
755  return r;
756 }
757 
758 template <typename DTOS>
759 Box get_dst_subbox (DTOS const& dtos, std::pair<Box,Box> const& sdboxes,
760  Box const& srcsubbox)
761 {
762  Box const& srcbox = sdboxes.first;
763  Box const& dstbox = sdboxes.second;
764  if (srcbox == srcsubbox) {
765  return dstbox;
766  } else {
767  auto sign = dtos.sign(amrex::lbound(dstbox));
768  auto perm = dtos.permutation(amrex::lbound(dstbox));
769  Box dstsubbox = dstbox;
770  for (int ddim = 0; ddim < AMREX_SPACEDIM; ++ddim) {
771  int sdim = perm[ddim];
772  if (sign[ddim] > 0) {
773  dstsubbox.growLo(ddim, srcbox.smallEnd(sdim)-srcsubbox.smallEnd(sdim));
774  dstsubbox.growHi(ddim, srcsubbox.bigEnd(sdim)-srcbox.bigEnd(sdim));
775  } else {
776  dstsubbox.growLo(ddim, srcsubbox.bigEnd(sdim)-srcbox.bigEnd(sdim));
777  dstsubbox.growHi(ddim, srcbox.smallEnd(sdim)-srcsubbox.smallEnd(sdim));
778  }
779  }
780  return dstsubbox;
781  }
782 }
783 
784 namespace detail {
785  void split_boxes (BoxList& bl, Box const& domain);
786 }
787 
788 template <typename FAB, typename DTOS>
789 std::enable_if_t<IsBaseFab<FAB>() && IsCallableR<Dim3,DTOS,Dim3>(),
790  FabArrayBase::CommMetaData>
791 makeFillBoundaryMetaData (FabArray<FAB>& mf, IntVect const& nghost,
792  Geometry const& geom, DTOS const& dtos)
793 {
794  FabArrayBase::CommMetaData cmd;
795  cmd.m_LocTags = std::make_unique<FabArrayBase::CopyComTagsContainer>();
796  cmd.m_SndTags = std::make_unique<FabArrayBase::MapOfCopyComTagContainers>();
797  cmd.m_RcvTags = std::make_unique<FabArrayBase::MapOfCopyComTagContainers>();
798 
799  // Normal FillBoundary part
800  mf.define_fb_metadata(cmd, nghost, false, geom.periodicity(), false);
801 
802  BoxArray const& ba = mf.boxArray();
803  DistributionMapping const& dm = mf.DistributionMap();
804  Box dombox = amrex::convert(geom.Domain(), ba.ixType());
805  Box pdombox = amrex::convert(geom.growPeriodicDomain(nghost), ba.ixType());
806 
807  const int myproc = ParallelDescriptor::MyProc();
808  const auto nboxes = static_cast<int>(ba.size());
809  std::vector<std::pair<int,Box> > isects;
810 
811  for (int i = 0; i < nboxes; ++i) {
812  Box const& gbx = amrex::grow(ba[i], nghost);
813  BoxList bl = amrex::boxDiff(gbx, pdombox);
814  if (bl.isEmpty()) { continue; }
815 
816  detail::split_boxes(bl, dombox);
817 
818  const int dst_owner = dm[i];
819  for (auto const& dst_box : bl) {
820  auto const& src_dst_boxes = get_src_dst_boxes(dtos, dst_box, dombox);
821  for (auto const& sd_box_pair : src_dst_boxes) {
822  ba.intersections(sd_box_pair.first, isects);
823  for (auto const& is : isects) {
824  int const k = is.first;
825  Box const src_b = is.second;
826  int const src_owner = dm[k];
827  if (dst_owner == myproc || src_owner == myproc) {
828  Box const& dst_b = get_dst_subbox(dtos, sd_box_pair, src_b);
829  if (src_owner == dst_owner) {
830  cmd.m_LocTags->emplace_back(dst_b, src_b, i, k);
831  } else {
832  auto& tags = (dst_owner == myproc) ?
833  (*cmd.m_RcvTags)[src_owner] :
834  (*cmd.m_SndTags)[dst_owner];
835  tags.emplace_back(dst_b, src_b, i, k);
836  }
837  }
838  }
839  }
840  }
841  }
842 
843  return cmd;
844 }
845 
846 struct SphThetaPhiRIndexMapping
847 {
848  SphThetaPhiRIndexMapping (Box const& a_domain)
849  : nx(a_domain.length(0)),
850  ny(a_domain.length(1)),
851  nz(a_domain.length(2))
852  {
853  AMREX_ASSERT(a_domain.smallEnd() == 0);
854  }
855 
856  [[nodiscard]] AMREX_GPU_HOST_DEVICE
857  Dim3 operator() (Dim3 const& ijk) const noexcept
858  {
859  const int i = ijk.x;
860  const int j = ijk.y;
861  const int k = ijk.z;
862  bool ilo = i < 0;
863  bool ihi = i >= nx;
864  bool imd = i >= 0 && i < nx;
865  bool jlo = j < 0;
866  bool jhi = j >= ny;
867  bool jmd = j >= 0 && j < ny;
868  bool klo = k < 0;
869  bool kmd = k >= 0 && k < nz;
870  // We do not need to do anything at the theta-lo/r-lo edge,
871  // theta-hi/r-lo edge, and r > r-hi.
872  if (ilo && jmd && kmd)
873  {
874  return {-1-i, (j+ny/2)%ny, k};
875  }
876  else if (ihi && jmd && kmd)
877  {
878  return {2*nx-1-i, (j+ny/2)%ny, k};
879  }
880  else if (imd && jlo && kmd)
881  {
882  return {i, j+ny, k};
883  }
884  else if (imd && jhi && kmd)
885  {
886  return {i, j-ny, k};
887  }
888  else if (imd && jmd & klo)
889  {
890  return {nx-1-i, (j+ny/2)%ny, -1-k};
891  }
892  else if (ilo && jlo && kmd)
893  {
894  return {-1-i, (j+ny/2)%ny, k};
895  }
896  else if (ihi && jlo && kmd)
897  {
898  return {2*nx-1-i, (j+ny/2)%ny, k};
899  }
900  else if (ilo && jhi && kmd)
901  {
902  return {-1-i, (j+ny/2)%ny, k};
903  }
904  else if (ihi && jhi && kmd)
905  {
906  return {2*nx-1-i, (j+ny/2)%ny, k};
907  }
908  else if (imd && jlo && klo)
909  {
910  return {nx-1-i, (j+ny/2)%ny, -1-k};
911  }
912  else if (imd && jhi && klo)
913  {
914  return {nx-1-i, (j+ny/2)%ny, -1-k};
915  }
916  else
917  {
918  return ijk;
919  }
920  }
921 
922  [[nodiscard]] IntVect sign (Dim3 const& ijk) const noexcept
923  {
924  if (ijk.z < 0) {
925  return IntVect{AMREX_D_DECL(-1, 1,-1)};
926  } else if (ijk.z >=0 && ijk.z < nz &&
927  (ijk.x < 0 || ijk.x >= nx)) {
928  return IntVect{AMREX_D_DECL(-1, 1, 1)};
929  } else {
930  return IntVect{AMREX_D_DECL( 1, 1, 1)};
931  }
932  }
933 
934  [[nodiscard]] IntVect permutation (Dim3 const& /*ijk*/) const noexcept // NOLINT(readability-convert-member-functions-to-static)
935  {
936  return IntVect(AMREX_D_DECL(0,1,2));
937  }
938 
939 private:
940  int nx, ny, nz;
941 };
942 
943 struct SphThetaPhiRComponentMapping
944 {
945  SphThetaPhiRComponentMapping (Box const& a_domain, int a_start_index)
946  : nx(a_domain.length(0)),
947  ny(a_domain.length(1)),
948  nz(a_domain.length(2)),
949  scomp(a_start_index) {}
950 
951  template <typename T>
952  [[nodiscard]] AMREX_GPU_HOST_DEVICE
953  T operator()(Array4<const T> const& a, Dim3 const& ijk, int n) const noexcept
954  {
955  const int i = ijk.x;
956  const int j = ijk.y;
957  const int k = ijk.z;
958  auto r = a(i,j,k,n);
959  if (n == scomp) {
960  if ((i >= 0 && i < nx) &&
961  (j < 0 || j >= ny) &&
962  (k >= 0 && k < nz)) {
963  return r;
964  } else {
965  // We do not need to worry about the theta-lo/r-lo edge,
966  // theta-hi/r-lo edge, and r > r-hi.
967  return -r;
968  }
969  } else if (n == scomp+2) {
970  if (k < 0) {
971  return -r;
972  } else {
973  return r;
974  }
975  } else {
976  return r;
977  }
978  }
979 private:
980  int nx, ny, nz;
981  int scomp;
982 };
983 
984 extern template MultiBlockCommMetaData ParallelCopy(FabArray<FArrayBox>& dest, const Box& destbox,
985  const FabArray<FArrayBox>& src, int destcomp,
986  int srccomp, int numcomp, const IntVect& ngrow,
987  MultiBlockIndexMapping const&, Identity const&);
988 }
989 
990 #endif
#define BL_PROFILE(a)
Definition: AMReX_BLProfiler.H:551
#define BL_ASSERT(EX)
Definition: AMReX_BLassert.H:39
#define AMREX_ASSERT(EX)
Definition: AMReX_BLassert.H:38
#define AMREX_NODISCARD
Definition: AMReX_Extension.H:251
#define AMREX_FORCE_INLINE
Definition: AMReX_Extension.H:119
#define AMREX_HOST_DEVICE_PARALLEL_FOR_4D(...)
Definition: AMReX_GpuLaunch.nolint.H:55
#define AMREX_GPU_DEVICE
Definition: AMReX_GpuQualifiers.H:18
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
Array4< int const > offset
Definition: AMReX_HypreMLABecLap.cpp:1089
#define AMREX_D_TERM(a, b, c)
Definition: AMReX_SPACE.H:129
#define AMREX_D_DECL(a, b, c)
Definition: AMReX_SPACE.H:104
virtual void free(void *pt)=0
A pure virtual function for deleting the arena pointed to by pt.
virtual void * alloc(std::size_t sz)=0
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
CopyComTag::CopyComTagsContainer CopyComTagsContainer
Definition: AMReX_FabArrayBase.H:219
void copyAsync(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:233
static constexpr HostToDevice hostToDevice
Definition: AMReX_GpuContainers.H:98
void streamSynchronize() noexcept
Definition: AMReX_GpuDevice.H:237
bool inLaunchRegion() noexcept
Definition: AMReX_GpuControl.H:86
int MyProc()
Definition: AMReX_MPMD.cpp:117
void split_boxes(BoxList &bl, Box const &domain)
Definition: AMReX_NonLocalBC.cpp:4
Definition: AMReX_NonLocalBC.cpp:3
std::enable_if_t< IsBaseFab< FAB >::value > FillPolar(FabArray< FAB > &mf, int scomp, int ncomp, IntVect const &nghost, Box const &domain)
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >), FabArrayBase::CommMetaData > makeFillBoundaryMetaData(FabArray< FAB > &mf, IntVect const &nghost, Geometry const &geom, DTOS const &dtos)
Make metadata for FillBoundary.
std::enable_if_t< IsBaseFab< FAB >::value > Rotate90(FabArray< FAB > &mf, int scomp, int ncomp, IntVect const &nghost, Box const &domain)
std::enable_if_t< IsBaseFab< FAB >::value > Rotate180(FabArray< FAB > &mf, int scomp, int ncomp, IntVect const &nghost, Box const &domain)
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >) &&IsFabProjection< Proj, FAB >)> local_copy_cpu(FabArray< FAB > &dest, const FabArray< FAB > &src, int dcomp, int scomp, int ncomp, FabArrayBase::CopyComTagsContainer const &local_tags, DTOS const &dtos=DTOS{}, Proj const &proj=Proj{}) noexcept
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >) &&IsFabProjection< Proj, FAB >)> unpack_recv_buffer_cpu(FabArray< FAB > &mf, int dcomp, int ncomp, Vector< char * > const &recv_data, Vector< std::size_t > const &recv_size, Vector< FabArrayBase::CopyComTagsContainer const * > const &recv_cctc, DTOS const &dtos=DTOS{}, Proj const &proj=Proj{}) noexcept
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >) &&IsFabProjection< Proj, FAB >)> local_copy_gpu(FabArray< FAB > &dest, const FabArray< FAB > &src, int dcomp, int scomp, int ncomp, FabArrayBase::CopyComTagsContainer const &local_tags, DTOS const &dtos=DTOS{}, Proj const &proj=Proj{}) noexcept
std::enable_if_t< IsBaseFab< FAB >::value > PrepareSendBuffers(const PackComponents &components, FabArray< FAB > &dest, const FabArray< FAB > &src, CommData &comm, const FabArrayBase::MapOfCopyComTagContainers &cctc)
Calls PrepareComBuffers.
Definition: AMReX_NonLocalBC.H:555
std::enable_if_t< IsCallableR< Dim3, DTOS, Dim3 >::value &&!IsCallableR< IndexType, DTOS, IndexType >::value, Box > Image(DTOS const &dtos, const Box &box)
Applies the Dim3 to Dim3 mapping onto Boxes but does not change the index type.
Definition: AMReX_NonLocalBC.H:118
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >) &&IsFabProjection< Proj, FAB >), CommHandler > FillBoundary_nowait(FabArray< FAB > &mf, const FabArrayBase::CommMetaData &cmd, int scomp, int ncomp, DTOS const &dtos, Proj const &proj=Proj{})
Start communication to fill boundary.
std::enable_if_t< HasInverseMemFn< DTOS >::value &&!IsCallableR< IndexType, DTOS, IndexType >::value, Box > InverseImage(DTOS const &dtos, const Box &box)
Applies the inverse Dim3 to Dim3 mapping onto Boxes without changing the index type.
Definition: AMReX_NonLocalBC.H:178
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >) &&IsFabProjection< Proj, FAB >)> unpack_recv_buffer_gpu(FabArray< FAB > &mf, int scomp, int ncomp, Vector< char * > const &recv_data, Vector< std::size_t > const &recv_size, Vector< FabArrayBase::CopyComTagsContainer const * > const &recv_cctc, DTOS const &dtos=DTOS{}, Proj const &proj=Proj{})
std::enable_if_t< IsBaseFab< FAB >) &&IsCallableR< Dim3, DTOS, Dim3 >) &&IsFabProjection< Proj, FAB >)> FillBoundary_finish(CommHandler handler, FabArray< FAB > &mf, const FabArrayBase::CommMetaData &cmd, int scomp, int ncomp, DTOS const &dtos, Proj const &proj=Proj{})
Finish communication started by FillBoundary_nowait.
int NProcsSub() noexcept
number of ranks in current frame
Definition: AMReX_ParallelContext.H:74
bool UseGpuAwareMpi()
Definition: AMReX_ParallelDescriptor.H:111
void Waitall(Vector< MPI_Request > &, Vector< MPI_Status > &)
Definition: AMReX_ParallelDescriptor.cpp:1295
int SeqNum() noexcept
Returns sequential message sequence numbers, usually used as tags for send/recv.
Definition: AMReX_ParallelDescriptor.H:613
double second() noexcept
Returns wall-clock seconds since start of execution.
Definition: AMReX_ParallelDescriptor.cpp:1285
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
DistributionMapping const & DistributionMap(FabArrayBase const &fa)
BoxND< AMREX_SPACEDIM > Box
Definition: AMReX_BaseFwd.H:27
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Array4< T > makeArray4(T *p, Box const &bx, int ncomp) noexcept
Definition: AMReX_BaseFab.H:87
AMREX_ATTRIBUTE_FLATTEN_FOR void LoopConcurrentOnCpu(Dim3 lo, Dim3 hi, F const &f) noexcept
Definition: AMReX_Loop.H:378
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > grow(const BoxND< dim > &b, int i) noexcept
Grow BoxND in all directions by given amount.
Definition: AMReX_Box.H:1211
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > convert(const BoxND< dim > &b, const IntVectND< dim > &typ) noexcept
Returns a BoxND with different type.
Definition: AMReX_Box.H:1435
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 ubound(Array4< T > const &a) noexcept
Definition: AMReX_Array4.H:315
BoxArray const & boxArray(FabArrayBase const &fa)
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 lbound(Array4< T > const &a) noexcept
Definition: AMReX_Array4.H:308
BoxList boxDiff(const Box &b1in, const Box &b2)
Returns BoxList defining the compliment of b2 in b1in.
IntVectND< AMREX_SPACEDIM > IntVect
Definition: AMReX_BaseFwd.H:30
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void ignore_unused(const Ts &...)
This shuts up the compiler about unused variables.
Definition: AMReX.H:111
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 length(Array4< T > const &a) noexcept
Definition: AMReX_Array4.H:322
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
Arena * The_Arena()
Definition: AMReX_Arena.cpp:609
Definition: AMReX_FabArrayCommI.H:896