Block-Structured AMR Software Framework
AMReX_WriteBinaryParticleData.H
Go to the documentation of this file.
1 #ifndef AMREX_WRITE_BINARY_PARTICLE_DATA_H
2 #define AMREX_WRITE_BINARY_PARTICLE_DATA_H
3 #include <AMReX_Config.H>
4 
5 #include <AMReX_TypeTraits.H>
6 #include <AMReX_ParticleUtil.H>
7 #include <AMReX_GpuDevice.H>
8 
10 {
11  template <typename SrcData>
13  int operator() (const SrcData& src, int i) const noexcept
14  {
15  return (src.id(i) > 0);
16  }
17 };
18 
19 namespace particle_detail {
20 
21 template <typename ParticleReal>
22 std::size_t PSizeInFile (const Vector<int>& wrc, const Vector<int>& wic)
23 {
24  std::size_t rsize = sizeof(ParticleReal)*std::accumulate(wrc.begin(), wrc.end(), 0);
25  std::size_t isize = sizeof(int)*std::accumulate(wic.begin(), wic.end(), 0);
26  return rsize + isize + AMREX_SPACEDIM*sizeof(ParticleReal) + 2*sizeof(int);
27 }
28 
29 template <template <class, class> class Container,
30  class Allocator,
31  class PTile,
32  class F>
33 std::enable_if_t<RunOnGpu<typename Container<int, Allocator>::allocator_type>::value>
34 fillFlags (Container<int, Allocator>& pflags, const PTile& ptile, F const& f)
35 {
36  const auto& ptd = ptile.getConstParticleTileData();
37  const auto np = ptile.numParticles();
38  pflags.resize(np, 0);
39  auto flag_ptr = pflags.data();
41  [=] AMREX_GPU_DEVICE (int k, amrex::RandomEngine const& engine) noexcept
42  {
43  const auto p = ptd.getSuperParticle(k);
44  amrex::ignore_unused(flag_ptr, f, engine);
45  if constexpr (IsCallable<F,decltype(p),RandomEngine>::value) {
46  flag_ptr[k] = f(p,engine);
47  } else {
48  flag_ptr[k] = f(p);
49  }
50  });
51 }
52 
53 template <template <class, class> class Container,
54  class Allocator,
55  class PTile,
56  class F>
57 std::enable_if_t<!RunOnGpu<typename Container<int, Allocator>::allocator_type>::value>
58 fillFlags (Container<int, Allocator>& pflags, const PTile& ptile, F const& f)
59 {
60  const auto& ptd = ptile.getConstParticleTileData();
61  const auto np = ptile.numParticles();
62  pflags.resize(np, 0);
63  auto flag_ptr = pflags.data();
64  for (int k = 0; k < np; ++k) {
65  const auto p = ptd.getSuperParticle(k);
66  if constexpr (IsCallable<F,decltype(p),RandomEngine>::value) {
67  flag_ptr[k] = f(p,RandomEngine{});
68  } else {
69  flag_ptr[k] = f(p);
70  }
71  }
72 }
73 
74 template <template <class, class> class Container, class Allocator, class PC>
75 std::enable_if_t<RunOnGpu<typename Container<int, Allocator>::allocator_type>::value, amrex::Long>
76 countFlags (const Vector<std::map<std::pair<int,int>,Container<int,Allocator>>>& particle_io_flags, const PC& pc)
77 {
78  ReduceOps<ReduceOpSum> reduce_op;
79  ReduceData<Long> reduce_data(reduce_op);
80  using ReduceTuple = typename decltype(reduce_data)::Type;
81 
82  for (int lev = 0; lev < pc.GetParticles().size(); lev++)
83  {
84  const auto& pmap = pc.GetParticles(lev);
85  for (const auto& kv : pmap)
86  {
87  const auto& pflags = particle_io_flags[lev].at(kv.first);
88  const auto flag_ptr = pflags.data();
89  reduce_op.eval(pflags.size(), reduce_data,
90  [=] AMREX_GPU_DEVICE (const int i) -> ReduceTuple
91  {
92  return flag_ptr[i] ? 1 : 0;
93  });
94  }
95  }
96  ReduceTuple hv = reduce_data.value(reduce_op);
97  return amrex::get<0>(hv);
98 }
99 
100 template <template <class, class> class Container, class Allocator>
101 std::enable_if_t<RunOnGpu<typename Container<int, Allocator>::allocator_type>::value, int>
102 countFlags (const Container<int,Allocator>& pflags)
103 {
104  ReduceOps<ReduceOpSum> reduce_op;
105  ReduceData<Long> reduce_data(reduce_op);
106  using ReduceTuple = typename decltype(reduce_data)::Type;
107 
108  const auto flag_ptr = pflags.data();
109  reduce_op.eval(pflags.size(), reduce_data,
110  [=] AMREX_GPU_DEVICE (const int i) -> ReduceTuple
111  {
112  return flag_ptr[i] ? 1 : 0;
113  });
114  ReduceTuple hv = reduce_data.value(reduce_op);
115  return amrex::get<0>(hv);
116 }
117 
118 template <template <class, class> class Container, class Allocator, class PC>
119 std::enable_if_t<!RunOnGpu<typename Container<int, Allocator>::allocator_type>::value, amrex::Long>
120 countFlags (const Vector<std::map<std::pair<int,int>,Container<int,Allocator>>>& particle_io_flags, const PC& pc)
121 {
122  amrex::Long nparticles = 0;
123  for (int lev = 0; lev < pc.GetParticles().size(); lev++)
124  {
125  const auto& pmap = pc.GetParticles(lev);
126  for (const auto& kv : pmap)
127  {
128  const auto& pflags = particle_io_flags[lev].at(kv.first);
129  for (int k = 0; k < kv.second.numParticles(); ++k)
130  {
131  if (pflags[k]) { nparticles++; }
132  }
133  }
134  }
135  return nparticles;
136 }
137 
138 template <template <class, class> class Container, class Allocator>
139 std::enable_if_t<!RunOnGpu<typename Container<int, Allocator>::allocator_type>::value, int>
140 countFlags (const Container<int,Allocator>& pflags)
141 {
142  int nparticles = 0;
143  for (std::size_t k = 0; k < pflags.size(); ++k)
144  {
145  if (pflags[k]) { nparticles++; }
146  }
147  return nparticles;
148 }
149 
150 template <typename P, typename I>
152 void packParticleIDs (I* idata, const P& p, bool is_checkpoint) noexcept
153 {
154  if (is_checkpoint) {
155  std::int32_t xi, yi;
156  std::uint32_t xu, yu;
157  xu = (std::uint32_t)((p.m_idcpu & 0xFFFFFFFF00000000LL) >> 32);
158  yu = (std::uint32_t)( p.m_idcpu & 0xFFFFFFFFLL);
159  amrex::Gpu::memcpy(&xi, &xu, sizeof(xu));
160  amrex::Gpu::memcpy(&yi, &yu, sizeof(yu));
161  idata[0] = xi;
162  idata[1] = yi;
163  } else {
164  idata[0] = p.id();
165  idata[1] = p.cpu();
166  }
167 }
168 
169 template <class PC>
170 std::enable_if_t<RunOnGpu<typename PC::template AllocatorType<int>>::value>
171 packIOData (Vector<int>& idata, Vector<ParticleReal>& rdata, const PC& pc, int lev, int grid,
172  const Vector<int>& write_real_comp, const Vector<int>& write_int_comp,
173  const Vector<std::map<std::pair<int, int>, typename PC::IntVector>>& particle_io_flags,
174  const Vector<int>& tiles, int np, bool is_checkpoint)
175 {
176  int num_output_int = 0;
177  for (int i = 0; i < pc.NumIntComps() + PC::NStructInt; ++i) {
178  if (write_int_comp[i]) { ++num_output_int; }
179  }
180 
181  const Long iChunkSize = 2 + num_output_int;
182  idata.resize(np*iChunkSize);
183 
184  int num_output_real = 0;
185  for (int i : write_real_comp) {
186  if (i) { ++num_output_real; }
187  }
188 
189  const Long rChunkSize = AMREX_SPACEDIM + num_output_real;
190  rdata.resize(np*rChunkSize);
191 
192  typename PC::IntVector write_int_comp_d(write_int_comp.size());
193  typename PC::IntVector write_real_comp_d(write_real_comp.size());
194  Gpu::copyAsync(Gpu::hostToDevice, write_int_comp.begin(), write_int_comp.end(),
195  write_int_comp_d.begin());
196  Gpu::copyAsync(Gpu::hostToDevice, write_real_comp.begin(), write_real_comp.end(),
197  write_real_comp_d.begin());
198 
199  const auto write_int_comp_d_ptr = write_int_comp_d.data();
200  const auto write_real_comp_d_ptr = write_real_comp_d.data();
201 
202  std::size_t poffset = 0;
203  for (int tile : tiles) {
204  const auto& ptile = pc.ParticlesAt(lev, grid, tile);
205  const auto& pflags = particle_io_flags[lev].at(std::make_pair(grid, tile));
206  int np_tile = ptile.numParticles();
207  typename PC::IntVector offsets(np_tile);
208  int num_copies = Scan::ExclusiveSum(np_tile, pflags.begin(), offsets.begin(), Scan::retSum);
209 
210  typename PC::IntVector idata_d(num_copies*iChunkSize);
211  typename PC::RealVector rdata_d(num_copies*rChunkSize);
212 
213  const auto flag_ptr = pflags.data();
214 
215  auto idata_d_ptr = idata_d.data();
216  auto rdata_d_ptr = rdata_d.data();
217 
218  const auto& ptd = ptile.getConstParticleTileData();
219  amrex::ParallelFor(num_copies,
220  [=] AMREX_GPU_DEVICE (int pindex) noexcept
221  {
222  // might be worth using shared memory here
223  const auto p = ptd.getSuperParticle(pindex);
224 
225  if (flag_ptr[pindex]) {
226  std::size_t iout_index = pindex*iChunkSize;
227  packParticleIDs(&idata_d_ptr[iout_index], p, is_checkpoint);
228  iout_index += 2;
229 
230  std::size_t rout_index = pindex*rChunkSize;
231  for (int j = 0; j < AMREX_SPACEDIM; j++) {
232  rdata_d_ptr[rout_index] = p.pos(j);
233  rout_index++;
234  }
235 
236  for (int j = 0; j < PC::SuperParticleType::NInt; j++) {
237  if (write_int_comp_d_ptr[j]) {
238  idata_d_ptr[iout_index] = p.idata(j);
239  iout_index++;
240  }
241  }
242 
243  for (int j = 0; j < ptd.m_num_runtime_int; j++) {
244  if (write_int_comp_d_ptr[PC::SuperParticleType::NInt + j]) {
245  idata_d_ptr[iout_index] = ptd.m_runtime_idata[j][pindex];
246  iout_index++;
247  }
248  }
249 
250  // extra SoA Real components
251  const int real_start_offset = PC::ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0; // pure SoA: skip positions
252  for (int j = real_start_offset; j < PC::SuperParticleType::NReal; j++) {
253  const int write_comp_index = j-real_start_offset;
254  if (write_real_comp_d_ptr[write_comp_index]) {
255  rdata_d_ptr[rout_index] = p.rdata(j);
256  rout_index++;
257  }
258  }
259 
260  for (int j = 0; j < ptd.m_num_runtime_real; j++) {
261  if (write_real_comp_d_ptr[PC::SuperParticleType::NReal+j-real_start_offset]) {
262  rdata_d_ptr[rout_index] = ptd.m_runtime_rdata[j][pindex];
263  rout_index++;
264  }
265  }
266  }
267  });
268 
269  Gpu::copyAsync(Gpu::deviceToHost, idata_d.begin(), idata_d.end(),
270  idata.begin() + typename PC::IntVector::difference_type(poffset));
271  Gpu::copyAsync(Gpu::deviceToHost, rdata_d.begin(), rdata_d.end(),
272  rdata.begin() + typename PC::RealVector::difference_type(poffset));
274 
275  poffset += num_copies;
276  }
277 }
278 
279 template <class PC>
280 std::enable_if_t<!RunOnGpu<typename PC::template AllocatorType<int>>::value>
281 packIOData (Vector<int>& idata, Vector<ParticleReal>& rdata, const PC& pc, int lev, int grid,
282  const Vector<int>& write_real_comp, const Vector<int>& write_int_comp,
283  const Vector<std::map<std::pair<int, int>, typename PC::IntVector>>& particle_io_flags,
284  const Vector<int>& tiles, int np, bool is_checkpoint)
285 {
286  int num_output_int = 0;
287  for (int i = 0; i < pc.NumIntComps() + PC::NStructInt; ++i) {
288  if (write_int_comp[i]) { ++num_output_int; }
289  }
290 
291  const Long iChunkSize = 2 + num_output_int;
292  idata.resize(np*iChunkSize);
293 
294  int num_output_real = 0;
295  for (int i : write_real_comp) {
296  if (i) { ++num_output_real; }
297  }
298 
299  const Long rChunkSize = AMREX_SPACEDIM + num_output_real;
300  rdata.resize(np*rChunkSize);
301 
302  int* iptr = idata.dataPtr();
303  ParticleReal* rptr = rdata.dataPtr();
304  for (int tile : tiles) {
305  const auto& ptile = pc.ParticlesAt(lev, grid, tile);
306  const auto& pflags = particle_io_flags[lev].at(std::make_pair(grid, tile));
307  for (int pindex = 0; pindex < ptile.numParticles(); ++pindex) {
308  if (pflags[pindex]) {
309  const auto& soa = ptile.GetStructOfArrays();
310 
311  // note: for pure SoA particle layouts, we do write the id, cpu and positions as a struct
312  // for backwards compatibility with readers
313  if constexpr(!PC::ParticleType::is_soa_particle)
314  {
315  const auto& aos = ptile.GetArrayOfStructs();
316  const auto& p = aos[pindex];
317 
318  // Int: id, cpu
319  packParticleIDs(iptr, p, is_checkpoint);
320  iptr += 2;
321 
322  // Real: positions
323  for (int j = 0; j < AMREX_SPACEDIM; j++) { rptr[j] = p.pos(j); }
324  rptr += AMREX_SPACEDIM;
325 
326  // extra AoS Int components
327  for (int j = 0; j < PC::NStructInt; j++) {
328  if (write_int_comp[j]) {
329  *iptr = p.idata(j);
330  ++iptr;
331  }
332  }
333  // extra AoS Real components
334  for (int j = 0; j < PC::NStructReal; j++) {
335  if (write_real_comp[j]) {
336  *rptr = p.rdata(j);
337  ++rptr;
338  }
339  }
340  }
341  else {
342  uint64_t idcpu = soa.GetIdCPUData()[pindex];
343  if (is_checkpoint) {
344  std::int32_t xi, yi;
345  std::uint32_t xu, yu;
346  xu = (std::uint32_t)((idcpu & 0xFFFFFFFF00000000LL) >> 32);
347  yu = (std::uint32_t)( idcpu & 0xFFFFFFFFLL);
348  std::memcpy(&xi, &xu, sizeof(xu));
349  std::memcpy(&yi, &yu, sizeof(yu));
350  *iptr = xi;
351  iptr += 1;
352  *iptr = yi;
353  iptr += 1;
354  } else {
355  // Int: id, cpu
356  *iptr = (int) ParticleIDWrapper(idcpu);
357  iptr += 1;
358  *iptr = (int) ParticleCPUWrapper(idcpu);
359  iptr += 1;
360  }
361 
362  // Real: position
363  for (int j = 0; j < AMREX_SPACEDIM; j++) { rptr[j] = soa.GetRealData(j)[pindex]; }
364  rptr += AMREX_SPACEDIM;
365  }
366 
367  // SoA int data
368  const int int_start_offset = 0;
369  for (int j = int_start_offset; j < pc.NumIntComps(); j++) {
370  if (write_int_comp[PC::NStructInt+j]) {
371  *iptr = soa.GetIntData(j)[pindex];
372  ++iptr;
373  }
374  }
375 
376  // extra SoA Real components
377  const int real_start_offset = PC::ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0; // pure SoA: skip positions
378  for (int j = real_start_offset; j < pc.NumRealComps(); j++) {
379  const int write_comp_index = PC::NStructReal+j-real_start_offset;
380  if (write_real_comp[write_comp_index]) {
381  *rptr = (ParticleReal) soa.GetRealData(j)[pindex];
382  ++rptr;
383  }
384  }
385  }
386  }
387  }
388 }
389 }
390 
391 template <class PC, class F, std::enable_if_t<IsParticleContainer<PC>::value, int> foo = 0>
392 void WriteBinaryParticleDataSync (PC const& pc,
393  const std::string& dir, const std::string& name,
394  const Vector<int>& write_real_comp,
395  const Vector<int>& write_int_comp,
396  const Vector<std::string>& real_comp_names,
397  const Vector<std::string>& int_comp_names,
398  F const& f, bool is_checkpoint)
399 {
400  BL_PROFILE("WriteBinaryParticleData()");
401  AMREX_ASSERT(pc.OK());
402 
403  AMREX_ASSERT(sizeof(typename PC::ParticleType::RealType) == 4 ||
404  sizeof(typename PC::ParticleType::RealType) == 8);
405 
406  constexpr int NStructReal = PC::NStructReal;
407  constexpr int NStructInt = PC::NStructInt;
408 
409  const int NProcs = ParallelDescriptor::NProcs();
410  const int IOProcNumber = ParallelDescriptor::IOProcessorNumber();
411 
412  if constexpr(PC::ParticleType::is_soa_particle) {
413  AMREX_ALWAYS_ASSERT(real_comp_names.size() == pc.NumRealComps() + NStructReal - AMREX_SPACEDIM); // pure SoA: skip positions
414  } else {
415  AMREX_ALWAYS_ASSERT(real_comp_names.size() == pc.NumRealComps() + NStructReal);
416  }
417  AMREX_ALWAYS_ASSERT( int_comp_names.size() == pc.NumIntComps() + NStructInt);
418 
419  std::string pdir = dir;
420  if ( ! pdir.empty() && pdir[pdir.size()-1] != '/') { pdir += '/'; }
421  pdir += name;
422 
423  if ( ! pc.GetLevelDirectoriesCreated()) {
425  {
426  if ( ! amrex::UtilCreateDirectory(pdir, 0755))
427  {
429  }
430  }
432  }
433 
434  std::ofstream HdrFile;
435 
436  Long nparticles = 0;
437  Long maxnextid;
438 
439  // evaluate f for every particle to determine which ones to output
440  Vector<std::map<std::pair<int, int>, typename PC::IntVector > >
441  particle_io_flags(pc.GetParticles().size());
442  for (int lev = 0; lev < pc.GetParticles().size(); lev++)
443  {
444  const auto& pmap = pc.GetParticles(lev);
445  for (const auto& kv : pmap)
446  {
447  auto& flags = particle_io_flags[lev][kv.first];
448  particle_detail::fillFlags(flags, kv.second, f);
449  }
450  }
451 
453 
454  if(pc.GetUsePrePost())
455  {
456  nparticles = pc.GetNParticlesPrePost();
457  maxnextid = pc.GetMaxNextIDPrePost();
458  }
459  else
460  {
461  nparticles = particle_detail::countFlags(particle_io_flags, pc);
462  maxnextid = PC::ParticleType::NextID();
463  ParallelDescriptor::ReduceLongSum(nparticles, IOProcNumber);
464  PC::ParticleType::NextID(maxnextid);
465  ParallelDescriptor::ReduceLongMax(maxnextid, IOProcNumber);
466  }
467 
469  {
470  std::string HdrFileName = pdir;
471 
472  if ( ! HdrFileName.empty() && HdrFileName[HdrFileName.size()-1] != '/') {
473  HdrFileName += '/';
474  }
475 
476  HdrFileName += "Header";
477  pc.HdrFileNamePrePost = HdrFileName;
478 
479  HdrFile.open(HdrFileName.c_str(), std::ios::out|std::ios::trunc);
480 
481  if ( ! HdrFile.good()) { amrex::FileOpenFailed(HdrFileName); }
482 
483  //
484  // First thing written is our version string.
485  // We append "_single" or "_double" to the version string indicating
486  // whether we're using "float" or "double" floating point data.
487  //
488  std::string version_string = is_checkpoint ? PC::CheckpointVersion() : PC::PlotfileVersion();
489  if (sizeof(typename PC::ParticleType::RealType) == 4)
490  {
491  HdrFile << version_string << "_single" << '\n';
492  }
493  else
494  {
495  HdrFile << version_string << "_double" << '\n';
496  }
497 
498  int num_output_real = 0;
499  for (int i : write_real_comp) {
500  if (i) { ++num_output_real; }
501  }
502 
503  int num_output_int = 0;
504  for (int i = 0; i < pc.NumIntComps() + NStructInt; ++i) {
505  if (write_int_comp[i]) { ++num_output_int; }
506  }
507 
508  // AMREX_SPACEDIM and N for sanity checking.
509  HdrFile << AMREX_SPACEDIM << '\n';
510 
511  // The number of extra real parameters
512  HdrFile << num_output_real << '\n';
513 
514  // Real component names
515  for (int i = 0; i < (int) real_comp_names.size(); ++i ) {
516  if (write_real_comp[i]) { HdrFile << real_comp_names[i] << '\n'; }
517  }
518 
519  // The number of extra int parameters
520  HdrFile << num_output_int << '\n';
521 
522  // int component names
523  for (int i = 0; i < NStructInt + pc.NumIntComps(); ++i ) {
524  if (write_int_comp[i]) { HdrFile << int_comp_names[i] << '\n'; }
525  }
526 
527  bool is_checkpoint_legacy = true; // legacy
528  HdrFile << is_checkpoint_legacy << '\n';
529 
530  // The total number of particles.
531  HdrFile << nparticles << '\n';
532 
533  // The value of nextid that we need to restore on restart.
534  HdrFile << maxnextid << '\n';
535 
536  // Then the finest level of the AMR hierarchy.
537  HdrFile << pc.finestLevel() << '\n';
538 
539  // Then the number of grids at each level.
540  for (int lev = 0; lev <= pc.finestLevel(); lev++) {
541  HdrFile << pc.ParticleBoxArray(lev).size() << '\n';
542  }
543  }
544 
545  // We want to write the data out in parallel.
546  // We'll allow up to nOutFiles active writers at a time.
547  int nOutFiles(256);
548 
549  ParmParse pp("particles");
550  pp.queryAdd("particles_nfiles",nOutFiles);
551  if(nOutFiles == -1) { nOutFiles = NProcs; }
552  nOutFiles = std::max(1, std::min(nOutFiles,NProcs));
553  pc.nOutFilesPrePost = nOutFiles;
554 
555  for (int lev = 0; lev <= pc.finestLevel(); lev++)
556  {
557  bool gotsome;
558  if(pc.usePrePost)
559  {
560  gotsome = (pc.nParticlesAtLevelPrePost[lev] > 0);
561  }
562  else
563  {
564  gotsome = (pc.NumberOfParticlesAtLevel(lev) > 0);
565  }
566 
567  // We store the particles at each level in their own subdirectory.
568  std::string LevelDir = pdir;
569 
570  if (gotsome)
571  {
572  if ( ! LevelDir.empty() && LevelDir[LevelDir.size()-1] != '/') { LevelDir += '/'; }
573 
574  LevelDir = amrex::Concatenate(LevelDir.append("Level_"), lev, 1);
575 
576  if ( ! pc.GetLevelDirectoriesCreated())
577  {
579  if ( ! amrex::UtilCreateDirectory(LevelDir, 0755)) {
581  }
582  }
584  }
585  }
586 
587  // Write out the header for each particle
588  if (gotsome && ParallelDescriptor::IOProcessor()) {
589  std::string HeaderFileName = LevelDir;
590  HeaderFileName += "/Particle_H";
591  std::ofstream ParticleHeader(HeaderFileName);
592 
593  pc.ParticleBoxArray(lev).writeOn(ParticleHeader);
594  ParticleHeader << '\n';
595 
596  ParticleHeader.flush();
597  ParticleHeader.close();
598  }
599 
600  MFInfo info;
601  info.SetAlloc(false);
602  MultiFab state(pc.ParticleBoxArray(lev),
603  pc.ParticleDistributionMap(lev),
604  1,0,info);
605 
606  // We eventually want to write out the file name and the offset
607  // into that file into which each grid of particles is written.
608  Vector<int> which(state.size(),0);
609  Vector<int > count(state.size(),0);
610  Vector<Long> where(state.size(),0);
611 
612  std::string filePrefix(LevelDir);
613  filePrefix += '/';
614  filePrefix += PC::DataPrefix();
615  if(pc.usePrePost) {
616  pc.filePrefixPrePost[lev] = filePrefix;
617  }
618  bool groupSets(false), setBuf(true);
619 
620  if (gotsome)
621  {
622  for(NFilesIter nfi(nOutFiles, filePrefix, groupSets, setBuf); nfi.ReadyToWrite(); ++nfi)
623  {
624  auto& myStream = (std::ofstream&) nfi.Stream();
625  pc.WriteParticles(lev, myStream, nfi.FileNumber(), which, count, where,
626  write_real_comp, write_int_comp, particle_io_flags, is_checkpoint);
627  }
628 
629  if(pc.usePrePost) {
630  pc.whichPrePost[lev] = which;
631  pc.countPrePost[lev] = count;
632  pc.wherePrePost[lev] = where;
633  } else {
634  ParallelDescriptor::ReduceIntSum (which.dataPtr(), static_cast<int>(which.size()), IOProcNumber);
635  ParallelDescriptor::ReduceIntSum (count.dataPtr(), static_cast<int>(count.size()), IOProcNumber);
636  ParallelDescriptor::ReduceLongSum(where.dataPtr(), static_cast<int>(where.size()), IOProcNumber);
637  }
638  }
639 
641  {
642  if(pc.GetUsePrePost()) {
643  // ---- write to the header and unlink in CheckpointPost
644  } else {
645  for (int j = 0; j < state.size(); j++)
646  {
647  HdrFile << which[j] << ' ' << count[j] << ' ' << where[j] << '\n';
648  }
649 
650  if (gotsome && pc.doUnlink)
651  {
652  // Unlink any zero-length data files.
653  Vector<Long> cnt(nOutFiles,0);
654 
655  for (int i = 0, N=static_cast<int>(count.size()); i < N; i++) {
656  cnt[which[i]] += count[i];
657  }
658 
659  for (int i = 0, N=static_cast<int>(cnt.size()); i < N; i++)
660  {
661  if (cnt[i] == 0)
662  {
663  std::string FullFileName = NFilesIter::FileName(i, filePrefix);
664  FileSystem::Remove(FullFileName);
665  }
666  }
667  }
668  }
669  }
670  }
671 
673  {
674  HdrFile.flush();
675  HdrFile.close();
676  if ( ! HdrFile.good())
677  {
678  amrex::Abort("amrex::WriteBinaryParticleDataSync(): problem writing HdrFile");
679  }
680  }
681 }
682 
683 template <class PC, std::enable_if_t<IsParticleContainer<PC>::value, int> foo = 0>
684 void WriteBinaryParticleDataAsync (PC const& pc,
685  const std::string& dir, const std::string& name,
686  const Vector<int>& write_real_comp,
687  const Vector<int>& write_int_comp,
688  const Vector<std::string>& real_comp_names,
689  const Vector<std::string>& int_comp_names, bool is_checkpoint)
690 {
691  BL_PROFILE("WriteBinaryParticleDataAsync");
692  AMREX_ASSERT(pc.OK());
693 
694  AMREX_ASSERT(sizeof(typename PC::ParticleType::RealType) == 4 ||
695  sizeof(typename PC::ParticleType::RealType) == 8);
696 
697  constexpr int NStructReal = PC::NStructReal;
698  constexpr int NStructInt = PC::NStructInt;
699  constexpr int NArrayReal = PC::NArrayReal;
700  constexpr int NArrayInt = PC::NArrayInt;
701 
702  const int MyProc = ParallelDescriptor::MyProc();
703  const int NProcs = ParallelDescriptor::NProcs();
704  const int IOProcNumber = NProcs - 1;
705 
706  if constexpr(PC::ParticleType::is_soa_particle) {
707  AMREX_ALWAYS_ASSERT(real_comp_names.size() == pc.NumRealComps() + NStructReal - AMREX_SPACEDIM); // pure SoA: skip positions
708  } else {
709  AMREX_ALWAYS_ASSERT(real_comp_names.size() == pc.NumRealComps() + NStructReal);
710  }
711  AMREX_ALWAYS_ASSERT( int_comp_names.size() == pc.NumIntComps() + NStructInt);
712 
713  Vector<LayoutData<Long> > np_per_grid_local(pc.finestLevel()+1);
714  for (int lev = 0; lev <= pc.finestLevel(); lev++)
715  {
716  np_per_grid_local[lev].define(pc.ParticleBoxArray(lev), pc.ParticleDistributionMap(lev));
717  using ParIter = typename PC::ParConstIterType;
718  for (ParIter pti(pc, lev); pti.isValid(); ++pti)
719  {
720  int gid = pti.index();
721  const auto& ptile = pc.ParticlesAt(lev, pti);
722  const auto& ptd = ptile.getConstParticleTileData();
723  const int np = ptile.numParticles();
724 
725  ReduceOps<ReduceOpSum> reduce_op;
726  ReduceData<int> reduce_data(reduce_op);
727  using ReduceTuple = typename decltype(reduce_data)::Type;
728 
729  reduce_op.eval(np, reduce_data,
730  [=] AMREX_GPU_DEVICE (int i) -> ReduceTuple
731  {
732  return (ptd.id(i) > 0) ? 1 : 0;
733  });
734 
735  int np_valid = amrex::get<0>(reduce_data.value(reduce_op));
736  np_per_grid_local[lev][gid] += np_valid;
737  }
738  }
739 
740  Vector<Vector<Long> > np_per_grid_global(pc.finestLevel()+1);
741  Long total_np = 0;
742  Vector<Long> np_per_level(pc.finestLevel()+1);
743  for (int lev = 0; lev <= pc.finestLevel(); lev++)
744  {
745  np_per_grid_global[lev].resize(np_per_grid_local[lev].size());
746  ParallelDescriptor::GatherLayoutDataToVector(np_per_grid_local[lev],
747  np_per_grid_global[lev],
748  IOProcNumber);
749  np_per_level[lev] = std::accumulate(np_per_grid_global[lev].begin(),
750  np_per_grid_global[lev].end(), 0L);
751  total_np += np_per_level[lev];
752  }
753 
754  std::string pdir = dir;
755  if ( ! pdir.empty() && pdir[pdir.size()-1] != '/') { pdir += '/'; }
756  pdir += name;
757 
758  if (MyProc == IOProcNumber)
759  {
760  if ( ! pc.GetLevelDirectoriesCreated())
761  {
762  if ( ! amrex::UtilCreateDirectory(pdir, 0755))
763  {
765  }
766  }
767 
768  for (int lev = 0; lev <= pc.finestLevel(); lev++)
769  {
770  std::string LevelDir = pdir;
771  bool gotsome = np_per_level[lev];
772 
773  if (gotsome)
774  {
775  if ( ! LevelDir.empty() && LevelDir[LevelDir.size()-1] != '/') { LevelDir += '/'; }
776 
777  LevelDir = amrex::Concatenate(LevelDir.append("Level_"), lev, 1);
778 
779  if ( ! pc.GetLevelDirectoriesCreated())
780  {
781  if ( ! amrex::UtilCreateDirectory(LevelDir, 0755))
782  {
784  }
785  }
786 
787  std::string HeaderFileName = LevelDir;
788  HeaderFileName += "/Particle_H";
789  std::ofstream ParticleHeader(HeaderFileName);
790 
791  pc.ParticleBoxArray(lev).writeOn(ParticleHeader);
792  ParticleHeader << '\n';
793 
794  ParticleHeader.flush();
795  ParticleHeader.close();
796  }
797  }
798  }
800 
801  Long maxnextid = PC::ParticleType::NextID();
802  ParallelDescriptor::ReduceLongMax(maxnextid, IOProcNumber);
803 
804  Vector<Long> np_on_rank(NProcs, 0L);
805  std::size_t psize = particle_detail::PSizeInFile<ParticleReal>(write_real_comp, write_int_comp);
806  Vector<int64_t> rank_start_offset(NProcs);
807  if (MyProc == IOProcNumber)
808  {
809  for (int lev = 0; lev <= pc.finestLevel(); lev++)
810  {
811  for (int k = 0; k < pc.ParticleBoxArray(lev).size(); ++k)
812  {
813  int rank = pc.ParticleDistributionMap(lev)[k];
814  np_on_rank[rank] += np_per_grid_global[lev][k];
815  }
816  }
817 
818  for (int ip = 0; ip < NProcs; ++ip)
819  {
820  auto info = AsyncOut::GetWriteInfo(ip);
821  rank_start_offset[ip] = (info.ispot == 0) ? 0 : static_cast<int64_t>(rank_start_offset[ip-1] + np_on_rank[ip-1]*psize);
822  }
823  }
824 
825  // make tmp particle tiles in pinned memory to write
826  using PinnedPTile = ParticleTile<typename PC::ParticleType, NArrayReal, NArrayInt,
827  PinnedArenaAllocator>;
828  auto myptiles = std::make_shared<Vector<std::map<std::pair<int, int>,PinnedPTile> > >();
829  myptiles->resize(pc.finestLevel()+1);
830  for (int lev = 0; lev <= pc.finestLevel(); lev++)
831  {
832  for (MFIter mfi = pc.MakeMFIter(lev); mfi.isValid(); ++mfi)
833  {
834  auto& new_ptile = (*myptiles)[lev][std::make_pair(mfi.index(),
835  mfi.LocalTileIndex())];
836 
837  if (np_per_grid_local[lev][mfi.index()] > 0)
838  {
839  const auto& ptile = pc.ParticlesAt(lev, mfi);
840 
841  const auto np = np_per_grid_local[lev][mfi.index()];
842 
843  new_ptile.resize(np);
844 
845  const auto runtime_real_comps = ptile.NumRuntimeRealComps();
846  const auto runtime_int_comps = ptile.NumRuntimeIntComps();
847 
848  new_ptile.define(runtime_real_comps, runtime_int_comps);
849 
850  for (auto comp(0); comp < runtime_real_comps; ++comp) {
851  new_ptile.push_back_real(NArrayReal+comp, np, 0.);
852  }
853 
854  for (auto comp(0); comp < runtime_int_comps; ++comp) {
855  new_ptile.push_back_int(NArrayInt+comp, np, 0);
856  }
857 
858  amrex::filterParticles(new_ptile, ptile, KeepValidFilter());
859  }
860  }
861  }
862 
863  int finest_level = pc.finestLevel();
864  Vector<BoxArray> bas;
865  Vector<DistributionMapping> dms;
866  for (int lev = 0; lev <= pc.finestLevel(); lev++)
867  {
868  bas.push_back(pc.ParticleBoxArray(lev));
869  dms.push_back(pc.ParticleDistributionMap(lev));
870  }
871 
872  int nrc = pc.NumRealComps();
873  int nic = pc.NumIntComps();
874  int rnames_size = (int) real_comp_names.size();
875 
876  auto RD = pc.ParticleRealDescriptor;
877 
878  AsyncOut::Submit([=] ()
879 #if defined(__GNUC__) && (__GNUC__ == 8) && (__GNUC_MINOR__ == 1)
880  mutable // workaround for bug in gcc 8.1
881 #endif
882  {
883  if (MyProc == IOProcNumber)
884  {
885  std::string HdrFileName = pdir;
886  std::ofstream HdrFile;
887 
888  if ( ! HdrFileName.empty() && HdrFileName[HdrFileName.size()-1] != '/') {
889  HdrFileName += '/';
890  }
891 
892  HdrFileName += "Header";
893 
894  HdrFile.open(HdrFileName.c_str(), std::ios::out|std::ios::trunc);
895 
896  if ( ! HdrFile.good()) { amrex::FileOpenFailed(HdrFileName); }
897 
898  std::string version_string = is_checkpoint ? PC::CheckpointVersion() : PC::PlotfileVersion();
899  if (sizeof(typename PC::ParticleType::RealType) == 4)
900  {
901  HdrFile << version_string << "_single" << '\n';
902  }
903  else
904  {
905  HdrFile << version_string << "_double" << '\n';
906  }
907 
908  int num_output_real = 0;
909  for (int i = 0; i < rnames_size; ++i) {
910  if (write_real_comp[i]) { ++num_output_real; }
911  }
912 
913  int num_output_int = 0;
914  for (int i = 0; i < nic + NStructInt; ++i) {
915  if (write_int_comp[i]) { ++num_output_int; }
916  }
917 
918  // AMREX_SPACEDIM and N for sanity checking.
919  HdrFile << AMREX_SPACEDIM << '\n';
920 
921  // The number of extra real parameters
922  HdrFile << num_output_real << '\n';
923 
924  // Real component names
925  for (int i = 0; i < rnames_size; ++i ) {
926  if (write_real_comp[i]) { HdrFile << real_comp_names[i] << '\n'; }
927  }
928 
929  // The number of extra int parameters
930  HdrFile << num_output_int << '\n';
931 
932  // int component names
933  for (int i = 0; i < NStructInt + nic; ++i ) {
934  if (write_int_comp[i]) { HdrFile << int_comp_names[i] << '\n'; }
935  }
936 
937  bool is_checkpoint_legacy = true; // legacy
938  HdrFile << is_checkpoint_legacy << '\n';
939 
940  // The total number of particles.
941  HdrFile << total_np << '\n';
942 
943  // The value of nextid that we need to restore on restart.
944  HdrFile << maxnextid << '\n';
945 
946  // Then the finest level of the AMR hierarchy.
947  HdrFile << finest_level << '\n';
948 
949  // Then the number of grids at each level.
950  for (int lev = 0; lev <= finest_level; lev++) {
951  HdrFile << dms[lev].size() << '\n';
952  }
953 
954  for (int lev = 0; lev <= finest_level; lev++)
955  {
956  Vector<int64_t> grid_offset(NProcs, 0);
957  for (int k = 0; k < bas[lev].size(); ++k)
958  {
959  int rank = dms[lev][k];
960  auto info = AsyncOut::GetWriteInfo(rank);
961  HdrFile << info.ifile << ' '
962  << np_per_grid_global[lev][k] << ' '
963  << grid_offset[rank] + rank_start_offset[rank] << '\n';
964  grid_offset[rank] += static_cast<int64_t>(np_per_grid_global[lev][k]*psize);
965  }
966  }
967 
968  HdrFile.flush();
969  HdrFile.close();
970  if ( ! HdrFile.good())
971  {
972  amrex::Abort("amrex::WriteBinaryParticleDataAsync(): problem writing HdrFile");
973  }
974  }
975 
976  AsyncOut::Wait(); // Wait for my turn
977 
978  for (int lev = 0; lev <= finest_level; lev++)
979  {
980  // For a each grid, the tiles it contains
981  std::map<int, Vector<int> > tile_map;
982 
983  for (const auto& kv : (*myptiles)[lev])
984  {
985  const int grid = kv.first.first;
986  const int tile = kv.first.second;
987  tile_map[grid].push_back(tile);
988  }
989 
990  std::string LevelDir = pdir;
991  if ( ! LevelDir.empty() && LevelDir[LevelDir.size()-1] != '/') { LevelDir += '/'; }
992  LevelDir = amrex::Concatenate(LevelDir.append("Level_"), lev, 1);
993  std::string filePrefix(LevelDir);
994  filePrefix += '/';
995  filePrefix += PC::DataPrefix();
996  auto info = AsyncOut::GetWriteInfo(MyProc);
997  std::string file_name = amrex::Concatenate(filePrefix, info.ifile, 5);
998  std::ofstream ofs;
999  ofs.open(file_name.c_str(), (info.ispot == 0) ? (std::ios::binary | std::ios::trunc)
1000  : (std::ios::binary | std::ios::app));
1001 
1002  for (int k = 0; k < bas[lev].size(); ++k)
1003  {
1004  int rank = dms[lev][k];
1005  if (rank != MyProc) { continue; }
1006  const int grid = k;
1007  if (np_per_grid_local[lev][grid] == 0) { continue; }
1008 
1009  // First write out the integer data in binary.
1010  int num_output_int = 0;
1011  for (int i = 0; i < nic + NStructInt; ++i) {
1012  if (write_int_comp[i]) { ++num_output_int; }
1013  }
1014 
1015  const Long iChunkSize = 2 + num_output_int;
1016  Vector<int> istuff(np_per_grid_local[lev][grid]*iChunkSize);
1017  int* iptr = istuff.dataPtr();
1018 
1019  for (unsigned i = 0; i < tile_map[grid].size(); i++) {
1020  auto ptile_index = std::make_pair(grid, tile_map[grid][i]);
1021  const auto& pbox = (*myptiles)[lev][ptile_index];
1022  const auto& ptd = pbox.getConstParticleTileData();
1023  for (int pindex = 0; pindex < pbox.numParticles(); ++pindex)
1024  {
1025  const auto& soa = pbox.GetStructOfArrays();
1026 
1028  const auto& p = mp(ptd, pindex);
1029  if (p.id() <= 0) { continue; }
1030 
1031  // note: for pure SoA particle layouts, we do write the id, cpu and positions as a struct
1032  // for backwards compatibility with readers
1033  if constexpr(!PC::ParticleType::is_soa_particle)
1034  {
1035  // Ints: id, cpu
1036  particle_detail::packParticleIDs(iptr, p, is_checkpoint);
1037  iptr += 2;
1038 
1039  // extra AoS Int components
1040  for (int j = 0; j < NStructInt; j++)
1041  {
1042  if (write_int_comp[j])
1043  {
1044  *iptr = p.idata(j);
1045  ++iptr;
1046  }
1047  }
1048  }
1049  else {
1050  uint64_t idcpu = soa.GetIdCPUData()[pindex];
1051  if (is_checkpoint) {
1052  std::int32_t xi, yi;
1053  std::uint32_t xu, yu;
1054  xu = (std::uint32_t)((idcpu & 0xFFFFFFFF00000000LL) >> 32);
1055  yu = (std::uint32_t)( idcpu & 0xFFFFFFFFLL);
1056  std::memcpy(&xi, &xu, sizeof(xu));
1057  std::memcpy(&yi, &yu, sizeof(yu));
1058  *iptr = xi;
1059  iptr += 1;
1060  *iptr = yi;
1061  iptr += 1;
1062  } else {
1063  // Int: id, cpu
1064  *iptr = (int) ParticleIDWrapper(idcpu);
1065  iptr += 1;
1066  *iptr = (int) ParticleCPUWrapper(idcpu);
1067  iptr += 1;
1068  }
1069  }
1070 
1071  // extra SoA Ints
1072  const int int_start_offset = 0;
1073  for (int j = int_start_offset; j < nic; j++)
1074  {
1075  if (write_int_comp[NStructInt+j])
1076  {
1077  *iptr = soa.GetIntData(j)[pindex];
1078  ++iptr;
1079  }
1080  }
1081  }
1082  }
1083 
1084  writeIntData(istuff.dataPtr(), istuff.size(), ofs);
1085  ofs.flush(); // Some systems require this flush() (probably due to a bug)
1086 
1087  // Write the Real data in binary.
1088  int num_output_real = 0;
1089  for (int i = 0; i < rnames_size; ++i) {
1090  if (write_real_comp[i]) { ++num_output_real; }
1091  }
1092 
1093  const Long rChunkSize = AMREX_SPACEDIM + num_output_real;
1094  Vector<typename PC::ParticleType::RealType> rstuff(np_per_grid_local[lev][grid]*rChunkSize);
1095  typename PC::ParticleType::RealType* rptr = rstuff.dataPtr();
1096 
1097  for (unsigned i = 0; i < tile_map[grid].size(); i++) {
1098  auto ptile_index = std::make_pair(grid, tile_map[grid][i]);
1099  const auto& pbox = (*myptiles)[lev][ptile_index];
1100  const auto& ptd = pbox.getConstParticleTileData();
1101  for (int pindex = 0;
1102  pindex < pbox.numParticles(); ++pindex)
1103  {
1104  const auto& soa = pbox.GetStructOfArrays();
1106  const auto& p = mp(ptd, pindex);
1107 
1108  if (p.id() <= 0) { continue; }
1109 
1110  if constexpr(!PC::ParticleType::is_soa_particle)
1111  {
1112  // Real: position
1113  for (int j = 0; j < AMREX_SPACEDIM; j++) { rptr[j] = p.pos(j); }
1114  rptr += AMREX_SPACEDIM;
1115 
1116  // extra AoS real
1117  for (int j = 0; j < NStructReal; j++)
1118  {
1119  if (write_real_comp[j])
1120  {
1121  *rptr = p.rdata(j);
1122  ++rptr;
1123  }
1124  }
1125  }
1126  else {
1127  // Real: position
1128  for (int j = 0; j < AMREX_SPACEDIM; j++) { rptr[j] = soa.GetRealData(j)[pindex]; }
1129  rptr += AMREX_SPACEDIM;
1130  }
1131 
1132  // extra SoA real
1133  const int real_start_offset = PC::ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0; // pure SoA: positions
1134  for (int j = real_start_offset; j < nrc; j++)
1135  {
1136  const int write_comp_offset = PC::ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0; // pure SoA: skip positions
1137  const int write_comp_index = PC::NStructReal+j-write_comp_offset;
1138  if (write_real_comp[write_comp_index])
1139  {
1140  *rptr = (typename PC::ParticleType::RealType) soa.GetRealData(j)[pindex];
1141  ++rptr;
1142  }
1143  }
1144  }
1145  }
1146 
1147  if (sizeof(typename PC::ParticleType::RealType) == 4) {
1148  writeFloatData((float*) rstuff.dataPtr(), rstuff.size(), ofs, RD);
1149  }
1150  else if (sizeof(typename PC::ParticleType::RealType) == 8) {
1151  writeDoubleData((double*) rstuff.dataPtr(), rstuff.size(), ofs, RD);
1152  }
1153 
1154  ofs.flush(); // Some systems require this flush() (probably due to a bug)
1155  }
1156  }
1157  AsyncOut::Notify(); // Notify others I am done
1158  });
1159 }
1160 
1161 #ifdef AMREX_USE_HDF5
1163 #endif
1164 
1165 #endif /*AMREX_WRITE_BINARY_PARTICLE_DATA_H*/
#define BL_PROFILE(a)
Definition: AMReX_BLProfiler.H:551
#define AMREX_ASSERT(EX)
Definition: AMReX_BLassert.H:38
#define AMREX_ALWAYS_ASSERT(EX)
Definition: AMReX_BLassert.H:50
#define AMREX_GPU_DEVICE
Definition: AMReX_GpuQualifiers.H:18
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
amrex::ParmParse pp
Input file parser instance for the given namespace.
Definition: AMReX_HypreIJIface.cpp:15
void WriteBinaryParticleDataAsync(PC const &pc, const std::string &dir, const std::string &name, const Vector< int > &write_real_comp, const Vector< int > &write_int_comp, const Vector< std::string > &real_comp_names, const Vector< std::string > &int_comp_names, bool is_checkpoint)
Definition: AMReX_WriteBinaryParticleData.H:684
void WriteBinaryParticleDataSync(PC const &pc, const std::string &dir, const std::string &name, const Vector< int > &write_real_comp, const Vector< int > &write_int_comp, const Vector< std::string > &real_comp_names, const Vector< std::string > &int_comp_names, F const &f, bool is_checkpoint)
Definition: AMReX_WriteBinaryParticleData.H:392
if(!(yy_init))
Definition: amrex_iparser.lex.nolint.H:935
int queryAdd(const char *name, T &ref)
If name is found, the value in the ParmParse database will be stored in the ref argument....
Definition: AMReX_ParmParse.H:993
WriteInfo GetWriteInfo(int rank)
Definition: AMReX_AsyncOut.cpp:72
void Wait()
Definition: AMReX_AsyncOut.cpp:112
void Notify()
Definition: AMReX_AsyncOut.cpp:127
void Submit(std::function< void()> &&a_f)
Definition: AMReX_AsyncOut.cpp:95
bool Remove(std::string const &filename)
Definition: AMReX_FileSystem.cpp:190
AMREX_GPU_HOST_DEVICE Long size(T const &b) noexcept
integer version
Definition: AMReX_GpuRange.H:26
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 DeviceToHost deviceToHost
Definition: AMReX_GpuContainers.H:99
static constexpr HostToDevice hostToDevice
Definition: AMReX_GpuContainers.H:98
void streamSynchronize() noexcept
Definition: AMReX_GpuDevice.H:237
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void * memcpy(void *dest, const void *src, std::size_t count)
Definition: AMReX_GpuUtility.H:214
int NProcs()
Process ID in MPI_COMM_WORLD.
Definition: AMReX_MPMD.cpp:122
int MyProc()
Definition: AMReX_MPMD.cpp:117
void ReduceIntSum(int &)
Integer sum reduction.
Definition: AMReX_ParallelDescriptor.cpp:1252
void ReduceLongMax(Long &)
Long max reduction.
Definition: AMReX_ParallelDescriptor.cpp:1224
void ReduceLongSum(Long &)
Long sum reduction.
Definition: AMReX_ParallelDescriptor.cpp:1223
void GatherLayoutDataToVector(const LayoutData< T > &sendbuf, Vector< T > &recvbuf, int root)
Gather LayoutData values to a vector on root.
Definition: AMReX_ParallelDescriptor.H:1211
int IOProcessorNumber() noexcept
Definition: AMReX_ParallelDescriptor.H:266
void Barrier(const std::string &)
Definition: AMReX_ParallelDescriptor.cpp:1202
bool IOProcessor() noexcept
Is this CPU the I/O Processor? To get the rank number, call IOProcessorNumber()
Definition: AMReX_ParallelDescriptor.H:275
T ExclusiveSum(N n, T const *in, T *out, RetSum a_ret_sum=retSum)
Definition: AMReX_Scan.H:1229
static constexpr RetSum retSum
Definition: AMReX_Scan.H:29
static int f(amrex::Real t, N_Vector y_data, N_Vector y_rhs, void *user_data)
Definition: AMReX_SundialsIntegrator.H:44
@ min
Definition: AMReX_ParallelReduce.H:18
@ max
Definition: AMReX_ParallelReduce.H:17
static constexpr int P
Definition: AMReX_OpenBC.H:14
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
void writeIntData(const From *data, std::size_t size, std::ostream &os, const amrex::IntDescriptor &id)
Definition: AMReX_IntConv.H:23
void FileOpenFailed(const std::string &file)
Output a message and abort when couldn't open the file.
Definition: AMReX_Utility.cpp:131
ParIter_impl< Particle< T_NStructReal, T_NStructInt >, T_NArrayReal, T_NArrayInt, Allocator, CellAssignor > ParIter
Definition: AMReX_ParIter.H:229
void writeFloatData(const float *data, std::size_t size, std::ostream &os, const RealDescriptor &rd=FPC::Native32RealDescriptor())
Definition: AMReX_VectorIO.cpp:114
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 end(BoxND< dim > const &box) noexcept
Definition: AMReX_Box.H:1890
std::string Concatenate(const std::string &root, int num, int mindigits)
Returns rootNNNN where NNNN == num.
Definition: AMReX_String.cpp:34
Index filterParticles(DstTile &dst, const SrcTile &src, const Index *mask) noexcept
Conditionally copy particles from src to dst based on the value of mask.
Definition: AMReX_ParticleTransformation.H:325
void CreateDirectoryFailed(const std::string &dir)
Output a message and abort when couldn't create the directory.
Definition: AMReX_Utility.cpp:123
bool UtilCreateDirectory(const std::string &path, mode_t mode, bool verbose=false)
Creates the specified directories. path may be either a full pathname or a relative pathname....
Definition: AMReX_Utility.cpp:116
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void ignore_unused(const Ts &...)
This shuts up the compiler about unused variables.
Definition: AMReX.H:111
void writeDoubleData(const double *data, std::size_t size, std::ostream &os, const RealDescriptor &rd=FPC::Native64RealDescriptor())
Definition: AMReX_VectorIO.cpp:126
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 begin(BoxND< dim > const &box) noexcept
Definition: AMReX_Box.H:1881
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
AMREX_ATTRIBUTE_FLATTEN_FOR void ParallelForRNG(T n, L const &f) noexcept
Definition: AMReX_GpuLaunchFunctsC.H:1221
Definition: AMReX_WriteBinaryParticleData.H:19
std::enable_if_t< RunOnGpu< typename Container< int, Allocator >::allocator_type >::value > fillFlags(Container< int, Allocator > &pflags, const PTile &ptile, F const &f)
Definition: AMReX_WriteBinaryParticleData.H:34
AMREX_GPU_HOST_DEVICE void packParticleIDs(I *idata, const P &p, bool is_checkpoint) noexcept
Definition: AMReX_WriteBinaryParticleData.H:152
std::enable_if_t< RunOnGpu< typename PC::template AllocatorType< int > >::value > packIOData(Vector< int > &idata, Vector< ParticleReal > &rdata, const PC &pc, int lev, int grid, const Vector< int > &write_real_comp, const Vector< int > &write_int_comp, const Vector< std::map< std::pair< int, int >, typename PC::IntVector >> &particle_io_flags, const Vector< int > &tiles, int np, bool is_checkpoint)
Definition: AMReX_WriteBinaryParticleData.H:171
std::enable_if_t< RunOnGpu< typename Container< int, Allocator >::allocator_type >::value, amrex::Long > countFlags(const Vector< std::map< std::pair< int, int >, Container< int, Allocator >>> &particle_io_flags, const PC &pc)
Definition: AMReX_WriteBinaryParticleData.H:76
std::size_t PSizeInFile(const Vector< int > &wrc, const Vector< int > &wic)
Definition: AMReX_WriteBinaryParticleData.H:22
Definition: AMReX_WriteBinaryParticleData.H:10
AMREX_GPU_HOST_DEVICE int operator()(const SrcData &src, int i) const noexcept
Definition: AMReX_WriteBinaryParticleData.H:13
Definition: AMReX_RandomEngine.H:57
Definition: AMReX_MakeParticle.H:16