Block-Structured AMR Software Framework
AMReX_ParticleIO.H
Go to the documentation of this file.
1 #ifndef AMREX_PARTICLEIO_H
2 #define AMREX_PARTICLEIO_H
3 #include <AMReX_Config.H>
4 
6 
7 template <typename ParticleType, int NArrayReal, int NArrayInt,
8  template<class> class Allocator, class CellAssignor>
9 void
10 ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt, Allocator, CellAssignor>
11 ::WriteParticleRealData (void* data, size_t size, std::ostream& os) const
12 {
13  if (sizeof(typename ParticleType::RealType) == 4) {
14  writeFloatData((float*) data, size, os, ParticleRealDescriptor);
15  }
16  else if (sizeof(typename ParticleType::RealType) == 8) {
17  writeDoubleData((double*) data, size, os, ParticleRealDescriptor);
18  }
19 }
20 
21 template <typename ParticleType, int NArrayReal, int NArrayInt,
22  template<class> class Allocator, class CellAssignor>
23 void
25 ::ReadParticleRealData (void* data, size_t size, std::istream& is)
26 {
27  if (sizeof(typename ParticleType::RealType) == 4) {
28  readFloatData((float*) data, size, is, ParticleRealDescriptor);
29  }
30  else if (sizeof(typename ParticleType::RealType) == 8) {
31  readDoubleData((double*) data, size, is, ParticleRealDescriptor);
32  }
33 }
34 
36  template <typename P>
38  int operator() (const P& p) const {
39  return p.id() > 0;
40  }
41 };
42 
43 template <typename ParticleType, int NArrayReal, int NArrayInt,
44  template<class> class Allocator, class CellAssignor>
45 void
46 ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt, Allocator, CellAssignor>
47 ::Checkpoint (const std::string& dir,
48  const std::string& name, bool /*is_checkpoint*/,
49  const Vector<std::string>& real_comp_names,
50  const Vector<std::string>& int_comp_names) const
51 {
52  Vector<int> write_real_comp;
53  Vector<std::string> tmp_real_comp_names;
54 
55  int first_rcomp = ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0;
56  for (int i = first_rcomp; i < NStructReal + NumRealComps(); ++i )
57  {
58  write_real_comp.push_back(1);
59  if (real_comp_names.empty())
60  {
61  tmp_real_comp_names.push_back(getDefaultCompNameReal<ParticleType>(i));
62  }
63  else
64  {
65  tmp_real_comp_names.push_back(real_comp_names[i-first_rcomp]);
66  }
67  }
68 
69  Vector<int> write_int_comp;
70  Vector<std::string> tmp_int_comp_names;
71  for (int i = 0; i < NStructInt + NumIntComps(); ++i )
72  {
73  write_int_comp.push_back(1);
74  if (int_comp_names.empty())
75  {
76  tmp_int_comp_names.push_back(getDefaultCompNameInt<ParticleType>(i));
77  }
78  else
79  {
80  tmp_int_comp_names.push_back(int_comp_names[i]);
81  }
82  }
83 
84  WriteBinaryParticleData(dir, name, write_real_comp, write_int_comp,
85  tmp_real_comp_names, tmp_int_comp_names,
86  FilterPositiveID{}, true);
87 }
88 
89 template <typename ParticleType, int NArrayReal, int NArrayInt,
90  template<class> class Allocator, class CellAssignor>
91 void
93 ::WritePlotFile (const std::string& dir, const std::string& name) const
94 {
95  Vector<int> write_real_comp;
96  Vector<std::string> real_comp_names;
97 
98  int first_rcomp = ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0;
99  for (int i = first_rcomp; i < NStructReal + NumRealComps(); ++i )
100  {
101  write_real_comp.push_back(1);
102  real_comp_names.push_back(getDefaultCompNameReal<ParticleType>(i));
103  }
104 
105  Vector<int> write_int_comp;
106  Vector<std::string> int_comp_names;
107  for (int i = 0; i < NStructInt + NumIntComps(); ++i )
108  {
109  write_int_comp.push_back(1);
110  int_comp_names.push_back(getDefaultCompNameInt<ParticleType>(i));
111  }
112 
113  WriteBinaryParticleData(dir, name, write_real_comp, write_int_comp,
114  real_comp_names, int_comp_names,
115  FilterPositiveID{});
116 }
117 
118 template <typename ParticleType, int NArrayReal, int NArrayInt,
119  template<class> class Allocator, class CellAssignor>
120 void
122 ::WritePlotFile (const std::string& dir, const std::string& name,
123  const Vector<std::string>& real_comp_names,
124  const Vector<std::string>& int_comp_names) const
125 {
126  if constexpr(ParticleType::is_soa_particle) {
127  AMREX_ALWAYS_ASSERT(real_comp_names.size() == NumRealComps() + NStructReal - AMREX_SPACEDIM); // pure SoA: skip positions
128  } else {
129  AMREX_ALWAYS_ASSERT(real_comp_names.size() == NumRealComps() + NStructReal);
130  }
131  AMREX_ASSERT( int_comp_names.size() == NStructInt + NumIntComps() );
132 
133  Vector<int> write_real_comp;
134  int nrc = ParticleType::is_soa_particle ? NStructReal + NumRealComps() - AMREX_SPACEDIM : NStructReal + NumRealComps();
135  for (int i = 0; i < nrc; ++i) {
136  write_real_comp.push_back(1);
137  }
138 
139  Vector<int> write_int_comp;
140  for (int i = 0; i < NStructInt + NumIntComps(); ++i) {
141  write_int_comp.push_back(1);
142  }
143 
144  WriteBinaryParticleData(dir, name,
145  write_real_comp, write_int_comp,
146  real_comp_names, int_comp_names,
147  FilterPositiveID{});
148 }
149 
150 template <typename ParticleType, int NArrayReal, int NArrayInt,
151  template<class> class Allocator, class CellAssignor>
152 void
154 ::WritePlotFile (const std::string& dir, const std::string& name,
155  const Vector<std::string>& real_comp_names) const
156 {
157  if constexpr(ParticleType::is_soa_particle) {
158  AMREX_ALWAYS_ASSERT(real_comp_names.size() == NumRealComps() + NStructReal - AMREX_SPACEDIM); // pure SoA: skip positions
159  } else {
160  AMREX_ALWAYS_ASSERT(real_comp_names.size() == NumRealComps() + NStructReal);
161  }
162 
163  Vector<int> write_real_comp;
164  int nrc = ParticleType::is_soa_particle ? NStructReal + NumRealComps() - AMREX_SPACEDIM : NStructReal + NumRealComps();
165  for (int i = 0; i < nrc; ++i) {
166  write_real_comp.push_back(1);
167  }
168 
169  Vector<int> write_int_comp;
170  for (int i = 0; i < NStructInt + NumIntComps(); ++i) {
171  write_int_comp.push_back(1);
172  }
173 
174  Vector<std::string> int_comp_names;
175  for (int i = 0; i < NStructInt + NumIntComps(); ++i )
176  {
177  int_comp_names.push_back(getDefaultCompNameInt<ParticleType>(i));
178  }
179 
180  WriteBinaryParticleData(dir, name,
181  write_real_comp, write_int_comp,
182  real_comp_names, int_comp_names,
183  FilterPositiveID{});
184 }
185 
186 template <typename ParticleType, int NArrayReal, int NArrayInt,
187  template<class> class Allocator, class CellAssignor>
188 void
190 ::WritePlotFile (const std::string& dir,
191  const std::string& name,
192  const Vector<int>& write_real_comp,
193  const Vector<int>& write_int_comp) const
194 {
195 
196  if constexpr(ParticleType::is_soa_particle) {
197  AMREX_ALWAYS_ASSERT(write_real_comp.size() == NumRealComps() + NStructReal - AMREX_SPACEDIM); // pure SoA: skip positions
198  } else {
199  AMREX_ALWAYS_ASSERT(write_real_comp.size() == NumRealComps() + NStructReal);
200  }
201  AMREX_ASSERT(write_int_comp.size() == NStructInt + NArrayInt );
202 
203  Vector<std::string> real_comp_names;
204  int first_rcomp = ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0;
205  for (int i = first_rcomp; i < NStructReal + NumRealComps(); ++i )
206  {
207  real_comp_names.push_back(getDefaultCompNameReal<ParticleType>(i));
208  }
209 
210  Vector<std::string> int_comp_names;
211  for (int i = 0; i < NStructInt + NumIntComps(); ++i )
212  {
213  int_comp_names.push_back(getDefaultCompNameInt<ParticleType>(i));
214  }
215 
216  WriteBinaryParticleData(dir, name, write_real_comp, write_int_comp,
217  real_comp_names, int_comp_names,
218  FilterPositiveID{});
219 }
220 
221 template <typename ParticleType, int NArrayReal, int NArrayInt,
222  template<class> class Allocator, class CellAssignor>
223 void
225 WritePlotFile (const std::string& dir, const std::string& name,
226  const Vector<int>& write_real_comp,
227  const Vector<int>& write_int_comp,
228  const Vector<std::string>& real_comp_names,
229  const Vector<std::string>& int_comp_names) const
230 {
231  BL_PROFILE("ParticleContainer::WritePlotFile()");
232 
233  WriteBinaryParticleData(dir, name,
234  write_real_comp, write_int_comp,
235  real_comp_names, int_comp_names,
236  FilterPositiveID{});
237 }
238 
239 template <typename ParticleType, int NArrayReal, int NArrayInt,
240  template<class> class Allocator, class CellAssignor>
241 template <class F, std::enable_if_t<!std::is_same_v<std::decay_t<F>, Vector<std::string>>>*>
242 void
244 ::WritePlotFile (const std::string& dir, const std::string& name, F&& f) const
245 {
246  Vector<int> write_real_comp;
247  Vector<std::string> real_comp_names;
248 
249  int first_rcomp = ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0;
250  for (int i = first_rcomp; i < NStructReal + NumRealComps(); ++i )
251  {
252  write_real_comp.push_back(1);
253  real_comp_names.push_back(getDefaultCompNameReal<ParticleType>(i));
254  }
255 
256  Vector<int> write_int_comp;
257  Vector<std::string> int_comp_names;
258  for (int i = 0; i < NStructInt + NumIntComps(); ++i )
259  {
260  write_int_comp.push_back(1);
261  int_comp_names.push_back(getDefaultCompNameInt<ParticleType>(i));
262  }
263 
264  WriteBinaryParticleData(dir, name, write_real_comp, write_int_comp,
265  real_comp_names, int_comp_names,
266  std::forward<F>(f));
267 }
268 
269 template <typename ParticleType, int NArrayReal, int NArrayInt,
270  template<class> class Allocator, class CellAssignor>
271 template <class F>
272 void
274 ::WritePlotFile (const std::string& dir, const std::string& name,
275  const Vector<std::string>& real_comp_names,
276  const Vector<std::string>& int_comp_names, F&& f) const
277 {
278  if constexpr(ParticleType::is_soa_particle) {
279  AMREX_ALWAYS_ASSERT(real_comp_names.size() == NumRealComps() + NStructReal - AMREX_SPACEDIM); // pure SoA: skip positions
280  } else {
281  AMREX_ALWAYS_ASSERT(real_comp_names.size() == NumRealComps() + NStructReal);
282  }
283  AMREX_ASSERT( int_comp_names.size() == NStructInt + NArrayInt );
284 
285  Vector<int> write_real_comp;
286  int nrc = ParticleType::is_soa_particle ? NStructReal + NumRealComps() - AMREX_SPACEDIM : NStructReal + NumRealComps();
287  for (int i = 0; i < nrc; ++i) {
288  write_real_comp.push_back(1);
289  }
290 
291  Vector<int> write_int_comp;
292  for (int i = 0; i < NStructInt + NumIntComps(); ++i) {
293  write_int_comp.push_back(1);
294  }
295 
296  WriteBinaryParticleData(dir, name,
297  write_real_comp, write_int_comp,
298  real_comp_names, int_comp_names,
299  std::forward<F>(f));
300 }
301 
302 template <typename ParticleType, int NArrayReal, int NArrayInt,
303  template<class> class Allocator, class CellAssignor>
304 template <class F, std::enable_if_t<!std::is_same_v<std::decay_t<F>, Vector<std::string>>>*>
305 void
307 ::WritePlotFile (const std::string& dir, const std::string& name,
308  const Vector<std::string>& real_comp_names, F&& f) const
309 {
310  if constexpr(ParticleType::is_soa_particle) {
311  AMREX_ALWAYS_ASSERT(real_comp_names.size() == NumRealComps() + NStructReal - AMREX_SPACEDIM); // pure SoA: skip positions
312  } else {
313  AMREX_ALWAYS_ASSERT(real_comp_names.size() == NumRealComps() + NStructReal);
314  }
315 
316  Vector<int> write_real_comp;
317  int nrc = ParticleType::is_soa_particle ? NStructReal + NumRealComps() - AMREX_SPACEDIM : NStructReal + NumRealComps();
318  for (int i = 0; i < nrc; ++i) {
319  write_real_comp.push_back(1);
320  }
321 
322  Vector<int> write_int_comp;
323  for (int i = 0; i < NStructInt + NumIntComps(); ++i) {
324  write_int_comp.push_back(1);
325  }
326 
327  Vector<std::string> int_comp_names;
328  for (int i = 0; i < NStructInt + NumIntComps(); ++i )
329  {
330  int_comp_names.push_back(getDefaultCompNameInt<ParticleType>(i));
331  }
332 
333  WriteBinaryParticleData(dir, name,
334  write_real_comp, write_int_comp,
335  real_comp_names, int_comp_names,
336  std::forward<F>(f));
337 }
338 
339 template <typename ParticleType, int NArrayReal, int NArrayInt,
340  template<class> class Allocator, class CellAssignor>
341 template <class F>
342 void
344 ::WritePlotFile (const std::string& dir,
345  const std::string& name,
346  const Vector<int>& write_real_comp,
347  const Vector<int>& write_int_comp, F&& f) const
348 {
349  if constexpr(ParticleType::is_soa_particle) {
350  AMREX_ALWAYS_ASSERT(write_real_comp.size() == NumRealComps() + NStructReal - AMREX_SPACEDIM); // pure SoA: skip positions
351  } else {
352  AMREX_ALWAYS_ASSERT(write_real_comp.size() == NumRealComps() + NStructReal);
353  }
354  AMREX_ASSERT(write_int_comp.size() == NStructInt + NumIntComps() );
355 
356  Vector<std::string> real_comp_names;
357  int first_rcomp = ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0;
358  for (int i = first_rcomp; i < NStructReal + NumRealComps(); ++i )
359  {
360  real_comp_names.push_back(getDefaultCompNameReal<ParticleType>(i));
361  }
362 
363  Vector<std::string> int_comp_names;
364  for (int i = 0; i < NStructInt + NumIntComps(); ++i )
365  {
366  int_comp_names.push_back(getDefaultCompNameInt<ParticleType>(i));
367  }
368 
369  WriteBinaryParticleData(dir, name, write_real_comp, write_int_comp,
370  real_comp_names, int_comp_names,
371  std::forward<F>(f));
372 }
373 
374 template <typename ParticleType, int NArrayReal, int NArrayInt,
375  template<class> class Allocator, class CellAssignor>
376 template <class F>
377 void
379 WritePlotFile (const std::string& dir, const std::string& name,
380  const Vector<int>& write_real_comp,
381  const Vector<int>& write_int_comp,
382  const Vector<std::string>& real_comp_names,
383  const Vector<std::string>& int_comp_names,
384  F&& f) const
385 {
386  BL_PROFILE("ParticleContainer::WritePlotFile()");
387 
388  WriteBinaryParticleData(dir, name,
389  write_real_comp, write_int_comp,
390  real_comp_names, int_comp_names,
391  std::forward<F>(f));
392 }
393 
394 template <typename ParticleType, int NArrayReal, int NArrayInt,
395  template<class> class Allocator, class CellAssignor>
396 template <class F>
397 void
399 ::WriteBinaryParticleData (const std::string& dir, const std::string& name,
400  const Vector<int>& write_real_comp,
401  const Vector<int>& write_int_comp,
402  const Vector<std::string>& real_comp_names,
403  const Vector<std::string>& int_comp_names,
404  F&& f, bool is_checkpoint) const
405 {
406  if (AsyncOut::UseAsyncOut()) {
407  WriteBinaryParticleDataAsync(*this, dir, name,
408  write_real_comp, write_int_comp,
409  real_comp_names, int_comp_names, is_checkpoint);
410  } else
411  {
412  WriteBinaryParticleDataSync(*this, dir, name,
413  write_real_comp, write_int_comp,
414  real_comp_names, int_comp_names,
415  std::forward<F>(f), is_checkpoint);
416  }
417 }
418 
419 template <typename ParticleType, int NArrayReal, int NArrayInt,
420  template<class> class Allocator, class CellAssignor>
421 void
424 {
425  if( ! usePrePost) {
426  return;
427  }
428 
429  BL_PROFILE("ParticleContainer::CheckpointPre()");
430 
431  const int IOProcNumber = ParallelDescriptor::IOProcessorNumber();
432  Long nparticles = 0;
433  Long maxnextid = ParticleType::NextID();
434 
435  for (int lev = 0; lev < m_particles.size(); lev++) {
436  const auto& pmap = m_particles[lev];
437  for (const auto& kv : pmap) {
438  const auto& aos = kv.second.GetArrayOfStructs();
439  for (int k = 0; k < aos.numParticles(); ++k) {
440  const ParticleType& p = aos[k];
441  if (p.id() > 0) {
442  //
443  // Only count (and checkpoint) valid particles.
444  //
445  nparticles++;
446  }
447  }
448  }
449  }
450  ParallelDescriptor::ReduceLongSum(nparticles, IOProcNumber);
451 
452  ParticleType::NextID(maxnextid);
453  ParallelDescriptor::ReduceLongMax(maxnextid, IOProcNumber);
454 
455  nparticlesPrePost = nparticles;
456  maxnextidPrePost = int(maxnextid);
457 
458  nParticlesAtLevelPrePost.clear();
459  nParticlesAtLevelPrePost.resize(finestLevel() + 1, 0);
460  for(int lev(0); lev <= finestLevel(); ++lev) {
461  nParticlesAtLevelPrePost[lev] = NumberOfParticlesAtLevel(lev);
462  }
463 
464  whichPrePost.clear();
465  whichPrePost.resize(finestLevel() + 1);
466  countPrePost.clear();
467  countPrePost.resize(finestLevel() + 1);
468  wherePrePost.clear();
469  wherePrePost.resize(finestLevel() + 1);
470 
471  filePrefixPrePost.clear();
472  filePrefixPrePost.resize(finestLevel() + 1);
473 }
474 
475 
476 template <typename ParticleType, int NArrayReal, int NArrayInt,
477  template<class> class Allocator, class CellAssignor>
478 void
481 {
482  if( ! usePrePost) {
483  return;
484  }
485 
486  BL_PROFILE("ParticleContainer::CheckpointPost()");
487 
488  const int IOProcNumber = ParallelDescriptor::IOProcessorNumber();
489  std::ofstream HdrFile;
490  HdrFile.open(HdrFileNamePrePost.c_str(), std::ios::out | std::ios::app);
491 
492  for(int lev(0); lev <= finestLevel(); ++lev) {
493  ParallelDescriptor::ReduceIntSum (whichPrePost[lev].dataPtr(), whichPrePost[lev].size(), IOProcNumber);
494  ParallelDescriptor::ReduceIntSum (countPrePost[lev].dataPtr(), countPrePost[lev].size(), IOProcNumber);
495  ParallelDescriptor::ReduceLongSum(wherePrePost[lev].dataPtr(), wherePrePost[lev].size(), IOProcNumber);
496 
497 
499  for(int j(0); j < whichPrePost[lev].size(); ++j) {
500  HdrFile << whichPrePost[lev][j] << ' ' << countPrePost[lev][j] << ' ' << wherePrePost[lev][j] << '\n';
501  }
502 
503  const bool gotsome = (nParticlesAtLevelPrePost[lev] > 0);
504  if(gotsome && doUnlink) {
505 // BL_PROFILE_VAR("PC<NNNN>::Checkpoint:unlink", unlink_post);
506  // Unlink any zero-length data files.
507  Vector<Long> cnt(nOutFilesPrePost,0);
508 
509  for(int i(0), N = countPrePost[lev].size(); i < N; ++i) {
510  cnt[whichPrePost[lev][i]] += countPrePost[lev][i];
511  }
512 
513  for(int i(0), N = int(cnt.size()); i < N; ++i) {
514  if(cnt[i] == 0) {
515  std::string FullFileName = NFilesIter::FileName(i, filePrefixPrePost[lev]);
516  FileSystem::Remove(FullFileName);
517  }
518  }
519  }
520  }
521  }
522 
524  HdrFile.flush();
525  HdrFile.close();
526  if( ! HdrFile.good()) {
527  amrex::Abort("ParticleContainer::CheckpointPost(): problem writing HdrFile");
528  }
529  }
530 }
531 
532 template <typename ParticleType, int NArrayReal, int NArrayInt,
533  template<class> class Allocator, class CellAssignor>
534 void
537 {
538  CheckpointPre();
539 }
540 
541 
542 template <typename ParticleType, int NArrayReal, int NArrayInt,
543  template<class> class Allocator, class CellAssignor>
544 void
547 {
548  CheckpointPost();
549 }
550 
551 
552 template <typename ParticleType, int NArrayReal, int NArrayInt,
553  template<class> class Allocator, class CellAssignor>
554 void
556 ::WriteParticles (int lev, std::ofstream& ofs, int fnum,
557  Vector<int>& which, Vector<int>& count, Vector<Long>& where,
558  const Vector<int>& write_real_comp,
559  const Vector<int>& write_int_comp,
560  const Vector<std::map<std::pair<int, int>, IntVector>>& particle_io_flags,
561  bool is_checkpoint) const
562 {
563  BL_PROFILE("ParticleContainer::WriteParticles()");
564 
565  // For a each grid, the tiles it contains
566  std::map<int, Vector<int> > tile_map;
567 
568  for (const auto& kv : m_particles[lev])
569  {
570  const int grid = kv.first.first;
571  const int tile = kv.first.second;
572  tile_map[grid].push_back(tile);
573  const auto& pflags = particle_io_flags[lev].at(kv.first);
574 
575  // Only write out valid particles.
576  count[grid] += particle_detail::countFlags(pflags);
577  }
578 
579  MFInfo info;
580  info.SetAlloc(false);
581  MultiFab state(ParticleBoxArray(lev), ParticleDistributionMap(lev), 1,0,info);
582 
583  for (MFIter mfi(state); mfi.isValid(); ++mfi)
584  {
585  const int grid = mfi.index();
586 
587  which[grid] = fnum;
588  where[grid] = VisMF::FileOffset(ofs);
589 
590  if (count[grid] == 0) { continue; }
591 
592  Vector<int> istuff;
593  Vector<ParticleReal> rstuff;
594  particle_detail::packIOData(istuff, rstuff, *this, lev, grid,
595  write_real_comp, write_int_comp,
596  particle_io_flags, tile_map[grid], count[grid], is_checkpoint);
597 
598  writeIntData(istuff.dataPtr(), istuff.size(), ofs);
599  ofs.flush(); // Some systems require this flush() (probably due to a bug)
600 
601  WriteParticleRealData(rstuff.dataPtr(), rstuff.size(), ofs);
602  ofs.flush(); // Some systems require this flush() (probably due to a bug)
603  }
604 }
605 
606 
607 template <typename ParticleType, int NArrayReal, int NArrayInt,
608  template<class> class Allocator, class CellAssignor>
609 void
611 ::Restart (const std::string& dir, const std::string& file, bool /*is_checkpoint*/)
612 {
613  Restart(dir, file);
614 }
615 
616 template <typename ParticleType, int NArrayReal, int NArrayInt,
617  template<class> class Allocator, class CellAssignor>
618 void
619 ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt, Allocator, CellAssignor>
620 ::Restart (const std::string& dir, const std::string& file)
621 {
622  BL_PROFILE("ParticleContainer::Restart()");
623  AMREX_ASSERT(!dir.empty());
624  AMREX_ASSERT(!file.empty());
625 
626  const auto strttime = amrex::second();
627 
628  int DATA_Digits_Read(5);
629  ParmParse pp("particles");
630  pp.query("datadigits_read",DATA_Digits_Read);
631 
632  std::string fullname = dir;
633  if (!fullname.empty() && fullname[fullname.size()-1] != '/') {
634  fullname += '/';
635  }
636  fullname += file;
637  std::string HdrFileName = fullname;
638  if (!HdrFileName.empty() && HdrFileName[HdrFileName.size()-1] != '/') {
639  HdrFileName += '/';
640  }
641  HdrFileName += "Header";
642 
643  Vector<char> fileCharPtr;
644  ParallelDescriptor::ReadAndBcastFile(HdrFileName, fileCharPtr);
645  std::string fileCharPtrString(fileCharPtr.dataPtr());
646  std::istringstream HdrFile(fileCharPtrString, std::istringstream::in);
647 
648  std::string version;
649  HdrFile >> version;
650  AMREX_ASSERT(!version.empty());
651 
652  // What do our version strings mean?
653  // "Version_One_Dot_Zero" -- hard-wired to write out in double precision.
654  // "Version_One_Dot_One" -- can write out either as either single or double precision.
655  // Appended to the latter version string are either "_single" or "_double" to
656  // indicate how the particles were written.
657  // "Version_Two_Dot_Zero" -- this is the AMReX particle file format
658  // "Version_Two_Dot_One" -- expanded particle ids to allow for 2**39-1 per proc
659  std::string how;
660  bool convert_ids = false;
661  if (version.find("Version_Two_Dot_One") != std::string::npos) {
662  convert_ids = true;
663  }
664  if (version.find("Version_One_Dot_Zero") != std::string::npos) {
665  how = "double";
666  }
667  else if (version.find("Version_One_Dot_One") != std::string::npos ||
668  version.find("Version_Two_Dot_Zero") != std::string::npos ||
669  version.find("Version_Two_Dot_One") != std::string::npos) {
670  if (version.find("_single") != std::string::npos) {
671  how = "single";
672  }
673  else if (version.find("_double") != std::string::npos) {
674  how = "double";
675  }
676  else {
677  std::string msg("ParticleContainer::Restart(): bad version string: ");
678  msg += version;
679  amrex::Error(version.c_str());
680  }
681  }
682  else {
683  std::string msg("ParticleContainer::Restart(): unknown version string: ");
684  msg += version;
685  amrex::Abort(msg.c_str());
686  }
687 
688  int dm;
689  HdrFile >> dm;
690  if (dm != AMREX_SPACEDIM) {
691  amrex::Abort("ParticleContainer::Restart(): dm != AMREX_SPACEDIM");
692  }
693 
694  int nr;
695  HdrFile >> nr;
696  int nrc = ParticleType::is_soa_particle ? NStructReal + NumRealComps() - AMREX_SPACEDIM : NStructReal + NumRealComps();
697  if (nr != nrc) {
698  amrex::Abort("ParticleContainer::Restart(): nr not the expected value");
699  }
700 
701  std::string comp_name;
702  for (int i = 0; i < nr; ++i) {
703  HdrFile >> comp_name;
704  }
705 
706  int ni;
707  HdrFile >> ni;
708  if (ni != NStructInt + NumIntComps()) {
709  amrex::Abort("ParticleContainer::Restart(): ni != NStructInt");
710  }
711 
712  for (int i = 0; i < ni; ++i) {
713  HdrFile >> comp_name;
714  }
715 
716  bool checkpoint;
717  HdrFile >> checkpoint;
718 
719  Long nparticles;
720  HdrFile >> nparticles;
721  AMREX_ASSERT(nparticles >= 0);
722 
723  Long maxnextid;
724  HdrFile >> maxnextid;
725  AMREX_ASSERT(maxnextid > 0);
726  ParticleType::NextID(maxnextid);
727 
728  int finest_level_in_file;
729  HdrFile >> finest_level_in_file;
730  AMREX_ASSERT(finest_level_in_file >= 0);
731 
732  // Determine whether this is a dual-grid restart or not.
733  Vector<DistributionMapping> old_dms(finest_level_in_file + 1);
734  Vector<BoxArray> old_bas(finest_level_in_file + 1);
735  Vector<BoxArray> particle_box_arrays(finest_level_in_file + 1);
736  bool dual_grid = false;
737 
738  bool have_pheaders = false;
739  for (int lev = 0; lev <= finest_level_in_file; lev++)
740  {
741  std::string phdr_name = fullname;
742  phdr_name += "/Level_";
743  phdr_name = amrex::Concatenate(phdr_name, lev, 1);
744  phdr_name += "/Particle_H";
745 
746  if (amrex::FileExists(phdr_name)) {
747  have_pheaders = true;
748  break;
749  }
750  }
751 
752  if (have_pheaders)
753  {
754  for (int lev = 0; lev <= finestLevel(); lev++)
755  {
756  old_dms[lev] = ParticleDistributionMap(lev);
757  old_bas[lev] = ParticleBoxArray(lev);
758  std::string phdr_name = fullname;
759  phdr_name += "/Level_";
760  phdr_name = amrex::Concatenate(phdr_name, lev, 1);
761  phdr_name += "/Particle_H";
762 
763  if (! amrex::FileExists(phdr_name)) { continue; }
764 
765  Vector<char> phdr_chars;
766  ParallelDescriptor::ReadAndBcastFile(phdr_name, phdr_chars);
767  std::string phdr_string(phdr_chars.dataPtr());
768  std::istringstream phdr_file(phdr_string, std::istringstream::in);
769 
770  particle_box_arrays[lev].readFrom(phdr_file);
771  if (! particle_box_arrays[lev].CellEqual(ParticleBoxArray(lev))) { dual_grid = true; }
772  }
773  } else // if no particle box array information exists in the file, we assume a single grid restart
774  {
775  dual_grid = false;
776  }
777 
778  if (dual_grid) {
779  for (int lev = 0; lev <= finestLevel(); lev++) {
780  // this can happen if there are no particles at a given level in the checkpoint
781  if (particle_box_arrays[lev].empty()) {
782  particle_box_arrays[lev] = BoxArray(Geom(lev).Domain());
783  }
784  SetParticleBoxArray(lev, particle_box_arrays[lev]);
785  DistributionMapping pdm(particle_box_arrays[lev]);
786  SetParticleDistributionMap(lev, pdm);
787  }
788  }
789 
790  Vector<int> ngrids(finest_level_in_file+1);
791  for (int lev = 0; lev <= finest_level_in_file; lev++) {
792  HdrFile >> ngrids[lev];
793  AMREX_ASSERT(ngrids[lev] > 0);
794  }
795 
796  resizeData();
797 
798  if (finest_level_in_file > finestLevel()) {
799  m_particles.resize(finest_level_in_file+1);
800  }
801 
802  for (int lev = 0; lev <= finest_level_in_file; lev++) {
803  Vector<int> which(ngrids[lev]);
804  Vector<int> count(ngrids[lev]);
805  Vector<Long> where(ngrids[lev]);
806  for (int i = 0; i < ngrids[lev]; i++) {
807  HdrFile >> which[i] >> count[i] >> where[i];
808  }
809 
810  Vector<int> grids_to_read;
811  if (lev <= finestLevel()) {
812  for (MFIter mfi(*m_dummy_mf[lev]); mfi.isValid(); ++mfi) {
813  grids_to_read.push_back(mfi.index());
814  }
815  } else {
816 
817  // we lost a level on restart. we still need to read in particles
818  // on finer levels, and put them in the right place via Redistribute()
819 
820  const int rank = ParallelDescriptor::MyProc();
821  const int NReaders = MaxReaders();
822  if (rank >= NReaders) { return; }
823 
824  const int Navg = ngrids[lev] / NReaders;
825  const int Nleft = ngrids[lev] - Navg * NReaders;
826 
827  int lo, hi;
828  if (rank < Nleft) {
829  lo = rank*(Navg + 1);
830  hi = lo + Navg + 1;
831  }
832  else {
833  lo = rank * Navg + Nleft;
834  hi = lo + Navg;
835  }
836 
837  for (int i = lo; i < hi; ++i) {
838  grids_to_read.push_back(i);
839  }
840  }
841 
842  for(int grid : grids_to_read) {
843  if (count[grid] <= 0) { continue; }
844 
845  // The file names in the header file are relative.
846  std::string name = fullname;
847 
848  if (!name.empty() && name[name.size()-1] != '/') {
849  name += '/';
850  }
851 
852  name += "Level_";
853  name += amrex::Concatenate("", lev, 1);
854  name += '/';
855  name += DataPrefix();
856  name += amrex::Concatenate("", which[grid], DATA_Digits_Read);
857 
858  std::ifstream ParticleFile;
859 
860  ParticleFile.open(name.c_str(), std::ios::in | std::ios::binary);
861 
862  if (!ParticleFile.good()) {
863  amrex::FileOpenFailed(name);
864  }
865 
866  ParticleFile.seekg(where[grid], std::ios::beg);
867 
868  // Use if constexpr to avoid instantiating the mis-matched
869  // type case and triggering the static_assert on the
870  // underlying copy calls
871  if (how == "single") {
872  if constexpr (std::is_same_v<ParticleReal, float>) {
873  ReadParticles<float>(count[grid], grid, lev, ParticleFile, finest_level_in_file, convert_ids);
874  } else {
875  amrex::Error("File contains single-precision data, while AMReX is compiled with ParticleReal==double");
876  }
877  }
878  else if (how == "double") {
879  if constexpr (std::is_same_v<ParticleReal, double>) {
880  ReadParticles<double>(count[grid], grid, lev, ParticleFile, finest_level_in_file, convert_ids);
881  } else {
882  amrex::Error("File contains double-precision data, while AMReX is compiled with ParticleReal==float");
883  }
884  }
885  else {
886  std::string msg("ParticleContainer::Restart(): bad parameter: ");
887  msg += how;
888  amrex::Error(msg.c_str());
889  }
890 
891  ParticleFile.close();
892 
893  if (!ParticleFile.good()) {
894  amrex::Abort("ParticleContainer::Restart(): problem reading particles");
895  }
896  }
897  }
898 
899  if (dual_grid) {
900  for (int lev = 0; lev <= finestLevel(); lev++) {
901  SetParticleBoxArray(lev, old_bas[lev]);
902  SetParticleDistributionMap(lev, old_dms[lev]);
903  }
904  }
905 
906  Redistribute();
907 
908  AMREX_ASSERT(OK());
909 
910  if (m_verbose > 1) {
911  auto stoptime = amrex::second() - strttime;
913  amrex::Print() << "ParticleContainer::Restart() time: " << stoptime << '\n';
914  }
915 }
916 
917 // Read a batch of particles from the checkpoint file
918 template <typename ParticleType, int NArrayReal, int NArrayInt,
919  template<class> class Allocator, class CellAssignor>
920 template <class RTYPE>
921 void
922 ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt, Allocator, CellAssignor>
923 ::ReadParticles (int cnt, int grd, int lev, std::ifstream& ifs,
924  int finest_level_in_file, bool convert_ids)
925 {
926  BL_PROFILE("ParticleContainer::ReadParticles()");
927  AMREX_ASSERT(cnt > 0);
928  AMREX_ASSERT(lev < int(m_particles.size()));
929 
930  // First read in the integer data in binary. We do not store
931  // the m_lev and m_grid data on disk. We can easily recreate
932  // that given the structure of the checkpoint file.
933  const int iChunkSize = 2 + NStructInt + NumIntComps();
934  Vector<int> istuff(std::size_t(cnt)*iChunkSize);
935  readIntData(istuff.dataPtr(), istuff.size(), ifs, FPC::NativeIntDescriptor());
936 
937  // Then the real data in binary.
938  const int rChunkSize = ParticleType::is_soa_particle ? NStructReal + NumRealComps() : AMREX_SPACEDIM + NStructReal + NumRealComps();
939  Vector<RTYPE> rstuff(std::size_t(cnt)*rChunkSize);
940  ReadParticleRealData(rstuff.dataPtr(), rstuff.size(), ifs);
941 
942  // Now reassemble the particles.
943  int* iptr = istuff.dataPtr();
944  RTYPE* rptr = rstuff.dataPtr();
945 
947  ParticleLocData pld;
948 
950  host_particles.reserve(15);
951  host_particles.resize(finest_level_in_file+1);
952 
954  std::vector<Gpu::HostVector<RTYPE> > > > host_real_attribs;
955  host_real_attribs.reserve(15);
956  host_real_attribs.resize(finest_level_in_file+1);
957 
959  std::vector<Gpu::HostVector<int> > > > host_int_attribs;
960  host_int_attribs.reserve(15);
961  host_int_attribs.resize(finest_level_in_file+1);
962 
964  host_idcpu.reserve(15);
965  host_idcpu.resize(finestLevel()+1);
966 
967  for (int i = 0; i < cnt; i++) {
968  // note: for pure SoA particle layouts, we do write the id, cpu and positions as a struct
969  // for backwards compatibility with readers
970  if (convert_ids) {
971  std::int32_t xi, yi;
972  std::uint32_t xu, yu;
973  xi = iptr[0];
974  yi = iptr[1];
975  std::memcpy(&xu, &xi, sizeof(xi));
976  std::memcpy(&yu, &yi, sizeof(yi));
977  ptemp.m_idcpu = ((std::uint64_t)xu) << 32 | yu;
978  } else {
979  ptemp.id() = iptr[0];
980  ptemp.cpu() = iptr[1];
981  }
982  iptr += 2;
983 
984  for (int j = 0; j < NStructInt; j++)
985  {
986  ptemp.idata(j) = *iptr;
987  ++iptr;
988  }
989 
990  AMREX_ASSERT(ptemp.id() > 0);
991 
992  AMREX_D_TERM(ptemp.pos(0) = ParticleReal(rptr[0]);,
993  ptemp.pos(1) = ParticleReal(rptr[1]);,
994  ptemp.pos(2) = ParticleReal(rptr[2]););
995 
996  rptr += AMREX_SPACEDIM;
997 
998  for (int j = 0; j < NStructReal; j++)
999  {
1000  ptemp.rdata(j) = ParticleReal(*rptr);
1001  ++rptr;
1002  }
1003 
1004  locateParticle(ptemp, pld, 0, finestLevel(), 0);
1005 
1006  std::pair<int, int> ind(grd, pld.m_tile);
1007 
1008  host_real_attribs[lev][ind].resize(NumRealComps());
1009  host_int_attribs[lev][ind].resize(NumIntComps());
1010 
1011  // add the struct
1012  if constexpr(!ParticleType::is_soa_particle)
1013  {
1014  host_particles[lev][ind].push_back(ptemp);
1015 
1016  // add the real...
1017  for (int icomp = 0; icomp < NumRealComps(); icomp++) {
1018  host_real_attribs[lev][ind][icomp].push_back(*rptr);
1019  ++rptr;
1020  }
1021 
1022  // ... and int array data
1023  for (int icomp = 0; icomp < NumIntComps(); icomp++) {
1024  host_int_attribs[lev][ind][icomp].push_back(*iptr);
1025  ++iptr;
1026  }
1027  } else {
1028  host_particles[lev][ind];
1029 
1030  for (int j = 0; j < AMREX_SPACEDIM; j++) {
1031  host_real_attribs[pld.m_lev][ind][j].push_back(ptemp.pos(j));
1032  }
1033 
1034  host_idcpu[pld.m_lev][ind].push_back(ptemp.m_idcpu);
1035 
1036  // read all other SoA
1037  // add the real...
1038  for (int icomp = AMREX_SPACEDIM; icomp < NumRealComps(); icomp++) {
1039  host_real_attribs[lev][ind][icomp].push_back(*rptr);
1040  ++rptr;
1041  }
1042 
1043  // ... and int array data
1044  for (int icomp = 0; icomp < NumIntComps(); icomp++) {
1045  host_int_attribs[lev][ind][icomp].push_back(*iptr);
1046  ++iptr;
1047  }
1048  }
1049  }
1050 
1051  for (int host_lev = 0; host_lev < static_cast<int>(host_particles.size()); ++host_lev)
1052  {
1053  for (auto& kv : host_particles[host_lev]) {
1054  auto grid = kv.first.first;
1055  auto tile = kv.first.second;
1056  const auto& src_tile = kv.second;
1057 
1058  auto& dst_tile = DefineAndReturnParticleTile(host_lev, grid, tile);
1059  auto old_size = dst_tile.size();
1060  auto new_size = old_size;
1061  if constexpr(!ParticleType::is_soa_particle)
1062  {
1063  new_size += src_tile.size();
1064  } else {
1065  amrex::ignore_unused(src_tile);
1066  new_size += host_real_attribs[host_lev][std::make_pair(grid,tile)][0].size();
1067  }
1068  dst_tile.resize(new_size);
1069 
1070  if constexpr(!ParticleType::is_soa_particle)
1071  {
1072  Gpu::copyAsync(Gpu::hostToDevice, src_tile.begin(), src_tile.end(),
1073  dst_tile.GetArrayOfStructs().begin() + old_size);
1074  } else {
1076  host_idcpu[host_lev][std::make_pair(grid,tile)].begin(),
1077  host_idcpu[host_lev][std::make_pair(grid,tile)].end(),
1078  dst_tile.GetStructOfArrays().GetIdCPUData().begin() + old_size);
1079  }
1080 
1081  for (int i = 0; i < NumRealComps(); ++i) { // NOLINT(readability-misleading-indentation)
1083  host_real_attribs[host_lev][std::make_pair(grid,tile)][i].begin(),
1084  host_real_attribs[host_lev][std::make_pair(grid,tile)][i].end(),
1085  dst_tile.GetStructOfArrays().GetRealData(i).begin() + old_size);
1086  }
1087 
1088  for (int i = 0; i < NumIntComps(); ++i) {
1090  host_int_attribs[host_lev][std::make_pair(grid,tile)][i].begin(),
1091  host_int_attribs[host_lev][std::make_pair(grid,tile)][i].end(),
1092  dst_tile.GetStructOfArrays().GetIntData(i).begin() + old_size);
1093  }
1094  }
1095  }
1096 
1098 }
1099 
1100 template <typename ParticleType, int NArrayReal, int NArrayInt,
1101  template<class> class Allocator, class CellAssignor>
1102 void
1103 ParticleContainer_impl<ParticleType, NArrayReal, NArrayInt, Allocator, CellAssignor>
1104 ::WriteAsciiFile (const std::string& filename)
1105 {
1106  BL_PROFILE("ParticleContainer::WriteAsciiFile()");
1107  AMREX_ASSERT(!filename.empty());
1108 
1109  const auto strttime = amrex::second();
1110  //
1111  // Count # of valid particles.
1112  //
1113  Long nparticles = 0;
1114 
1115  for (int lev = 0; lev < m_particles.size(); lev++) {
1116  auto& pmap = m_particles[lev];
1117  for (const auto& kv : pmap) {
1118  const auto& aos = kv.second.GetArrayOfStructs();
1119  auto np = aos.numParticles();
1120  Gpu::HostVector<ParticleType> host_aos(np);
1121  Gpu::copy(Gpu::deviceToHost, aos.begin(), aos.begin() + np, host_aos.begin());
1122  for (int k = 0; k < np; ++k) {
1123  const ParticleType& p = host_aos[k];
1124  if (p.id() > 0) {
1125  //
1126  // Only count (and checkpoint) valid particles.
1127  //
1128  nparticles++;
1129  }
1130  }
1131  }
1132  }
1133 
1134  //
1135  // And send count to I/O processor.
1136  //
1138 
1140  {
1141  //
1142  // Have I/O processor open file and write out particle metadata.
1143  //
1144  std::ofstream File;
1145 
1146  File.open(filename.c_str(), std::ios::out|std::ios::trunc);
1147 
1148  if (!File.good()) {
1149  amrex::FileOpenFailed(filename);
1150  }
1151 
1152  File << nparticles << '\n';
1153  File << NStructReal << '\n';
1154  File << NStructInt << '\n';
1155  File << NumRealComps() << '\n';
1156  File << NumIntComps() << '\n';
1157 
1158  File.flush();
1159 
1160  File.close();
1161 
1162  if (!File.good()) {
1163  amrex::Abort("ParticleContainer::WriteAsciiFile(): problem writing file");
1164  }
1165  }
1166 
1168 
1169  const int MyProc = ParallelDescriptor::MyProc();
1170 
1171  for (int proc = 0; proc < ParallelDescriptor::NProcs(); proc++)
1172  {
1173  if (MyProc == proc)
1174  {
1175  //
1176  // Each CPU opens the file for appending and adds its particles.
1177  //
1178  VisMF::IO_Buffer io_buffer(VisMF::IO_Buffer_Size);
1179 
1180  std::ofstream File;
1181 
1182  File.rdbuf()->pubsetbuf(io_buffer.dataPtr(), io_buffer.size());
1183 
1184  File.open(filename.c_str(), std::ios::out|std::ios::app);
1185 
1186  File.precision(15);
1187 
1188  if (!File.good()) {
1189  amrex::FileOpenFailed(filename);
1190  }
1191 
1192  for (int lev = 0; lev < m_particles.size(); lev++) {
1193  auto& pmap = m_particles[lev];
1194  for (const auto& kv : pmap) {
1195  ParticleTile<ParticleType, NArrayReal, NArrayInt,
1196  amrex::PinnedArenaAllocator> pinned_ptile;
1197  pinned_ptile.define(NumRuntimeRealComps(), NumRuntimeIntComps());
1198  pinned_ptile.resize(kv.second.numParticles());
1199  amrex::copyParticles(pinned_ptile, kv.second);
1200  const auto& host_aos = pinned_ptile.GetArrayOfStructs();
1201  const auto& host_soa = pinned_ptile.GetStructOfArrays();
1202 
1203  auto np = host_aos.numParticles();
1204  for (int index = 0; index < np; ++index) {
1205  const ParticleType* it = &host_aos[index];
1206  if (it->id() > 0) {
1207 
1208  // write out the particle struct first...
1209  AMREX_D_TERM(File << it->pos(0) << ' ',
1210  << it->pos(1) << ' ',
1211  << it->pos(2) << ' ');
1212 
1213  for (int i = 0; i < NStructReal; i++) {
1214  File << it->rdata(i) << ' ';
1215  }
1216 
1217  File << it->id() << ' ';
1218  File << it->cpu() << ' ';
1219 
1220  for (int i = 0; i < NStructInt; i++) {
1221  File << it->idata(i) << ' ';
1222  }
1223 
1224  // then the particle attributes.
1225  for (int i = 0; i < NumRealComps(); i++) {
1226  File << host_soa.GetRealData(i)[index] << ' ';
1227  }
1228 
1229  for (int i = 0; i < NumIntComps(); i++) {
1230  File << host_soa.GetIntData(i)[index] << ' ';
1231  }
1232 
1233  File << '\n';
1234  }
1235  }
1236  }
1237  }
1238 
1239  File.flush();
1240 
1241  File.close();
1242 
1243  if (!File.good()) {
1244  amrex::Abort("ParticleContainer::WriteAsciiFile(): problem writing file");
1245  }
1246 
1247  }
1248 
1250  }
1251 
1252  if (m_verbose > 1)
1253  {
1254  auto stoptime = amrex::second() - strttime;
1255 
1257 
1258  amrex::Print() << "ParticleContainer::WriteAsciiFile() time: " << stoptime << '\n';
1259  }
1260 }
1261 
1262 #endif /*AMREX_PARTICLEIO_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_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
amrex::ParmParse pp
Input file parser instance for the given namespace.
Definition: AMReX_HypreIJIface.cpp:15
#define AMREX_D_TERM(a, b, c)
Definition: AMReX_SPACE.H:129
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
Definition: AMReX_MFIter.H:57
bool isValid() const noexcept
Is the iterator valid i.e. is it associated with a FAB?
Definition: AMReX_MFIter.H:141
A collection (stored as an array) of FArrayBox objects.
Definition: AMReX_MultiFab.H:38
Definition: AMReX_PODVector.H:246
iterator begin() noexcept
Definition: AMReX_PODVector.H:601
Parse Parameters From Command Line and Input Files.
Definition: AMReX_ParmParse.H:320
int query(const char *name, bool &ref, int ival=FIRST) const
Same as querykth() but searches for the last occurrence of name.
Definition: AMReX_ParmParse.cpp:1309
A distributed container for Particles sorted onto the levels, grids, and tiles of a block-structured ...
Definition: AMReX_ParticleContainer.H:145
typename SoA::IntVector IntVector
Definition: AMReX_ParticleContainer.H:191
T_ParticleType ParticleType
Definition: AMReX_ParticleContainer.H:147
Definition: AMReX_GpuAllocators.H:126
This class provides the user with a few print options.
Definition: AMReX_Print.H:35
Long size() const noexcept
Definition: AMReX_Vector.H:50
T * dataPtr() noexcept
get access to the underlying data pointer
Definition: AMReX_Vector.H:46
bool UseAsyncOut()
Definition: AMReX_AsyncOut.cpp:70
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 copy(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:121
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 ReadAndBcastFile(const std::string &filename, Vector< char > &charBuf, bool bExitOnError, const MPI_Comm &comm)
Definition: AMReX_ParallelDescriptor.cpp:1459
void ReduceLongMax(Long &)
Long max reduction.
Definition: AMReX_ParallelDescriptor.cpp:1224
void ReduceLongSum(Long &)
Long sum reduction.
Definition: AMReX_ParallelDescriptor.cpp:1223
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
void ReduceRealMax(Vector< std::reference_wrapper< Real > > const &)
Definition: AMReX_ParallelDescriptor.cpp:1215
static int f(amrex::Real t, N_Vector y_data, N_Vector y_rhs, void *user_data)
Definition: AMReX_SundialsIntegrator.H:44
static constexpr int P
Definition: AMReX_OpenBC.H:14
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
void copyParticles(DstTile &dst, const SrcTile &src) noexcept
Copy particles from src to dst. This version copies all the particles, writing them to the beginning ...
Definition: AMReX_ParticleTransformation.H:158
void readFloatData(float *data, std::size_t size, std::istream &is, const RealDescriptor &rd)
Definition: AMReX_VectorIO.cpp:120
void WritePlotFile(const Vector< MultiFab * > &mfa, const Vector< Box > &probDomain, AmrData &amrdToMimic, const std::string &oFile, bool verbose, const Vector< std::string > &varNames)
Definition: AMReX_WritePlotFile.cpp:292
void writeFloatData(const float *data, std::size_t size, std::ostream &os, const RealDescriptor &rd=FPC::Native32RealDescriptor())
Definition: AMReX_VectorIO.cpp:114
bool FileExists(const std::string &filename)
Check if a file already exists. Return true if the filename is an existing file, directory,...
Definition: AMReX_Utility.cpp:139
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
double second() noexcept
Definition: AMReX_Utility.cpp:922
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
void Error(const std::string &msg)
Print out message to cerr and exit via amrex::Abort().
Definition: AMReX.cpp:219
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 begin(BoxND< dim > const &box) noexcept
Definition: AMReX_Box.H:1881
void readDoubleData(double *data, std::size_t size, std::istream &is, const RealDescriptor &rd)
Definition: AMReX_VectorIO.cpp:132
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
void readIntData(To *data, std::size_t size, std::istream &is, const amrex::IntDescriptor &id)
Definition: AMReX_IntConv.H:36
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
Definition: AMReX_ParticleIO.H:35
AMREX_GPU_HOST_DEVICE int operator()(const P &p) const
Definition: AMReX_ParticleIO.H:38
FabArray memory allocation information.
Definition: AMReX_FabArray.H:66
MFInfo & SetAlloc(bool a) noexcept
Definition: AMReX_FabArray.H:73
uint64_t m_idcpu
Definition: AMReX_Particle.H:252
A struct used for storing a particle's position in the AMR hierarchy.
Definition: AMReX_ParticleContainer.H:91
int m_tile
Definition: AMReX_ParticleContainer.H:94
int m_lev
Definition: AMReX_ParticleContainer.H:92
Definition: AMReX_ParticleTile.H:693
The struct used to store particles.
Definition: AMReX_Particle.H:295
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ParticleIDWrapper id() &
Definition: AMReX_Particle.H:315
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE RealVect pos() const &
Definition: AMReX_Particle.H:338
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE RealType & rdata(int index) &
Definition: AMReX_Particle.H:356
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE int & idata(int index) &
Definition: AMReX_Particle.H:427
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ParticleCPUWrapper cpu() &
Definition: AMReX_Particle.H:312