1 #ifndef AMREX_WRITE_BINARY_PARTICLE_DATA_H
2 #define AMREX_WRITE_BINARY_PARTICLE_DATA_H
3 #include <AMReX_Config.H>
11 template <
typename SrcData>
13 int operator() (
const SrcData& src,
int i)
const noexcept
15 return (src.id(i) > 0);
21 template <
typename ParticleReal>
22 std::size_t
PSizeInFile (
const Vector<int>& wrc,
const Vector<int>& wic)
24 std::size_t rsize =
sizeof(ParticleReal)*std::accumulate(wrc.begin(), wrc.end(), 0);
25 std::size_t isize =
sizeof(
int)*std::accumulate(wic.begin(), wic.end(), 0);
26 return rsize + isize + AMREX_SPACEDIM*
sizeof(ParticleReal) + 2*
sizeof(
int);
29 template <
template <
class,
class>
class Container,
33 std::enable_if_t<RunOnGpu<typename Container<int, Allocator>::allocator_type>::value>
34 fillFlags (Container<int, Allocator>& pflags,
const PTile& ptile, F
const&
f)
36 const auto& ptd = ptile.getConstParticleTileData();
37 const auto np = ptile.numParticles();
39 auto flag_ptr = pflags.data();
43 const auto p = ptd.getSuperParticle(k);
45 if constexpr (IsCallable<F,decltype(p),RandomEngine>::value) {
46 flag_ptr[k] =
f(p,engine);
53 template <
template <
class,
class>
class Container,
57 std::enable_if_t<!RunOnGpu<typename Container<int, Allocator>::allocator_type>::value>
58 fillFlags (Container<int, Allocator>& pflags,
const PTile& ptile, F
const&
f)
60 const auto& ptd = ptile.getConstParticleTileData();
61 const auto np = ptile.numParticles();
63 auto flag_ptr = pflags.data();
64 for (
int k = 0; k < np; ++k) {
65 const auto p = ptd.getSuperParticle(k);
66 if constexpr (IsCallable<F,decltype(p),RandomEngine>::value) {
67 flag_ptr[k] =
f(p,RandomEngine{});
74 template <
template <
class,
class>
class Container,
class Allocator,
class PC>
75 std::enable_if_t<RunOnGpu<typename Container<int, Allocator>::allocator_type>::value, amrex::Long>
76 countFlags (
const Vector<std::map<std::pair<int,int>,Container<int,Allocator>>>& particle_io_flags,
const PC& pc)
78 ReduceOps<ReduceOpSum> reduce_op;
79 ReduceData<Long> reduce_data(reduce_op);
80 using ReduceTuple =
typename decltype(reduce_data)::Type;
82 for (
int lev = 0; lev < pc.GetParticles().
size(); lev++)
84 const auto& pmap = pc.GetParticles(lev);
85 for (
const auto& kv : pmap)
87 const auto& pflags = particle_io_flags[lev].at(kv.first);
88 const auto flag_ptr = pflags.data();
89 reduce_op.eval(pflags.size(), reduce_data,
92 return flag_ptr[i] ? 1 : 0;
96 ReduceTuple hv = reduce_data.value(reduce_op);
97 return amrex::get<0>(hv);
100 template <
template <
class,
class>
class Container,
class Allocator>
101 std::enable_if_t<RunOnGpu<typename Container<int, Allocator>::allocator_type>::value,
int>
104 ReduceOps<ReduceOpSum> reduce_op;
105 ReduceData<Long> reduce_data(reduce_op);
106 using ReduceTuple =
typename decltype(reduce_data)::Type;
108 const auto flag_ptr = pflags.data();
109 reduce_op.eval(pflags.size(), reduce_data,
112 return flag_ptr[i] ? 1 : 0;
114 ReduceTuple hv = reduce_data.value(reduce_op);
115 return amrex::get<0>(hv);
118 template <
template <
class,
class>
class Container,
class Allocator,
class PC>
119 std::enable_if_t<!RunOnGpu<typename Container<int, Allocator>::allocator_type>::value, amrex::Long>
120 countFlags (
const Vector<std::map<std::pair<int,int>,Container<int,Allocator>>>& particle_io_flags,
const PC& pc)
122 amrex::Long nparticles = 0;
123 for (
int lev = 0; lev < pc.GetParticles().
size(); lev++)
125 const auto& pmap = pc.GetParticles(lev);
126 for (
const auto& kv : pmap)
128 const auto& pflags = particle_io_flags[lev].at(kv.first);
129 for (
int k = 0; k < kv.second.numParticles(); ++k)
131 if (pflags[k]) { nparticles++; }
138 template <
template <
class,
class>
class Container,
class Allocator>
139 std::enable_if_t<!RunOnGpu<typename Container<int, Allocator>::allocator_type>::value,
int>
143 for (std::size_t k = 0; k < pflags.size(); ++k)
145 if (pflags[k]) { nparticles++; }
150 template <
typename P,
typename I>
156 std::uint32_t xu, yu;
157 xu = (std::uint32_t)((p.m_idcpu & 0xFFFFFFFF00000000LL) >> 32);
158 yu = (std::uint32_t)( p.m_idcpu & 0xFFFFFFFFLL);
170 std::enable_if_t<RunOnGpu<typename PC::template AllocatorType<int>>::value>
171 packIOData (Vector<int>& idata, Vector<ParticleReal>& rdata,
const PC& pc,
int lev,
int grid,
172 const Vector<int>& write_real_comp,
const Vector<int>& write_int_comp,
173 const Vector<std::map<std::pair<int, int>,
typename PC::IntVector>>& particle_io_flags,
174 const Vector<int>& tiles,
int np,
bool is_checkpoint)
176 int num_output_int = 0;
177 for (
int i = 0; i < pc.NumIntComps() + PC::NStructInt; ++i) {
178 if (write_int_comp[i]) { ++num_output_int; }
181 const Long iChunkSize = 2 + num_output_int;
182 idata.resize(np*iChunkSize);
184 int num_output_real = 0;
185 for (
int i : write_real_comp) {
186 if (i) { ++num_output_real; }
189 const Long rChunkSize = AMREX_SPACEDIM + num_output_real;
190 rdata.resize(np*rChunkSize);
192 typename PC::IntVector write_int_comp_d(write_int_comp.size());
193 typename PC::IntVector write_real_comp_d(write_real_comp.size());
195 write_int_comp_d.begin());
197 write_real_comp_d.begin());
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();
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);
210 typename PC::IntVector idata_d(num_copies*iChunkSize);
211 typename PC::RealVector rdata_d(num_copies*rChunkSize);
213 const auto flag_ptr = pflags.data();
215 auto idata_d_ptr = idata_d.data();
216 auto rdata_d_ptr = rdata_d.data();
218 const auto& ptd = ptile.getConstParticleTileData();
223 const auto p = ptd.getSuperParticle(pindex);
225 if (flag_ptr[pindex]) {
226 std::size_t iout_index = pindex*iChunkSize;
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);
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);
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];
251 const int real_start_offset = PC::ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0;
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);
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];
270 idata.begin() +
typename PC::IntVector::difference_type(poffset));
272 rdata.begin() +
typename PC::RealVector::difference_type(poffset));
275 poffset += num_copies;
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)
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; }
291 const Long iChunkSize = 2 + num_output_int;
292 idata.resize(np*iChunkSize);
294 int num_output_real = 0;
295 for (
int i : write_real_comp) {
296 if (i) { ++num_output_real; }
299 const Long rChunkSize = AMREX_SPACEDIM + num_output_real;
300 rdata.resize(np*rChunkSize);
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();
313 if constexpr(!PC::ParticleType::is_soa_particle)
315 const auto& aos = ptile.GetArrayOfStructs();
316 const auto& p = aos[pindex];
323 for (
int j = 0; j < AMREX_SPACEDIM; j++) { rptr[j] = p.pos(j); }
324 rptr += AMREX_SPACEDIM;
327 for (
int j = 0; j < PC::NStructInt; j++) {
328 if (write_int_comp[j]) {
334 for (
int j = 0; j < PC::NStructReal; j++) {
335 if (write_real_comp[j]) {
342 uint64_t idcpu = soa.GetIdCPUData()[pindex];
345 std::uint32_t xu, yu;
346 xu = (std::uint32_t)((idcpu & 0xFFFFFFFF00000000LL) >> 32);
347 yu = (std::uint32_t)( idcpu & 0xFFFFFFFFLL);
356 *iptr = (
int) ParticleIDWrapper(idcpu);
358 *iptr = (
int) ParticleCPUWrapper(idcpu);
363 for (
int j = 0; j < AMREX_SPACEDIM; j++) { rptr[j] = soa.GetRealData(j)[pindex]; }
364 rptr += AMREX_SPACEDIM;
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];
377 const int real_start_offset = PC::ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0;
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];
391 template <class PC, class F, std::enable_if_t<IsParticleContainer<PC>::value,
int> foo = 0>
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)
403 AMREX_ASSERT(
sizeof(
typename PC::ParticleType::RealType) == 4 ||
404 sizeof(
typename PC::ParticleType::RealType) == 8);
406 constexpr
int NStructReal = PC::NStructReal;
407 constexpr
int NStructInt = PC::NStructInt;
412 if constexpr(PC::ParticleType::is_soa_particle) {
413 AMREX_ALWAYS_ASSERT(real_comp_names.size() == pc.NumRealComps() + NStructReal - AMREX_SPACEDIM);
419 std::string pdir = dir;
420 if ( ! pdir.empty() && pdir[pdir.size()-1] !=
'/') { pdir +=
'/'; }
423 if ( ! pc.GetLevelDirectoriesCreated()) {
434 std::ofstream HdrFile;
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++)
444 const auto& pmap = pc.GetParticles(lev);
445 for (
const auto& kv : pmap)
447 auto& flags = particle_io_flags[lev][kv.first];
454 if(pc.GetUsePrePost())
456 nparticles = pc.GetNParticlesPrePost();
457 maxnextid = pc.GetMaxNextIDPrePost();
462 maxnextid = PC::ParticleType::NextID();
464 PC::ParticleType::NextID(maxnextid);
470 std::string HdrFileName = pdir;
472 if ( ! HdrFileName.empty() && HdrFileName[HdrFileName.size()-1] !=
'/') {
476 HdrFileName +=
"Header";
477 pc.HdrFileNamePrePost = HdrFileName;
479 HdrFile.open(HdrFileName.c_str(), std::ios::out|std::ios::trunc);
488 std::string version_string = is_checkpoint ? PC::CheckpointVersion() : PC::PlotfileVersion();
489 if (
sizeof(
typename PC::ParticleType::RealType) == 4)
491 HdrFile << version_string <<
"_single" <<
'\n';
495 HdrFile << version_string <<
"_double" <<
'\n';
498 int num_output_real = 0;
499 for (
int i : write_real_comp) {
500 if (i) { ++num_output_real; }
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; }
509 HdrFile << AMREX_SPACEDIM <<
'\n';
512 HdrFile << num_output_real <<
'\n';
515 for (
int i = 0; i < (
int) real_comp_names.size(); ++i ) {
516 if (write_real_comp[i]) { HdrFile << real_comp_names[i] <<
'\n'; }
520 HdrFile << num_output_int <<
'\n';
523 for (
int i = 0; i < NStructInt + pc.NumIntComps(); ++i ) {
524 if (write_int_comp[i]) { HdrFile << int_comp_names[i] <<
'\n'; }
527 bool is_checkpoint_legacy =
true;
528 HdrFile << is_checkpoint_legacy <<
'\n';
531 HdrFile << nparticles <<
'\n';
534 HdrFile << maxnextid <<
'\n';
537 HdrFile << pc.finestLevel() <<
'\n';
540 for (
int lev = 0; lev <= pc.finestLevel(); lev++) {
541 HdrFile << pc.ParticleBoxArray(lev).size() <<
'\n';
549 ParmParse
pp(
"particles");
551 if(nOutFiles == -1) { nOutFiles =
NProcs; }
553 pc.nOutFilesPrePost = nOutFiles;
555 for (
int lev = 0; lev <= pc.finestLevel(); lev++)
560 gotsome = (pc.nParticlesAtLevelPrePost[lev] > 0);
564 gotsome = (pc.NumberOfParticlesAtLevel(lev) > 0);
568 std::string LevelDir = pdir;
572 if ( ! LevelDir.empty() && LevelDir[LevelDir.size()-1] !=
'/') { LevelDir +=
'/'; }
576 if ( ! pc.GetLevelDirectoriesCreated())
589 std::string HeaderFileName = LevelDir;
590 HeaderFileName +=
"/Particle_H";
591 std::ofstream ParticleHeader(HeaderFileName);
593 pc.ParticleBoxArray(lev).writeOn(ParticleHeader);
594 ParticleHeader <<
'\n';
596 ParticleHeader.flush();
597 ParticleHeader.close();
601 info.SetAlloc(
false);
602 MultiFab state(pc.ParticleBoxArray(lev),
603 pc.ParticleDistributionMap(lev),
608 Vector<int> which(state.size(),0);
609 Vector<int > count(state.size(),0);
610 Vector<Long> where(state.size(),0);
612 std::string filePrefix(LevelDir);
614 filePrefix += PC::DataPrefix();
616 pc.filePrefixPrePost[lev] = filePrefix;
618 bool groupSets(
false), setBuf(
true);
622 for(NFilesIter nfi(nOutFiles, filePrefix, groupSets, setBuf); nfi.ReadyToWrite(); ++nfi)
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);
630 pc.whichPrePost[lev] = which;
631 pc.countPrePost[lev] = count;
632 pc.wherePrePost[lev] = where;
642 if(pc.GetUsePrePost()) {
645 for (
int j = 0; j < state.size(); j++)
647 HdrFile << which[j] <<
' ' << count[j] <<
' ' << where[j] <<
'\n';
650 if (gotsome && pc.doUnlink)
653 Vector<Long> cnt(nOutFiles,0);
655 for (
int i = 0, N=
static_cast<int>(count.size()); i < N; i++) {
656 cnt[which[i]] += count[i];
659 for (
int i = 0, N=
static_cast<int>(cnt.size()); i < N; i++)
663 std::string FullFileName = NFilesIter::FileName(i, filePrefix);
676 if ( ! HdrFile.good())
678 amrex::Abort(
"amrex::WriteBinaryParticleDataSync(): problem writing HdrFile");
683 template <class PC, std::enable_if_t<IsParticleContainer<PC>::value,
int> foo = 0>
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)
694 AMREX_ASSERT(
sizeof(
typename PC::ParticleType::RealType) == 4 ||
695 sizeof(
typename PC::ParticleType::RealType) == 8);
697 constexpr
int NStructReal = PC::NStructReal;
698 constexpr
int NStructInt = PC::NStructInt;
699 constexpr
int NArrayReal = PC::NArrayReal;
700 constexpr
int NArrayInt = PC::NArrayInt;
704 const int IOProcNumber =
NProcs - 1;
706 if constexpr(PC::ParticleType::is_soa_particle) {
707 AMREX_ALWAYS_ASSERT(real_comp_names.size() == pc.NumRealComps() + NStructReal - AMREX_SPACEDIM);
713 Vector<LayoutData<Long> > np_per_grid_local(pc.finestLevel()+1);
714 for (
int lev = 0; lev <= pc.finestLevel(); lev++)
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)
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();
725 ReduceOps<ReduceOpSum> reduce_op;
726 ReduceData<int> reduce_data(reduce_op);
727 using ReduceTuple =
typename decltype(reduce_data)::Type;
729 reduce_op.eval(np, reduce_data,
732 return (ptd.id(i) > 0) ? 1 : 0;
735 int np_valid = amrex::get<0>(reduce_data.value(reduce_op));
736 np_per_grid_local[lev][gid] += np_valid;
740 Vector<Vector<Long> > np_per_grid_global(pc.finestLevel()+1);
742 Vector<Long> np_per_level(pc.finestLevel()+1);
743 for (
int lev = 0; lev <= pc.finestLevel(); lev++)
745 np_per_grid_global[lev].resize(np_per_grid_local[lev].
size());
747 np_per_grid_global[lev],
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];
754 std::string pdir = dir;
755 if ( ! pdir.empty() && pdir[pdir.size()-1] !=
'/') { pdir +=
'/'; }
758 if (
MyProc == IOProcNumber)
760 if ( ! pc.GetLevelDirectoriesCreated())
768 for (
int lev = 0; lev <= pc.finestLevel(); lev++)
770 std::string LevelDir = pdir;
771 bool gotsome = np_per_level[lev];
775 if ( ! LevelDir.empty() && LevelDir[LevelDir.size()-1] !=
'/') { LevelDir +=
'/'; }
779 if ( ! pc.GetLevelDirectoriesCreated())
787 std::string HeaderFileName = LevelDir;
788 HeaderFileName +=
"/Particle_H";
789 std::ofstream ParticleHeader(HeaderFileName);
791 pc.ParticleBoxArray(lev).writeOn(ParticleHeader);
792 ParticleHeader <<
'\n';
794 ParticleHeader.flush();
795 ParticleHeader.close();
801 Long maxnextid = PC::ParticleType::NextID();
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)
809 for (
int lev = 0; lev <= pc.finestLevel(); lev++)
811 for (
int k = 0; k < pc.ParticleBoxArray(lev).
size(); ++k)
813 int rank = pc.ParticleDistributionMap(lev)[k];
814 np_on_rank[rank] += np_per_grid_global[lev][k];
818 for (
int ip = 0; ip <
NProcs; ++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);
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++)
832 for (MFIter mfi = pc.MakeMFIter(lev); mfi.isValid(); ++mfi)
834 auto& new_ptile = (*myptiles)[lev][std::make_pair(mfi.index(),
835 mfi.LocalTileIndex())];
837 if (np_per_grid_local[lev][mfi.index()] > 0)
839 const auto& ptile = pc.ParticlesAt(lev, mfi);
841 const auto np = np_per_grid_local[lev][mfi.index()];
843 new_ptile.resize(np);
845 const auto runtime_real_comps = ptile.NumRuntimeRealComps();
846 const auto runtime_int_comps = ptile.NumRuntimeIntComps();
848 new_ptile.define(runtime_real_comps, runtime_int_comps);
850 for (
auto comp(0); comp < runtime_real_comps; ++comp) {
851 new_ptile.push_back_real(NArrayReal+comp, np, 0.);
854 for (
auto comp(0); comp < runtime_int_comps; ++comp) {
855 new_ptile.push_back_int(NArrayInt+comp, np, 0);
863 int finest_level = pc.finestLevel();
864 Vector<BoxArray> bas;
865 Vector<DistributionMapping> dms;
866 for (
int lev = 0; lev <= pc.finestLevel(); lev++)
868 bas.push_back(pc.ParticleBoxArray(lev));
869 dms.push_back(pc.ParticleDistributionMap(lev));
872 int nrc = pc.NumRealComps();
873 int nic = pc.NumIntComps();
874 int rnames_size = (
int) real_comp_names.size();
876 auto RD = pc.ParticleRealDescriptor;
879 #
if defined(__GNUC__) && (__GNUC__ == 8) && (__GNUC_MINOR__ == 1)
883 if (
MyProc == IOProcNumber)
885 std::string HdrFileName = pdir;
886 std::ofstream HdrFile;
888 if ( ! HdrFileName.empty() && HdrFileName[HdrFileName.size()-1] !=
'/') {
892 HdrFileName +=
"Header";
894 HdrFile.open(HdrFileName.c_str(), std::ios::out|std::ios::trunc);
896 if ( ! HdrFile.good()) { amrex::FileOpenFailed(HdrFileName); }
898 std::string version_string = is_checkpoint ? PC::CheckpointVersion() : PC::PlotfileVersion();
899 if (
sizeof(
typename PC::ParticleType::RealType) == 4)
901 HdrFile << version_string <<
"_single" <<
'\n';
905 HdrFile << version_string <<
"_double" <<
'\n';
908 int num_output_real = 0;
909 for (
int i = 0; i < rnames_size; ++i) {
910 if (write_real_comp[i]) { ++num_output_real; }
913 int num_output_int = 0;
914 for (
int i = 0; i < nic + NStructInt; ++i) {
915 if (write_int_comp[i]) { ++num_output_int; }
919 HdrFile << AMREX_SPACEDIM <<
'\n';
922 HdrFile << num_output_real <<
'\n';
925 for (
int i = 0; i < rnames_size; ++i ) {
926 if (write_real_comp[i]) { HdrFile << real_comp_names[i] <<
'\n'; }
930 HdrFile << num_output_int <<
'\n';
933 for (
int i = 0; i < NStructInt + nic; ++i ) {
934 if (write_int_comp[i]) { HdrFile << int_comp_names[i] <<
'\n'; }
937 bool is_checkpoint_legacy =
true;
938 HdrFile << is_checkpoint_legacy <<
'\n';
941 HdrFile << total_np <<
'\n';
944 HdrFile << maxnextid <<
'\n';
947 HdrFile << finest_level <<
'\n';
950 for (
int lev = 0; lev <= finest_level; lev++) {
951 HdrFile << dms[lev].size() <<
'\n';
954 for (
int lev = 0; lev <= finest_level; lev++)
956 Vector<int64_t> grid_offset(
NProcs, 0);
957 for (
int k = 0; k < bas[lev].size(); ++k)
959 int rank = dms[lev][k];
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);
970 if ( ! HdrFile.good())
972 amrex::Abort(
"amrex::WriteBinaryParticleDataAsync(): problem writing HdrFile");
978 for (
int lev = 0; lev <= finest_level; lev++)
981 std::map<int, Vector<int> > tile_map;
983 for (
const auto& kv : (*myptiles)[lev])
985 const int grid = kv.first.first;
986 const int tile = kv.first.second;
987 tile_map[grid].push_back(tile);
990 std::string LevelDir = pdir;
991 if ( ! LevelDir.empty() && LevelDir[LevelDir.size()-1] !=
'/') { LevelDir +=
'/'; }
993 std::string filePrefix(LevelDir);
995 filePrefix += PC::DataPrefix();
999 ofs.open(file_name.c_str(), (info.ispot == 0) ? (std::ios::binary | std::ios::trunc)
1000 : (std::ios::binary | std::ios::app));
1002 for (
int k = 0; k < bas[lev].size(); ++k)
1004 int rank = dms[lev][k];
1005 if (rank !=
MyProc) {
continue; }
1007 if (np_per_grid_local[lev][grid] == 0) {
continue; }
1010 int num_output_int = 0;
1011 for (
int i = 0; i < nic + NStructInt; ++i) {
1012 if (write_int_comp[i]) { ++num_output_int; }
1015 const Long iChunkSize = 2 + num_output_int;
1016 Vector<int> istuff(np_per_grid_local[lev][grid]*iChunkSize);
1017 int* iptr = istuff.dataPtr();
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)
1025 const auto& soa = pbox.GetStructOfArrays();
1028 const auto& p = mp(ptd, pindex);
1029 if (p.id() <= 0) {
continue; }
1033 if constexpr(!PC::ParticleType::is_soa_particle)
1040 for (
int j = 0; j < NStructInt; j++)
1042 if (write_int_comp[j])
1050 uint64_t idcpu = soa.GetIdCPUData()[pindex];
1051 if (is_checkpoint) {
1052 std::int32_t xi, yi;
1053 std::uint32_t xu, yu;
1054 xu = (std::uint32_t)((idcpu & 0xFFFFFFFF00000000LL) >> 32);
1055 yu = (std::uint32_t)( idcpu & 0xFFFFFFFFLL);
1064 *iptr = (
int) ParticleIDWrapper(idcpu);
1066 *iptr = (
int) ParticleCPUWrapper(idcpu);
1072 const int int_start_offset = 0;
1073 for (
int j = int_start_offset; j < nic; j++)
1075 if (write_int_comp[NStructInt+j])
1077 *iptr = soa.GetIntData(j)[pindex];
1088 int num_output_real = 0;
1089 for (
int i = 0; i < rnames_size; ++i) {
1090 if (write_real_comp[i]) { ++num_output_real; }
1093 const Long rChunkSize = AMREX_SPACEDIM + num_output_real;
1094 Vector<typename PC::ParticleType::RealType> rstuff(np_per_grid_local[lev][grid]*rChunkSize);
1095 typename PC::ParticleType::RealType* rptr = rstuff.dataPtr();
1097 for (
unsigned i = 0; i < tile_map[grid].size(); i++) {
1098 auto ptile_index = std::make_pair(grid, tile_map[grid][i]);
1099 const auto& pbox = (*myptiles)[lev][ptile_index];
1100 const auto& ptd = pbox.getConstParticleTileData();
1101 for (
int pindex = 0;
1102 pindex < pbox.numParticles(); ++pindex)
1104 const auto& soa = pbox.GetStructOfArrays();
1106 const auto& p = mp(ptd, pindex);
1108 if (p.id() <= 0) {
continue; }
1110 if constexpr(!PC::ParticleType::is_soa_particle)
1113 for (
int j = 0; j < AMREX_SPACEDIM; j++) { rptr[j] = p.pos(j); }
1114 rptr += AMREX_SPACEDIM;
1117 for (
int j = 0; j < NStructReal; j++)
1119 if (write_real_comp[j])
1128 for (
int j = 0; j < AMREX_SPACEDIM; j++) { rptr[j] = soa.GetRealData(j)[pindex]; }
1129 rptr += AMREX_SPACEDIM;
1133 const int real_start_offset = PC::ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0;
1134 for (
int j = real_start_offset; j < nrc; j++)
1136 const int write_comp_offset = PC::ParticleType::is_soa_particle ? AMREX_SPACEDIM : 0;
1137 const int write_comp_index = PC::NStructReal+j-write_comp_offset;
1138 if (write_real_comp[write_comp_index])
1140 *rptr = (
typename PC::ParticleType::RealType) soa.GetRealData(j)[pindex];
1147 if (
sizeof(
typename PC::ParticleType::RealType) == 4) {
1148 writeFloatData((
float*) rstuff.dataPtr(), rstuff.size(), ofs, RD);
1150 else if (
sizeof(
typename PC::ParticleType::RealType) == 8) {
1161 #ifdef AMREX_USE_HDF5
#define BL_PROFILE(a)
Definition: AMReX_BLProfiler.H:551
#define AMREX_ASSERT(EX)
Definition: AMReX_BLassert.H:38
#define AMREX_ALWAYS_ASSERT(EX)
Definition: AMReX_BLassert.H:50
#define AMREX_GPU_DEVICE
Definition: AMReX_GpuQualifiers.H:18
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
amrex::ParmParse pp
Input file parser instance for the given namespace.
Definition: AMReX_HypreIJIface.cpp:15
void WriteBinaryParticleDataAsync(PC const &pc, const std::string &dir, const std::string &name, const Vector< int > &write_real_comp, const Vector< int > &write_int_comp, const Vector< std::string > &real_comp_names, const Vector< std::string > &int_comp_names, bool is_checkpoint)
Definition: AMReX_WriteBinaryParticleData.H:684
void WriteBinaryParticleDataSync(PC const &pc, const std::string &dir, const std::string &name, const Vector< int > &write_real_comp, const Vector< int > &write_int_comp, const Vector< std::string > &real_comp_names, const Vector< std::string > &int_comp_names, F const &f, bool is_checkpoint)
Definition: AMReX_WriteBinaryParticleData.H:392
if(!(yy_init))
Definition: amrex_iparser.lex.nolint.H:935
int queryAdd(const char *name, T &ref)
If name is found, the value in the ParmParse database will be stored in the ref argument....
Definition: AMReX_ParmParse.H:993
WriteInfo GetWriteInfo(int rank)
Definition: AMReX_AsyncOut.cpp:72
void Wait()
Definition: AMReX_AsyncOut.cpp:112
void Notify()
Definition: AMReX_AsyncOut.cpp:127
void Submit(std::function< void()> &&a_f)
Definition: AMReX_AsyncOut.cpp:95
bool Remove(std::string const &filename)
Definition: AMReX_FileSystem.cpp:190
AMREX_GPU_HOST_DEVICE Long size(T const &b) noexcept
integer version
Definition: AMReX_GpuRange.H:26
void copyAsync(HostToDevice, InIter begin, InIter end, OutIter result) noexcept
A host-to-device copy routine. Note this is just a wrapper around memcpy, so it assumes contiguous st...
Definition: AMReX_GpuContainers.H:233
static constexpr DeviceToHost deviceToHost
Definition: AMReX_GpuContainers.H:99
static constexpr HostToDevice hostToDevice
Definition: AMReX_GpuContainers.H:98
void streamSynchronize() noexcept
Definition: AMReX_GpuDevice.H:237
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void * memcpy(void *dest, const void *src, std::size_t count)
Definition: AMReX_GpuUtility.H:214
int NProcs()
Process ID in MPI_COMM_WORLD.
Definition: AMReX_MPMD.cpp:122
int MyProc()
Definition: AMReX_MPMD.cpp:117
void ReduceIntSum(int &)
Integer sum reduction.
Definition: AMReX_ParallelDescriptor.cpp:1252
void ReduceLongMax(Long &)
Long max reduction.
Definition: AMReX_ParallelDescriptor.cpp:1224
void ReduceLongSum(Long &)
Long sum reduction.
Definition: AMReX_ParallelDescriptor.cpp:1223
void GatherLayoutDataToVector(const LayoutData< T > &sendbuf, Vector< T > &recvbuf, int root)
Gather LayoutData values to a vector on root.
Definition: AMReX_ParallelDescriptor.H:1211
int IOProcessorNumber() noexcept
Definition: AMReX_ParallelDescriptor.H:266
void Barrier(const std::string &)
Definition: AMReX_ParallelDescriptor.cpp:1202
bool IOProcessor() noexcept
Is this CPU the I/O Processor? To get the rank number, call IOProcessorNumber()
Definition: AMReX_ParallelDescriptor.H:275
T ExclusiveSum(N n, T const *in, T *out, RetSum a_ret_sum=retSum)
Definition: AMReX_Scan.H:1229
static constexpr RetSum retSum
Definition: AMReX_Scan.H:29
static int f(amrex::Real t, N_Vector y_data, N_Vector y_rhs, void *user_data)
Definition: AMReX_SundialsIntegrator.H:44
@ min
Definition: AMReX_ParallelReduce.H:18
@ max
Definition: AMReX_ParallelReduce.H:17
static constexpr int P
Definition: AMReX_OpenBC.H:14
std::enable_if_t< std::is_integral_v< T > > ParallelFor(TypeList< CTOs... > ctos, std::array< int, sizeof...(CTOs)> const &runtime_options, T N, F &&f)
Definition: AMReX_CTOParallelForImpl.H:200
void writeIntData(const From *data, std::size_t size, std::ostream &os, const amrex::IntDescriptor &id)
Definition: AMReX_IntConv.H:23
void FileOpenFailed(const std::string &file)
Output a message and abort when couldn't open the file.
Definition: AMReX_Utility.cpp:131
ParIter_impl< Particle< T_NStructReal, T_NStructInt >, T_NArrayReal, T_NArrayInt, Allocator, CellAssignor > ParIter
Definition: AMReX_ParIter.H:229
void writeFloatData(const float *data, std::size_t size, std::ostream &os, const RealDescriptor &rd=FPC::Native32RealDescriptor())
Definition: AMReX_VectorIO.cpp:114
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 end(BoxND< dim > const &box) noexcept
Definition: AMReX_Box.H:1890
std::string Concatenate(const std::string &root, int num, int mindigits)
Returns rootNNNN where NNNN == num.
Definition: AMReX_String.cpp:34
Index filterParticles(DstTile &dst, const SrcTile &src, const Index *mask) noexcept
Conditionally copy particles from src to dst based on the value of mask.
Definition: AMReX_ParticleTransformation.H:325
void CreateDirectoryFailed(const std::string &dir)
Output a message and abort when couldn't create the directory.
Definition: AMReX_Utility.cpp:123
bool UtilCreateDirectory(const std::string &path, mode_t mode, bool verbose=false)
Creates the specified directories. path may be either a full pathname or a relative pathname....
Definition: AMReX_Utility.cpp:116
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void ignore_unused(const Ts &...)
This shuts up the compiler about unused variables.
Definition: AMReX.H:111
void writeDoubleData(const double *data, std::size_t size, std::ostream &os, const RealDescriptor &rd=FPC::Native64RealDescriptor())
Definition: AMReX_VectorIO.cpp:126
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Dim3 begin(BoxND< dim > const &box) noexcept
Definition: AMReX_Box.H:1881
void Abort(const std::string &msg)
Print out message to cerr and exit via abort().
Definition: AMReX.cpp:225
const int[]
Definition: AMReX_BLProfiler.cpp:1664
AMREX_ATTRIBUTE_FLATTEN_FOR void ParallelForRNG(T n, L const &f) noexcept
Definition: AMReX_GpuLaunchFunctsC.H:1221
Definition: AMReX_WriteBinaryParticleData.H:19
std::enable_if_t< RunOnGpu< typename Container< int, Allocator >::allocator_type >::value > fillFlags(Container< int, Allocator > &pflags, const PTile &ptile, F const &f)
Definition: AMReX_WriteBinaryParticleData.H:34
AMREX_GPU_HOST_DEVICE void packParticleIDs(I *idata, const P &p, bool is_checkpoint) noexcept
Definition: AMReX_WriteBinaryParticleData.H:152
std::enable_if_t< RunOnGpu< typename PC::template AllocatorType< int > >::value > packIOData(Vector< int > &idata, Vector< ParticleReal > &rdata, const PC &pc, int lev, int grid, const Vector< int > &write_real_comp, const Vector< int > &write_int_comp, const Vector< std::map< std::pair< int, int >, typename PC::IntVector >> &particle_io_flags, const Vector< int > &tiles, int np, bool is_checkpoint)
Definition: AMReX_WriteBinaryParticleData.H:171
std::enable_if_t< RunOnGpu< typename Container< int, Allocator >::allocator_type >::value, amrex::Long > countFlags(const Vector< std::map< std::pair< int, int >, Container< int, Allocator >>> &particle_io_flags, const PC &pc)
Definition: AMReX_WriteBinaryParticleData.H:76
std::size_t PSizeInFile(const Vector< int > &wrc, const Vector< int > &wic)
Definition: AMReX_WriteBinaryParticleData.H:22
Definition: AMReX_WriteBinaryParticleData.H:10
AMREX_GPU_HOST_DEVICE int operator()(const SrcData &src, int i) const noexcept
Definition: AMReX_WriteBinaryParticleData.H:13
Definition: AMReX_RandomEngine.H:57
Definition: AMReX_MakeParticle.H:16