Block-Structured AMR Software Framework
AMReX_NeighborParticlesI.H
Go to the documentation of this file.
1 
2 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
3 bool NeighborParticleContainer<NStructReal, NStructInt, NArrayReal, NArrayInt>::use_mask = false;
4 
5 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
6 bool NeighborParticleContainer<NStructReal, NStructInt, NArrayReal, NArrayInt>::enable_inverse = false;
7 
8 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
9 NeighborParticleContainer<NStructReal, NStructInt, NArrayReal, NArrayInt>
10 ::NeighborParticleContainer (ParGDBBase* gdb, int ncells)
11  : ParticleContainer<NStructReal, NStructInt, NArrayReal, NArrayInt> (gdb),
12  m_num_neighbor_cells(ncells)
13 {
15 }
16 
17 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
20  const DistributionMapping & dmap,
21  const BoxArray & ba,
22  int nneighbor)
23  : ParticleContainer<NStructReal, NStructInt, NArrayReal, NArrayInt> (geom, dmap, ba),
24  m_num_neighbor_cells(nneighbor)
25 {
26  initializeCommComps();
27 }
28 
29 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
32  const Vector<DistributionMapping> & dmap,
33  const Vector<BoxArray> & ba,
34  const Vector<int> & rr,
35  int nneighbor)
36  : ParticleContainer<NStructReal, NStructInt, NArrayReal, NArrayInt> (geom, dmap, ba, rr),
37  m_num_neighbor_cells(nneighbor)
38 {
39  initializeCommComps();
40 }
41 
42 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
43 void
46  for (int ii = 0; ii < AMREX_SPACEDIM + NStructReal + this->NumRealComps(); ++ii) {
47  ghost_real_comp.push_back(1);
48  }
49  for (int ii = 0; ii < 2 + NStructInt + this->NumIntComps(); ++ii) {
50  ghost_int_comp.push_back(1);
51  }
52  calcCommSize();
53 }
54 
55 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
56 void
58 ::setRealCommComp (int i, bool value) {
59  ghost_real_comp[i] = value;
60  calcCommSize();
61 }
62 
63 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
64 void
66 ::setIntCommComp (int i, bool value) {
67  ghost_int_comp[i] = value;
68  calcCommSize();
69 }
70 
71 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
72 void
75  size_t comm_size = 0;
76  for (int ii = 0; ii < AMREX_SPACEDIM + NStructReal + this->NumRealComps(); ++ii) {
77  if (ghost_real_comp[ii]) {
78  comm_size += sizeof(typename ParticleType::RealType);
79  }
80  }
81  for (int ii = 0; ii < 2 + NStructInt + this->NumIntComps(); ++ii) {
82  if (ghost_int_comp[ii]) {
83  comm_size += sizeof(int);
84  }
85  }
86  if ( enableInverse() ) { comm_size += 4*sizeof(int); }
87  cdata_size = comm_size;
88 }
89 
90 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
91 void
93 ::Regrid (const DistributionMapping &dmap, const BoxArray &ba ) {
94  const int lev = 0;
95  AMREX_ASSERT(this->finestLevel() == 0);
96  this->SetParticleBoxArray(lev, ba);
97  this->SetParticleDistributionMap(lev, dmap);
98  this->Redistribute();
99 }
100 
101 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
102 void
104 ::Regrid (const DistributionMapping &dmap, const BoxArray &ba, int lev) {
105  AMREX_ASSERT(lev <= this->finestLevel());
106  this->SetParticleBoxArray(lev, ba);
107  this->SetParticleDistributionMap(lev, dmap);
108  this->Redistribute();
109 }
110 
111 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
112 void
115  AMREX_ASSERT(ba.size() == this->finestLevel()+1);
116  for (int lev = 0; lev < this->numLevels(); ++lev)
117  {
118  this->SetParticleBoxArray(lev, ba[lev]);
119  this->SetParticleDistributionMap(lev, dmap[lev]);
120  }
121  this->Redistribute();
122 }
123 
124 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
125 bool
128 
129  BL_PROFILE("NeighborParticleContainer::areMasksValid");
130 
131  resizeContainers(this->numLevels());
132 
133  for (int lev = 0; lev < this->numLevels(); ++lev)
134  {
135  BoxArray ba = this->ParticleBoxArray(lev);
136  const DistributionMapping& dmap = this->ParticleDistributionMap(lev);
137 
138  if (mask_ptr[lev] == nullptr ||
139  ! BoxArray::SameRefs(mask_ptr[lev]->boxArray(), ba) ||
140  ! DistributionMapping::SameRefs(mask_ptr[lev]->DistributionMap(), dmap))
141  {
142  return false;
143  }
144  }
145  return true;
146 }
147 
148 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
149 void
152 
153  BL_PROFILE("NeighborParticleContainer::BuildMasks");
154 
155  if (this->numLevels() == 1) { use_mask = true; }
156  else { use_mask = false; }
157 
158  resizeContainers(this->numLevels());
159 
160  for (int lev = 0; lev < this->numLevels(); ++lev)
161  {
162  BoxArray ba = this->ParticleBoxArray(lev);
163  const DistributionMapping& dmap = this->ParticleDistributionMap(lev);
164 
165  const Geometry& geom = this->Geom(lev);
166 
167  mask_ptr[lev] = std::make_unique<iMultiFab>(ba, dmap, int(num_mask_comps), m_num_neighbor_cells);
168  mask_ptr[lev]->setVal(-1, m_num_neighbor_cells);
169 
170 #ifdef AMREX_USE_OMP
171 #pragma omp parallel
172 #endif
173  for (MFIter mfi(*mask_ptr[lev],this->do_tiling ? this->tile_size : IntVect::TheZeroVector());
174  mfi.isValid(); ++mfi) {
175  const Box& box = mfi.tilebox();
176  const int grid_id = mfi.index();
177  const int tile_id = mfi.LocalTileIndex();
178  (*mask_ptr[lev])[mfi].template setVal<RunOn::Host>(grid_id, box, MaskComps::grid, 1);
179  (*mask_ptr[lev])[mfi].template setVal<RunOn::Host>(tile_id, box, MaskComps::tile, 1);
180  (*mask_ptr[lev])[mfi].template setVal<RunOn::Host>(lev , box, MaskComps::level, 1);
181  }
182 
183  mask_ptr[lev]->FillBoundary(geom.periodicity());
184  }
185 }
186 
187 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
188 void
191 {
192  BL_PROFILE("NeighborParticleContainer::GetNeighborCommTags");
193 
194  local_neighbors.clear();
195  neighbor_procs.clear();
196 
197  if (use_mask)
198  {
199  AMREX_ASSERT(this->finestLevel() == 0);
200  const int lev = 0;
201  for (MFIter mfi(*mask_ptr[lev],this->do_tiling ? this->tile_size : IntVect::TheZeroVector());
202  mfi.isValid(); ++mfi) {
203  const Box& box = mfi.growntilebox();
204  for (IntVect iv = box.smallEnd(); iv <= box.bigEnd(); box.next(iv)) {
205  const int grid = (*mask_ptr[lev])[mfi](iv, MaskComps::grid);
206  if (grid >= 0) {
207  const int tile = (*mask_ptr[lev])[mfi](iv, MaskComps::tile);
208  const int level = (*mask_ptr[lev])[mfi](iv, MaskComps::level);
209  const int global_proc = this->ParticleDistributionMap(level)[grid];
210  const int proc = ParallelContext::global_to_local_rank(global_proc);
211  NeighborCommTag comm_tag(proc, level, grid, tile);
212  local_neighbors.push_back(comm_tag);
213  if (proc != ParallelContext::MyProcSub()) {
214  neighbor_procs.push_back(proc);
215  }
216  }
217  }
218  }
219  }
220  else
221  {
222  for (int lev = 0; lev < this->numLevels(); ++lev)
223  {
224  for (MFIter mfi(*mask_ptr[lev],this->do_tiling ? this->tile_size : IntVect::TheZeroVector());
225  mfi.isValid(); ++mfi) {
226  const Box& box = mfi.validbox();
227  Vector<NeighborCommTag> comm_tags;
228  GetCommTagsBox(comm_tags, lev, box);
229  for (auto const& tag : comm_tags) {
230  local_neighbors.push_back(tag);
231  if (tag.proc_id != ParallelContext::MyProcSub()) {
232  neighbor_procs.push_back(tag.proc_id);
233  }
234  }
235  }
236  }
237  }
238 
239  RemoveDuplicates(local_neighbors);
240  RemoveDuplicates(neighbor_procs);
241 }
242 
243 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
244 IntVect
246 ::computeRefFac (int src_lev, int lev)
247 {
248  IntVect ref_fac(1);
249  if (src_lev < lev) {
250  for (int l = src_lev; l < lev; ++l) {
251  ref_fac *= this->GetParGDB()->refRatio(l);
252  }
253  } else if (src_lev > lev) {
254  for (int l = src_lev; l > lev; --l) {
255  ref_fac *= this->GetParGDB()->refRatio(l-1);
256  }
257  ref_fac *= -1;
258  }
259  return ref_fac;
260 }
261 
262 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
263 void
265 ::GetCommTagsBox (Vector<NeighborCommTag>& tags, int src_lev, const Box& in_box)
266 {
267  std::vector< std::pair<int, Box> > isects;
268  Box tbx;
269 
270  for (int lev = 0; lev < this->numLevels(); ++lev) {
271  Box box = in_box;
272  const IntVect& ref_fac = computeRefFac(src_lev, lev);
273  if (ref_fac < IntVect::TheZeroVector())
274  {
275  box.coarsen(-1*ref_fac);
276  }
277  else if (ref_fac > IntVect::TheZeroVector())
278  {
279  box.refine(ref_fac);
280  }
281  box.grow(computeRefFac(0, lev)*m_num_neighbor_cells);
282  const Periodicity& periodicity = this->Geom(lev).periodicity();
283  const std::vector<IntVect>& pshifts = periodicity.shiftIntVect();
284  const BoxArray& ba = this->ParticleBoxArray(lev);
285 
286  for (auto const& pshift : pshifts)
287  {
288  const Box& pbox = box + pshift;
289  bool first_only = false;
290  ba.intersections(pbox, isects, first_only, 0);
291  for (const auto& isec : isects) {
292  const int grid = isec.first;
293  const int global_proc = this->ParticleDistributionMap(lev)[grid];
294  const int proc = ParallelContext::global_to_local_rank(global_proc);
295  for (IntVect iv = pbox.smallEnd(); iv <= pbox.bigEnd(); pbox.next(iv))
296  {
297  if (ba[grid].contains(iv))
298  {
299  int tile = getTileIndex(iv, ba[grid],
300  this->do_tiling, this->tile_size, tbx);
301  tags.push_back(NeighborCommTag(proc, lev, grid, tile));
302  }
303  }
304  }
305  }
306  }
307 }
308 
309 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
310 void
313 
314  BL_PROFILE("NeighborParticleContainer::cacheNeighborInfo");
315 
316  AMREX_ASSERT(this->OK());
317 
318  resizeContainers(this->numLevels());
319 
320  clearNeighbors();
321 
322  AMREX_ASSERT(hasNeighbors() == false);
323 
324  const int MyProc = ParallelContext::MyProcSub();
325 
327  std::map<NeighborCommTag, Vector<NeighborIndexMap> > remote_map;
328 
329  // tmp data structures used for OMP reduction
331  std::map<NeighborCommTag, Vector<Vector<NeighborIndexMap> > > tmp_remote_map;
332 
333  local_map.resize(this->numLevels());
334  tmp_local_map.resize(this->numLevels());
335 
336  int num_threads = OpenMP::get_max_threads();
337 
338  for (int lev = 0; lev < this->numLevels(); ++lev) {
339  // resize our temporaries in serial
340  for (int i = 0; i < static_cast<int>(local_neighbors.size()); ++i) {
341  const NeighborCommTag& comm_tag = local_neighbors[i];
342  tmp_remote_map[comm_tag].resize(num_threads);
343  remote_map[comm_tag];
344  PairIndex index(comm_tag.grid_id, comm_tag.tile_id);
345  tmp_local_map[lev][index].resize(num_threads);
346  local_map[lev][index];
347  buffer_tag_cache[lev][index].resize(num_threads);
348  }
349  }
350 
351  for (int lev = 0; lev < this->numLevels(); ++lev) {
352  // First pass - each thread collects the NeighborIndexMaps it owes to other
353  // grids / tiles / procs
354 #ifdef AMREX_USE_OMP
355 #pragma omp parallel
356 #endif
357  {
359  tags.reserve(AMREX_D_TERM(3, *3, *3));
360  for (MyParIter pti(*this, lev); pti.isValid(); ++pti) {
361  int thread_num = OpenMP::get_thread_num();
362  const int& grid = pti.index();
363  const int& tile = pti.LocalTileIndex();
364  PairIndex src_index(grid, tile);
365 
366  NeighborCopyTag src_tag(lev, grid, tile);
367 
368  auto& cache = buffer_tag_cache[lev][src_index][thread_num];
369 
370  auto& particles = pti.GetArrayOfStructs();
371  for (int i = 0; i < pti.numParticles(); ++i) {
372  const ParticleType& p = particles[i];
373 
374  getNeighborTags(tags, p, m_num_neighbor_cells, src_tag, pti);
375 
376  // Add neighbors to buffers
377  for (int j = 0; j < static_cast<int>(tags.size()); ++j) {
378  NeighborCopyTag& tag = tags[j];
379  PairIndex dst_index(tag.grid, tag.tile);
380  if (tag.grid < 0) { continue; }
381 
382  tag.src_index = i;
383  const int cache_index = cache.size();
384  cache.push_back(tag);
385 
386  const int global_who = this->ParticleDistributionMap(tag.level)[tag.grid];
387  const int who = ParallelContext::global_to_local_rank(global_who);
388  NeighborIndexMap nim(tag.level, dst_index.first, dst_index.second, -1,
389  lev, src_index.first, src_index.second,
390  cache_index, thread_num);
391  if (who == MyProc) {
392  auto& tmp = tmp_local_map[tag.level][dst_index];
393  Vector<NeighborIndexMap>& buffer = tmp[thread_num];
394  buffer.push_back(nim);
395  } else {
396  NeighborCommTag comm_tag(who, tag.level, tag.grid, tag.tile);
397  Vector<NeighborIndexMap>& buffer = tmp_remote_map[comm_tag][thread_num];
398  buffer.push_back(nim);
399  }
400  }
401  tags.clear();
402  }
403  }
404  }
405  }
406 
407  for (int lev = 0; lev < this->numLevels(); ++lev) {
408  // second pass - for each tile, collect the neighbors owed from all threads
409 #ifdef AMREX_USE_OMP
410 #pragma omp parallel
411 #endif
412  for (MFIter mfi = this->MakeMFIter(lev); mfi.isValid(); ++mfi) {
413  const int grid = mfi.index();
414  const int tile = mfi.LocalTileIndex();
415  PairIndex index(grid, tile);
416  for (int i = 0; i < num_threads; ++i) {
417  local_map[lev][index].insert(local_map[lev][index].end(),
418  tmp_local_map[lev][index][i].begin(),
419  tmp_local_map[lev][index][i].end());
420  tmp_local_map[lev][index][i].erase(tmp_local_map[lev][index][i].begin(),
421  tmp_local_map[lev][index][i].end());
422  }
423  }
424  }
425 
426  // do the same for the remote neighbors
427  typename std::map<NeighborCommTag, Vector<Vector<NeighborIndexMap> > >::iterator it;
428 #ifdef AMREX_USE_OMP
429 #pragma omp parallel
430 #pragma omp single nowait
431 #endif
432  for (it=tmp_remote_map.begin(); it != tmp_remote_map.end(); it++) {
433 #ifdef AMREX_USE_OMP
434 #pragma omp task firstprivate(it)
435 #endif
436  {
437  const NeighborCommTag& tag = it->first;
438  Vector<Vector<NeighborIndexMap> >& tmp = it->second;
439  for (int i = 0; i < num_threads; ++i) {
440  remote_map[tag].insert(remote_map[tag].end(), tmp[i].begin(), tmp[i].end());
441  tmp[i].erase(tmp[i].begin(), tmp[i].end());
442  }
443  }
444  }
445 
446  for (int lev = 0; lev < this->numLevels(); ++lev) {
447  // now for the local neighbors, allocate buffers and cache
448  for (MFIter mfi = this->MakeMFIter(lev); mfi.isValid(); ++mfi) {
449  const int grid = mfi.index();
450  const int tile = mfi.LocalTileIndex();
451  PairIndex dst_index(grid, tile);
452  const Vector<NeighborIndexMap>& map = local_map[lev][dst_index];
453  const int num_ghosts = map.size();
454  neighbors[lev][dst_index].define(this->NumRuntimeRealComps(),
455  this->NumRuntimeIntComps());
456  neighbors[lev][dst_index].resize(num_ghosts);
457  local_neighbor_sizes[lev][dst_index] = neighbors[lev][dst_index].size();
458  }
459  }
460 
461  for (int lev = 0; lev < this->numLevels(); ++lev) {
462  for (MFIter mfi = this->MakeMFIter(lev); mfi.isValid(); ++mfi) {
463  const int grid = mfi.index();
464  const int tile = mfi.LocalTileIndex();
465  PairIndex dst_index(grid, tile);
466  const Vector<NeighborIndexMap>& map = local_map[lev][dst_index];
467  const int num_ghosts = map.size();
468 #ifdef AMREX_USE_OMP
469 #pragma omp parallel for
470 #endif
471  for (int i = 0; i < num_ghosts; ++i) {
472  const NeighborIndexMap& nim = map[i];
473  PairIndex src_index(nim.src_grid, nim.src_tile);
474  Vector<NeighborCopyTag>& tags = buffer_tag_cache[nim.src_level][src_index][nim.thread_num];
475  AMREX_ASSERT(nim.src_index < tags.size());
476  tags[nim.src_index].dst_index = i;
477  AMREX_ASSERT(size_t(tags[nim.src_index].dst_index) < neighbors[nim.dst_level][dst_index].size());
478  }
479  }
480  }
481 
482  // now we allocate the send buffers and cache the remotes
483  std::map<int, int> tile_counts;
484  for (const auto& kv: remote_map) {
485  tile_counts[kv.first.proc_id] += 1;
486  }
487 
488  for (const auto& kv: remote_map) {
489  if (kv.first.proc_id == MyProc) { continue; }
490  Vector<char>& buffer = send_data[kv.first.proc_id];
491  buffer.resize(sizeof(int));
492  std::memcpy(buffer.data(), &tile_counts[kv.first.proc_id], sizeof(int));
493  }
494 
495  for (auto& kv : remote_map) {
496  if (kv.first.proc_id == MyProc) { continue; }
497  int np = kv.second.size();
498  int data_size = np * cdata_size;
499  Vector<char>& buffer = send_data[kv.first.proc_id];
500  size_t old_size = buffer.size();
501  size_t new_size = buffer.size() + 4*sizeof(int) + data_size;
502  buffer.resize(new_size);
503  char* dst = &buffer[old_size];
504  std::memcpy(dst, &(kv.first.level_id), sizeof(int)); dst += sizeof(int);
505  std::memcpy(dst, &(kv.first.grid_id ), sizeof(int)); dst += sizeof(int);
506  std::memcpy(dst, &(kv.first.tile_id ), sizeof(int)); dst += sizeof(int);
507  std::memcpy(dst, &data_size, sizeof(int)); dst += sizeof(int);
508  size_t buffer_offset = old_size + 4*sizeof(int);
509 #ifdef AMREX_USE_OMP
510 #pragma omp parallel for
511 #endif
512  for (int i = 0; i < np; ++i) {
513  const NeighborIndexMap& nim = kv.second[i];
514  PairIndex src_index(nim.src_grid, nim.src_tile);
515  Vector<NeighborCopyTag>& tags = buffer_tag_cache[nim.src_level][src_index][nim.thread_num];
516  tags[nim.src_index].dst_index = buffer_offset + i*cdata_size;
517  }
518  }
519 
520  if ( enableInverse() )
521  {
522  for (int lev = 0; lev < this->numLevels(); ++lev)
523  {
524  for (const auto& kv : neighbors[lev])
525  {
526  inverse_tags[lev][kv.first].resize(kv.second.size());
527  }
528  }
529  }
530 }
531 
532 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
533 void
536  int nGrow, const NeighborCopyTag& src_tag, const MyParIter& pti)
537 {
538  getNeighborTags(tags, p, IntVect(AMREX_D_DECL(nGrow, nGrow, nGrow)), src_tag, pti);
539 }
540 
541 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
542 void
545  const IntVect& nGrow, const NeighborCopyTag& src_tag, const MyParIter& pti)
546 {
547  Box shrink_box = pti.tilebox();
548  shrink_box.grow(-nGrow);
549 
550  if (use_mask) {
551  const BaseFab<int>& mask = (*mask_ptr[src_tag.level])[src_tag.grid];
552  AMREX_ASSERT(this->finestLevel() == 0);
553  AMREX_ASSERT(src_tag.level == 0);
554 
555  const int lev = 0;
556  const IntVect& iv = this->Index(p, lev);
557  if (shrink_box.contains(iv)) { return; }
558 
559  const Periodicity& periodicity = this->Geom(lev).periodicity();
560  const Box& domain = this->Geom(lev).Domain();
561  const IntVect& lo = domain.smallEnd();
562  const IntVect& hi = domain.bigEnd();
563 
564  // Figure out all our neighbors, removing duplicates
565  AMREX_D_TERM(
566  for (int ii = -nGrow[0]; ii < nGrow[0] + 1; ii += nGrow[0]) {,
567  for (int jj = -nGrow[1]; jj < nGrow[1] + 1; jj += nGrow[1]) {,
568  for (int kk = -nGrow[2]; kk < nGrow[2] + 1; kk += nGrow[2]) {)
569  if (AMREX_D_TERM((ii == 0), && (jj == 0), && (kk == 0))) { continue; }
570  IntVect shift(AMREX_D_DECL(ii, jj, kk));
571  IntVect neighbor_cell = iv + shift;
572 
573  NeighborCopyTag tag;
574  tag.grid = mask(neighbor_cell, MaskComps::grid);
575  tag.tile = mask(neighbor_cell, MaskComps::tile);
576  tag.level = mask(neighbor_cell, MaskComps::level);
577  if (periodicity.isAnyPeriodic()) {
578  for (int dir = 0; dir < AMREX_SPACEDIM; ++dir) {
579  if (! periodicity.isPeriodic(dir)) { continue; }
580  if (neighbor_cell[dir] < lo[dir]) {
581  tag.periodic_shift[dir] = -1;
582  } else if (neighbor_cell[dir] > hi[dir]) {
583  tag.periodic_shift[dir] = 1;
584  }
585  }
586  }
587 
588  if (tag != src_tag) { tags.push_back(tag); }
589 
590  AMREX_D_TERM(
591  },
592  },
593  })
594 
595  RemoveDuplicates(tags);
596  return;
597  }
598  else
599  {
600  std::vector< std::pair<int, Box> > isects;
601  Box tbx;
602  for (int lev = 0; lev < this->numLevels(); ++lev)
603  {
604  IntVect ref_fac = computeRefFac(0, lev);
605  const Periodicity& periodicity = this->Geom(lev).periodicity();
606  const std::vector<IntVect>& pshifts = periodicity.shiftIntVect();
607  const BoxArray& ba = this->ParticleBoxArray(lev);
608  const IntVect& iv = this->Index(p, lev);
609  for (auto const& pshift : pshifts)
610  {
611  Box pbox = amrex::grow(Box(iv, iv), ref_fac*nGrow) + pshift;
612  bool first_only = false;
613  ba.intersections(pbox, isects, first_only, 0);
614  for (const auto& isec : isects)
615  {
616  const Box& grid_box = ba[isec.first];
617  for (IntVect cell = pbox.smallEnd(); cell <= pbox.bigEnd(); pbox.next(cell)) {
618  if ( !grid_box.contains(cell) ) { continue; }
619  int tile = getTileIndex(cell, grid_box,
620  this->do_tiling, this->tile_size, tbx);
621  auto nbor = NeighborCopyTag(lev, isec.first, tile);
622  nbor.periodic_shift = -pshift;
623  if (src_tag != nbor) { tags.push_back(nbor); }
624  }
625  }
626  }
627  }
628 
629  RemoveDuplicates(tags);
630  return;
631  }
632 }
633 
634 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
635 void
638 #ifdef AMREX_USE_GPU
639  fillNeighborsGPU();
640 #else
641  fillNeighborsCPU();
642 #endif
643  m_has_neighbors = true;
644 }
645 
646 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
647 void
649 ::sumNeighbors (int real_start_comp, int real_num_comp,
650  int int_start_comp, int int_num_comp) {
651 #ifdef AMREX_USE_GPU
652  amrex::ignore_unused(real_start_comp,real_num_comp,int_start_comp,int_num_comp);
653  amrex::Abort("Not implemented.");
654 #else
655  sumNeighborsCPU(real_start_comp, real_num_comp, int_start_comp, int_num_comp);
656 #endif
657 }
658 
659 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
660 void
662 ::updateNeighbors (bool boundary_neighbors_only)
663 {
664 
665  AMREX_ASSERT(hasNeighbors());
666 
667 #ifdef AMREX_USE_GPU
668  updateNeighborsGPU(boundary_neighbors_only);
669 #else
670  amrex::ignore_unused(boundary_neighbors_only);
671  updateNeighborsCPU(true);
672 #endif
673  m_has_neighbors = true;
674 }
675 
676 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
677 void
680 {
681 #ifdef AMREX_USE_GPU
682  clearNeighborsGPU();
683 #else
684  clearNeighborsCPU();
685 #endif
686  m_has_neighbors = false;
687 }
688 
689 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
690 template <class CheckPair>
691 void
693 buildNeighborList (CheckPair const& check_pair, bool /*sort*/)
694 {
695  AMREX_ASSERT(numParticlesOutOfRange(*this, m_num_neighbor_cells) == 0);
696 
697  BL_PROFILE("NeighborParticleContainer::buildNeighborList");
698 
699  resizeContainers(this->numLevels());
700 
701  for (int lev = 0; lev < this->numLevels(); ++lev)
702  {
703  m_neighbor_list[lev].clear();
704 
705  for (MyParIter pti(*this, lev); pti.isValid(); ++pti) {
706  PairIndex index(pti.index(), pti.LocalTileIndex());
707  m_neighbor_list[lev][index];
708  }
709 
710 #ifndef AMREX_USE_GPU
711  neighbor_list[lev].clear();
712  for (MyParIter pti(*this, lev); pti.isValid(); ++pti) {
713  PairIndex index(pti.index(), pti.LocalTileIndex());
714  neighbor_list[lev][index];
715  }
716 #endif
717 
718  auto& plev = this->GetParticles(lev);
719  const auto& geom = this->Geom(lev);
720 
721 #ifdef AMREX_USE_OMP
722 #pragma omp parallel if (Gpu::notInLaunchRegion())
723 #endif
724  for (MyParIter pti(*this, lev); pti.isValid(); ++pti)
725  {
726  int gid = pti.index();
727  int tid = pti.LocalTileIndex();
728  auto index = std::make_pair(gid, tid);
729 
730  auto& ptile = plev[index];
731 
732  if (ptile.numParticles() == 0) { continue; }
733 
734  Box bx = pti.tilebox();
735  int ng = computeRefFac(0, lev).max()*m_num_neighbor_cells;
736  bx.grow(ng);
737 
738  Gpu::DeviceVector<int> off_bins_v;
743 
744  off_bins_v.push_back(0);
745  off_bins_v.push_back(int(bx.numPts()));
746  lo_v.push_back(lbound(bx));
747  hi_v.push_back(ubound(bx));
748  dxi_v.push_back(geom.InvCellSizeArray());
749  plo_v.push_back(geom.ProbLoArray());
750 
751  m_neighbor_list[lev][index].build(ptile,
752  check_pair,
753  off_bins_v, dxi_v, plo_v, lo_v, hi_v, ng);
754 
755 #ifndef AMREX_USE_GPU
756  const auto& counts = m_neighbor_list[lev][index].GetCounts();
757  const auto& list = m_neighbor_list[lev][index].GetList();
758 
759  int li = 0;
760  for (int i = 0; i < ptile.numParticles(); ++i)
761  {
762  auto cnt = counts[i];
763  neighbor_list[lev][index].push_back(cnt);
764  for (size_t j = 0; j < cnt; ++j)
765  {
766  neighbor_list[lev][index].push_back(list[li++]+1);
767  }
768  }
769 #endif
770  }
771  }
772 }
773 
774 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
775 template <class CheckPair, class OtherPCType>
776 void
778 buildNeighborList (CheckPair const& check_pair, OtherPCType& other,
779  Vector<std::map<std::pair<int, int>, amrex::NeighborList<typename OtherPCType::ParticleType> > >& neighbor_lists,
780  bool /*sort*/)
781 {
782  BL_PROFILE("NeighborParticleContainer::buildNeighborList");
783 
784  AMREX_ASSERT(numParticlesOutOfRange(*this, m_num_neighbor_cells) == 0);
785  AMREX_ASSERT(numParticlesOutOfRange(other, m_num_neighbor_cells) == 0);
786  AMREX_ASSERT(SameIteratorsOK(*this, other));
787 
788  EnsureThreadSafeTiles(*this);
789  EnsureThreadSafeTiles(other);
790 
791  resizeContainers(this->numLevels());
792  neighbor_lists.resize(this->numLevels());
793 
794  for (int lev = 0; lev < this->numLevels(); ++lev)
795  {
796  neighbor_lists[lev].clear();
797 
798  for (MyParIter pti(*this, lev); pti.isValid(); ++pti) {
799  PairIndex index(pti.index(), pti.LocalTileIndex());
800  neighbor_lists[lev][index];
801  }
802 
803  auto& plev = this->GetParticles(lev);
804  const auto& geom = this->Geom(lev);
805 
806 #ifdef AMREX_USE_OMP
807 #pragma omp parallel if (Gpu::notInLaunchRegion())
808 #endif
809  for (MyParIter pti(*this, lev); pti.isValid(); ++pti)
810  {
811  int gid = pti.index();
812  int tid = pti.LocalTileIndex();
813  auto index = std::make_pair(gid, tid);
814 
815  const auto& ptile = plev[index];
816  auto& other_ptile = other.ParticlesAt(lev, pti);
817  if (ptile.numParticles() == 0) { continue; }
818 
819  Box bx = pti.tilebox();
820  int ng = computeRefFac(0, lev).max()*m_num_neighbor_cells;
821  bx.grow(ng);
822 
823  Gpu::DeviceVector<int> off_bins_v;
828 
829  off_bins_v.push_back(0);
830  off_bins_v.push_back(int(bx.numPts()));
831  lo_v.push_back(lbound(bx));
832  hi_v.push_back(ubound(bx));
833  dxi_v.push_back(geom.InvCellSizeArray());
834  plo_v.push_back(geom.ProbLoArray());
835 
836  neighbor_lists[lev][index].build(ptile, other_ptile,
837  check_pair,
838  off_bins_v, dxi_v, plo_v, lo_v, hi_v, ng);
839  }
840  }
841 }
842 
843 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
844 template <class CheckPair>
845 void
847 buildNeighborList (CheckPair const& check_pair, int type_ind, int* ref_ratio,
848  int num_bin_types, bool /*sort*/)
849 {
850  AMREX_ASSERT(numParticlesOutOfRange(*this, m_num_neighbor_cells) == 0);
851 
852  if (num_bin_types == 1) { AMREX_ASSERT(ref_ratio[0] == 1); }
853 
854  BL_PROFILE("NeighborParticleContainer::buildNeighborList");
855 
856  resizeContainers(this->numLevels());
857 
858  for (int lev = 0; lev < this->numLevels(); ++lev)
859  {
860  m_neighbor_list[lev].clear();
861 
862  for (MyParIter pti(*this, lev); pti.isValid(); ++pti) {
863  PairIndex index(pti.index(), pti.LocalTileIndex());
864  m_neighbor_list[lev][index];
865  }
866 
867 #ifndef AMREX_USE_GPU
868  neighbor_list[lev].clear();
869  for (MyParIter pti(*this, lev); pti.isValid(); ++pti) {
870  PairIndex index(pti.index(), pti.LocalTileIndex());
871  neighbor_list[lev][index];
872  }
873 #endif
874 
875  auto& plev = this->GetParticles(lev);
876  const auto& geom = this->Geom(lev);
877 
878 #ifdef AMREX_USE_OMP
879 #pragma omp parallel if (Gpu::notInLaunchRegion())
880 #endif
881  for (MyParIter pti(*this, lev); pti.isValid(); ++pti)
882  {
883  int gid = pti.index();
884  int tid = pti.LocalTileIndex();
885  auto index = std::make_pair(gid, tid);
886  auto& ptile = plev[index];
887 
888  if (ptile.numParticles() == 0) { continue; }
889 
890  Box bx = pti.tilebox();
891  int ng = 1;
892 
893  auto& soa = pti.GetStructOfArrays();
894  auto TypeVec = soa.GetIntData(type_ind);
895  int* bin_type_array = TypeVec.data();
896 
897  Gpu::DeviceVector<int> off_bins_v(num_bin_types+1,0);
898  Gpu::DeviceVector<int> nbins_v(num_bin_types+1,0);
899  Gpu::DeviceVector<Dim3> lo_v(num_bin_types);
900  Gpu::DeviceVector<Dim3> hi_v(num_bin_types);
903 
904  for (int type(0); type<num_bin_types; ++type) {
905  // Domain, RB, Coord, Per
906  Box dom = geom.Domain();
907  const Real* plo = geom.ProbLo();
908  const Real* phi = geom.ProbHi();
909  auto lcoord = geom.Coord();
910  Array<int,AMREX_SPACEDIM> lper = geom.isPeriodic();
911 
912  // Refined tile box and domain
913  Box lbx( bx.smallEnd(), bx.bigEnd(), bx.ixType() );
914  Box ldom( dom.smallEnd(),dom.bigEnd(),dom.ixType() );
915  lbx.refine( ref_ratio[type] );
916  ldom.refine( ref_ratio[type] );
917 
918  // Local copy of RB for refined geom
919  RealBox lrb(plo,phi);
920 
921  // New geometry with refined domain
922  Geometry lgeom(ldom,lrb,lcoord,lper);
923 
924  // Grow for ghost cells
925  int NGhost = ref_ratio[type]*m_num_neighbor_cells;
926  lbx.grow(NGhost);
927 
928  // Store for memcpy
929  auto nbins = int(lbx.numPts());
930  Dim3 lo = lbound( lbx );
931  Dim3 hi = ubound( lbx );
932 
933  auto dxInv = lgeom.InvCellSizeArray();
934  auto ploa = lgeom.ProbLoArray();
935 
936 #ifdef AMREX_USE_GPU
937  Gpu::htod_memcpy_async( dxi_v.data() + type, dxInv.data(), sizeof(dxInv) );
938  Gpu::htod_memcpy_async( plo_v.data() + type, ploa.data() , sizeof(ploa) );
939  Gpu::htod_memcpy_async( lo_v.data() + type, &lo , sizeof(lo) );
940  Gpu::htod_memcpy_async( hi_v.data() + type, &hi , sizeof(hi) );
941  Gpu::htod_memcpy_async( nbins_v.data() + type, &nbins , sizeof(nbins) );
942 #else
943  std::memcpy( dxi_v.data() + type, dxInv.data(), sizeof(dxInv) );
944  std::memcpy( plo_v.data() + type, ploa.data() , sizeof(ploa) );
945  std::memcpy( lo_v.data() + type, &lo , sizeof(lo) );
946  std::memcpy( hi_v.data() + type, &hi , sizeof(hi) );
947  std::memcpy( nbins_v.data() + type, &nbins , sizeof(nbins) );
948 #endif
949  }
950 
951  Gpu::exclusive_scan(nbins_v.begin(), nbins_v.end(), off_bins_v.begin());
952 
953  m_neighbor_list[lev][index].build(ptile,
954  check_pair,
955  off_bins_v, dxi_v, plo_v, lo_v, hi_v,
956  ng, num_bin_types, bin_type_array);
957 
958 #ifndef AMREX_USE_GPU
959  BL_PROFILE_VAR("CPU_CopyNeighborList()",CPUCNL);
960 
961  const auto& counts = m_neighbor_list[lev][index].GetCounts();
962  const auto& list = m_neighbor_list[lev][index].GetList();
963 
964  int li = 0;
965  for (int i = 0; i < ptile.numParticles(); ++i) {
966  auto cnt = counts[i];
967  neighbor_list[lev][index].push_back(cnt);
968  for (size_t j = 0; j < cnt; ++j) {
969  neighbor_list[lev][index].push_back(list[li++]+1);
970  }
971  }
972 
973  BL_PROFILE_VAR_STOP(CPUCNL);
974 #endif
975  } //ParIter
976  } //Lev
977 }
978 
979 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
980 template <class CheckPair>
981 void
983 selectActualNeighbors (CheckPair const& check_pair, int num_cells)
984 {
985  BL_PROFILE("NeighborParticleContainer::selectActualNeighbors");
986  const auto& geom_fine = this->Geom(0);
987  const auto& ba_fine = this->ParticleBoxArray(0);
988  if (ba_fine.size() == 1 && !geom_fine.isAnyPeriodic()) {
989  return;
990  }
991 
992  for (int lev = 0; lev < this->numLevels(); ++lev)
993  {
994  // clear previous neighbor particle ids
995  if (!m_boundary_particle_ids.empty()) {
996  for (auto& keyval: m_boundary_particle_ids[lev]) {
997  keyval.second.clear();
998  }
999  }
1000 
1001  for (MyParIter pti(*this, lev); pti.isValid(); ++pti) {
1002  PairIndex index(pti.index(), pti.LocalTileIndex());
1003 
1004  // id of actual particles that need to be sent
1005  m_boundary_particle_ids[lev][index];
1006  m_boundary_particle_ids[lev][index].resize(pti.numNeighborParticles());
1007  auto* p_boundary_particle_ids = m_boundary_particle_ids[lev][index].dataPtr();
1008 
1009  const auto& aos = pti.GetArrayOfStructs();
1010  const auto* pstruct = aos().dataPtr();
1011  const auto ptile_data = this->ParticlesAt(lev, pti).getConstParticleTileData();
1012 
1013  Box box = pti.validbox();
1014  Box grownBox = pti.tilebox();
1015  grownBox.grow(computeRefFac(0, lev).max()*m_num_neighbor_cells);
1016  const auto lo = lbound(grownBox);
1017  const auto hi = ubound(grownBox);
1018 
1019  const auto& geom = this->Geom(lev);
1020  const auto domain = geom.Domain();
1021  const auto dxi = geom.InvCellSizeArray();
1022  const auto plo = geom.ProbLoArray();
1023 
1024  const size_t np_real = pti.numRealParticles();
1025  const size_t np_total = aos().size();
1026 
1028  bins.build(np_total, pstruct, grownBox,
1029  [=] AMREX_GPU_DEVICE (const ParticleType& p) noexcept -> IntVect
1030  {
1031  AMREX_D_TERM(int i = static_cast<int>(amrex::Math::floor((p.pos(0)-plo[0])*dxi[0]) - lo.x);,
1032  int j = static_cast<int>(amrex::Math::floor((p.pos(1)-plo[1])*dxi[1]) - lo.y);,
1033  int k = static_cast<int>(amrex::Math::floor((p.pos(2)-plo[2])*dxi[2]) - lo.z));
1034  AMREX_D_TERM(AMREX_ASSERT(i >= 0);, AMREX_ASSERT(j >= 0);, AMREX_ASSERT(k >= 0));
1035 
1036  return IntVect(AMREX_D_DECL(i, j, k));
1037  });
1038 
1039  auto pperm = bins.permutationPtr();
1040  auto poffset = bins.offsetsPtr();
1041 
1042  Gpu::Buffer<unsigned int> np_boundary({0});
1043  unsigned int* p_np_boundary = np_boundary.data();
1044 
1045  AMREX_FOR_1D ( np_real, i,
1046  {
1047  IntVect iv(AMREX_D_DECL(
1048  static_cast<int>(amrex::Math::floor((pstruct[i].pos(0)-plo[0])*dxi[0])) - lo.x,
1049  static_cast<int>(amrex::Math::floor((pstruct[i].pos(1)-plo[1])*dxi[1])) - lo.y,
1050  static_cast<int>(amrex::Math::floor((pstruct[i].pos(2)-plo[2])*dxi[2])) - lo.z));
1051  auto iv3 = iv.dim3();
1052 
1053  int ix = iv3.x;
1054  int iy = iv3.y;
1055  int iz = iv3.z;
1056 
1057  int nx = hi.x-lo.x+1;
1058  int ny = hi.y-lo.y+1;
1059  int nz = hi.z-lo.z+1;
1060 
1061  bool isActualNeighbor = false;
1062  for (int ii = amrex::max(ix-num_cells, 0); ii <= amrex::min(ix+num_cells, nx); ++ii) {
1063  for (int jj = amrex::max(iy-num_cells, 0); jj <= amrex::min(iy+num_cells, ny); ++jj) {
1064  for (int kk = amrex::max(iz-num_cells, 0); kk <= amrex::min(iz+num_cells, nz); ++kk) {
1065  if (isActualNeighbor) { break; }
1066  int nbr_cell_id = (ii * ny + jj) * nz + kk;
1067  for (auto p = poffset[nbr_cell_id]; p < poffset[nbr_cell_id+1]; ++p) {
1068  if (pperm[p] == int(i)) { continue; }
1069  if (detail::call_check_pair(check_pair, ptile_data, ptile_data, i, pperm[p])) {
1070  IntVect cell_ijk = getParticleCell(pstruct[pperm[p]], plo, dxi, domain);
1071  if (!box.contains(cell_ijk)) {
1072  unsigned int loc = Gpu::Atomic::Add(p_np_boundary, 1U);
1073  p_boundary_particle_ids[loc] = i;
1074  isActualNeighbor = true;
1075  break;
1076  }
1077  }// end if check_pair
1078  }
1079  }
1080  }
1081  }
1082  });// end amrex_for_1d
1083 
1084  unsigned int* p_np_boundary_h = np_boundary.copyToHost();
1085  m_boundary_particle_ids[lev][index].resize(*p_np_boundary_h);
1086 
1087  }// end mypariter
1088  }// end lev
1089 }
1090 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
1091 void
1094 {
1095  BL_PROFILE("NeighborParticleContainer::printNeighborList");
1096 
1097  for (int lev = 0; lev < this->numLevels(); ++lev)
1098  {
1099  for(MFIter mfi = this->MakeMFIter(lev); mfi.isValid(); ++mfi)
1100  {
1101  int gid = mfi.index();
1102  int tid = mfi.LocalTileIndex();
1103  auto index = std::make_pair(gid, tid);
1104  m_neighbor_list[lev][index].print();
1105  }
1106  }
1107 }
1108 
1109 template <int NStructReal, int NStructInt, int NArrayReal, int NArrayInt>
1110 void
1112 resizeContainers (int num_levels)
1113 {
1114  this->reserveData();
1115  this->resizeData();
1116  if ( static_cast<int>(neighbors.size()) <= num_levels )
1117  {
1118  neighbors.resize(num_levels);
1119  m_neighbor_list.resize(num_levels);
1120  neighbor_list.resize(num_levels);
1121  mask_ptr.resize(num_levels);
1122  buffer_tag_cache.resize(num_levels);
1123  local_neighbor_sizes.resize(num_levels);
1124  if ( enableInverse() ) { inverse_tags.resize(num_levels); }
1125  }
1126 
1127  AMREX_ASSERT((neighbors.size() == m_neighbor_list.size()) &&
1128  (neighbors.size() == mask_ptr.size() ) );
1129 }
#define BL_PROFILE(a)
Definition: AMReX_BLProfiler.H:551
#define BL_PROFILE_VAR_STOP(vname)
Definition: AMReX_BLProfiler.H:563
#define BL_PROFILE_VAR(fname, vname)
Definition: AMReX_BLProfiler.H:560
#define AMREX_ASSERT(EX)
Definition: AMReX_BLassert.H:38
#define AMREX_FOR_1D(...)
Definition: AMReX_GpuLaunch.nolint.H:41
#define AMREX_GPU_DEVICE
Definition: AMReX_GpuQualifiers.H:18
Array4< int const > mask
Definition: AMReX_InterpFaceRegister.cpp:93
#define AMREX_D_TERM(a, b, c)
Definition: AMReX_SPACE.H:129
#define AMREX_D_DECL(a, b, c)
Definition: AMReX_SPACE.H:104
A collection of Boxes stored in an Array.
Definition: AMReX_BoxArray.H:550
std::vector< std::pair< int, Box > > intersections(const Box &bx) const
Return intersections of Box and BoxArray.
static bool SameRefs(const BoxArray &lhs, const BoxArray &rhs)
whether two BoxArrays share the same data
Definition: AMReX_BoxArray.H:820
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 BoxND & refine(int ref_ratio) noexcept
Refine BoxND by given (positive) refinement ratio. NOTE: if type(dir) = CELL centered: lo <- lo*ratio...
Definition: AMReX_Box.H:684
AMREX_GPU_HOST_DEVICE void next(IntVectND< dim > &) const noexcept
Step through the rectangle. It is a runtime error to give a point not inside rectangle....
Definition: AMReX_Box.H:1059
AMREX_GPU_HOST_DEVICE BoxND & grow(int i) noexcept
Definition: AMReX_Box.H:627
AMREX_GPU_HOST_DEVICE IndexTypeND< dim > ixType() const noexcept
Returns the indexing type.
Definition: AMReX_Box.H:127
AMREX_GPU_HOST_DEVICE Long index(const IntVectND< dim > &v) const noexcept
Returns offset of point from smallend; i.e. index(smallend) -> 0, bigend would return numPts()-1....
Definition: AMReX_Box.H:993
AMREX_GPU_HOST_DEVICE const IntVectND< dim > & bigEnd() const &noexcept
Get the bigend.
Definition: AMReX_Box.H:116
AMREX_GPU_HOST_DEVICE BoxND & coarsen(int ref_ratio) noexcept
Coarsen BoxND by given (positive) refinement ratio. NOTE: if type(dir) = CELL centered: lo <- lo/rati...
Definition: AMReX_Box.H:708
AMREX_GPU_HOST_DEVICE Long numPts() const noexcept
Returns the number of points contained in the BoxND.
Definition: AMReX_Box.H:346
AMREX_GPU_HOST_DEVICE bool contains(const IntVectND< dim > &p) const noexcept
Returns true if argument is contained within BoxND.
Definition: AMReX_Box.H:204
GpuArray< Real, AMREX_SPACEDIM > InvCellSizeArray() const noexcept
Definition: AMReX_CoordSys.H:87
index_type * offsetsPtr() noexcept
returns the pointer to the offsets array
Definition: AMReX_DenseBins.H:510
void build(N nitems, const_pointer_input_type v, const Box &bx, F &&f)
Populate the bins with a set of items.
Definition: AMReX_DenseBins.H:130
index_type * permutationPtr() noexcept
returns the pointer to the permutation array
Definition: AMReX_DenseBins.H:507
Calculates the distribution of FABs to MPI processes.
Definition: AMReX_DistributionMapping.H:41
static bool SameRefs(const DistributionMapping &lhs, const DistributionMapping &rhs)
Definition: AMReX_DistributionMapping.H:187
Rectangular problem domain geometry.
Definition: AMReX_Geometry.H:73
GpuArray< Real, AMREX_SPACEDIM > ProbLoArray() const noexcept
Definition: AMReX_Geometry.H:186
Periodicity periodicity() const noexcept
Definition: AMReX_Geometry.H:355
Definition: AMReX_GpuBuffer.H:17
T const * data() const noexcept
Definition: AMReX_GpuBuffer.H:65
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 dim3() const noexcept
Definition: AMReX_IntVect.H:163
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE int max() const noexcept
maximum (no absolute values) value
Definition: AMReX_IntVect.H:214
AMREX_GPU_HOST_DEVICE static constexpr AMREX_FORCE_INLINE IntVectND< dim > TheZeroVector() noexcept
This static member function returns a reference to a constant IntVectND object, all of whose dim argu...
Definition: AMReX_IntVect.H:672
Definition: AMReX_MFIter.H:57
Box tilebox() const noexcept
Return the tile Box at the current index.
Definition: AMReX_MFIter.cpp:385
bool isValid() const noexcept
Is the iterator valid i.e. is it associated with a FAB?
Definition: AMReX_MFIter.H:141
Definition: AMReX_NeighborList.H:247
void GetNeighborCommTags()
Definition: AMReX_NeighborParticlesI.H:190
std::pair< int, int > PairIndex
Definition: AMReX_NeighborParticles.H:196
void selectActualNeighbors(CheckPair const &check_pair, int num_cells=1)
bool areMasksValid()
Definition: AMReX_NeighborParticlesI.H:127
void clearNeighbors()
Definition: AMReX_NeighborParticlesI.H:679
void fillNeighbors()
Definition: AMReX_NeighborParticlesI.H:637
void cacheNeighborInfo()
Definition: AMReX_NeighborParticlesI.H:312
void resizeContainers(int num_levels)
Definition: AMReX_NeighborParticlesI.H:1112
void setIntCommComp(int i, bool value)
Definition: AMReX_NeighborParticlesI.H:66
void updateNeighbors(bool boundary_neighbors_only=false)
Definition: AMReX_NeighborParticlesI.H:662
void buildNeighborList(CheckPair const &check_pair, bool sort=false)
void sumNeighbors(int real_start_comp, int real_num_comp, int int_start_comp, int int_num_comp)
Definition: AMReX_NeighborParticlesI.H:649
void getNeighborTags(Vector< NeighborCopyTag > &tags, const ParticleType &p, int nGrow, const NeighborCopyTag &src_tag, const MyParIter &pti)
Definition: AMReX_NeighborParticlesI.H:535
void setRealCommComp(int i, bool value)
Definition: AMReX_NeighborParticlesI.H:58
void initializeCommComps()
Definition: AMReX_NeighborParticlesI.H:45
void calcCommSize()
Definition: AMReX_NeighborParticlesI.H:74
void printNeighborList()
Definition: AMReX_NeighborParticlesI.H:1093
typename ParticleContainerType::ParticleType ParticleType
Definition: AMReX_NeighborParticles.H:38
void GetCommTagsBox(Vector< NeighborCommTag > &tags, int lev, const Box &in_box)
Definition: AMReX_NeighborParticlesI.H:265
void Regrid(const DistributionMapping &dmap, const BoxArray &ba)
Definition: AMReX_NeighborParticlesI.H:93
NeighborParticleContainer(ParGDBBase *gdb, int ncells)
Definition: AMReX_NeighborParticlesI.H:10
IntVect computeRefFac(int src_lev, int lev)
Definition: AMReX_NeighborParticlesI.H:246
void BuildMasks()
Definition: AMReX_NeighborParticlesI.H:151
Definition: AMReX_PODVector.H:246
T * data() noexcept
Definition: AMReX_PODVector.H:593
iterator begin() noexcept
Definition: AMReX_PODVector.H:601
iterator end() noexcept
Definition: AMReX_PODVector.H:605
void push_back(const T &a_value)
Definition: AMReX_PODVector.H:556
Definition: AMReX_ParGDB.H:13
Definition: AMReX_ParIter.H:113
A distributed container for Particles sorted onto the levels, grids, and tiles of a block-structured ...
Definition: AMReX_ParticleContainer.H:145
This provides length of period for periodic domains. 0 means it is not periodic in that direction....
Definition: AMReX_Periodicity.H:17
std::vector< IntVect > shiftIntVect(IntVect const &nghost=IntVect(0)) const
Definition: AMReX_Periodicity.cpp:8
bool isAnyPeriodic() const noexcept
Definition: AMReX_Periodicity.H:22
bool isPeriodic(int dir) const noexcept
Definition: AMReX_Periodicity.H:26
A Box with real dimensions. A RealBox is OK iff volume >= 0.
Definition: AMReX_RealBox.H:21
This class is a thin wrapper around std::vector. Unlike vector, Vector::operator[] provides bound che...
Definition: AMReX_Vector.H:27
Long size() const noexcept
Definition: AMReX_Vector.H:50
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE T Add(T *sum, T value) noexcept
Definition: AMReX_GpuAtomic.H:198
OutIter exclusive_scan(InIter begin, InIter end, OutIter result)
Definition: AMReX_Scan.H:1377
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void * memcpy(void *dest, const void *src, std::size_t count)
Definition: AMReX_GpuUtility.H:214
void htod_memcpy_async(void *p_d, const void *p_h, const std::size_t sz) noexcept
Definition: AMReX_GpuDevice.H:251
int MyProc()
Definition: AMReX_MPMD.cpp:117
constexpr int get_thread_num()
Definition: AMReX_OpenMP.H:37
constexpr int get_max_threads()
Definition: AMReX_OpenMP.H:36
int MyProcSub() noexcept
my sub-rank in current frame
Definition: AMReX_ParallelContext.H:76
int global_to_local_rank(int rank) noexcept
Definition: AMReX_ParallelContext.H:98
AMREX_GPU_HOST_DEVICE auto call_check_pair(F const &check_pair, const SrcData &src_tile, const DstData &dst_tile, N1 i, N2 j) noexcept -> decltype(check_pair(src_tile.m_aos[i], dst_tile.m_aos[j]))
Definition: AMReX_NeighborList.H:19
constexpr int iz
Definition: AMReX_Interp_3D_C.H:37
constexpr int iy
Definition: AMReX_Interp_2D_C.H:33
constexpr int ix
Definition: AMReX_Interp_2D_C.H:32
DistributionMapping const & DistributionMap(FabArrayBase const &fa)
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE BoxND< dim > shift(const BoxND< dim > &b, int dir, int nzones) noexcept
Return a BoxND with indices shifted by nzones in dir direction.
Definition: AMReX_Box.H:1372
BoxND< AMREX_SPACEDIM > Box
Definition: AMReX_BaseFwd.H:27
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE IntVect getParticleCell(P const &p, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &plo, amrex::GpuArray< amrex::Real, AMREX_SPACEDIM > const &dxi) noexcept
Returns the cell index for a given particle using the provided lower bounds and cell sizes.
Definition: AMReX_ParticleUtil.H:374
void EnsureThreadSafeTiles(PC &pc)
Definition: AMReX_ParticleUtil.H:705
AMREX_GPU_HOST_DEVICE constexpr AMREX_FORCE_INLINE const T & max(const T &a, const T &b) noexcept
Definition: AMReX_Algorithm.H:35
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 constexpr AMREX_FORCE_INLINE const T & min(const T &a, const T &b) noexcept
Definition: AMReX_Algorithm.H:21
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 end(BoxND< dim > const &box) noexcept
Definition: AMReX_Box.H:1890
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
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 begin(BoxND< dim > const &box) noexcept
Definition: AMReX_Box.H:1881
bool SameIteratorsOK(const PC1 &pc1, const PC2 &pc2)
Definition: AMReX_ParticleUtil.H:693
int numParticlesOutOfRange(Iterator const &pti, int nGrow)
Returns the number of particles that are more than nGrow cells from the box correspond to the input i...
Definition: AMReX_ParticleUtil.H:34
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition: AMReX.cpp:225
const int[]
Definition: AMReX_BLProfiler.cpp:1664
IntVect computeRefFac(const ParGDBBase *a_gdb, int src_lev, int lev)
Definition: AMReX_ParticleUtil.cpp:6
void RemoveDuplicates(Vector< T > &vec)
Definition: AMReX_Vector.H:190
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE int getTileIndex(const IntVect &iv, const Box &box, const bool a_do_tiling, const IntVect &a_tile_size, Box &tbx)
Definition: AMReX_ParticleUtil.H:222
std::array< T, N > Array
Definition: AMReX_Array.H:24
Definition: AMReX_Dim3.H:12
Definition: AMReX_NeighborParticles.H:152
int tile_id
Definition: AMReX_NeighborParticles.H:161
int grid_id
Definition: AMReX_NeighborParticles.H:160
Definition: AMReX_NeighborParticles.H:82
int src_index
Definition: AMReX_NeighborParticles.H:86
int tile
Definition: AMReX_NeighborParticles.H:85
int grid
Definition: AMReX_NeighborParticles.H:84
int level
Definition: AMReX_NeighborParticles.H:83
IntVect periodic_shift
Definition: AMReX_NeighborParticles.H:88
Definition: AMReX_NeighborParticles.H:51
int thread_num
Definition: AMReX_NeighborParticles.H:60
int src_level
Definition: AMReX_NeighborParticles.H:56
int src_index
Definition: AMReX_NeighborParticles.H:59
int src_grid
Definition: AMReX_NeighborParticles.H:57
int dst_level
Definition: AMReX_NeighborParticles.H:52
int src_tile
Definition: AMReX_NeighborParticles.H:58