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