Block-Structured AMR Software Framework
AMReX_WriteBinaryParticleDataHDF5.H
Go to the documentation of this file.
1 #ifndef AMREX_WRITE_BINARY_PARTICLE_DATA_HDF5_H
2 #define AMREX_WRITE_BINARY_PARTICLE_DATA_HDF5_H
3 
4 static int CreateWriteHDF5Attr(hid_t loc, const char *name, hsize_t n, void *data, hid_t dtype)
5 {
6  herr_t ret;
7  hid_t attr, attr_space;
8  hsize_t dims = n;
9 
10  attr_space = H5Screate_simple(1, &dims, NULL);
11 
12  attr = H5Acreate(loc, name, dtype, attr_space, H5P_DEFAULT, H5P_DEFAULT);
13  if (attr < 0) {
14  amrex::Print() << " Error with H5Acreate [" << name << "]\n";
15  return -1;
16  }
17 
18  ret = H5Awrite(attr, dtype, (void*)data);
19  if (ret < 0) {
20  amrex::Print() << " Error with H5Awrite [" << name << "]\n";
21  return -1;
22  }
23  H5Sclose(attr_space);
24  H5Aclose(attr);
25  return 1;
26 }
27 
28 static int CreateWriteHDF5AttrString(hid_t loc, const char *name, const char* str)
29 {
30  hid_t attr, atype, space;
31  herr_t ret;
32 
33  AMREX_ASSERT(name);
34  AMREX_ASSERT(str);
35 
36  space = H5Screate(H5S_SCALAR);
37  atype = H5Tcopy(H5T_C_S1);
38  H5Tset_size(atype, strlen(str)+1);
39  H5Tset_strpad(atype,H5T_STR_NULLTERM);
40  attr = H5Acreate(loc, name, atype, space, H5P_DEFAULT, H5P_DEFAULT);
41  if (attr < 0) {
42  amrex::Print() << " Error with H5Acreate [" << name << "]\n";
43  return -1;
44  }
45 
46  ret = H5Awrite(attr, atype, str);
47  if (ret < 0) {
48  amrex::Print() << " Error with H5Awrite [" << name << "]\n";
49  return -1;
50  }
51 
52  H5Tclose(atype);
53  H5Sclose(space);
54  H5Aclose(attr);
55 
56  return 1;
57 }
58 
59 static int ReadHDF5Attr(hid_t loc, const char *name, void *data, hid_t dtype)
60 {
61  herr_t ret;
62  hid_t attr;
63 
64  attr = H5Aopen(loc, name, H5P_DEFAULT);
65  if (attr < 0) {
66  amrex::Print() << " Error with H5Aopen [" << name << "]\n";
67  return -1;
68  }
69 
70  ret = H5Aread(attr, dtype, data);
71  if (ret < 0) {
72  amrex::Print() << " Error with H5Aread [" << name << "]\n";
73  return -1;
74  }
75  H5Aclose(attr);
76  return 1;
77 }
78 
79 #ifdef AMREX_USE_MPI
80 static void SetHDF5fapl(hid_t fapl, MPI_Comm comm)
81 #else
82 static void SetHDF5fapl(hid_t fapl)
83 #endif
84 {
85 #ifdef AMREX_USE_MPI
86  H5Pset_fapl_mpio(fapl, comm, MPI_INFO_NULL);
87 
88  // Alignment and metadata block size
89  int alignment = 16 * 1024 * 1024;
90  int blocksize = 4 * 1024 * 1024;
91  H5Pset_alignment(fapl, alignment, alignment);
92  H5Pset_meta_block_size(fapl, blocksize);
93 
94  // Collective metadata ops
95  H5Pset_coll_metadata_write(fapl, true);
96  H5Pset_all_coll_metadata_ops(fapl, true);
97 
98  // Defer cache flush
99  H5AC_cache_config_t cache_config;
100  cache_config.version = H5AC__CURR_CACHE_CONFIG_VERSION;
101  H5Pget_mdc_config(fapl, &cache_config);
102  cache_config.set_initial_size = 1;
103  cache_config.initial_size = 16 * 1024 * 1024;
104  cache_config.evictions_enabled = 0;
105  cache_config.incr_mode = H5C_incr__off;
106  cache_config.flash_incr_mode = H5C_flash_incr__off;
107  cache_config.decr_mode = H5C_decr__off;
108  H5Pset_mdc_config (fapl, &cache_config);
109 #else
110  H5Pset_fapl_sec2(fapl);
111 #endif
112 
113 }
114 template <class PC, class F, std::enable_if_t<IsParticleContainer<PC>::value, int> foo = 0>
115 void WriteHDF5ParticleDataSync (PC const& pc,
116  const std::string& dir, const std::string& name,
117  const Vector<int>& write_real_comp,
118  const Vector<int>& write_int_comp,
119  const Vector<std::string>& real_comp_names,
120  const Vector<std::string>& int_comp_names,
121  const std::string& compression,
122  F&& f, bool is_checkpoint)
123 {
124  BL_PROFILE("WriteHDF5ParticleDataSync()");
125  AMREX_ASSERT(pc.OK());
126 
127  AMREX_ASSERT(sizeof(typename PC::ParticleType::RealType) == 4 ||
128  sizeof(typename PC::ParticleType::RealType) == 8);
129 
130  constexpr int NStructReal = PC::NStructReal;
131  constexpr int NStructInt = PC::NStructInt;
132 
133  const int NProcs = ParallelDescriptor::NProcs();
134  const int IOProcNumber = ParallelDescriptor::IOProcessorNumber();
135 
136  AMREX_ALWAYS_ASSERT(real_comp_names.size() == pc.NumRealComps() + NStructReal);
137  AMREX_ALWAYS_ASSERT( int_comp_names.size() == pc.NumIntComps() + NStructInt);
138 
139 #ifdef AMREX_USE_HDF5_ASYNC
140  // For HDF5 async VOL, block and wait previous tasks have all completed
141  if (es_par_g != 0)
142  async_vol_es_wait_particle();
143  else {
144  ExecOnFinalize(async_vol_es_wait_close_particle);
145  es_par_g = H5EScreate();
146  }
147 #endif
148 
149  std::string pdir = dir;
150  if ( ! pdir.empty() && pdir[pdir.size()-1] != '/') { pdir += '/'; }
151  pdir += name;
152 
153  if ( ! pc.GetLevelDirectoriesCreated()) {
155  {
156  if ( ! amrex::UtilCreateDirectory(pdir, 0755))
157  {
159  }
160  }
162  }
163 
164  /* std::ofstream HdrFile; */
165  hid_t fid, grp, fapl, comp_dtype;
166  int status;
167 
168  Long nparticles = 0;
169  Long maxnextid;
170 
171  comp_dtype = H5Tcreate (H5T_COMPOUND, 2 * AMREX_SPACEDIM * sizeof(int));
172  if (1 == AMREX_SPACEDIM) {
173  H5Tinsert (comp_dtype, "lo_i", 0 * sizeof(int), H5T_NATIVE_INT);
174  H5Tinsert (comp_dtype, "hi_i", 1 * sizeof(int), H5T_NATIVE_INT);
175  }
176  else if (2 == AMREX_SPACEDIM) {
177  H5Tinsert (comp_dtype, "lo_i", 0 * sizeof(int), H5T_NATIVE_INT);
178  H5Tinsert (comp_dtype, "lo_j", 1 * sizeof(int), H5T_NATIVE_INT);
179  H5Tinsert (comp_dtype, "hi_i", 2 * sizeof(int), H5T_NATIVE_INT);
180  H5Tinsert (comp_dtype, "hi_j", 3 * sizeof(int), H5T_NATIVE_INT);
181  }
182  else if (3 == AMREX_SPACEDIM) {
183  H5Tinsert (comp_dtype, "lo_i", 0 * sizeof(int), H5T_NATIVE_INT);
184  H5Tinsert (comp_dtype, "lo_j", 1 * sizeof(int), H5T_NATIVE_INT);
185  H5Tinsert (comp_dtype, "lo_k", 2 * sizeof(int), H5T_NATIVE_INT);
186  H5Tinsert (comp_dtype, "hi_i", 3 * sizeof(int), H5T_NATIVE_INT);
187  H5Tinsert (comp_dtype, "hi_j", 4 * sizeof(int), H5T_NATIVE_INT);
188  H5Tinsert (comp_dtype, "hi_k", 5 * sizeof(int), H5T_NATIVE_INT);
189  }
190 
191  // evaluate f for every particle to determine which ones to output
192  Vector<std::map<std::pair<int, int>, typename PC::IntVector > >
193  particle_io_flags(pc.GetParticles().size());
194  for (int lev = 0; lev < pc.GetParticles().size(); lev++)
195  {
196  const auto& pmap = pc.GetParticles(lev);
197  for (const auto& kv : pmap)
198  {
199  auto& flags = particle_io_flags[lev][kv.first];
200  particle_detail::fillFlags(flags, kv.second, f);
201  }
202  }
203 
205 
206  if(pc.GetUsePrePost())
207  {
208  nparticles = pc.GetNParticlesPrePost();
209  maxnextid = pc.GetMaxNextIDPrePost();
210  }
211  else
212  {
213  nparticles = particle_detail::countFlags(particle_io_flags, pc);
214  maxnextid = PC::ParticleType::NextID();
215  ParallelDescriptor::ReduceLongSum(nparticles, IOProcNumber);
216  PC::ParticleType::NextID(maxnextid);
217  ParallelDescriptor::ReduceLongMax(maxnextid, IOProcNumber);
218  }
219 
220  std::string HDF5FileName = pdir;
221  if ( ! HDF5FileName.empty() && HDF5FileName[HDF5FileName.size()-1] != '/')
222  HDF5FileName += '/';
223 
224  HDF5FileName += name;
225  HDF5FileName += ".h5";
226  pc.HdrFileNamePrePost = HDF5FileName;
227 
228  BL_PROFILE_VAR("H5FileCreate", h5fg);
229 
231  {
232  int set_stripe = 0;
233  char setstripe[1024];
234  int stripe_count = 128;
235  int stripe_size = 1;
236  char *stripe_count_str = getenv("HDF5_STRIPE_COUNT");
237  char *stripe_size_str = getenv("HDF5_STRIPE_SIZE");
238  if (stripe_count_str) {
239  stripe_count = atoi(stripe_count_str);
240  set_stripe = 1;
241  }
242  if (stripe_size_str) {
243  stripe_size = atoi(stripe_size_str);
244  set_stripe = 1;
245  }
246  if (set_stripe == 1) {
247  snprintf(setstripe, sizeof setstripe, "lfs setstripe -c %d -S %dm %s", stripe_count, stripe_size, pdir.c_str());
248  std::cout << "Setting stripe parameters for HDF5 output: " << setstripe << std::endl;
249  amrex::ignore_unused(std::system(setstripe));
250  }
251 
252  fid = H5Fcreate(HDF5FileName.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
253  if (fid < 0) { amrex::FileOpenFailed(HDF5FileName.c_str()); }
254 
255  //
256  // First thing written is our Checkpoint/Restart version string.
257  // We append "_single" or "_double" to the version string indicating
258  // whether we're using "float" or "double" floating point data in the
259  // particles so that we can Restart from the checkpoint files.
260  //
261  std::string versionName = is_checkpoint ? PC::CheckpointVersion() : PC::PlotfileVersion();
262  if (sizeof(typename PC::ParticleType::RealType) == 4) {
263  versionName += "_single";
264  } else {
265  versionName += "_double";
266  }
267 
268  CreateWriteHDF5AttrString(fid, "version_name", versionName.c_str());
269 
270  int num_output_real = 0;
271  for (int i = 0; i < pc.NumRealComps() + NStructReal; ++i) {
272  if (write_real_comp[i]) { ++num_output_real; }
273  }
274 
275  int num_output_int = 0;
276  for (int i = 0; i < pc.NumIntComps() + NStructInt; ++i) {
277  if (write_int_comp[i]) { ++num_output_int; }
278  }
279 
280  // AMREX_SPACEDIM and N for sanity checking.
281  int ndim = AMREX_SPACEDIM;
282  grp = H5Gcreate(fid, "Chombo_global", H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
283  CreateWriteHDF5Attr(grp, "SpaceDim", 1, &ndim, H5T_NATIVE_INT);
284  H5Gclose(grp);
285 
286  // The number of extra real parameters
287  int ncomp = num_output_real + num_output_int;
288  CreateWriteHDF5Attr(fid, "num_component", 1, &ncomp, H5T_NATIVE_INT);
289  CreateWriteHDF5Attr(fid, "num_component_real", 1, &num_output_real, H5T_NATIVE_INT);
290 
291  char comp_name[128];
292 
293  // Real component names
294  int real_comp_count = AMREX_SPACEDIM; // position values
295 
296  for (int i = 0; i < NStructReal + pc.NumRealComps(); ++i ) {
297  if (write_real_comp[i]) {
298  /* HdrFile << real_comp_names[i] << '\n'; */
299  snprintf(comp_name, sizeof comp_name, "real_component_%d", real_comp_count);
300  CreateWriteHDF5AttrString(fid, comp_name, real_comp_names[i].c_str());
301  ++real_comp_count;
302  }
303  }
304 
305  // The number of extra int parameters
306  CreateWriteHDF5Attr(fid, "num_component_int", 1, &num_output_int, H5T_NATIVE_INT);
307 
308  // int component names
309  int int_comp_count = 2;
310 
311  for (int i = 0; i < NStructInt + pc.NumIntComps(); ++i ) {
312  if (write_int_comp[i]) {
313  snprintf(comp_name, sizeof comp_name, "int_component_%d", int_comp_count);
314  CreateWriteHDF5AttrString(fid, comp_name, int_comp_names[i].c_str());
315  ++int_comp_count;
316  }
317  }
318 
319  // The total number of particles.
320  CreateWriteHDF5Attr(fid, "nparticles", 1, &nparticles, H5T_NATIVE_LLONG);
321 
322  // The value of nextid that we need to restore on restart.
323  CreateWriteHDF5Attr(fid, "maxnextid", 1, &maxnextid, H5T_NATIVE_LLONG);
324 
325  // Then the finest level of the AMR hierarchy.
326  int finest_level = pc.finestLevel();
327  CreateWriteHDF5Attr(fid, "finest_level", 1, &finest_level, H5T_NATIVE_INT);
328 
329  char level_name[128];
330  int ngrids;
331 
332  hid_t dcpl_id = H5Pcreate(H5P_DATASET_CREATE);
333  H5Pset_fill_time(dcpl_id, H5D_FILL_TIME_NEVER);
334 
335  for (int lev = 0; lev <= finest_level; ++lev) {
336  snprintf(level_name, sizeof level_name, "level_%d", lev);
337 
338  grp = H5Gcreate(fid, level_name, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
339  if (grp < 0) {
340  std::cout << "H5Gcreate [" << level_name << "] failed!" << std::endl;
341  continue;
342  }
343 
344  // Then the number of grids at each level.
345  ngrids = pc.ParticleBoxArray(lev).size();
346  CreateWriteHDF5Attr(grp, "ngrids", 1, &ngrids, H5T_NATIVE_INT);
347 
348  /* int ratio = 1; */
349  /* CreateWriteHDF5Attr(grp, "ref_ratio", 1, &ratio); */
350 
351  int mfs_size = 2 * AMREX_SPACEDIM;
352  hsize_t mfs_dim = (hsize_t)ngrids;
353 
354  hid_t mfs_dset_space = H5Screate_simple(1, &mfs_dim, NULL);
355  hid_t mfs_dset= H5Dcreate(grp, "boxes", comp_dtype, mfs_dset_space, H5P_DEFAULT, dcpl_id, H5P_DEFAULT);
356 
357  Vector<int> vbox(ngrids * mfs_size);
358  for(int j = 0; j < pc.ParticleBoxArray(lev).size(); ++j) {
359  for(int i = 0; i < AMREX_SPACEDIM; ++i) {
360  vbox[(j * mfs_size) + i] = pc.ParticleBoxArray(lev)[j].smallEnd(i);
361  vbox[(j * mfs_size) + i + AMREX_SPACEDIM] = pc.ParticleBoxArray(lev)[j].bigEnd(i);
362  }
363  }
364 
365  status = H5Dwrite(mfs_dset, comp_dtype, H5S_ALL, H5S_ALL, H5P_DEFAULT, &(vbox[0]));
366  if (status < 0) {
367  std::string msg("ParticleContainer::WriteHDF5ParticleDataSync(): unable to write boxes dataset");
368  amrex::Abort(msg.c_str());
369  }
370 
371  H5Sclose(mfs_dset_space);
372  H5Dclose(mfs_dset);
373 
374  H5Gclose(grp);
375  }
376  H5Pclose(dcpl_id);
377  H5Fclose(fid);
378  }
379 
381  BL_PROFILE_VAR_STOP(h5fg);
382 
383  // We want to write the data out in parallel.
384  // We'll allow up to nOutFiles active writers at a time.
385  int nOutFiles(-1);
386 
387  ParmParse pp("particles");
388  pp.query("particles_nfiles",nOutFiles);
389  if(nOutFiles == -1) nOutFiles = NProcs;
390  /* nOutFiles = std::max(1, std::min(nOutFiles,NProcs)); */
391  pc.nOutFilesPrePost = nOutFiles;
392 
393  fapl = H5Pcreate (H5P_FILE_ACCESS);
394 #ifdef AMREX_USE_MPI
396 #else
397  SetHDF5fapl(fapl);
398 #endif
399 
400 #ifdef AMREX_USE_HDF5_ASYNC
401  fid = H5Fopen_async(HDF5FileName.c_str(), H5F_ACC_RDWR, fapl, es_par_g);
402 #else
403  fid = H5Fopen(HDF5FileName.c_str(), H5F_ACC_RDWR, fapl);
404 #endif
405  if (fid < 0) {
406  FileOpenFailed(HDF5FileName.c_str());
407  }
408 
409  char level_name[64];
410  for (int lev = 0; lev <= pc.finestLevel(); lev++)
411  {
412  snprintf(level_name, sizeof level_name, "level_%d", lev);
413 #ifdef AMREX_USE_HDF5_ASYNC
414  grp = H5Gopen_async(fid, level_name, H5P_DEFAULT, es_par_g);
415 #else
416  grp = H5Gopen(fid, level_name, H5P_DEFAULT);
417 #endif
418 
419  bool gotsome;
420  if(pc.usePrePost) {
421  gotsome = (pc.nParticlesAtLevelPrePost[lev] > 0);
422  } else {
423  gotsome = (pc.NumberOfParticlesAtLevel(lev) > 0);
424  }
425 
426  MFInfo info;
427  info.SetAlloc(false);
428  MultiFab state(pc.ParticleBoxArray(lev),
429  pc.ParticleDistributionMap(lev),
430  1,0,info);
431 
432  // We eventually want to write out the file name and the offset
433  // into that file into which each grid of particles is written.
434  Vector<int> which(state.size(),0);
435  Vector<int > count(state.size(),0);
436  Vector<Long> where(state.size(),0);
437 
438  if(pc.usePrePost) {
439  pc.filePrefixPrePost[lev] = HDF5FileName;
440  }
441 
442  if (gotsome)
443  {
444  pc.WriteParticlesHDF5(lev, grp, which, count, where,
445  write_real_comp, write_int_comp,
446  compression,
447  particle_io_flags, is_checkpoint);
448 
449  if(pc.usePrePost) {
450  pc.whichPrePost[lev] = which;
451  pc.countPrePost[lev] = count;
452  pc.wherePrePost[lev] = where;
453  } else {
454  ParallelDescriptor::ReduceIntSum (which.dataPtr(), which.size(), IOProcNumber);
455  ParallelDescriptor::ReduceIntSum (count.dataPtr(), count.size(), IOProcNumber);
456  ParallelDescriptor::ReduceLongSum(where.dataPtr(), where.size(), IOProcNumber);
457  }
458  }
459 
460 #ifdef AMREX_USE_HDF5_ASYNC
461  H5Gclose_async(grp, es_par_g);
462 #else
463  H5Gclose(grp);
464 #endif
465  } // end for (lev...)
466 
467  H5Tclose(comp_dtype);
468  H5Pclose(fapl);
469 
470 #ifdef AMREX_USE_HDF5_ASYNC
471  H5Fclose_async(fid, es_par_g);
472 #else
473  H5Fclose(fid);
474 #endif
475 
476  return;
477 } // End WriteHDF5ParticleDataSync
478 
479 #endif /*AMREX_WRITE_BINARY_PARTICLE_DATA_HDF5_H*/
#define BL_PROFILE(a)
Definition: AMReX_BLProfiler.H:551
#define BL_PROFILE_VAR_STOP(vname)
Definition: AMReX_BLProfiler.H:563
#define BL_PROFILE_VAR(fname, vname)
Definition: AMReX_BLProfiler.H:560
#define AMREX_ASSERT(EX)
Definition: AMReX_BLassert.H:38
#define AMREX_ALWAYS_ASSERT(EX)
Definition: AMReX_BLassert.H:50
amrex::ParmParse pp
Input file parser instance for the given namespace.
Definition: AMReX_HypreIJIface.cpp:15
static int CreateWriteHDF5Attr(hid_t loc, const char *name, hsize_t n, void *data, hid_t dtype)
Definition: AMReX_WriteBinaryParticleDataHDF5.H:4
static int ReadHDF5Attr(hid_t loc, const char *name, void *data, hid_t dtype)
Definition: AMReX_WriteBinaryParticleDataHDF5.H:59
static int CreateWriteHDF5AttrString(hid_t loc, const char *name, const char *str)
Definition: AMReX_WriteBinaryParticleDataHDF5.H:28
static void SetHDF5fapl(hid_t fapl, MPI_Comm comm)
Definition: AMReX_WriteBinaryParticleDataHDF5.H:80
void WriteHDF5ParticleDataSync(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, const std::string &compression, F &&f, bool is_checkpoint)
Definition: AMReX_WriteBinaryParticleDataHDF5.H:115
int MPI_Comm
Definition: AMReX_ccse-mpi.H:47
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
This class provides the user with a few print options.
Definition: AMReX_Print.H:35
AMREX_GPU_HOST_DEVICE Long size(T const &b) noexcept
integer version
Definition: AMReX_GpuRange.H:26
void streamSynchronize() noexcept
Definition: AMReX_GpuDevice.H:237
int NProcs()
Process ID in MPI_COMM_WORLD.
Definition: AMReX_MPMD.cpp:122
MPI_Comm Communicator() noexcept
Definition: AMReX_ParallelDescriptor.H:210
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
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
static int f(amrex::Real t, N_Vector y_data, N_Vector y_rhs, void *user_data)
Definition: AMReX_SundialsIntegrator.H:44
void FileOpenFailed(const std::string &file)
Output a message and abort when couldn't open the file.
Definition: AMReX_Utility.cpp:131
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 Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition: AMReX.cpp:225
void ExecOnFinalize(std::function< void()>)
We maintain a stack of functions that need to be called in Finalize(). The functions are called in LI...
Definition: AMReX.cpp:314
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
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