diff --git a/amrex/docs_html/_downloads/008eb6dbfab802633dff40122ece848c/amrex.pdf b/amrex/docs_html/_downloads/008eb6dbfab802633dff40122ece848c/amrex.pdf index 442f3a284c..3e7c81a2d8 100644 Binary files a/amrex/docs_html/_downloads/008eb6dbfab802633dff40122ece848c/amrex.pdf and b/amrex/docs_html/_downloads/008eb6dbfab802633dff40122ece848c/amrex.pdf differ diff --git a/amrex/docs_html/doxygen/AMReX__ParticleIO_8H_source.html b/amrex/docs_html/doxygen/AMReX__ParticleIO_8H_source.html index 1e7a2a6877..ef224c9554 100644 --- a/amrex/docs_html/doxygen/AMReX__ParticleIO_8H_source.html +++ b/amrex/docs_html/doxygen/AMReX__ParticleIO_8H_source.html @@ -1398,8 +1398,8 @@
pp
amrex::ParmParse pp
Input file parser instance for the given namespace.
Definition: AMReX_HypreIJIface.cpp:15
AMREX_D_TERM
#define AMREX_D_TERM(a, b, c)
Definition: AMReX_SPACE.H:129
AMReX_WriteBinaryParticleData.H
-
WriteBinaryParticleDataAsync
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:683
-
WriteBinaryParticleDataSync
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:391
+
WriteBinaryParticleDataAsync
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
+
WriteBinaryParticleDataSync
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
amrex::BoxArray
A collection of Boxes stored in an Array.
Definition: AMReX_BoxArray.H:529
amrex::DistributionMapping
Calculates the distribution of FABs to MPI processes.
Definition: AMReX_DistributionMapping.H:41
amrex::MFIter
Definition: AMReX_MFIter.H:57
diff --git a/amrex/docs_html/doxygen/AMReX__WriteBinaryParticleData_8H_source.html b/amrex/docs_html/doxygen/AMReX__WriteBinaryParticleData_8H_source.html index 941912f298..3ee1b6348a 100644 --- a/amrex/docs_html/doxygen/AMReX__WriteBinaryParticleData_8H_source.html +++ b/amrex/docs_html/doxygen/AMReX__WriteBinaryParticleData_8H_source.html @@ -290,977 +290,978 @@
189  const Long rChunkSize = AMREX_SPACEDIM + num_output_real;
190  rdata.resize(np*rChunkSize);
191 
-
192  typename PC::IntVector idata_d(idata.size());
-
193  typename PC::RealVector rdata_d(rdata.size());
-
194 
-
195  typename PC::IntVector write_int_comp_d(write_int_comp.size());
-
196  typename PC::IntVector write_real_comp_d(write_real_comp.size());
-
197  Gpu::copyAsync(Gpu::hostToDevice, write_int_comp.begin(), write_int_comp.end(),
-
198  write_int_comp_d.begin());
-
199  Gpu::copyAsync(Gpu::hostToDevice, write_real_comp.begin(), write_real_comp.end(),
-
200  write_real_comp_d.begin());
-
201  Gpu::Device::streamSynchronize();
-
202 
-
203  const auto write_int_comp_d_ptr = write_int_comp_d.data();
-
204  const auto write_real_comp_d_ptr = write_real_comp_d.data();
-
205 
-
206  std::size_t poffset = 0;
-
207  for (int tile : tiles) {
-
208  const auto& ptile = pc.ParticlesAt(lev, grid, tile);
-
209  const auto& pflags = particle_io_flags[lev].at(std::make_pair(grid, tile));
-
210  int np_tile = ptile.numParticles();
-
211  typename PC::IntVector offsets(np_tile);
-
212  int num_copies = Scan::ExclusiveSum(np_tile, pflags.begin(), offsets.begin(), Scan::retSum);
-
213 
-
214  const auto flag_ptr = pflags.data();
-
215 
-
216  auto idata_d_ptr = idata_d.data();
-
217  auto rdata_d_ptr = rdata_d.data();
-
218 
-
219  const auto ptd = ptile.getConstParticleTileData();
-
220  amrex::ParallelFor(num_copies,
-
221  [=] AMREX_GPU_DEVICE (int pindex) noexcept
-
222  {
-
223  // might be worth using shared memory here
-
224  const auto p = ptd.getSuperParticle(pindex);
-
225 
-
226  if (flag_ptr[pindex]) {
-
227  std::size_t iout_index = (pindex+poffset)*iChunkSize;
-
228  packParticleIDs(&idata_d_ptr[iout_index], p, is_checkpoint);
-
229  iout_index += 2;
-
230 
-
231  std::size_t rout_index = (pindex+poffset)*rChunkSize;
-
232  for (int j = 0; j < AMREX_SPACEDIM; j++) {
-
233  rdata_d_ptr[rout_index] = p.pos(j);
-
234  rout_index++;
-
235  }
-
236 
-
237  for (int j = 0; j < PC::SuperParticleType::NInt; j++) {
-
238  if (write_int_comp_d_ptr[j]) {
-
239  idata_d_ptr[iout_index] = p.idata(j);
-
240  iout_index++;
-
241  }
-
242  }
-
243 
-
244  for (int j = 0; j < ptd.m_num_runtime_int; j++) {
-
245  if (write_int_comp_d_ptr[PC::SuperParticleType::NInt + j]) {
-
246  idata_d_ptr[iout_index] = ptd.m_runtime_idata[j][pindex];
-
247  iout_index++;
-
248  }
-
249  }
-
250 
-
251  // extra SoA Real components
-
252  const int real_start_offset = PC::ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0; // pure SoA: skip positions
-
253  for (int j = real_start_offset; j < PC::SuperParticleType::NReal; j++) {
-
254  const int write_comp_index = j-real_start_offset;
-
255  if (write_real_comp_d_ptr[write_comp_index]) {
-
256  rdata_d_ptr[rout_index] = p.rdata(j);
-
257  rout_index++;
-
258  }
-
259  }
-
260 
-
261  for (int j = 0; j < ptd.m_num_runtime_real; j++) {
-
262  if (write_real_comp_d_ptr[PC::SuperParticleType::NReal+j-real_start_offset]) {
-
263  rdata_d_ptr[rout_index] = ptd.m_runtime_rdata[j][pindex];
-
264  rout_index++;
-
265  }
-
266  }
-
267  }
-
268  });
-
269 
-
270  poffset += num_copies;
-
271  }
-
272 
-
273  Gpu::copyAsync(Gpu::deviceToHost, idata_d.begin(), idata_d.end(), idata.begin());
-
274  Gpu::copyAsync(Gpu::deviceToHost, rdata_d.begin(), rdata_d.end(), rdata.begin());
-
275  Gpu::Device::streamSynchronize();
-
276 }
-
277 
-
278 template <class PC>
-
279 std::enable_if_t<!RunOnGpu<typename PC::template AllocatorType<int>>::value>
-
280 packIOData (Vector<int>& idata, Vector<ParticleReal>& rdata, const PC& pc, int lev, int grid,
-
281  const Vector<int>& write_real_comp, const Vector<int>& write_int_comp,
-
282  const Vector<std::map<std::pair<int, int>, typename PC::IntVector>>& particle_io_flags,
-
283  const Vector<int>& tiles, int np, bool is_checkpoint)
-
284 {
-
285  int num_output_int = 0;
-
286  for (int i = 0; i < pc.NumIntComps() + PC::NStructInt; ++i) {
-
287  if (write_int_comp[i]) { ++num_output_int; }
-
288  }
-
289 
-
290  const Long iChunkSize = 2 + num_output_int;
-
291  idata.resize(np*iChunkSize);
-
292 
-
293  int num_output_real = 0;
-
294  for (int i : write_real_comp) {
-
295  if (i) { ++num_output_real; }
-
296  }
-
297 
-
298  const Long rChunkSize = AMREX_SPACEDIM + num_output_real;
-
299  rdata.resize(np*rChunkSize);
-
300 
-
301  int* iptr = idata.dataPtr();
-
302  ParticleReal* rptr = rdata.dataPtr();
-
303  for (int tile : tiles) {
-
304  const auto& ptile = pc.ParticlesAt(lev, grid, tile);
-
305  const auto& pflags = particle_io_flags[lev].at(std::make_pair(grid, tile));
-
306  for (int pindex = 0; pindex < ptile.numParticles(); ++pindex) {
-
307  if (pflags[pindex]) {
-
308  const auto& soa = ptile.GetStructOfArrays();
-
309 
-
310  // note: for pure SoA particle layouts, we do write the id, cpu and positions as a struct
-
311  // for backwards compatibility with readers
-
312  if constexpr(!PC::ParticleType::is_soa_particle)
-
313  {
-
314  const auto& aos = ptile.GetArrayOfStructs();
-
315  const auto& p = aos[pindex];
-
316 
-
317  // Int: id, cpu
-
318  packParticleIDs(iptr, p, is_checkpoint);
-
319  iptr += 2;
-
320 
-
321  // Real: positions
-
322  for (int j = 0; j < AMREX_SPACEDIM; j++) { rptr[j] = p.pos(j); }
-
323  rptr += AMREX_SPACEDIM;
-
324 
-
325  // extra AoS Int components
-
326  for (int j = 0; j < PC::NStructInt; j++) {
-
327  if (write_int_comp[j]) {
-
328  *iptr = p.idata(j);
-
329  ++iptr;
-
330  }
-
331  }
-
332  // extra AoS Real components
-
333  for (int j = 0; j < PC::NStructReal; j++) {
-
334  if (write_real_comp[j]) {
-
335  *rptr = p.rdata(j);
-
336  ++rptr;
-
337  }
-
338  }
-
339  }
-
340  else {
-
341  uint64_t idcpu = soa.GetIdCPUData()[pindex];
-
342  if (is_checkpoint) {
-
343  std::int32_t xi, yi;
-
344  std::uint32_t xu, yu;
-
345  xu = (std::uint32_t)((idcpu & 0xFFFFFFFF00000000LL) >> 32);
-
346  yu = (std::uint32_t)( idcpu & 0xFFFFFFFFLL);
-
347  std::memcpy(&xi, &xu, sizeof(xu));
-
348  std::memcpy(&yi, &yu, sizeof(yu));
-
349  *iptr = xi;
-
350  iptr += 1;
-
351  *iptr = yi;
-
352  iptr += 1;
-
353  } else {
-
354  // Int: id, cpu
-
355  *iptr = (int) ParticleIDWrapper(idcpu);
-
356  iptr += 1;
-
357  *iptr = (int) ParticleCPUWrapper(idcpu);
-
358  iptr += 1;
-
359  }
-
360 
-
361  // Real: position
-
362  for (int j = 0; j < AMREX_SPACEDIM; j++) { rptr[j] = soa.GetRealData(j)[pindex]; }
-
363  rptr += AMREX_SPACEDIM;
-
364  }
-
365 
-
366  // SoA int data
-
367  const int int_start_offset = 0;
-
368  for (int j = int_start_offset; j < pc.NumIntComps(); j++) {
-
369  if (write_int_comp[PC::NStructInt+j]) {
-
370  *iptr = soa.GetIntData(j)[pindex];
-
371  ++iptr;
-
372  }
-
373  }
-
374 
-
375  // extra SoA Real components
-
376  const int real_start_offset = PC::ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0; // pure SoA: skip positions
-
377  for (int j = real_start_offset; j < pc.NumRealComps(); j++) {
-
378  const int write_comp_index = PC::NStructReal+j-real_start_offset;
-
379  if (write_real_comp[write_comp_index]) {
-
380  *rptr = (ParticleReal) soa.GetRealData(j)[pindex];
-
381  ++rptr;
-
382  }
-
383  }
-
384  }
-
385  }
-
386  }
-
387 }
+
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));
+
273  Gpu::Device::streamSynchronize();
+
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 template <class PC, class F, std::enable_if_t<IsParticleContainer<PC>::value, int> foo = 0>
-
391 void WriteBinaryParticleDataSync (PC const& pc,
-
392  const std::string& dir, const std::string& name,
-
393  const Vector<int>& write_real_comp,
-
394  const Vector<int>& write_int_comp,
-
395  const Vector<std::string>& real_comp_names,
-
396  const Vector<std::string>& int_comp_names,
-
397  F const& f, bool is_checkpoint)
-
398 {
-
399  BL_PROFILE("WriteBinaryParticleData()");
-
400  AMREX_ASSERT(pc.OK());
-
401 
-
402  AMREX_ASSERT(sizeof(typename PC::ParticleType::RealType) == 4 ||
-
403  sizeof(typename PC::ParticleType::RealType) == 8);
-
404 
-
405  constexpr int NStructReal = PC::NStructReal;
-
406  constexpr int NStructInt = PC::NStructInt;
-
407 
-
408  const int NProcs = ParallelDescriptor::NProcs();
-
409  const int IOProcNumber = ParallelDescriptor::IOProcessorNumber();
-
410 
-
411  if constexpr(PC::ParticleType::is_soa_particle) {
-
412  AMREX_ALWAYS_ASSERT(real_comp_names.size() == pc.NumRealComps() + NStructReal - AMREX_SPACEDIM); // pure SoA: skip positions
-
413  } else {
-
414  AMREX_ALWAYS_ASSERT(real_comp_names.size() == pc.NumRealComps() + NStructReal);
-
415  }
-
416  AMREX_ALWAYS_ASSERT( int_comp_names.size() == pc.NumIntComps() + NStructInt);
-
417 
-
418  std::string pdir = dir;
-
419  if ( ! pdir.empty() && pdir[pdir.size()-1] != '/') { pdir += '/'; }
-
420  pdir += name;
-
421 
-
422  if ( ! pc.GetLevelDirectoriesCreated()) {
-
423  if (ParallelDescriptor::IOProcessor())
-
424  {
-
425  if ( ! amrex::UtilCreateDirectory(pdir, 0755))
-
426  {
-
427  amrex::CreateDirectoryFailed(pdir);
-
428  }
-
429  }
-
430  ParallelDescriptor::Barrier();
-
431  }
-
432 
-
433  std::ofstream HdrFile;
-
434 
-
435  Long nparticles = 0;
-
436  Long maxnextid;
-
437 
-
438  // evaluate f for every particle to determine which ones to output
-
439  Vector<std::map<std::pair<int, int>, typename PC::IntVector > >
-
440  particle_io_flags(pc.GetParticles().size());
-
441  for (int lev = 0; lev < pc.GetParticles().size(); lev++)
-
442  {
-
443  const auto& pmap = pc.GetParticles(lev);
-
444  for (const auto& kv : pmap)
-
445  {
-
446  auto& flags = particle_io_flags[lev][kv.first];
-
447  particle_detail::fillFlags(flags, kv.second, f);
-
448  }
-
449  }
-
450 
-
451  Gpu::Device::streamSynchronize();
-
452 
-
453  if(pc.GetUsePrePost())
-
454  {
-
455  nparticles = pc.GetNParticlesPrePost();
-
456  maxnextid = pc.GetMaxNextIDPrePost();
-
457  }
-
458  else
-
459  {
-
460  nparticles = particle_detail::countFlags(particle_io_flags, pc);
-
461  maxnextid = PC::ParticleType::NextID();
-
462  ParallelDescriptor::ReduceLongSum(nparticles, IOProcNumber);
-
463  PC::ParticleType::NextID(maxnextid);
-
464  ParallelDescriptor::ReduceLongMax(maxnextid, IOProcNumber);
-
465  }
-
466 
-
467  if (ParallelDescriptor::IOProcessor())
-
468  {
-
469  std::string HdrFileName = pdir;
-
470 
-
471  if ( ! HdrFileName.empty() && HdrFileName[HdrFileName.size()-1] != '/') {
-
472  HdrFileName += '/';
-
473  }
-
474 
-
475  HdrFileName += "Header";
-
476  pc.HdrFileNamePrePost = HdrFileName;
-
477 
-
478  HdrFile.open(HdrFileName.c_str(), std::ios::out|std::ios::trunc);
-
479 
-
480  if ( ! HdrFile.good()) { amrex::FileOpenFailed(HdrFileName); }
-
481 
-
482  //
-
483  // First thing written is our version string.
-
484  // We append "_single" or "_double" to the version string indicating
-
485  // whether we're using "float" or "double" floating point data.
-
486  //
-
487  std::string version_string = is_checkpoint ? PC::CheckpointVersion() : PC::PlotfileVersion();
-
488  if (sizeof(typename PC::ParticleType::RealType) == 4)
-
489  {
-
490  HdrFile << version_string << "_single" << '\n';
-
491  }
-
492  else
-
493  {
-
494  HdrFile << version_string << "_double" << '\n';
-
495  }
-
496 
-
497  int num_output_real = 0;
-
498  for (int i : write_real_comp) {
-
499  if (i) { ++num_output_real; }
-
500  }
-
501 
-
502  int num_output_int = 0;
-
503  for (int i = 0; i < pc.NumIntComps() + NStructInt; ++i) {
-
504  if (write_int_comp[i]) { ++num_output_int; }
-
505  }
-
506 
-
507  // AMREX_SPACEDIM and N for sanity checking.
-
508  HdrFile << AMREX_SPACEDIM << '\n';
-
509 
-
510  // The number of extra real parameters
-
511  HdrFile << num_output_real << '\n';
-
512 
-
513  // Real component names
-
514  for (int i = 0; i < (int) real_comp_names.size(); ++i ) {
-
515  if (write_real_comp[i]) { HdrFile << real_comp_names[i] << '\n'; }
-
516  }
-
517 
-
518  // The number of extra int parameters
-
519  HdrFile << num_output_int << '\n';
-
520 
-
521  // int component names
-
522  for (int i = 0; i < NStructInt + pc.NumIntComps(); ++i ) {
-
523  if (write_int_comp[i]) { HdrFile << int_comp_names[i] << '\n'; }
-
524  }
-
525 
-
526  bool is_checkpoint_legacy = true; // legacy
-
527  HdrFile << is_checkpoint_legacy << '\n';
-
528 
-
529  // The total number of particles.
-
530  HdrFile << nparticles << '\n';
-
531 
-
532  // The value of nextid that we need to restore on restart.
-
533  HdrFile << maxnextid << '\n';
-
534 
-
535  // Then the finest level of the AMR hierarchy.
-
536  HdrFile << pc.finestLevel() << '\n';
-
537 
-
538  // Then the number of grids at each level.
-
539  for (int lev = 0; lev <= pc.finestLevel(); lev++) {
-
540  HdrFile << pc.ParticleBoxArray(lev).size() << '\n';
-
541  }
-
542  }
-
543 
-
544  // We want to write the data out in parallel.
-
545  // We'll allow up to nOutFiles active writers at a time.
-
546  int nOutFiles(256);
-
547 
-
548  ParmParse pp("particles");
-
549  pp.queryAdd("particles_nfiles",nOutFiles);
-
550  if(nOutFiles == -1) { nOutFiles = NProcs; }
-
551  nOutFiles = std::max(1, std::min(nOutFiles,NProcs));
-
552  pc.nOutFilesPrePost = nOutFiles;
-
553 
-
554  for (int lev = 0; lev <= pc.finestLevel(); lev++)
-
555  {
-
556  bool gotsome;
-
557  if(pc.usePrePost)
-
558  {
-
559  gotsome = (pc.nParticlesAtLevelPrePost[lev] > 0);
-
560  }
-
561  else
-
562  {
-
563  gotsome = (pc.NumberOfParticlesAtLevel(lev) > 0);
-
564  }
-
565 
-
566  // We store the particles at each level in their own subdirectory.
-
567  std::string LevelDir = pdir;
-
568 
-
569  if (gotsome)
-
570  {
-
571  if ( ! LevelDir.empty() && LevelDir[LevelDir.size()-1] != '/') { LevelDir += '/'; }
-
572 
-
573  LevelDir = amrex::Concatenate(LevelDir.append("Level_"), lev, 1);
-
574 
-
575  if ( ! pc.GetLevelDirectoriesCreated())
-
576  {
-
577  if (ParallelDescriptor::IOProcessor()) {
-
578  if ( ! amrex::UtilCreateDirectory(LevelDir, 0755)) {
-
579  amrex::CreateDirectoryFailed(LevelDir);
-
580  }
-
581  }
-
582  ParallelDescriptor::Barrier();
-
583  }
-
584  }
-
585 
-
586  // Write out the header for each particle
-
587  if (gotsome && ParallelDescriptor::IOProcessor()) {
-
588  std::string HeaderFileName = LevelDir;
-
589  HeaderFileName += "/Particle_H";
-
590  std::ofstream ParticleHeader(HeaderFileName);
-
591 
-
592  pc.ParticleBoxArray(lev).writeOn(ParticleHeader);
-
593  ParticleHeader << '\n';
-
594 
-
595  ParticleHeader.flush();
-
596  ParticleHeader.close();
-
597  }
-
598 
-
599  MFInfo info;
-
600  info.SetAlloc(false);
-
601  MultiFab state(pc.ParticleBoxArray(lev),
-
602  pc.ParticleDistributionMap(lev),
-
603  1,0,info);
-
604 
-
605  // We eventually want to write out the file name and the offset
-
606  // into that file into which each grid of particles is written.
-
607  Vector<int> which(state.size(),0);
-
608  Vector<int > count(state.size(),0);
-
609  Vector<Long> where(state.size(),0);
-
610 
-
611  std::string filePrefix(LevelDir);
-
612  filePrefix += '/';
-
613  filePrefix += PC::DataPrefix();
-
614  if(pc.usePrePost) {
-
615  pc.filePrefixPrePost[lev] = filePrefix;
-
616  }
-
617  bool groupSets(false), setBuf(true);
-
618 
-
619  if (gotsome)
-
620  {
-
621  for(NFilesIter nfi(nOutFiles, filePrefix, groupSets, setBuf); nfi.ReadyToWrite(); ++nfi)
-
622  {
-
623  auto& myStream = (std::ofstream&) nfi.Stream();
-
624  pc.WriteParticles(lev, myStream, nfi.FileNumber(), which, count, where,
-
625  write_real_comp, write_int_comp, particle_io_flags, is_checkpoint);
-
626  }
-
627 
-
628  if(pc.usePrePost) {
-
629  pc.whichPrePost[lev] = which;
-
630  pc.countPrePost[lev] = count;
-
631  pc.wherePrePost[lev] = where;
-
632  } else {
-
633  ParallelDescriptor::ReduceIntSum (which.dataPtr(), static_cast<int>(which.size()), IOProcNumber);
-
634  ParallelDescriptor::ReduceIntSum (count.dataPtr(), static_cast<int>(count.size()), IOProcNumber);
-
635  ParallelDescriptor::ReduceLongSum(where.dataPtr(), static_cast<int>(where.size()), IOProcNumber);
-
636  }
-
637  }
-
638 
-
639  if (ParallelDescriptor::IOProcessor())
-
640  {
-
641  if(pc.GetUsePrePost()) {
-
642  // ---- write to the header and unlink in CheckpointPost
-
643  } else {
-
644  for (int j = 0; j < state.size(); j++)
-
645  {
-
646  HdrFile << which[j] << ' ' << count[j] << ' ' << where[j] << '\n';
-
647  }
-
648 
-
649  if (gotsome && pc.doUnlink)
-
650  {
-
651  // Unlink any zero-length data files.
-
652  Vector<Long> cnt(nOutFiles,0);
-
653 
-
654  for (int i = 0, N=static_cast<int>(count.size()); i < N; i++) {
-
655  cnt[which[i]] += count[i];
-
656  }
-
657 
-
658  for (int i = 0, N=static_cast<int>(cnt.size()); i < N; i++)
-
659  {
-
660  if (cnt[i] == 0)
-
661  {
-
662  std::string FullFileName = NFilesIter::FileName(i, filePrefix);
-
663  FileSystem::Remove(FullFileName);
-
664  }
-
665  }
-
666  }
-
667  }
-
668  }
-
669  }
-
670 
-
671  if (ParallelDescriptor::IOProcessor())
-
672  {
-
673  HdrFile.flush();
-
674  HdrFile.close();
-
675  if ( ! HdrFile.good())
-
676  {
-
677  amrex::Abort("ParticleContainer::Checkpoint(): problem writing HdrFile");
-
678  }
-
679  }
-
680 }
-
681 
-
682 template <class PC, std::enable_if_t<IsParticleContainer<PC>::value, int> foo = 0>
-
683 void WriteBinaryParticleDataAsync (PC const& pc,
-
684  const std::string& dir, const std::string& name,
-
685  const Vector<int>& write_real_comp,
-
686  const Vector<int>& write_int_comp,
-
687  const Vector<std::string>& real_comp_names,
-
688  const Vector<std::string>& int_comp_names, bool is_checkpoint)
-
689 {
-
690  BL_PROFILE("WriteBinaryParticleDataAsync");
-
691  AMREX_ASSERT(pc.OK());
-
692 
-
693  AMREX_ASSERT(sizeof(typename PC::ParticleType::RealType) == 4 ||
-
694  sizeof(typename PC::ParticleType::RealType) == 8);
-
695 
-
696  constexpr int NStructReal = PC::NStructReal;
-
697  constexpr int NStructInt = PC::NStructInt;
-
698  constexpr int NArrayReal = PC::NArrayReal;
-
699  constexpr int NArrayInt = PC::NArrayInt;
-
700 
-
701  const int MyProc = ParallelDescriptor::MyProc();
-
702  const int NProcs = ParallelDescriptor::NProcs();
-
703  const int IOProcNumber = NProcs - 1;
-
704 
-
705  if constexpr(PC::ParticleType::is_soa_particle) {
-
706  AMREX_ALWAYS_ASSERT(real_comp_names.size() == pc.NumRealComps() + NStructReal - AMREX_SPACEDIM); // pure SoA: skip positions
-
707  } else {
-
708  AMREX_ALWAYS_ASSERT(real_comp_names.size() == pc.NumRealComps() + NStructReal);
-
709  }
-
710  AMREX_ALWAYS_ASSERT( int_comp_names.size() == pc.NumIntComps() + NStructInt);
-
711 
-
712  Vector<LayoutData<Long> > np_per_grid_local(pc.finestLevel()+1);
-
713  for (int lev = 0; lev <= pc.finestLevel(); lev++)
-
714  {
-
715  np_per_grid_local[lev].define(pc.ParticleBoxArray(lev), pc.ParticleDistributionMap(lev));
-
716  using ParIter = typename PC::ParConstIterType;
-
717  for (ParIter pti(pc, lev); pti.isValid(); ++pti)
-
718  {
-
719  int gid = pti.index();
-
720  const auto& ptile = pc.ParticlesAt(lev, pti);
-
721  const auto& ptd = ptile.getConstParticleTileData();
-
722  const int np = ptile.numParticles();
-
723 
-
724  ReduceOps<ReduceOpSum> reduce_op;
-
725  ReduceData<int> reduce_data(reduce_op);
-
726  using ReduceTuple = typename decltype(reduce_data)::Type;
-
727 
-
728  reduce_op.eval(np, reduce_data,
-
729  [=] AMREX_GPU_DEVICE (int i) -> ReduceTuple
-
730  {
-
731  return (ptd.id(i) > 0) ? 1 : 0;
-
732  });
-
733 
-
734  int np_valid = amrex::get<0>(reduce_data.value(reduce_op));
-
735  np_per_grid_local[lev][gid] += np_valid;
-
736  }
-
737  }
-
738 
-
739  Vector<Vector<Long> > np_per_grid_global(pc.finestLevel()+1);
-
740  Long total_np = 0;
-
741  Vector<Long> np_per_level(pc.finestLevel()+1);
-
742  for (int lev = 0; lev <= pc.finestLevel(); lev++)
-
743  {
-
744  np_per_grid_global[lev].resize(np_per_grid_local[lev].size());
-
745  ParallelDescriptor::GatherLayoutDataToVector(np_per_grid_local[lev],
-
746  np_per_grid_global[lev],
-
747  IOProcNumber);
-
748  np_per_level[lev] = std::accumulate(np_per_grid_global[lev].begin(),
-
749  np_per_grid_global[lev].end(), 0L);
-
750  total_np += np_per_level[lev];
-
751  }
-
752 
-
753  std::string pdir = dir;
-
754  if ( ! pdir.empty() && pdir[pdir.size()-1] != '/') { pdir += '/'; }
-
755  pdir += name;
-
756 
-
757  if (MyProc == IOProcNumber)
-
758  {
-
759  if ( ! pc.GetLevelDirectoriesCreated())
-
760  {
-
761  if ( ! amrex::UtilCreateDirectory(pdir, 0755))
-
762  {
-
763  amrex::CreateDirectoryFailed(pdir);
-
764  }
-
765  }
-
766 
-
767  for (int lev = 0; lev <= pc.finestLevel(); lev++)
-
768  {
-
769  std::string LevelDir = pdir;
-
770  bool gotsome = np_per_level[lev];
-
771 
-
772  if (gotsome)
-
773  {
-
774  if ( ! LevelDir.empty() && LevelDir[LevelDir.size()-1] != '/') { LevelDir += '/'; }
-
775 
-
776  LevelDir = amrex::Concatenate(LevelDir.append("Level_"), lev, 1);
-
777 
-
778  if ( ! pc.GetLevelDirectoriesCreated())
-
779  {
-
780  if ( ! amrex::UtilCreateDirectory(LevelDir, 0755))
-
781  {
-
782  amrex::CreateDirectoryFailed(LevelDir);
-
783  }
-
784  }
-
785 
-
786  std::string HeaderFileName = LevelDir;
-
787  HeaderFileName += "/Particle_H";
-
788  std::ofstream ParticleHeader(HeaderFileName);
-
789 
-
790  pc.ParticleBoxArray(lev).writeOn(ParticleHeader);
-
791  ParticleHeader << '\n';
-
792 
-
793  ParticleHeader.flush();
-
794  ParticleHeader.close();
-
795  }
-
796  }
-
797  }
-
798  ParallelDescriptor::Barrier();
-
799 
-
800  Long maxnextid = PC::ParticleType::NextID();
-
801  ParallelDescriptor::ReduceLongMax(maxnextid, IOProcNumber);
-
802 
-
803  Vector<Long> np_on_rank(NProcs, 0L);
-
804  std::size_t psize = particle_detail::PSizeInFile<ParticleReal>(write_real_comp, write_int_comp);
-
805  Vector<int64_t> rank_start_offset(NProcs);
-
806  if (MyProc == IOProcNumber)
-
807  {
-
808  for (int lev = 0; lev <= pc.finestLevel(); lev++)
-
809  {
-
810  for (int k = 0; k < pc.ParticleBoxArray(lev).size(); ++k)
-
811  {
-
812  int rank = pc.ParticleDistributionMap(lev)[k];
-
813  np_on_rank[rank] += np_per_grid_global[lev][k];
-
814  }
-
815  }
-
816 
-
817  for (int ip = 0; ip < NProcs; ++ip)
-
818  {
-
819  auto info = AsyncOut::GetWriteInfo(ip);
-
820  rank_start_offset[ip] = (info.ispot == 0) ? 0 : static_cast<int64_t>(rank_start_offset[ip-1] + np_on_rank[ip-1]*psize);
-
821  }
-
822  }
-
823 
-
824  // make tmp particle tiles in pinned memory to write
-
825  using PinnedPTile = ParticleTile<typename PC::ParticleType, NArrayReal, NArrayInt,
-
826  PinnedArenaAllocator>;
-
827  auto myptiles = std::make_shared<Vector<std::map<std::pair<int, int>,PinnedPTile> > >();
-
828  myptiles->resize(pc.finestLevel()+1);
-
829  for (int lev = 0; lev <= pc.finestLevel(); lev++)
-
830  {
-
831  for (MFIter mfi = pc.MakeMFIter(lev); mfi.isValid(); ++mfi)
-
832  {
-
833  auto& new_ptile = (*myptiles)[lev][std::make_pair(mfi.index(),
-
834  mfi.LocalTileIndex())];
-
835 
-
836  if (np_per_grid_local[lev][mfi.index()] > 0)
-
837  {
-
838  const auto& ptile = pc.ParticlesAt(lev, mfi);
-
839 
-
840  const auto np = np_per_grid_local[lev][mfi.index()];
-
841 
-
842  new_ptile.resize(np);
-
843 
-
844  const auto runtime_real_comps = ptile.NumRuntimeRealComps();
-
845  const auto runtime_int_comps = ptile.NumRuntimeIntComps();
-
846 
-
847  new_ptile.define(runtime_real_comps, runtime_int_comps);
-
848 
-
849  for (auto comp(0); comp < runtime_real_comps; ++comp) {
-
850  new_ptile.push_back_real(NArrayReal+comp, np, 0.);
-
851  }
-
852 
-
853  for (auto comp(0); comp < runtime_int_comps; ++comp) {
-
854  new_ptile.push_back_int(NArrayInt+comp, np, 0);
-
855  }
-
856 
-
857  amrex::filterParticles(new_ptile, ptile, KeepValidFilter());
-
858  }
-
859  }
-
860  }
-
861 
-
862  int finest_level = pc.finestLevel();
-
863  Vector<BoxArray> bas;
-
864  Vector<DistributionMapping> dms;
-
865  for (int lev = 0; lev <= pc.finestLevel(); lev++)
-
866  {
-
867  bas.push_back(pc.ParticleBoxArray(lev));
-
868  dms.push_back(pc.ParticleDistributionMap(lev));
-
869  }
-
870 
-
871  int nrc = pc.NumRealComps();
-
872  int nic = pc.NumIntComps();
-
873  int rnames_size = (int) real_comp_names.size();
-
874 
-
875  auto RD = pc.ParticleRealDescriptor;
-
876 
-
877  AsyncOut::Submit([=] ()
-
878 #if defined(__GNUC__) && (__GNUC__ == 8) && (__GNUC_MINOR__ == 1)
-
879  mutable // workaround for bug in gcc 8.1
-
880 #endif
-
881  {
-
882  if (MyProc == IOProcNumber)
-
883  {
-
884  std::string HdrFileName = pdir;
-
885  std::ofstream HdrFile;
-
886 
-
887  if ( ! HdrFileName.empty() && HdrFileName[HdrFileName.size()-1] != '/') {
-
888  HdrFileName += '/';
-
889  }
-
890 
-
891  HdrFileName += "Header";
-
892 
-
893  HdrFile.open(HdrFileName.c_str(), std::ios::out|std::ios::trunc);
-
894 
-
895  if ( ! HdrFile.good()) { amrex::FileOpenFailed(HdrFileName); }
-
896 
-
897  std::string version_string = is_checkpoint ? PC::CheckpointVersion() : PC::PlotfileVersion();
-
898  if (sizeof(typename PC::ParticleType::RealType) == 4)
-
899  {
-
900  HdrFile << version_string << "_single" << '\n';
-
901  }
-
902  else
-
903  {
-
904  HdrFile << version_string << "_double" << '\n';
-
905  }
-
906 
-
907  int num_output_real = 0;
-
908  for (int i = 0; i < rnames_size; ++i) {
-
909  if (write_real_comp[i]) { ++num_output_real; }
-
910  }
-
911 
-
912  int num_output_int = 0;
-
913  for (int i = 0; i < nic + NStructInt; ++i) {
-
914  if (write_int_comp[i]) { ++num_output_int; }
-
915  }
-
916 
-
917  // AMREX_SPACEDIM and N for sanity checking.
-
918  HdrFile << AMREX_SPACEDIM << '\n';
-
919 
-
920  // The number of extra real parameters
-
921  HdrFile << num_output_real << '\n';
-
922 
-
923  // Real component names
-
924  for (int i = 0; i < rnames_size; ++i ) {
-
925  if (write_real_comp[i]) { HdrFile << real_comp_names[i] << '\n'; }
-
926  }
-
927 
-
928  // The number of extra int parameters
-
929  HdrFile << num_output_int << '\n';
-
930 
-
931  // int component names
-
932  for (int i = 0; i < NStructInt + nic; ++i ) {
-
933  if (write_int_comp[i]) { HdrFile << int_comp_names[i] << '\n'; }
-
934  }
-
935 
-
936  bool is_checkpoint_legacy = true; // legacy
-
937  HdrFile << is_checkpoint_legacy << '\n';
-
938 
-
939  // The total number of particles.
-
940  HdrFile << total_np << '\n';
-
941 
-
942  // The value of nextid that we need to restore on restart.
-
943  HdrFile << maxnextid << '\n';
-
944 
-
945  // Then the finest level of the AMR hierarchy.
-
946  HdrFile << finest_level << '\n';
-
947 
-
948  // Then the number of grids at each level.
-
949  for (int lev = 0; lev <= finest_level; lev++) {
-
950  HdrFile << dms[lev].size() << '\n';
-
951  }
-
952 
-
953  for (int lev = 0; lev <= finest_level; lev++)
-
954  {
-
955  Vector<int64_t> grid_offset(NProcs, 0);
-
956  for (int k = 0; k < bas[lev].size(); ++k)
-
957  {
-
958  int rank = dms[lev][k];
-
959  auto info = AsyncOut::GetWriteInfo(rank);
-
960  HdrFile << info.ifile << ' '
-
961  << np_per_grid_global[lev][k] << ' '
-
962  << grid_offset[rank] + rank_start_offset[rank] << '\n';
-
963  grid_offset[rank] += static_cast<int64_t>(np_per_grid_global[lev][k]*psize);
-
964  }
-
965  }
-
966 
-
967  HdrFile.flush();
-
968  HdrFile.close();
-
969  if ( ! HdrFile.good())
-
970  {
-
971  amrex::Abort("ParticleContainer::Checkpoint(): problem writing HdrFile");
-
972  }
-
973  }
-
974 
-
975  AsyncOut::Wait(); // Wait for my turn
-
976 
-
977  for (int lev = 0; lev <= finest_level; lev++)
-
978  {
-
979  // For a each grid, the tiles it contains
-
980  std::map<int, Vector<int> > tile_map;
-
981 
-
982  for (const auto& kv : (*myptiles)[lev])
-
983  {
-
984  const int grid = kv.first.first;
-
985  const int tile = kv.first.second;
-
986  tile_map[grid].push_back(tile);
-
987  }
-
988 
-
989  std::string LevelDir = pdir;
-
990  if ( ! LevelDir.empty() && LevelDir[LevelDir.size()-1] != '/') { LevelDir += '/'; }
-
991  LevelDir = amrex::Concatenate(LevelDir.append("Level_"), lev, 1);
-
992  std::string filePrefix(LevelDir);
-
993  filePrefix += '/';
-
994  filePrefix += PC::DataPrefix();
-
995  auto info = AsyncOut::GetWriteInfo(MyProc);
-
996  std::string file_name = amrex::Concatenate(filePrefix, info.ifile, 5);
-
997  std::ofstream ofs;
-
998  ofs.open(file_name.c_str(), (info.ispot == 0) ? (std::ios::binary | std::ios::trunc)
-
999  : (std::ios::binary | std::ios::app));
-
1000 
-
1001  for (int k = 0; k < bas[lev].size(); ++k)
-
1002  {
-
1003  int rank = dms[lev][k];
-
1004  if (rank != MyProc) { continue; }
-
1005  const int grid = k;
-
1006  if (np_per_grid_local[lev][grid] == 0) { continue; }
-
1007 
-
1008  // First write out the integer data in binary.
-
1009  int num_output_int = 0;
-
1010  for (int i = 0; i < nic + NStructInt; ++i) {
-
1011  if (write_int_comp[i]) { ++num_output_int; }
-
1012  }
-
1013 
-
1014  const Long iChunkSize = 2 + num_output_int;
-
1015  Vector<int> istuff(np_per_grid_local[lev][grid]*iChunkSize);
-
1016  int* iptr = istuff.dataPtr();
-
1017 
-
1018  for (unsigned i = 0; i < tile_map[grid].size(); i++) {
-
1019  auto ptile_index = std::make_pair(grid, tile_map[grid][i]);
-
1020  const auto& pbox = (*myptiles)[lev][ptile_index];
-
1021  const auto ptd = pbox.getConstParticleTileData();
-
1022  for (int pindex = 0; pindex < pbox.numParticles(); ++pindex)
-
1023  {
-
1024  const auto& soa = pbox.GetStructOfArrays();
-
1025 
-
1026  const auto& p = make_particle<typename PC::ConstParticleType>{}(ptd, pindex);
-
1027  if (p.id() <= 0) { continue; }
-
1028 
-
1029  // note: for pure SoA particle layouts, we do write the id, cpu and positions as a struct
-
1030  // for backwards compatibility with readers
-
1031  if constexpr(!PC::ParticleType::is_soa_particle)
-
1032  {
-
1033  // Ints: id, cpu
-
1034  particle_detail::packParticleIDs(iptr, p, is_checkpoint);
-
1035  iptr += 2;
-
1036 
-
1037  // extra AoS Int components
-
1038  for (int j = 0; j < NStructInt; j++)
-
1039  {
-
1040  if (write_int_comp[j])
-
1041  {
-
1042  *iptr = p.idata(j);
-
1043  ++iptr;
-
1044  }
-
1045  }
-
1046  }
-
1047  else {
-
1048  uint64_t idcpu = soa.GetIdCPUData()[pindex];
-
1049  if (is_checkpoint) {
-
1050  std::int32_t xi, yi;
-
1051  std::uint32_t xu, yu;
-
1052  xu = (std::uint32_t)((idcpu & 0xFFFFFFFF00000000LL) >> 32);
-
1053  yu = (std::uint32_t)( idcpu & 0xFFFFFFFFLL);
-
1054  std::memcpy(&xi, &xu, sizeof(xu));
-
1055  std::memcpy(&yi, &yu, sizeof(yu));
-
1056  *iptr = xi;
-
1057  iptr += 1;
-
1058  *iptr = yi;
-
1059  iptr += 1;
-
1060  } else {
-
1061  // Int: id, cpu
-
1062  *iptr = (int) ParticleIDWrapper(idcpu);
-
1063  iptr += 1;
-
1064  *iptr = (int) ParticleCPUWrapper(idcpu);
-
1065  iptr += 1;
-
1066  }
-
1067  }
-
1068 
-
1069  // extra SoA Ints
-
1070  const int int_start_offset = 0;
-
1071  for (int j = int_start_offset; j < nic; j++)
-
1072  {
-
1073  if (write_int_comp[NStructInt+j])
-
1074  {
-
1075  *iptr = soa.GetIntData(j)[pindex];
-
1076  ++iptr;
-
1077  }
-
1078  }
-
1079  }
-
1080  }
-
1081 
-
1082  writeIntData(istuff.dataPtr(), istuff.size(), ofs);
-
1083  ofs.flush(); // Some systems require this flush() (probably due to a bug)
-
1084 
-
1085  // Write the Real data in binary.
-
1086  int num_output_real = 0;
-
1087  for (int i = 0; i < rnames_size; ++i) {
-
1088  if (write_real_comp[i]) { ++num_output_real; }
-
1089  }
-
1090 
-
1091  const Long rChunkSize = AMREX_SPACEDIM + num_output_real;
-
1092  Vector<typename PC::ParticleType::RealType> rstuff(np_per_grid_local[lev][grid]*rChunkSize);
-
1093  typename PC::ParticleType::RealType* rptr = rstuff.dataPtr();
-
1094 
-
1095  for (unsigned i = 0; i < tile_map[grid].size(); i++) {
-
1096  auto ptile_index = std::make_pair(grid, tile_map[grid][i]);
-
1097  const auto& pbox = (*myptiles)[lev][ptile_index];
-
1098  const auto ptd = pbox.getConstParticleTileData();
-
1099  for (int pindex = 0;
-
1100  pindex < pbox.numParticles(); ++pindex)
-
1101  {
-
1102  const auto& soa = pbox.GetStructOfArrays();
-
1103  const auto& p = make_particle<typename PC::ConstParticleType>{}(ptd, pindex);
-
1104 
-
1105  if (p.id() <= 0) { continue; }
-
1106 
-
1107  if constexpr(!PC::ParticleType::is_soa_particle)
-
1108  {
-
1109  // Real: position
-
1110  for (int j = 0; j < AMREX_SPACEDIM; j++) { rptr[j] = p.pos(j); }
-
1111  rptr += AMREX_SPACEDIM;
-
1112 
-
1113  // extra AoS real
-
1114  for (int j = 0; j < NStructReal; j++)
-
1115  {
-
1116  if (write_real_comp[j])
-
1117  {
-
1118  *rptr = p.rdata(j);
-
1119  ++rptr;
-
1120  }
-
1121  }
-
1122  }
-
1123  else {
-
1124  // Real: position
-
1125  for (int j = 0; j < AMREX_SPACEDIM; j++) { rptr[j] = soa.GetRealData(j)[pindex]; }
-
1126  rptr += AMREX_SPACEDIM;
-
1127  }
-
1128 
-
1129  // extra SoA real
-
1130  const int real_start_offset = PC::ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0; // pure SoA: positions
-
1131  for (int j = real_start_offset; j < nrc; j++)
-
1132  {
-
1133  const int write_comp_offset = PC::ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0; // pure SoA: skip positions
-
1134  const int write_comp_index = PC::NStructReal+j-write_comp_offset;
-
1135  if (write_real_comp[write_comp_index])
-
1136  {
-
1137  *rptr = (typename PC::ParticleType::RealType) soa.GetRealData(j)[pindex];
-
1138  ++rptr;
-
1139  }
-
1140  }
-
1141  }
-
1142  }
-
1143 
-
1144  if (sizeof(typename PC::ParticleType::RealType) == 4) {
-
1145  writeFloatData((float*) rstuff.dataPtr(), rstuff.size(), ofs, RD);
-
1146  }
-
1147  else if (sizeof(typename PC::ParticleType::RealType) == 8) {
-
1148  writeDoubleData((double*) rstuff.dataPtr(), rstuff.size(), ofs, RD);
-
1149  }
-
1150 
-
1151  ofs.flush(); // Some systems require this flush() (probably due to a bug)
-
1152  }
-
1153  }
-
1154  AsyncOut::Notify(); // Notify others I am done
-
1155  });
-
1156 }
-
1157 
-
1158 #ifdef AMREX_USE_HDF5
-
1159 #include <AMReX_WriteBinaryParticleDataHDF5.H>
-
1160 #endif
-
1161 
-
1162 #endif /*AMREX_WRITE_BINARY_PARTICLE_DATA_H*/
+
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()) {
+
424  if (ParallelDescriptor::IOProcessor())
+
425  {
+
426  if ( ! amrex::UtilCreateDirectory(pdir, 0755))
+
427  {
+
428  amrex::CreateDirectoryFailed(pdir);
+
429  }
+
430  }
+
431  ParallelDescriptor::Barrier();
+
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 
+
452  Gpu::Device::streamSynchronize();
+
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 
+
468  if (ParallelDescriptor::IOProcessor())
+
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  {
+
578  if (ParallelDescriptor::IOProcessor()) {
+
579  if ( ! amrex::UtilCreateDirectory(LevelDir, 0755)) {
+
580  amrex::CreateDirectoryFailed(LevelDir);
+
581  }
+
582  }
+
583  ParallelDescriptor::Barrier();
+
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 
+
640  if (ParallelDescriptor::IOProcessor())
+
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 
+
672  if (ParallelDescriptor::IOProcessor())
+
673  {
+
674  HdrFile.flush();
+
675  HdrFile.close();
+
676  if ( ! HdrFile.good())
+
677  {
+
678  amrex::Abort("ParticleContainer::Checkpoint(): 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  {
+
764  amrex::CreateDirectoryFailed(pdir);
+
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  {
+
783  amrex::CreateDirectoryFailed(LevelDir);
+
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  }
+
799  ParallelDescriptor::Barrier();
+
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("ParticleContainer::Checkpoint(): 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 
+
1027  const auto& p = make_particle<typename PC::ConstParticleType>{}(ptd, pindex);
+
1028  if (p.id() <= 0) { continue; }
+
1029 
+
1030  // note: for pure SoA particle layouts, we do write the id, cpu and positions as a struct
+
1031  // for backwards compatibility with readers
+
1032  if constexpr(!PC::ParticleType::is_soa_particle)
+
1033  {
+
1034  // Ints: id, cpu
+
1035  particle_detail::packParticleIDs(iptr, p, is_checkpoint);
+
1036  iptr += 2;
+
1037 
+
1038  // extra AoS Int components
+
1039  for (int j = 0; j < NStructInt; j++)
+
1040  {
+
1041  if (write_int_comp[j])
+
1042  {
+
1043  *iptr = p.idata(j);
+
1044  ++iptr;
+
1045  }
+
1046  }
+
1047  }
+
1048  else {
+
1049  uint64_t idcpu = soa.GetIdCPUData()[pindex];
+
1050  if (is_checkpoint) {
+
1051  std::int32_t xi, yi;
+
1052  std::uint32_t xu, yu;
+
1053  xu = (std::uint32_t)((idcpu & 0xFFFFFFFF00000000LL) >> 32);
+
1054  yu = (std::uint32_t)( idcpu & 0xFFFFFFFFLL);
+
1055  std::memcpy(&xi, &xu, sizeof(xu));
+
1056  std::memcpy(&yi, &yu, sizeof(yu));
+
1057  *iptr = xi;
+
1058  iptr += 1;
+
1059  *iptr = yi;
+
1060  iptr += 1;
+
1061  } else {
+
1062  // Int: id, cpu
+
1063  *iptr = (int) ParticleIDWrapper(idcpu);
+
1064  iptr += 1;
+
1065  *iptr = (int) ParticleCPUWrapper(idcpu);
+
1066  iptr += 1;
+
1067  }
+
1068  }
+
1069 
+
1070  // extra SoA Ints
+
1071  const int int_start_offset = 0;
+
1072  for (int j = int_start_offset; j < nic; j++)
+
1073  {
+
1074  if (write_int_comp[NStructInt+j])
+
1075  {
+
1076  *iptr = soa.GetIntData(j)[pindex];
+
1077  ++iptr;
+
1078  }
+
1079  }
+
1080  }
+
1081  }
+
1082 
+
1083  writeIntData(istuff.dataPtr(), istuff.size(), ofs);
+
1084  ofs.flush(); // Some systems require this flush() (probably due to a bug)
+
1085 
+
1086  // Write the Real data in binary.
+
1087  int num_output_real = 0;
+
1088  for (int i = 0; i < rnames_size; ++i) {
+
1089  if (write_real_comp[i]) { ++num_output_real; }
+
1090  }
+
1091 
+
1092  const Long rChunkSize = AMREX_SPACEDIM + num_output_real;
+
1093  Vector<typename PC::ParticleType::RealType> rstuff(np_per_grid_local[lev][grid]*rChunkSize);
+
1094  typename PC::ParticleType::RealType* rptr = rstuff.dataPtr();
+
1095 
+
1096  for (unsigned i = 0; i < tile_map[grid].size(); i++) {
+
1097  auto ptile_index = std::make_pair(grid, tile_map[grid][i]);
+
1098  const auto& pbox = (*myptiles)[lev][ptile_index];
+
1099  const auto ptd = pbox.getConstParticleTileData();
+
1100  for (int pindex = 0;
+
1101  pindex < pbox.numParticles(); ++pindex)
+
1102  {
+
1103  const auto& soa = pbox.GetStructOfArrays();
+
1104  const auto& p = make_particle<typename PC::ConstParticleType>{}(ptd, pindex);
+
1105 
+
1106  if (p.id() <= 0) { continue; }
+
1107 
+
1108  if constexpr(!PC::ParticleType::is_soa_particle)
+
1109  {
+
1110  // Real: position
+
1111  for (int j = 0; j < AMREX_SPACEDIM; j++) { rptr[j] = p.pos(j); }
+
1112  rptr += AMREX_SPACEDIM;
+
1113 
+
1114  // extra AoS real
+
1115  for (int j = 0; j < NStructReal; j++)
+
1116  {
+
1117  if (write_real_comp[j])
+
1118  {
+
1119  *rptr = p.rdata(j);
+
1120  ++rptr;
+
1121  }
+
1122  }
+
1123  }
+
1124  else {
+
1125  // Real: position
+
1126  for (int j = 0; j < AMREX_SPACEDIM; j++) { rptr[j] = soa.GetRealData(j)[pindex]; }
+
1127  rptr += AMREX_SPACEDIM;
+
1128  }
+
1129 
+
1130  // extra SoA real
+
1131  const int real_start_offset = PC::ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0; // pure SoA: positions
+
1132  for (int j = real_start_offset; j < nrc; j++)
+
1133  {
+
1134  const int write_comp_offset = PC::ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0; // pure SoA: skip positions
+
1135  const int write_comp_index = PC::NStructReal+j-write_comp_offset;
+
1136  if (write_real_comp[write_comp_index])
+
1137  {
+
1138  *rptr = (typename PC::ParticleType::RealType) soa.GetRealData(j)[pindex];
+
1139  ++rptr;
+
1140  }
+
1141  }
+
1142  }
+
1143  }
+
1144 
+
1145  if (sizeof(typename PC::ParticleType::RealType) == 4) {
+
1146  writeFloatData((float*) rstuff.dataPtr(), rstuff.size(), ofs, RD);
+
1147  }
+
1148  else if (sizeof(typename PC::ParticleType::RealType) == 8) {
+
1149  writeDoubleData((double*) rstuff.dataPtr(), rstuff.size(), ofs, RD);
+
1150  }
+
1151 
+
1152  ofs.flush(); // Some systems require this flush() (probably due to a bug)
+
1153  }
+
1154  }
+
1155  AsyncOut::Notify(); // Notify others I am done
+
1156  });
+
1157 }
+
1158 
+
1159 #ifdef AMREX_USE_HDF5
+
1160 #include <AMReX_WriteBinaryParticleDataHDF5.H>
+
1161 #endif
+
1162 
+
1163 #endif /*AMREX_WRITE_BINARY_PARTICLE_DATA_H*/
BL_PROFILE
#define BL_PROFILE(a)
Definition: AMReX_BLProfiler.H:558
AMREX_ASSERT
#define AMREX_ASSERT(EX)
Definition: AMReX_BLassert.H:38
AMREX_ALWAYS_ASSERT
#define AMREX_ALWAYS_ASSERT(EX)
Definition: AMReX_BLassert.H:50
@@ -1271,8 +1272,8 @@
AMReX_ParticleUtil.H
AMReX_TypeTraits.H
AMReX_WriteBinaryParticleDataHDF5.H
-
WriteBinaryParticleDataAsync
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:683
-
WriteBinaryParticleDataSync
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:391
+
WriteBinaryParticleDataAsync
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
+
WriteBinaryParticleDataSync
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
if(!(yy_init))
Definition: amrex_iparser.lex.nolint.H:924
amrex::ParmParse::queryAdd
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:954
amrex::AsyncOut::GetWriteInfo
WriteInfo GetWriteInfo(int rank)
Definition: AMReX_AsyncOut.cpp:72
diff --git a/amrex/docs_xml/doxygen/AMReX__WriteBinaryParticleData_8H.xml b/amrex/docs_xml/doxygen/AMReX__WriteBinaryParticleData_8H.xml index 05c61a4458..fdbce5e4ca 100644 --- a/amrex/docs_xml/doxygen/AMReX__WriteBinaryParticleData_8H.xml +++ b/amrex/docs_xml/doxygen/AMReX__WriteBinaryParticleData_8H.xml @@ -2284,7 +2284,7 @@ - + @@ -2340,7 +2340,7 @@ - + @@ -2539,977 +2539,978 @@ constLongrChunkSize=AMREX_SPACEDIM+num_output_real; rdata.resize(np*rChunkSize); -typenamePC::IntVectoridata_d(idata.size()); -typenamePC::RealVectorrdata_d(rdata.size()); - -typenamePC::IntVectorwrite_int_comp_d(write_int_comp.size()); -typenamePC::IntVectorwrite_real_comp_d(write_real_comp.size()); -Gpu::copyAsync(Gpu::hostToDevice,write_int_comp.begin(),write_int_comp.end(), -write_int_comp_d.begin()); -Gpu::copyAsync(Gpu::hostToDevice,write_real_comp.begin(),write_real_comp.end(), -write_real_comp_d.begin()); -Gpu::Device::streamSynchronize(); - -constautowrite_int_comp_d_ptr=write_int_comp_d.data(); -constautowrite_real_comp_d_ptr=write_real_comp_d.data(); - -std::size_tpoffset=0; -for(inttile:tiles){ -constauto&ptile=pc.ParticlesAt(lev,grid,tile); -constauto&pflags=particle_io_flags[lev].at(std::make_pair(grid,tile)); -intnp_tile=ptile.numParticles(); -typenamePC::IntVectoroffsets(np_tile); -intnum_copies=Scan::ExclusiveSum(np_tile,pflags.begin(),offsets.begin(),Scan::retSum); - -constautoflag_ptr=pflags.data(); - -autoidata_d_ptr=idata_d.data(); -autordata_d_ptr=rdata_d.data(); - -constautoptd=ptile.getConstParticleTileData(); -amrex::ParallelFor(num_copies, -[=]AMREX_GPU_DEVICE(intpindex)noexcept -{ -//mightbeworthusingsharedmemoryhere -constautop=ptd.getSuperParticle(pindex); - -if(flag_ptr[pindex]){ -std::size_tiout_index=(pindex+poffset)*iChunkSize; -packParticleIDs(&idata_d_ptr[iout_index],p,is_checkpoint); -iout_index+=2; - -std::size_trout_index=(pindex+poffset)*rChunkSize; -for(intj=0;j<AMREX_SPACEDIM;j++){ -rdata_d_ptr[rout_index]=p.pos(j); -rout_index++; -} - -for(intj=0;j<PC::SuperParticleType::NInt;j++){ -if(write_int_comp_d_ptr[j]){ -idata_d_ptr[iout_index]=p.idata(j); -iout_index++; -} -} - -for(intj=0;j<ptd.m_num_runtime_int;j++){ -if(write_int_comp_d_ptr[PC::SuperParticleType::NInt+j]){ -idata_d_ptr[iout_index]=ptd.m_runtime_idata[j][pindex]; -iout_index++; -} -} - -//extraSoARealcomponents -constintreal_start_offset=PC::ParticleType::is_soa_particle?AMREX_SPACEDIM:0;//pureSoA:skippositions -for(intj=real_start_offset;j<PC::SuperParticleType::NReal;j++){ -constintwrite_comp_index=j-real_start_offset; -if(write_real_comp_d_ptr[write_comp_index]){ -rdata_d_ptr[rout_index]=p.rdata(j); -rout_index++; -} -} - -for(intj=0;j<ptd.m_num_runtime_real;j++){ -if(write_real_comp_d_ptr[PC::SuperParticleType::NReal+j-real_start_offset]){ -rdata_d_ptr[rout_index]=ptd.m_runtime_rdata[j][pindex]; -rout_index++; -} -} -} -}); - -poffset+=num_copies; -} - -Gpu::copyAsync(Gpu::deviceToHost,idata_d.begin(),idata_d.end(),idata.begin()); -Gpu::copyAsync(Gpu::deviceToHost,rdata_d.begin(),rdata_d.end(),rdata.begin()); -Gpu::Device::streamSynchronize(); -} - -template<classPC> -std::enable_if_t<!RunOnGpu<typenamePC::templateAllocatorType<int>>::value> -packIOData(Vector<int>&idata,Vector<ParticleReal>&rdata,constPC&pc,intlev,intgrid, -constVector<int>&write_real_comp,constVector<int>&write_int_comp, -constVector<std::map<std::pair<int,int>,typenamePC::IntVector>>&particle_io_flags, -constVector<int>&tiles,intnp,boolis_checkpoint) -{ -intnum_output_int=0; -for(inti=0;i<pc.NumIntComps()+PC::NStructInt;++i){ -if(write_int_comp[i]){++num_output_int;} -} - -constLongiChunkSize=2+num_output_int; -idata.resize(np*iChunkSize); - -intnum_output_real=0; -for(inti:write_real_comp){ -if(i){++num_output_real;} -} - -constLongrChunkSize=AMREX_SPACEDIM+num_output_real; -rdata.resize(np*rChunkSize); - -int*iptr=idata.dataPtr(); -ParticleReal*rptr=rdata.dataPtr(); -for(inttile:tiles){ -constauto&ptile=pc.ParticlesAt(lev,grid,tile); -constauto&pflags=particle_io_flags[lev].at(std::make_pair(grid,tile)); -for(intpindex=0;pindex<ptile.numParticles();++pindex){ -if(pflags[pindex]){ -constauto&soa=ptile.GetStructOfArrays(); - -//note:forpureSoAparticlelayouts,wedowritetheid,cpuandpositionsasastruct -//forbackwardscompatibilitywithreaders -ifconstexpr(!PC::ParticleType::is_soa_particle) -{ -constauto&aos=ptile.GetArrayOfStructs(); -constauto&p=aos[pindex]; - -//Int:id,cpu -packParticleIDs(iptr,p,is_checkpoint); -iptr+=2; - -//Real:positions -for(intj=0;j<AMREX_SPACEDIM;j++){rptr[j]=p.pos(j);} -rptr+=AMREX_SPACEDIM; - -//extraAoSIntcomponents -for(intj=0;j<PC::NStructInt;j++){ -if(write_int_comp[j]){ -*iptr=p.idata(j); -++iptr; -} -} -//extraAoSRealcomponents -for(intj=0;j<PC::NStructReal;j++){ -if(write_real_comp[j]){ -*rptr=p.rdata(j); -++rptr; -} -} -} -else{ -uint64_tidcpu=soa.GetIdCPUData()[pindex]; -if(is_checkpoint){ -std::int32_txi,yi; -std::uint32_txu,yu; -xu=(std::uint32_t)((idcpu&0xFFFFFFFF00000000LL)>>32); -yu=(std::uint32_t)(idcpu&0xFFFFFFFFLL); -std::memcpy(&xi,&xu,sizeof(xu)); -std::memcpy(&yi,&yu,sizeof(yu)); -*iptr=xi; -iptr+=1; -*iptr=yi; -iptr+=1; -}else{ -//Int:id,cpu -*iptr=(int)ParticleIDWrapper(idcpu); -iptr+=1; -*iptr=(int)ParticleCPUWrapper(idcpu); -iptr+=1; -} - -//Real:position -for(intj=0;j<AMREX_SPACEDIM;j++){rptr[j]=soa.GetRealData(j)[pindex];} -rptr+=AMREX_SPACEDIM; -} - -//SoAintdata -constintint_start_offset=0; -for(intj=int_start_offset;j<pc.NumIntComps();j++){ -if(write_int_comp[PC::NStructInt+j]){ -*iptr=soa.GetIntData(j)[pindex]; -++iptr; -} -} - -//extraSoARealcomponents -constintreal_start_offset=PC::ParticleType::is_soa_particle?AMREX_SPACEDIM:0;//pureSoA:skippositions -for(intj=real_start_offset;j<pc.NumRealComps();j++){ -constintwrite_comp_index=PC::NStructReal+j-real_start_offset; -if(write_real_comp[write_comp_index]){ -*rptr=(ParticleReal)soa.GetRealData(j)[pindex]; -++rptr; -} -} -} -} -} -} +typenamePC::IntVectorwrite_int_comp_d(write_int_comp.size()); +typenamePC::IntVectorwrite_real_comp_d(write_real_comp.size()); +Gpu::copyAsync(Gpu::hostToDevice,write_int_comp.begin(),write_int_comp.end(), +write_int_comp_d.begin()); +Gpu::copyAsync(Gpu::hostToDevice,write_real_comp.begin(),write_real_comp.end(), +write_real_comp_d.begin()); + +constautowrite_int_comp_d_ptr=write_int_comp_d.data(); +constautowrite_real_comp_d_ptr=write_real_comp_d.data(); + +std::size_tpoffset=0; +for(inttile:tiles){ +constauto&ptile=pc.ParticlesAt(lev,grid,tile); +constauto&pflags=particle_io_flags[lev].at(std::make_pair(grid,tile)); +intnp_tile=ptile.numParticles(); +typenamePC::IntVectoroffsets(np_tile); +intnum_copies=Scan::ExclusiveSum(np_tile,pflags.begin(),offsets.begin(),Scan::retSum); + +typenamePC::IntVectoridata_d(num_copies*iChunkSize); +typenamePC::RealVectorrdata_d(num_copies*rChunkSize); + +constautoflag_ptr=pflags.data(); + +autoidata_d_ptr=idata_d.data(); +autordata_d_ptr=rdata_d.data(); + +constautoptd=ptile.getConstParticleTileData(); +amrex::ParallelFor(num_copies, +[=]AMREX_GPU_DEVICE(intpindex)noexcept +{ +//mightbeworthusingsharedmemoryhere +constautop=ptd.getSuperParticle(pindex); + +if(flag_ptr[pindex]){ +std::size_tiout_index=pindex*iChunkSize; +packParticleIDs(&idata_d_ptr[iout_index],p,is_checkpoint); +iout_index+=2; + +std::size_trout_index=pindex*rChunkSize; +for(intj=0;j<AMREX_SPACEDIM;j++){ +rdata_d_ptr[rout_index]=p.pos(j); +rout_index++; +} + +for(intj=0;j<PC::SuperParticleType::NInt;j++){ +if(write_int_comp_d_ptr[j]){ +idata_d_ptr[iout_index]=p.idata(j); +iout_index++; +} +} + +for(intj=0;j<ptd.m_num_runtime_int;j++){ +if(write_int_comp_d_ptr[PC::SuperParticleType::NInt+j]){ +idata_d_ptr[iout_index]=ptd.m_runtime_idata[j][pindex]; +iout_index++; +} +} + +//extraSoARealcomponents +constintreal_start_offset=PC::ParticleType::is_soa_particle?AMREX_SPACEDIM:0;//pureSoA:skippositions +for(intj=real_start_offset;j<PC::SuperParticleType::NReal;j++){ +constintwrite_comp_index=j-real_start_offset; +if(write_real_comp_d_ptr[write_comp_index]){ +rdata_d_ptr[rout_index]=p.rdata(j); +rout_index++; +} +} + +for(intj=0;j<ptd.m_num_runtime_real;j++){ +if(write_real_comp_d_ptr[PC::SuperParticleType::NReal+j-real_start_offset]){ +rdata_d_ptr[rout_index]=ptd.m_runtime_rdata[j][pindex]; +rout_index++; +} +} +} +}); + +Gpu::copyAsync(Gpu::deviceToHost,idata_d.begin(),idata_d.end(), +idata.begin()+typenamePC::IntVector::difference_type(poffset)); +Gpu::copyAsync(Gpu::deviceToHost,rdata_d.begin(),rdata_d.end(), +rdata.begin()+typenamePC::RealVector::difference_type(poffset)); +Gpu::Device::streamSynchronize(); + +poffset+=num_copies; +} +} + +template<classPC> +std::enable_if_t<!RunOnGpu<typenamePC::templateAllocatorType<int>>::value> +packIOData(Vector<int>&idata,Vector<ParticleReal>&rdata,constPC&pc,intlev,intgrid, +constVector<int>&write_real_comp,constVector<int>&write_int_comp, +constVector<std::map<std::pair<int,int>,typenamePC::IntVector>>&particle_io_flags, +constVector<int>&tiles,intnp,boolis_checkpoint) +{ +intnum_output_int=0; +for(inti=0;i<pc.NumIntComps()+PC::NStructInt;++i){ +if(write_int_comp[i]){++num_output_int;} +} + +constLongiChunkSize=2+num_output_int; +idata.resize(np*iChunkSize); + +intnum_output_real=0; +for(inti:write_real_comp){ +if(i){++num_output_real;} +} + +constLongrChunkSize=AMREX_SPACEDIM+num_output_real; +rdata.resize(np*rChunkSize); + +int*iptr=idata.dataPtr(); +ParticleReal*rptr=rdata.dataPtr(); +for(inttile:tiles){ +constauto&ptile=pc.ParticlesAt(lev,grid,tile); +constauto&pflags=particle_io_flags[lev].at(std::make_pair(grid,tile)); +for(intpindex=0;pindex<ptile.numParticles();++pindex){ +if(pflags[pindex]){ +constauto&soa=ptile.GetStructOfArrays(); + +//note:forpureSoAparticlelayouts,wedowritetheid,cpuandpositionsasastruct +//forbackwardscompatibilitywithreaders +ifconstexpr(!PC::ParticleType::is_soa_particle) +{ +constauto&aos=ptile.GetArrayOfStructs(); +constauto&p=aos[pindex]; + +//Int:id,cpu +packParticleIDs(iptr,p,is_checkpoint); +iptr+=2; + +//Real:positions +for(intj=0;j<AMREX_SPACEDIM;j++){rptr[j]=p.pos(j);} +rptr+=AMREX_SPACEDIM; + +//extraAoSIntcomponents +for(intj=0;j<PC::NStructInt;j++){ +if(write_int_comp[j]){ +*iptr=p.idata(j); +++iptr; +} +} +//extraAoSRealcomponents +for(intj=0;j<PC::NStructReal;j++){ +if(write_real_comp[j]){ +*rptr=p.rdata(j); +++rptr; +} +} +} +else{ +uint64_tidcpu=soa.GetIdCPUData()[pindex]; +if(is_checkpoint){ +std::int32_txi,yi; +std::uint32_txu,yu; +xu=(std::uint32_t)((idcpu&0xFFFFFFFF00000000LL)>>32); +yu=(std::uint32_t)(idcpu&0xFFFFFFFFLL); +std::memcpy(&xi,&xu,sizeof(xu)); +std::memcpy(&yi,&yu,sizeof(yu)); +*iptr=xi; +iptr+=1; +*iptr=yi; +iptr+=1; +}else{ +//Int:id,cpu +*iptr=(int)ParticleIDWrapper(idcpu); +iptr+=1; +*iptr=(int)ParticleCPUWrapper(idcpu); +iptr+=1; +} + +//Real:position +for(intj=0;j<AMREX_SPACEDIM;j++){rptr[j]=soa.GetRealData(j)[pindex];} +rptr+=AMREX_SPACEDIM; +} + +//SoAintdata +constintint_start_offset=0; +for(intj=int_start_offset;j<pc.NumIntComps();j++){ +if(write_int_comp[PC::NStructInt+j]){ +*iptr=soa.GetIntData(j)[pindex]; +++iptr; +} +} + +//extraSoARealcomponents +constintreal_start_offset=PC::ParticleType::is_soa_particle?AMREX_SPACEDIM:0;//pureSoA:skippositions +for(intj=real_start_offset;j<pc.NumRealComps();j++){ +constintwrite_comp_index=PC::NStructReal+j-real_start_offset; +if(write_real_comp[write_comp_index]){ +*rptr=(ParticleReal)soa.GetRealData(j)[pindex]; +++rptr; +} +} +} +} +} } - -template<classPC,classF,std::enable_if_t<IsParticleContainer<PC>::value,int>foo=0> -voidWriteBinaryParticleDataSync(PCconst&pc, -conststd::string&dir,conststd::string&name, -constVector<int>&write_real_comp, -constVector<int>&write_int_comp, -constVector<std::string>&real_comp_names, -constVector<std::string>&int_comp_names, -Fconst&f,boolis_checkpoint) -{ -BL_PROFILE("WriteBinaryParticleData()"); -AMREX_ASSERT(pc.OK()); - -AMREX_ASSERT(sizeof(typenamePC::ParticleType::RealType)==4|| -sizeof(typenamePC::ParticleType::RealType)==8); - -constexprintNStructReal=PC::NStructReal; -constexprintNStructInt=PC::NStructInt; - -constintNProcs=ParallelDescriptor::NProcs(); -constintIOProcNumber=ParallelDescriptor::IOProcessorNumber(); - -ifconstexpr(PC::ParticleType::is_soa_particle){ -AMREX_ALWAYS_ASSERT(real_comp_names.size()==pc.NumRealComps()+NStructReal-AMREX_SPACEDIM);//pureSoA:skippositions -}else{ -AMREX_ALWAYS_ASSERT(real_comp_names.size()==pc.NumRealComps()+NStructReal); -} -AMREX_ALWAYS_ASSERT(int_comp_names.size()==pc.NumIntComps()+NStructInt); - -std::stringpdir=dir; -if(!pdir.empty()&&pdir[pdir.size()-1]!='/'){pdir+='/';} -pdir+=name; - -if(!pc.GetLevelDirectoriesCreated()){ -if(ParallelDescriptor::IOProcessor()) -{ -if(!amrex::UtilCreateDirectory(pdir,0755)) -{ -amrex::CreateDirectoryFailed(pdir); -} -} -ParallelDescriptor::Barrier(); -} - -std::ofstreamHdrFile; - -Longnparticles=0; -Longmaxnextid; - -//evaluatefforeveryparticletodeterminewhichonestooutput -Vector<std::map<std::pair<int,int>,typenamePC::IntVector>> -particle_io_flags(pc.GetParticles().size()); -for(intlev=0;lev<pc.GetParticles().size();lev++) -{ -constauto&pmap=pc.GetParticles(lev); -for(constauto&kv:pmap) -{ -auto&flags=particle_io_flags[lev][kv.first]; -particle_detail::fillFlags(flags,kv.second,f); -} -} - -Gpu::Device::streamSynchronize(); - -if(pc.GetUsePrePost()) -{ -nparticles=pc.GetNParticlesPrePost(); -maxnextid=pc.GetMaxNextIDPrePost(); -} -else -{ -nparticles=particle_detail::countFlags(particle_io_flags,pc); -maxnextid=PC::ParticleType::NextID(); -ParallelDescriptor::ReduceLongSum(nparticles,IOProcNumber); -PC::ParticleType::NextID(maxnextid); -ParallelDescriptor::ReduceLongMax(maxnextid,IOProcNumber); -} - -if(ParallelDescriptor::IOProcessor()) -{ -std::stringHdrFileName=pdir; - -if(!HdrFileName.empty()&&HdrFileName[HdrFileName.size()-1]!='/'){ -HdrFileName+='/'; -} - -HdrFileName+="Header"; -pc.HdrFileNamePrePost=HdrFileName; - -HdrFile.open(HdrFileName.c_str(),std::ios::out|std::ios::trunc); - -if(!HdrFile.good()){amrex::FileOpenFailed(HdrFileName);} - -// -//Firstthingwrittenisourversionstring. -//Weappend"_single"or"_double"totheversionstringindicating -//whetherwe'reusing"float"or"double"floatingpointdata. -// -std::stringversion_string=is_checkpoint?PC::CheckpointVersion():PC::PlotfileVersion(); -if(sizeof(typenamePC::ParticleType::RealType)==4) -{ -HdrFile<<version_string<<"_single"<<'\n'; -} -else -{ -HdrFile<<version_string<<"_double"<<'\n'; -} - -intnum_output_real=0; -for(inti:write_real_comp){ -if(i){++num_output_real;} -} - -intnum_output_int=0; -for(inti=0;i<pc.NumIntComps()+NStructInt;++i){ -if(write_int_comp[i]){++num_output_int;} -} - -//AMREX_SPACEDIMandNforsanitychecking. -HdrFile<<AMREX_SPACEDIM<<'\n'; - -//Thenumberofextrarealparameters -HdrFile<<num_output_real<<'\n'; - -//Realcomponentnames -for(inti=0;i<(int)real_comp_names.size();++i){ -if(write_real_comp[i]){HdrFile<<real_comp_names[i]<<'\n';} -} - -//Thenumberofextraintparameters -HdrFile<<num_output_int<<'\n'; - -//intcomponentnames -for(inti=0;i<NStructInt+pc.NumIntComps();++i){ -if(write_int_comp[i]){HdrFile<<int_comp_names[i]<<'\n';} -} - -boolis_checkpoint_legacy=true;//legacy -HdrFile<<is_checkpoint_legacy<<'\n'; - -//Thetotalnumberofparticles. -HdrFile<<nparticles<<'\n'; - -//Thevalueofnextidthatweneedtorestoreonrestart. -HdrFile<<maxnextid<<'\n'; - -//ThenthefinestleveloftheAMRhierarchy. -HdrFile<<pc.finestLevel()<<'\n'; - -//Thenthenumberofgridsateachlevel. -for(intlev=0;lev<=pc.finestLevel();lev++){ -HdrFile<<pc.ParticleBoxArray(lev).size()<<'\n'; -} -} - -//Wewanttowritethedataoutinparallel. -//We'llallowuptonOutFilesactivewritersatatime. -intnOutFiles(256); - -ParmParsepp("particles"); -pp.queryAdd("particles_nfiles",nOutFiles); -if(nOutFiles==-1){nOutFiles=NProcs;} -nOutFiles=std::max(1,std::min(nOutFiles,NProcs)); -pc.nOutFilesPrePost=nOutFiles; - -for(intlev=0;lev<=pc.finestLevel();lev++) -{ -boolgotsome; -if(pc.usePrePost) -{ -gotsome=(pc.nParticlesAtLevelPrePost[lev]>0); -} -else -{ -gotsome=(pc.NumberOfParticlesAtLevel(lev)>0); -} - -//Westoretheparticlesateachlevelintheirownsubdirectory. -std::stringLevelDir=pdir; - -if(gotsome) -{ -if(!LevelDir.empty()&&LevelDir[LevelDir.size()-1]!='/'){LevelDir+='/';} - -LevelDir=amrex::Concatenate(LevelDir.append("Level_"),lev,1); - -if(!pc.GetLevelDirectoriesCreated()) -{ -if(ParallelDescriptor::IOProcessor()){ -if(!amrex::UtilCreateDirectory(LevelDir,0755)){ -amrex::CreateDirectoryFailed(LevelDir); -} -} -ParallelDescriptor::Barrier(); -} -} - -//Writeouttheheaderforeachparticle -if(gotsome&&ParallelDescriptor::IOProcessor()){ -std::stringHeaderFileName=LevelDir; -HeaderFileName+="/Particle_H"; -std::ofstreamParticleHeader(HeaderFileName); - -pc.ParticleBoxArray(lev).writeOn(ParticleHeader); -ParticleHeader<<'\n'; - -ParticleHeader.flush(); -ParticleHeader.close(); -} - -MFInfoinfo; -info.SetAlloc(false); -MultiFabstate(pc.ParticleBoxArray(lev), -pc.ParticleDistributionMap(lev), -1,0,info); - -//Weeventuallywanttowriteoutthefilenameandtheoffset -//intothatfileintowhicheachgridofparticlesiswritten. -Vector<int>which(state.size(),0); -Vector<int>count(state.size(),0); -Vector<Long>where(state.size(),0); - -std::stringfilePrefix(LevelDir); -filePrefix+='/'; -filePrefix+=PC::DataPrefix(); -if(pc.usePrePost){ -pc.filePrefixPrePost[lev]=filePrefix; -} -boolgroupSets(false),setBuf(true); - -if(gotsome) -{ -for(NFilesIternfi(nOutFiles,filePrefix,groupSets,setBuf);nfi.ReadyToWrite();++nfi) -{ -auto&myStream=(std::ofstream&)nfi.Stream(); -pc.WriteParticles(lev,myStream,nfi.FileNumber(),which,count,where, -write_real_comp,write_int_comp,particle_io_flags,is_checkpoint); -} - -if(pc.usePrePost){ -pc.whichPrePost[lev]=which; -pc.countPrePost[lev]=count; -pc.wherePrePost[lev]=where; -}else{ -ParallelDescriptor::ReduceIntSum(which.dataPtr(),static_cast<int>(which.size()),IOProcNumber); -ParallelDescriptor::ReduceIntSum(count.dataPtr(),static_cast<int>(count.size()),IOProcNumber); -ParallelDescriptor::ReduceLongSum(where.dataPtr(),static_cast<int>(where.size()),IOProcNumber); -} -} - -if(ParallelDescriptor::IOProcessor()) -{ -if(pc.GetUsePrePost()){ -//----writetotheheaderandunlinkinCheckpointPost -}else{ -for(intj=0;j<state.size();j++) -{ -HdrFile<<which[j]<<''<<count[j]<<''<<where[j]<<'\n'; -} - -if(gotsome&&pc.doUnlink) -{ -//Unlinkanyzero-lengthdatafiles. -Vector<Long>cnt(nOutFiles,0); - -for(inti=0,N=static_cast<int>(count.size());i<N;i++){ -cnt[which[i]]+=count[i]; -} - -for(inti=0,N=static_cast<int>(cnt.size());i<N;i++) -{ -if(cnt[i]==0) -{ -std::stringFullFileName=NFilesIter::FileName(i,filePrefix); -FileSystem::Remove(FullFileName); -} -} -} -} -} -} - -if(ParallelDescriptor::IOProcessor()) -{ -HdrFile.flush(); -HdrFile.close(); -if(!HdrFile.good()) -{ -amrex::Abort("ParticleContainer::Checkpoint():problemwritingHdrFile"); -} -} -} - -template<classPC,std::enable_if_t<IsParticleContainer<PC>::value,int>foo=0> -voidWriteBinaryParticleDataAsync(PCconst&pc, -conststd::string&dir,conststd::string&name, -constVector<int>&write_real_comp, -constVector<int>&write_int_comp, -constVector<std::string>&real_comp_names, -constVector<std::string>&int_comp_names,boolis_checkpoint) -{ -BL_PROFILE("WriteBinaryParticleDataAsync"); -AMREX_ASSERT(pc.OK()); - -AMREX_ASSERT(sizeof(typenamePC::ParticleType::RealType)==4|| -sizeof(typenamePC::ParticleType::RealType)==8); - -constexprintNStructReal=PC::NStructReal; -constexprintNStructInt=PC::NStructInt; -constexprintNArrayReal=PC::NArrayReal; -constexprintNArrayInt=PC::NArrayInt; - -constintMyProc=ParallelDescriptor::MyProc(); -constintNProcs=ParallelDescriptor::NProcs(); -constintIOProcNumber=NProcs-1; - -ifconstexpr(PC::ParticleType::is_soa_particle){ -AMREX_ALWAYS_ASSERT(real_comp_names.size()==pc.NumRealComps()+NStructReal-AMREX_SPACEDIM);//pureSoA:skippositions -}else{ -AMREX_ALWAYS_ASSERT(real_comp_names.size()==pc.NumRealComps()+NStructReal); -} -AMREX_ALWAYS_ASSERT(int_comp_names.size()==pc.NumIntComps()+NStructInt); - -Vector<LayoutData<Long>>np_per_grid_local(pc.finestLevel()+1); -for(intlev=0;lev<=pc.finestLevel();lev++) -{ -np_per_grid_local[lev].define(pc.ParticleBoxArray(lev),pc.ParticleDistributionMap(lev)); -usingParIter=typenamePC::ParConstIterType; -for(ParIterpti(pc,lev);pti.isValid();++pti) -{ -intgid=pti.index(); -constauto&ptile=pc.ParticlesAt(lev,pti); -constauto&ptd=ptile.getConstParticleTileData(); -constintnp=ptile.numParticles(); - -ReduceOps<ReduceOpSum>reduce_op; -ReduceData<int>reduce_data(reduce_op); -usingReduceTuple=typenamedecltype(reduce_data)::Type; - -reduce_op.eval(np,reduce_data, -[=]AMREX_GPU_DEVICE(inti)->ReduceTuple -{ -return(ptd.id(i)>0)?1:0; -}); - -intnp_valid=amrex::get<0>(reduce_data.value(reduce_op)); -np_per_grid_local[lev][gid]+=np_valid; -} -} - -Vector<Vector<Long>>np_per_grid_global(pc.finestLevel()+1); -Longtotal_np=0; -Vector<Long>np_per_level(pc.finestLevel()+1); -for(intlev=0;lev<=pc.finestLevel();lev++) -{ -np_per_grid_global[lev].resize(np_per_grid_local[lev].size()); -ParallelDescriptor::GatherLayoutDataToVector(np_per_grid_local[lev], -np_per_grid_global[lev], -IOProcNumber); -np_per_level[lev]=std::accumulate(np_per_grid_global[lev].begin(), -np_per_grid_global[lev].end(),0L); -total_np+=np_per_level[lev]; -} - -std::stringpdir=dir; -if(!pdir.empty()&&pdir[pdir.size()-1]!='/'){pdir+='/';} -pdir+=name; - -if(MyProc==IOProcNumber) -{ -if(!pc.GetLevelDirectoriesCreated()) -{ -if(!amrex::UtilCreateDirectory(pdir,0755)) -{ -amrex::CreateDirectoryFailed(pdir); -} -} - -for(intlev=0;lev<=pc.finestLevel();lev++) -{ -std::stringLevelDir=pdir; -boolgotsome=np_per_level[lev]; - -if(gotsome) -{ -if(!LevelDir.empty()&&LevelDir[LevelDir.size()-1]!='/'){LevelDir+='/';} - -LevelDir=amrex::Concatenate(LevelDir.append("Level_"),lev,1); - -if(!pc.GetLevelDirectoriesCreated()) -{ -if(!amrex::UtilCreateDirectory(LevelDir,0755)) -{ -amrex::CreateDirectoryFailed(LevelDir); -} -} - -std::stringHeaderFileName=LevelDir; -HeaderFileName+="/Particle_H"; -std::ofstreamParticleHeader(HeaderFileName); - -pc.ParticleBoxArray(lev).writeOn(ParticleHeader); -ParticleHeader<<'\n'; - -ParticleHeader.flush(); -ParticleHeader.close(); -} -} -} -ParallelDescriptor::Barrier(); - -Longmaxnextid=PC::ParticleType::NextID(); -ParallelDescriptor::ReduceLongMax(maxnextid,IOProcNumber); - -Vector<Long>np_on_rank(NProcs,0L); -std::size_tpsize=particle_detail::PSizeInFile<ParticleReal>(write_real_comp,write_int_comp); -Vector<int64_t>rank_start_offset(NProcs); -if(MyProc==IOProcNumber) -{ -for(intlev=0;lev<=pc.finestLevel();lev++) -{ -for(intk=0;k<pc.ParticleBoxArray(lev).size();++k) -{ -intrank=pc.ParticleDistributionMap(lev)[k]; -np_on_rank[rank]+=np_per_grid_global[lev][k]; -} -} - -for(intip=0;ip<NProcs;++ip) -{ -autoinfo=AsyncOut::GetWriteInfo(ip); -rank_start_offset[ip]=(info.ispot==0)?0:static_cast<int64_t>(rank_start_offset[ip-1]+np_on_rank[ip-1]*psize); -} -} - -//maketmpparticletilesinpinnedmemorytowrite -usingPinnedPTile=ParticleTile<typenamePC::ParticleType,NArrayReal,NArrayInt, -PinnedArenaAllocator>; -automyptiles=std::make_shared<Vector<std::map<std::pair<int,int>,PinnedPTile>>>(); -myptiles->resize(pc.finestLevel()+1); -for(intlev=0;lev<=pc.finestLevel();lev++) -{ -for(MFItermfi=pc.MakeMFIter(lev);mfi.isValid();++mfi) -{ -auto&new_ptile=(*myptiles)[lev][std::make_pair(mfi.index(), -mfi.LocalTileIndex())]; - -if(np_per_grid_local[lev][mfi.index()]>0) -{ -constauto&ptile=pc.ParticlesAt(lev,mfi); - -constautonp=np_per_grid_local[lev][mfi.index()]; - -new_ptile.resize(np); - -constautoruntime_real_comps=ptile.NumRuntimeRealComps(); -constautoruntime_int_comps=ptile.NumRuntimeIntComps(); - -new_ptile.define(runtime_real_comps,runtime_int_comps); - -for(autocomp(0);comp<runtime_real_comps;++comp){ -new_ptile.push_back_real(NArrayReal+comp,np,0.); -} - -for(autocomp(0);comp<runtime_int_comps;++comp){ -new_ptile.push_back_int(NArrayInt+comp,np,0); -} - -amrex::filterParticles(new_ptile,ptile,KeepValidFilter()); -} -} -} - -intfinest_level=pc.finestLevel(); -Vector<BoxArray>bas; -Vector<DistributionMapping>dms; -for(intlev=0;lev<=pc.finestLevel();lev++) -{ -bas.push_back(pc.ParticleBoxArray(lev)); -dms.push_back(pc.ParticleDistributionMap(lev)); -} - -intnrc=pc.NumRealComps(); -intnic=pc.NumIntComps(); -intrnames_size=(int)real_comp_names.size(); - -autoRD=pc.ParticleRealDescriptor; - -AsyncOut::Submit([=]() -#ifdefined(__GNUC__)&&(__GNUC__==8)&&(__GNUC_MINOR__==1) -mutable//workaroundforbugingcc8.1 -#endif -{ -if(MyProc==IOProcNumber) -{ -std::stringHdrFileName=pdir; -std::ofstreamHdrFile; - -if(!HdrFileName.empty()&&HdrFileName[HdrFileName.size()-1]!='/'){ -HdrFileName+='/'; -} - -HdrFileName+="Header"; - -HdrFile.open(HdrFileName.c_str(),std::ios::out|std::ios::trunc); - -if(!HdrFile.good()){amrex::FileOpenFailed(HdrFileName);} - -std::stringversion_string=is_checkpoint?PC::CheckpointVersion():PC::PlotfileVersion(); -if(sizeof(typenamePC::ParticleType::RealType)==4) -{ -HdrFile<<version_string<<"_single"<<'\n'; -} -else -{ -HdrFile<<version_string<<"_double"<<'\n'; -} - -intnum_output_real=0; -for(inti=0;i<rnames_size;++i){ -if(write_real_comp[i]){++num_output_real;} -} - -intnum_output_int=0; -for(inti=0;i<nic+NStructInt;++i){ -if(write_int_comp[i]){++num_output_int;} -} - -//AMREX_SPACEDIMandNforsanitychecking. -HdrFile<<AMREX_SPACEDIM<<'\n'; - -//Thenumberofextrarealparameters -HdrFile<<num_output_real<<'\n'; - -//Realcomponentnames -for(inti=0;i<rnames_size;++i){ -if(write_real_comp[i]){HdrFile<<real_comp_names[i]<<'\n';} -} - -//Thenumberofextraintparameters -HdrFile<<num_output_int<<'\n'; - -//intcomponentnames -for(inti=0;i<NStructInt+nic;++i){ -if(write_int_comp[i]){HdrFile<<int_comp_names[i]<<'\n';} -} - -boolis_checkpoint_legacy=true;//legacy -HdrFile<<is_checkpoint_legacy<<'\n'; - -//Thetotalnumberofparticles. -HdrFile<<total_np<<'\n'; - -//Thevalueofnextidthatweneedtorestoreonrestart. -HdrFile<<maxnextid<<'\n'; - -//ThenthefinestleveloftheAMRhierarchy. -HdrFile<<finest_level<<'\n'; - -//Thenthenumberofgridsateachlevel. -for(intlev=0;lev<=finest_level;lev++){ -HdrFile<<dms[lev].size()<<'\n'; -} - -for(intlev=0;lev<=finest_level;lev++) -{ -Vector<int64_t>grid_offset(NProcs,0); -for(intk=0;k<bas[lev].size();++k) -{ -intrank=dms[lev][k]; -autoinfo=AsyncOut::GetWriteInfo(rank); -HdrFile<<info.ifile<<'' -<<np_per_grid_global[lev][k]<<'' -<<grid_offset[rank]+rank_start_offset[rank]<<'\n'; -grid_offset[rank]+=static_cast<int64_t>(np_per_grid_global[lev][k]*psize); -} -} - -HdrFile.flush(); -HdrFile.close(); -if(!HdrFile.good()) -{ -amrex::Abort("ParticleContainer::Checkpoint():problemwritingHdrFile"); -} -} - -AsyncOut::Wait();//Waitformyturn - -for(intlev=0;lev<=finest_level;lev++) -{ -//Foraeachgrid,thetilesitcontains -std::map<int,Vector<int>>tile_map; - -for(constauto&kv:(*myptiles)[lev]) -{ -constintgrid=kv.first.first; -constinttile=kv.first.second; -tile_map[grid].push_back(tile); -} - -std::stringLevelDir=pdir; -if(!LevelDir.empty()&&LevelDir[LevelDir.size()-1]!='/'){LevelDir+='/';} -LevelDir=amrex::Concatenate(LevelDir.append("Level_"),lev,1); -std::stringfilePrefix(LevelDir); -filePrefix+='/'; -filePrefix+=PC::DataPrefix(); -autoinfo=AsyncOut::GetWriteInfo(MyProc); -std::stringfile_name=amrex::Concatenate(filePrefix,info.ifile,5); -std::ofstreamofs; -ofs.open(file_name.c_str(),(info.ispot==0)?(std::ios::binary|std::ios::trunc) -:(std::ios::binary|std::ios::app)); - -for(intk=0;k<bas[lev].size();++k) -{ -intrank=dms[lev][k]; -if(rank!=MyProc){continue;} -constintgrid=k; -if(np_per_grid_local[lev][grid]==0){continue;} - -//Firstwriteouttheintegerdatainbinary. -intnum_output_int=0; -for(inti=0;i<nic+NStructInt;++i){ -if(write_int_comp[i]){++num_output_int;} -} - -constLongiChunkSize=2+num_output_int; -Vector<int>istuff(np_per_grid_local[lev][grid]*iChunkSize); -int*iptr=istuff.dataPtr(); - -for(unsignedi=0;i<tile_map[grid].size();i++){ -autoptile_index=std::make_pair(grid,tile_map[grid][i]); -constauto&pbox=(*myptiles)[lev][ptile_index]; -constautoptd=pbox.getConstParticleTileData(); -for(intpindex=0;pindex<pbox.numParticles();++pindex) -{ -constauto&soa=pbox.GetStructOfArrays(); - -constauto&p=make_particle<typename PC::ConstParticleType>{}(ptd,pindex); -if(p.id()<=0){continue;} - -//note:forpureSoAparticlelayouts,wedowritetheid,cpuandpositionsasastruct -//forbackwardscompatibilitywithreaders -ifconstexpr(!PC::ParticleType::is_soa_particle) -{ -//Ints:id,cpu -particle_detail::packParticleIDs(iptr,p,is_checkpoint); -iptr+=2; - -//extraAoSIntcomponents -for(intj=0;j<NStructInt;j++) -{ -if(write_int_comp[j]) -{ -*iptr=p.idata(j); -++iptr; -} -} -} -else{ -uint64_tidcpu=soa.GetIdCPUData()[pindex]; -if(is_checkpoint){ -std::int32_txi,yi; -std::uint32_txu,yu; -xu=(std::uint32_t)((idcpu&0xFFFFFFFF00000000LL)>>32); -yu=(std::uint32_t)(idcpu&0xFFFFFFFFLL); -std::memcpy(&xi,&xu,sizeof(xu)); -std::memcpy(&yi,&yu,sizeof(yu)); -*iptr=xi; -iptr+=1; -*iptr=yi; -iptr+=1; -}else{ -//Int:id,cpu -*iptr=(int)ParticleIDWrapper(idcpu); -iptr+=1; -*iptr=(int)ParticleCPUWrapper(idcpu); -iptr+=1; -} -} - -//extraSoAInts -constintint_start_offset=0; -for(intj=int_start_offset;j<nic;j++) -{ -if(write_int_comp[NStructInt+j]) -{ -*iptr=soa.GetIntData(j)[pindex]; -++iptr; -} -} -} -} - -writeIntData(istuff.dataPtr(),istuff.size(),ofs); -ofs.flush();//Somesystemsrequirethisflush()(probablyduetoabug) - -//WritetheRealdatainbinary. -intnum_output_real=0; -for(inti=0;i<rnames_size;++i){ -if(write_real_comp[i]){++num_output_real;} -} - -constLongrChunkSize=AMREX_SPACEDIM+num_output_real; -Vector<typenamePC::ParticleType::RealType>rstuff(np_per_grid_local[lev][grid]*rChunkSize); -typenamePC::ParticleType::RealType*rptr=rstuff.dataPtr(); - -for(unsignedi=0;i<tile_map[grid].size();i++){ -autoptile_index=std::make_pair(grid,tile_map[grid][i]); -constauto&pbox=(*myptiles)[lev][ptile_index]; -constautoptd=pbox.getConstParticleTileData(); -for(intpindex=0; -pindex<pbox.numParticles();++pindex) -{ -constauto&soa=pbox.GetStructOfArrays(); -constauto&p=make_particle<typename PC::ConstParticleType>{}(ptd,pindex); - -if(p.id()<=0){continue;} - -ifconstexpr(!PC::ParticleType::is_soa_particle) -{ -//Real:position -for(intj=0;j<AMREX_SPACEDIM;j++){rptr[j]=p.pos(j);} -rptr+=AMREX_SPACEDIM; - -//extraAoSreal -for(intj=0;j<NStructReal;j++) -{ -if(write_real_comp[j]) -{ -*rptr=p.rdata(j); -++rptr; -} -} -} -else{ -//Real:position -for(intj=0;j<AMREX_SPACEDIM;j++){rptr[j]=soa.GetRealData(j)[pindex];} -rptr+=AMREX_SPACEDIM; -} - -//extraSoAreal -constintreal_start_offset=PC::ParticleType::is_soa_particle?AMREX_SPACEDIM:0;//pureSoA:positions -for(intj=real_start_offset;j<nrc;j++) -{ -constintwrite_comp_offset=PC::ParticleType::is_soa_particle?AMREX_SPACEDIM:0;//pureSoA:skippositions -constintwrite_comp_index=PC::NStructReal+j-write_comp_offset; -if(write_real_comp[write_comp_index]) -{ -*rptr=(typenamePC::ParticleType::RealType)soa.GetRealData(j)[pindex]; -++rptr; -} -} -} -} - -if(sizeof(typenamePC::ParticleType::RealType)==4){ -writeFloatData((float*)rstuff.dataPtr(),rstuff.size(),ofs,RD); -} -elseif(sizeof(typenamePC::ParticleType::RealType)==8){ -writeDoubleData((double*)rstuff.dataPtr(),rstuff.size(),ofs,RD); -} - -ofs.flush();//Somesystemsrequirethisflush()(probablyduetoabug) -} -} -AsyncOut::Notify();//NotifyothersIamdone -}); -} - -#ifdefAMREX_USE_HDF5 -#include<AMReX_WriteBinaryParticleDataHDF5.H> -#endif - -#endif/*AMREX_WRITE_BINARY_PARTICLE_DATA_H*/ +} + +template<classPC,classF,std::enable_if_t<IsParticleContainer<PC>::value,int>foo=0> +voidWriteBinaryParticleDataSync(PCconst&pc, +conststd::string&dir,conststd::string&name, +constVector<int>&write_real_comp, +constVector<int>&write_int_comp, +constVector<std::string>&real_comp_names, +constVector<std::string>&int_comp_names, +Fconst&f,boolis_checkpoint) +{ +BL_PROFILE("WriteBinaryParticleData()"); +AMREX_ASSERT(pc.OK()); + +AMREX_ASSERT(sizeof(typenamePC::ParticleType::RealType)==4|| +sizeof(typenamePC::ParticleType::RealType)==8); + +constexprintNStructReal=PC::NStructReal; +constexprintNStructInt=PC::NStructInt; + +constintNProcs=ParallelDescriptor::NProcs(); +constintIOProcNumber=ParallelDescriptor::IOProcessorNumber(); + +ifconstexpr(PC::ParticleType::is_soa_particle){ +AMREX_ALWAYS_ASSERT(real_comp_names.size()==pc.NumRealComps()+NStructReal-AMREX_SPACEDIM);//pureSoA:skippositions +}else{ +AMREX_ALWAYS_ASSERT(real_comp_names.size()==pc.NumRealComps()+NStructReal); +} +AMREX_ALWAYS_ASSERT(int_comp_names.size()==pc.NumIntComps()+NStructInt); + +std::stringpdir=dir; +if(!pdir.empty()&&pdir[pdir.size()-1]!='/'){pdir+='/';} +pdir+=name; + +if(!pc.GetLevelDirectoriesCreated()){ +if(ParallelDescriptor::IOProcessor()) +{ +if(!amrex::UtilCreateDirectory(pdir,0755)) +{ +amrex::CreateDirectoryFailed(pdir); +} +} +ParallelDescriptor::Barrier(); +} + +std::ofstreamHdrFile; + +Longnparticles=0; +Longmaxnextid; + +//evaluatefforeveryparticletodeterminewhichonestooutput +Vector<std::map<std::pair<int,int>,typenamePC::IntVector>> +particle_io_flags(pc.GetParticles().size()); +for(intlev=0;lev<pc.GetParticles().size();lev++) +{ +constauto&pmap=pc.GetParticles(lev); +for(constauto&kv:pmap) +{ +auto&flags=particle_io_flags[lev][kv.first]; +particle_detail::fillFlags(flags,kv.second,f); +} +} + +Gpu::Device::streamSynchronize(); + +if(pc.GetUsePrePost()) +{ +nparticles=pc.GetNParticlesPrePost(); +maxnextid=pc.GetMaxNextIDPrePost(); +} +else +{ +nparticles=particle_detail::countFlags(particle_io_flags,pc); +maxnextid=PC::ParticleType::NextID(); +ParallelDescriptor::ReduceLongSum(nparticles,IOProcNumber); +PC::ParticleType::NextID(maxnextid); +ParallelDescriptor::ReduceLongMax(maxnextid,IOProcNumber); +} + +if(ParallelDescriptor::IOProcessor()) +{ +std::stringHdrFileName=pdir; + +if(!HdrFileName.empty()&&HdrFileName[HdrFileName.size()-1]!='/'){ +HdrFileName+='/'; +} + +HdrFileName+="Header"; +pc.HdrFileNamePrePost=HdrFileName; + +HdrFile.open(HdrFileName.c_str(),std::ios::out|std::ios::trunc); + +if(!HdrFile.good()){amrex::FileOpenFailed(HdrFileName);} + +// +//Firstthingwrittenisourversionstring. +//Weappend"_single"or"_double"totheversionstringindicating +//whetherwe'reusing"float"or"double"floatingpointdata. +// +std::stringversion_string=is_checkpoint?PC::CheckpointVersion():PC::PlotfileVersion(); +if(sizeof(typenamePC::ParticleType::RealType)==4) +{ +HdrFile<<version_string<<"_single"<<'\n'; +} +else +{ +HdrFile<<version_string<<"_double"<<'\n'; +} + +intnum_output_real=0; +for(inti:write_real_comp){ +if(i){++num_output_real;} +} + +intnum_output_int=0; +for(inti=0;i<pc.NumIntComps()+NStructInt;++i){ +if(write_int_comp[i]){++num_output_int;} +} + +//AMREX_SPACEDIMandNforsanitychecking. +HdrFile<<AMREX_SPACEDIM<<'\n'; + +//Thenumberofextrarealparameters +HdrFile<<num_output_real<<'\n'; + +//Realcomponentnames +for(inti=0;i<(int)real_comp_names.size();++i){ +if(write_real_comp[i]){HdrFile<<real_comp_names[i]<<'\n';} +} + +//Thenumberofextraintparameters +HdrFile<<num_output_int<<'\n'; + +//intcomponentnames +for(inti=0;i<NStructInt+pc.NumIntComps();++i){ +if(write_int_comp[i]){HdrFile<<int_comp_names[i]<<'\n';} +} + +boolis_checkpoint_legacy=true;//legacy +HdrFile<<is_checkpoint_legacy<<'\n'; + +//Thetotalnumberofparticles. +HdrFile<<nparticles<<'\n'; + +//Thevalueofnextidthatweneedtorestoreonrestart. +HdrFile<<maxnextid<<'\n'; + +//ThenthefinestleveloftheAMRhierarchy. +HdrFile<<pc.finestLevel()<<'\n'; + +//Thenthenumberofgridsateachlevel. +for(intlev=0;lev<=pc.finestLevel();lev++){ +HdrFile<<pc.ParticleBoxArray(lev).size()<<'\n'; +} +} + +//Wewanttowritethedataoutinparallel. +//We'llallowuptonOutFilesactivewritersatatime. +intnOutFiles(256); + +ParmParsepp("particles"); +pp.queryAdd("particles_nfiles",nOutFiles); +if(nOutFiles==-1){nOutFiles=NProcs;} +nOutFiles=std::max(1,std::min(nOutFiles,NProcs)); +pc.nOutFilesPrePost=nOutFiles; + +for(intlev=0;lev<=pc.finestLevel();lev++) +{ +boolgotsome; +if(pc.usePrePost) +{ +gotsome=(pc.nParticlesAtLevelPrePost[lev]>0); +} +else +{ +gotsome=(pc.NumberOfParticlesAtLevel(lev)>0); +} + +//Westoretheparticlesateachlevelintheirownsubdirectory. +std::stringLevelDir=pdir; + +if(gotsome) +{ +if(!LevelDir.empty()&&LevelDir[LevelDir.size()-1]!='/'){LevelDir+='/';} + +LevelDir=amrex::Concatenate(LevelDir.append("Level_"),lev,1); + +if(!pc.GetLevelDirectoriesCreated()) +{ +if(ParallelDescriptor::IOProcessor()){ +if(!amrex::UtilCreateDirectory(LevelDir,0755)){ +amrex::CreateDirectoryFailed(LevelDir); +} +} +ParallelDescriptor::Barrier(); +} +} + +//Writeouttheheaderforeachparticle +if(gotsome&&ParallelDescriptor::IOProcessor()){ +std::stringHeaderFileName=LevelDir; +HeaderFileName+="/Particle_H"; +std::ofstreamParticleHeader(HeaderFileName); + +pc.ParticleBoxArray(lev).writeOn(ParticleHeader); +ParticleHeader<<'\n'; + +ParticleHeader.flush(); +ParticleHeader.close(); +} + +MFInfoinfo; +info.SetAlloc(false); +MultiFabstate(pc.ParticleBoxArray(lev), +pc.ParticleDistributionMap(lev), +1,0,info); + +//Weeventuallywanttowriteoutthefilenameandtheoffset +//intothatfileintowhicheachgridofparticlesiswritten. +Vector<int>which(state.size(),0); +Vector<int>count(state.size(),0); +Vector<Long>where(state.size(),0); + +std::stringfilePrefix(LevelDir); +filePrefix+='/'; +filePrefix+=PC::DataPrefix(); +if(pc.usePrePost){ +pc.filePrefixPrePost[lev]=filePrefix; +} +boolgroupSets(false),setBuf(true); + +if(gotsome) +{ +for(NFilesIternfi(nOutFiles,filePrefix,groupSets,setBuf);nfi.ReadyToWrite();++nfi) +{ +auto&myStream=(std::ofstream&)nfi.Stream(); +pc.WriteParticles(lev,myStream,nfi.FileNumber(),which,count,where, +write_real_comp,write_int_comp,particle_io_flags,is_checkpoint); +} + +if(pc.usePrePost){ +pc.whichPrePost[lev]=which; +pc.countPrePost[lev]=count; +pc.wherePrePost[lev]=where; +}else{ +ParallelDescriptor::ReduceIntSum(which.dataPtr(),static_cast<int>(which.size()),IOProcNumber); +ParallelDescriptor::ReduceIntSum(count.dataPtr(),static_cast<int>(count.size()),IOProcNumber); +ParallelDescriptor::ReduceLongSum(where.dataPtr(),static_cast<int>(where.size()),IOProcNumber); +} +} + +if(ParallelDescriptor::IOProcessor()) +{ +if(pc.GetUsePrePost()){ +//----writetotheheaderandunlinkinCheckpointPost +}else{ +for(intj=0;j<state.size();j++) +{ +HdrFile<<which[j]<<''<<count[j]<<''<<where[j]<<'\n'; +} + +if(gotsome&&pc.doUnlink) +{ +//Unlinkanyzero-lengthdatafiles. +Vector<Long>cnt(nOutFiles,0); + +for(inti=0,N=static_cast<int>(count.size());i<N;i++){ +cnt[which[i]]+=count[i]; +} + +for(inti=0,N=static_cast<int>(cnt.size());i<N;i++) +{ +if(cnt[i]==0) +{ +std::stringFullFileName=NFilesIter::FileName(i,filePrefix); +FileSystem::Remove(FullFileName); +} +} +} +} +} +} + +if(ParallelDescriptor::IOProcessor()) +{ +HdrFile.flush(); +HdrFile.close(); +if(!HdrFile.good()) +{ +amrex::Abort("ParticleContainer::Checkpoint():problemwritingHdrFile"); +} +} +} + +template<classPC,std::enable_if_t<IsParticleContainer<PC>::value,int>foo=0> +voidWriteBinaryParticleDataAsync(PCconst&pc, +conststd::string&dir,conststd::string&name, +constVector<int>&write_real_comp, +constVector<int>&write_int_comp, +constVector<std::string>&real_comp_names, +constVector<std::string>&int_comp_names,boolis_checkpoint) +{ +BL_PROFILE("WriteBinaryParticleDataAsync"); +AMREX_ASSERT(pc.OK()); + +AMREX_ASSERT(sizeof(typenamePC::ParticleType::RealType)==4|| +sizeof(typenamePC::ParticleType::RealType)==8); + +constexprintNStructReal=PC::NStructReal; +constexprintNStructInt=PC::NStructInt; +constexprintNArrayReal=PC::NArrayReal; +constexprintNArrayInt=PC::NArrayInt; + +constintMyProc=ParallelDescriptor::MyProc(); +constintNProcs=ParallelDescriptor::NProcs(); +constintIOProcNumber=NProcs-1; + +ifconstexpr(PC::ParticleType::is_soa_particle){ +AMREX_ALWAYS_ASSERT(real_comp_names.size()==pc.NumRealComps()+NStructReal-AMREX_SPACEDIM);//pureSoA:skippositions +}else{ +AMREX_ALWAYS_ASSERT(real_comp_names.size()==pc.NumRealComps()+NStructReal); +} +AMREX_ALWAYS_ASSERT(int_comp_names.size()==pc.NumIntComps()+NStructInt); + +Vector<LayoutData<Long>>np_per_grid_local(pc.finestLevel()+1); +for(intlev=0;lev<=pc.finestLevel();lev++) +{ +np_per_grid_local[lev].define(pc.ParticleBoxArray(lev),pc.ParticleDistributionMap(lev)); +usingParIter=typenamePC::ParConstIterType; +for(ParIterpti(pc,lev);pti.isValid();++pti) +{ +intgid=pti.index(); +constauto&ptile=pc.ParticlesAt(lev,pti); +constauto&ptd=ptile.getConstParticleTileData(); +constintnp=ptile.numParticles(); + +ReduceOps<ReduceOpSum>reduce_op; +ReduceData<int>reduce_data(reduce_op); +usingReduceTuple=typenamedecltype(reduce_data)::Type; + +reduce_op.eval(np,reduce_data, +[=]AMREX_GPU_DEVICE(inti)->ReduceTuple +{ +return(ptd.id(i)>0)?1:0; +}); + +intnp_valid=amrex::get<0>(reduce_data.value(reduce_op)); +np_per_grid_local[lev][gid]+=np_valid; +} +} + +Vector<Vector<Long>>np_per_grid_global(pc.finestLevel()+1); +Longtotal_np=0; +Vector<Long>np_per_level(pc.finestLevel()+1); +for(intlev=0;lev<=pc.finestLevel();lev++) +{ +np_per_grid_global[lev].resize(np_per_grid_local[lev].size()); +ParallelDescriptor::GatherLayoutDataToVector(np_per_grid_local[lev], +np_per_grid_global[lev], +IOProcNumber); +np_per_level[lev]=std::accumulate(np_per_grid_global[lev].begin(), +np_per_grid_global[lev].end(),0L); +total_np+=np_per_level[lev]; +} + +std::stringpdir=dir; +if(!pdir.empty()&&pdir[pdir.size()-1]!='/'){pdir+='/';} +pdir+=name; + +if(MyProc==IOProcNumber) +{ +if(!pc.GetLevelDirectoriesCreated()) +{ +if(!amrex::UtilCreateDirectory(pdir,0755)) +{ +amrex::CreateDirectoryFailed(pdir); +} +} + +for(intlev=0;lev<=pc.finestLevel();lev++) +{ +std::stringLevelDir=pdir; +boolgotsome=np_per_level[lev]; + +if(gotsome) +{ +if(!LevelDir.empty()&&LevelDir[LevelDir.size()-1]!='/'){LevelDir+='/';} + +LevelDir=amrex::Concatenate(LevelDir.append("Level_"),lev,1); + +if(!pc.GetLevelDirectoriesCreated()) +{ +if(!amrex::UtilCreateDirectory(LevelDir,0755)) +{ +amrex::CreateDirectoryFailed(LevelDir); +} +} + +std::stringHeaderFileName=LevelDir; +HeaderFileName+="/Particle_H"; +std::ofstreamParticleHeader(HeaderFileName); + +pc.ParticleBoxArray(lev).writeOn(ParticleHeader); +ParticleHeader<<'\n'; + +ParticleHeader.flush(); +ParticleHeader.close(); +} +} +} +ParallelDescriptor::Barrier(); + +Longmaxnextid=PC::ParticleType::NextID(); +ParallelDescriptor::ReduceLongMax(maxnextid,IOProcNumber); + +Vector<Long>np_on_rank(NProcs,0L); +std::size_tpsize=particle_detail::PSizeInFile<ParticleReal>(write_real_comp,write_int_comp); +Vector<int64_t>rank_start_offset(NProcs); +if(MyProc==IOProcNumber) +{ +for(intlev=0;lev<=pc.finestLevel();lev++) +{ +for(intk=0;k<pc.ParticleBoxArray(lev).size();++k) +{ +intrank=pc.ParticleDistributionMap(lev)[k]; +np_on_rank[rank]+=np_per_grid_global[lev][k]; +} +} + +for(intip=0;ip<NProcs;++ip) +{ +autoinfo=AsyncOut::GetWriteInfo(ip); +rank_start_offset[ip]=(info.ispot==0)?0:static_cast<int64_t>(rank_start_offset[ip-1]+np_on_rank[ip-1]*psize); +} +} + +//maketmpparticletilesinpinnedmemorytowrite +usingPinnedPTile=ParticleTile<typenamePC::ParticleType,NArrayReal,NArrayInt, +PinnedArenaAllocator>; +automyptiles=std::make_shared<Vector<std::map<std::pair<int,int>,PinnedPTile>>>(); +myptiles->resize(pc.finestLevel()+1); +for(intlev=0;lev<=pc.finestLevel();lev++) +{ +for(MFItermfi=pc.MakeMFIter(lev);mfi.isValid();++mfi) +{ +auto&new_ptile=(*myptiles)[lev][std::make_pair(mfi.index(), +mfi.LocalTileIndex())]; + +if(np_per_grid_local[lev][mfi.index()]>0) +{ +constauto&ptile=pc.ParticlesAt(lev,mfi); + +constautonp=np_per_grid_local[lev][mfi.index()]; + +new_ptile.resize(np); + +constautoruntime_real_comps=ptile.NumRuntimeRealComps(); +constautoruntime_int_comps=ptile.NumRuntimeIntComps(); + +new_ptile.define(runtime_real_comps,runtime_int_comps); + +for(autocomp(0);comp<runtime_real_comps;++comp){ +new_ptile.push_back_real(NArrayReal+comp,np,0.); +} + +for(autocomp(0);comp<runtime_int_comps;++comp){ +new_ptile.push_back_int(NArrayInt+comp,np,0); +} + +amrex::filterParticles(new_ptile,ptile,KeepValidFilter()); +} +} +} + +intfinest_level=pc.finestLevel(); +Vector<BoxArray>bas; +Vector<DistributionMapping>dms; +for(intlev=0;lev<=pc.finestLevel();lev++) +{ +bas.push_back(pc.ParticleBoxArray(lev)); +dms.push_back(pc.ParticleDistributionMap(lev)); +} + +intnrc=pc.NumRealComps(); +intnic=pc.NumIntComps(); +intrnames_size=(int)real_comp_names.size(); + +autoRD=pc.ParticleRealDescriptor; + +AsyncOut::Submit([=]() +#ifdefined(__GNUC__)&&(__GNUC__==8)&&(__GNUC_MINOR__==1) +mutable//workaroundforbugingcc8.1 +#endif +{ +if(MyProc==IOProcNumber) +{ +std::stringHdrFileName=pdir; +std::ofstreamHdrFile; + +if(!HdrFileName.empty()&&HdrFileName[HdrFileName.size()-1]!='/'){ +HdrFileName+='/'; +} + +HdrFileName+="Header"; + +HdrFile.open(HdrFileName.c_str(),std::ios::out|std::ios::trunc); + +if(!HdrFile.good()){amrex::FileOpenFailed(HdrFileName);} + +std::stringversion_string=is_checkpoint?PC::CheckpointVersion():PC::PlotfileVersion(); +if(sizeof(typenamePC::ParticleType::RealType)==4) +{ +HdrFile<<version_string<<"_single"<<'\n'; +} +else +{ +HdrFile<<version_string<<"_double"<<'\n'; +} + +intnum_output_real=0; +for(inti=0;i<rnames_size;++i){ +if(write_real_comp[i]){++num_output_real;} +} + +intnum_output_int=0; +for(inti=0;i<nic+NStructInt;++i){ +if(write_int_comp[i]){++num_output_int;} +} + +//AMREX_SPACEDIMandNforsanitychecking. +HdrFile<<AMREX_SPACEDIM<<'\n'; + +//Thenumberofextrarealparameters +HdrFile<<num_output_real<<'\n'; + +//Realcomponentnames +for(inti=0;i<rnames_size;++i){ +if(write_real_comp[i]){HdrFile<<real_comp_names[i]<<'\n';} +} + +//Thenumberofextraintparameters +HdrFile<<num_output_int<<'\n'; + +//intcomponentnames +for(inti=0;i<NStructInt+nic;++i){ +if(write_int_comp[i]){HdrFile<<int_comp_names[i]<<'\n';} +} + +boolis_checkpoint_legacy=true;//legacy +HdrFile<<is_checkpoint_legacy<<'\n'; + +//Thetotalnumberofparticles. +HdrFile<<total_np<<'\n'; + +//Thevalueofnextidthatweneedtorestoreonrestart. +HdrFile<<maxnextid<<'\n'; + +//ThenthefinestleveloftheAMRhierarchy. +HdrFile<<finest_level<<'\n'; + +//Thenthenumberofgridsateachlevel. +for(intlev=0;lev<=finest_level;lev++){ +HdrFile<<dms[lev].size()<<'\n'; +} + +for(intlev=0;lev<=finest_level;lev++) +{ +Vector<int64_t>grid_offset(NProcs,0); +for(intk=0;k<bas[lev].size();++k) +{ +intrank=dms[lev][k]; +autoinfo=AsyncOut::GetWriteInfo(rank); +HdrFile<<info.ifile<<'' +<<np_per_grid_global[lev][k]<<'' +<<grid_offset[rank]+rank_start_offset[rank]<<'\n'; +grid_offset[rank]+=static_cast<int64_t>(np_per_grid_global[lev][k]*psize); +} +} + +HdrFile.flush(); +HdrFile.close(); +if(!HdrFile.good()) +{ +amrex::Abort("ParticleContainer::Checkpoint():problemwritingHdrFile"); +} +} + +AsyncOut::Wait();//Waitformyturn + +for(intlev=0;lev<=finest_level;lev++) +{ +//Foraeachgrid,thetilesitcontains +std::map<int,Vector<int>>tile_map; + +for(constauto&kv:(*myptiles)[lev]) +{ +constintgrid=kv.first.first; +constinttile=kv.first.second; +tile_map[grid].push_back(tile); +} + +std::stringLevelDir=pdir; +if(!LevelDir.empty()&&LevelDir[LevelDir.size()-1]!='/'){LevelDir+='/';} +LevelDir=amrex::Concatenate(LevelDir.append("Level_"),lev,1); +std::stringfilePrefix(LevelDir); +filePrefix+='/'; +filePrefix+=PC::DataPrefix(); +autoinfo=AsyncOut::GetWriteInfo(MyProc); +std::stringfile_name=amrex::Concatenate(filePrefix,info.ifile,5); +std::ofstreamofs; +ofs.open(file_name.c_str(),(info.ispot==0)?(std::ios::binary|std::ios::trunc) +:(std::ios::binary|std::ios::app)); + +for(intk=0;k<bas[lev].size();++k) +{ +intrank=dms[lev][k]; +if(rank!=MyProc){continue;} +constintgrid=k; +if(np_per_grid_local[lev][grid]==0){continue;} + +//Firstwriteouttheintegerdatainbinary. +intnum_output_int=0; +for(inti=0;i<nic+NStructInt;++i){ +if(write_int_comp[i]){++num_output_int;} +} + +constLongiChunkSize=2+num_output_int; +Vector<int>istuff(np_per_grid_local[lev][grid]*iChunkSize); +int*iptr=istuff.dataPtr(); + +for(unsignedi=0;i<tile_map[grid].size();i++){ +autoptile_index=std::make_pair(grid,tile_map[grid][i]); +constauto&pbox=(*myptiles)[lev][ptile_index]; +constautoptd=pbox.getConstParticleTileData(); +for(intpindex=0;pindex<pbox.numParticles();++pindex) +{ +constauto&soa=pbox.GetStructOfArrays(); + +constauto&p=make_particle<typename PC::ConstParticleType>{}(ptd,pindex); +if(p.id()<=0){continue;} + +//note:forpureSoAparticlelayouts,wedowritetheid,cpuandpositionsasastruct +//forbackwardscompatibilitywithreaders +ifconstexpr(!PC::ParticleType::is_soa_particle) +{ +//Ints:id,cpu +particle_detail::packParticleIDs(iptr,p,is_checkpoint); +iptr+=2; + +//extraAoSIntcomponents +for(intj=0;j<NStructInt;j++) +{ +if(write_int_comp[j]) +{ +*iptr=p.idata(j); +++iptr; +} +} +} +else{ +uint64_tidcpu=soa.GetIdCPUData()[pindex]; +if(is_checkpoint){ +std::int32_txi,yi; +std::uint32_txu,yu; +xu=(std::uint32_t)((idcpu&0xFFFFFFFF00000000LL)>>32); +yu=(std::uint32_t)(idcpu&0xFFFFFFFFLL); +std::memcpy(&xi,&xu,sizeof(xu)); +std::memcpy(&yi,&yu,sizeof(yu)); +*iptr=xi; +iptr+=1; +*iptr=yi; +iptr+=1; +}else{ +//Int:id,cpu +*iptr=(int)ParticleIDWrapper(idcpu); +iptr+=1; +*iptr=(int)ParticleCPUWrapper(idcpu); +iptr+=1; +} +} + +//extraSoAInts +constintint_start_offset=0; +for(intj=int_start_offset;j<nic;j++) +{ +if(write_int_comp[NStructInt+j]) +{ +*iptr=soa.GetIntData(j)[pindex]; +++iptr; +} +} +} +} + +writeIntData(istuff.dataPtr(),istuff.size(),ofs); +ofs.flush();//Somesystemsrequirethisflush()(probablyduetoabug) + +//WritetheRealdatainbinary. +intnum_output_real=0; +for(inti=0;i<rnames_size;++i){ +if(write_real_comp[i]){++num_output_real;} +} + +constLongrChunkSize=AMREX_SPACEDIM+num_output_real; +Vector<typenamePC::ParticleType::RealType>rstuff(np_per_grid_local[lev][grid]*rChunkSize); +typenamePC::ParticleType::RealType*rptr=rstuff.dataPtr(); + +for(unsignedi=0;i<tile_map[grid].size();i++){ +autoptile_index=std::make_pair(grid,tile_map[grid][i]); +constauto&pbox=(*myptiles)[lev][ptile_index]; +constautoptd=pbox.getConstParticleTileData(); +for(intpindex=0; +pindex<pbox.numParticles();++pindex) +{ +constauto&soa=pbox.GetStructOfArrays(); +constauto&p=make_particle<typename PC::ConstParticleType>{}(ptd,pindex); + +if(p.id()<=0){continue;} + +ifconstexpr(!PC::ParticleType::is_soa_particle) +{ +//Real:position +for(intj=0;j<AMREX_SPACEDIM;j++){rptr[j]=p.pos(j);} +rptr+=AMREX_SPACEDIM; + +//extraAoSreal +for(intj=0;j<NStructReal;j++) +{ +if(write_real_comp[j]) +{ +*rptr=p.rdata(j); +++rptr; +} +} +} +else{ +//Real:position +for(intj=0;j<AMREX_SPACEDIM;j++){rptr[j]=soa.GetRealData(j)[pindex];} +rptr+=AMREX_SPACEDIM; +} + +//extraSoAreal +constintreal_start_offset=PC::ParticleType::is_soa_particle?AMREX_SPACEDIM:0;//pureSoA:positions +for(intj=real_start_offset;j<nrc;j++) +{ +constintwrite_comp_offset=PC::ParticleType::is_soa_particle?AMREX_SPACEDIM:0;//pureSoA:skippositions +constintwrite_comp_index=PC::NStructReal+j-write_comp_offset; +if(write_real_comp[write_comp_index]) +{ +*rptr=(typenamePC::ParticleType::RealType)soa.GetRealData(j)[pindex]; +++rptr; +} +} +} +} + +if(sizeof(typenamePC::ParticleType::RealType)==4){ +writeFloatData((float*)rstuff.dataPtr(),rstuff.size(),ofs,RD); +} +elseif(sizeof(typenamePC::ParticleType::RealType)==8){ +writeDoubleData((double*)rstuff.dataPtr(),rstuff.size(),ofs,RD); +} + +ofs.flush();//Somesystemsrequirethisflush()(probablyduetoabug) +} +} +AsyncOut::Notify();//NotifyothersIamdone +}); +} + +#ifdefAMREX_USE_HDF5 +#include<AMReX_WriteBinaryParticleDataHDF5.H> +#endif + +#endif/*AMREX_WRITE_BINARY_PARTICLE_DATA_H*/ diff --git a/amrex/docs_xml/doxygen/namespaceparticle__detail.xml b/amrex/docs_xml/doxygen/namespaceparticle__detail.xml index 8bdfe0422a..f401706b3d 100644 --- a/amrex/docs_xml/doxygen/namespaceparticle__detail.xml +++ b/amrex/docs_xml/doxygen/namespaceparticle__detail.xml @@ -326,7 +326,7 @@ - + @@ -388,7 +388,7 @@ - +